From 819453a39bef93fb058f64a10f4c45b97c5b0169 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 13 Feb 2024 10:10:35 -0500 Subject: [PATCH 01/15] Adding volume_mass as unit and renaming molar_volume --- idaes/core/base/property_meta.py | 15 ++++++++++++++- idaes/core/base/property_set.py | 9 +++++++-- idaes/core/base/tests/test_property_meta.py | 2 ++ 3 files changed, 23 insertions(+), 3 deletions(-) diff --git a/idaes/core/base/property_meta.py b/idaes/core/base/property_meta.py index b12aee115c..6c44031f66 100644 --- a/idaes/core/base/property_meta.py +++ b/idaes/core/base/property_meta.py @@ -257,9 +257,22 @@ def VOLUME(self): return self._length**3 @property - def MOLAR_VOLUME(self): + def VOLUME_MASS(self): + return self._length**3 * self._mass**-1 + + @property + def VOLUME_MOL(self): return self._length**3 * self._amount**-1 + # Backward compatibility name + @property + def MOLAR_VOLUME(self): + msg = ( + f"The unit name MOLAR_VOLUME is being deprecated in " "favor of VOLUME_MOL." + ) + deprecation_warning(msg=msg, logger=_log, version="2.3.0", remove_in="3.0.0") + return self.VOLUME_MOL + # Flows @property def FLOW_MASS(self): diff --git a/idaes/core/base/property_set.py b/idaes/core/base/property_set.py index 70af3f584f..a98167b472 100644 --- a/idaes/core/base/property_set.py +++ b/idaes/core/base/property_set.py @@ -627,6 +627,11 @@ class StandardPropertySet(PropertySetBase): doc="Compressiblity Factor", units=pyunits.dimensionless, ) + compress_fact_crit = PropertyMetadata( + name="compress_fact_crit", + doc="Compressiblity Factor at Critical Point", + units=pyunits.dimensionless, + ) conc_mass = PropertyMetadata( name="conc_mass", doc="Concentration on a Mass Basis", @@ -890,12 +895,12 @@ class StandardPropertySet(PropertySetBase): vol_mol = PropertyMetadata( name="vol_mol", doc="Molar Volume", - units="MOLAR_VOLUME", + units="VOLUME_MOL", ) vol_mol_crit = PropertyMetadata( name="vol_mol", doc="Molar Volume at Critical Point", - units="MOLAR_VOLUME", + units="VOLUME_MOL", ) # Log terms log_act = PropertyMetadata( diff --git a/idaes/core/base/tests/test_property_meta.py b/idaes/core/base/tests/test_property_meta.py index 4e4e30b6d2..e67053b474 100644 --- a/idaes/core/base/tests/test_property_meta.py +++ b/idaes/core/base/tests/test_property_meta.py @@ -175,6 +175,8 @@ def test_luminous_intensity(unit_set): derived_quantities = { "area": units.m**2, "volume": units.m**3, + "volume_mass": units.m**3 / units.kg, + "volume_mol": units.m**3 / units.mol, "flow_mass": units.kg * units.s**-1, "flow_mole": units.mol * units.s**-1, "flow_vol": units.m**3 * units.s**-1, From 9e3e0983e6fcb0b88d45c07ec7c4bf5232e714ff Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 13 Feb 2024 10:54:23 -0500 Subject: [PATCH 02/15] Hooks for mixture critical properties --- idaes/core/base/components.py | 1 + idaes/core/base/tests/test_components.py | 8 + .../base/generic_property.py | 64 ++++++++ .../base/tests/test_generic_property.py | 147 ++++++++++++++++++ 4 files changed, 220 insertions(+) diff --git a/idaes/core/base/components.py b/idaes/core/base/components.py index f0e8baa507..5f9bcbe091 100644 --- a/idaes/core/base/components.py +++ b/idaes/core/base/components.py @@ -246,6 +246,7 @@ def build(self): # Create Vars for common parameters var_dict = { + "compress_fact_crit": pyunits.dimensionless, "dens_mol_crit": base_units.DENSITY_MOLE, "omega": pyunits.dimensionless, "pressure_crit": base_units.PRESSURE, diff --git a/idaes/core/base/tests/test_components.py b/idaes/core/base/tests/test_components.py index a0fda6ce6a..31d378d24f 100644 --- a/idaes/core/base/tests/test_components.py +++ b/idaes/core/base/tests/test_components.py @@ -185,6 +185,8 @@ def get_metadata(self): m.comp = Component( parameter_data={ "mw": 10, + "compress_fact_crit": 1, + "dens_mol_crit": 55, "pressure_crit": 1e5, "temperature_crit": 500, } @@ -193,6 +195,12 @@ def get_metadata(self): assert isinstance(m.comp.mw, Param) assert m.comp.mw.value == 10 + assert isinstance(m.comp.compress_fact_crit, Var) + assert m.comp.compress_fact_crit.value == 1 + + assert isinstance(m.comp.dens_mol_crit, Var) + assert m.comp.dens_mol_crit.value == 55 + assert isinstance(m.comp.pressure_crit, Var) assert m.comp.pressure_crit.value == 1e5 diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index f90609b0fb..6ea12736ad 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -1073,6 +1073,7 @@ def define_metadata(cls, obj): "method": "_act_coeff_phase_comp_apparent" }, "compress_fact_phase": {"method": "_compress_fact_phase"}, + "compress_fact_crit": {"method": "_critical_props"}, "conc_mol_comp": {"method": "_conc_mol_comp"}, "conc_mol_phase_comp": {"method": "_conc_mol_phase_comp"}, "conc_mol_phase_comp_apparent": { @@ -1094,6 +1095,7 @@ def define_metadata(cls, obj): "dens_mass": {"method": "_dens_mass"}, "dens_mass_phase": {"method": "_dens_mass_phase"}, "dens_mol": {"method": "_dens_mol"}, + "dens_mol_crit": {"method": "_critical_props"}, "dens_mol_phase": {"method": "_dens_mol_phase"}, "energy_internal_mol": {"method": "_energy_internal_mol"}, "energy_internal_mol_phase": {"method": "_energy_internal_mol_phase"}, @@ -1132,6 +1134,7 @@ def define_metadata(cls, obj): "mw_comp": {"method": "_mw_comp"}, "mw_phase": {"method": "_mw_phase"}, "prandtl_number_phase": {"method": "_prandtl_number_phase"}, + "pressure_crit": {"method": "_critical_props"}, "pressure_phase_comp": {"method": "_pressure_phase_comp"}, "pressure_phase_comp_true": {"method": "_pressure_phase_comp_true"}, "pressure_phase_comp_apparent": { @@ -1142,6 +1145,7 @@ def define_metadata(cls, obj): "pressure_osm_phase": {"method": "_pressure_osm_phase"}, "pressure_sat_comp": {"method": "_pressure_sat_comp"}, "surf_tens_phase": {"method": "_surf_tens_phase"}, + "temperature_crit": {"method": "_critical_props"}, "temperature_bubble": {"method": "_temperature_bubble"}, "temperature_dew": {"method": "_temperature_dew"}, "therm_cond_phase": {"method": "_therm_cond_phase"}, @@ -2920,6 +2924,66 @@ def get_mole_frac(self, phase=None): return self.mole_frac_phase_comp_apparent return self.mole_frac_phase_comp_true + # ------------------------------------------------------------------------- + # Mixture critical point + # Critical properties will be based off liquid phase if present, as we assume + # supercritical fluid is liquid-like. Otherwise, the vapor phase is used. + # Get first liquid phase we find + def _get_critical_ref_phase(self): + ref_phase = None + vap_phase = None + for p in self.params.phase_list: + if self.params.get_phase(p).is_liquid_phase(): + return p + elif self.params.get_phase(p).is_vapor_phase(): + vap_phase = p + + if ref_phase is None and vap_phase is not None: + # Use vapor phase for reference phase + return vap_phase + if ref_phase is None: + # If still no reference phase, raise an exception + raise PropertyPackageError( + "No liquid or vapor phase found to use as reference phase " + "for calculating critical properties." + ) + + def _critical_props(self): + ref_phase = self._get_critical_ref_phase() + + try: + base_units = self.params.get_metadata().default_units + + self.compress_fact_crit = Var( + doc="Critical compressibility factor of mixture", + units=pyunits.dimensionless, + ) + + self.dens_mol_crit = Var( + doc="Critical molar density of mixture", + units=base_units.DENS_MOL, + ) + + self.pressure_crit = Var( + doc="Critical pressure of mixture", + units=base_units.PRESSURE, + ) + + self.temperature_crit = Var( + doc="Critical temperature of mixture", + units=base_units.TEMPERATURE, + ) + + p_config = self.params.get_phase(ref_phase).config + p_config.equation_of_state.build_critical_properties(self, ref_phase) + + except AttributeError: + self.del_component(self.compress_fact_crit) + self.del_component(self.dens_mol_crit) + self.del_component(self.pressure_crit) + self.del_component(self.temperature_crit) + raise + # ------------------------------------------------------------------------- # Bubble and Dew Points diff --git a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py index 85b4950afb..0c383fb276 100644 --- a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py +++ b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py @@ -1461,3 +1461,150 @@ def test_visc_d_phase(self, frame): for p in frame.props[1].phase_list: assert value(frame.props[1].visc_d_phase[p]) == 4 + + +class TestCriticalProps: + @pytest.mark.unit + def test_get_critical_ref_phase_exception(self): + # No vapor or liquid phases + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": {}, + }, + phases={ + "p1": { + "equation_of_state": DummyEoS, + }, + "p2": { + "equation_of_state": DummyEoS, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + with pytest.raises( + AttributeError, + match="No liquid or vapor phase found to use as reference phase " + "for calculating critical properties.", + ): + m.props[1]._get_critical_ref_phase() + + @pytest.mark.unit + def test_get_critical_ref_phase_exception(self): + # No vapor or liquid phases + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": {}, + }, + phases={ + "p1": { + "equation_of_state": DummyEoS, + }, + "p2": { + "equation_of_state": DummyEoS, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + with pytest.raises( + PropertyPackageError, + match="No liquid or vapor phase found to use as reference phase " + "for calculating critical properties.", + ): + m.props[1]._get_critical_ref_phase() + + @pytest.mark.unit + def test_get_critical_ref_phase_VL(self): + # No vapor or liquid phases + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": {}, + }, + phases={ + "p1": { + "type": LiquidPhase, + "equation_of_state": DummyEoS, + }, + "p2": { + "type": VaporPhase, + "equation_of_state": DummyEoS, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + assert m.props[1]._get_critical_ref_phase() == "p1" + + @pytest.mark.unit + def test_get_critical_ref_phase_LL(self): + # No vapor or liquid phases + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": {}, + }, + phases={ + "p1": { + "type": LiquidPhase, + "equation_of_state": DummyEoS, + }, + "p2": { + "type": LiquidPhase, + "equation_of_state": DummyEoS, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + assert m.props[1]._get_critical_ref_phase() == "p1" + + @pytest.mark.unit + def test_get_critical_ref_phase_VX(self): + # No vapor or liquid phases + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": {}, + }, + phases={ + "p1": { + "type": VaporPhase, + "equation_of_state": DummyEoS, + }, + "p2": { + "equation_of_state": DummyEoS, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + assert m.props[1]._get_critical_ref_phase() == "p1" From 027ad8325c7be7bfa4fdf78ce4efe1690a778a99 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 13 Feb 2024 13:56:39 -0500 Subject: [PATCH 03/15] Hook for link ot EoS for critical properties --- idaes/core/base/property_meta.py | 4 +- idaes/core/base/property_set.py | 4 +- idaes/core/base/tests/test_property_meta.py | 2 +- .../base/generic_property.py | 2 +- .../base/tests/test_generic_property.py | 103 ++++++++++++++++++ .../modular_properties/eos/eos_base.py | 8 ++ 6 files changed, 117 insertions(+), 6 deletions(-) diff --git a/idaes/core/base/property_meta.py b/idaes/core/base/property_meta.py index 6c44031f66..e06b26eb7d 100644 --- a/idaes/core/base/property_meta.py +++ b/idaes/core/base/property_meta.py @@ -261,7 +261,7 @@ def VOLUME_MASS(self): return self._length**3 * self._mass**-1 @property - def VOLUME_MOL(self): + def VOLUME_MOLE(self): return self._length**3 * self._amount**-1 # Backward compatibility name @@ -271,7 +271,7 @@ def MOLAR_VOLUME(self): f"The unit name MOLAR_VOLUME is being deprecated in " "favor of VOLUME_MOL." ) deprecation_warning(msg=msg, logger=_log, version="2.3.0", remove_in="3.0.0") - return self.VOLUME_MOL + return self.VOLUME_MOLE # Flows @property diff --git a/idaes/core/base/property_set.py b/idaes/core/base/property_set.py index a98167b472..e9b8a41629 100644 --- a/idaes/core/base/property_set.py +++ b/idaes/core/base/property_set.py @@ -895,12 +895,12 @@ class StandardPropertySet(PropertySetBase): vol_mol = PropertyMetadata( name="vol_mol", doc="Molar Volume", - units="VOLUME_MOL", + units="VOLUME_MOLE", ) vol_mol_crit = PropertyMetadata( name="vol_mol", doc="Molar Volume at Critical Point", - units="VOLUME_MOL", + units="VOLUME_MOLE", ) # Log terms log_act = PropertyMetadata( diff --git a/idaes/core/base/tests/test_property_meta.py b/idaes/core/base/tests/test_property_meta.py index e67053b474..3ed90a2086 100644 --- a/idaes/core/base/tests/test_property_meta.py +++ b/idaes/core/base/tests/test_property_meta.py @@ -176,7 +176,7 @@ def test_luminous_intensity(unit_set): "area": units.m**2, "volume": units.m**3, "volume_mass": units.m**3 / units.kg, - "volume_mol": units.m**3 / units.mol, + "volume_mole": units.m**3 / units.mol, "flow_mass": units.kg * units.s**-1, "flow_mole": units.mol * units.s**-1, "flow_vol": units.m**3 * units.s**-1, diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index 6ea12736ad..442a3d1124 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -2961,7 +2961,7 @@ def _critical_props(self): self.dens_mol_crit = Var( doc="Critical molar density of mixture", - units=base_units.DENS_MOL, + units=base_units.DENSITY_MOLE, ) self.pressure_crit = Var( diff --git a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py index 0c383fb276..f64f355c99 100644 --- a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py +++ b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py @@ -1276,6 +1276,13 @@ def test_build_all(self, frame): frame.params.get_metadata().properties[p.name].method, ) continue + elif p.name.endswith(("_crit")): + # Critical properties will be tested elsewhere + assert hasattr( + frame.props[1], + frame.params.get_metadata().properties[p.name].method, + ) + continue elif p.name.endswith(("bubble", "bub", "dew")): # Bubble and dew properties require phase equilibria, which are # not tested here @@ -1608,3 +1615,99 @@ def test_get_critical_ref_phase_VX(self): m.props = m.params.build_state_block([1], defined_state=False) assert m.props[1]._get_critical_ref_phase() == "p1" + + @pytest.mark.unit + def test_critical_props_not_implemented(self): + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": {}, + }, + phases={ + "p1": { + "type": LiquidPhase, + "equation_of_state": DummyEoS, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + with pytest.raises( + NotImplementedError, + match="props\[1\] Equation of State module has not implemented a method for " + "build_critical_properties. Please contact the EoS developer or use a " + "different module.", + ): + m.props[1]._critical_props() + + @pytest.mark.unit + def test_critical_props_dummy_method(self): + class DummyEoS2(DummyEoS): + @staticmethod + def build_critical_properties(b, *args, **kwargs): + b._dummy_crit_executed = True + + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": {}, + }, + phases={ + "p1": { + "type": LiquidPhase, + "equation_of_state": DummyEoS2, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + m.props[1]._critical_props() + + assert m.props[1]._dummy_crit_executed + + assert isinstance(m.props[1].compress_fact_crit, Var) + assert isinstance(m.props[1].dens_mol_crit, Var) + assert isinstance(m.props[1].pressure_crit, Var) + assert isinstance(m.props[1].temperature_crit, Var) + + @pytest.mark.unit + def test_critical_props_attribute_errord(self): + class DummyEoS2(DummyEoS): + @staticmethod + def build_critical_properties(b, *args, **kwargs): + raise AttributeError() + + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": {}, + }, + phases={ + "p1": { + "type": LiquidPhase, + "equation_of_state": DummyEoS2, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + try: + m.props[1]._critical_props() + except AttributeError: + assert not hasattr(m.props[1], "compress_fact_crit") + assert not hasattr(m.props[1], "dens_mol_crit") + assert not hasattr(m.props[1], "pressure_crit") + assert not hasattr(m.props[1], "temperature_crit") diff --git a/idaes/models/properties/modular_properties/eos/eos_base.py b/idaes/models/properties/modular_properties/eos/eos_base.py index 7cf4cd8b04..3689b601b0 100644 --- a/idaes/models/properties/modular_properties/eos/eos_base.py +++ b/idaes/models/properties/modular_properties/eos/eos_base.py @@ -50,6 +50,14 @@ def calculate_scaling_factors(b, pobj): def build_parameters(b): raise NotImplementedError(_msg(b, "build_parameters")) + @staticmethod + def build_critical_properties(b, p): + raise NotImplementedError(_msg(b, "build_critical_properties")) + + @staticmethod + def initialize_critical_properties(b, p): + raise NotImplementedError(_msg(b, "initialize_critical_properties")) + @staticmethod def get_vol_mol_pure(b, phase, comp, temperature): try: From 4af1dda3ecba5865c6dba1a3cc6ced9c0ce9c564 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 13 Feb 2024 16:03:40 -0500 Subject: [PATCH 04/15] Adding critical properties to ideal EoS --- .../base/generic_property.py | 2 +- .../modular_properties/eos/eos_base.py | 6 +- .../modular_properties/eos/ideal.py | 38 +++++++++- .../eos/tests/test_ideal.py | 73 ++++++++++++++++++- 4 files changed, 111 insertions(+), 8 deletions(-) diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index 442a3d1124..5eadfcd6d7 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -2975,7 +2975,7 @@ def _critical_props(self): ) p_config = self.params.get_phase(ref_phase).config - p_config.equation_of_state.build_critical_properties(self, ref_phase) + p_config.equation_of_state.build_critical_properties(self) except AttributeError: self.del_component(self.compress_fact_crit) diff --git a/idaes/models/properties/modular_properties/eos/eos_base.py b/idaes/models/properties/modular_properties/eos/eos_base.py index 3689b601b0..a8f0d3f53d 100644 --- a/idaes/models/properties/modular_properties/eos/eos_base.py +++ b/idaes/models/properties/modular_properties/eos/eos_base.py @@ -51,13 +51,9 @@ def build_parameters(b): raise NotImplementedError(_msg(b, "build_parameters")) @staticmethod - def build_critical_properties(b, p): + def build_critical_properties(b): raise NotImplementedError(_msg(b, "build_critical_properties")) - @staticmethod - def initialize_critical_properties(b, p): - raise NotImplementedError(_msg(b, "initialize_critical_properties")) - @staticmethod def get_vol_mol_pure(b, phase, comp, temperature): try: diff --git a/idaes/models/properties/modular_properties/eos/ideal.py b/idaes/models/properties/modular_properties/eos/ideal.py index 8986e7286a..c791263b04 100644 --- a/idaes/models/properties/modular_properties/eos/ideal.py +++ b/idaes/models/properties/modular_properties/eos/ideal.py @@ -21,7 +21,7 @@ # TODO: Look into protected access issues # pylint: disable=protected-access -from pyomo.environ import Expression, log +from pyomo.environ import Constraint, Expression, log, units from idaes.core import Apparent from idaes.core.util.exceptions import ConfigurationError, PropertyNotSupportedError @@ -34,6 +34,7 @@ log_henry_pressure, ) from .eos_base import EoSBase +from idaes.core.util.constants import Constants as CONST # TODO: Add support for ideal solids @@ -48,6 +49,41 @@ def common(b, pobj): # No common components required for ideal property calculations pass + @staticmethod + def build_critical_properties(b): + base_units = b.params.get_metadata().default_units + + # Use Kay's method for critical pressure and temperature + def p_crit_rule(b): + return b.pressure_crit == sum( + b.mole_frac_comp[j] * b.params.get_component(j).pressure_crit + for j in b.component_list + ) + + b.pressure_crit_constraint = Constraint(rule=p_crit_rule) + + def t_crit_rule(b): + return b.temperature_crit == sum( + b.mole_frac_comp[j] * b.params.get_component(j).temperature_crit + for j in b.component_list + ) + + b.temperature_crit_constraint = Constraint(rule=t_crit_rule) + + # Fix critical compressibility factor to 1 - user can override if necessary + b.compress_fact_crit.fix(1) + + def rho_crit_rule(b): + return b.pressure_crit == units.convert( + b.compress_fact_crit + * CONST.gas_constant + * b.temperature_crit + * b.dens_mol_crit, + to_units=base_units.PRESSURE, + ) + + b.dens_mol_crit_constraint = Constraint(rule=rho_crit_rule) + @staticmethod def calculate_scaling_factors(b, pobj): pass diff --git a/idaes/models/properties/modular_properties/eos/tests/test_ideal.py b/idaes/models/properties/modular_properties/eos/tests/test_ideal.py index 61dd2082fa..c9a272b4f3 100644 --- a/idaes/models/properties/modular_properties/eos/tests/test_ideal.py +++ b/idaes/models/properties/modular_properties/eos/tests/test_ideal.py @@ -18,7 +18,7 @@ import pytest from sys import modules -from pyomo.environ import ConcreteModel, Var, units as pyunits, value +from pyomo.environ import ConcreteModel, Constraint, Var, units as pyunits, value from pyomo.util.check_units import assert_units_equivalent from idaes.core import ( @@ -686,3 +686,74 @@ def test_vol_mol_phase(): + 1 / 42.0 * m.props[1].mole_frac_phase_comp["Liq", "b"] + 42.0 * m.props[1].mole_frac_phase_comp["Liq", "c"] ) + + +@pytest.mark.unit +def test_critical_props(): + m = ConcreteModel() + + # Dummy params block + m.params = DummyParameterBlock( + components={ + "a": {"parameter_data": {"pressure_crit": 1e5, "temperature_crit": 100}}, + "b": {"parameter_data": {"pressure_crit": 2e5, "temperature_crit": 200}}, + "c": {"parameter_data": {"pressure_crit": 3e5, "temperature_crit": 300}}, + }, + phases={ + "Vap": {"type": VaporPhase, "equation_of_state": Ideal}, + "Liq": {"type": LiquidPhase, "equation_of_state": Ideal}, + }, + base_units={ + "time": pyunits.s, + "length": pyunits.m, + "mass": pyunits.kg, + "amount": pyunits.mol, + "temperature": pyunits.K, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + ) + + m.props = m.params.state_block_class([1], defined_state=False, parameters=m.params) + + # Add common variables + m.props[1].mole_frac_comp = Var(m.params.component_list, initialize=0.5) + + m.props[1]._critical_props() + + assert isinstance(m.props[1].compress_fact_crit, Var) + assert isinstance(m.props[1].dens_mol_crit, Var) + assert isinstance(m.props[1].pressure_crit, Var) + assert isinstance(m.props[1].temperature_crit, Var) + + assert isinstance(m.props[1].pressure_crit_constraint, Constraint) + assert str(m.props[1].pressure_crit_constraint.expr) == str( + m.props[1].pressure_crit + == m.props[1].mole_frac_comp["a"] * m.params.a.pressure_crit + + m.props[1].mole_frac_comp["b"] * m.params.b.pressure_crit + + m.props[1].mole_frac_comp["c"] * m.params.c.pressure_crit + ) + + assert isinstance(m.props[1].temperature_crit_constraint, Constraint) + assert str(m.props[1].temperature_crit_constraint.expr) == str( + m.props[1].temperature_crit + == m.props[1].mole_frac_comp["a"] * m.params.a.temperature_crit + + m.props[1].mole_frac_comp["b"] * m.params.b.temperature_crit + + m.props[1].mole_frac_comp["c"] * m.params.c.temperature_crit + ) + + assert m.props[1].compress_fact_crit.fixed + assert value(m.props[1].compress_fact_crit) == 1 + + assert isinstance(m.props[1].dens_mol_crit_constraint, Constraint) + assert str(m.props[1].dens_mol_crit_constraint.expr) == str( + m.props[1].pressure_crit + == pyunits.convert( + m.props[1].compress_fact_crit + * const.gas_constant + * m.props[1].temperature_crit + * m.props[1].dens_mol_crit, + to_units=pyunits.kg / pyunits.m / pyunits.s**2, + ) + ) From d795063627ff4443a8654127ea8ba4c57d471691 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Thu, 15 Feb 2024 13:06:44 -0500 Subject: [PATCH 05/15] Initial working cubic critical properties --- .../base/generic_property.py | 2 +- .../properties/modular_properties/eos/ceos.py | 131 +++++++- .../modular_properties/eos/eos_base.py | 2 +- .../modular_properties/eos/ideal.py | 2 +- .../eos/tests/test_ceos_PR.py | 282 ++++++++++++++++++ .../eos/tests/test_ideal.py | 1 - 6 files changed, 415 insertions(+), 5 deletions(-) diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index 5eadfcd6d7..442a3d1124 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -2975,7 +2975,7 @@ def _critical_props(self): ) p_config = self.params.get_phase(ref_phase).config - p_config.equation_of_state.build_critical_properties(self) + p_config.equation_of_state.build_critical_properties(self, ref_phase) except AttributeError: self.del_component(self.compress_fact_crit) diff --git a/idaes/models/properties/modular_properties/eos/ceos.py b/idaes/models/properties/modular_properties/eos/ceos.py index 2157d78db5..002ff4f729 100644 --- a/idaes/models/properties/modular_properties/eos/ceos.py +++ b/idaes/models/properties/modular_properties/eos/ceos.py @@ -25,6 +25,7 @@ from copy import deepcopy from pyomo.environ import ( + Constraint, exp, Expression, log, @@ -350,7 +351,7 @@ def rule_bm(m, p): else: raise ConfigurationError( "{} Unrecognized option for Equation of State " - "mixing_rule_a: {}. Must be an instance of MixingRuleB " + "mixing_rule_b: {}. Must be an instance of MixingRuleB " "Enum.".format(b.name, mixing_rule_b) ) @@ -901,6 +902,116 @@ def vol_mol_phase_comp(b, p, j): * (b.compress_fact_phase[p] + _N_dZ_dNj(b, p, j)) ) + @staticmethod + def build_critical_properties(m, ref_phase): + pobj = m.params.get_phase(ref_phase) + cname = pobj.config.equation_of_state_options["type"].name + ctype = pobj._cubic_type + + # Add components for calculation of mixture critical properties + def func_a_crit(m, j): + cobj = m.params.get_component(j) + fw = getattr(m, cname + "_fw") + return ( + EoS_param[ctype]["omegaA"] + * ( + (Cubic.gas_constant(m) * cobj.temperature_crit) ** 2 + / cobj.pressure_crit + ) + * ( + (1 + fw[j] * (1 - sqrt(m.temperature_crit / cobj.temperature_crit))) + ** 2 + ) + ) + + m.add_component( + "a_crit", + Expression( + m.component_list, rule=func_a_crit, doc="Component a coefficient" + ), + ) + + def rule_am_crit(m): + try: + rule = m.params.get_phase(ref_phase).config.equation_of_state_options[ + "mixing_rule_a" + ] + except (KeyError, TypeError): + rule = MixingRuleA.default + + a_crit = getattr(m, "a_crit") + if rule == MixingRuleA.default: + return rule_am_crit_default(m, cname, a_crit) + else: + raise ConfigurationError( + "{} Unrecognized option for Equation of State " + "mixing_rule_a: {}. Must be an instance of MixingRuleA " + "Enum.".format(m.name, rule) + ) + + m.add_component("am_crit", Expression(rule=rule_am_crit)) + + def rule_bm_crit(m): + try: + rule = m.params.get_phase(ref_phase).config.equation_of_state_options[ + "mixing_rule_b" + ] + except KeyError: + rule = MixingRuleB.default + + b = getattr(m, cname + "_b") + if rule == MixingRuleB.default: + return rule_bm_crit_default(m, b) + else: + raise ConfigurationError( + "{} Unrecognized option for Equation of State " + "mixing_rule_b: {}. Must be an instance of MixingRuleB " + "Enum.".format(m.name, rule) + ) + + m.add_component("bm_crit", Expression(rule=rule_bm_crit)) + + def rule_A_crit(m): + am_crit = getattr(m, "am_crit") + return EoS_param[ctype]["omegaA"] * ( + Cubic.gas_constant(m) * m.temperature_crit + ) ** 2 == (am_crit * m.pressure_crit) + + m.add_component("A_crit", Constraint(rule=rule_A_crit)) + + def rule_B_crit(m): + bm_crit = getattr(m, "bm_crit") + return EoS_param[ctype]["coeff_b"] * Cubic.gas_constant( + m + ) * m.temperature_crit == (bm_crit * m.pressure_crit) + + m.add_component("B_crit", Constraint(rule=rule_B_crit)) + + Acrit = ( + m.am_crit + * m.pressure_crit + / (Cubic.gas_constant(m) * m.temperature_crit) ** 2 + ) + Bcrit = ( + m.bm_crit * m.pressure_crit / (Cubic.gas_constant(m) * m.temperature_crit) + ) + + expr_write = CubicThermoExpressions(m) + if pobj.is_vapor_phase(): + Z_crit = expr_write.z_vap(eos=pobj._cubic_type, A=Acrit, B=Bcrit) + elif pobj.is_liquid_phase(): + Z_crit = expr_write.z_liq(eos=pobj._cubic_type, A=Acrit, B=Bcrit) + + m.compress_fact_crit_eq = Constraint(rule=m.compress_fact_crit == Z_crit) + + m.dens_mol_crit_eq = Constraint( + rule=m.pressure_crit + == m.compress_fact_crit + * Cubic.gas_constant(m) + * m.temperature_crit + * m.dens_mol_crit + ) + def _dZ_dT(blk, p): pobj = blk.params.get_phase(p) @@ -1277,5 +1388,23 @@ def rule_am_default(m, cname, a, p, pp=()): ) +def rule_am_crit_default(m, cname, a_crit): + k = getattr(m.params, cname + "_kappa") + return sum( + sum( + m.mole_frac_comp[i] + * m.mole_frac_comp[j] + * sqrt(a_crit[i] * a_crit[j]) + * (1 - k[i, j]) + for j in m.component_list + ) + for i in m.component_list + ) + + def rule_bm_default(m, b, p): return sum(m.mole_frac_phase_comp[p, i] * b[i] for i in m.components_in_phase(p)) + + +def rule_bm_crit_default(m, b): + return sum(m.mole_frac_comp[i] * b[i] for i in m.component_list) diff --git a/idaes/models/properties/modular_properties/eos/eos_base.py b/idaes/models/properties/modular_properties/eos/eos_base.py index a8f0d3f53d..054a7ad2d1 100644 --- a/idaes/models/properties/modular_properties/eos/eos_base.py +++ b/idaes/models/properties/modular_properties/eos/eos_base.py @@ -51,7 +51,7 @@ def build_parameters(b): raise NotImplementedError(_msg(b, "build_parameters")) @staticmethod - def build_critical_properties(b): + def build_critical_properties(b, ref_phase): raise NotImplementedError(_msg(b, "build_critical_properties")) @staticmethod diff --git a/idaes/models/properties/modular_properties/eos/ideal.py b/idaes/models/properties/modular_properties/eos/ideal.py index c791263b04..75ae312641 100644 --- a/idaes/models/properties/modular_properties/eos/ideal.py +++ b/idaes/models/properties/modular_properties/eos/ideal.py @@ -50,7 +50,7 @@ def common(b, pobj): pass @staticmethod - def build_critical_properties(b): + def build_critical_properties(b, ref_phase): base_units = b.params.get_metadata().default_units # Use Kay's method for critical pressure and temperature diff --git a/idaes/models/properties/modular_properties/eos/tests/test_ceos_PR.py b/idaes/models/properties/modular_properties/eos/tests/test_ceos_PR.py index 7fbaffa440..6b734b0f46 100644 --- a/idaes/models/properties/modular_properties/eos/tests/test_ceos_PR.py +++ b/idaes/models/properties/modular_properties/eos/tests/test_ceos_PR.py @@ -19,7 +19,9 @@ from sys import modules from pyomo.environ import ( + assert_optimal_termination, ConcreteModel, + Constraint, Expression, log, sqrt, @@ -28,6 +30,7 @@ units as pyunits, ) from pyomo.core.expr.numeric_expr import ExternalFunctionExpression +from pyomo.util.check_units import assert_units_consistent from idaes.core import declare_process_block_class, LiquidPhase, VaporPhase, SolidPhase from idaes.models.properties.modular_properties.eos.ceos import Cubic, CubicType @@ -37,6 +40,8 @@ from idaes.core.util.exceptions import PropertyNotSupportedError, ConfigurationError from idaes.core.util.constants import Constants as const from idaes.models.properties.modular_properties.eos.ceos import cubic_roots_available +from idaes.core.solvers import get_solver + # Dummy method for property method calls def dummy_call(b, j, T): @@ -1065,3 +1070,280 @@ def test_vol_mol_phase_comp(m): for j in m.params.component_list ) ) + + +class TestCEOSCriticalProps: + @pytest.fixture(scope="class") + def model(self): + m = ConcreteModel() + + # Dummy params block + m.params = DummyParameterBlock( + components={ + "a": { + "parameter_data": { + "omega": 0.1, + "pressure_crit": 1e5, + "temperature_crit": 100, + } + }, + "b": { + "parameter_data": { + "omega": 0.2, + "pressure_crit": 2e5, + "temperature_crit": 200, + } + }, + "c": { + "parameter_data": { + "omega": 0.3, + "pressure_crit": 3e5, + "temperature_crit": 300, + } + }, + }, + phases={ + "Vap": { + "type": LiquidPhase, + "equation_of_state": Cubic, + "equation_of_state_options": {"type": CubicType.PR}, + }, + "Liq": { + "type": LiquidPhase, + "equation_of_state": Cubic, + "equation_of_state_options": {"type": CubicType.PR}, + }, + }, + base_units={ + "time": pyunits.s, + "length": pyunits.m, + "mass": pyunits.kg, + "amount": pyunits.mol, + "temperature": pyunits.K, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + parameter_data={ + "PR_kappa": { + ("a", "a"): 0.0, + ("a", "b"): 0.0, + ("a", "c"): 0.0, + ("b", "a"): 0.0, + ("b", "b"): 0.0, + ("b", "c"): 0.0, + ("c", "a"): 0.0, + ("c", "b"): 0.0, + ("c", "c"): 0.0, + }, + }, + ) + + m.props = m.params.state_block_class( + [1], defined_state=False, parameters=m.params + ) + + # Add common variables + m.props[1].mole_frac_comp = Var( + m.params.component_list, initialize=0.5, bounds=(1e-12, 1) + ) + m.props[1].mole_frac_comp["a"].fix(0.4) + m.props[1].mole_frac_comp["b"].fix(0.6) + m.props[1].mole_frac_comp["c"].fix(1e-8) + + # Build critical properties + m.props[1]._critical_props() + + return m + + @pytest.mark.unit + def test_build_critical_props(self, model): + assert isinstance(model.props[1].a_crit, Expression) + assert isinstance(model.props[1].am_crit, Expression) + assert isinstance(model.props[1].bm_crit, Expression) + + base_units = model.params.get_metadata().default_units + + assert str(model.props[1].a_crit["a"].expr) == str( + 0.45724 + * ( + ( + pyunits.convert( + const.gas_constant, to_units=base_units.GAS_CONSTANT + ) + * model.params.a.temperature_crit + ) + ** 2 + / model.params.a.pressure_crit + ) + * ( + 1 + + model.props[1].PR_fw["a"] + * ( + 1 + - sqrt( + model.props[1].temperature_crit + / model.params.a.temperature_crit + ) + ) + ) + ** 2 + ) + assert str(model.props[1].a_crit["b"].expr) == str( + 0.45724 + * ( + ( + pyunits.convert( + const.gas_constant, to_units=base_units.GAS_CONSTANT + ) + * model.params.b.temperature_crit + ) + ** 2 + / model.params.b.pressure_crit + ) + * ( + 1 + + model.props[1].PR_fw["b"] + * ( + 1 + - sqrt( + model.props[1].temperature_crit + / model.params.b.temperature_crit + ) + ) + ) + ** 2 + ) + assert str(model.props[1].a_crit["c"].expr) == str( + 0.45724 + * ( + ( + pyunits.convert( + const.gas_constant, to_units=base_units.GAS_CONSTANT + ) + * model.params.c.temperature_crit + ) + ** 2 + / model.params.c.pressure_crit + ) + * ( + 1 + + model.props[1].PR_fw["c"] + * ( + 1 + - sqrt( + model.props[1].temperature_crit + / model.params.c.temperature_crit + ) + ) + ) + ** 2 + ) + + assert str(model.props[1].am_crit.expr) == str( + ( + model.props[1].mole_frac_comp["a"] + * model.props[1].mole_frac_comp["a"] + * sqrt(model.props[1].a_crit["a"] * model.props[1].a_crit["a"]) + * (1 - model.params.PR_kappa["a", "a"]) + ) + + ( + model.props[1].mole_frac_comp["a"] + * model.props[1].mole_frac_comp["b"] + * sqrt(model.props[1].a_crit["a"] * model.props[1].a_crit["b"]) + * (1 - model.params.PR_kappa["a", "b"]) + ) + + ( + model.props[1].mole_frac_comp["a"] + * model.props[1].mole_frac_comp["c"] + * sqrt(model.props[1].a_crit["a"] * model.props[1].a_crit["c"]) + * (1 - model.params.PR_kappa["a", "c"]) + ) + + ( + model.props[1].mole_frac_comp["b"] + * model.props[1].mole_frac_comp["a"] + * sqrt(model.props[1].a_crit["b"] * model.props[1].a_crit["a"]) + * (1 - model.params.PR_kappa["b", "a"]) + ) + + ( + model.props[1].mole_frac_comp["b"] + * model.props[1].mole_frac_comp["b"] + * sqrt(model.props[1].a_crit["b"] * model.props[1].a_crit["b"]) + * (1 - model.params.PR_kappa["b", "b"]) + ) + + ( + model.props[1].mole_frac_comp["b"] + * model.props[1].mole_frac_comp["c"] + * sqrt(model.props[1].a_crit["b"] * model.props[1].a_crit["c"]) + * (1 - model.params.PR_kappa["b", "c"]) + ) + + ( + model.props[1].mole_frac_comp["c"] + * model.props[1].mole_frac_comp["a"] + * sqrt(model.props[1].a_crit["c"] * model.props[1].a_crit["a"]) + * (1 - model.params.PR_kappa["c", "a"]) + ) + + ( + model.props[1].mole_frac_comp["c"] + * model.props[1].mole_frac_comp["b"] + * sqrt(model.props[1].a_crit["c"] * model.props[1].a_crit["b"]) + * (1 - model.params.PR_kappa["c", "b"]) + ) + + ( + model.props[1].mole_frac_comp["c"] + * model.props[1].mole_frac_comp["c"] + * sqrt(model.props[1].a_crit["c"] * model.props[1].a_crit["c"]) + * (1 - model.params.PR_kappa["c", "c"]) + ) + ) + + assert str(model.props[1].bm_crit.expr) == str( + model.props[1].mole_frac_comp["a"] * model.props[1].PR_b["a"] + + model.props[1].mole_frac_comp["b"] * model.props[1].PR_b["b"] + + model.props[1].mole_frac_comp["c"] * model.props[1].PR_b["c"] + ) + + assert isinstance(model.props[1].A_crit, Constraint) + assert isinstance(model.props[1].B_crit, Constraint) + + assert str(model.props[1].A_crit.expr) == str( + 0.45724 + * ( + pyunits.convert(const.gas_constant, to_units=base_units.GAS_CONSTANT) + * model.props[1].temperature_crit + ) + ** 2 + == model.props[1].am_crit * model.props[1].pressure_crit + ) + + assert str(model.props[1].B_crit.expr) == str( + 0.07780 + * pyunits.convert(const.gas_constant, to_units=base_units.GAS_CONSTANT) + * model.props[1].temperature_crit + == model.props[1].bm_crit * model.props[1].pressure_crit + ) + + assert_units_consistent(model) + + @pytest.mark.component + @pytest.mark.solver + def test_build_critical_props(self, model): + # Initialize model + model.props[1].compress_fact_crit.set_value(1) + model.props[1].pressure_crit.set_value(1.6e5) + model.props[1].temperature_crit.set_value(160) + + # Solve model + solver = get_solver() + res = solver.solve(model, tee=True) + assert_optimal_termination(res) + + # Confirm results + assert value(model.props[1].compress_fact_crit) == pytest.approx( + 0.321379, rel=1e-5 + ) + assert value(model.props[1].pressure_crit) == pytest.approx(1.58138e5, rel=1e-5) + assert value(model.props[1].temperature_crit) == pytest.approx( + 158.138, rel=1e-5 + ) + assert value(model.props[1].dens_mol_crit) == pytest.approx(374.238, rel=1e-5) diff --git a/idaes/models/properties/modular_properties/eos/tests/test_ideal.py b/idaes/models/properties/modular_properties/eos/tests/test_ideal.py index c9a272b4f3..76eb169765 100644 --- a/idaes/models/properties/modular_properties/eos/tests/test_ideal.py +++ b/idaes/models/properties/modular_properties/eos/tests/test_ideal.py @@ -217,7 +217,6 @@ def test_cv_mol_phase_comp(m): @pytest.mark.unit def test_heat_capacity_ratio_phase(m): - m.props[1].cp_mol_phase = Var(m.params.phase_list) m.props[1].cv_mol_phase = Var(m.params.phase_list) From 59525d6b1ab862344b67f20887d0245ba4f747cb Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Fri, 16 Feb 2024 12:31:25 -0500 Subject: [PATCH 06/15] Adding basic initialization --- .../base/generic_property.py | 36 +++++++ .../base/tests/test_generic_property.py | 100 +++++++++++++++++- 2 files changed, 135 insertions(+), 1 deletion(-) diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index 442a3d1124..fe9157fe6d 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -1295,6 +1295,13 @@ def initialization_routine( ): k.inherent_equilibrium_constraint.deactivate() + # --------------------------------------------------------------------- + # If present, initialize critical properties + for k in model.values(): + # We only need to look for one critical property as it is all or nothing + if hasattr(k, "pressure_crit"): + self.initialize_critical_props(k) + # --------------------------------------------------------------------- # If present, initialize bubble and dew point calculations for k in model.values(): @@ -1561,6 +1568,35 @@ def initialization_routine( return None + def initialize_critical_props(self, state_data): + params = state_data.params + # Use mole weighted sum of component critical properties + state_data.compress_fact_crit.set_value( + sum( + state_data.mole_frac_comp[j] + * params.get_component(j).compress_fact_crit + for j in state_data.component_list + ) + ) + state_data.dens_mol_crit.set_value( + sum( + state_data.mole_frac_comp[j] * params.get_component(j).dens_mol_crit + for j in state_data.component_list + ) + ) + state_data.pressure_crit.set_value( + sum( + state_data.mole_frac_comp[j] * params.get_component(j).pressure_crit + for j in state_data.component_list + ) + ) + state_data.temperature_crit.set_value( + sum( + state_data.mole_frac_comp[j] * params.get_component(j).temperature_crit + for j in state_data.component_list + ) + ) + class _GenericStateBlock(StateBlock): """ diff --git a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py index f64f355c99..24975dc849 100644 --- a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py +++ b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py @@ -25,6 +25,7 @@ from idaes.models.properties.modular_properties.base.generic_property import ( GenericParameterData, GenericStateBlock, + ModularPropertiesInitializer, ) from idaes.models.properties.modular_properties.base.tests.dummy_eos import DummyEoS @@ -42,6 +43,7 @@ from idaes.models.properties.modular_properties.phase_equil.henry import HenryType from idaes.core.base.property_meta import UnitSet from idaes.core.initialization import BlockTriangularizationInitializer +from idaes.models.properties.modular_properties.eos.ceos import Cubic, CubicType import idaes.logger as idaeslog @@ -1680,7 +1682,7 @@ def build_critical_properties(b, *args, **kwargs): assert isinstance(m.props[1].temperature_crit, Var) @pytest.mark.unit - def test_critical_props_attribute_errord(self): + def test_critical_props_attribute_error(self): class DummyEoS2(DummyEoS): @staticmethod def build_critical_properties(b, *args, **kwargs): @@ -1711,3 +1713,99 @@ def build_critical_properties(b, *args, **kwargs): assert not hasattr(m.props[1], "dens_mol_crit") assert not hasattr(m.props[1], "pressure_crit") assert not hasattr(m.props[1], "temperature_crit") + + @pytest.mark.unit + def test_initialize_critical_props(self): + m = ConcreteModel() + + class DummyEoS2(DummyEoS): + @staticmethod + def build_critical_properties(b, *args, **kwargs): + b._dummy_crit_executed = True + + # Dummy params block + m.params = DummyParameterBlock( + components={ + "a": { + "parameter_data": { + "compress_fact_crit": 0.1, + "dens_mol_crit": 1, + "pressure_crit": 1e5, + "temperature_crit": 100, + } + }, + "b": { + "parameter_data": { + "compress_fact_crit": 0.2, + "dens_mol_crit": 2, + "pressure_crit": 2e5, + "temperature_crit": 200, + } + }, + "c": { + "parameter_data": { + "compress_fact_crit": 0.3, + "dens_mol_crit": 3, + "pressure_crit": 3e5, + "temperature_crit": 300, + } + }, + }, + phases={ + "Vap": { + "type": LiquidPhase, + "equation_of_state": DummyEoS2, + }, + "Liq": { + "type": LiquidPhase, + "equation_of_state": DummyEoS2, + }, + }, + base_units={ + "time": pyunits.s, + "length": pyunits.m, + "mass": pyunits.kg, + "amount": pyunits.mol, + "temperature": pyunits.K, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + parameter_data={ + "PR_kappa": { + ("a", "a"): 0.0, + ("a", "b"): 0.0, + ("a", "c"): 0.0, + ("b", "a"): 0.0, + ("b", "b"): 0.0, + ("b", "c"): 0.0, + ("c", "a"): 0.0, + ("c", "b"): 0.0, + ("c", "c"): 0.0, + }, + }, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + # Add common variables + m.props[1].mole_frac_comp = Var( + m.params.component_list, initialize=0.5, bounds=(1e-12, 1) + ) + m.props[1].mole_frac_comp["a"].fix(0.4) + m.props[1].mole_frac_comp["b"].fix(0.6) + m.props[1].mole_frac_comp["c"].fix(1e-8) + + # Build critical props + m.props[1]._critical_props() + + # Initialize critical props + initializer = ModularPropertiesInitializer() + initializer.initialize_critical_props(m.props[1]) + + assert ( + value(m.props[1].compress_fact_crit) == 0.4 * 0.1 + 0.6 * 0.2 + 1e-8 * 0.3 + ) + assert value(m.props[1].dens_mol_crit) == 0.4 * 1 + 0.6 * 2 + 1e-8 * 3 + assert value(m.props[1].pressure_crit) == 0.4 * 1e5 + 0.6 * 2e5 + 1e-8 * 3e5 + assert value(m.props[1].temperature_crit) == 0.4 * 1e2 + 0.6 * 2e2 + 1e-8 * 3e2 From 9d06c2efd9e17888d0880b02da66cc805101911d Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Fri, 16 Feb 2024 16:22:43 -0500 Subject: [PATCH 07/15] Refining crit prop initialization --- .../base/generic_property.py | 87 ++++++++++++------- .../base/tests/test_generic_property.py | 51 ++++++++++- 2 files changed, 101 insertions(+), 37 deletions(-) diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index fe9157fe6d..79bfebbf61 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -1299,8 +1299,9 @@ def initialization_routine( # If present, initialize critical properties for k in model.values(): # We only need to look for one critical property as it is all or nothing - if hasattr(k, "pressure_crit"): - self.initialize_critical_props(k) + with k.lock_attribute_creation_context(): + if hasattr(k, "pressure_crit"): + _initialize_critical_props(k) # --------------------------------------------------------------------- # If present, initialize bubble and dew point calculations @@ -1568,35 +1569,6 @@ def initialization_routine( return None - def initialize_critical_props(self, state_data): - params = state_data.params - # Use mole weighted sum of component critical properties - state_data.compress_fact_crit.set_value( - sum( - state_data.mole_frac_comp[j] - * params.get_component(j).compress_fact_crit - for j in state_data.component_list - ) - ) - state_data.dens_mol_crit.set_value( - sum( - state_data.mole_frac_comp[j] * params.get_component(j).dens_mol_crit - for j in state_data.component_list - ) - ) - state_data.pressure_crit.set_value( - sum( - state_data.mole_frac_comp[j] * params.get_component(j).pressure_crit - for j in state_data.component_list - ) - ) - state_data.temperature_crit.set_value( - sum( - state_data.mole_frac_comp[j] * params.get_component(j).temperature_crit - for j in state_data.component_list - ) - ) - class _GenericStateBlock(StateBlock): """ @@ -1768,6 +1740,14 @@ def initialize( # Create solver opt = get_solver(solver, optarg) + # --------------------------------------------------------------------- + # If present, initialize critical properties + for k in blk.values(): + # We only need to look for one critical property as it is all or nothing + with k.lock_attribute_creation_context(): + if hasattr(k, "pressure_crit"): + _initialize_critical_props(k) + # --------------------------------------------------------------------- # If present, initialize bubble and dew point calculations for k in blk.values(): @@ -2963,11 +2943,23 @@ def get_mole_frac(self, phase=None): # ------------------------------------------------------------------------- # Mixture critical point # Critical properties will be based off liquid phase if present, as we assume - # supercritical fluid is liquid-like. Otherwise, the vapor phase is used. - # Get first liquid phase we find + # supercritical fluid is liquid-like. def _get_critical_ref_phase(self): ref_phase = None vap_phase = None + # First, look for VLE pair and use liquid phase if present + if self.params.config.phases_in_equilibrium is not None: + for pp in self.params.config.phases_in_equilibrium: + p1 = self.params.get_phase(pp[0]) + p2 = self.params.get_phase(pp[1]) + + if p1.is_liquid_phase() and p2.is_vapor_phase(): + return pp[0] + elif p2.is_liquid_phase() and p1.is_vapor_phase(): + return pp[1] + + # Next, iterate all phases and return either the first liquid phase + # Also collect vapor phase for final fall back for p in self.params.phase_list: if self.params.get_phase(p).is_liquid_phase(): return p @@ -5145,3 +5137,32 @@ def rule_log_mole_frac(b, p1, p2, j): eqn = getattr(b, "log_mole_frac_" + abbrv + "_eqn") b.del_component(eqn) raise + + +def _initialize_critical_props(state_data): + params = state_data.params + # Use mole weighted sum of component critical properties + state_data.compress_fact_crit.set_value( + sum( + state_data.mole_frac_comp[j] * params.get_component(j).compress_fact_crit + for j in state_data.component_list + ) + ) + state_data.dens_mol_crit.set_value( + sum( + state_data.mole_frac_comp[j] * params.get_component(j).dens_mol_crit + for j in state_data.component_list + ) + ) + state_data.pressure_crit.set_value( + sum( + state_data.mole_frac_comp[j] * params.get_component(j).pressure_crit + for j in state_data.component_list + ) + ) + state_data.temperature_crit.set_value( + sum( + state_data.mole_frac_comp[j] * params.get_component(j).temperature_crit + for j in state_data.component_list + ) + ) diff --git a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py index 24975dc849..cc7fbef1ea 100644 --- a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py +++ b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py @@ -25,7 +25,7 @@ from idaes.models.properties.modular_properties.base.generic_property import ( GenericParameterData, GenericStateBlock, - ModularPropertiesInitializer, + _initialize_critical_props, ) from idaes.models.properties.modular_properties.base.tests.dummy_eos import DummyEoS @@ -43,7 +43,6 @@ from idaes.models.properties.modular_properties.phase_equil.henry import HenryType from idaes.core.base.property_meta import UnitSet from idaes.core.initialization import BlockTriangularizationInitializer -from idaes.models.properties.modular_properties.eos.ceos import Cubic, CubicType import idaes.logger as idaeslog @@ -1127,6 +1126,7 @@ def test_default_scaling_convert(self): # Dummy methods for testing build calls to sub-modules def define_state(b): b.state_defined = True + b.temperature = Var(initialize=100) def get_material_flow_basis(self, *args, **kwargs): @@ -1591,6 +1591,50 @@ def test_get_critical_ref_phase_LL(self): assert m.props[1]._get_critical_ref_phase() == "p1" + @pytest.mark.unit + def test_get_critical_ref_phase_VLL(self): + class DummyPE: + def phase_equil(self, *args, **kwargs): + pass + + @staticmethod + def return_expression(b, *args, **kwargs): + return b.temperature == 42 + + # No vapor or liquid phases + m = ConcreteModel() + m.params = DummyParameterBlock( + components={ + "a": { + "phase_equilibrium_form": {("p2", "p3"): DummyPE}, + }, + }, + phases={ + "p1": { + "type": LiquidPhase, + "equation_of_state": DummyEoS, + }, + "p2": { + "type": VaporPhase, + "equation_of_state": DummyEoS, + }, + "p3": { + "type": LiquidPhase, + "equation_of_state": DummyEoS, + }, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + base_units=base_units, + phases_in_equilibrium=[("p2", "p3")], + phase_equilibrium_state={("p2", "p3"): DummyPE}, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + assert m.props[1]._get_critical_ref_phase() == "p3" + @pytest.mark.unit def test_get_critical_ref_phase_VX(self): # No vapor or liquid phases @@ -1800,8 +1844,7 @@ def build_critical_properties(b, *args, **kwargs): m.props[1]._critical_props() # Initialize critical props - initializer = ModularPropertiesInitializer() - initializer.initialize_critical_props(m.props[1]) + _initialize_critical_props(m.props[1]) assert ( value(m.props[1].compress_fact_crit) == 0.4 * 0.1 + 0.6 * 0.2 + 1e-8 * 0.3 From a7a0e92471b486e33e6053df51f85f7cff0e67f3 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Sun, 18 Feb 2024 18:20:08 -0500 Subject: [PATCH 08/15] Finishing initialization --- .../base/generic_property.py | 119 ++++++++++-------- .../properties/modular_properties/eos/ceos.py | 9 ++ .../modular_properties/eos/eos_base.py | 16 +++ .../eos/tests/test_ceos_PR.py | 68 +++++++--- 4 files changed, 143 insertions(+), 69 deletions(-) diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index 79bfebbf61..3d7004454b 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -1296,17 +1296,41 @@ def initialization_routine( k.inherent_equilibrium_constraint.deactivate() # --------------------------------------------------------------------- - # If present, initialize critical properties + # If present, initialize bubble and dew point and critical property calculations for k in model.values(): - # We only need to look for one critical property as it is all or nothing + T_units = k.params.get_metadata().default_units.TEMPERATURE + + # List of bubble and dew point constraints + cons_list = [ + "eq_pressure_dew", + "eq_pressure_bubble", + "eq_temperature_dew", + "eq_temperature_bubble", + "eq_mole_frac_tbub", + "eq_mole_frac_tdew", + "eq_mole_frac_pbub", + "eq_mole_frac_pdew", + "log_mole_frac_tbub_eqn", + "log_mole_frac_tdew_eqn", + "log_mole_frac_pbub_eqn", + "log_mole_frac_pdew_eqn", + "mole_frac_comp_eq", + "log_mole_frac_comp_eqn", + ] + + # Critical point with k.lock_attribute_creation_context(): + # Only need to look for one, as it is all-or-nothing if hasattr(k, "pressure_crit"): + # Initialize critical point properties _initialize_critical_props(k) + # Add critical point constraints to cons_list + ref_phase = k._get_critical_ref_phase() + p_config = k.params.get_phase(ref_phase).config + cons_list += ( + p_config.equation_of_state.list_critical_property_constraint_names() + ) - # --------------------------------------------------------------------- - # If present, initialize bubble and dew point calculations - for k in model.values(): - T_units = k.params.get_metadata().default_units.TEMPERATURE # Bubble temperature initialization if hasattr(k, "_mole_frac_tbub"): model._init_Tbub(k, T_units) @@ -1323,26 +1347,11 @@ def initialization_routine( if hasattr(k, "_mole_frac_pdew"): model._init_Pdew(k, T_units) - # Solve bubble and dew point constraints + # Solve bubble, dew and critical point constraints for c in k.component_objects(Constraint): # Deactivate all constraints not associated with bubble and dew - # points - if c.local_name not in ( - "eq_pressure_dew", - "eq_pressure_bubble", - "eq_temperature_dew", - "eq_temperature_bubble", - "eq_mole_frac_tbub", - "eq_mole_frac_tdew", - "eq_mole_frac_pbub", - "eq_mole_frac_pdew", - "log_mole_frac_tbub_eqn", - "log_mole_frac_tdew_eqn", - "log_mole_frac_pbub_eqn", - "log_mole_frac_pdew_eqn", - "mole_frac_comp_eq", - "log_mole_frac_comp_eqn", - ): + # points or critical points + if c.local_name not in cons_list: c.deactivate() # If StateBlock has active constraints (i.e. has bubble and/or dew @@ -1352,7 +1361,8 @@ def initialization_routine( if not degrees_of_freedom(b) == 0: raise InitializationError( f"{b.name} Unexpected degrees of freedom during " - f"initialization at bubble and dew point step: {degrees_of_freedom(b)}." + f"initialization at bubble and dew and critical point step: " + f"{degrees_of_freedom(b)}." ) with idaeslog.solver_log(solve_log, idaeslog.DEBUG) as slc: solve_strongly_connected_components( @@ -1361,7 +1371,7 @@ def initialization_routine( solve_kwds={"tee": slc.tee}, calc_var_kwds=self.config.calculate_variable_options, ) - init_log.info("Dew and bubble point initialization completed.") + init_log.info("Bubble, dew and critical point initialization completed.") # --------------------------------------------------------------------- # Calculate _teq if required @@ -1741,17 +1751,37 @@ def initialize( opt = get_solver(solver, optarg) # --------------------------------------------------------------------- - # If present, initialize critical properties + # If present, initialize bubble and dew point, and critical property calculations for k in blk.values(): - # We only need to look for one critical property as it is all or nothing + T_units = k.params.get_metadata().default_units.TEMPERATURE + + # List of bubble and dew point constraints + cons_list = [ + "eq_pressure_dew", + "eq_pressure_bubble", + "eq_temperature_dew", + "eq_temperature_bubble", + "eq_mole_frac_tbub", + "eq_mole_frac_tdew", + "eq_mole_frac_pbub", + "eq_mole_frac_pdew", + "log_mole_frac_tbub_eqn", + "log_mole_frac_tdew_eqn", + "log_mole_frac_pbub_eqn", + "log_mole_frac_pdew_eqn", + "mole_frac_comp_eq", + "log_mole_frac_comp_eqn", + ] + + # Critical point with k.lock_attribute_creation_context(): + # Only need to look for one, as it is all-or-nothing if hasattr(k, "pressure_crit"): + # Initialize critical point properties _initialize_critical_props(k) + # Add critical point constraints to cons_list + cons_list += k.list_critical_property_constraint_names() - # --------------------------------------------------------------------- - # If present, initialize bubble and dew point calculations - for k in blk.values(): - T_units = k.params.get_metadata().default_units.TEMPERATURE # Bubble temperature initialization if hasattr(k, "_mole_frac_tbub"): blk._init_Tbub(k, T_units) @@ -1768,26 +1798,11 @@ def initialize( if hasattr(k, "_mole_frac_pdew"): blk._init_Pdew(k, T_units) - # Solve bubble and dew point constraints + # Solve bubble, dew and critical point constraints for c in k.component_objects(Constraint): # Deactivate all constraints not associated with bubble and dew - # points - if c.local_name not in ( - "eq_pressure_dew", - "eq_pressure_bubble", - "eq_temperature_dew", - "eq_temperature_bubble", - "eq_mole_frac_tbub", - "eq_mole_frac_tdew", - "eq_mole_frac_pbub", - "eq_mole_frac_pdew", - "log_mole_frac_tbub_eqn", - "log_mole_frac_tdew_eqn", - "log_mole_frac_pbub_eqn", - "log_mole_frac_pdew_eqn", - "mole_frac_comp_eq", - "log_mole_frac_comp_eqn", - ): + # points or critical points + if c.local_name not in cons_list: c.deactivate() # If StateBlock has active constraints (i.e. has bubble and/or dew @@ -1801,12 +1816,12 @@ def initialize( if dof > 0: raise InitializationError( f"{blk.name} Unexpected degrees of freedom during " - f"initialization at bubble and dew point step: {dof}." + f"initialization at bubble, dew and critical point step: {dof}." ) with idaeslog.solver_log(solve_log, idaeslog.DEBUG) as slc: res = solve_indexed_blocks(opt, [blk], tee=slc.tee) init_log.info( - "Dew and bubble point initialization: {}.".format( + "Bubble, dew and critical point initialization: {}.".format( idaeslog.condition(res) ) ) diff --git a/idaes/models/properties/modular_properties/eos/ceos.py b/idaes/models/properties/modular_properties/eos/ceos.py index 002ff4f729..d419ef149f 100644 --- a/idaes/models/properties/modular_properties/eos/ceos.py +++ b/idaes/models/properties/modular_properties/eos/ceos.py @@ -1012,6 +1012,15 @@ def rule_B_crit(m): * m.dens_mol_crit ) + @staticmethod + def list_critical_property_constraint_names(): + return [ + "dens_mol_crit_eq", + "compress_fact_crit_eq", + "A_crit", + "B_crit", + ] + def _dZ_dT(blk, p): pobj = blk.params.get_phase(p) diff --git a/idaes/models/properties/modular_properties/eos/eos_base.py b/idaes/models/properties/modular_properties/eos/eos_base.py index 054a7ad2d1..52944a8314 100644 --- a/idaes/models/properties/modular_properties/eos/eos_base.py +++ b/idaes/models/properties/modular_properties/eos/eos_base.py @@ -52,8 +52,24 @@ def build_parameters(b): @staticmethod def build_critical_properties(b, ref_phase): + """ + This method is used to define the constraints used to calculate mixture + critical properties. + """ raise NotImplementedError(_msg(b, "build_critical_properties")) + @staticmethod + def list_critical_property_constraint_names(): + """ + If critical properties are supported, this method is required during + initialization to identify the constraints that are associated + with calculating the critical properties. + + Returns: + list of constraint names + """ + raise NotImplementedError(_msg(b, "list_critical_property_constraint_names")) + @staticmethod def get_vol_mol_pure(b, phase, comp, temperature): try: diff --git a/idaes/models/properties/modular_properties/eos/tests/test_ceos_PR.py b/idaes/models/properties/modular_properties/eos/tests/test_ceos_PR.py index 6b734b0f46..2240c0160e 100644 --- a/idaes/models/properties/modular_properties/eos/tests/test_ceos_PR.py +++ b/idaes/models/properties/modular_properties/eos/tests/test_ceos_PR.py @@ -35,12 +35,16 @@ from idaes.core import declare_process_block_class, LiquidPhase, VaporPhase, SolidPhase from idaes.models.properties.modular_properties.eos.ceos import Cubic, CubicType from idaes.models.properties.modular_properties.base.generic_property import ( + GenericParameterBlock, GenericParameterData, + ModularPropertiesInitializer, ) from idaes.core.util.exceptions import PropertyNotSupportedError, ConfigurationError from idaes.core.util.constants import Constants as const from idaes.models.properties.modular_properties.eos.ceos import cubic_roots_available from idaes.core.solvers import get_solver +from idaes.core.initialization.initializer_base import InitializationStatus +from idaes.models.properties.modular_properties.state_definitions import FTPx # Dummy method for property method calls @@ -1078,11 +1082,13 @@ def model(self): m = ConcreteModel() # Dummy params block - m.params = DummyParameterBlock( + m.params = GenericParameterBlock( components={ "a": { "parameter_data": { "omega": 0.1, + "compress_fact_crit": 1, + "dens_mol_crit": 1, "pressure_crit": 1e5, "temperature_crit": 100, } @@ -1090,6 +1096,8 @@ def model(self): "b": { "parameter_data": { "omega": 0.2, + "compress_fact_crit": 1, + "dens_mol_crit": 1, "pressure_crit": 2e5, "temperature_crit": 200, } @@ -1097,17 +1105,14 @@ def model(self): "c": { "parameter_data": { "omega": 0.3, + "compress_fact_crit": 1, + "dens_mol_crit": 1, "pressure_crit": 3e5, "temperature_crit": 300, } }, }, phases={ - "Vap": { - "type": LiquidPhase, - "equation_of_state": Cubic, - "equation_of_state_options": {"type": CubicType.PR}, - }, "Liq": { "type": LiquidPhase, "equation_of_state": Cubic, @@ -1121,7 +1126,7 @@ def model(self): "amount": pyunits.mol, "temperature": pyunits.K, }, - state_definition=modules[__name__], + state_definition=FTPx, pressure_ref=100000.0, temperature_ref=300, parameter_data={ @@ -1140,13 +1145,14 @@ def model(self): ) m.props = m.params.state_block_class( - [1], defined_state=False, parameters=m.params + [1], defined_state=True, parameters=m.params ) - # Add common variables - m.props[1].mole_frac_comp = Var( - m.params.component_list, initialize=0.5, bounds=(1e-12, 1) - ) + # Fix state variables + m.props[1].flow_mol.fix(1) + m.props[1].pressure.fix(1e5) + m.props[1].temperature.fix(300) + m.props[1].mole_frac_comp["a"].fix(0.4) m.props[1].mole_frac_comp["b"].fix(0.6) m.props[1].mole_frac_comp["c"].fix(1e-8) @@ -1325,14 +1331,42 @@ def test_build_critical_props(self, model): assert_units_consistent(model) + # Check list_critical_property_constraint_names + con_list = Cubic.list_critical_property_constraint_names() + + assert con_list == [ + "dens_mol_crit_eq", + "compress_fact_crit_eq", + "A_crit", + "B_crit", + ] + + # Check that names match constraints that were built + for c in con_list: + assert isinstance(getattr(model.props[1], c), Constraint) + @pytest.mark.component @pytest.mark.solver - def test_build_critical_props(self, model): - # Initialize model - model.props[1].compress_fact_crit.set_value(1) - model.props[1].pressure_crit.set_value(1.6e5) - model.props[1].temperature_crit.set_value(160) + def test_initialize_with_critical_props(self, model): + initializer = ModularPropertiesInitializer() + + status = initializer.initialize(model.props) + assert status == InitializationStatus.Ok + + # Confirm results + assert value(model.props[1].compress_fact_crit) == pytest.approx( + 0.321379, rel=1e-5 + ) + assert value(model.props[1].pressure_crit) == pytest.approx(1.58138e5, rel=1e-5) + assert value(model.props[1].temperature_crit) == pytest.approx( + 158.138, rel=1e-5 + ) + assert value(model.props[1].dens_mol_crit) == pytest.approx(374.238, rel=1e-5) + + @pytest.mark.component + @pytest.mark.solver + def test_solve_critical_props(self, model): # Solve model solver = get_solver() res = solver.solve(model, tee=True) From ea8ff6e3a0dc2b05fd9a14d1810715bd7984a1e6 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 20 Feb 2024 09:06:35 -0500 Subject: [PATCH 09/15] Fixing pylint and removing ideal crit props --- idaes/core/base/property_meta.py | 4 +- .../modular_properties/eos/ideal.py | 38 +--------- .../eos/tests/test_ideal.py | 73 +------------------ 3 files changed, 3 insertions(+), 112 deletions(-) diff --git a/idaes/core/base/property_meta.py b/idaes/core/base/property_meta.py index e06b26eb7d..502e11b15d 100644 --- a/idaes/core/base/property_meta.py +++ b/idaes/core/base/property_meta.py @@ -267,9 +267,7 @@ def VOLUME_MOLE(self): # Backward compatibility name @property def MOLAR_VOLUME(self): - msg = ( - f"The unit name MOLAR_VOLUME is being deprecated in " "favor of VOLUME_MOL." - ) + msg = "The unit name MOLAR_VOLUME is being deprecated in favor of VOLUME_MOL." deprecation_warning(msg=msg, logger=_log, version="2.3.0", remove_in="3.0.0") return self.VOLUME_MOLE diff --git a/idaes/models/properties/modular_properties/eos/ideal.py b/idaes/models/properties/modular_properties/eos/ideal.py index 75ae312641..8986e7286a 100644 --- a/idaes/models/properties/modular_properties/eos/ideal.py +++ b/idaes/models/properties/modular_properties/eos/ideal.py @@ -21,7 +21,7 @@ # TODO: Look into protected access issues # pylint: disable=protected-access -from pyomo.environ import Constraint, Expression, log, units +from pyomo.environ import Expression, log from idaes.core import Apparent from idaes.core.util.exceptions import ConfigurationError, PropertyNotSupportedError @@ -34,7 +34,6 @@ log_henry_pressure, ) from .eos_base import EoSBase -from idaes.core.util.constants import Constants as CONST # TODO: Add support for ideal solids @@ -49,41 +48,6 @@ def common(b, pobj): # No common components required for ideal property calculations pass - @staticmethod - def build_critical_properties(b, ref_phase): - base_units = b.params.get_metadata().default_units - - # Use Kay's method for critical pressure and temperature - def p_crit_rule(b): - return b.pressure_crit == sum( - b.mole_frac_comp[j] * b.params.get_component(j).pressure_crit - for j in b.component_list - ) - - b.pressure_crit_constraint = Constraint(rule=p_crit_rule) - - def t_crit_rule(b): - return b.temperature_crit == sum( - b.mole_frac_comp[j] * b.params.get_component(j).temperature_crit - for j in b.component_list - ) - - b.temperature_crit_constraint = Constraint(rule=t_crit_rule) - - # Fix critical compressibility factor to 1 - user can override if necessary - b.compress_fact_crit.fix(1) - - def rho_crit_rule(b): - return b.pressure_crit == units.convert( - b.compress_fact_crit - * CONST.gas_constant - * b.temperature_crit - * b.dens_mol_crit, - to_units=base_units.PRESSURE, - ) - - b.dens_mol_crit_constraint = Constraint(rule=rho_crit_rule) - @staticmethod def calculate_scaling_factors(b, pobj): pass diff --git a/idaes/models/properties/modular_properties/eos/tests/test_ideal.py b/idaes/models/properties/modular_properties/eos/tests/test_ideal.py index 76eb169765..087f8a2b3b 100644 --- a/idaes/models/properties/modular_properties/eos/tests/test_ideal.py +++ b/idaes/models/properties/modular_properties/eos/tests/test_ideal.py @@ -18,7 +18,7 @@ import pytest from sys import modules -from pyomo.environ import ConcreteModel, Constraint, Var, units as pyunits, value +from pyomo.environ import ConcreteModel, Var, units as pyunits, value from pyomo.util.check_units import assert_units_equivalent from idaes.core import ( @@ -685,74 +685,3 @@ def test_vol_mol_phase(): + 1 / 42.0 * m.props[1].mole_frac_phase_comp["Liq", "b"] + 42.0 * m.props[1].mole_frac_phase_comp["Liq", "c"] ) - - -@pytest.mark.unit -def test_critical_props(): - m = ConcreteModel() - - # Dummy params block - m.params = DummyParameterBlock( - components={ - "a": {"parameter_data": {"pressure_crit": 1e5, "temperature_crit": 100}}, - "b": {"parameter_data": {"pressure_crit": 2e5, "temperature_crit": 200}}, - "c": {"parameter_data": {"pressure_crit": 3e5, "temperature_crit": 300}}, - }, - phases={ - "Vap": {"type": VaporPhase, "equation_of_state": Ideal}, - "Liq": {"type": LiquidPhase, "equation_of_state": Ideal}, - }, - base_units={ - "time": pyunits.s, - "length": pyunits.m, - "mass": pyunits.kg, - "amount": pyunits.mol, - "temperature": pyunits.K, - }, - state_definition=modules[__name__], - pressure_ref=100000.0, - temperature_ref=300, - ) - - m.props = m.params.state_block_class([1], defined_state=False, parameters=m.params) - - # Add common variables - m.props[1].mole_frac_comp = Var(m.params.component_list, initialize=0.5) - - m.props[1]._critical_props() - - assert isinstance(m.props[1].compress_fact_crit, Var) - assert isinstance(m.props[1].dens_mol_crit, Var) - assert isinstance(m.props[1].pressure_crit, Var) - assert isinstance(m.props[1].temperature_crit, Var) - - assert isinstance(m.props[1].pressure_crit_constraint, Constraint) - assert str(m.props[1].pressure_crit_constraint.expr) == str( - m.props[1].pressure_crit - == m.props[1].mole_frac_comp["a"] * m.params.a.pressure_crit - + m.props[1].mole_frac_comp["b"] * m.params.b.pressure_crit - + m.props[1].mole_frac_comp["c"] * m.params.c.pressure_crit - ) - - assert isinstance(m.props[1].temperature_crit_constraint, Constraint) - assert str(m.props[1].temperature_crit_constraint.expr) == str( - m.props[1].temperature_crit - == m.props[1].mole_frac_comp["a"] * m.params.a.temperature_crit - + m.props[1].mole_frac_comp["b"] * m.params.b.temperature_crit - + m.props[1].mole_frac_comp["c"] * m.params.c.temperature_crit - ) - - assert m.props[1].compress_fact_crit.fixed - assert value(m.props[1].compress_fact_crit) == 1 - - assert isinstance(m.props[1].dens_mol_crit_constraint, Constraint) - assert str(m.props[1].dens_mol_crit_constraint.expr) == str( - m.props[1].pressure_crit - == pyunits.convert( - m.props[1].compress_fact_crit - * const.gas_constant - * m.props[1].temperature_crit - * m.props[1].dens_mol_crit, - to_units=pyunits.kg / pyunits.m / pyunits.s**2, - ) - ) From 66db51f3c6bcaf1d17212aabe32f7311e0247e7f Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 20 Feb 2024 12:23:41 -0500 Subject: [PATCH 10/15] Docs and extra error message --- .../property_package/general/eos/cubic.rst | 27 +++++++++-- .../property_package/general/eos/ideal.rst | 5 ++ .../property_package/general/phase_def.rst | 18 +++++++- .../base/generic_property.py | 46 +++++++++---------- 4 files changed, 67 insertions(+), 29 deletions(-) diff --git a/docs/explanations/components/property_package/general/eos/cubic.rst b/docs/explanations/components/property_package/general/eos/cubic.rst index 1893289ff0..8123d8d801 100644 --- a/docs/explanations/components/property_package/general/eos/cubic.rst +++ b/docs/explanations/components/property_package/general/eos/cubic.rst @@ -67,6 +67,7 @@ Cubic equations of state require the following parameters to be defined: 1. `omega` (Pitzer acentricity factor) needs to be defined for each component (in the `parameter_data` for each component). 2. `kappa` (binary interaction parameters) needs to be defined for each component pair in the system. This parameter needs to be defined in the general `parameter_data` argument for the overall property package (as it can be used in multiple phases). +3. Component critical properties (`compresss_fact_crit`, `dens_mol_crit`, `pressure_crit`, and `temperature_crit`) are often required as well. Calculation of Properties ------------------------- @@ -77,6 +78,26 @@ Many thermophysical properties are calculated using an ideal and residual term, The residual term is derived from the partial derivatives of the cubic equation of state, whilst the ideal term is determined using pure component properties for the ideal gas phase defined for each component. +Critical Properties of Mixtures +------------------------------- + +The critical properties of the mixture are calculated by solving the cubic equation of state at the critical point using the following constraints. + +.. math:: \Omega_A R^2 T_c^2 = a_{c,m}P_c +.. math:: \Omega_BRT_c = b_{c,m}P_c +.. math:: 0 = Z_c^3 - (1+B_c-uB_c)Z_c^2 + (A_c-uB_c-(u-w)B_c^2)Z_c - A_cB_c-wB_c^2-wB_c^3 +.. math:: P_c = Z_cRT_c\rho_{c} + + +with the following expressions: + +.. math:: a_{c,j} = \frac{\Omega_AR^2T_{c,j}^2}{P_{c, j}}\alpha_j(T_c) +.. math:: a_{c,m} = \sum_i{\sum_j{y_iy_j(a_{c,i}a_{c,j})^{1/2}(1-\kappa_{ij})}} +.. math:: b_{c,m} = \sum_i{y_ib_i} +.. math:: A_c = \frac{a_{c,m}P_c}{R^2T_c^2} +.. math:: B_c = \frac{b_{c,m}P_c}{RT_c} + + Mass Density by Phase --------------------- @@ -109,14 +130,14 @@ Component Molar Enthalpy by Phase Component molar enthalpies by phase are calculated using the pure component method provided by the users in the property package configuration arguments. -Molar Isobaric Heat Capcity (:math:`C_p`) +Molar Isobaric Heat Capacity (:math:`C_p`) ----------------------------------------- -The ideal molar isobaric heat capcity term is calculated from the weighted sum of the (ideal) component molar isobaric heat capacity: +The ideal molar isobaric heat capacity term is calculated from the weighted sum of the (ideal) component molar isobaric heat capacity: .. math:: C_{p, ig}^0 = \sum_j y_j C_{p, ig, j} -The residual molar isobaric heat capcity term is given by: +The residual molar isobaric heat capacity term is given by: .. math:: C_p^r = R \left[ T \left(\frac{\partial Z}{\partial T}\right)_P + Z - 1 \right] + \frac{ T \frac{d^2a_m}{dT^2}}{\sqrt{u^2 - 4w} \cdot b_m} \ln \left[ \frac{2Z + uB + \sqrt{u^2 - 4w} B}{2Z + uB - \sqrt{u^2 - 4w} B} \right] .. math:: + \left(a_m - T \frac{da_m}{dT}\right) \cdot \frac{B}{b_m} \cdot \frac{\left(\frac{\partial Z}{\partial T}\right)_P + \frac{Z}{T}}{Z^2 + Z uB + wB^2} diff --git a/docs/explanations/components/property_package/general/eos/ideal.rst b/docs/explanations/components/property_package/general/eos/ideal.rst index f4b1a37d36..38fc9cd32b 100644 --- a/docs/explanations/components/property_package/general/eos/ideal.rst +++ b/docs/explanations/components/property_package/general/eos/ideal.rst @@ -9,6 +9,11 @@ Introduction Ideal behavior represents the simplest possible equation of state that ensures thermodynamic consistency between different properties. +Critical Properties of Mixtures +------------------------------- + +The Ideal Equation of State module does not support the calculation of mixture critical properties as the ideal equation of state cannot represent the critical point. + Mass Density by Phase --------------------- diff --git a/docs/explanations/components/property_package/general/phase_def.rst b/docs/explanations/components/property_package/general/phase_def.rst index 22d64e9f17..b80e0a2a4e 100644 --- a/docs/explanations/components/property_package/general/phase_def.rst +++ b/docs/explanations/components/property_package/general/phase_def.rst @@ -34,7 +34,7 @@ Equations of state (or equivalent methods) describe the relationship between dif A wide range of equations of states are available in literature for different applications and levels of rigor, and the IDAES Generic Property Package Framework provides a number of prebuilt modules for users, which are listed below. -Equation of state packages may allow for user options (e.g. choosing a specific type of cubic equation of state). The options are set using the `equation_of_state_options` argument, and the options available are described in the documentation of each equation of state module. +Equation of state packages may allow for user options (e.g., choosing a specific type of cubic equation of state). The options are set using the `equation_of_state_options` argument, and the options available are described in the documentation of each equation of state module. Equation of State Libraries """"""""""""""""""""""""""" @@ -44,6 +44,20 @@ Equation of State Libraries eos/ideal eos/cubic + +Critical Properties of Mixtures +""""""""""""""""""""""""""""""" + +Calculation of the critical properties of mixtures depends on the equation of state being used. As the Modular Property Package framework allows users to specify different equations of state for each phase, the following logic is used to determine which equation of state to use for calculating critical properties. + +1. If a vapor-liquid equilibrium pair is defined in the `"phases_in_equilibrium"` configuration argument, then the liquid phase from this pair is used (IDAES generally assumed supercritical fluids are liquid like). +2. If no vapor-liquid equilibrium pair is defined, then the first liquid phase define is used. +3. If no liquid phases are defined then the vapor phase is used (it is assumed there will only be one vapor phase). +4. If no vapor phase is defined, then a `PropertyPackageError` is returned as there is no suitable phase for calculating critical properties. + +Note that not all Equations of State are suitable for calculating critical properties as many cannot represent the critical conditions. The Ideal equation of State is one common example of an equation of state that DOES NOT support critical properties. + +During initialization, the critical properties of the mixture are approximated using the mole fraction weighted sum of the component critical properties (note that users must provide values for all component critical properties (i.e., `compress_fact_crit`, `dens_mol_crit`, `pressure_crit` and `temperature_crit`) for initialization if critical properties are to be calculated. Phase-Specific Parameter ^^^^^^^^^^^^^^^^^^^^^^^^ @@ -53,6 +67,6 @@ In some cases, a property package may include parameters which are specific to a Phases with Partial Component Lists ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ -In many applications a mixture will contain species that only appear in a single phase (either by nature or assumption). Common examples include crystalline solids and non-condensable gases. The IDAES Generic Property Package Framework provides support for these behaviors and allows users to specify phase-specific component lists (i.e. a list of components which appear in a given phase). +In many applications a mixture will contain species that only appear in a single phase (either by nature or assumption). Common examples include crystalline solids and non-condensable gases. The IDAES Generic Property Package Framework provides support for these behaviors and allows users to specify phase-specific component lists (i.e., a list of components which appear in a given phase). This is done by providing a phase with a `component_list` argument, which provides a `list` of component names which appear in the phase. The framework automatically validates the `component_list` argument to ensure that it is a sub-set of the master component list for the property package, and will inform the user if an unrecognized component is included. If a phase is not provided with a `component_list` argument it is assumed that all components defined in the master component list may be present in the phase. diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index 3d7004454b..350fc4c5dc 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -5157,27 +5157,25 @@ def rule_log_mole_frac(b, p1, p2, j): def _initialize_critical_props(state_data): params = state_data.params # Use mole weighted sum of component critical properties - state_data.compress_fact_crit.set_value( - sum( - state_data.mole_frac_comp[j] * params.get_component(j).compress_fact_crit - for j in state_data.component_list - ) - ) - state_data.dens_mol_crit.set_value( - sum( - state_data.mole_frac_comp[j] * params.get_component(j).dens_mol_crit - for j in state_data.component_list - ) - ) - state_data.pressure_crit.set_value( - sum( - state_data.mole_frac_comp[j] * params.get_component(j).pressure_crit - for j in state_data.component_list - ) - ) - state_data.temperature_crit.set_value( - sum( - state_data.mole_frac_comp[j] * params.get_component(j).temperature_crit - for j in state_data.component_list - ) - ) + crit_props = [ + "compress_fact_crit", + "dens_mol_crit", + "pressure_crit", + "temperature_crit", + ] + + for prop in crit_props: + try: + getattr(state_data, prop).set_value( + sum( + state_data.mole_frac_comp[j] + * getattr(params.get_component(j), prop) + for j in state_data.component_list + ) + ) + except AttributeError: + raise AttributeError( + f"Missing attribute found when initializing {prop}. " + f"Make sure you have provided values for {prop} in all Component " + "declarations." + ) From 46253e380e1283f183acd59158b2fa9a7b8f2f8d Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 20 Feb 2024 12:26:23 -0500 Subject: [PATCH 11/15] Test for new error message --- .../base/tests/test_generic_property.py | 93 +++++++++++++++++++ 1 file changed, 93 insertions(+) diff --git a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py index cc7fbef1ea..8f00f5c26b 100644 --- a/idaes/models/properties/modular_properties/base/tests/test_generic_property.py +++ b/idaes/models/properties/modular_properties/base/tests/test_generic_property.py @@ -1852,3 +1852,96 @@ def build_critical_properties(b, *args, **kwargs): assert value(m.props[1].dens_mol_crit) == 0.4 * 1 + 0.6 * 2 + 1e-8 * 3 assert value(m.props[1].pressure_crit) == 0.4 * 1e5 + 0.6 * 2e5 + 1e-8 * 3e5 assert value(m.props[1].temperature_crit) == 0.4 * 1e2 + 0.6 * 2e2 + 1e-8 * 3e2 + + @pytest.mark.unit + def test_initialize_critical_props_missing_value(self): + m = ConcreteModel() + + class DummyEoS2(DummyEoS): + @staticmethod + def build_critical_properties(b, *args, **kwargs): + b._dummy_crit_executed = True + + # Dummy params block + m.params = DummyParameterBlock( + components={ + "a": { + "parameter_data": { + "dens_mol_crit": 1, # Missing Z_crit + "pressure_crit": 1e5, + "temperature_crit": 100, + } + }, + "b": { + "parameter_data": { + "compress_fact_crit": 0.2, + "dens_mol_crit": 2, + "pressure_crit": 2e5, + "temperature_crit": 200, + } + }, + "c": { + "parameter_data": { + "compress_fact_crit": 0.3, + "dens_mol_crit": 3, + "pressure_crit": 3e5, + "temperature_crit": 300, + } + }, + }, + phases={ + "Vap": { + "type": LiquidPhase, + "equation_of_state": DummyEoS2, + }, + "Liq": { + "type": LiquidPhase, + "equation_of_state": DummyEoS2, + }, + }, + base_units={ + "time": pyunits.s, + "length": pyunits.m, + "mass": pyunits.kg, + "amount": pyunits.mol, + "temperature": pyunits.K, + }, + state_definition=modules[__name__], + pressure_ref=100000.0, + temperature_ref=300, + parameter_data={ + "PR_kappa": { + ("a", "a"): 0.0, + ("a", "b"): 0.0, + ("a", "c"): 0.0, + ("b", "a"): 0.0, + ("b", "b"): 0.0, + ("b", "c"): 0.0, + ("c", "a"): 0.0, + ("c", "b"): 0.0, + ("c", "c"): 0.0, + }, + }, + ) + + m.props = m.params.build_state_block([1], defined_state=False) + + # Add common variables + m.props[1].mole_frac_comp = Var( + m.params.component_list, initialize=0.5, bounds=(1e-12, 1) + ) + m.props[1].mole_frac_comp["a"].fix(0.4) + m.props[1].mole_frac_comp["b"].fix(0.6) + m.props[1].mole_frac_comp["c"].fix(1e-8) + + # Build critical props + m.props[1]._critical_props() + + # Initialize critical props + with pytest.raises( + AttributeError, + match="Missing attribute found when initializing compress_fact_crit. " + "Make sure you have provided values for compress_fact_crit in all " + "Component declarations.", + ): + _initialize_critical_props(m.props[1]) From 16520bdf54c7422322ea69e19c54067c8d407463 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 20 Feb 2024 12:39:42 -0500 Subject: [PATCH 12/15] Fixing typo --- .../components/property_package/general/eos/cubic.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/explanations/components/property_package/general/eos/cubic.rst b/docs/explanations/components/property_package/general/eos/cubic.rst index 8123d8d801..773c3904a1 100644 --- a/docs/explanations/components/property_package/general/eos/cubic.rst +++ b/docs/explanations/components/property_package/general/eos/cubic.rst @@ -67,7 +67,7 @@ Cubic equations of state require the following parameters to be defined: 1. `omega` (Pitzer acentricity factor) needs to be defined for each component (in the `parameter_data` for each component). 2. `kappa` (binary interaction parameters) needs to be defined for each component pair in the system. This parameter needs to be defined in the general `parameter_data` argument for the overall property package (as it can be used in multiple phases). -3. Component critical properties (`compresss_fact_crit`, `dens_mol_crit`, `pressure_crit`, and `temperature_crit`) are often required as well. +3. Component critical properties (`compress_fact_crit`, `dens_mol_crit`, `pressure_crit`, and `temperature_crit`) are often required as well. Calculation of Properties ------------------------- From 7bdff7000bea8743cf5b202759b118b9d2157d94 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Tue, 20 Feb 2024 12:53:34 -0500 Subject: [PATCH 13/15] Fixing Sphinx and pylint errors --- .../components/property_package/general/eos/cubic.rst | 2 +- idaes/models/properties/modular_properties/eos/eos_base.py | 6 +++++- 2 files changed, 6 insertions(+), 2 deletions(-) diff --git a/docs/explanations/components/property_package/general/eos/cubic.rst b/docs/explanations/components/property_package/general/eos/cubic.rst index 773c3904a1..115896bb98 100644 --- a/docs/explanations/components/property_package/general/eos/cubic.rst +++ b/docs/explanations/components/property_package/general/eos/cubic.rst @@ -131,7 +131,7 @@ Component Molar Enthalpy by Phase Component molar enthalpies by phase are calculated using the pure component method provided by the users in the property package configuration arguments. Molar Isobaric Heat Capacity (:math:`C_p`) ------------------------------------------ +------------------------------------------ The ideal molar isobaric heat capacity term is calculated from the weighted sum of the (ideal) component molar isobaric heat capacity: diff --git a/idaes/models/properties/modular_properties/eos/eos_base.py b/idaes/models/properties/modular_properties/eos/eos_base.py index 52944a8314..37db8a1f40 100644 --- a/idaes/models/properties/modular_properties/eos/eos_base.py +++ b/idaes/models/properties/modular_properties/eos/eos_base.py @@ -68,7 +68,11 @@ def list_critical_property_constraint_names(): Returns: list of constraint names """ - raise NotImplementedError(_msg(b, "list_critical_property_constraint_names")) + raise NotImplementedError( + "Equation of State module has not implemented a method for " + "list_critical_property_constraint_names. Please contact the " + "EoS developer or use a different module." + ) @staticmethod def get_vol_mol_pure(b, phase, comp, temperature): From c1a1f7522485e70526955d985b12a0dbfad5b2f2 Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Wed, 21 Feb 2024 09:27:28 -0500 Subject: [PATCH 14/15] Address PR review --- .../property_package/general/phase_def.rst | 6 ++--- .../base/generic_property.py | 24 +++++++++---------- 2 files changed, 15 insertions(+), 15 deletions(-) diff --git a/docs/explanations/components/property_package/general/phase_def.rst b/docs/explanations/components/property_package/general/phase_def.rst index b80e0a2a4e..e6814ad839 100644 --- a/docs/explanations/components/property_package/general/phase_def.rst +++ b/docs/explanations/components/property_package/general/phase_def.rst @@ -50,14 +50,14 @@ Critical Properties of Mixtures Calculation of the critical properties of mixtures depends on the equation of state being used. As the Modular Property Package framework allows users to specify different equations of state for each phase, the following logic is used to determine which equation of state to use for calculating critical properties. -1. If a vapor-liquid equilibrium pair is defined in the `"phases_in_equilibrium"` configuration argument, then the liquid phase from this pair is used (IDAES generally assumed supercritical fluids are liquid like). -2. If no vapor-liquid equilibrium pair is defined, then the first liquid phase define is used. +1. If a vapor-liquid equilibrium pair is defined in the `"phases_in_equilibrium"` configuration argument, then the liquid phase from this pair is used (IDAES generally assumed supercritical fluids are liquid-like). +2. If no vapor-liquid equilibrium pair is defined, then the first liquid phase defined is used. 3. If no liquid phases are defined then the vapor phase is used (it is assumed there will only be one vapor phase). 4. If no vapor phase is defined, then a `PropertyPackageError` is returned as there is no suitable phase for calculating critical properties. Note that not all Equations of State are suitable for calculating critical properties as many cannot represent the critical conditions. The Ideal equation of State is one common example of an equation of state that DOES NOT support critical properties. -During initialization, the critical properties of the mixture are approximated using the mole fraction weighted sum of the component critical properties (note that users must provide values for all component critical properties (i.e., `compress_fact_crit`, `dens_mol_crit`, `pressure_crit` and `temperature_crit`) for initialization if critical properties are to be calculated. +During initialization, the critical properties of the mixture are approximated using the mole fraction weighted sum of the component critical properties (note that users must provide values for all component critical properties (i.e., `compress_fact_crit`, `dens_mol_crit`, `pressure_crit` and `temperature_crit`) for initialization if critical properties are to be calculated. Note that it is up to the user to ensure these values are consistent (if necessary). Phase-Specific Parameter ^^^^^^^^^^^^^^^^^^^^^^^^ diff --git a/idaes/models/properties/modular_properties/base/generic_property.py b/idaes/models/properties/modular_properties/base/generic_property.py index 350fc4c5dc..1a21f239c2 100644 --- a/idaes/models/properties/modular_properties/base/generic_property.py +++ b/idaes/models/properties/modular_properties/base/generic_property.py @@ -1296,7 +1296,7 @@ def initialization_routine( k.inherent_equilibrium_constraint.deactivate() # --------------------------------------------------------------------- - # If present, initialize bubble and dew point and critical property calculations + # If present, initialize bubble, dew, and critical point calculations for k in model.values(): T_units = k.params.get_metadata().default_units.TEMPERATURE @@ -1347,21 +1347,21 @@ def initialization_routine( if hasattr(k, "_mole_frac_pdew"): model._init_Pdew(k, T_units) - # Solve bubble, dew and critical point constraints + # Solve bubble, dew, and critical point constraints for c in k.component_objects(Constraint): # Deactivate all constraints not associated with bubble and dew # points or critical points if c.local_name not in cons_list: c.deactivate() - # If StateBlock has active constraints (i.e. has bubble and/or dew + # If StateBlock has active constraints (i.e. has bubble, dew, or critical # point calculations), solve the block to converge these for b in model.values(): if number_activated_constraints(b) > 0: if not degrees_of_freedom(b) == 0: raise InitializationError( f"{b.name} Unexpected degrees of freedom during " - f"initialization at bubble and dew and critical point step: " + f"initialization at bubble, dew, and critical point step: " f"{degrees_of_freedom(b)}." ) with idaeslog.solver_log(solve_log, idaeslog.DEBUG) as slc: @@ -1371,7 +1371,7 @@ def initialization_routine( solve_kwds={"tee": slc.tee}, calc_var_kwds=self.config.calculate_variable_options, ) - init_log.info("Bubble, dew and critical point initialization completed.") + init_log.info("Bubble, dew, and critical point initialization completed.") # --------------------------------------------------------------------- # Calculate _teq if required @@ -1751,7 +1751,7 @@ def initialize( opt = get_solver(solver, optarg) # --------------------------------------------------------------------- - # If present, initialize bubble and dew point, and critical property calculations + # If present, initialize bubble, dew , and critical point calculations for k in blk.values(): T_units = k.params.get_metadata().default_units.TEMPERATURE @@ -1798,14 +1798,14 @@ def initialize( if hasattr(k, "_mole_frac_pdew"): blk._init_Pdew(k, T_units) - # Solve bubble, dew and critical point constraints + # Solve bubble, dew, and critical point constraints for c in k.component_objects(Constraint): - # Deactivate all constraints not associated with bubble and dew - # points or critical points + # Deactivate all constraints not associated with bubble, dew, + # or critical points if c.local_name not in cons_list: c.deactivate() - # If StateBlock has active constraints (i.e. has bubble and/or dew + # If StateBlock has active constraints (i.e. has bubble, dew, or critical # point calculations), solve the block to converge these n_cons = 0 dof = 0 @@ -1816,12 +1816,12 @@ def initialize( if dof > 0: raise InitializationError( f"{blk.name} Unexpected degrees of freedom during " - f"initialization at bubble, dew and critical point step: {dof}." + f"initialization at bubble, dew, and critical point step: {dof}." ) with idaeslog.solver_log(solve_log, idaeslog.DEBUG) as slc: res = solve_indexed_blocks(opt, [blk], tee=slc.tee) init_log.info( - "Bubble, dew and critical point initialization: {}.".format( + "Bubble, dew, and critical point initialization: {}.".format( idaeslog.condition(res) ) ) From c4fff4f80140991972c26bc8dc783eaf897d7dfd Mon Sep 17 00:00:00 2001 From: Andrew Lee Date: Wed, 21 Feb 2024 10:41:08 -0500 Subject: [PATCH 15/15] Fixing typo --- .../components/property_package/general/phase_def.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/explanations/components/property_package/general/phase_def.rst b/docs/explanations/components/property_package/general/phase_def.rst index e6814ad839..d7df89afb3 100644 --- a/docs/explanations/components/property_package/general/phase_def.rst +++ b/docs/explanations/components/property_package/general/phase_def.rst @@ -50,7 +50,7 @@ Critical Properties of Mixtures Calculation of the critical properties of mixtures depends on the equation of state being used. As the Modular Property Package framework allows users to specify different equations of state for each phase, the following logic is used to determine which equation of state to use for calculating critical properties. -1. If a vapor-liquid equilibrium pair is defined in the `"phases_in_equilibrium"` configuration argument, then the liquid phase from this pair is used (IDAES generally assumed supercritical fluids are liquid-like). +1. If a vapor-liquid equilibrium pair is defined in the `"phases_in_equilibrium"` configuration argument, then the liquid phase from this pair is used (IDAES generally assumes supercritical fluids are liquid-like). 2. If no vapor-liquid equilibrium pair is defined, then the first liquid phase defined is used. 3. If no liquid phases are defined then the vapor phase is used (it is assumed there will only be one vapor phase). 4. If no vapor phase is defined, then a `PropertyPackageError` is returned as there is no suitable phase for calculating critical properties.