Clean up of some of the test code

This commit is contained in:
Ray Speth 2012-10-24 15:44:33 +00:00
parent 1d19803885
commit fcb46460cc
4 changed files with 25 additions and 91 deletions

View file

@ -1,11 +1,9 @@
#include <cstdio> #include <cstdio>
#include <cstdlib> #include <cstdlib>
#include <cmath> #include <cmath>
using namespace std; using namespace std;
/* /*
* Values of A_J/R : tabular form * Values of A_J/R : tabular form
* units sqrt(kg/gmol), * units sqrt(kg/gmol),
@ -95,24 +93,13 @@ double Cphi(double temp, int ifunc)
return retn; return retn;
} }
double calc(double temp, double Iionic) void calc(double temp, double Iionic)
{ {
/* /*
* Gas Constant in J gmol-1 K-1 * Gas Constant in J gmol-1 K-1
*/ */
double GasConst = 8.314472; double GasConst = 8.314472;
double Aphi = 0.0;
if (temp == 323.15) {
Aphi = 0.4102995331359;
} else if (temp == 473.15) {
Aphi = 0.622777;
} else {
printf("ERROR: unknown temp\n");
exit(-1);
}
//printf(" Aphi = %g\n", Aphi);
/* /*
* Calculate A_H in J gmol-1 sqrt(kg/gmol) * Calculate A_H in J gmol-1 sqrt(kg/gmol)
*/ */
@ -174,13 +161,10 @@ double calc(double temp, double Iionic)
double xo = 1.0 / (molecWeight/1000. * 2 * m + 1.0); double xo = 1.0 / (molecWeight/1000. * 2 * m + 1.0);
printf(" no = %g\n", xo); printf(" no = %g\n", xo);
return phiJ;
} }
int main() int main()
{ {
printf("Standalone test of the apparent relative molal enthalpy, phiL:\n"); printf("Standalone test of the apparent relative molal enthalpy, phiL:\n");
printf(" (Check against simple formula in Silvester&Pitzer, J. Phys. Chem. 81, 1822 (1977)\n"); printf(" (Check against simple formula in Silvester&Pitzer, J. Phys. Chem. 81, 1822 (1977)\n");
@ -188,10 +172,10 @@ int main()
double Iionic = 6.146; double Iionic = 6.146;
printf("Ionic Strength = %g\n", Iionic); printf("Ionic Strength = %g\n", Iionic);
double res = calc(273.15 + 50., Iionic); calc(273.15 + 50., Iionic);
printf("T = 200C\n"); printf("T = 200C\n");
printf("Ionic Strength = %g\n", Iionic); printf("Ionic Strength = %g\n", Iionic);
res = calc(273.15 + 200., Iionic); calc(273.15 + 200., Iionic);
return 0; return 0;
} }

View file

@ -1,60 +1,45 @@
#include <cstdio> #include <cstdio>
#include <cstdlib> #include <cstdlib>
#include <cmath> #include <cmath>
using namespace std; using namespace std;
double Beta0(double temp)
double Beta0(double temp, int ifunc)
{ {
double q1 = 0.0765; double q1 = 0.0765;
double q2 = -777.03; double q2 = -777.03;
double q3 = -4.4706; double q3 = -4.4706;
double q4 = 0.008946; double q4 = 0.008946;
double q5 = -3.3158E-6; double q5 = -3.3158E-6;
double retn;
double tref = 298.15; double tref = 298.15;
if (ifunc == 0) { return q1 + q2 * (1.0/temp - 1.0/tref)
retn = q1 + q2 * (1.0/temp - 1.0/tref) + q3 * (log(temp/tref)) + q4 * (temp - tref)
+ q3 * (log(temp/tref)) + q4 * (temp - tref) + q5 * (temp * temp - tref * tref);
+ q5 * (temp * temp - tref * tref);
}
return retn;
} }
double Beta1(double temp, int ifunc) double Beta1(double temp)
{ {
double q6 = 0.2664; double q6 = 0.2664;
double q9 = 6.1608E-5; double q9 = 6.1608E-5;
double q10 = 1.0715E-6; double q10 = 1.0715E-6;
double retn;
double tref = 298.15; double tref = 298.15;
if (ifunc == 0) { return q6 + q9 * (temp - tref)
retn = q6 + q9 * (temp - tref) + q10 * (temp * temp - tref * tref);
+ q10 * (temp * temp - tref * tref);
}
return retn;
} }
double Cphi(double temp, int ifunc) double Cphi(double temp)
{ {
double q11 = 0.00127; double q11 = 0.00127;
double q12 = 33.317; double q12 = 33.317;
double q13 = 0.09421; double q13 = 0.09421;
double q14 = -4.655E-5; double q14 = -4.655E-5;
double retn;
double tref = 298.15; double tref = 298.15;
if (ifunc == 0) { return q11 + q12 * (1.0/temp - 1.0/tref)
retn = q11 + q12 * (1.0/temp - 1.0/tref) + q13 * (log(temp/tref)) + q14 * (temp - tref);
+ q13 * (log(temp/tref)) + q14 * (temp - tref);
}
return retn;
} }
double calc(double temp, double Iionic) void calc(double temp, double Iionic)
{ {
double Aphi = 0.0; double Aphi = 0.0;
if (temp == 323.15) { if (temp == 323.15) {
Aphi = 0.4102995331359; Aphi = 0.4102995331359;
@ -67,13 +52,13 @@ double calc(double temp, double Iionic)
printf(" Aphi = %g\n", Aphi); printf(" Aphi = %g\n", Aphi);
double beta0 = Beta0(temp, 0); double beta0 = Beta0(temp);
printf(" beta0 = %g\n", beta0); printf(" beta0 = %g\n", beta0);
double beta1 = Beta1(temp, 0); double beta1 = Beta1(temp);
printf(" beta1 = %g\n", beta1); printf(" beta1 = %g\n", beta1);
double cphi = Cphi(temp, 0); double cphi = Cphi(temp);
printf(" Cphi = %g\n", cphi); printf(" Cphi = %g\n", cphi);
double vm = 1.0; double vm = 1.0;
@ -95,7 +80,6 @@ double calc(double temp, double Iionic)
double os = osm1 + 1.0; double os = osm1 + 1.0;
double a2 = alpha * alpha; double a2 = alpha * alpha;
printf("osmotic coeff = %20.13g\n", os); printf("osmotic coeff = %20.13g\n", os);
double lnmeanAct = - zm * zx * Aphi * double lnmeanAct = - zm * zx * Aphi *
@ -138,20 +122,18 @@ double calc(double temp, double Iionic)
double tmp = diff / (RT); double tmp = diff / (RT);
double actCoefWater = exp(tmp) / xo; double actCoefWater = exp(tmp) / xo;
printf("actCoefWater = %g\n", actCoefWater); printf("actCoefWater = %g\n", actCoefWater);
return gex;
} }
int main() int main()
{ {
printf("standalone test of Gibbs excess free energy:\n"); printf("standalone test of Gibbs excess free energy:\n");
printf("T = 50C\n"); printf("T = 50C\n");
double Iionic = 6.146; double Iionic = 6.146;
printf("Ionic Strength = %g\n", Iionic); printf("Ionic Strength = %g\n", Iionic);
double res = calc(273.15 + 50., Iionic); calc(273.15 + 50., Iionic);
printf("T = 200C\n"); printf("T = 200C\n");
res = calc(273.15 + 200., Iionic); calc(273.15 + 200., Iionic);
return 0; return 0;
} }

View file

@ -1,11 +1,9 @@
#include <cstdio> #include <cstdio>
#include <cstdlib> #include <cstdlib>
#include <cmath> #include <cmath>
using namespace std; using namespace std;
/* /*
* Values of A_L/RT : tabular form * Values of A_L/RT : tabular form
* units sqrt(kg/gmol), * units sqrt(kg/gmol),
@ -86,24 +84,13 @@ double Cphi(double temp, int ifunc)
return retn; return retn;
} }
double calc(double temp, double Iionic) void calc(double temp, double Iionic)
{ {
/* /*
* Gas Constant in J gmol-1 K-1 * Gas Constant in J gmol-1 K-1
*/ */
double GasConst = 8.314472; double GasConst = 8.314472;
double Aphi = 0.0;
if (temp == 323.15) {
Aphi = 0.4102995331359;
} else if (temp == 473.15) {
Aphi = 0.622777;
} else {
printf("ERROR: unknown temp\n");
exit(-1);
}
//printf(" Aphi = %g\n", Aphi);
/* /*
* Calculate A_H in J gmol-1 sqrt(kg/gmol) * Calculate A_H in J gmol-1 sqrt(kg/gmol)
*/ */
@ -151,13 +138,10 @@ double calc(double temp, double Iionic)
double xo = 1.0 / (molecWeight/1000. * 2 * m + 1.0); double xo = 1.0 / (molecWeight/1000. * 2 * m + 1.0);
printf(" no = %g\n", xo); printf(" no = %g\n", xo);
return phiL;
} }
int main() int main()
{ {
printf("Standalone test of the apparent relative molal enthalpy, phiL:\n"); printf("Standalone test of the apparent relative molal enthalpy, phiL:\n");
printf(" (Check against simple formula in Silvester&Pitzer, J. Phys. Chem. 81, 1822 (1977)\n"); printf(" (Check against simple formula in Silvester&Pitzer, J. Phys. Chem. 81, 1822 (1977)\n");
@ -165,10 +149,10 @@ int main()
double Iionic = 6.146; double Iionic = 6.146;
printf("Ionic Strength = %g\n", Iionic); printf("Ionic Strength = %g\n", Iionic);
double res = calc(273.15 + 50., Iionic); calc(273.15 + 50., Iionic);
printf("T = 200C\n"); printf("T = 200C\n");
printf("Ionic Strength = %g\n", Iionic); printf("Ionic Strength = %g\n", Iionic);
res = calc(273.15 + 200., Iionic); calc(273.15 + 200., Iionic);
return 0; return 0;
} }

View file

@ -1,12 +1,9 @@
#include <cstdio> #include <cstdio>
#include <cstdlib> #include <cstdlib>
#include <cmath> #include <cmath>
using namespace std; using namespace std;
/* /*
* Values of A_V : tabular form * Values of A_V : tabular form
* units sqrt(kg/gmol)cm3/gmol * units sqrt(kg/gmol)cm3/gmol
@ -102,24 +99,13 @@ double Cphi(double temp, int ifunc)
return retn; return retn;
} }
double calc(double temp, double Iionic) void calc(double temp, double Iionic)
{ {
/* /*
* Gas Constant in J gmol-1 K-1 * Gas Constant in J gmol-1 K-1
*/ */
double GasConst = 8.314472; double GasConst = 8.314472;
double Aphi = 0.0;
if (temp == 323.15) {
Aphi = 0.4102995331359;
} else if (temp == 473.15) {
Aphi = 0.622777;
} else {
printf("ERROR: unknown temp\n");
exit(-1);
}
/* /*
* Calculate A_V in sqrt(kg/gmol)cm3/gmol * Calculate A_V in sqrt(kg/gmol)cm3/gmol
*/ */
@ -162,12 +148,10 @@ double calc(double temp, double Iionic)
//double RT = GasConst * temp * 1.0E-3; //double RT = GasConst * temp * 1.0E-3;
double xo = 1.0 / (molecWeight/1000. * 2 * m + 1.0); double xo = 1.0 / (molecWeight/1000. * 2 * m + 1.0);
printf(" no = %g\n", xo); printf(" no = %g\n", xo);
return phiV;
} }
int main() int main()
{ {
printf("Standalone test of the apparent relative molal excess volume, phiV:\n"); printf("Standalone test of the apparent relative molal excess volume, phiV:\n");
printf(" (Check against simple formula in \n"); printf(" (Check against simple formula in \n");
printf(" Activity Coefficients in Eletrolyte Solutions, 2nd Ed K. S. Pitzer, " printf(" Activity Coefficients in Eletrolyte Solutions, 2nd Ed K. S. Pitzer, "
@ -176,10 +160,10 @@ int main()
double Iionic = 6.146; double Iionic = 6.146;
printf("Ionic Strength = %g\n", Iionic); printf("Ionic Strength = %g\n", Iionic);
double res = calc(273.15 + 50., Iionic); calc(273.15 + 50., Iionic);
printf("T = 200C\n"); printf("T = 200C\n");
printf("Ionic Strength = %g\n", Iionic); printf("Ionic Strength = %g\n", Iionic);
res = calc(273.15 + 200., Iionic); calc(273.15 + 200., Iionic);
return 0; return 0;
} }