配平化學方程式的C++代碼實現
紀念一下我今天寫過了 20171006。
(去年的這個時候我就有了這個大膽的想法， 當時的思路是：字符串處理->暴力搜系數，可是太年輕寫不對，我那會還是個只會模擬的孩子啊，(現在也是))
主要思路：
先做字符串處理，把每個物質的的每種原子數都找出來，
然后利用每種原子的守恆 關於系數 列出方程組 進行求解 （化合價好像不太現實，我化學不好）
先說方程的解法，
解線性方程組當然是要用高斯消元了。
#include<bits/stdc++.h>
using namespace std; double M[105][105]; int N; inline bool Gauss() { for(int k=1;k<=N;k++){ double maxm=-1;int maxi; for(int i=k;i<=N;i++) if(maxm<fabs(M[i][k])) maxm=fabs(M[i][k]),maxi=i; if(fabs(maxm)<1e-7) return false; if(maxi-k) for(int j=1;j<=N+1;j++) swap(M[maxi][j],M[k][j]); double tmp=M[k][k]; for(int j=1;j<=N+1;j++) M[k][j]/=tmp; for(int i=k-1?1:2;i<=N;i++){ if(i==k)continue; double tmp=M[i][k]; for(int j=1;j<=N+1;j++) M[i][j]-=tmp*M[k][j]; } } return true; } int main() { scanf("%d",&N); for(int i=1;i<=N;i++) for(int j=1;j<=N+1;j++) scanf("%lf",&M[i][j]); if(Gauss()) for(int i=1;i<=N;i++) printf("%.2lf\n",M[i][N+1]); else printf("No Solution"); return 0; }
我們在處理字符串的過程中， 直接把對應系數 放到矩陣中 ，
考慮一個問題 ， 化學方程式的系數是可以按比例變化的，所以這個方程組應該是無窮解，
我們的方法是 設其中的一個系數為1 先解方程 ，然后把解出的分數通分 ，就可以得到最簡整數解。
舉個例子
C16H18O9+O2=CO2+H2O 設系數分別是 x1,x2,x3,x4 ,三個原子守恆方程 :
C : x1*16+x2*0-x3*1-x4*0=0
H : x1*18+x2*0-x3*0-x4*2=0
O : x1*9+x2*2-x3*2-x4*1=0
設最后一個也就是x4為1 ，那么方程就變成了
16 0 -1 0
18 0 0 2
9 2 -2 1
現在就可以直接套高斯消元模板了。
對了 還有 ，因為要出分數，我們需要一個分數類模板
struct frac{ //分數類
int a,b; void reduce(){ int x=gcd(a,b); a/=x,b/=x; }; frac operator = (int x){ a=x,b=1; return *this; }; frac operator = (const frac x){ a=x.a,b=x.b; reduce(); return *this; }; frac operator + (const frac x){ return (frac){b*x.a+a*x.b,b*x.b}; }; frac operator - (const frac x){ return (frac){a*x.b-b*x.a,b*x.b}; }; frac operator * (const frac x){ return (frac){a*x.a,b*x.b}; }; frac operator / (const frac x){ return (frac){a*x.b,b*x.a}; }; bool operator < (const frac x){ return a*x.b<b*x.a; }; bool operator == (const frac x){ return a*x.b==b*x.a; }; void print(){ if(b==1)printf("%d\n",a); else printf("%d/%d\n",a,b); }; }; inline frac Abs(frac x){ int p=x.a>0?x.a:-x.a,q=x.b>0?x.b:-x.b; return (frac){p,q}; }
實現過程小問題 :
1.要記好每個物質的名字輸出的時候用
2.元素的名字拿Map這種東西判重，我是開了一個26*26數組手動判重
3.轉類型賦值可能會出現問題，我可是調了好久，心態爆炸
4.把系數放進矩陣的時候有很多小細節 ：比如說存儲的位置、括號的倍數、正負什么的
5.高斯消元里面那個交換操作很關鍵。而且這個矩陣不一定是正方形，注意寫在循環里的長和寬邊界。
下面是完全版代碼
/* Chemical Equation Balancer HiJ1m 2017.10.6 */ #include<bits/stdc++.h>
using namespace std; inline int gcd(int x,int y){ return x%y==0?y:gcd(y,x%y); } inline int lcm(int x,int y){ return x*y/gcd(x,y); } struct frac{ //分數類
int a,b; void reduce(){ int x=gcd(a,b); a/=x,b/=x; }; frac operator = (int x){ a=x,b=1; return *this; }; frac operator = (const frac x){ a=x.a,b=x.b; reduce(); return *this; }; frac operator + (const frac x){ return (frac){b*x.a+a*x.b,b*x.b}; }; frac operator - (const frac x){ return (frac){a*x.b-b*x.a,b*x.b}; }; frac operator * (const frac x){ return (frac){a*x.a,b*x.b}; }; frac operator / (const frac x){ return (frac){a*x.b,b*x.a}; }; bool operator < (const frac x){ return a*x.b<b*x.a; }; bool operator == (const frac x){ return a*x.b==b*x.a; }; void print(){ if(b==1)printf("%d\n",a); else printf("%d/%d\n",a,b); }; }; inline frac Abs(frac x){ int p=x.a>0?x.a:-x.a,q=x.b>0?x.b:-x.b; return (frac){p,q}; } char s[55]; int fun[55][55]; int Map[27][27]; //手動MAP
frac M[55][55]; //求解矩陣
frac ans[55]; //解
int Ans[55]; //整數解
int cnt,c1,c2,flag=1,N,K; //cnt數元素，c1數反應物，c2總數 （未知數的數量）
char mat[55][55]; //存儲物質的名稱
void print(){ printf("%d %d\n",N,K); for(int i=1;i<=K;i++){ for(int j=1;j<=N+1;j++) printf("%d ",M[i][j].a); printf("\n"); } printf("\n"); } inline int getint(int pos){ //讀數
pos++; if(s[pos]>='a'&&s[pos]<='z')pos++; if(s[pos]<'0'||s[pos]>'9')return 1; //沒數就是1
else { int x=0; while(s[pos]>='0'&&s[pos]<='9')x=x*10+s[pos]-'0',pos++; //讀元素后面的數字
return x; } } inline void scan(int l,int r){ //處理物質
c2++; for(int i=0;i<=r-l;i++)mat[c2][i]=s[l+i]; //存下元素的名字
if(flag==1)c1++; //統計一下反應物數量
int tmp=1; //tmp是小括號倍數
for(int i=l;i<=r;i++){ if(s[i]==')')tmp=1; if(s[i]=='('){ int j=i+1;while(s[j]!=')')j++; //找這個括號的范圍
tmp=getint(j); //讀")"右邊的數字
} if(s[i]>='A'&&s[i]<='Z'){ //發現元素
int x=s[i]-'A'+1,y=0; if(s[i+1]>='a'&&s[i]<='z') //看一眼是一個字母的還是兩個的
y=s[i+1]-'a'+1; if(!Map[x][y])Map[x][y]=++cnt; //判重
fun[Map[x][y]][c2]+=flag*getint(i)*tmp; //把這個物質里的這種元素數量放進矩陣里，坐標（map[x][y]，c2）
} } } inline bool Solve(){ //解方程 （矩陣 高cnt，寬c2+1,c2+1列常數全0）
ans[c2]=1; //令最后一個解為1
for(int i=1;i<=cnt;i++){ for(int j=1;j<=c2;j++) M[i][j]=fun[i][j]; } for(int i=1;i<=cnt;i++) M[i][c2].a=-M[i][c2].a; //移到常數 //高斯消元過程
N=c2-1,K=cnt; for(int k=1;k<=N;k++){ frac maxm=(frac){-1,1};int maxi; for(int i=k;i<=K;i++) if(maxm<Abs(M[i][k])) maxm=Abs(M[i][k]),maxi=i; if(maxm==(frac){0,1}) return false; if(maxi!=k) for(int j=1;j<=N+1;j++){ swap(M[k][j],M[maxi][j]); } frac tmp=M[k][k]; for(int j=1;j<=N+1;j++) M[k][j]=M[k][j]/tmp; for(int i=k-1?1:2;i<=K;i++){ if(i==k)continue; frac tmp=M[i][k]; for(int j=1;j<=N+1;j++) M[i][j]=M[i][j]-tmp*M[k][j]; } } return true; } int main() { // printf("Chemical Equation Balancer\n"); // printf("\nEnter the chemical equation:\n");
scanf("%s",s); int lst=0; for(int i=1;i<strlen(s);i++){ if(i==strlen(s)-1)scan(lst,i); if(s[i]=='+'||s[i]=='=')scan(lst,i-1),lst=i+1; if(s[i]=='=')flag=-1; //等號后面的系數變負
} if(Solve()) for(int i=1;i<=c2-1;i++) ans[i]=M[i][N+1]; else printf("No Solution"); int tmp=lcm(ans[1].b,ans[2].b); for(int i=3;i<=c2;i++)tmp=lcm(tmp,ans[i].b); for(int i=1;i<=c2;i++)Ans[i]=ans[i].a*tmp/ans[i].b; //取分母Lcm，把分數變整數
for(int i=1;i<=c2;i++) { if(Ans[i]>1)printf("%d",Ans[i]); for(int j=0;j<strlen(mat[i]);j++) printf("%c",mat[i][j]); if(i==c2)return 0; else if(i==c1)printf("="); else printf("+"); } }
可能會被有的方程式卡住，歡迎Hack，歡迎Debug