高斯消元法可以用來找出一個可逆矩陣的逆矩陣。設A 為一個N * N的矩陣,其逆矩陣可被兩個分塊矩陣表示出來。將一個N * N單位矩陣 放在A 的右手邊,形成一個N * 2N的分塊矩陣B = [A,I] 。經過高斯消元法的計算程序後,矩陣B 的左手邊會變成一個單位矩陣I ,而逆矩陣A ^(-1) 會出現在B 的右手邊。假如高斯消元法不能將A 化為三角形的格式,那就代表A 是一個不可逆的矩陣。應用上,高斯消元法極少被用來求出逆矩陣。高斯消元法通常只為線性方程組求解。
//******************************** //*** 求任何一個實矩陣的逆*** //******************************** #include "stdafx.h" #include實驗結果:#include #include #include using namespace std; #define N 10 //定義方陣的最大階數為10 //函數的聲明部分 float MatDet(float *p, int n); //求矩陣的行列式 float Creat_M(float *p, int m, int n, int k); //求矩陣元素A(m, n)的代數余之式 void print(float *p, int n); //輸出矩陣n*n bool Gauss(float A[][N], float B[][N], int n); //采用部分主元的高斯消去法求方陣A的逆矩陣B int main() { float *buffer, *p; //定義數組首地址指針變量 int row, num; //定義矩陣的行數和矩陣元素個數 int i, j; float determ; //定義矩陣的行列式 float a[N][N], b[N][N]; int n; cout << "采用逆矩陣的定義法求矩陣的逆矩陣!\n"; cout << "請輸入矩陣的行數: "; cin >> row; num = 2 * row * row; buffer = (float *)calloc(num, sizeof(float)); //分配內存單元 p = buffer; if (NULL != p) { for (i = 0; i < row; i++) { cout << "Please input the number of " << i+1 << " row: "; for (j = 0; j < row; j++) { cin >> *p++; } } } else { cout << "Can't distribute memory\n"; } cout << "The original matrix : \n"; print(buffer, row); //打印該矩陣 determ = MatDet(buffer, row); //求整個矩陣的行列式 p = buffer + row * row; if (determ != 0) { cout << "The determinant of the matrix is " << determ << endl; for (i = 0; i < row; i++) //求逆矩陣 { for (j = 0; j < row; j++) { *(p+j*row+i) = Creat_M(buffer, i, j, row)/determ; } } cout << "The inverse matrix is: " << endl; print(p, row); //打印該矩陣 } else { cout << "The determinant is 0, and there is no inverse matrix!\n"; } free(buffer); //釋放內存空間 cout << "采用部分主元的高斯消去法求方陣的逆矩陣!\n"; cout << "請輸入方陣的階數: "; cin >> n; cout << "請輸入" << n << "階方陣: \n"; //輸入一個n階方陣 for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { cin >> a[i][j]; } } //運用高斯消去法求該矩陣的逆矩陣並輸出 if (Gauss(a, b, n)) { cout << "該方陣的逆矩陣為: \n"; for (i = 0; i < n; i++) { cout << setw(4); for (j = 0; j < n; j++) { cout << b[i][j] << setw(10); } cout << endl; } } getchar(); return 0; } //----------------------------------------------- //功能: 求矩陣(n*n)的行列式 //入口參數: 矩陣的首地址,矩陣的行數 //返回值: 矩陣的行列式值 //---------------------------------------------- float MatDet(float *p, int n) { int r, c, m; int lop = 0; float result = 0; float mid = 1; if (n != 1) { lop = (n == 2) ? 1 : n; //控制求和循環次數,若為2階,則循環1次,否則為n次 for (m = 0; m < lop; m++) { mid = 1; //順序求和, 主對角線元素相乘之和 for (r = 0, c = m; r < n; r++, c++) { mid = mid * (*(p+r*n+c%n)); } result += mid; } for (m = 0; m < lop; m++) { mid = 1; //逆序相減, 減去次對角線元素乘積 for (r = 0, c = n-1-m+n; r < n; r++, c--) { mid = mid * (*(p+r*n+c%n)); } result -= mid; } } else result = *p; return result; } //---------------------------------------------------------------------------- //功能: 求k*k矩陣中元素A(m, n)的代數余之式 //入口參數: k*k矩陣的首地址,矩陣元素A的下標m,n,矩陣行數k //返回值: k*k矩陣中元素A(m, n)的代數余之式 //---------------------------------------------------------------------------- float Creat_M(float *p, int m, int n, int k) { int len; int i, j; float mid_result = 0; int sign = 1; float *p_creat, *p_mid; len = (k-1)*(k-1); //k階矩陣的代數余之式為k-1階矩陣 p_creat = (float*)calloc(len, sizeof(float)); //分配內存單元 p_mid = p_creat; for (i = 0; i < k; i++) { for (j = 0; j < k; j++) { if (i != m && j != n) //將除第i行和第j列外的所有元素存儲到以p_mid為首地址的內存單元 { *p_mid++ = *(p+i*k+j); } } } sign = (m+n)%2 == 0 ? 1 : -1; //代數余之式前面的正、負號 mid_result = (float)sign*MatDet(p_creat, k-1); free(p_creat); return mid_result; } //----------------------------------------------------- //功能: 打印n*n矩陣 //入口參數: n*n矩陣的首地址,矩陣的行數n //返回值: 無返回值 //----------------------------------------------------- void print(float *p, int n) { int i, j; for (i = 0; i < n; i++) { cout << setw(4); for (j = 0; j < n; j++) { cout << setiosflags(ios::right) << *p++ << setw(10); } cout << endl; } } //------------------------------------------------------------------ //功能: 采用部分主元的高斯消去法求方陣A的逆矩陣B //入口參數: 輸入方陣,輸出方陣,方陣階數 //返回值: true or false //------------------------------------------------------------------- bool Gauss(float A[][N], float B[][N], int n) { int i, j, k; float max, temp; float t[N][N]; //臨時矩陣 //將A矩陣存放在臨時矩陣t[n][n]中 for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { t[i][j] = A[i][j]; } } //初始化B矩陣為單位陣 for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { B[i][j] = (i == j) ? (float)1 : 0; } } for (i = 0; i < n; i++) { //尋找主元 max = t[i][i]; k = i; for (j = i+1; j < n; j++) { if (fabs(t[j][i]) > fabs(max)) { max = t[j][i]; k = j; } } //如果主元所在行不是第i行,進行行交換 if (k != i) { for (j = 0; j < n; j++) { temp = t[i][j]; t[i][j] = t[k][j]; t[k][j] = temp; //B伴隨交換 temp = B[i][j]; B[i][j] = B[k][j]; B[k][j] = temp; } } //判斷主元是否為0, 若是, 則矩陣A不是滿秩矩陣,不存在逆矩陣 if (t[i][i] == 0) { cout << "There is no inverse matrix!"; return false; } //消去A的第i列除去i行以外的各行元素 temp = t[i][i]; for (j = 0; j < n; j++) { t[i][j] = t[i][j] / temp; //主對角線上的元素變為1 B[i][j] = B[i][j] / temp; //伴隨計算 } for (j = 0; j < n; j++) //第0行->第n行 { if (j != i) //不是第i行 { temp = t[j][i]; for (k = 0; k < n; k++) //第j行元素 - i行元素*j列i行元素 { t[j][k] = t[j][k] - t[i][k]*temp; B[j][k] = B[j][k] - B[i][k]*temp; } } } } getchar(); return true; }