Skip to content

Commit

Permalink
WIP db refresh
Browse files Browse the repository at this point in the history
  • Loading branch information
rkingsbury committed Aug 10, 2023
1 parent a682aa6 commit 29f999d
Show file tree
Hide file tree
Showing 6 changed files with 6,951 additions and 62 deletions.
1 change: 1 addition & 0 deletions pyproject.toml
Original file line number Diff line number Diff line change
Expand Up @@ -86,3 +86,4 @@ no_implicit_optional = false
[tool.codespell]
ignore-words-list = "nd"
skip = 'tests/test_files/*'
exclude = 'src/pyEQL/database/pyeql_db.json'
2 changes: 1 addition & 1 deletion src/pyEQL/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -47,7 +47,7 @@
unit.default_format = "P~"

# this must be imported after instantiating the UnitRegistry()
from pyEQL.database import Paramsdb # noqa
# from pyEQL.database import Paramsdb

# initialize the parameters database
# paramsDB = Paramsdb()
Expand Down
6,843 changes: 6,842 additions & 1 deletion src/pyEQL/database/pyeql_db.json

Large diffs are not rendered by default.

16 changes: 8 additions & 8 deletions src/pyEQL/engines.py
Original file line number Diff line number Diff line change
Expand Up @@ -222,7 +222,7 @@ def get_activity_coefficient(self, solution, solute):
if verbose is True:
print("Calculating activity coefficient based on parent salt %s" % Salt.formula)

param = db.get_parameter(Salt.formula, "pitzer_parameters_activity")
param = db.get_property(Salt.formula, "pitzer_parameters_activity")

# determine alpha1 and alpha2 based on the type of salt
# see the May reference for the rules used to determine
Expand Down Expand Up @@ -438,7 +438,7 @@ def get_osmotic_coefficient(self, solution):
db.search_parameters(item.formula)

if db.has_parameter(item.formula, "pitzer_parameters_activity"):
param = db.get_parameter(item.formula, "pitzer_parameters_activity")
param = db.get_property(item.formula, "pitzer_parameters_activity")

osmotic_coefficient = ac.get_osmotic_coefficient_pitzer(
ionic_strength,
Expand Down Expand Up @@ -478,17 +478,16 @@ def get_solute_volume(self, solution):
Salt = solution.get_salt()

# search the database for pitzer parameters for 'salt'
db.search_parameters(Salt.formula)
# db.search_parameters(Salt.formula)

solute_vol = 0 * unit("L")

# use the pitzer approach if parameters are available

pitzer_calc = False

if db.has_parameter(Salt.formula, "pitzer_parameters_volume"):
param = db.get_parameter(Salt.formula, "pitzer_parameters_volume")

param = solution.get_property(Salt.formula, "pitzer_parameters_volume")
if param is not None:
# determine the average molality of the salt
# this is necessary for solutions inside e.g. an ion exchange
# membrane, where the cation and anion concentrations may be
Expand Down Expand Up @@ -550,8 +549,9 @@ def get_solute_volume(self, solution):
if pitzer_calc is True and solute in [Salt.anion, Salt.cation]:
continue

if db.has_parameter(solute, "partial_molar_volume"):
solute_vol += solution.get_parameter(solute, "partial_molar_volume") * mol
part_vol = solution.get_property(solute, "partial_molar_volume")
if part_vol is not None:
solute_vol += part_vol * mol
logger.info("Updated solution volume using direct partial molar volume for solute %s" % solute)

else:
Expand Down
111 changes: 59 additions & 52 deletions src/pyEQL/solution.py
Original file line number Diff line number Diff line change
Expand Up @@ -127,7 +127,7 @@ def __init__(
json = Path(database_dir) / "pyeql_db.json"
else:
json = database if isinstance(database, str) else str(database)
db_store = JSONStore(json)
db_store = JSONStore(json, key="formula")
logger.info(f"Created maggma JSONStore from .json file {database}")
else:
db_store = database
Expand Down Expand Up @@ -207,7 +207,7 @@ def add_solute(self, formula, amount):

# add the new solute
quantity = unit.Quantity(amount)
mw = unit.Quantity(self.get_property(formula, "molecular_weight"))
mw = self.get_property(formula, "molecular_weight") # returns a quantity
target_mol = quantity.to("moles", "chem", mw=mw, volume=self.volume, solvent_mass=self.get_solvent_mass())
self.components[formula] = target_mol.to("moles").magnitude

Expand All @@ -222,7 +222,10 @@ def add_solute(self, formula, amount):
target_mass = target_vol.to("L").magnitude * self.water_substance.rho * unit.Quantity("1 g")
# mw = unit.Quantity(self.get_property(self.solvent_name, "molecular_weight"))
mw = self.get_property(self.solvent_name, "molecular_weight")
print(target_mass, mw)
if mw is None:
raise ValueError(
f"Molecular weight for solvent {self.solvent_name} not found in database. This is required to proceed."
)
target_mol = target_mass.to("g") / mw.to("g/mol")
self.components[self.solvent_name] = target_mol.magnitude

Expand Down Expand Up @@ -454,15 +457,16 @@ def dielectric_constant(self) -> Quantity:
# ignore water
if item != "H2O":
# skip over solutes that don't have parameters
try:
fraction = self.get_amount(item, "fraction")
coefficient = self.get_property(item, "dielectric_parameter_water")
# try:
fraction = self.get_amount(item, "fraction")
coefficient = self.get_property(item, "model_parameters.dielectric_zuber")
if coefficient is not None:
denominator += coefficient * fraction
except TypeError:
logger.warning("No dielectric parameters found for species %s." % item)
continue
# except TypeError:
# logger.warning("No dielectric parameters found for species %s." % item)
# continue

return di_water / denominator
return unit.Quantity(di_water / denominator, "dimensionless")

# TODO - need tests for viscosity
@property
Expand Down Expand Up @@ -1040,7 +1044,6 @@ def get_amount(self, solute, units):
# retrieve the number of moles of solute and its molecular weight
try:
moles = unit.Quantity(self.components[solute], "mol")
mw = unit.Quantity(self.get_property(solute, "molecular_weight"), "g/mol")
# if the solute is not present in the solution, we'll get a KeyError
# In that case, the amount is zero
except KeyError:
Expand All @@ -1056,6 +1059,7 @@ def get_amount(self, solute, units):
# function calls.
if units == "fraction":
return moles / (self.get_moles_solvent() + self.get_total_moles_solute())
mw = self.get_property(solute, "molecular_weight").to("g/mol")
if units == "%":
return moles.to("kg", "chem", mw=mw) / self.mass.to("kg") * 100
if unit.Quantity(units).dimensionality in (
Expand Down Expand Up @@ -1324,15 +1328,15 @@ def get_osmolality(self, activity_correction=False):
factor = self.get_osmotic_coefficient() if activity_correction is True else 1
return factor * self.get_total_moles_solute() / self.get_solvent_mass().to("kg")

def get_total_moles_solute(self):
def get_total_moles_solute(self) -> Quantity:
"""Return the total moles of all solute in the solution."""
tot_mol = 0
for item in self.components:
if item != self.solvent_name:
tot_mol += self.components[item].moles
return tot_mol
tot_mol += self.components[item]
return unit.Quantity(tot_mol, "mol")

def get_moles_solvent(self):
def get_moles_solvent(self) -> Quantity:
"""
Return the moles of solvent present in the solution.
Expand Down Expand Up @@ -1812,7 +1816,7 @@ def get_mobility(self, solute):

return mobility.to("m**2/V/s")

def _get_property(self, solute, name):
def _get_property(self, solute: str, name: str) -> Optional[Quantity]:
"""Retrieve a thermodynamic property (such as diffusion coefficient)
for solute, and adjust it from the reference conditions to the conditions
of the solution.
Expand All @@ -1827,45 +1831,47 @@ def _get_property(self, solute, name):
Returns:
-------
Quantity: The desired parameter
Quantity: The desired parameter or None if not found
"""
# retrieve the base value and the conditions of measurement from the
# database

# TODO - replace with a Store query to the database.
base_value = None
rform = Ion.from_formula(solute).reduced_formula
if self.database.get(rform):
base_value = self.database[rform].get(name)
# base_value = self.get_property(solute, name) if db.has_parameter(solute, name) else None

base_temperature = unit.Quantity("25 degC")
# base_pressure = unit.Quantity("1 atm")

# query the database using the sanitized formula
rform = Ion.from_formula(solute).reduced_formula
data = list(self.database.query({"formula": rform, name: {"$exists": True}}, properties=["formula", name]))
# formulas should always be unique in the database. len==0 indicates no
# data. len>1 indicates duplicate data.
if len(data) != 1:
print(len(data))
logger.warning(f"Property {name} for solute {solute} not found in database. Returning None.")
return None
data = data[0]

# perform temperature-corrections or other adjustments for certain
# parameter types
if name == "diffusion_coefficient":
if base_value is not None:
# correct for temperature and viscosity
# .. math:: D_1 \over D_2 = T_1 \over T_2 * \mu_2 \over \mu_1
# where :math:`\mu` is the dynamic viscosity
# assume that the base viscosity is that of pure water
return (
base_value
* self.temperature
/ base_temperature
* self.water_substance.mu
* unit.Quantity("1 Pa*s")
/ self.get_viscosity_dynamic()
)
if name == "transport.diffusion_coefficient":
base_value = data.get(name)["value"]

# correct for temperature and viscosity
# .. math:: D_1 \over D_2 = T_1 \over T_2 * \mu_2 \over \mu_1
# where :math:`\mu` is the dynamic viscosity
# assume that the base viscosity is that of pure water
return (
base_value
* self.temperature
/ base_temperature
* self.water_substance.mu
* unit.Quantity("1 Pa*s")
/ self.get_viscosity_dynamic()
)

logger.warning("Diffusion coefficient not found for species %s. Assuming zero." % (solute))
return unit.Quantity("0 m**2/s")
# logger.warning("Diffusion coefficient not found for species %s. Assuming zero." % (solute))
# return unit.Quantity("0 m**2/s")

# just return the base-value molar volume for now; find a way to adjust for
# concentration later
if name == "partial_molar_volume":
if name == "size.molar_volume":
# calculate the partial molar volume for water since it isn't in the database
if solute == "H2O":
vol = (
Expand All @@ -1876,17 +1882,18 @@ def _get_property(self, solute, name):

return vol.to("cm **3 / mol")

if base_value is not None:
if self.temperature != base_temperature:
logger.warning("Partial molar volume for species %s not corrected for temperature" % solute)
return base_value
base_value = unit.Quantity(data["size"]["molar_volume"]["value"])
if self.temperature != base_temperature:
logger.warning("Partial molar volume for species %s not corrected for temperature" % solute)
return base_value

logger.warning("Partial molar volume not found for species %s. Assuming zero." % solute)
return unit.Quantity("0 cm **3 / mol")
if name == "model_parameters.dielectric_zuber":
return unit.Quantity(data["model_parameters"]["dielectric_zuber"]["value"])

# for parameters not named above, just return the base value
logger.warning("%s has not been corrected for solution conditions" % name)
return base_value
val = data.get(name) if isinstance(data.get(name), str) else data[name].get("value")
# logger.warning("%s has not been corrected for solution conditions" % name)
return unit.Quantity(val)

def get_chemical_potential_energy(self, activity_correction=True):
"""
Expand Down
40 changes: 40 additions & 0 deletions tests/test_solution.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,40 @@
"""
pyEQL volume and concentration methods test suite
=================================================
This file contains tests for the volume and concentration-related methods
used by pyEQL's Solution class
"""

import pytest
from pyEQL import Solution


@pytest.fixture()
def s1():
return Solution(volume="2 L")


@pytest.fixture()
def s2():
return Solution([["Na+", "4 mol/L"], ["Cl-", "4 mol/L"]], volume="2 L")


@pytest.fixture()
def s3():
return Solution([["Na+", "4 mol/kg"], ["Cl-", "4 mol/kg"]], volume="2 L")


@pytest.fixture()
def s4():
return Solution([["Na+", "8 mol"], ["Cl-", "8 mol"]], volume="2 L")


def test_db_connect():
Solution(database=None)


def test_serialization(s1, s2, s3):
assert isinstance(s1.as_dict(), dict)
s1_new = Solution.from_dict(s1.as_dict())
assert s1_new.volume.magnitude == 2

0 comments on commit 29f999d

Please sign in to comment.