diff --git a/MaterialsLib/Adsorption/Adsorption.cpp b/MaterialsLib/Adsorption/Adsorption.cpp
index bba67ad6c34b06b25b7ce9868647b389180774fe..0871724a840050d3d54ee36ba3b7d0e8a03b62a6 100644
--- a/MaterialsLib/Adsorption/Adsorption.cpp
+++ b/MaterialsLib/Adsorption/Adsorption.cpp
@@ -11,13 +11,13 @@
#include "Adsorption.h"
namespace {
- const double k_rate = 6.0e-3; //to be specified
+ const double k_rate = 6.0e-3; //to be specified
- template <typename T>
- T square(const T& v)
- {
- return v * v;
- }
+ template <typename T>
+ T square(const T& v)
+ {
+ return v * v;
+ }
}
namespace Adsorption
@@ -26,168 +26,168 @@ namespace Adsorption
//Saturation pressure for water used in Nunez
double Adsorption::get_equilibrium_vapour_pressure(const double T_Ads)
{
- //critical T and p
- const double Tc = 647.3; //K
- const double pc = 221.2e5; //Pa
- //dimensionless T
- const double Tr = T_Ads/Tc;
- const double theta = 1. - Tr;
- //empirical constants
- const double c[] = {-7.69123,-26.08023,-168.17065,64.23285,-118.96462,4.16717,20.97506,1.0e9,6.0};
- const double K[] = {c[0]*theta + c[1]*pow(theta,2) + c[2]*pow(theta,3) + c[3]*pow(theta,4) + c[4]*pow(theta,5),
- 1. + c[5]*theta + c[6]*pow(theta,2)};
-
- const double exponent = K[0]/(K[1]*Tr) - theta/(c[7]*pow(theta,2) + c[8]);
- return pc * exp(exponent); //in Pa
+ //critical T and p
+ const double Tc = 647.3; //K
+ const double pc = 221.2e5; //Pa
+ //dimensionless T
+ const double Tr = T_Ads/Tc;
+ const double theta = 1. - Tr;
+ //empirical constants
+ const double c[] = {-7.69123,-26.08023,-168.17065,64.23285,-118.96462,4.16717,20.97506,1.0e9,6.0};
+ const double K[] = {c[0]*theta + c[1]*pow(theta,2) + c[2]*pow(theta,3) + c[3]*pow(theta,4) + c[4]*pow(theta,5),
+ 1. + c[5]*theta + c[6]*pow(theta,2)};
+
+ const double exponent = K[0]/(K[1]*Tr) - theta/(c[7]*pow(theta,2) + c[8]);
+ return pc * exp(exponent); //in Pa
}
//Evaporation enthalpy of water from Nunez
double Adsorption::get_evaporation_enthalpy(double T_Ads) //in kJ/kg
{
- T_Ads -= 273.15;
- if (T_Ads <= 10.){
- const double c[] = {2.50052e3,-2.1068,-3.57500e-1,1.905843e-1,-5.11041e-2,7.52511e-3,-6.14313e-4,2.59674e-5,-4.421e-7};
- double hv = 0.;
- for (size_t i=0; i< sizeof(c)/sizeof(c[0]);i++)
- hv += c[i] * pow(T_Ads,i);
- return hv;
- } else if (T_Ads <= 300.){
- const double c[] = {2.50043e3,-2.35209,1.91685e-4,-1.94824e-5,2.89539e-7,-3.51199e-9,2.06926e-11,-6.4067e-14,8.518e-17,1.558e-20,-1.122e-22};
- double hv = 0.;
- for (size_t i=0; i< sizeof(c)/sizeof(c[0]);i++)
- hv += c[i] * pow(T_Ads,i);
- return hv;
- } else {
- const double c[] = {2.99866e3,-3.1837e-3,-1.566964e1,-2.514e-6,2.045933e-2,1.0389e-8};
- return ((c[0] + c[2]*T_Ads + c[4]*pow(T_Ads,2))/(1. + c[1]*T_Ads + c[3]*pow(T_Ads,2) + c[5]*pow(T_Ads,3)));
- }
+ T_Ads -= 273.15;
+ if (T_Ads <= 10.){
+ const double c[] = {2.50052e3,-2.1068,-3.57500e-1,1.905843e-1,-5.11041e-2,7.52511e-3,-6.14313e-4,2.59674e-5,-4.421e-7};
+ double hv = 0.;
+ for (size_t i=0; i< sizeof(c)/sizeof(c[0]);i++)
+ hv += c[i] * pow(T_Ads,i);
+ return hv;
+ } else if (T_Ads <= 300.){
+ const double c[] = {2.50043e3,-2.35209,1.91685e-4,-1.94824e-5,2.89539e-7,-3.51199e-9,2.06926e-11,-6.4067e-14,8.518e-17,1.558e-20,-1.122e-22};
+ double hv = 0.;
+ for (size_t i=0; i< sizeof(c)/sizeof(c[0]);i++)
+ hv += c[i] * pow(T_Ads,i);
+ return hv;
+ } else {
+ const double c[] = {2.99866e3,-3.1837e-3,-1.566964e1,-2.514e-6,2.045933e-2,1.0389e-8};
+ return ((c[0] + c[2]*T_Ads + c[4]*pow(T_Ads,2))/(1. + c[1]*T_Ads + c[3]*pow(T_Ads,2) + c[5]*pow(T_Ads,3)));
+ }
}
//evaluate specific heat capacity of adsorbate follwing Nunez
double Adsorption::get_specific_heat_capacity(const double T_Ads)
{
- const double c[] = {4.224,-3.716e-3,9.351e-5,-7.1786e-7,-9.1266e-9,2.69247e-10,-2.773104e-12,1.553177e-14,-4.982795e-17,8.578e-20,-6.12423e-23};
- double cp = 0.;
- for (unsigned i=0; i< sizeof(c)/sizeof(c[0]);i++)
- cp += c[i] * pow(T_Ads,i);
- return cp; //kJ/(kg*K)
+ const double c[] = {4.224,-3.716e-3,9.351e-5,-7.1786e-7,-9.1266e-9,2.69247e-10,-2.773104e-12,1.553177e-14,-4.982795e-17,8.578e-20,-6.12423e-23};
+ double cp = 0.;
+ for (unsigned i=0; i< sizeof(c)/sizeof(c[0]);i++)
+ cp += c[i] * pow(T_Ads,i);
+ return cp; //kJ/(kg*K)
}
double Adsorption::get_molar_fraction(double xm, double M_this, double M_other)
{
- return M_other*xm/(M_other*xm + M_this*(1.0-xm));
+ return M_other*xm/(M_other*xm + M_this*(1.0-xm));
}
double Adsorption::get_mass_fraction(double xn, double M_this, double M_other)
{
- return M_this*xn/(M_this*xn + M_other*(1.0-xn));
+ return M_this*xn/(M_this*xn + M_other*(1.0-xn));
}
double Adsorption::d_molar_fraction(double xm, double M_this, double M_other)
{
- return M_other * M_this
- / square(M_other * xm + M_this * (1.0 - xm));
+ return M_other * M_this
+ / square(M_other * xm + M_this * (1.0 - xm));
}
double Adsorption::get_reaction_rate(const double p_Ads, const double T_Ads,
- const double M_Ads, const double loading) const
+ const double M_Ads, const double loading) const
{
- // const double k_rate = 3.0e-3; //to be specified
+ // const double k_rate = 3.0e-3; //to be specified
- const double A = get_potential(p_Ads, T_Ads, M_Ads);
- double C_eq = get_adsorbate_density(T_Ads) * characteristic_curve(A);
- if (C_eq < 0.0) C_eq = 0.0;
+ const double A = get_potential(p_Ads, T_Ads, M_Ads);
+ double C_eq = get_adsorbate_density(T_Ads) * characteristic_curve(A);
+ if (C_eq < 0.0) C_eq = 0.0;
- // return 0.0; // TODO [CL] for testing only
+ // return 0.0; // TODO [CL] for testing only
- return k_rate * (C_eq - loading); //scaled with mass fraction
- // this the rate in terms of loading!
+ return k_rate * (C_eq - loading); //scaled with mass fraction
+ // this the rate in terms of loading!
}
void Adsorption::get_d_reaction_rate(const double p_Ads, const double T_Ads,
- const double M_Ads, const double /*loading*/,
- std::array<double, 3> &dqdr) const
+ const double M_Ads, const double /*loading*/,
+ std::array<double, 3> &dqdr) const
{
- const double A = get_potential(p_Ads, T_Ads, M_Ads);
- const double p_S = get_equilibrium_vapour_pressure(T_Ads);
- const double dAdT = GAS_CONST * log(p_S/p_Ads) / (M_Ads*1.e3);
- const double dAdp = - GAS_CONST * T_Ads / M_Ads / p_Ads;
-
- const double W = characteristic_curve(A);
- const double dWdA = d_characteristic_curve(A);
-
- const double rho_Ads = get_adsorbate_density(T_Ads);
- const double drhodT = - rho_Ads * get_alphaT(T_Ads);
-
- dqdr = std::array<double, 3>{{
- rho_Ads*dWdA*dAdp,
- drhodT*W + rho_Ads*dWdA*dAdT,
- -k_rate
- }};
+ const double A = get_potential(p_Ads, T_Ads, M_Ads);
+ const double p_S = get_equilibrium_vapour_pressure(T_Ads);
+ const double dAdT = GAS_CONST * log(p_S/p_Ads) / (M_Ads*1.e3);
+ const double dAdp = - GAS_CONST * T_Ads / M_Ads / p_Ads;
+
+ const double W = characteristic_curve(A);
+ const double dWdA = d_characteristic_curve(A);
+
+ const double rho_Ads = get_adsorbate_density(T_Ads);
+ const double drhodT = - rho_Ads * get_alphaT(T_Ads);
+
+ dqdr = std::array<double, 3>{{
+ rho_Ads*dWdA*dAdp,
+ drhodT*W + rho_Ads*dWdA*dAdT,
+ -k_rate
+ }};
}
//Evaluate adsorbtion potential A
double Adsorption::get_potential(const double p_Ads, double T_Ads, const double M_Ads) const
{
- double A = GAS_CONST * T_Ads * log(get_equilibrium_vapour_pressure(T_Ads)/p_Ads) / (M_Ads*1.e3); //in kJ/kg = J/g
- if (A < 0.0) {
- // vapour partial pressure > saturation pressure
- // A = 0.0; // TODO [CL] debug output
- }
- return A;
+ double A = GAS_CONST * T_Ads * log(get_equilibrium_vapour_pressure(T_Ads)/p_Ads) / (M_Ads*1.e3); //in kJ/kg = J/g
+ if (A < 0.0) {
+ // vapour partial pressure > saturation pressure
+ // A = 0.0; // TODO [CL] debug output
+ }
+ return A;
}
double Adsorption::get_loading(const double rho_curr, const double rho_dry)
{
- return rho_curr / rho_dry - 1.0;
+ return rho_curr / rho_dry - 1.0;
}
//Calculate sorption entropy
double Adsorption::get_entropy(const double T_Ads, const double A) const
{
- const double epsilon = 1.0e-8;
-
- //* // This change will also change simulation results.
- const double W_p = characteristic_curve(A+epsilon);
- const double W_m = characteristic_curve(A-epsilon);
- const double dAdlnW = 2.0*epsilon/(log(W_p/W_m));
- // */
-
- if (W_p <= 0.0 || W_m <= 0.0)
- {
- ERR("characteristic curve in negative region (W-, W+): %g, %g", W_m, W_p);
- return 0.0;
- }
-
- // const double dAdlnW = 2.0*epsilon/(log(characteristic_curve(A+epsilon)) - log(characteristic_curve(A-epsilon)));
- return dAdlnW * get_alphaT(T_Ads);
+ const double epsilon = 1.0e-8;
+
+ //* // This change will also change simulation results.
+ const double W_p = characteristic_curve(A+epsilon);
+ const double W_m = characteristic_curve(A-epsilon);
+ const double dAdlnW = 2.0*epsilon/(log(W_p/W_m));
+ // */
+
+ if (W_p <= 0.0 || W_m <= 0.0)
+ {
+ ERR("characteristic curve in negative region (W-, W+): %g, %g", W_m, W_p);
+ return 0.0;
+ }
+
+ // const double dAdlnW = 2.0*epsilon/(log(characteristic_curve(A+epsilon)) - log(characteristic_curve(A-epsilon)));
+ return dAdlnW * get_alphaT(T_Ads);
}
//Calculate sorption enthalpy
double Adsorption::get_enthalpy(const double p_Ads, const double T_Ads, const double M_Ads) const
{
- // TODO [CL] consider using A as obtained from current loading (needs inverse CC A(W)) instead of p_Vapour, T_Vapour
- const double A = get_potential(p_Ads, T_Ads, M_Ads);
+ // TODO [CL] consider using A as obtained from current loading (needs inverse CC A(W)) instead of p_Vapour, T_Vapour
+ const double A = get_potential(p_Ads, T_Ads, M_Ads);
- // return 0.0; // TODO [CL] for testing only
- return (get_evaporation_enthalpy(T_Ads) + A - T_Ads * get_entropy(T_Ads,A))*1000.0; //in J/kg
+ // return 0.0; // TODO [CL] for testing only
+ return (get_evaporation_enthalpy(T_Ads) + A - T_Ads * get_entropy(T_Ads,A))*1000.0; //in J/kg
}
double Adsorption::get_equilibrium_loading(const double p_Ads, const double T_Ads,
- const double M_Ads) const
+ const double M_Ads) const
{
- const double A = get_potential(p_Ads, T_Ads, M_Ads);
- return get_adsorbate_density(T_Ads) * characteristic_curve(A);
+ const double A = get_potential(p_Ads, T_Ads, M_Ads);
+ return get_adsorbate_density(T_Ads) * characteristic_curve(A);
}
diff --git a/MaterialsLib/Adsorption/Adsorption.h b/MaterialsLib/Adsorption/Adsorption.h
index c2030cd4db1cecc1499cdd9d380e8a0c5a4bb7ab..7ea6df0a4c9dd26d16d92002acb8d32775576394 100644
--- a/MaterialsLib/Adsorption/Adsorption.h
+++ b/MaterialsLib/Adsorption/Adsorption.h
@@ -17,80 +17,80 @@ const double M_H2O = 0.018;
class Adsorption : public Reaction
{
public:
- // TODO [CL] move those three methods to water properties class
- static double get_evaporation_enthalpy(const double T_Ads);
- static double get_equilibrium_vapour_pressure(const double T_Ads);
- static double get_specific_heat_capacity(const double T_Ads); // TODO [CL] why unused?
+ // TODO [CL] move those three methods to water properties class
+ static double get_evaporation_enthalpy(const double T_Ads);
+ static double get_equilibrium_vapour_pressure(const double T_Ads);
+ static double get_specific_heat_capacity(const double T_Ads); // TODO [CL] why unused?
- static double get_molar_fraction(double xm, double M_this, double M_other);
- static double get_mass_fraction(double xn, double M_this, double M_other);
- static double d_molar_fraction(double xm, double M_this, double M_other);
+ static double get_molar_fraction(double xm, double M_this, double M_other);
+ static double get_mass_fraction(double xn, double M_this, double M_other);
+ static double d_molar_fraction(double xm, double M_this, double M_other);
- static double get_loading(const double rho_curr, const double rho_dry);
+ static double get_loading(const double rho_curr, const double rho_dry);
// non-virtual members
- double get_equilibrium_loading(const double p_Ads, const double T_Ads, const double M_Ads)
- const override;
+ double get_equilibrium_loading(const double p_Ads, const double T_Ads, const double M_Ads)
+ const override;
// virtual members:
- virtual ~Adsorption() = default;
-
- virtual double get_enthalpy(const double p_Ads, const double T_Ads, const double M_Ads) const override;
- virtual double get_reaction_rate(const double p_Ads, const double T_Ads,
- const double M_Ads, const double loading) const override;
- /**
- * @brief get_d_reaction_rate
- * @param p_Ads
- * @param T_ads
- * @param M_Ads
- * @param loading
- * @param dqdr array containing the differentials wrt: p, T, C
- */
- virtual void get_d_reaction_rate(const double p_Ads, const double T_Ads,
- const double M_Ads, const double loading,
- std::array<double, 3>& dqdr) const;
+ virtual ~Adsorption() = default;
+
+ virtual double get_enthalpy(const double p_Ads, const double T_Ads, const double M_Ads) const override;
+ virtual double get_reaction_rate(const double p_Ads, const double T_Ads,
+ const double M_Ads, const double loading) const override;
+ /**
+ * @brief get_d_reaction_rate
+ * @param p_Ads
+ * @param T_ads
+ * @param M_Ads
+ * @param loading
+ * @param dqdr array containing the differentials wrt: p, T, C
+ */
+ virtual void get_d_reaction_rate(const double p_Ads, const double T_Ads,
+ const double M_Ads, const double loading,
+ std::array<double, 3>& dqdr) const;
protected:
- virtual double get_adsorbate_density(const double T_Ads) const = 0;
- virtual double get_alphaT(const double T_Ads) const = 0;
- virtual double characteristic_curve(const double A) const = 0;
- virtual double d_characteristic_curve(const double A) const = 0;
+ virtual double get_adsorbate_density(const double T_Ads) const = 0;
+ virtual double get_alphaT(const double T_Ads) const = 0;
+ virtual double characteristic_curve(const double A) const = 0;
+ virtual double d_characteristic_curve(const double A) const = 0;
private:
// non-virtual members
- double get_potential(const double p_Ads, const double T_Ads, const double M_Ads) const;
- double get_entropy(const double T_Ads, const double A) const;
+ double get_potential(const double p_Ads, const double T_Ads, const double M_Ads) const;
+ double get_entropy(const double T_Ads, const double A) const;
};
inline double curve_polyfrac(const double* coeffs, const double x)
{
- return ( coeffs[0] + coeffs[2] * x + coeffs[4] * pow(x,2) + coeffs[6] * pow(x,3) )
- / ( 1.0 + coeffs[1] * x + coeffs[3] * pow(x,2) + coeffs[5] * pow(x,3) );
+ return ( coeffs[0] + coeffs[2] * x + coeffs[4] * pow(x,2) + coeffs[6] * pow(x,3) )
+ / ( 1.0 + coeffs[1] * x + coeffs[3] * pow(x,2) + coeffs[5] * pow(x,3) );
- // Apparently, even pow and std::pow are different
- // return ( coeffs[0] + coeffs[2] * x + coeffs[4] * std::pow(x,2) + coeffs[6] * std::pow(x,3) )
- // / ( 1.0 + coeffs[1] * x + coeffs[3] * std::pow(x,2) + coeffs[5] * std::pow(x,3) );
+ // Apparently, even pow and std::pow are different
+ // return ( coeffs[0] + coeffs[2] * x + coeffs[4] * std::pow(x,2) + coeffs[6] * std::pow(x,3) )
+ // / ( 1.0 + coeffs[1] * x + coeffs[3] * std::pow(x,2) + coeffs[5] * std::pow(x,3) );
- // Analytically the same, but numerically quite different
- // return ( coeffs[0] + x * ( coeffs[2] + x * (coeffs[4] + x * coeffs[6] ) ) )
- // / ( 1.0 + x * ( coeffs[1] + x * (coeffs[3] + x * coeffs[5] ) ) );
+ // Analytically the same, but numerically quite different
+ // return ( coeffs[0] + x * ( coeffs[2] + x * (coeffs[4] + x * coeffs[6] ) ) )
+ // / ( 1.0 + x * ( coeffs[1] + x * (coeffs[3] + x * coeffs[5] ) ) );
- // Analytically the same, but numerically quite different
- // return ( coeffs[0] + x * coeffs[2] + x*x * coeffs[4] + x*x*x * coeffs[6] )
- // / ( 1.0 + x * coeffs[1] + x*x * coeffs[3] + x*x*x * coeffs[5] );
+ // Analytically the same, but numerically quite different
+ // return ( coeffs[0] + x * coeffs[2] + x*x * coeffs[4] + x*x*x * coeffs[6] )
+ // / ( 1.0 + x * coeffs[1] + x*x * coeffs[3] + x*x*x * coeffs[5] );
}
inline double d_curve_polyfrac(const double* coeffs, const double x)
{
- const double x2 = x*x;
- const double x3 = x2*x;
- const double u = coeffs[0] + coeffs[2] * x + coeffs[4] * x2 + coeffs[6] * x3;
- const double du = coeffs[2] + 2.0*coeffs[4] * x + 3.0*coeffs[6] * x2;
- const double v = 1.0 + coeffs[1] * x + coeffs[3] * x2 + coeffs[5] * x3;
- const double dv = coeffs[1] + 2.0*coeffs[3] * x + 3.0*coeffs[5] * x2;
-
- return (du*v - u*dv) / v / v;
+ const double x2 = x*x;
+ const double x3 = x2*x;
+ const double u = coeffs[0] + coeffs[2] * x + coeffs[4] * x2 + coeffs[6] * x3;
+ const double du = coeffs[2] + 2.0*coeffs[4] * x + 3.0*coeffs[6] * x2;
+ const double v = 1.0 + coeffs[1] * x + coeffs[3] * x2 + coeffs[5] * x3;
+ const double dv = coeffs[1] + 2.0*coeffs[3] * x + 3.0*coeffs[5] * x2;
+
+ return (du*v - u*dv) / v / v;
}
}
diff --git a/MaterialsLib/Adsorption/Density100MPa.cpp b/MaterialsLib/Adsorption/Density100MPa.cpp
index 92f214e953263f2b7321377fcdf919ae6218f4f8..b9a4161c2490bace092380782ef08a32ac8f2791 100644
--- a/MaterialsLib/Adsorption/Density100MPa.cpp
+++ b/MaterialsLib/Adsorption/Density100MPa.cpp
@@ -6,13 +6,13 @@ namespace
// NaX_HighP_polyfrac_CC.pickle
// date extracted 2015-06-23 15:38:35 file mtime 2015-06-23 15:19:57
const double c[] = {
- 0.3490302932983226, /* a0 */
- -0.0014061345691831226, /* a1 */
- -0.0007399303393402753, /* a2 */
- 5.129318840267485e-09, /* a3 */
- 5.243619689772646e-07, /* a4 */
- 6.347011955956523e-10, /* a5 */
- -9.919599580166727e-11 /* a6 */
+ 0.3490302932983226, /* a0 */
+ -0.0014061345691831226, /* a1 */
+ -0.0007399303393402753, /* a2 */
+ 5.129318840267485e-09, /* a3 */
+ 5.243619689772646e-07, /* a4 */
+ 6.347011955956523e-10, /* a5 */
+ -9.919599580166727e-11 /* a6 */
};
}
@@ -22,35 +22,35 @@ namespace Adsorption
double Density100MPa::get_adsorbate_density(const double T_Ads) const
{
- return -0.0013*T_Ads*T_Ads + 0.3529*T_Ads + 1049.2;
+ return -0.0013*T_Ads*T_Ads + 0.3529*T_Ads + 1049.2;
}
//Thermal expansivity model for water found in the works of Hauer
double Density100MPa::get_alphaT(const double T_Ads) const
{
- const double rho = -0.0013*T_Ads*T_Ads+0.3529*T_Ads+1049.2;
- const double drhodT = -0.0026*T_Ads + 0.3529;
+ const double rho = -0.0013*T_Ads*T_Ads+0.3529*T_Ads+1049.2;
+ const double drhodT = -0.0026*T_Ads + 0.3529;
- return - drhodT / rho;
+ return - drhodT / rho;
}
//Characteristic curve. Return W (A)
double Density100MPa::characteristic_curve(const double A) const
{
- double W = curve_polyfrac(c, A); //cm^3/g
+ double W = curve_polyfrac(c, A); //cm^3/g
- if (W < 0.0) {
- W = 0.0; // TODO [CL] debug output
- }
+ if (W < 0.0) {
+ W = 0.0; // TODO [CL] debug output
+ }
- return W/1.e3; //m^3/kg
+ return W/1.e3; //m^3/kg
}
double Density100MPa::d_characteristic_curve(const double A) const
{
- return d_curve_polyfrac(c, A);
+ return d_curve_polyfrac(c, A);
}
}
diff --git a/MaterialsLib/Adsorption/Density100MPa.h b/MaterialsLib/Adsorption/Density100MPa.h
index 78bf2c0e08bcae674385a18597a43913c738e059..b3268f4f740a6ea9f7b392dfd803b702ff36ecdd 100644
--- a/MaterialsLib/Adsorption/Density100MPa.h
+++ b/MaterialsLib/Adsorption/Density100MPa.h
@@ -8,10 +8,10 @@ namespace Adsorption
class Density100MPa : public Adsorption
{
public:
- double get_adsorbate_density(const double T_Ads) const;
- double get_alphaT(const double T_Ads) const;
- double characteristic_curve(const double A) const;
- double d_characteristic_curve(const double A) const;
+ double get_adsorbate_density(const double T_Ads) const;
+ double get_alphaT(const double T_Ads) const;
+ double characteristic_curve(const double A) const;
+ double d_characteristic_curve(const double A) const;
};
}
diff --git a/MaterialsLib/Adsorption/DensityConst.cpp b/MaterialsLib/Adsorption/DensityConst.cpp
index 5af46d852819ccc0724b01ddfb329ae0054bced9..01f7e06f8d275901e4c929c1098ec4d1bb9948bf 100644
--- a/MaterialsLib/Adsorption/DensityConst.cpp
+++ b/MaterialsLib/Adsorption/DensityConst.cpp
@@ -7,13 +7,13 @@ namespace
// NaX_Constant_polyfrac_CC.pickle
// date extracted 2015-06-23 15:38:35 file mtime 2015-06-23 15:20:05
const double c[] = {
- 0.3824098506898007, /* a0 */
- -0.001316857559708455, /* a1 */
- -0.0007935756090263691, /* a2 */
- -1.1600036977157845e-07, /* a3 */
- 5.610354459181838e-07, /* a4 */
- 7.113664938298873e-10, /* a5 */
- -1.0668790477629686e-10 /* a6 */
+ 0.3824098506898007, /* a0 */
+ -0.001316857559708455, /* a1 */
+ -0.0007935756090263691, /* a2 */
+ -1.1600036977157845e-07, /* a3 */
+ 5.610354459181838e-07, /* a4 */
+ 7.113664938298873e-10, /* a5 */
+ -1.0668790477629686e-10 /* a6 */
};
}
@@ -23,32 +23,32 @@ namespace Adsorption
double DensityConst::get_adsorbate_density(const double /*T_Ads*/) const
{
- return rho_water_Hauer(150.0+273.15);
+ return rho_water_Hauer(150.0+273.15);
}
//Thermal expansivity model for water found in the works of Hauer
double DensityConst::get_alphaT(const double /*T_Ads*/) const
{
- return 0.0;
+ return 0.0;
}
//Characteristic curve. Return W (A)
double DensityConst::characteristic_curve(const double A) const
{
- double W = curve_polyfrac(c, A); //cm^3/g
+ double W = curve_polyfrac(c, A); //cm^3/g
- if (W < 0.0) {
- W = 0.0; // TODO [CL] debug output
- }
+ if (W < 0.0) {
+ W = 0.0; // TODO [CL] debug output
+ }
- return W/1.e3; //m^3/kg
+ return W/1.e3; //m^3/kg
}
double DensityConst::d_characteristic_curve(const double A) const
{
- return d_curve_polyfrac(c, A);
+ return d_curve_polyfrac(c, A);
}
}
diff --git a/MaterialsLib/Adsorption/DensityConst.h b/MaterialsLib/Adsorption/DensityConst.h
index 4d9d97603dd8080eb0bd906d6a99fde227f3804e..a181b630623baaccb84ae48ed51431841f4b2cdc 100644
--- a/MaterialsLib/Adsorption/DensityConst.h
+++ b/MaterialsLib/Adsorption/DensityConst.h
@@ -8,10 +8,10 @@ namespace Adsorption
class DensityConst : public Adsorption
{
public:
- double get_adsorbate_density(const double T_Ads) const;
- double get_alphaT(const double T_Ads) const;
- double characteristic_curve(const double A) const;
- double d_characteristic_curve(const double A) const;
+ double get_adsorbate_density(const double T_Ads) const;
+ double get_alphaT(const double T_Ads) const;
+ double characteristic_curve(const double A) const;
+ double d_characteristic_curve(const double A) const;
};
}
diff --git a/MaterialsLib/Adsorption/DensityCook.cpp b/MaterialsLib/Adsorption/DensityCook.cpp
index 5e2f905f78697e922165e24a21a55e45ba3d36a5..52e3de22973123bb13336756789bc04d382b20ba 100644
--- a/MaterialsLib/Adsorption/DensityCook.cpp
+++ b/MaterialsLib/Adsorption/DensityCook.cpp
@@ -6,13 +6,13 @@ namespace
// NaX_Dean_polyfrac_CC.pickle
// date extracted 2015-06-23 15:38:35 file mtime 2015-06-23 15:19:42
const double c[] = {
- 0.3632627555646154, /* a0 */
- -0.0014090624975800715, /* a1 */
- -0.0007717609035743321, /* a2 */
- 5.03634836561135e-09, /* a3 */
- 5.478509959282738e-07, /* a4 */
- 6.36458510620815e-10, /* a5 */
- -1.037977321231462e-10 /* a6 */
+ 0.3632627555646154, /* a0 */
+ -0.0014090624975800715, /* a1 */
+ -0.0007717609035743321, /* a2 */
+ 5.03634836561135e-09, /* a3 */
+ 5.478509959282738e-07, /* a4 */
+ 6.36458510620815e-10, /* a5 */
+ -1.037977321231462e-10 /* a6 */
};
}
@@ -22,32 +22,32 @@ namespace Adsorption
double DensityCook::get_adsorbate_density(const double T_Ads) const
{
- return rho_water_Dean(T_Ads);
+ return rho_water_Dean(T_Ads);
}
//Thermal expansivity model for water found in the works of Hauer
double DensityCook::get_alphaT(const double T_Ads) const
{
- return alphaT_water_Dean(T_Ads);
+ return alphaT_water_Dean(T_Ads);
}
//Characteristic curve. Return W (A)
double DensityCook::characteristic_curve(const double A) const
{
- double W = curve_polyfrac(c, A); //cm^3/g
+ double W = curve_polyfrac(c, A); //cm^3/g
- if (W < 0.0) {
- W = 0.0; // TODO [CL] debug output
- }
+ if (W < 0.0) {
+ W = 0.0; // TODO [CL] debug output
+ }
- return W/1.e3; //m^3/kg
+ return W/1.e3; //m^3/kg
}
double DensityCook::d_characteristic_curve(const double A) const
{
- return d_curve_polyfrac(c, A);
+ return d_curve_polyfrac(c, A);
}
}
diff --git a/MaterialsLib/Adsorption/DensityCook.h b/MaterialsLib/Adsorption/DensityCook.h
index 6f6bc5fdf59ccc2c8fff50e66d550a83e73ee596..8876eba6d06be142da1d770b75a2bba569c6eaf5 100644
--- a/MaterialsLib/Adsorption/DensityCook.h
+++ b/MaterialsLib/Adsorption/DensityCook.h
@@ -8,40 +8,40 @@ namespace Adsorption
class DensityCook : public Adsorption
{
public:
- double get_adsorbate_density(const double T_Ads) const;
- double get_alphaT(const double T_Ads) const;
- double characteristic_curve(const double A) const;
- double d_characteristic_curve(const double A) const;
+ double get_adsorbate_density(const double T_Ads) const;
+ double get_alphaT(const double T_Ads) const;
+ double characteristic_curve(const double A) const;
+ double d_characteristic_curve(const double A) const;
};
inline double rho_water_Dean(const double T_Ads)
{
- const double Tcel = T_Ads - 273.15;
- const double b[] = { 999.9,2.03E-02,-6.16E-03,2.26E-05,-4.68E-08 };
- if (Tcel <= 100.) {
- return b[0] + Tcel * (b[1] + Tcel * (b[2] + Tcel * (b[3] + Tcel * b[4]) ) );
- }
- else {
- const double rho_100 = b[0] + b[1]*1.e2 + b[2]*1.e4 + b[3]*1.e6 + b[4]*1.e8;
- const double aT_100 = -1./rho_100 * (b[1] + 2.*b[2]*1.e2 + 3.*b[3]*1.e4 + 4.*b[4]*1.e6);
- return rho_100 * (1. - aT_100*(Tcel-100.));
- }
+ const double Tcel = T_Ads - 273.15;
+ const double b[] = { 999.9,2.03E-02,-6.16E-03,2.26E-05,-4.68E-08 };
+ if (Tcel <= 100.) {
+ return b[0] + Tcel * (b[1] + Tcel * (b[2] + Tcel * (b[3] + Tcel * b[4]) ) );
+ }
+ else {
+ const double rho_100 = b[0] + b[1]*1.e2 + b[2]*1.e4 + b[3]*1.e6 + b[4]*1.e8;
+ const double aT_100 = -1./rho_100 * (b[1] + 2.*b[2]*1.e2 + 3.*b[3]*1.e4 + 4.*b[4]*1.e6);
+ return rho_100 * (1. - aT_100*(Tcel-100.));
+ }
}
inline double alphaT_water_Dean(const double T_Ads)
{
- const double Tcel = T_Ads - 273.15;
- const double b[] = { 999.9,2.03E-02,-6.16E-03,2.26E-05,-4.68E-08 };
- if (Tcel <= 100.) {
- const double r = b[0] + Tcel * (b[1] + Tcel * (b[2] + Tcel * (b[3] + Tcel * b[4]) ) );
- return -1.0/r * ( b[1] + Tcel * (2.0*b[2] + Tcel * (3.0*b[3] + Tcel * 4.0*b[4]) ) );
- }
- else {
- const double rho_100 = b[0] + b[1]*1.e2 + b[2]*1.e4 + b[3]*1.e6 + b[4]*1.e8;
- const double aT_100 = -1./rho_100 * (b[1] + 2.*b[2]*1.e2 + 3.*b[3]*1.e4 + 4.*b[4]*1.e6);
- return aT_100 / (1. - aT_100*(Tcel-100.));
- }
+ const double Tcel = T_Ads - 273.15;
+ const double b[] = { 999.9,2.03E-02,-6.16E-03,2.26E-05,-4.68E-08 };
+ if (Tcel <= 100.) {
+ const double r = b[0] + Tcel * (b[1] + Tcel * (b[2] + Tcel * (b[3] + Tcel * b[4]) ) );
+ return -1.0/r * ( b[1] + Tcel * (2.0*b[2] + Tcel * (3.0*b[3] + Tcel * 4.0*b[4]) ) );
+ }
+ else {
+ const double rho_100 = b[0] + b[1]*1.e2 + b[2]*1.e4 + b[3]*1.e6 + b[4]*1.e8;
+ const double aT_100 = -1./rho_100 * (b[1] + 2.*b[2]*1.e2 + 3.*b[3]*1.e4 + 4.*b[4]*1.e6);
+ return aT_100 / (1. - aT_100*(Tcel-100.));
+ }
}
}
diff --git a/MaterialsLib/Adsorption/DensityDubinin.cpp b/MaterialsLib/Adsorption/DensityDubinin.cpp
index 3b3138de4350e2fa928a55be92f62d99abbf1790..71dd07cfa77ca65544af29aa994b978eabdb5959 100644
--- a/MaterialsLib/Adsorption/DensityDubinin.cpp
+++ b/MaterialsLib/Adsorption/DensityDubinin.cpp
@@ -8,13 +8,13 @@ namespace
// NaX_Dubinin_polyfrac_CC.pickle
// date extracted 2015-06-23 16:47:50 file mtime 2015-06-23 16:47:23
const double c[] = {
- 0.3635538371322433, /* a0 */
- -0.0014521033261199435, /* a1 */
- -0.0007855160157616825, /* a2 */
- 4.385666000850253e-08, /* a3 */
- 5.567776459188524e-07, /* a4 */
- 6.026002134230559e-10, /* a5 */
- -1.0477401124006098e-10 /* a6 */
+ 0.3635538371322433, /* a0 */
+ -0.0014521033261199435, /* a1 */
+ -0.0007855160157616825, /* a2 */
+ 4.385666000850253e-08, /* a3 */
+ 5.567776459188524e-07, /* a4 */
+ 6.026002134230559e-10, /* a5 */
+ -1.0477401124006098e-10 /* a6 */
};
}
@@ -24,66 +24,66 @@ namespace Adsorption
double DensityDubinin::get_adsorbate_density(const double T_Ads) const
{
- const double Tb = 373.1;
+ const double Tb = 373.1;
- if (T_Ads < Tb) {
- return rho_water_Dean(T_Ads);
- } else {
- const double Tc = 647.3; //K
- // const double rhoc = 322.; //kg/m^3
- const double pc = 221.2e5; //Pa
- //boiling point density
- const double rhob = rho_water_Dean(Tb);
- //state values
- const double R = GAS_CONST;
- const double M = M_H2O;
- const double b = R * Tc/(8. * pc); //m^3/mol
- const double rhom = M/b; //kg/m^3
- const double rho = rhob - (rhob-rhom)/(Tc-Tb)*(T_Ads-Tb);
- return rho;
- }
+ if (T_Ads < Tb) {
+ return rho_water_Dean(T_Ads);
+ } else {
+ const double Tc = 647.3; //K
+ // const double rhoc = 322.; //kg/m^3
+ const double pc = 221.2e5; //Pa
+ //boiling point density
+ const double rhob = rho_water_Dean(Tb);
+ //state values
+ const double R = GAS_CONST;
+ const double M = M_H2O;
+ const double b = R * Tc/(8. * pc); //m^3/mol
+ const double rhom = M/b; //kg/m^3
+ const double rho = rhob - (rhob-rhom)/(Tc-Tb)*(T_Ads-Tb);
+ return rho;
+ }
}
//Thermal expansivity model for water found in the works of Hauer
double DensityDubinin::get_alphaT(const double T_Ads) const
{
- const double Tb = 373.1;
- if (T_Ads <= Tb) {
- return alphaT_water_Dean(T_Ads);
- } else {
- //critical T and p
- const double Tc = 647.3; //K
- // const double rhoc = 322.; //kg/m^3
- const double pc = 221.2e5; //Pa
- //boiling point density
- const double rhob = rho_water_Dean(Tb);
- //state values
- const double R = GAS_CONST;
- const double M = M_H2O;
- const double b = R * Tc/(8. * pc); //m^3/mol
- const double rhom = M/(b); //kg/m^3
- const double rho = rhob - (rhob-rhom)/(Tc-Tb)*(T_Ads-Tb);
- return ((rhob-rhom)/(Tc-Tb)*1./rho);
- }
+ const double Tb = 373.1;
+ if (T_Ads <= Tb) {
+ return alphaT_water_Dean(T_Ads);
+ } else {
+ //critical T and p
+ const double Tc = 647.3; //K
+ // const double rhoc = 322.; //kg/m^3
+ const double pc = 221.2e5; //Pa
+ //boiling point density
+ const double rhob = rho_water_Dean(Tb);
+ //state values
+ const double R = GAS_CONST;
+ const double M = M_H2O;
+ const double b = R * Tc/(8. * pc); //m^3/mol
+ const double rhom = M/(b); //kg/m^3
+ const double rho = rhob - (rhob-rhom)/(Tc-Tb)*(T_Ads-Tb);
+ return ((rhob-rhom)/(Tc-Tb)*1./rho);
+ }
}
//Characteristic curve. Return W (A)
double DensityDubinin::characteristic_curve(const double A) const
{
- double W = curve_polyfrac(c, A); //cm^3/g
+ double W = curve_polyfrac(c, A); //cm^3/g
- if (W < 0.0) {
- W = 0.0; // TODO [CL] debug output
- }
+ if (W < 0.0) {
+ W = 0.0; // TODO [CL] debug output
+ }
- return W/1.e3; //m^3/kg
+ return W/1.e3; //m^3/kg
}
double DensityDubinin::d_characteristic_curve(const double A) const
{
- return d_curve_polyfrac(c, A);
+ return d_curve_polyfrac(c, A);
}
}
diff --git a/MaterialsLib/Adsorption/DensityDubinin.h b/MaterialsLib/Adsorption/DensityDubinin.h
index a82738d54eaee80a776a00afda88ff004643b914..2e1bb4b8c0c30d588307863caec0fbffd20141bb 100644
--- a/MaterialsLib/Adsorption/DensityDubinin.h
+++ b/MaterialsLib/Adsorption/DensityDubinin.h
@@ -8,10 +8,10 @@ namespace Adsorption
class DensityDubinin : public Adsorption
{
public:
- double get_adsorbate_density(const double T_Ads) const;
- double get_alphaT(const double T_Ads) const;
- double characteristic_curve(const double A) const;
- double d_characteristic_curve(const double A) const;
+ double get_adsorbate_density(const double T_Ads) const;
+ double get_alphaT(const double T_Ads) const;
+ double characteristic_curve(const double A) const;
+ double d_characteristic_curve(const double A) const;
};
}
diff --git a/MaterialsLib/Adsorption/DensityHauer.cpp b/MaterialsLib/Adsorption/DensityHauer.cpp
index 4e441676ecc4441d2c66d8c02c04e8cd9de6562a..be7bee83c026fc4f6c818afc42ea33501fc04b8d 100644
--- a/MaterialsLib/Adsorption/DensityHauer.cpp
+++ b/MaterialsLib/Adsorption/DensityHauer.cpp
@@ -6,13 +6,13 @@ namespace
// NaX_Hauer_polyfrac_CC.pickle
// date extracted 2015-06-23 15:38:35 file mtime 2015-06-23 15:19:19
const double c[] = {
- 0.36490158988356747, /* a0 */
- -0.0013723270478333963, /* a1 */
- -0.0007655780628099964, /* a2 */
- -3.353324854315774e-08, /* a3 */
- 5.424357157710913e-07, /* a4 */
- 6.613430586648678e-10, /* a5 */
- -1.0300151379421499e-10 /* a6 */
+ 0.36490158988356747, /* a0 */
+ -0.0013723270478333963, /* a1 */
+ -0.0007655780628099964, /* a2 */
+ -3.353324854315774e-08, /* a3 */
+ 5.424357157710913e-07, /* a4 */
+ 6.613430586648678e-10, /* a5 */
+ -1.0300151379421499e-10 /* a6 */
};
}
@@ -22,35 +22,35 @@ namespace Adsorption
double DensityHauer::get_adsorbate_density(const double T_Ads) const
{
- return rho_water_Hauer(T_Ads);
+ return rho_water_Hauer(T_Ads);
}
//Thermal expansivity model for water found in the works of Hauer
double DensityHauer::get_alphaT(const double T_Ads) const
{
- // data like in python script
- const double T0 = 283.15, alpha0 = 3.781e-4; //K; 1/K
+ // data like in python script
+ const double T0 = 283.15, alpha0 = 3.781e-4; //K; 1/K
- return alpha0/(1. - alpha0 * (T_Ads-T0)); //in 1/K
+ return alpha0/(1. - alpha0 * (T_Ads-T0)); //in 1/K
}
//Characteristic curve. Return W (A)
double DensityHauer::characteristic_curve(const double A) const
{
- double W = curve_polyfrac(c, A); //cm^3/g
+ double W = curve_polyfrac(c, A); //cm^3/g
- if (W < 0.0) {
- W = 0.0; // TODO [CL] debug output
- }
+ if (W < 0.0) {
+ W = 0.0; // TODO [CL] debug output
+ }
- return W/1.e3; //m^3/kg
+ return W/1.e3; //m^3/kg
}
double DensityHauer::d_characteristic_curve(const double A) const
{
- return d_curve_polyfrac(c, A);
+ return d_curve_polyfrac(c, A);
}
}
diff --git a/MaterialsLib/Adsorption/DensityHauer.h b/MaterialsLib/Adsorption/DensityHauer.h
index 8180a7756ca050867de1aac472c1b2258722869e..6c8a3c521b8fdd6aa22cd9fcaba805e7d1e0b7a1 100644
--- a/MaterialsLib/Adsorption/DensityHauer.h
+++ b/MaterialsLib/Adsorption/DensityHauer.h
@@ -9,18 +9,18 @@ namespace Adsorption
class DensityHauer : public Adsorption
{
public:
- double get_adsorbate_density(const double T_Ads) const;
- double get_alphaT(const double T_Ads) const;
- double characteristic_curve(const double A) const;
- double d_characteristic_curve(const double A) const;
+ double get_adsorbate_density(const double T_Ads) const;
+ double get_alphaT(const double T_Ads) const;
+ double characteristic_curve(const double A) const;
+ double d_characteristic_curve(const double A) const;
};
inline double rho_water_Hauer(const double T_Ads)
{
- // data like in python script
- const double T0 = 283.15, rho0 = rho_water_Dean(T0), alpha0 = 3.781e-4; //K; kg/m^3; 1/K
+ // data like in python script
+ const double T0 = 283.15, rho0 = rho_water_Dean(T0), alpha0 = 3.781e-4; //K; kg/m^3; 1/K
- return rho0 * (1. - alpha0 * (T_Ads-T0)); //in kg/m^3
+ return rho0 * (1. - alpha0 * (T_Ads-T0)); //in kg/m^3
}
}
diff --git a/MaterialsLib/Adsorption/DensityLegacy.cpp b/MaterialsLib/Adsorption/DensityLegacy.cpp
index 6969aa1448ac8c03d02492637f4741deb11a82e7..298cbeeb453d9dad5e77ae8b82f898abb207065a 100644
--- a/MaterialsLib/Adsorption/DensityLegacy.cpp
+++ b/MaterialsLib/Adsorption/DensityLegacy.cpp
@@ -5,13 +5,13 @@ namespace
//parameters from least squares fit (experimental data)
const double c[] = { 0.34102920966608297,
- -0.0013106032830951296,
- -0.00060754147575378876,
- 3.7843404172683339e-07,
- 4.0107503869519016e-07,
- 3.1274595098338057e-10,
- -7.610441241719489e-11
- };
+ -0.0013106032830951296,
+ -0.00060754147575378876,
+ 3.7843404172683339e-07,
+ 4.0107503869519016e-07,
+ 3.1274595098338057e-10,
+ -7.610441241719489e-11
+ };
}
@@ -20,40 +20,40 @@ namespace Adsorption
double DensityLegacy::get_adsorbate_density(const double T_Ads) const
{
- //set reference state for adsorbate EOS in Hauer
- const double T0 = 293.15, rho0 = 998.084, alpha0 = 2.06508e-4; //K; kg/m^3; 1/K
+ //set reference state for adsorbate EOS in Hauer
+ const double T0 = 293.15, rho0 = 998.084, alpha0 = 2.06508e-4; //K; kg/m^3; 1/K
- return (rho0 * (1. - alpha0 * (T_Ads-T0))); //in kg/m^3
+ return (rho0 * (1. - alpha0 * (T_Ads-T0))); //in kg/m^3
}
//Thermal expansivity model for water found in the works of Hauer
double DensityLegacy::get_alphaT(const double T_Ads) const
{
- //set reference state for adsorbate EOS in Hauer
- const double T0 = 293.15, /*rho0 = 998.084,*/ alpha0 = 2.06508e-4; //K; kg/m^3; 1/K
+ //set reference state for adsorbate EOS in Hauer
+ const double T0 = 293.15, /*rho0 = 998.084,*/ alpha0 = 2.06508e-4; //K; kg/m^3; 1/K
- return (alpha0/(1. - alpha0 * (T_Ads-T0))); //in 1/K
+ return (alpha0/(1. - alpha0 * (T_Ads-T0))); //in 1/K
}
//Characteristic curve. Return W (A)
double DensityLegacy::characteristic_curve(const double A) const
{
- double W = curve_polyfrac(c, A); //cm^3/g
+ double W = curve_polyfrac(c, A); //cm^3/g
- /*
- if (W < 0.0) {
- W = 0.0; // TODO [CL] debug output
- }
- */
+ /*
+ if (W < 0.0) {
+ W = 0.0; // TODO [CL] debug output
+ }
+ */
- return W/1.e3; //m^3/kg
+ return W/1.e3; //m^3/kg
}
double DensityLegacy::d_characteristic_curve(const double A) const
{
- return d_curve_polyfrac(c, A);
+ return d_curve_polyfrac(c, A);
}
}
diff --git a/MaterialsLib/Adsorption/DensityLegacy.h b/MaterialsLib/Adsorption/DensityLegacy.h
index 721155c7b51bc4d507489b01deb5d0861281c2d7..2028a87eb5fe9deae2652b3099335d3978a045ba 100644
--- a/MaterialsLib/Adsorption/DensityLegacy.h
+++ b/MaterialsLib/Adsorption/DensityLegacy.h
@@ -8,10 +8,10 @@ namespace Adsorption
class DensityLegacy : public Adsorption
{
public:
- double get_adsorbate_density(const double T_Ads) const;
- double get_alphaT(const double T_Ads) const;
- double characteristic_curve(const double A) const;
- double d_characteristic_curve(const double A) const;
+ double get_adsorbate_density(const double T_Ads) const;
+ double get_alphaT(const double T_Ads) const;
+ double characteristic_curve(const double A) const;
+ double d_characteristic_curve(const double A) const;
};
}
diff --git a/MaterialsLib/Adsorption/DensityMette.cpp b/MaterialsLib/Adsorption/DensityMette.cpp
index f2f14c0cd0a0d4473d667f800e8dfa8bee9c47e6..eea5b9e6f96a7dccf2e0349c1db5c281e26e5867 100644
--- a/MaterialsLib/Adsorption/DensityMette.cpp
+++ b/MaterialsLib/Adsorption/DensityMette.cpp
@@ -6,13 +6,13 @@ namespace
// NaX_Mette_polyfrac_CC.pickle
// date extracted 2015-06-23 15:38:35 file mtime 2015-06-23 15:19:26
const double c[] = {
- 0.36340572890087813, /* a0 */
- -0.0013449597038375108, /* a1 */
- -0.0007581210111121073, /* a2 */
- -7.331279615575401e-08, /* a3 */
- 5.365656973806218e-07, /* a4 */
- 6.854673678427112e-10, /* a5 */
- -1.0197050219481966e-10 /* a6 */
+ 0.36340572890087813, /* a0 */
+ -0.0013449597038375108, /* a1 */
+ -0.0007581210111121073, /* a2 */
+ -7.331279615575401e-08, /* a3 */
+ 5.365656973806218e-07, /* a4 */
+ 6.854673678427112e-10, /* a5 */
+ -1.0197050219481966e-10 /* a6 */
};
}
@@ -21,37 +21,37 @@ namespace Adsorption
double DensityMette::get_adsorbate_density(const double T_Ads) const
{
- const double T0 = 293.15;
- const double rho0 = rho_water_Dean(T0);
- const double alpha20 = alphaT_water_Dean(T0);
- return rho0 / (1. + alpha20*(T_Ads-T0));
+ const double T0 = 293.15;
+ const double rho0 = rho_water_Dean(T0);
+ const double alpha20 = alphaT_water_Dean(T0);
+ return rho0 / (1. + alpha20*(T_Ads-T0));
}
//Thermal expansivity model for water found in the works of Hauer
double DensityMette::get_alphaT(const double T_Ads) const
{
- const double T0 = 293.15;
- const double alpha20 = alphaT_water_Dean(T0);
- return alpha20 / (1. + alpha20 * (T_Ads-T0));
+ const double T0 = 293.15;
+ const double alpha20 = alphaT_water_Dean(T0);
+ return alpha20 / (1. + alpha20 * (T_Ads-T0));
}
//Characteristic curve. Return W (A)
double DensityMette::characteristic_curve(const double A) const
{
- double W = curve_polyfrac(c, A); //cm^3/g
+ double W = curve_polyfrac(c, A); //cm^3/g
- if (W < 0.0) {
- W = 0.0; // TODO [CL] debug output
- }
+ if (W < 0.0) {
+ W = 0.0; // TODO [CL] debug output
+ }
- return W/1.e3; //m^3/kg
+ return W/1.e3; //m^3/kg
}
double DensityMette::d_characteristic_curve(const double A) const
{
- return d_curve_polyfrac(c, A);
+ return d_curve_polyfrac(c, A);
}
}
diff --git a/MaterialsLib/Adsorption/DensityMette.h b/MaterialsLib/Adsorption/DensityMette.h
index 42ddea7289a32a5077132365f4aeeab3bf408e00..e087bbcaaea5c80e4f8c2b9afed76baa8550c2a6 100644
--- a/MaterialsLib/Adsorption/DensityMette.h
+++ b/MaterialsLib/Adsorption/DensityMette.h
@@ -8,10 +8,10 @@ namespace Adsorption
class DensityMette : public Adsorption
{
public:
- double get_adsorbate_density(const double T_Ads) const;
- double get_alphaT(const double T_Ads) const;
- double characteristic_curve(const double A) const;
- double d_characteristic_curve(const double A) const;
+ double get_adsorbate_density(const double T_Ads) const;
+ double get_alphaT(const double T_Ads) const;
+ double characteristic_curve(const double A) const;
+ double d_characteristic_curve(const double A) const;
};
}
diff --git a/MaterialsLib/Adsorption/DensityNunez.cpp b/MaterialsLib/Adsorption/DensityNunez.cpp
index ab560501cb671ef43af5f4df236e34805dfaa44d..61f4bd55ef2cbc3d35a67e92f45c976c2c5ecef3 100644
--- a/MaterialsLib/Adsorption/DensityNunez.cpp
+++ b/MaterialsLib/Adsorption/DensityNunez.cpp
@@ -6,13 +6,13 @@ namespace
// NaX_Nunez_polyfrac_CC.pickle
// date extracted 2015-06-23 15:38:35 file mtime 2015-06-23 15:19:34
const double c[] = {
- 0.3631900485031771, /* a0 */
- -0.0014242280940080726, /* a1 */
- -0.0007751726942386291, /* a2 */
- 2.1775655036811842e-08, /* a3 */
- 5.488166913667265e-07, /* a4 */
- 6.204064716725214e-10, /* a5 */
- -1.0345385018952998e-10 /* a6 */
+ 0.3631900485031771, /* a0 */
+ -0.0014242280940080726, /* a1 */
+ -0.0007751726942386291, /* a2 */
+ 2.1775655036811842e-08, /* a3 */
+ 5.488166913667265e-07, /* a4 */
+ 6.204064716725214e-10, /* a5 */
+ -1.0345385018952998e-10 /* a6 */
};
}
@@ -22,54 +22,54 @@ namespace Adsorption
double DensityNunez::get_adsorbate_density(const double T_Ads) const
{
- if (T_Ads < 273.16 or T_Ads > 633.15) {
- // print('Value outside admissible range for rho.');
- // return -1;
- }
+ if (T_Ads < 273.16 or T_Ads > 633.15) {
+ // print('Value outside admissible range for rho.');
+ // return -1;
+ }
- const double a[] = { 1.0644e3,-8.01905,1.445348e-2,-4.19589e-6,-4.5294e-9 };
- const double b[] = { -8.039e-3,1.8698e-5,-2.3015e-8,2.3809e-11,-1.388e-14 };
- const double u = a[0] + T_Ads * (a[1] + T_Ads * (a[2] + T_Ads * (a[3] + T_Ads * a[4]) ) );
- const double v = 1.0 + T_Ads * (b[0] + T_Ads * (b[1] + T_Ads * (b[2] + T_Ads * (b[3] + T_Ads * b[4]) ) ) );
- return u/v;
+ const double a[] = { 1.0644e3,-8.01905,1.445348e-2,-4.19589e-6,-4.5294e-9 };
+ const double b[] = { -8.039e-3,1.8698e-5,-2.3015e-8,2.3809e-11,-1.388e-14 };
+ const double u = a[0] + T_Ads * (a[1] + T_Ads * (a[2] + T_Ads * (a[3] + T_Ads * a[4]) ) );
+ const double v = 1.0 + T_Ads * (b[0] + T_Ads * (b[1] + T_Ads * (b[2] + T_Ads * (b[3] + T_Ads * b[4]) ) ) );
+ return u/v;
}
//Thermal expansivity model for water found in the works of Hauer
double DensityNunez::get_alphaT(const double T_Ads) const
{
- if (T_Ads < 273.16 or T_Ads > 633.15) {
- // print('Value outside admissible range for rho.');
- // return -1;
- }
+ if (T_Ads < 273.16 or T_Ads > 633.15) {
+ // print('Value outside admissible range for rho.');
+ // return -1;
+ }
- const double a[] = { 1.0644e3,-8.01905,1.445348e-2,-4.19589e-6,-4.5294e-9 };
- const double b[] = { -8.039e-3,1.8698e-5,-2.3015e-8,2.3809e-11,-1.388e-14 };
- const double u = a[0] + T_Ads * (a[1] + T_Ads * (a[2] + T_Ads * (a[3] + T_Ads * a[4]) ) );
- const double v = 1.0 + T_Ads * (b[0] + T_Ads * (b[1] + T_Ads * (b[2] + T_Ads * (b[3] + T_Ads * b[4]) ) ) );
- const double du = a[1] + T_Ads * (2.0*a[2] + T_Ads * (3.0*a[3] + T_Ads * 4.0*a[4]) );
- const double dv = b[0] + T_Ads * (2.0*b[1] + T_Ads * (3.0*b[2] + T_Ads * (4.0*b[3] + T_Ads * 5.0*b[4]) ) );
- // return - (du*v - dv*u) / u / v;
- // return (dv*u - du*v) / u / v;
- return dv/v - du/u;
+ const double a[] = { 1.0644e3,-8.01905,1.445348e-2,-4.19589e-6,-4.5294e-9 };
+ const double b[] = { -8.039e-3,1.8698e-5,-2.3015e-8,2.3809e-11,-1.388e-14 };
+ const double u = a[0] + T_Ads * (a[1] + T_Ads * (a[2] + T_Ads * (a[3] + T_Ads * a[4]) ) );
+ const double v = 1.0 + T_Ads * (b[0] + T_Ads * (b[1] + T_Ads * (b[2] + T_Ads * (b[3] + T_Ads * b[4]) ) ) );
+ const double du = a[1] + T_Ads * (2.0*a[2] + T_Ads * (3.0*a[3] + T_Ads * 4.0*a[4]) );
+ const double dv = b[0] + T_Ads * (2.0*b[1] + T_Ads * (3.0*b[2] + T_Ads * (4.0*b[3] + T_Ads * 5.0*b[4]) ) );
+ // return - (du*v - dv*u) / u / v;
+ // return (dv*u - du*v) / u / v;
+ return dv/v - du/u;
}
//Characteristic curve. Return W (A)
double DensityNunez::characteristic_curve(const double A) const
{
- double W = curve_polyfrac(c, A); //cm^3/g
+ double W = curve_polyfrac(c, A); //cm^3/g
- if (W < 0.0) {
- W = 0.0; // TODO [CL] debug output
- }
+ if (W < 0.0) {
+ W = 0.0; // TODO [CL] debug output
+ }
- return W/1.e3; //m^3/kg
+ return W/1.e3; //m^3/kg
}
double DensityNunez::d_characteristic_curve(const double A) const
{
- return d_curve_polyfrac(c, A);
+ return d_curve_polyfrac(c, A);
}
}
diff --git a/MaterialsLib/Adsorption/DensityNunez.h b/MaterialsLib/Adsorption/DensityNunez.h
index e3eea37d9142f08ad13770f0e6eeb2f0a5ef84f2..b93aad4a99ca4f5f79560c292db6e3db25093c55 100644
--- a/MaterialsLib/Adsorption/DensityNunez.h
+++ b/MaterialsLib/Adsorption/DensityNunez.h
@@ -8,10 +8,10 @@ namespace Adsorption
class DensityNunez : public Adsorption
{
public:
- double get_adsorbate_density(const double T_Ads) const;
- double get_alphaT(const double T_Ads) const;
- double characteristic_curve(const double A) const;
- double d_characteristic_curve(const double A) const;
+ double get_adsorbate_density(const double T_Ads) const;
+ double get_alphaT(const double T_Ads) const;
+ double characteristic_curve(const double A) const;
+ double d_characteristic_curve(const double A) const;
};
}
diff --git a/MaterialsLib/Adsorption/Reaction.cpp b/MaterialsLib/Adsorption/Reaction.cpp
index e79a19ee81fcfc5a5e3adbc4d62e79f4b2d0b755..6b0587b92cb809eddf6de30b1c22f63e12050384 100644
--- a/MaterialsLib/Adsorption/Reaction.cpp
+++ b/MaterialsLib/Adsorption/Reaction.cpp
@@ -32,35 +32,35 @@ std::unique_ptr<Reaction>
Reaction::
newInstance(BaseLib::ConfigTree const& conf)
{
- auto const type = conf.getConfParam<std::string>("type");
+ auto const type = conf.getConfParam<std::string>("type");
- if (type == "Z13XBF")
- return std::unique_ptr<Reaction>(new DensityLegacy);
- else if (type == "Z13XBF_100MPa")
- return std::unique_ptr<Reaction>(new Density100MPa);
- else if (type == "Z13XBF_Const")
- return std::unique_ptr<Reaction>(new DensityConst);
- else if (type == "Z13XBF_Cook")
- return std::unique_ptr<Reaction>(new DensityCook);
- else if (type == "Z13XBF_Dubinin")
- return std::unique_ptr<Reaction>(new DensityDubinin);
- else if (type == "Z13XBF_Hauer")
- return std::unique_ptr<Reaction>(new DensityHauer);
- else if (type == "Z13XBF_Mette")
- return std::unique_ptr<Reaction>(new DensityMette);
- else if (type == "Z13XBF_Nunez")
- return std::unique_ptr<Reaction>(new DensityNunez);
- else if (type == "Inert")
- return std::unique_ptr<Reaction>(new ReactionInert);
- else if (type == "Sinusoidal")
- return std::unique_ptr<Reaction>(new ReactionSinusoidal(conf));
- else if (type == "CaOH2")
- return std::unique_ptr<Reaction>(new ReactionCaOH2(conf));
+ if (type == "Z13XBF")
+ return std::unique_ptr<Reaction>(new DensityLegacy);
+ else if (type == "Z13XBF_100MPa")
+ return std::unique_ptr<Reaction>(new Density100MPa);
+ else if (type == "Z13XBF_Const")
+ return std::unique_ptr<Reaction>(new DensityConst);
+ else if (type == "Z13XBF_Cook")
+ return std::unique_ptr<Reaction>(new DensityCook);
+ else if (type == "Z13XBF_Dubinin")
+ return std::unique_ptr<Reaction>(new DensityDubinin);
+ else if (type == "Z13XBF_Hauer")
+ return std::unique_ptr<Reaction>(new DensityHauer);
+ else if (type == "Z13XBF_Mette")
+ return std::unique_ptr<Reaction>(new DensityMette);
+ else if (type == "Z13XBF_Nunez")
+ return std::unique_ptr<Reaction>(new DensityNunez);
+ else if (type == "Inert")
+ return std::unique_ptr<Reaction>(new ReactionInert);
+ else if (type == "Sinusoidal")
+ return std::unique_ptr<Reaction>(new ReactionSinusoidal(conf));
+ else if (type == "CaOH2")
+ return std::unique_ptr<Reaction>(new ReactionCaOH2(conf));
- ERR("Unknown reactive system: %s.", type.c_str());
- std::abort();
+ ERR("Unknown reactive system: %s.", type.c_str());
+ std::abort();
- return nullptr;
+ return nullptr;
}
} // namespace Ads
diff --git a/MaterialsLib/Adsorption/Reaction.h b/MaterialsLib/Adsorption/Reaction.h
index 29eba1038a2f2b3c422f02c958b1e6923d4fffd9..0dd998b9d027fec09f101901396f1bd9190593d9 100644
--- a/MaterialsLib/Adsorption/Reaction.h
+++ b/MaterialsLib/Adsorption/Reaction.h
@@ -19,16 +19,16 @@ namespace Adsorption
class Reaction
{
public:
- static std::unique_ptr<Reaction> newInstance(BaseLib::ConfigTree const& rsys);
+ static std::unique_ptr<Reaction> newInstance(BaseLib::ConfigTree const& rsys);
- virtual double get_enthalpy(const double p_Ads, const double T_Ads, const double M_Ads) const = 0;
- virtual double get_reaction_rate(const double p_Ads, const double T_Ads,
- const double M_Ads, const double loading) const = 0;
+ virtual double get_enthalpy(const double p_Ads, const double T_Ads, const double M_Ads) const = 0;
+ virtual double get_reaction_rate(const double p_Ads, const double T_Ads,
+ const double M_Ads, const double loading) const = 0;
- // TODO: get rid of
- virtual double get_equilibrium_loading(const double, const double, const double) const {
- return -1.0;
- }
+ // TODO: get rid of
+ virtual double get_equilibrium_loading(const double, const double, const double) const {
+ return -1.0;
+ }
virtual ~Reaction() = default;
};
diff --git a/MaterialsLib/Adsorption/ReactionCaOH2.cpp b/MaterialsLib/Adsorption/ReactionCaOH2.cpp
index 0ba026ba826f0931f5c3f47b249ff1bc6bb32313..ec3dbf9a06deff7a557ddeee217502c9eff7c293 100644
--- a/MaterialsLib/Adsorption/ReactionCaOH2.cpp
+++ b/MaterialsLib/Adsorption/ReactionCaOH2.cpp
@@ -47,88 +47,88 @@ ReactionCaOH2::get_reaction_rate(const double, const double, const double, const
double ReactionCaOH2::getReactionRate(double const solid_density)
{
- rho_s = solid_density;
- calculate_qR();
- return qR;
+ rho_s = solid_density;
+ calculate_qR();
+ return qR;
}
void ReactionCaOH2::update_param(
- double T_solid,
- double p_gas,
- double x_react,
- double rho_s_initial)
+ double T_solid,
+ double p_gas,
+ double x_react,
+ double rho_s_initial)
{
- T_s = T_solid;
- this->p_gas = p_gas / 1e5; // convert Pa to bar
- this->x_react = x_react;
- rho_s = rho_s_initial;
+ T_s = T_solid;
+ this->p_gas = p_gas / 1e5; // convert Pa to bar
+ this->x_react = x_react;
+ rho_s = rho_s_initial;
}
void ReactionCaOH2::calculate_qR()
{
- //Convert mass fraction into mole fraction
- // const double mol_frac_react = get_mole_fraction(x_react);
- const double mol_frac_react = Adsorption::Adsorption::get_molar_fraction(x_react, M_react, M_carrier);
-
- p_r_g = std::max(mol_frac_react * p_gas, 1.0e-3); //avoid illdefined log
- set_chemical_equilibrium();
- const double dXdt = Ca_hydration();
- qR = (rho_up - rho_low) * dXdt;
+ //Convert mass fraction into mole fraction
+ // const double mol_frac_react = get_mole_fraction(x_react);
+ const double mol_frac_react = Adsorption::Adsorption::get_molar_fraction(x_react, M_react, M_carrier);
+
+ p_r_g = std::max(mol_frac_react * p_gas, 1.0e-3); //avoid illdefined log
+ set_chemical_equilibrium();
+ const double dXdt = Ca_hydration();
+ qR = (rho_up - rho_low) * dXdt;
}
//determine equilibrium temperature and pressure according to van't Hoff
void ReactionCaOH2::set_chemical_equilibrium()
{
- X_D = (rho_s - rho_up - tol_rho)/(rho_low - rho_up - 2.0*tol_rho) ;
- X_D = (X_D < 0.5) ? std::max(tol_l,X_D) : std::min(X_D,tol_u); //constrain to interval [tol_l;tol_u]
+ X_D = (rho_s - rho_up - tol_rho)/(rho_low - rho_up - 2.0*tol_rho) ;
+ X_D = (X_D < 0.5) ? std::max(tol_l,X_D) : std::min(X_D,tol_u); //constrain to interval [tol_l;tol_u]
- X_H = 1.0 - X_D;
+ X_H = 1.0 - X_D;
- //calculate equilibrium
- // using the p_eq to calculate the T_eq - Clausius-Clapeyron
- T_eq = (reaction_enthalpy/R) / ((reaction_entropy/R) + std::log(p_r_g)); // unit of p in bar
- //Alternative: Use T_s as T_eq and calculate p_eq - for Schaube kinetics
- p_eq = exp((reaction_enthalpy/R)/T_s - (reaction_entropy/R));
+ //calculate equilibrium
+ // using the p_eq to calculate the T_eq - Clausius-Clapeyron
+ T_eq = (reaction_enthalpy/R) / ((reaction_entropy/R) + std::log(p_r_g)); // unit of p in bar
+ //Alternative: Use T_s as T_eq and calculate p_eq - for Schaube kinetics
+ p_eq = exp((reaction_enthalpy/R)/T_s - (reaction_entropy/R));
}
double ReactionCaOH2::Ca_hydration()
{
- double dXdt;
- // step 3, calculate dX/dt
+ double dXdt;
+ // step 3, calculate dX/dt
#ifdef SIMPLE_KINETICS
- if ( T_s < T_eq ) // hydration - simple model
+ if ( T_s < T_eq ) // hydration - simple model
#else
- if ( p_r_g > p_eq ) // hydration - Schaube model
+ if ( p_r_g > p_eq ) // hydration - Schaube model
#endif
- {
- //X_H = max(tol_l,X_H); //lower tolerance to avoid oscillations at onset of hydration reaction. Set here so that no residual reaction rate occurs at end of hydration.
+ {
+ //X_H = max(tol_l,X_H); //lower tolerance to avoid oscillations at onset of hydration reaction. Set here so that no residual reaction rate occurs at end of hydration.
#ifdef SIMPLE_KINETICS // this is from P. Schmidt
- dXdt = -1.0*(1.0-X_H) * (T_s - T_eq) / T_eq * 0.2 * conversion_rate::x_react;
+ dXdt = -1.0*(1.0-X_H) * (T_s - T_eq) / T_eq * 0.2 * conversion_rate::x_react;
#else //this is from Schaube
- if (X_H == tol_u || rho_s == rho_up)
- dXdt = 0.0;
- else if ( (T_eq-T_s) >= 50.0)
- dXdt = 13945.0 * exp(-89486.0/R/T_s) * std::pow(p_r_g/p_eq - 1.0,0.83) * 3.0 * (X_D) * std::pow(-1.0*log(X_D),0.666);
- else
- dXdt = 1.0004e-34 * exp(5.3332e4/T_s) * std::pow(p_r_g, 6.0) * (X_D);
+ if (X_H == tol_u || rho_s == rho_up)
+ dXdt = 0.0;
+ else if ( (T_eq-T_s) >= 50.0)
+ dXdt = 13945.0 * exp(-89486.0/R/T_s) * std::pow(p_r_g/p_eq - 1.0,0.83) * 3.0 * (X_D) * std::pow(-1.0*log(X_D),0.666);
+ else
+ dXdt = 1.0004e-34 * exp(5.3332e4/T_s) * std::pow(p_r_g, 6.0) * (X_D);
#endif
- }
- else // dehydration
- {
- //X_D = max(tol_l,X_D); //lower tolerance to avoid oscillations at onset of dehydration reaction. Set here so that no residual reaction rate occurs at end of dehydration.
+ }
+ else // dehydration
+ {
+ //X_D = max(tol_l,X_D); //lower tolerance to avoid oscillations at onset of dehydration reaction. Set here so that no residual reaction rate occurs at end of dehydration.
#ifdef SIMPLE_KINETICS // this is from P. Schmidt
- dXdt = -1.0* (1.0-X_D) * (T_s - T_eq) / T_eq * 0.05;
+ dXdt = -1.0* (1.0-X_D) * (T_s - T_eq) / T_eq * 0.05;
#else
- if (X_D == tol_u || rho_s == rho_low)
- dXdt = 0.0;
- else if (X_D < 0.2)
- dXdt = -1.9425e12 * exp( -1.8788e5/R/T_s ) * std::pow(1.0-p_r_g/p_eq,3.0)*(X_H);
- else
- dXdt = -8.9588e9 * exp( -1.6262e5/R/T_s ) * std::pow(1.0-p_r_g/p_eq,3.0)*2.0 * std::pow(X_H, 0.5);
+ if (X_D == tol_u || rho_s == rho_low)
+ dXdt = 0.0;
+ else if (X_D < 0.2)
+ dXdt = -1.9425e12 * exp( -1.8788e5/R/T_s ) * std::pow(1.0-p_r_g/p_eq,3.0)*(X_H);
+ else
+ dXdt = -8.9588e9 * exp( -1.6262e5/R/T_s ) * std::pow(1.0-p_r_g/p_eq,3.0)*2.0 * std::pow(X_H, 0.5);
#endif
- }
- return dXdt;
+ }
+ return dXdt;
}