Common matrix calculation functions and matrix calculation functions in C Language
1. Matrix transpose Function
Void matrix_t (double ** a_matrix, const double ** B _matrix, int krow, int kline) //////////////////////////////////////// /// // //////// // a_matrix: transpose matrix // B _matrix: matrix before transpose // krow: number of rows // kline: number of columns /////////////////////////////////////// /// // {int k, k2; for (k = 0; k <krow; k ++) {for (k2 = 0; k2 <kline; k2 ++) {a_matrix [k2] [k] = B _matrix [k] [k2] ;}}
2. matrix addition (subtraction) Function
Void matrix_a (double ** a_matrix, const double ** B _matrix, const double ** c_matrix, int krow, int kline, int ktrl) //////////////////////////////////////// /// // a_matrix = B _matrix + c_matrix // krow: number of rows // kline: Number of columns // ktrl: greater than 0: addition not greater than 0: subtraction /////////////////////////////////////// /// // {int k, k2; for (k = 0; k <krow; k ++) {for (k2 = 0; k2 <k Line; k2 ++) {a_matrix [k] [k2] = B _matrix [k] [k2] + (ktrl> 0 )? C_matrix [k] [k2]:-c_matrix [k] [k2]) ;}}
3. matrix multiplication function
Void matrix_m (double ** a_matrix, const double ** B _matrix, const double ** c_matrix, int krow, int kline, int kmiddle, int ktrl) //////////////////////////////////////// /// // a_matrix = B _matrix * c_matrix // krow: number of rows // kline: Number of columns // ktrl: greater than 0: two positive number matrices are multiplied not greater than 0: positive Matrix multiplied by negative matrix /////////////////////////////////// //////////////////////////////////////// /{int k, k2, k4; double stmp; for (k = 0; k <krow; k ++) {for (k2 = 0; k2 <kline; k2 ++) {stmp = 0.0; for (k4 = 0; k4 <kmiddle; k4 ++) {stmp + = B _matrix [k] [k4] * c_matrix [k4] [k2];} a_matrix [k] [k2] = stmp ;}} if (ktrl <= 0) {for (k = 0; k <krow; k ++) {for (k2 = 0; k2 <kline; k2 ++) {a_matrix [k] [k2] =-a_matrix [k] [k2] ;}}
4. Inverse Matrix Function
Int matrix_inv (double ** a_matrix, int ndimen) //////////////////////////////////////// /// // //////// // a_matrix: matrix // ndimen: dimension /////////////////////////////////////// /// // {double tmp, tmp2, B _tmp [20], c_tmp [20]; int k, k1, k2, k3, j, I, j2, i2, kme [20], kmf [20]; i2 = j2 = 0; for (k = 0; k <ndimen; k ++) {tmp2 = 0.0; for (I = k; I <ndimen; I ++) {For (j = k; j <ndimen; j ++) {if (fabs (a_matrix [I] [j]) <= fabs (tmp2) continue; tmp2 = a_matrix [I] [j]; i2 = I; j2 = j ;}} if (i2! = K) {for (j = 0; j <ndimen; j ++) {tmp = a_matrix [i2] [j]; a_matrix [i2] [j] = a_matrix [k] [j]; a_matrix [k] [j] = tmp ;}} if (j2! = K) {for (I = 0; I <ndimen; I ++) {tmp = a_matrix [I] [j2]; a_matrix [I] [j2] = a_matrix [I] [k]; a_matrix [I] [k] = tmp ;}} kme [k] = i2; kmf [k] = j2; for (j = 0; j <ndimen; j ++) {if (j = k) {B _tmp [j] = 1.0/tmp2; c_tmp [j] = 1.0;} else {B _tmp [j] =-a_matrix [k] [j]/tmp2; c_tmp [j] = a_matrix [j] [k];} a_matrix [k] [j] = 0.0; a_matrix [j] [k] = 0.0;} for (I = 0; I <ndimen; I ++) {for (j = 0; j <ndim En; j ++) {a_matrix [I] [j] = a_matrix [I] [j] + c_tmp [I] * B _tmp [j] ;}} for (k3 = 0; k3 <ndimen; k3 ++) {k = ndimen-k3-1; k1 = kme [k]; k2 = kmf [k]; if (k1! = K) {for (I = 0; I <ndimen; I ++) {tmp = a_matrix [I] [k1]; a_matrix [I] [k1] = a_matrix [I] [k]; a_matrix [I] [k] = tmp ;}} if (k2! = K) {for (j = 0; j <ndimen; j ++) {tmp = a_matrix [k2] [j]; a_matrix [k2] [j] = a_matrix [k] [j]; a_matrix [k] [j] = tmp ;}} return (0 );}
5. Matrix jorisky decomposition Functions
Void chol (double ** a_matrix, const double ** B _matrix, int ndimen) //////////////////////////////////////// /// // input parameters: // B _matrix: symmetric positive matrix ndimen: matrix dimension // return value: // a_matrix: bottom Triangle Matrix ///////////////////////////////////// ///////////////////////////////////////{ int I, j, r; double m = 0; static double ** c_matrix; static int flag = 0; if (flag = 0) {flag = 1; c_matrix = (doub Le **) malloc (ndimen * sizeof (double *); for (I = 0; I <ndimen; I ++) c_matrix [I] = (double *) malloc (ndimen * sizeof (double);} for (I = 0; I <ndimen; I ++) {for (j = 0; j <ndimen; j ++) c_matrix [I] [j] = 0;} c_matrix [0] [0] = sqrt (B _matrix [0] [0]); for (I = 1; I <ndimen; I ++) {if (c_matrix [0] [0]! = 0) c_matrix [I] [0] = B _matrix [I] [0]/c_matrix [0] [0];} for (I = 1; I <ndimen; I ++) {for (r = 0; r <I; r ++) m = m + c_matrix [I] [r] * c_matrix [I] [r]; c_matrix [I] [I] = sqrt (B _matrix [I] [I]-m); m = 0.0; for (j = I + 1; j <ndimen; j ++) {for (r = 0; r <I; r ++) m = m + c_matrix [I] [r] * c_matrix [j] [r]; c_matrix [j] [I] = (B _matrix [I] [j]-m)/c_matrix [I] [I]; m = 0 ;}} for (I = 0; I <ndimen; I ++) {for (j = 0; j <ndimen; j ++) a_matrix [I] [j] = c_matrix [I] [j];}