Solution
Class Reference¶
This page contains detailed information on each of the methods, attributes, and properties in Solution
. Use the sidebar on the right for easier navigation.
Solution¶
- class pyEQL.Solution(solutes: list[list[str]] | dict[str, str] | None = None, volume: str | None = None, temperature: str = '298.15 K', pressure: str = '1 atm', pH: float = 7, pE: float = 8.5, balance_charge: str | None = None, solvent: str | list = 'H2O', engine: EOS | Literal['native', 'ideal', 'phreeqc'] = 'native', database: str | Path | Store | None = None)[source]¶
Class representing the properties of a solution. Instances of this class contain information about the solutes, solvent, and bulk properties.
- __init__(solutes: list[list[str]] | dict[str, str] | None = None, volume: str | None = None, temperature: str = '298.15 K', pressure: str = '1 atm', pH: float = 7, pE: float = 8.5, balance_charge: str | None = None, solvent: str | list = 'H2O', engine: EOS | Literal['native', 'ideal', 'phreeqc'] = 'native', database: str | Path | Store | None = None)[source]¶
- Parameters:
solutes –
dict, optional. Keys must be the chemical formula, while values must be str Quantity representing the amount. For example:
{“Na+”: “0.1 mol/L”, “Cl-”: “0.1 mol/L”}
Note that an older “list of lists” syntax is also supported; however this will be deprecated in the future and is no longer recommended. The equivalent list syntax for the above example is
[[“Na+”, “0.1 mol/L”], [“Cl-”, “0.1 mol/L”]]
Defaults to empty (pure solvent) if omitted
volume – str, optional Volume of the solvent, including the unit. Defaults to ‘1 L’ if omitted. Note that the total solution volume will be computed using partial molar volumes of the respective solutes as they are added to the solution.
temperature – str, optional The solution temperature, including the ureg. Defaults to ‘25 degC’ if omitted.
pressure – Quantity, optional The ambient pressure of the solution, including the unit. Defaults to ‘1 atm’ if omitted.
pH – number, optional Negative log of H+ activity. If omitted, the solution will be initialized to pH 7 (neutral) with appropriate quantities of H+ and OH- ions
pe – the pE value (redox potential) of the solution. Lower values = more reducing, higher values = more oxidizing. At pH 7, water is stable between approximately -7 to +14. The default value corresponds to a pE value typical of natural waters in equilibrium with the atmosphere.
balance_charge – The strategy for balancing charge during init and equilibrium calculations. Valid options are ‘pH’, which will adjust the solution pH to balance charge, ‘pE’ which will adjust the redox equilibrium to balance charge, or the name of a dissolved species e.g. ‘Ca+2’ or ‘Cl-’ that will be added/subtracted to balance charge. If set to None, no charge balancing will be performed either on init or when equilibrate() is called.
solvent – Formula of the solvent. Solvents other than water are not supported at this time.
engine – Electrolyte modeling engine to use. See documentation for details on the available engines.
database – path to a .json file (str or Path) or maggma Store instance that contains serialized SoluteDocs. None (default) will use the built-in pyEQL database.
Examples
>>> s1 = pyEQL.Solution({'Na+': '1 mol/L','Cl-': '1 mol/L'},temperature='20 degC',volume='500 mL') >>> print(s1) Components: Volume: 0.500 l Pressure: 1.000 atm Temperature: 293.150 K Components: ['H2O(aq)', 'H[+1]', 'OH[-1]', 'Na[+1]', 'Cl[-1]']
- property mass: Quantity¶
Return the total mass of the solution.
The mass is calculated each time this method is called. :param None: :param Returns: :param ——-: :param Quantity: :type Quantity: the mass of the solution, in kg
- property solvent_mass: Quantity¶
Return the mass of the solvent.
This property is used whenever mol/kg (or similar) concentrations are requested by get_amount()
Returns:¶
Quantity: the mass of the solvent, in kg
See Also:¶
- property volume: Quantity¶
Return the volume of the solution.
Returns:¶
Quantity: the volume of the solution, in L
- property temperature: Quantity¶
Return the temperature of the solution in Kelvin.
- property pressure: Quantity¶
Return the hydrostatic pressure of the solution in atm.
- p(solute: str, activity=True) float | None [source]¶
Return the negative log of the activity of solute.
Generally used for expressing concentration of hydrogen ions (pH)
- Parameters:
solute (str) – String representing the formula of the solute
activity (bool, optional) – If False, the function will use the molar concentration rather than the activity to calculate p. Defaults to True.
Returns –
------- –
Quantity – The negative log10 of the activity (or molar concentration if activity = False) of the solute.
- property density: Quantity¶
Return the density of the solution.
Density is calculated from the mass and volume each time this method is called.
Returns:¶
Quantity: The density of the solution.
- property dielectric_constant: Quantity¶
Returns the dielectric constant of the solution.
- Parameters:
None –
Returns –
------- –
Quantity (the dielectric constant of the solution, dimensionless.) –
Notes –
----- –
al. (Implements the following equation as given by Zuber et) –
math: (..) – epsilon = epsilon_{solvent} over 1 + sum_i alpha_i x_i:
:param where \(\alpha_i\) is a coefficient specific to the solvent and ion: :param and \(x_i\): :param is the mole fraction of the ion in solution.:
References:¶
.A. Zuber, L. Cardozo-Filho, V.F. Cabral, R.F. Checoni, M. Castier, An empirical equation for the dielectric constant in aqueous and nonaqueous electrolyte mixtures, Fluid Phase Equilib. 376 (2014) 116-123. doi:10.1016/j.fluid.2014.05.037.
- property chemical_system: str¶
Return the chemical system of the Solution as a “-” separated list of elements, sorted alphabetically. For example, a solution containing CaCO3 would have a chemical system of “C-Ca-H-O”.
- property elements: list¶
Return a list of elements that are present in the solution. For example, a solution containing CaCO3 would return [“C”, “Ca”, “H”, “O”]
- property cations: dict[str, float]¶
moles} that are cations. The returned dict is sorted by amount in descending order.
- Type:
Returns the subset of components {formula
- property anions: dict[str, float]¶
moles} that are anions. The returned dict is sorted by amount in descending order.
- Type:
Returns the subset of components {formula
- property neutrals: dict[str, float]¶
moles} that are neutral (not charged). The returned dict is sorted by amount in descending order.
- Type:
Returns the subset of components {formula
- property viscosity_dynamic: Quantity¶
Return the dynamic (absolute) viscosity of the solution.
Calculated from the kinematic viscosity
See Also:¶
viscosity_kinematic
- property viscosity_kinematic: Quantity¶
Return the kinematic viscosity of the solution.
Notes:¶
The calculation is based on a model derived from the Eyring equation and presented in
\[\ln \nu = \ln {\nu_w MW_w \over \sum_i x_i MW_i } + 15 x_+^2 + x_+^3 \delta G^*_{123} + 3 x_+ \delta G^*_{23} (1-0.05x_+)\]Where:
\[\delta G^*_{123} = a_o + a_1 (T)^{0.75}\]\[\delta G^*_{23} = b_o + b_1 (T)^{0.5}\]In which \(\nu\) is the kinematic viscosity, MW is the molecular weight, \(x_{+}\) is the mole fraction of cations, and \(T\) is the temperature in degrees C.
The a and b fitting parameters for a variety of common salts are included in the database.
References:¶
Vásquez-Castillo, G.; Iglesias-Silva, G. a.; Hall, K. R. An extension of the McAllister model to correlate kinematic viscosity of electrolyte solutions. Fluid Phase Equilib. 2013, 358, 44-49.
See Also:¶
- property conductivity: Quantity¶
Compute the electrical conductivity of the solution.
- Parameters:
None –
Returns –
------- –
Quantity – The electrical conductivity of the solution in Siemens / meter.
Notes –
----- –
respective (Conductivity is calculated by summing the molar conductivities of the) –
solutes –
exponent. (but they are activity-corrected and adjusted using an empricial) –
hc (This approach is used in PHREEQC and Aqion models aq) –
math:: (..) – EC = {F^2 over R T} sum_i D_i z_i ^ 2 gamma_i ^ {alpha} m_i
Where –
math:: –
alpha = begin{cases}
{frac{0.6}{sqrt{| z_{i} | }}} & {I < 0.36 | z_{i} | } {frac{sqrt{I}}{| z_i |}} & otherwise
end{cases}
Note (PHREEQC uses the molal rather than molar concentration according to) –
http (.. [hc]) –
References –
---------- –
https (.. [aq]) –
http –
Also (See) –
-------- –
:param
ionic_strength
: :paramget_molar_conductivity()
: :paramget_activity_coefficient()
:
- property ionic_strength: Quantity¶
Return the ionic strength of the solution.
Return the ionic strength of the solution, calculated as 1/2 * sum ( molality * charge ^2) over all the ions.
Molal (mol/kg) scale concentrations are used for compatibility with the activity correction formulas.
Returns:¶
- Quantity :
The ionic strength of the parent solution, mol/kg.
See Also:¶
Notes:¶
The ionic strength is calculated according to:
\[I = \sum_i m_i z_i^2\]Where \(m_i\) is the molal concentration and \(z_i\) is the charge on species i.
Examples:¶
>>> s1 = pyEQL.Solution([['Na+','0.2 mol/kg'],['Cl-','0.2 mol/kg']]) >>> s1.ionic_strength <Quantity(0.20000010029672785, 'mole / kilogram')>
>>> s1 = pyEQL.Solution([['Mg+2','0.3 mol/kg'],['Na+','0.1 mol/kg'],['Cl-','0.7 mol/kg']],temperature='30 degC') >>> s1.ionic_strength <Quantity(1.0000001004383303, 'mole / kilogram')>
- property charge_balance: float¶
Return the charge balance of the solution.
Return the charge balance of the solution. The charge balance represents the net electric charge on the solution and SHOULD equal zero at all times, but due to numerical errors will usually have a small nonzero value. It is calculated according to:
\[CB = \sum_i n_i z_i\]where \(n_i\) is the number of moles, and \(z_i\) is the charge on species i.
Returns:¶
- float :
The charge balance of the solution, in equivalents (mol of charge) per L.
- property alkalinity: Quantity¶
Return the alkalinity or acid neutralizing capacity of a solution.
Returns:¶
- Quantity :
The alkalinity of the solution in mg/L as CaCO3
Notes:¶
The alkalinity is calculated according to [stm]
\[Alk = \sum_{i} z_{i} C_{B} + \sum_{i} z_{i} C_{A}\]Where \(C_{B}\) and \(C_{A}\) are conservative cations and anions, respectively (i.e. ions that do not participate in acid-base reactions), and \(z_{i}\) is their signed charge. In this method, the set of conservative cations is all Group I and Group II cations, and the conservative anions are all the anions of strong acids.
References:¶
[stm]Stumm, Werner and Morgan, James J. Aquatic Chemistry, 3rd ed, pp 165. Wiley Interscience, 1996.
- property hardness: Quantity¶
Return the hardness of a solution.
Hardness is defined as the sum of the equivalent concentrations of multivalent cations as calcium carbonate.
NOTE: at present pyEQL cannot distinguish between mg/L as CaCO3 and mg/L units. Use with caution.
- Parameters:
None –
Returns –
------- –
Quantity – The hardness of the solution in mg/L as CaCO3
- property total_dissolved_solids: Quantity¶
Total dissolved solids in mg/L (equivalent to ppm) including both charged and uncharged species.
The TDS is defined as the sum of the concentrations of all aqueous solutes (not including the solvent), except for H[+1] and OH[-1]].
- property TDS: Quantity¶
Alias of
total_dissolved_solids()
- property debye_length: Quantity¶
Return the Debye length of a solution.
Debye length is calculated as [wk3]
\[\kappa^{-1} = \sqrt({\epsilon_r \epsilon_o k_B T \over (2 N_A e^2 I)})\]where \(I\) is the ionic strength, \(\epsilon_r\) and \(\epsilon_r\) are the relative permittivity and vacuum permittivity, \(k_B\) is the Boltzmann constant, and \(T\) is the temperature, \(e\) is the elementary charge, and \(N_A\) is Avogadro’s number.
Returns The Debye length, in nanometers.
References: .. [wk3] https://en.wikipedia.org/wiki/Debye_length#Debye_length_in_an_electrolyte
See also
- property bjerrum_length: Quantity¶
Return the Bjerrum length of a solution.
Bjerrum length represents the distance at which electrostatic interactions between particles become comparable in magnitude to the thermal energy.:math:lambda_B is calculated as
\[\lambda_B = {e^2 \over (4 \pi \epsilon_r \epsilon_o k_B T)}\]where \(e\) is the fundamental charge, \(\epsilon_r\) and \(\epsilon_r\) are the relative permittivity and vacuum permittivity, \(k_B\) is the Boltzmann constant, and \(T\) is the temperature.
- Parameters:
None –
Returns –
------- –
Quantity – The Bjerrum length, in nanometers.
References –
---------- –
https (//en.wikipedia.org/wiki/Bjerrum_length) –
Examples –
-------- –
pyEQL.Solution() (>>> s1 =) –
s1.bjerrum_length (>>>) –
<Quantity(0.7152793009386953 –
'nanometer')> –
Also (See) –
-------- –
:param
dielectric_constant
:
- property osmotic_pressure: Quantity¶
Return the osmotic pressure of the solution relative to pure water.
- Returns:
The osmotic pressure of the solution relative to pure water in Pa
See also
get_water_activity get_osmotic_coefficient get_salt
Notes
Osmotic pressure is calculated based on the water activity [sata] [wk]
\[\Pi = -\frac{RT}{V_{w}} \ln a_{w}\]Where \(\Pi\) is the osmotic pressure, \(V_{w}\) is the partial molar volume of water (18.2 cm**3/mol), and \(a_{w}\) is the water activity.
References
[sata]Sata, Toshikatsu. Ion Exchange Membranes: Preparation, Characterization, and Modification. Royal Society of Chemistry, 2004, p. 10.
Examples
>>> s1=pyEQL.Solution() >>> s1.osmotic_pressure <Quantity(0.495791416, 'pascal')>
>>> s1 = pyEQL.Solution([['Na+','0.2 mol/kg'],['Cl-','0.2 mol/kg']]) >>> soln.osmotic_pressure <Quantity(906516.7318131207, 'pascal')>
- get_amount(solute: str, units: str = 'mol/L') Quantity [source]¶
Return the amount of ‘solute’ in the parent solution.
The amount of a solute can be given in a variety of unit types. 1. substance per volume (e.g., ‘mol/L’, ‘M’) 2. equivalents (i.e., moles of charge) per volume (e.g., ‘eq/L’, ‘meq/L’) 3. substance per mass of solvent (e.g., ‘mol/kg’, ‘m’) 4. mass of substance (e.g., ‘kg’) 5. moles of substance (‘mol’) 6. mole fraction (‘fraction’) 7. percent by weight (%) 8. number of molecules (‘count’) 9. “parts-per-x” units, where ppm = mg/L, ppb = ug/L ppt = ng/L
- Parameters:
solute (str) – String representing the name of the solute of interest
units (in the specified) – Units desired for the output. Examples of valid units are ‘mol/L’,’mol/kg’,’mol’, ‘kg’, and ‘g/L’ Use ‘fraction’ to return the mole fraction. Use ‘%’ to return the mass percent
Returns –
------- –
question (The amount of the solute in) –
units –
See Also:¶
add_amount set_amount get_total_amount get_osmolarity get_osmolality get_mass get_total_moles_solute
- get_components_by_element() dict[str, list] [source]¶
Return a list of all species associated with a given element. Elements (keys) are suffixed with their oxidation state in parentheses, e.g.,
{“Na(1)”:[“Na[+1]”, “NaOH(aq)”]}
Species associated with each element are sorted in descending order of the amount present (i.e., the first species listed is the most abundant).
- get_el_amt_dict()[source]¶
Return a dict of Element: amount in mol
Elements (keys) are suffixed with their oxidation state in parentheses, e.g. “Fe(2)”, “Cl(-1)”.
- get_total_amount(element: str, units) Quantity [source]¶
Return the total amount of ‘element’ (across all solutes) in the solution.
- Parameters:
element – The symbol of the element of interest. The symbol can optionally be followed by the oxidation state in parentheses, e.g., “Na(1)”, “Fe(2)”, or “O(0)”. If no oxidation state is given, the total concentration of the element (over all oxidation states) is returned.
- unitsstr
Units desired for the output. Examples of valid units are ‘mol/L’,’mol/kg’,’mol’, ‘kg’, and ‘g/L’
Returns:¶
The total amount of the element in the solution, in the specified units
See also
- add_solute(formula: str, amount: str)[source]¶
Primary method for adding substances to a pyEQL solution.
- Parameters:
formula (str) – Chemical formula for the solute. Charged species must contain a + or - and (for polyvalent solutes) a number representing the net charge (e.g. ‘SO4-2’).
amount (str) – The amount of substance in the specified unit system. The string should contain both a quantity and a pint-compatible representation of a ureg. e.g. ‘5 mol/kg’ or ‘0.1 g/L’
- add_solvent(formula: str, amount: str)[source]¶
Same as add_solute but omits the need to pass solvent mass to pint.
- add_amount(solute: str, amount: str)[source]¶
Add the amount of ‘solute’ to the parent solution.
- Parameters:
solute (str) – String representing the name of the solute of interest
amount (str quantity) – String representing the concentration desired, e.g. ‘1 mol/kg’ If the units are given on a per-volume basis, the solution volume is not recalculated If the units are given on a mass, substance, per-mass, or per-substance basis, then the solution volume is recalculated based on the new composition
Returns –
------- –
modified. (Nothing. The concentration of solute is) –
- set_amount(solute: str, amount: str)[source]¶
Set the amount of ‘solute’ in the parent solution.
- Parameters:
solute (str) – String representing the name of the solute of interest
amount (str Quantity) –
String representing the concentration desired, e.g. ‘1 mol/kg’ If the units are given on a per-volume basis, the solution volume is not recalculated and the molar concentrations of other components in the solution are not altered, while the molal concentrations are modified.
If the units are given on a mass, substance, per-mass, or per-substance basis, then the solution volume is recalculated based on the new composition and the molal concentrations of other components are not altered, while the molar concentrations are modified.
Returns –
------- –
modified. (Nothing. The concentration of solute is) –
- get_moles_solvent() Quantity [source]¶
Return the moles of solvent present in the solution.
- Returns:
The moles of solvent in the solution.
- get_osmolarity(activity_correction=False) Quantity [source]¶
Return the osmolarity of the solution in Osm/L.
- Parameters:
activity_correction (bool) – If TRUE, the osmotic coefficient is used to calculate the osmolarity. This correction is appropriate when trying to predict the osmolarity that would be measured from e.g. freezing point depression. Defaults to FALSE if omitted.
- get_osmolality(activity_correction=False) Quantity [source]¶
Return the osmolality of the solution in Osm/kg.
- Parameters:
activity_correction (bool) – If TRUE, the osmotic coefficient is used to calculate the osmolarity. This correction is appropriate when trying to predict the osmolarity that would be measured from e.g. freezing point depression. Defaults to FALSE if omitted.
- get_salt() Salt [source]¶
Determine the predominant salt in a solution of ions.
Many empirical equations for solution properties such as activity coefficient, partial molar volume, or viscosity are based on the concentration of single salts (e.g., NaCl). When multiple ions are present (e.g., a solution containing Na+, Cl-, and Mg+2), it is generally not possible to directly model these quantities. pyEQL works around this problem by treating such solutions as single salt solutions.
The get_salt() method examines the ionic composition of a solution and returns an object that identifies the single most predominant salt in the solution, defined by the cation and anion with the highest mole fraction. The Salt object contains information about the stoichiometry of the salt to enable its effective concentration to be calculated (e.g., if a solution contains 0.5 mol/kg of Na+ and Cl-, plus traces of H+ and OH-, the matched salt is 0.5 mol/kg NaCl).
- Returns:
Salt object containing information about the parent salt.
See Also:¶
get_activity()
get_activity_coefficient()
get_water_activity()
get_osmotic_coefficient()
get_osmotic_pressure()
get_viscosity_kinematic()
Examples
>>> s1 = Solution([['Na+','0.5 mol/kg'],['Cl-','0.5 mol/kg']]) >>> s1.get_salt() <pyEQL.salt_ion_match.Salt object at 0x7fe6d3542048> >>> s1.get_salt().formula 'NaCl' >>> s1.get_salt().nu_cation 1 >>> s1.get_salt().z_anion -1
>>> s2 = pyEQL.Solution([['Na+','0.1 mol/kg'],['Mg+2','0.2 mol/kg'],['Cl-','0.5 mol/kg']]) >>> s2.get_salt().formula 'MgCl2' >>> s2.get_salt().nu_anion 2 >>> s2.get_salt().z_cation 2
- get_salt_dict(cutoff: float = 0.01, use_totals: bool = True) dict[str, dict] [source]¶
Returns a dict of salts that approximates the composition of the Solution. Like components, the dict is keyed by formula and the values are the total moles present in the solution, e.g., {“NaCl(aq)”: 1}. If the Solution is pure water, the returned dict contains only ‘HOH’.
- Parameters:
cutoff – Lowest salt concentration to consider. Analysis will stop once the concentrations of Salts being analyzed goes below this value. Useful for excluding analysis of trace anions.
use_totals – Whether to base the analysis on total element concentrations or individual species concentrations.
Notes
Salts are identified by pairing the predominant cations and anions in the solution, in descending order of their respective equivalent amounts.
Many empirical equations for solution properties such as activity coefficient, partial molar volume, or viscosity are based on the concentration of single salts (e.g., NaCl). When multiple ions are present (e.g., a solution containing Na+, Cl-, and Mg+2), it is generally not possible to directly model these quantities.
The get_salt_dict() method examines the ionic composition of a solution and simplifies it into a list of salts. The method returns a dictionary of Salt objects where the keys are the salt formulas (e.g., ‘NaCl’). The Salt object contains information about the stoichiometry of the salt to enable its effective concentration to be calculated (e.g., 1 M MgCl2 yields 1 M Mg+2 and 2 M Cl-).
- Parameters:
None –
Returns –
------- –
dict – A dictionary of Salt objects, keyed to the salt formula
Also (See) –
-------- –
:param
get_activity()
: :paramget_activity_coefficient()
: :paramget_water_activity()
: :paramget_osmotic_coefficient()
: :paramget_osmotic_pressure()
: :paramget_viscosity_kinematic()
:
- equilibrate(**kwargs) None [source]¶
Update the composition of the Solution using the thermodynamic engine. Any kwargs specified are passed through to self.engine.equilibrate()
- Returns:
Nothing. The .components attribute of the Solution is updated.
- get_activity_coefficient(solute: str, scale: Literal['molal', 'molar', 'fugacity', 'rational'] = 'molal') Quantity [source]¶
Return the activity coefficient of a solute in solution.
The model used to calculate the activity coefficient is determined by the Solution’s equation of state engine.
- Parameters:
solute – The solute for which to retrieve the activity coefficient
scale – The activity coefficient concentration scale
verbose – If True, pyEQL will print a message indicating the parent salt that is being used for activity calculations. This option is useful when modeling multicomponent solutions. False by default.
- Returns:
the activity coefficient as a dimensionless pint Quantity
- Return type:
Quantity
- get_activity(solute: str, scale: Literal['molal', 'molar', 'rational'] = 'molal') Quantity [source]¶
Return the thermodynamic activity of the solute in solution on the chosen concentration scale.
- Parameters:
solute – String representing the name of the solute of interest
scale – The concentration scale for the returned activity. Valid options are “molal”, “molar”, and “rational” (i.e., mole fraction). By default, the molal scale activity is returned.
verbose – If True, pyEQL will print a message indicating the parent salt that is being used for activity calculations. This option is useful when modeling multicomponent solutions. False by default.
- Returns:
The thermodynamic activity of the solute in question (dimensionless Quantity)
Notes
The thermodynamic activity depends on the concentration scale used [rs] . By default, the ionic strength, activity coefficients, and activities are all calculated based on the molal (mol/kg) concentration scale.
References
[rs]Robinson, R. A.; Stokes, R. H. Electrolyte Solutions: Second Revised Edition; Butterworths: London, 1968, p.32.
- get_osmotic_coefficient(scale: Literal['molal', 'molar', 'rational'] = 'molal') Quantity [source]¶
Return the osmotic coefficient of an aqueous solution.
The method used depends on the Solution object’s equation of state engine.
- get_water_activity() Quantity [source]¶
Return the water activity.
Returns:¶
- Quantity :
The thermodynamic activity of water in the solution.
See Also:¶
Notes:¶
Water activity is related to the osmotic coefficient in a solution containing i solutes by:
\[\ln a_{w} = - \Phi M_{w} \sum_{i} m_{i}\]Where \(M_{w}\) is the molar mass of water (0.018015 kg/mol) and \(m_{i}\) is the molal concentration of each species.
If appropriate Pitzer model parameters are not available, the water activity is assumed equal to the mole fraction of water.
References:¶
Blandamer, Mike J., Engberts, Jan B. F. N., Gleeson, Peter T., Reis, Joao Carlos R., 2005. “Activity of water in aqueous systems: A frequently neglected property.” Chemical Society Review 34, 440-458.
Examples:¶
>>> s1 = pyEQL.Solution([['Na+','0.3 mol/kg'],['Cl-','0.3 mol/kg']]) >>> s1.get_water_activity() <Quantity(0.9900944932888518, 'dimensionless')>
- get_chemical_potential_energy(activity_correction: bool = True) Quantity [source]¶
Return the total chemical potential energy of a solution (not including pressure or electric effects).
- Parameters:
activity_correction (bool, optional) – If True, activities will be used to calculate the true chemical potential. If False, mole fraction will be used, resulting in a calculation of the ideal chemical potential.
Returns –
------- –
Quantity – The actual or ideal chemical potential energy of the solution, in Joules.
Notes –
----- –
is (The chemical potential energy (related to the Gibbs mixing energy)) –
follows (calculated as) –
math: (..) – E = R T sum_i n_i ln a_i:
or –
math: – E = R T sum_i n_i ln x_i:
:param Where \(n\) is the number of moles of substance: :param \(T\) is the temperature in kelvin: :param : :param \(R\) the ideal gas constant: :param \(x\) the mole fraction: :param and \(a\) the activity of: :param each component.: :param Note that dissociated ions must be counted as separate components: :param : :param so a simple salt dissolved in water is a three component solution (cation: :param : :param anion: :param and water).: :param References: :param ———-: :param .. [koga] Koga: :type .. [koga] Koga: A differential approach.* Elsevier, 2007, pp. 23-37. :param Yoshikata: :type Yoshikata: A differential approach.* Elsevier, 2007, pp. 23-37. :param 2007. Solution Thermodynamics and its Application to Aqueous Solutions: :type 2007. *Solution Thermodynamics and its Application to Aqueous Solutions: A differential approach. Elsevier, 2007, pp. 23-37.
- _get_property(solute: str, name: str) Any | None [source]¶
Retrieve a thermodynamic property (such as diffusion coefficient) for solute, and adjust it from the reference conditions to the conditions of the solution.
- get_transport_number(solute, activity_correction=False) Quantity [source]¶
Calculate the transport number of the solute in the solution.
- Parameters:
solute – String identifying the solute for which the transport number is to be calculated.
activity_correction – If True, the transport number will be corrected for activity following the same method used for solution conductivity. Defaults to False if omitted.
- Returns:
The transport number of solute
Notes
Transport number is calculated according to :
\[t_i = {D_i z_i^2 C_i \over \sum D_i z_i^2 C_i}\]Where \(C_i\) is the concentration in mol/L, \(D_i\) is the diffusion coefficient, and \(z_i\) is the charge, and the summation extends over all species in the solution.
If activity_correction is True, the contribution of each ion to the transport number is corrected with an activity factor. See the documentation for Solution.conductivity for an explanation of this correction.
References
Geise, G. M.; Cassady, H. J.; Paul, D. R.; Logan, E.; Hickner, M. A. “Specific ion effects on membrane potential and the permselectivity of ion exchange membranes.”” Phys. Chem. Chem. Phys. 2014, 16, 21673-21681.
- get_molar_conductivity(solute: str) Quantity [source]¶
Calculate the molar (equivalent) conductivity for a solute.
- Parameters:
solute – String identifying the solute for which the molar conductivity is to be calculated.
- Returns:
The molar or equivalent conductivity of the species in the solution. Zero if the solute is not charged.
Notes
Molar conductivity is calculated from the Nernst-Einstein relation [smed]
\[\kappa_i = {z_i^2 D_i F^2 \over RT}\]Note that the diffusion coefficient is strongly variable with temperature.
References
[smed]Smedley, Stuart. The Interpretation of Ionic Conductivity in Liquids, pp 1-9. Plenum Press, 1980.
- get_mobility(solute: str) Quantity [source]¶
Calculate the ionic mobility of the solute.
- Parameters:
solute (str) – String identifying the solute for which the mobility is to be calculated.
Returns –
------- –
float (The ionic mobility. Zero if the solute is not charged.) –
Notes:¶
This function uses the Einstein relation to convert a diffusion coefficient into an ionic mobility [smed]
\[\mu_i = {F |z_i| D_i \over RT}\]References:¶
[smed]Smedley, Stuart I. The Interpretation of Ionic Conductivity in Liquids. Plenum Press, 1980.
- get_lattice_distance(solute: str) Quantity [source]¶
Calculate the average distance between molecules.
Calculate the average distance between molecules of the given solute, assuming that the molecules are uniformly distributed throughout the solution.
- Parameters:
solute (str) – String representing the name of the solute of interest
Returns –
------- –
Quantity (The average distance between solute molecules) –
Examples –
-------- –
Solution([['Na+' (>>> soln =) –
mol/kg'] ('0.5) –
['Cl-' –
mol/kg']]) ('0.5) –
soln.get_lattice_distance('Na+') (>>>) –
nanometer (1.492964....) –
Notes –
----- –
follows (The lattice distance is related to the molar concentration as) –
math: (..) – d = ( C_i N_A ) ^ {-{1 over 3}}:
- as_dict() dict [source]¶
Convert the Solution into a dict representation that can be serialized to .json or other format.
- classmethod from_dict(d: dict) Solution [source]¶
Instantiate a Solution from a dictionary generated by as_dict().
- print(mode: Literal['all', 'ions', 'cations', 'anions', 'neutrals'] = 'all', units: Literal['ppm', 'mol', 'mol/kg', 'mol/L', '%', 'activity'] = 'mol', places=4)[source]¶
Print details about the Solution.
- Parameters:
mode – Whether to list the amounts of all solutes, or only anions, cations, any ion, or any neutral solute.
units – The units to list solute amounts in. “activity” will list dimensionless activities instead of concentrations.
places – The number of decimal places to round the solute amounts.
- unsafe_hash()¶
Returns an hash of the current object. This uses a generic but low performance method of converting the object to a dictionary, flattening any nested keys, and then performing a hash on the resulting object
- classmethod validate_monty_v1(_MSONable__input_value)¶
Pydantic validator with correct signature for pydantic v1.x
- classmethod validate_monty_v2(_MSONable__input_value, _)¶
Pydantic validator with correct signature for pydantic v2.x