parent
fd7abbc9d4
commit
e2f545283c
|
@ -0,0 +1,207 @@
|
|||
#include<math.h>
|
||||
#include<stdarg.h>
|
||||
#include<stdbool.h>
|
||||
#include<stdio.h>
|
||||
#include<sys/param.h>
|
||||
|
||||
#define MAX_DIMENSION 8
|
||||
|
||||
/* we want to build up matrices from vectors in order to make the functions compatible */
|
||||
|
||||
typedef struct {
|
||||
int len;
|
||||
double _[MAX_DIMENSION];
|
||||
_Bool valid;
|
||||
} vector_t;
|
||||
|
||||
typedef struct {
|
||||
int rows;
|
||||
int cols;
|
||||
double _[MAX_DIMENSION][MAX_DIMENSION];
|
||||
_Bool valid;
|
||||
} matrix_t;
|
||||
|
||||
/* missing braces around initializer? */
|
||||
|
||||
const matrix_t EMPTY_MATRIX = {
|
||||
.rows = 0,
|
||||
.cols = 0,
|
||||
.valid = false,
|
||||
._ = {0}
|
||||
};
|
||||
|
||||
const vector_t EMPTY_VECTOR = {
|
||||
.len = 0,
|
||||
.valid = false,
|
||||
._ = {0}
|
||||
};
|
||||
|
||||
/* function prototypes */
|
||||
matrix_t new_matrix(int, int, ...);
|
||||
|
||||
/* arithmetic */
|
||||
matrix_t matrix_add(matrix_t, matrix_t);
|
||||
matrix_t matrix_sub(matrix_t, matrix_t);
|
||||
matrix_t matrix_scaler_mul(double, matrix_t);
|
||||
matrix_t matrix_mul(matrix_t, matrix_t);
|
||||
|
||||
/* miscellaneous, what do we need? */
|
||||
matrix_t new_zero(int, int);
|
||||
matrix_t new_identity(int diag);
|
||||
matrix_t scaler_mul(double scale, matrix_t A);
|
||||
matrix_t matrix_add(matrix_t A, matrix_t B);
|
||||
|
||||
/* helper functions */
|
||||
void matrix_print(FILE *f, matrix_t m);
|
||||
vector_t get_row(matrix_t m, int row);
|
||||
vector_t get_col(matrix_t m, int col);
|
||||
vector_t zero_vec(int len);
|
||||
double dot_product(vector_t a, vector_t b);
|
||||
|
||||
|
||||
matrix_t new_matrix(int rows, int cols, ...) {
|
||||
matrix_t m;
|
||||
m.rows = rows;
|
||||
m.cols = cols;
|
||||
m.valid = true;
|
||||
|
||||
if(MIN(rows, cols) <= 0) return EMPTY_MATRIX;
|
||||
if(MAX_DIMENSION < MAX(rows, cols)) return EMPTY_MATRIX;
|
||||
|
||||
va_list ap;
|
||||
va_start(ap, cols);
|
||||
|
||||
for(int j = 0; j < rows; j++) {
|
||||
for(int i = 0; i < cols; i++) {
|
||||
// this may crash is arg list isn't long enough or invalid args
|
||||
m._[j][i] = va_arg(ap, double);
|
||||
}
|
||||
}
|
||||
|
||||
va_end(ap);
|
||||
return m;
|
||||
}
|
||||
|
||||
matrix_t new_zero(int rows, int cols) {
|
||||
matrix_t m = EMPTY_MATRIX;
|
||||
m.rows = rows;
|
||||
m.cols = cols;
|
||||
m.valid = true;
|
||||
|
||||
if(MIN(rows, cols) <= 0) return EMPTY_MATRIX;
|
||||
if(MAX_DIMENSION < MAX(rows, cols)) return EMPTY_MATRIX;
|
||||
|
||||
return m;
|
||||
}
|
||||
|
||||
matrix_t new_identity(int diag) {
|
||||
/* works only for mxm matrix */
|
||||
if(diag <= 0 || MAX_DIMENSION < diag) return EMPTY_MATRIX;
|
||||
|
||||
matrix_t m = new_zero(diag, diag);
|
||||
for(int i = 0; i < diag; i++) m._[i][i] = 1.0;
|
||||
return m;
|
||||
}
|
||||
|
||||
matrix_t matrix_add(matrix_t m, matrix_t n) {
|
||||
if (m.rows != n.rows || m.cols != n.cols) return m;
|
||||
|
||||
for (int i = 0; i < m.rows; ++i) {
|
||||
for (int j = 0; j < n.cols; ++j) {
|
||||
m._[i][j] += n._[i][j];
|
||||
}
|
||||
}
|
||||
return m;
|
||||
}
|
||||
|
||||
matrix_t matrix_sub(matrix_t m, matrix_t n) {
|
||||
if (m.rows != n.rows || m.cols != n.cols) return m;
|
||||
|
||||
for (int i = 0; i < m.rows; ++i) {
|
||||
for (int j = 0; j < n.cols; ++j) {
|
||||
m._[i][j] -= n._[i][j];
|
||||
}
|
||||
}
|
||||
return m;
|
||||
}
|
||||
|
||||
matrix_t matrix_scaler_mul(double s, matrix_t m) {
|
||||
for (int i = 0; i < m.rows; ++i) {
|
||||
for (int j = 0; j < m.cols; ++j) {
|
||||
m._[i][j] *= s;
|
||||
}
|
||||
}
|
||||
return m;
|
||||
}
|
||||
|
||||
matrix_t matrix_mul(matrix_t m, matrix_t n) {
|
||||
if (m.cols != n.rows) return EMPTY_MATRIX; /* required */
|
||||
|
||||
/* new matrix for new dimensions */
|
||||
matrix_t o = new_zero(m.rows, n.cols);
|
||||
|
||||
for (int j = 0; j < m.rows; ++j) {
|
||||
vector_t ar = get_row(m, j);
|
||||
for (int i = 0; i < m.cols; ++i) {
|
||||
vector_t Bc = get_col(n, i);
|
||||
m._[j][i] = dot_product(ar, Bc);
|
||||
}
|
||||
}
|
||||
return o;
|
||||
}
|
||||
|
||||
vector_t zero_vec(int len) {
|
||||
if(len <= 0 || MAX_DIMENSION < len) return EMPTY_VECTOR;
|
||||
vector_t v = EMPTY_VECTOR;
|
||||
v.len = len;
|
||||
v.valid = true;
|
||||
return v;
|
||||
}
|
||||
|
||||
vector_t get_row(matrix_t m, int row) {
|
||||
if(row < 0 || m.rows <= row) return EMPTY_VECTOR;
|
||||
vector_t v = zero_vec(m.cols);
|
||||
for(int i = 0; i < v.len; i++) v._[i] = m._[row][i];
|
||||
return v;
|
||||
}
|
||||
|
||||
vector_t get_col(matrix_t m, int col) {
|
||||
if(col < 0 || m.cols <= col) return EMPTY_VECTOR;
|
||||
vector_t v = zero_vec(m.rows);
|
||||
for(int i = 0; i < v.len; i++) v._[i] = m._[i][col];
|
||||
return v;
|
||||
}
|
||||
|
||||
double dot_product(vector_t a, vector_t b) {
|
||||
if(!a.valid || !b.valid || a.len != b.len || a.len == 0) return NAN;
|
||||
double acc = 0.0;
|
||||
for(int i = 0; i < a.len; i++) acc += a._[i]*b._[i];
|
||||
return acc;
|
||||
}
|
||||
|
||||
void matrix_print(FILE *f, matrix_t m)
|
||||
{
|
||||
char *bump = "";
|
||||
fprintf(f, "[");
|
||||
|
||||
for(int j = 0; j < m.rows; j++) {
|
||||
bump = "";
|
||||
char *split = "";
|
||||
|
||||
fprintf(f, "%s",
|
||||
bump);
|
||||
fprintf(f, "[");
|
||||
|
||||
for(int i = 0; i < m.cols; i++) {
|
||||
split = ", ";
|
||||
|
||||
fprintf(f, "%s%4.2lf",
|
||||
split, m._[j][i]);
|
||||
}
|
||||
|
||||
fprintf(f, "]");
|
||||
if (j < m.rows-1) fprintf(f, ",\n");
|
||||
}
|
||||
fprintf(f, "]\n\n");
|
||||
return;
|
||||
}
|
Loading…
Reference in New Issue