Matrix

矩阵乘法

bool matrix_product(const double* A,const double * B,double* C,int ARowNum,int AColNum,int BColNum){// memory errorif(A == NULL||B == NULL || C == NULL)return false;// column or row number errorif(!(ARowNum * AColNum *BColNum) )return false;// calculate matrix productfor(int rowNum = 0;rowNum < ARowNum;rowNum++){for(int colNum = 0;colNum < BColNum; colNum++){double sum = 0;for(int index = 0;index < AColNum;index++){sum+=A[rowNum*AColNum+index]*B[index*BColNum+colNum];}C[rowNum*BColNum+colNum]= sum;}}return true;}

矩阵转秩

bool matrix_transpose(const double* A,double* B,int ARowNum,int AColNum){// memory errorif(A == NULL||B == NULL || C == NULL)return false;// column or row number errorif(!(ARowNum * AColNum *BColNum) )return false;// swap data a[i][j] = a[j][i]for(int rowNum = 0;rowNum < AColNum;rowNum++){for(int colNum = 0;colNum < ARowNum; colNum++){B[rowNum*ARowNum+colNum] = A[colNum*AColNum+rowNum];}}return true;}

矩阵求逆

bool matrix_inverse(const double* A,double* A_inv,int rowNum){// memory errorif(A == NULL||B == NULL || C == NULL)return false;// column or row number errorif(!(ARowNum * AColNum *BColNum) )return false;// calculate determinant of Adouble d = calculateDeterminant(rowNum,A);// A is degradation.if(d == 0)   return false;for(int i= 0;i<rowNum;i++){for(int j=0;j<rowNum;j++){A_inv[i*rowNum+j] = CalculateAlgebraicComplement(i,j,A,rowNum)/d;}}return true;}double calculateDeterminant(int rowNum,const double* A){double* B = new double[rowNum*rowNum];memcpy(B,A,rowNum*rowNum*sizeof(double));int sign = 0;double result = 1.0;for(int i= 0;i<rowNum;i++){if(B[i*rowNum+i] == 0){if(!AdjustDeterminant(B,i,rowNum)){delete B;return 0;}sign++;}Elimination(B,i,rowNum);result*= B[i*rowNum+i];}delete B;return result*(sign%2?-1:1);}bool AdjustDeterminant(double* A,int i,int rowNum){for(int k= i+1;k<rowNum;k++){if(A[k*rowNum+i] != 0){SwapRow(A,i,k,rowNum);break;}}if(A[i*rowNum+i] == 0)return false;return true;}void Elimination(double*A,int i,int rowNum){for(int j= i+1;j<rowNum;j++){double scale = -1*(A[j*rowNum+i]/A[i*rowNum+i]);for(int k= i;k<rowNum;k++){A[j*rowNum+k] += A[i*rowNum+k]*scale;}}}void SwapRow(double* A,int i,int j,int rowNum ){double mid;for(int k=0;k<rowNum;k++){mid = A[i*rowNum+k];A[i*rowNum+k] = A[j*rowNum+k];A[j*rowNum+k] = mid;} }double CalculateAlgebraicComplement(int rowIndex,int colIndex,const double*A,int rowNum){double* adjointMatrix = GetAdjointMatrix(rowIndex,colIndex,A,rowNum);double result =  calculateDeterminant(rowNum-1,adjointMatrix);delete adjointMatrix;result *= (rowIndex+colIndex)%2?-1:1;return result;}double* GetAdjointMatrix(int rowIndex,int colIndex,const double*A,int rowNum){double* adjointMatrix = new double[(rowNum-1)*(rowNum-1)];int count = 0;for(int i= 0;i<rowNum;i++){if(i==rowIndex)continue;for(int j=0;j<rowNum;j++){if(j==colIndex)continue;adjointMatrix[count] = A[i*rowNum+j];count++;}}return adjointMatrix;}

矩阵数乘

bool matrix_scale(double *A,double scale,int rowNum,int colNum){if(A == NULL)return false;if(!(rowNum*colNum))return false;for(int i = 0;i<rowNum;i++)for(int j=0;j<colNum;j++)A[i*colNum+j]  *=scale;return true;}

矩阵加法

bool matrix_Sum(const double* A,const double*B,double* C,int ARowNum,int AColNum,int BRowNum,int BColNum){if(A == NULL||B == NULL || C == NULL)return false;if(ARowNum != BRowNum||AColNum != BColNum)return false;if(!(ARowNum*AColNum))return false;for(int i= 0;i<ARowNum;i++)for(int j=0;j<AColNum;j++)C[i*AColNum+j] = A[i*AColNum+j]+B[i*AColNum+j];return true;}