1 #include<iostream> 2 #include<cstdio> 3 #include<iomanip> 4 using namespace std; 5 #define e 0.00000001 6 #define maxn 50 7 8 int n;//规模nXn 9 double a[maxn][maxn];//系数矩阵 10 double b[maxn];//b矩阵 11 double m[maxn][maxn];//中间变量矩阵 12 double x[maxn];//最终解 13 int H=1;//扩大H被结算(优化) 14 /* 15 读取数据 16 */ 17 void read(){ 18 cout<<"请输入系数矩阵规模n:= "; 19 cin>>n; 20 cout<<"|-----------------------------\n"; 21 cout<<"|请输入系数矩阵,如:\n"; 22 cout<<"|1.1348 3.8326 1.1651 3.4017\n"; 23 cout<<"|0.5301 1.7875 2.5330 1.5435\n"; 24 cout<<"|3.4129 4.9317 8.7643 1.3142\n"; 25 cout<<"|1.2371 4.9998 10.6721 0.0147\n"; 26 cout<<"|-----------------------------\n"; 27 for(int i=1;i<=n;i++) 28 for(int j=1;j<=n;j++){ 29 cin>>a[i][j]; 30 a[i][j]*=H; 31 } 32 cout<<"|-----------------------------\n"; 33 cout<<"|请输入b矩阵,如:\n"; 34 cout<<"|9.5342 6.3941 18.4231 16.9237\n"; 35 cout<<"|-----------------------------\n"; 36 for(int i=1;i<=n;i++){ 37 cin>>b[i]; 38 b[i]*=H; 39 } 40 } 41 42 /* 43 中间矩阵输出 44 参数:消元次数 45 */ 46 void PrintProc(int cases){ 47 printf("--------第%d次消元结果如下:\n",cases); 48 for(int i=1;i<=n;i++){ 49 for(int j=1;j<=n;j++){ 50 cout<<setw(10)<<a[i][j]<<' '; 51 } 52 cout<<setw(10)<<b[i]<<'\n'; 53 } 54 cout<<"END THIS SHOW-------------\n"; 55 } 56 57 /* 58 显示结果 59 */ 60 void Print(){ 61 cout<<"|-----------------------------\n"; 62 cout<<"|结果为:\n"; 63 for(int i=1;i<=n;i++){ 64 printf("x[%d]= %lf\n",i,x[i]); 65 } 66 cout<<"|-----------------------------\n\n"; 67 } 68 69 /* 70 顺序消元法 71 */ 72 void ShunXuXiaoYuan(){ 73 //消元计算 74 for(int k=1;k<n;k++){ 75 for(int i=k+1;i<=n;i++){ 76 m[i][k]=a[i][k]/a[k][k]; 77 for(int j=k+1;j<=n;j++){ 78 a[i][j]-=m[i][k]*a[k][j]; 79 } 80 } 81 for(int i=k+1;i<=n;i++){ 82 b[i]-=m[i][k]*b[k]; 83 } 84 PrintProc(k);//输出中间计算过程 85 } 86 //回代求解 87 x[n]=b[n]/a[n][n]; 88 for(int i=n-1;i>0;i--){ 89 x[i]=b[i]; 90 for(int j=i+1;j<=n;j++) 91 x[i]-=a[i][j]*x[j]; 92 x[i]/=a[i][i]; 93 } 94 //输出结果 95 Print(); 96 } 97 98 /* 99 列主消元 100 */ 101 void LieZhuXiaoYuan(){ 102 for(int k=1;k<n;k++){ 103 //选主元[这一列的绝对值最大值] 104 double ab_max=-1; 105 int max_ik; 106 for(int i=k;i<=n;i++){ 107 if(abs(a[i][k])>ab_max){ 108 ab_max=abs(a[i][k]); 109 max_ik=i; 110 } 111 } 112 //交换行处理[先判断是否为0矩阵] 113 if(ab_max<e){//0矩阵情况 114 cout<<"det A=0\n"; 115 break; 116 }else if(max_ik!=k){//是否是当前行,不是交换 117 double temp; 118 for(int j=1;j<=n;j++){ 119 temp=a[max_ik][j]; 120 a[max_ik][j]=a[k][j]; 121 a[k][j]=temp; 122 } 123 temp=b[max_ik]; 124 b[max_ik]=b[k]; 125 b[k]=temp; 126 } 127 //消元计算 128 for(int i=k+1;i<=n;i++){ 129 a[i][k]/=a[k][k]; 130 for(int j=k+1;j<=n;j++){ 131 a[i][j]-=a[i][k]*a[k][j]; 132 } 133 b[i]-=a[i][k]*b[k]; 134 } 135 PrintProc(k);//输出中间计算过程 136 if(k<n-1)continue; 137 else{ 138 if(abs(a[n][n])<e){ 139 cout<<"det A=0\n"; 140 break; 141 }else{//回代求解 142 x[n]=b[n]/a[n][n]; 143 for(int i=n-1;i>0;i--){ 144 x[i]=b[i]; 145 for(int j=i+1;j<=n;j++) 146 x[i]-=a[i][j]*x[j]; 147 x[i]/=a[i][i]; 148 } 149 //输出结果 150 Print(); 151 } 152 } 153 } 154 } 155 156 /* 157 主函数 158 */ 159 int main(){ 160 while(1){ 161 read(); 162 LieZhuXiaoYuan(); 163 //ShunXuXiaoYuan(); 164 }return 0; 165 } 166 /* 167 书上高斯消元的例子: 168 1 1 1 169 1 3 -2 170 2 -2 1 171 172 6 1 1 173 */ 174 /* 175 书上列主消元的例子: 176 -0.002 2 2 177 1 0.78125 0 178 3.996 5.5625 4 179 180 0.4 1.3816 7.4178 181 */