Internal Reference Documentation

Activity Correction API

pyEQL activity correction library

This file contains functions for computing molal-scale activity coefficients of ions and salts in aqueous solution.

Individual functions for activity coefficients are defined here so that they can be used independently of a pyEQL solution object. Normally, these functions are called from within the get_activity_coefficient method of the Solution class.

copyright:2013-2016 by Ryan S. Kingsbury
license:LGPL, see LICENSE for more details.
pyEQL.activity_correction._debye_parameter_B(temperature='25 degC')

Return the constant B used in the extended Debye-Huckel equation

Parameters:

temperature : str Quantity, optional

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

Notes

The parameter B is equal to: [1]

\[B = ( {8 \pi N_A e^2 \over 1000 \epsilon k T} ) ^ {1 \over 2}\]
[1]Bockris and Reddy. /Modern Electrochemistry/, vol 1. Plenum/Rosetta, 1977, p.210.

Examples

>>> _debye_parameter_B() 
0.3291...
pyEQL.activity_correction._debye_parameter_activity(temperature='25 degC')

Return the constant A for use in the Debye-Huckel limiting law (base 10)

Parameters:

temperature : str Quantity, optional

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

Returns:

Quantity The parameter A for use in the Debye-Huckel limiting law (base e)

Notes

The parameter A is equal to: [2]

\[A^{\gamma} = {e^3 ( 2 \pi N_A {\rho})^{0.5} \over (4 \pi \epsilon_o \epsilon_r k T)^{1.5}}\]

Note that this equation returns the parameter value that can be used to calculate the natural logarithm of the activity coefficient. For base 10, divide the value returned by 2.303. The value is often given in base 10 terms (0.509 at 25 degC) in older textbooks.

References

[2]Archer, Donald G. and Wang, Peiming. “The Dielectric Constant of Water and Debye-Huckel Limiting Law Slopes.” /J. Phys. Chem. Ref. Data/ 19(2), 1990.

Examples

>>> _debye_parameter_activity() 
1.17499...
pyEQL.activity_correction._debye_parameter_osmotic(temperature='25 degC')

Return the constant A_phi for use in calculating the osmotic coefficient according to Debye-Huckel theory

Parameters:

temperature : str Quantity, optional

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

Notes

Not to be confused with the Debye-Huckel constant used for activity coefficients in the limiting law. Takes the value 0.392 at 25 C. This constant is calculated according to: [3] [4]

\[A^{\phi} = {1 \over 3} A^{\gamma}\]

References

[3]Kim, Hee-Talk and Frederick, William Jr, 1988. “Evaluation of Pitzer Ion Interaction Parameters of Aqueous Electrolytes at 25 C. 1. Single Salt Parameters,” J. Chemical Engineering Data 33, pp.177-184.
[4]Archer, Donald G. and Wang, Peiming. “The Dielectric Constant of Water and Debye-Huckel Limiting Law Slopes.” /J. Phys. Chem. Ref. Data/ 19(2), 1990.

Examples

>>> _debye_parameter_osmotic() 
0.3916... 
pyEQL.activity_correction._debye_parameter_volume(temperature='25 degC')

Return the constant A_V, the Debye-Huckel limiting slope for apparent molar volume.

Parameters:

temperature : str Quantity, optional

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

Notes

Takes the value 1.8305 cm ** 3 * kg ** 0.5 / mol ** 1.5 at 25 C. This constant is calculated according to: [5]

\[A_V = -2 A_{\phi} R T [ {3 \over \epsilon} {{\partial \epsilon \over \partial p} } - {{1 \over \rho}{\partial \rho \over \partial p} }]\]

NOTE: at this time, the term in brackets (containing the partial derivatives) is approximate. These approximations give the correct value of the slope at 25 degC and produce estimates with less than 10% error between 0 and 60 degC.

The derivative of epsilon with respect to pressure is assumed constant (for atmospheric pressure) at -0.01275 1/MPa. Note that the negative sign does not make sense in light of real data, but is required to give the correct result.

The second term is equivalent to the inverse of the bulk modulus of water, which is taken to be 2.2 GPa. [6]

References

[5]Archer, Donald G. and Wang, Peiming. “The Dielectric Constant of Water and Debye-Huckel Limiting Law Slopes.” /J. Phys. Chem. Ref. Data/ 19(2), 1990.
[6]http://hyperphysics.phy-astr.gsu.edu/hbase/permot3.html

Examples

TODO

pyEQL.activity_correction._pitzer_B_MX(ionic_strength, alpha1, alpha2, beta0, beta1, beta2)

Return the B_MX coefficient for the Pitzer ion interaction model.

\[B_MX = \beta_0 + \beta_1 f1(\alpha_1 I ^ {0.5}) + \beta_2 f2(\alpha_2 I ^ {0.5})\]
Parameters:

ionic_strength: number

The ionic strength of the parent solution, mol/kg

alpha1, alpha2: number

Coefficients for the Pitzer model, kg ** 0.5 / mol ** 0.5

beta0, beta1, beta2: number

Coefficients for the Pitzer model. These ion-interaction parameters are specific to each salt system.

Returns:

float

The B_MX parameter for the Pitzer ion interaction model.

See also

_pitzer_f1

References

Scharge, T., Munoz, A.G., and Moog, H.C. (2012). Activity Coefficients of Fission Products in Highly Salinary Solutions of Na+, K+, Mg2+, Ca2+, Cl-, and SO42- : Cs+. /Journal of Chemical& Engineering Data (57), p. 1637-1647.

Kim, H., & Jr, W. F. (1988). Evaluation of Pitzer ion interaction parameters of aqueous electrolytes at 25 degree C. 1. Single salt parameters. Journal of Chemical and Engineering Data, (2), 177–184.

pyEQL.activity_correction._pitzer_B_phi(ionic_strength, alpha1, alpha2, beta0, beta1, beta2)

Return the B^Phi coefficient for the Pitzer ion interaction model.

\[B^\Phi = \beta_0 + \beta1 \exp(-\alpha_1 I ^{0.5}) + \beta_2 \exp(-\alpha_2 I ^ {0.5})\]

or

\[B^\Phi = B^\gamma - B_{MX}\]
Parameters:

ionic_strength: number

The ionic strength of the parent solution, mol/kg

alpha1, alpha2: number

Coefficients for the Pitzer model, kg ** 0.5 / mol ** 0.5

beta0, beta1, beta2: number

Coefficients for the Pitzer model. These ion-interaction parameters are specific to each salt system.

Returns:

float

The B^Phi parameter for the Pitzer ion interaction model.

References

Scharge, T., Munoz, A.G., and Moog, H.C. (2012). Activity Coefficients of Fission Products in Highly Salinary Solutions of Na+, K+, Mg2+, Ca2+, Cl-, and SO42- : Cs+. /Journal of Chemical& Engineering Data (57), p. 1637-1647.

Kim, H., & Jr, W. F. (1988). Evaluation of Pitzer ion interaction parameters of aqueous electrolytes at 25 degree C. 1. Single salt parameters. Journal of Chemical and Engineering Data, (2), 177–184.

Beyer, R., & Steiger, M. (2010). Vapor Pressure Measurements of NaHCOO + H 2 O and KHCOO + H 2 O from 278 to 308 K and Representation with an Ion Interaction (Pitzer) Model. Journal of Chemical & Engineering Data, 55(2), 830–838. doi:10.1021/je900487a

pyEQL.activity_correction._pitzer_f1(x)

The function of ionic strength used to calculate eta_MX in the Pitzer ion intercation model.

\[f(x) = 2 [ 1- (1+x) \exp(-x)] / x ^ 2\]

References

Scharge, T., Munoz, A.G., and Moog, H.C. (2012). Activity Coefficients of Fission Products in Highly Salinary Solutions of Na+, K+, Mg2+, Ca2+, Cl-, and SO42- : Cs+. /Journal of Chemical& Engineering Data (57), p. 1637-1647.

Kim, H., & Jr, W. F. (1988). Evaluation of Pitzer ion interaction parameters of aqueous electrolytes at 25 degree C. 1. Single salt parameters. Journal of Chemical and Engineering Data, (2), 177–184.

pyEQL.activity_correction._pitzer_f2(x)

The function of ionic strength used to calculate eta_gamma in the Pitzer ion intercation model.

\[f(x) = -{2 \over x ^ 2} [ 1 - ({1+x+ x^2 \over 2}) \exp(-x)] \]

References

Scharge, T., Munoz, A.G., and Moog, H.C. (2012). Activity Coefficients of Fission Products in Highly Salinary Solutions of Na+, K+, Mg2+, Ca2+, Cl-, and SO42- : Cs+. /Journal of Chemical& Engineering Data (57), p. 1637-1647.

Kim, H., & Jr, W. F. (1988). Evaluation of Pitzer ion interaction parameters of aqueous electrolytes at 25 degree C. 1. Single salt parameters. Journal of Chemical and Engineering Data, (2), 177–184.

pyEQL.activity_correction._pitzer_log_gamma(ionic_strength, molality, B_MX, B_phi, C_phi, z_cation, z_anion, nu_cation, nu_anion, temperature='25 degC', b=<Quantity(1.2, 'kilogram ** 0.5 / mole ** 0.5')>)

Return the natural logarithm of the binary activity coefficient calculated by the Pitzer ion interaction model.

\[\ln \gamma_{MX} = -{|z_+ z_-| A^{Phi} ( I ^ {0.5} \over (1 + b I ^ {0.5})} + {2 \over b }\ln (1 + b I ^ {0.5}) )+ + {m (2 \nu_+ \nu_-) \over (\nu_+ + \nu_-)} (B_{MX} + B_{MX}^\Phi) + {m^2(3 (\nu_+ \nu_-)^{1.5} \over (\nu_+ + \nu_-))} C_{MX}^\Phi \]
Parameters:

ionic_strength: Quantity

The ionic strength of the parent solution, mol/kg

molality: Quantity

The concentration of the salt, mol/kg

B_MX,B_phi,C_phi: Quantity

Calculated paramters for the Pitzer ion interaction model.

z_cation, z_anion: int

The formal charge on the cation and anion, respectively

nu_cation, nu_anion: int

The stoichiometric coefficient of the cation and anion in the salt

temperature: str Quantity

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

b: number, optional

Coefficient. Usually set equal to 1.2 kg ** 0.5 / mol ** 0.5 and considered independent of temperature and pressure

Returns:

float

The natural logarithm of the binary activity coefficient calculated by the Pitzer ion interaction model.

References

Kim, H., & Jr, W. F. (1988). Evaluation of Pitzer ion interaction parameters of aqueous electrolytes at 25 degree C. 1. Single salt parameters. Journal of Chemical and Engineering Data, (2), 177–184.

May, P. M., Rowland, D., Hefter, G., & Königsberger, E. (2011). A Generic and Updatable Pitzer Characterization of Aqueous Binary Electrolyte Solutions at 1 bar and 25 °C. Journal of Chemical & Engineering Data, 56(12), 5066–5077. doi:10.1021/je2009329

pyEQL.activity_correction.get_activity_coefficient_davies(ionic_strength, formal_charge=1, temperature='25 degC')

Return the activity coefficient of solute in the parent solution according to the Davies equation.

Parameters:

formal_charge : int, optional

The charge on the solute, including sign. Defaults to +1 if not specified.

ionic_strength : Quantity

The ionic strength of the parent solution, mol/kg

temperature : str Quantity, optional

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

Returns:

Quantity

The mean molal (mol/kg) scale ionic activity coefficient of solute, dimensionless.

Notes

Activity coefficient is calculated according to: [7]

\[\ln \gamma = A^{\gamma} z_i^2 ({\sqrt I \over (1 + \sqrt I)} + 0.2 I)\]

Valid for 0.1 < I < 0.5

References

[7]Stumm, Werner and Morgan, James J. Aquatic Chemistry, 3rd ed, pp 103. Wiley Interscience, 1996.
pyEQL.activity_correction.get_activity_coefficient_debyehuckel(ionic_strength, formal_charge=1, temperature='25 degC')

Return the activity coefficient of solute in the parent solution according to the Debye-Huckel limiting law.

Parameters:

formal_charge : int, optional

The charge on the solute, including sign. Defaults to +1 if not specified.

ionic_strength : Quantity

The ionic strength of the parent solution, mol/kg

temperature : str Quantity, optional

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

Returns:

Quantity

The mean molal (mol/kg) scale ionic activity coefficient of solute, dimensionless.

Notes

Activity coefficient is calculated according to: [8]

\[\ln \gamma = A^{\gamma} z_i^2 \sqrt I\]

Valid only for I < 0.005

References

[8]Stumm, Werner and Morgan, James J. Aquatic Chemistry, 3rd ed, pp 103. Wiley Interscience, 1996.
pyEQL.activity_correction.get_activity_coefficient_guntelberg(ionic_strength, formal_charge=1, temperature='25 degC')

Return the activity coefficient of solute in the parent solution according to the Guntelberg approximation.

Parameters:

formal_charge : int, optional

The charge on the solute, including sign. Defaults to +1 if not specified.

ionic_strength : Quantity

The ionic strength of the parent solution, mol/kg

temperature : str Quantity, optional

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

Returns:

Quantity

The mean molal (mol/kg) scale ionic activity coefficient of solute, dimensionless.

Notes

Activity coefficient is calculated according to: [9]

\[\ln \gamma = A^{\gamma} z_i^2 {\sqrt I \over (1 + \sqrt I)}\]

Valid for I < 0.1

References

[9]Stumm, Werner and Morgan, James J. Aquatic Chemistry, 3rd ed, pp 103. Wiley Interscience, 1996.
pyEQL.activity_correction.get_activity_coefficient_pitzer(ionic_strength, molality, alpha1, alpha2, beta0, beta1, beta2, C_phi, z_cation, z_anion, nu_cation, nu_anion, temperature='25 degC', b=1.2)

Return the activity coefficient of solute in the parent solution according to the Pitzer model.

Parameters:

ionic_strength: Quantity

The ionic strength of the parent solution, mol/kg

molality: Quantity

The molal concentration of the parent salt, mol/kg

alpha1, alpha2: number

Coefficients for the Pitzer model. This function assigns the coefficients proper units of kg ** 0.5 / mol ** 0.5 after they are entered.

beta0, beta1, beta2, C_phi: number

Coefficients for the Pitzer model. These ion-interaction parameters are specific to each salt system.

z_cation, z_anion: int

The formal charge on the cation and anion, respectively

nu_cation, nu_anion: int

The stoichiometric coefficient of the cation and anion in the salt

temperature: str Quantity

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

b: number, optional

Coefficient. Usually set equal to 1.2 and considered independent of temperature and pressure. If provided, this coefficient is assigned proper units of kg ** 0.5 / mol ** 0.5 after entry.

Returns:

Quantity

The mean molal (mol/kg) scale ionic activity coefficient of solute, dimensionless

References

Scharge, T., Munoz, A.G., and Moog, H.C. (2012). Activity Coefficients of Fission Products in Highly Salinary Solutions of Na+, K+, Mg2+, Ca2+, Cl-, and SO42- : Cs+. /Journal of Chemical& Engineering Data (57), p. 1637-1647.

Kim, H., & Jr, W. F. (1988). Evaluation of Pitzer ion interaction parameters of aqueous electrolytes at 25 degree C. 1. Single salt parameters. Journal of Chemical and Engineering Data, (2), 177–184.

May, P. M., Rowland, D., Hefter, G., & Königsberger, E. (2011). A Generic and Updatable Pitzer Characterization of Aqueous Binary Electrolyte Solutions at 1 bar and 25 °C. Journal of Chemical & Engineering Data, 56(12), 5066–5077. doi:10.1021/je2009329

Beyer, R., & Steiger, M. (2010). Vapor Pressure Measurements of NaHCOO + H 2 O and KHCOO + H 2 O from 278 to 308 K and Representation with an Ion Interaction (Pitzer) Model. Journal of Chemical & Engineering Data, 55(2), 830–838. doi:10.1021/je900487a

Examples

>>> get_activity_coefficient_pitzer(0.5*unit('mol/kg'),0.5*unit('mol/kg'),1,0.5,-.0181191983,-.4625822071,.4682,.000246063,1,-1,1,1,b=1.2)
0.61915...   
>>> get_activity_coefficient_pitzer(5.6153*unit('mol/kg'),5.6153*unit('mol/kg'),3,0.5,0.0369993,0.354664,0.0997513,-0.00171868,1,-1,1,1,b=1.2)
0.76331...

NOTE: the examples below are for comparison with experimental and modeling data presented in the May et al reference below.

10 mol/kg ammonium nitrate. Estimated result (from graph) = 0.2725

>>> get_activity_coefficient_pitzer(10*unit('mol/kg'),10*unit('mol/kg'),2,0,-0.01709,0.09198,0,0.000419,1,-1,1,1,b=1.2)
0.22595 ...

5 mol/kg ammonium nitrate. Estimated result (from graph) = 0.3011

>>> get_activity_coefficient_pitzer(5*unit('mol/kg'),5*unit('mol/kg'),2,0,-0.01709,0.09198,0,0.000419,1,-1,1,1,b=1.2)
0.30249 ...

18 mol/kg ammonium nitrate. Estimated result (from graph) = 0.1653

>>> get_activity_coefficient_pitzer(18*unit('mol/kg'),18*unit('mol/kg'),2,0,-0.01709,0.09198,0,0.000419,1,-1,1,1,b=1.2)
0.16241 ...
pyEQL.activity_correction.get_apparent_volume_pitzer(ionic_strength, molality, alpha1, alpha2, beta0, beta1, beta2, C_phi, V_o, z_cation, z_anion, nu_cation, nu_anion, temperature='25 degC', b=1.2)

Return the apparent molar volume of solute in the parent solution according to the Pitzer model.

Parameters:

ionic_strength: Quantity

The ionic strength of the parent solution, mol/kg

molality: Quantity

The molal concentration of the parent salt, mol/kg

alpha1, alpha2: number

Coefficients for the Pitzer model. This function assigns the coefficients proper units of kg ** 0.5 / mol ** 0.5 after they are entered.

beta0, beta1, beta2, C_phi: number

Pitzer coefficients for the apparent molar volume. These ion-interaction parameters are specific to each salt system.

V_o: number

The V^o Pitzer coefficient for the apparent molar volume.

z_cation, z_anion: int

The formal charge on the cation and anion, respectively

nu_cation, nu_anion: int

The stoichiometric coefficient of the cation and anion in the salt

temperature: str Quantity

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

b: number, optional

Coefficient. Usually set equal to 1.2 and considered independent of temperature and pressure. If provided, this coefficient is assigned proper units of kg ** 0.5 / mol ** 0.5 after entry.

Returns:

Quantity

The apparent molar volume of the solute, cm ** 3 / mol

References

May, P. M., Rowland, D., Hefter, G., & Königsberger, E. (2011). A Generic and Updatable Pitzer Characterization of Aqueous Binary Electrolyte Solutions at 1 bar and 25 °C. Journal of Chemical & Engineering Data, 56(12), 5066–5077. doi:10.1021/je2009329

Krumgalz, Boris S., Pogorelsky, Rita (1996). Volumetric Properties of Single Aqueous Electrolytes from Zero to Saturation Concentration at 298.15 K Represented by Pitzer’s Ion-Interaction Equations. Journal of Physical Chemical Reference Data, 25(2), 663-689.

Examples

NOTE: the example below is for comparison with experimental and modeling data presented in the Krumgalz et al reference below.

0.25 mol/kg CuSO4. Expected result (from graph) = 0.5 cm ** 3 / mol

>>> get_apparent_volume_pitzer(1.0*unit('mol/kg'),0.25*unit('mol/kg'),1.4,12,0.001499,-0.008124,0.2203,-0.0002589,-6,2,-2,1,1,b=1.2)
0.404...    

1.0 mol/kg CuSO4. Expected result (from graph) = 4 cm ** 3 / mol

>>> get_apparent_volume_pitzer(4.0*unit('mol/kg'),1.0*unit('mol/kg'),1.4,12,0.001499,-0.008124,0.2203,-0.0002589,-6,2,-2,1,1,b=1.2)
4.424...

10.0 mol/kg ammonium nitrate. Expected result (from graph) = 50.3 cm ** 3 / mol

>>> get_apparent_volume_pitzer(10.0*unit('mol/kg'),10.0*unit('mol/kg'),2,0,0.000001742,0.0002926,0,0.000000424,46.9,1,-1,1,1,b=1.2)
50.286...

20.0 mol/kg ammonium nitrate. Expected result (from graph) = 51.2 cm ** 3 / mol

>>> get_apparent_volume_pitzer(20.0*unit('mol/kg'),20.0*unit('mol/kg'),2,0,0.000001742,0.0002926,0,0.000000424,46.9,1,-1,1,1,b=1.2)
51.145...

NOTE: the examples below are for comparison with experimental and modeling data presented in the Krumgalz et al reference below.

0.8 mol/kg NaF. Expected result = 0.03

>>> get_apparent_volume_pitzer(0.8*unit('mol/kg'),0.8*unit('mol/kg'),2,0,0.000024693,0.00003169,0,-0.000004068,-2.426,1,-1,1,1,b=1.2)
0.22595 ...
pyEQL.activity_correction.get_osmotic_coefficient_pitzer(ionic_strength, molality, alpha1, alpha2, beta0, beta1, beta2, C_phi, z_cation, z_anion, nu_cation, nu_anion, temperature='25 degC', b=1.2)

Return the osmotic coefficient of water in an electrolyte solution according to the Pitzer model.

Parameters:

ionic_strength: Quantity

The ionic strength of the parent solution, mol/kg

molality: Quantity

The molal concentration of the parent salt, mol/kg

alpha1, alpha2: number

Coefficients for the Pitzer model. This function assigns the coefficients proper units of kg ** 0.5 / mol ** 0.5 after they are entered.

beta0, beta1, beta2, C_phi

Coefficients for the Pitzer model. These ion-interaction parameters are specific to each salt system.

z_cation, z_anion: int

The formal charge on the cation and anion, respectively

nu_cation, nu_anion: int

The stoichiometric coefficient of the cation and anion in the salt

temperature: str Quantity

String representing the temperature of the solution. Defaults to ‘25 degC’ if not specified.

b: number, optional

Coefficient. Usually set equal to 1.2 and considered independent of temperature and pressure. If provided, this coefficient is assigned proper units of kg ** 0.5 / mol ** 0.5 after entry.

Returns:

Quantity

The osmotic coefficient of water, dimensionless

References

Scharge, T., Munoz, A.G., and Moog, H.C. (2012). Activity Coefficients of Fission Products in Highly Salinary Solutions of Na+, K+, Mg2+, Ca2+, Cl-, and SO42- : Cs+. /Journal of Chemical& Engineering Data (57), p. 1637-1647.

Kim, H., & Jr, W. F. (1988). Evaluation of Pitzer ion interaction parameters of aqueous electrolytes at 25 degree C. 1. Single salt parameters. Journal of Chemical and Engineering Data, (2), 177–184.

May, P. M., Rowland, D., Hefter, G., & Königsberger, E. (2011). A Generic and Updatable Pitzer Characterization of Aqueous Binary Electrolyte Solutions at 1 bar and 25 °C. Journal of Chemical & Engineering Data, 56(12), 5066–5077. doi:10.1021/je2009329

Beyer, R., & Steiger, M. (2010). Vapor Pressure Measurements of NaHCOO + H 2 O and KHCOO + H 2 O from 278 to 308 K and Representation with an Ion Interaction (Pitzer) Model. Journal of Chemical & Engineering Data, 55(2), 830–838. doi:10.1021/je900487a

Examples

Experimental value according to Beyer and Stieger reference is 1.3550

>>> get_osmotic_coefficient_pitzer(10.175*unit('mol/kg'),10.175*unit('mol/kg'),1,0.5,-.0181191983,-.4625822071,.4682,.000246063,1,-1,1,1,b=1.2)
1.3552 ...

Experimental value according to Beyer and Stieger reference is 1.084

>>> get_osmotic_coefficient_pitzer(5.6153*unit('mol/kg'),5.6153*unit('mol/kg'),3,0.5,0.0369993,0.354664,0.0997513,-0.00171868,1,-1,1,1,b=1.2)
1.0850 ...  

NOTE: the examples below are for comparison with experimental and modeling data presented in the May et al reference below.

10 mol/kg ammonium nitrate. Estimated result (from graph) = 0.62

>>> get_osmotic_coefficient_pitzer(10*unit('mol/kg'),10*unit('mol/kg'),2,0,-0.01709,0.09198,0,0.000419,1,-1,1,1,b=1.2)
0.6143 ...

5 mol/kg ammonium nitrate. Estimated result (from graph) = 0.7

>>> get_osmotic_coefficient_pitzer(5*unit('mol/kg'),5*unit('mol/kg'),2,0,-0.01709,0.09198,0,0.000419,1,-1,1,1,b=1.2)
0.6925 ...

18 mol/kg ammonium nitrate. Estimated result (from graph) = 0.555

>>> get_osmotic_coefficient_pitzer(18*unit('mol/kg'),18*unit('mol/kg'),2,0,-0.01709,0.09198,0,0.000419,1,-1,1,1,b=1.2)
0.5556 ...

Water Properties API

pyEQL water properties library

This file contains functions for retrieving various physical properties of water substance

copyright:2013-2016 by Ryan S. Kingsbury
license:LGPL, see LICENSE for more details.
pyEQL.water_properties.water_density(temperature=<Quantity(25, 'degC')>, pressure=<Quantity(1, 'atmosphere')>)

Return the density of water in kg/m3 at the specified temperature and pressure.

Parameters:

temperature : float or int, optional

The temperature in Celsius. Defaults to 25 degrees if not specified.

pressure : float or int, optional

The ambient pressure of the solution in Pascals (N/m2). Defaults to atmospheric pressure (101325 Pa) if not specified.

Returns:

float

The density of water in kg/m3.

Notes

Based on the following empirical equation reported in [10]

\[\rho_W = 999.65 + 0.20438 T - 6.1744e-2 T ^ {1.5}\]

Where \(T\) is the temperature in Celsius.

[10]Sohnel, O and Novotny, P. Densities of Aqueous Solutions of Inorganic Substances. Elsevier Science, Amsterdam, 1985.

Examples

>>> water_density(25*unit('degC')) 
<Quantity(997.0415, 'kilogram / meter ** 3')>
pyEQL.water_properties.water_dielectric_constant(temperature=<Quantity(25, 'degC')>)

Return the dielectric constant of water at the specified temperature.

Parameters:

temperature : Quantity, optional

The temperature. Defaults to 25 degC if omitted.

Returns:

float

The dielectric constant (or permittivity) of water relative to the permittivity of a vacuum. Dimensionless.

Notes

This function implements a quadratic fit of measured permittivity data as reported in the CRC Handbook [11]. The parameters given are valid over the range 273 K to 372 K. Permittivity should not be extrapolated beyond this range.

\[\epsilon(T) = a + b T + c T^2\]

References

[11]“Permittivity (Dielectric Constant) of Liquids.” CRC Handbook of Chemistry and Physics, 92nd ed, pp 6-187 - 6-208.

Examples

>>> water_dielectric_constant(unit('20 degC')) 
80.15060...

Display an error if ‘temperature’ is outside the valid range

>>> water_dielectric_constant(-5*unit('degC'))
pyEQL.water_properties.water_specific_weight(temperature=<Quantity(25, 'degC')>, pressure=<Quantity(1, 'atmosphere')>)

Return the specific weight of water in N/m3 at the specified temperature and pressure.

Parameters:

temperature : Quantity, optional

The temperature. Defaults to 25 degC if omitted.

pressure : Quantity, optional

The ambient pressure of the solution. Defaults to atmospheric pressure (1 atm) if omitted.

Returns:

Quantity

The specific weight of water in N/m3.

See also

water_density

Examples

>>> water_specific_weight() 
<Quantity(9777.637025975, 'newton / meter ** 3')>
pyEQL.water_properties.water_viscosity_dynamic(temperature=<Quantity(25, 'degC')>, pressure=<Quantity(1, 'atmosphere')>)

Return the dynamic (absolute) viscosity of water in N-s/m2 = Pa-s = kg/m-s at the specified temperature.

Parameters:

temperature : Quantity, optional

The temperature. Defaults to 25 degC if omitted.

pressure : Quantity, optional

The ambient pressure of the solution. Defaults to atmospheric pressure (1 atm) if omitted.

Returns:

Quantity

The dynamic (absolute) viscosity of water in N-s/m2 = Pa-s = kg/m-s

Notes

Implements the international equation for viscosity of water as specified by NIST [12]

Valid for 273 < temperature < 1073 K and 0 < pressure < 100,000,000 Pa

References

[12]Sengers, J.V. “Representative Equations for the Viscosity of Water Substance.” J. Phys. Chem. Ref. Data 13(1), 1984.http://www.nist.gov/data/PDFfiles/jpcrd243.pdf

Examples

>>> water_viscosity_dynamic(20*unit('degC')) 
<Quantity(0.000998588610804179, 'kilogram / meter / second')>
>>> water_viscosity_dynamic(unit('100 degC'),unit('25 MPa')) 
<Quantity(0.00028165034364318573, 'kilogram / meter / second')>
>>> water_viscosity_dynamic(25*unit('degC'),0.1*unit('MPa')) 
<Quantity(0.0008872817880143659, 'kilogram / meter / second')>

#TODO - check these again after I implement pressure-dependent density function

pyEQL.water_properties.water_viscosity_kinematic(temperature=<Quantity(25, 'degC')>, pressure=<Quantity(1, 'atmosphere')>)

Return the kinematic viscosity of water in m2/s = Stokes at the specified temperature.

Parameters:

temperature : Quantity, optional

The temperature. Defaults to 25 degC if omitted.

pressure : Quantity, optional

The ambient pressure of the solution. Defaults to atmospheric pressure (1 atm) if omitted.

Returns:

Quantity

The kinematic viscosity of water in Stokes (m2/s)

Examples

>>> water_viscosity_kinematic()  
<Quantity(8.899146003595295e-07, 'meter ** 2 / second')>