化学离子平衡作业偷懒神器
下一次加入沉淀平衡
这一段放到K.dat中，有需要自己添加
1e-14
4
1 CH3COOH 1.74e-5 CH3COO-
2 H2CO3 4.2e-7 HCO3- 5.6e-11 CO3--
2 H2S 1.3e-7 HS- 7.1e-15 S--
1 HCl 1e100 Cl-
3
1 NH3_H2O 1.78e-5 NH4+
3 Fe(OH)3 3.16e-9 Fe(OH)2- 5.01e-10 Fe(OH)-- 1.35e-12 Fe---
1 NaOH 1e100 Na+
以下为main.cpp
#include<bits/stdc++.h>
using namespace std;
const double eps = 1e-14;
double Kacid[110][5], Kbase[110][5];
int Na, Nb;
int Tacid[110], Tbase[110];
double Cacid[110][5], Cbase[110][5];
double CH, COH;
double Kw; // the ionic product of water
string nameAcid[110][5];
string nameBase[110][5];
void readData(){
ifstream in("K.dat");
in >> Kw;
in >> Na;
for (int i=1;i<=Na;i++){
in >> Tacid[i];
for (int j=1;j<=Tacid[i];j++){
in >> nameAcid[i][j];
in >> Kacid[i][j];
}
in >> nameAcid[i][Tacid[i]+1];
}
in >> Nb;
for (int i=1;i<=Nb;i++){
in >> Tbase[i];
for (int j=1;j<=Tbase[i];j++){
in >> nameBase[i][j];
in >> Kbase[i][j];
}
in >> nameBase[i][Tbase[i]+1];
}
in.close();
}
double Solve(double &C1, double &C2, double &C3, double &C4, double K){
if (K > 1e30){
double T = min(C1, C2);
C3 += T;
C4 += T;
C1 -= T;
C2 -= T;
return 0;
}
else if (K < 1e-30){
double T = min(C3, C4);
C3 -= T;
C4 -= T;
C1 += T;
C2 += T;
return 0;
}
double A = K-1;
double B = -K*(C1+C2)-(C3+C4);
double C = K*C1*C2 - C3*C4;
// cout << A << ' ' << B << ' ' << C << endl;
assert(B*B - 4 * A * C + eps >= 0);
double X1 = (-B + sqrt(B*B-4*A*C)) / 2 / A;
double X2 = (-B - sqrt(B*B-4*A*C)) / 2 / A;
if (abs(X1) < eps && abs(X2) < eps) return 0;
// cout << A << ' ' << B << ' ' << C << ' ' << X1 << ' ' << X2 << ' ' << K << endl;
if (X1 < min(C1,C2) + eps && -X1 < min(C3,C4) + eps) {
C1 -= X1; C2 -= X1; C3 += X1; C4 += X1;
return X1;
}
if (X2 < min(C1,C2) + eps && -X2 < min(C3,C4) + eps) {
C1 -= X2; C2 -= X2; C3 += X2; C4 += X2;
return X2;
}
assert(0);
return 0;
}
double Solve(double &C1, double &C2, double &C3, double K){
if (K > 1e30){
double T = C1;
C3 += T;
C1 -= T;
C2 += T;
return 0;
}
else if (K < 1e-30){
double T = min(C3, C2);
C1 += T;
C2 -= T;
C3 -= T;
return 0;
}
double A = 1;
double B = C2+C3+K*C1;
double C = C2*C3-K*C1;
assert(B*B - 4 * A * C + eps >= 0);
double X1 = (-B + sqrt(B*B-4*A*C)) / 2 / A;
double X2 = (-B - sqrt(B*B-4*A*C)) / 2 / A;
// cout << X1 << ' ' << X2 << ' ' << K << endl;
if (abs(X1) < eps && abs(X2) < eps) return 0;
if (X1 < C1 + eps && -X1 < min(C2,C3) + eps) {
C1 -= X1; C2 += X1; C3 += X1;
return X1;
}
if (X2 < C1 + eps && -X2 < min(C2,C3) + eps) {
C1 -= X2; C2 += X2; C3 += X2;
return X2;
}
assert(0);
return 0;
}
double SolveWater(){
double A = 1;
double B = -CH - COH;
double C = CH * COH - Kw;
// A = 1 B = -([H] + [OH]) C = [H][OH]-Kw;
assert(B*B - 4 * A * C >= 0);
double X1 = (-B + sqrt(B*B-4*A*C)) / 2 / A;
double X2 = (-B - sqrt(B*B-4*A*C)) / 2 / A;
if (X1 < min(CH,COH) + eps) {
CH -= X1; COH -= X1;
return X1;
}
if (X2 < min(CH,COH) + eps) {
CH -= X2; COH -= X2;
return X2;
}
assert(0);
return 0;
}
void calc(){
/* Find the reaction that can be do max*/
// int Ka = 1, Kal = 1, Kb = 1, Kbl = 1;
// double val = 0;
for (int i=1;i<=Na;i++){
for (int j=1;j<=Tacid[i];j++){
for (int k=1;k<=Nb;k++){
for (int l=1;l<=Tbase[k];l++){
double K = Kacid[i][j] * Kbase[k][l] / Kw;
// cout << Cacid[i][j] << ' ' << Cbase[k][l] << endl;
Solve(Cacid[i][j], Cbase[k][l], Cacid[i][j+1], Cbase[k][l+1], K);
}
}
}
}
for (int i=1;i<=Na;i++){
for (int j=1;j<=Tacid[i];j++){
double K = Kacid[i][j];
Solve(Cacid[i][j], Cacid[i][j+1], CH, K);
}
}
for (int i=1;i<=Nb;i++){
for (int j=1;j<=Tbase[i];j++){
double K = Kbase[i][j];
Solve(Cbase[i][j], Cbase[i][j+1], COH, K);
}
}
SolveWater();
}
/*
HA == A- + H+ Ka
BOH == B+ + OH- Kb
H+ + OH- == H2O 1/Kw
HA + BOH == AB + H2O
[A-][B+] [A-][H+][B+][OH-]
-------- = ----------------- = Ka * Kb / Kw
[AB] [AB][H2O]
Sol.
HA + BOH == A- + B+ + H2O K
C1 C2 C3 C4
(C3+x * C4+x) / (C1-x * C2-x) == K
K(C1C2 - (C1+C2)x + x^2) == C3C4 + (C3+C4)x + x^2
--> A = (K-1) B = -K*(C1+C2)-(C3+C4) C = KC1C2 - C3C4
--------------------------------------------------------
HA == A- + H+ Ka
C1 C2 C3
([A-]+x)([H+]+x)
------------------ = K
[HA]-x
--> A = 1 B = C2+C3+K*C1 C = C2*C3-K*C1
--------------------------------------------------------
[H] * [OH] = Kw;
([H]-x) * ([OH]-x) = Kw;
A = 1 B = -([H] + [OH]) C = [H][OH]-Kw;
*/
int main(){
readData();
int n;
int T = 10000;
Cacid[1][1] = 0.1;
Cacid[1][2] = 0.1;
Cbase[3][2] = 0.2;
while (T--) calc();
cout << "[H] = " << CH << endl;
cout << "[OH] = " << COH << endl;
cout << "pH = " << -log10(CH) << endl;
cout << "Acid:" << endl;
for (int i=1;i<=Na;i++){
for (int j=1;j<=Tacid[i]+1;j++)
cout << "[" << nameAcid[i][j] << "] = " << Cacid[i][j] << ' ';
cout << endl;
}
cout << "Base" << endl;
for (int i=1;i<=Nb;i++){
for (int j=1;j<=Tbase[i]+1;j++)
cout << "[" << nameBase[i][j] << "] = " << Cbase[i][j] << ' ';
cout << endl;
}
}
示例：
0.1mol/L CH3COOH 与 0.1mol/L CH3COONa混，最后pH值
Cacid[1][1] = 0.1;
Cacid[1][2] = 0.1;
Cbase[3][2] = 0.1;
Solution:
[H] = 1.73939e-05
[OH] = 5.74913e-10
pH = 4.7596
Acid:
[CH3COOH] = 0.0999826 [CH3COO-] = 0.100017
[H2CO3] = 0 [HCO3-] = 0 [CO3--] = 0
[H2S] = 0 [HS-] = 0 [S--] = 0
[HCl] = 0 [Cl-] = 0
Base
[NH3_H2O] = 0 [NH4+] = 0
[Fe(OH)3] = 0 [Fe(OH)2-] = 0 [Fe(OH)--] = 0 [Fe---] = 0
[NaOH] = 0 [Na+] = 0.1
示例：
0.01mol/L CH3COOH 与 0.40mol/L CH3COONa混，最后pH值
Cacid[1][1] = 0.01;
Cacid[1][2] = 0.4;
Cbase[3][2] = 0.4;
Solution:
[H] = 4.34982e-07
[OH] = 2.29895e-08
pH = 6.36153
Acid:
[CH3COOH] = 0.00999959 [CH3COO-] = 0.4
[H2CO3] = 0 [HCO3-] = 0 [CO3--] = 0
[H2S] = 0 [HS-] = 0 [S--] = 0
[HCl] = 0 [Cl-] = 0
Base
[NH3_H2O] = 0 [NH4+] = 0
[Fe(OH)3] = 0 [Fe(OH)2-] = 0 [Fe(OH)--] = 0 [Fe---] = 0
[NaOH] = 0 [Na+] = 0.4
Copyright (c) dgklr : Creative Commons - BY - NC
