# Include <math. h>
# Include <stdio. h> ////////////////////////////////////// //////////////////////////////////////// ////////////
// Matrix Structure
Struct Matrix
{
Int M, N; // M indicates the number of rows, and N indicates the number of columns.
Double ** PM; // pointer to a matrix two-dimensional array
};
// Initialize the matrix Mt. The behavior m of the parallel matrix is column N.
Void initmatrix (struct matrix * Mt, int M, int N)
{
Int I;
// Specify the rows and columns of the matrix
MT-> M = m;
MT-> N = N;
// Allocate memory for the matrix
MT-> PM = new double * [m];
For (I = 0; I <m; I ++)
{
MT-> PM [I] = new double [N];
}
}
// Use 0 to initialize the matrix Mt, and the behavior m of the parallel matrix is column N
Void initmatrix0 (struct matrix * Mt, int M, int N)
{
Int I, J;
Initmatrix (MT, m, n );
For (I = 0; I <m; I ++)
For (j = 0; j <n; j ++) MT-> PM [I] [J] = 0.0;
}
// Destroy the matrix Mt
Void destroymatrix (struct matrix * MT)
{
Int I;
// Release the matrix memory
For (I = 0; I <MT-> m; I ++)
{
Delete [] MT-> PM [I];
}
Delete [] MT-> PM;
}
// Multiply the matrix by mt3 = MT1 * MT2
// 1 is returned for success and 0 is returned for failure.
Int matrixmul (matrix * MT1, matrix * MT2, matrix * mt3)
{
Int I, J, K;
Double S;
// Check whether the row and column numbers match
If (MT1-> n! = MT2-> M | MT1-> M! = Mt3-> M | MT2-> n! = Mt3-> N) return 0;
For (I = 0; I <MT1-> m; I ++)
For (j = 0; j <MT2-> N; j ++)
{
S = 0.0;
For (k = 0; k <MT1-> N; k ++) S = S + MT1-> PM [I] [k] * MT2-> PM [k] [J];
Mt3-> PM [I] [J] = s;
}
Return 1;
}
// Matrix transpose MT2 = T (MT1)
// 1 is returned for success and 0 is returned for failure.
Int matrixtran (matrix * MT1, matrix * MT2)
{
Int I, J;
// Check whether the row and column numbers match
If (MT1-> M! = MT2-> N | MT1-> n! = MT2-> m) return 0;
For (I = 0; I <MT1-> m; I ++)
For (j = 0; j <MT1-> N; j ++) MT2-> PM [J] [I] = MT1-> PM [I] [J];
Return 1;
}
// Console display matrix Mt
Void consoleshowmatrix (matrix * MT)
{
Int I, J;
For (I = 0; I <MT-> m; I ++)
{
Printf ("/N ");
For (j = 0; j <MT-> N; j ++) printf ("% 2.4f", MT-> PM [I] [J]);
}
Printf ("/N ");
}
//////////////////////////////////////// //////////////////////////////////////// //////////
// Solve the equation Ax = B, A = (a, B) by using the principal element elimination method of the guass Column)
Int guassl (double ** A, double * X, int N)
{
Int I, J, K, numl, * h, T;
Double * l, Max;
L = new double [N];
H = new int [N];
For (I = 0; I <n; I ++) H [I] = I; // row mark
For (I = 1; I <n; I ++)
{
Max = FABS (A [H [I-1] [I-1]);
Numl = I-1;
// Maximum value of column element
For (j = I; j <n; j ++)
{
If (FABS (A [H [J] [I-1])> MAX)
{
Numl = H [J];
Max = FABS (A [H [J] [I-1]);
}
}
If (max <0.00000000001) return 0;
// Exchange the row number
If (numl> i-1)
{
T = H [I];
H [I] = H [numl];
H [numl] = T;
}
For (j = I; j <n; j ++) L [J] = A [H [J] [I-1]/A [H [I-1] [I-1];
For (j = I; j <n; j ++)
For (k = I; k <n + 1; k ++) A [H [J] [k] = A [H [J] [k]-l [J] * A [H [I-1] [k];
}
For (I = n-1; I> = 0; I --)
{
X [I] = A [H [I] [N];
For (j = I + 1; j <n; j ++) X [I] = x [I]-A [H [I] [J] * X [J];
X [I] = x [I]/A [H [I] [I];
}
// Clear the temporary Array Memory
Delete [] L;
Delete [] h;
Return 1;
}/////////////////////////////////////// //////////////////////////////////////// ////////////
// Least square algorithm for solving the matrix Ax = B
Int minmul (matrix A, matrix B, double * X)
{
Int I, J;
If (B. n! = 1) return 0;
If (A. M! = B. m) return 0;
Matrix trana; // defines the transpose of
Initmatrix0 (& trana, A. N, A. m );
Matrixtran (& A, & trana );
Matrix trana_a; // defines the product matrix of A's transpose and
Initmatrix0 (& trana_a, A. N, A. N );
Matrixmul (& trana, & A, & trana_a); // the product of a's transpose and
Matrix trana_ B; // defines the product matrix of the transpose of A and B.
Initmatrix0 (& trana_ B, A. N, 1 );
Matrixmul (& trana, & B, & trana_ B); // the product of the transpose of A and B
Destroymatrix (& trana); // release the transposed memory of
Matrix trana_a_ B; // defines the augmented matrix.
Initmatrix0 (& trana_a_ B, trana_a.m, trana_a.m + 1 );
// Assign values to the Augmented Matrix
For (I = 0; I <trana_a_ B .m; I ++)
{
For (j = 0; j <trana_a_ B .m; j ++) trana_a_ B .pm [I] [J] = trana_a.pm [I] [J];
Trana_a_ B .pm [I] [trana_a_ B .m] = trana_ B .pm [I] [0];
}
Destroymatrix (& trana_a );
Destroymatrix (& trana_ B );
// Solve the guass column primary Elimination Method
Guassl (trana_a_ B .pm, X, trana_a_a_ B .m );
Destroymatrix (& trana_a_ B );
Return 1;
}
Int minmul (double ** A, double * B, int M, int N, double * X)
{
Int R, I;
Matrix Al, BL;
Al. PM = new double * [m];
Al. M = m;
Al. n = N;
Initmatrix (& BL, M, 1 );
For (I = 0; I <m; I ++)
{
Al. PM [I] = A [I];
BL. PM [I] [0] = B [I];
}
R = minmul (Al, BL, X );
Delete [] Al. PM;
Destroymatrix (& BL );
Return R;
}