René's URL Explorer Experiment

{"@context":"https://schema.org","@type":"QAPage","mainEntity":{"@type":"Question","name":"ePCSAFT example","text":"
Hello everyone,
\n
First of all I want to thank you for providing this amazing eos package!
\nCurrently, I try to model gas solubilities in aqueous solutions. Recently your group published this paper https://pubs.acs.org/doi/10.1021/acs.jced.4c00347, where you extended and use the ePCSAFT model and refer to feos.
\nDo you have a minimal example of how to use the ePCSAFT in feos?
\n
Thanks and best regards,
\n
Sven
","upvoteCount":1,"answerCount":1,"acceptedAnswer":{"@type":"Answer","text":"
Hi Sven,
\n
here is an example for MIAC of NaCl (modeled as two components) in water. The paths used assume that you work inside the \"examples\" folder. I cut it into three parts.
\n
1. Imports and helper functions
\n
Imports for FeOs, SI units and plotting and helper functions to calculate the mean ionic activity coefficient:
\n
from feos.eos import *\nfrom feos.epcsaft import *\nimport si_units as si\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\nsns.set_context('notebook')\nsns.set_palette('Dark2')\nsns.set_style('ticks')\n\ndef molality_to_mole_fraction(\n    molality_solute: si.SIObject, \n    molarweight_solvent: si.SIObject\n) -> float:\n    \"\"\"Calculate mole fraction of solute from solute molality given molar weight of solvent.\n\n    Args:\n        molality_solute: Molality of solute (salt).\n        molarweight_solvent: Molar weight of solvent.\n\n    Returns:\n        Mole fraction of solute (salt)\n    \"\"\"\n    return molality_solute * molarweight_solvent / (1.0 + molality_solute * molarweight_solvent)\n\ndef mean_ionic_activity_coefficient(\n    eos: EquationOfState,\n    temperature: si.SIObject, \n    pressure: si.SIObject, \n    molality_solute: si.SIObject,\n    molarweight_solvent: si.SIObject,\n    stoichiometric_coefficients: np.ndarray\n) -> float:\n    \"\"\"Calculate mean ionic activity coefficient (MIAC). \n\n    1. calculates fugacity coefficients for given molality and inf. dilution\n    2. uses eq. 9 to calculate activity coefficients from fugacities\n    3. uses eq. 8 to get MIAC\n    \n    Args:\n        eos: Equation of state to use\n        temperature: Temperature\n        pressure: Pressure\n        molality_solute: Molality of solute (salt).\n        molarweight_solvent: Molar weight of solvent.\n        stoichiometric_coefficients: \n            Stoichiometric coefficients of salt ions.\n            Order has to match parameter order.\n\n    Returns:\n        Mean ionic activity coefficient.\n    \"\"\"\n    # mole fraction of salt from molality\n    x_ions = molality_to_mole_fraction(molality_solute=molality_solute, molarweight_solvent=molarweight_solvent)\n    # mole fractions from closure\n    molefracs = np.array([1 - x_ions, x_ions / 2, x_ions / 2])\n\n    # thermodynamic condition\n    state = State(eos, temperature=temperature, pressure=pressure, molefracs=molefracs)\n    # reference: infinite dilution (same for each ion so only calculated once)\n    infinite_dilution = State(eos, temperature=temperature, pressure=pressure, molefracs=np.array([1.0, 0.0, 0.0]))\n\n    # Equation 9\n    # only use ion activity coefficients (component 0 is solvent)\n    ac_ions = (np.exp(state.ln_phi()) / np.exp(infinite_dilution.ln_phi()))[1:]\n    # Equation 8\n    miac = np.prod(ac_ions**stoichiometric_coefficients)**(1 / stoichiometric_coefficients.sum())\n    return miac
\n
2. Initializing the substances, equation of state and thermodynamic conditions
\n
For this example and the above functions, we assume that the solvent is always component at index 0 and the two ions are 1 and 2.
\nIf you initialize your substances in a different order, you have to adjust the code accordingly.
\n
# parameters\nparameters = ElectrolytePcSaftParameters.from_json(\n    substances=['water', 'sodium ion', 'chloride ion'],\n    pure_path='../parameters/epcsaft/held2014_w_permittivity_added.json',\n    binary_path='../parameters/epcsaft/held2014_binary.json'\n)\nmolarweight_water = parameters.pure_records[0].molarweight * si.GRAM / si.MOL\n\n# equation of state\neos = EquationOfState.epcsaft(parameters, epcsaft_variant=ElectrolytePcSaftVariants.Advanced)\n\n# thermodynamic conditions\ntemperature = 298.15 * si.KELVIN\npressure = si.BAR\n# stoichiometric coefficients for [Na+] and [Cl-]\nnu = np.array([1.0, 1.0])
\n
3. Calculating the MIAC's and plotting the results
\n
# mean ionic activity coefficients\nmolality_salt = si.linspace(0 * si.MOL / si.KILOGRAM, 6.0 * si.MOL / si.KILOGRAM, 100) \ngammas = np.array([\n    mean_ionic_activity_coefficient(eos, temperature, pressure, b, molarweight_solvent=molarweight_water, stoichiometric_coefficients=nu) \n    for b in molality_salt\n])\n\n# plot results\nplt.plot(molality_salt / si.MOL * si.KILOGRAM, gammas, clip_on=False)\nplt.xlabel(r'Molality $m_\\text{salt}$ / mol / kg')\nplt.ylabel(r'$\\gamma^{*, x}_\\pm$')\nplt.xlim(0, 6)\nplt.ylim(0.7, 1.2)\nsns.despine(offset=10);
\n
I hope that helps - let us know if there are any more questions.
","upvoteCount":1,"url":"https://github.com/feos-org/feos/discussions/270#discussioncomment-12083231"}}}