这模版敲了我俩个小时+写注释,参考自kuangbin!
两百行的大模拟,累死了QAQ
下面附上模版!
1 #include <bits/stdc++.h> 2 using namespace std; 3 const int maxn=50; 4 typedef long long ll; 5 int a[maxn][maxn];///增广矩阵 6 int x[maxn];///解集 7 bool free_x[maxn];///标记是否是不确定的变元 8 /* 9 void debug(void) 10 { 11 int i,j; 12 for(i=0;i<equ;i++) 13 { 14 for(j=0;j<var+1;j++) 15 { 16 cout<<a[i][j]<<" "; 17 } 18 cout<<endl; 19 } 20 cout<<endl; 21 } 22 */ 23 inline int read() 24 { 25 int x=0,f=1; 26 char ch=getchar(); 27 while(ch<'0'||ch>'9') 28 { 29 if(ch=='-') 30 f=-1; 31 ch=getchar(); 32 } 33 while(ch>='0'&&ch<='9') 34 { 35 x=x*10+ch-'0'; 36 ch=getchar(); 37 } 38 return x*f; 39 } 40 inline void write(int x) 41 { 42 if(x<0) 43 { 44 putchar('-'); 45 x=-x; 46 } 47 if(x>9) 48 write(x/10); 49 putchar(x%10+'0'); 50 } 51 inline int gcd(int a,int b)///最大公约数 52 { 53 return b==0?a:gcd(b,a%b); 54 } 55 inline int lcm(int a,int b)///最小公倍数 56 { 57 return a/gcd(a,b)*b;///先除后乘防溢出 58 } 59 ///高斯消元法解方程组【-2表示有浮点型解,无整数解】 60 ///【-1表示无解,0表示唯一解,大于0表示无穷解,并返回自由元的个数】 61 ///【有equ方程,var个变元,增广矩阵行数为equ,分别为0->equ-1,列数为var+1,分别为0->var】 62 int Gauss(int equ,int var) 63 { 64 int i,j,k; 65 int max_r;///当前这列绝对值最大的行 66 int col;///当前处理的列 67 int ta,tb; 68 int LCM; 69 int temp; 70 int free_x_num; 71 int free_index; 72 for(i=0;i<=var;i++) 73 { 74 x[i]=0; 75 free_x[i]=true; 76 } 77 ///转化为阶梯型矩阵 78 col=0;///处理当前列 79 for(k=0;k<equ&&col<var;k++,col++)///枚举当前处理的行 80 { 81 ///找到该col列元素绝对值最大的那一行与第k行交换(为了在除法时减小误差) 82 max_r=k; 83 for(i=k+1;i<equ;i++) 84 { 85 if(abs(a[i][col])>abs(a[max_r][col])) 86 max_r=i; 87 } 88 if(max_r!=k) 89 { 90 ///与第k行交换 91 for(j=k;j<var+1;j++) 92 { 93 swap(a[k][j],a[max_r][j]); 94 } 95 } 96 if(a[k][col]==0) 97 { 98 ///说明该col列第k行以下全是0了,则处理当前行的下一列 99 k--; 100 continue; 101 } 102 for(i=k+1;i<equ;i++) 103 { 104 ///枚举要删除的行 105 if(a[i][col]!=0) 106 { 107 LCM=lcm(abs(a[i][col]),abs(a[k][col])); 108 ta=LCM/abs(a[i][col]); 109 tb=LCM/abs(a[k][col]); 110 if(a[i][col]*a[k][col]<0) 111 tb=-tb;///异号情况是相加 112 for(j=col;j<var+1;j++) 113 { 114 a[i][j]=a[i][j]*ta-a[k][j]*tb; 115 } 116 } 117 } 118 } 119 ///debug(); 120 ///无解的情况:化简的增广矩阵中存在(0,0,......a)这样的行(a!=0) 121 for(i=k;i<equ;i++) 122 { 123 if(a[i][col]!=0) 124 return -1; 125 } 126 ///对于无穷解来说,如果要判断哪些是自由元,那么初等变换中的交换就会影响,则要记录交换 127 ///无穷解的情况:在var*(var+1)的增广矩阵中,出现(0,0,......0)这样的行,即说明没有形成严格的上三角阵 128 ///且目前出现的行数即为自由元的个数 129 if(k<var) 130 { 131 ///首先,自由元有var-k个,即不确定的变元至少有var-k个 132 for(i=k-1;i>=0;i--) 133 { 134 ///第i行一定不会是(0,0,......0)的情况,因为这样的行是在第k行到第equ行 135 ///同样,第i行一定不会是(0,0,......a),(a!=0)这样的情况无解 136 free_x_num=0;///用于判断该行中的不确定的变元的个数,如果超过1就无法求解它们仍为不确定的变元 137 for(j=0;j<var;j++) 138 { 139 if(a[i][j]!=0&&free_x[j]!=0) 140 { 141 free_x_num++; 142 free_index=j; 143 } 144 } 145 if(free_x_num>1) 146 continue;///无法求出确定的变元 147 ///说明只有一个不确定的变元free_index,那么可以求解出该变元,且该变元是确定的, 148 temp=a[i][var]; 149 for(j=0;j<var;j++) 150 { 151 if(a[i][j]!=0&&j!=free_index) 152 temp-=a[i][j]*x[j]; 153 } 154 x[free_index]=temp/a[i][free_index];///求出该变元 155 free_x[free_index]=0;///该变元是确定的 156 } 157 return var-k; 158 } 159 ///唯一解的情况:在var*(var+1)的增广矩阵中形成严格的上三角形 160 ///求出Xn-1,Xn-2,......X0 161 for(i=var-1;i>=0;i--) 162 { 163 temp=a[i][var]; 164 for(j=i+1;j<var;j++) 165 { 166 if(a[i][j]!=0) 167 temp-=a[i][j]*x[j]; 168 } 169 if(temp%a[i][i]!=0) 170 return -2;///说明有浮点数解,但无整数解 171 x[i]=temp/a[i][i]; 172 } 173 return 0; 174 } 175 int main(void) 176 { 177 ///freopen("in.txt","r",stdin); 178 ///freopen("out.txt","w",stdout); 179 int i,j; 180 int equ,var; 181 while(equ=read(),var=read()) 182 { 183 memset(a,0,sizeof(a)); 184 for(i=0;i<equ;i++) 185 { 186 for(j=0;j<var+1;j++) 187 { 188 a[i][j]=read(); 189 } 190 } 191 ///debug(); 192 int free_num=Gauss(equ,var); 193 if(free_num==-1) 194 printf("No solution ");///无解 195 else if(free_num==-2) 196 printf("There are floating point numbers, no integer solutions ");///有浮点数解,无整数解 197 else if(free_num>0) 198 { 199 printf("The number of variables of infinite solution is:%d ",free_num);///无穷多解,自由变元个数为 200 for(i=0;i<var;i++) 201 { 202 if(free_x[i]!=0) 203 printf("x%d is not sure ",i+1); 204 else 205 printf("x%d:%d ",i+1,x[i]); 206 } 207 } 208 else 209 { 210 for(i=0;i<var;i++) 211 printf("x%d:%d ",i+1,x[i]); 212 } 213 printf(" "); 214 } 215 return 0; 216 }