diff --git a/docs/technical_reference/unit_models/index.rst b/docs/technical_reference/unit_models/index.rst index dd8508b40a..5794fd2051 100644 --- a/docs/technical_reference/unit_models/index.rst +++ b/docs/technical_reference/unit_models/index.rst @@ -26,6 +26,7 @@ Unit Models membrane_distillation_1D mvc nanofiltration_ZO + nanofiltration_0D nanofiltration_dspmde_0D osmotically_assisted_reverse_osmosis_0D osmotically_assisted_reverse_osmosis_1D diff --git a/docs/technical_reference/unit_models/nanofiltration_0D.rst b/docs/technical_reference/unit_models/nanofiltration_0D.rst new file mode 100644 index 0000000000..2262d22248 --- /dev/null +++ b/docs/technical_reference/unit_models/nanofiltration_0D.rst @@ -0,0 +1,70 @@ +Nanofiltration (0D) +==================== +This nanofiltration (NF) unit model is suitable for non-predictive case studies where the user wishes to specify, constant rejection fractions +whilst preserving electroneutrality in the permeate stream (optional). A single recovery fraction is assumed for all solvents, whilst +individual rejection fractions can be set for all solutes, with the exception of one ion specified for maintaining +electroneutrality. Solutes are assigned default rejection value by grouping them into either two categories, passing and excluded, with +default rejection values for each category. Solutes are categorised using either a user-provided list of species which freely pass the +membrane, or by assuming all neutral and monovalent species pass the membrane. + +Retentate pressure is assumed to be related to feed pressure with an optional pressure drop whilst permeate pressure is assumed to be +a degree of freedom, and temperature equality is assumed. Temperature and pressure constraints can be removed with configuration arguments. + +This model assumes supports a single liquid phase only and assumes steady-state. + +.. index:: + pair: watertap.unit_models.nanofiltration_0D;nanofiltration_0D + +.. currentmodule:: watertap.unit_models.nanofiltration_0D + +Degrees of Freedom +------------------ +The ``Nanofiltration0D`` model has the following degrees of freedom + + * inlet state + * permeate pressure + * solvent recovery fraction (``solvent_recovery``) + * solute rejection fractions (``rejection_comp``) for all solutes EXCEPT the one identified as the electroneutrality species. + +The following additional degrees of freedom may exist depending on configuration options + + * retentate pressure drop (``deltaP``) + +Model Structure +--------------- +The Nanofiltration0D model consists of three state blocks (``properties_in``, ``properties_retentate`` and ``properties_permeate``). + +.. _NF0D_variables: + +Variables +--------- + +.. csv-table:: + :header: "Description", "Symbol", "Variable Name", "Index", "Units" + + "Solvent recovery", ":math:`Q_{solvent}`", "solvent_recovery", None, ":math:`\text{dimensionless}`" + "Solute rejection", ":math:`R_j`", "rejection_comp", [j], ":math:`\text{dimensionless}`" + "Retentate pressure drop", ":math:`\deltaP`", "deltaP", [t], ":math:`\text{pressure}`" + +.. _NF0D_equations: + +Equations +--------- + +Here :math:`F` represents component flowrate and :math:`Z_j` is the charge on species j. + +.. csv-table:: + :header: "Description", "Equation" + + "Component material balances", ":math:`F_{in, j} = F_{retentate, j} + F_{permeate, j}`" + "Solvent recovery",":math:`Q_{solvent} \times F_{in, j} = F_{permeate_j}`" + "Solute rejection",":math:`R_j = 1 - \frac{C_{permeate, j}}{C_{feed_j}}`" + "Permeate electroneutrality",":math:`0 = F_{permeate, j} \times Z_{j}`" + "Retentate pressure balance",":math:`P_{in} = P_{retentate} - \deltaP`" + "Retentate temperature equality", ":math:`T_{in} = T_{retentate}`" + "Permeate temperature equality", ":math:`T_{in} = T_{permeate}`" + +Class Documentation +------------------- + +* :mod:`watertap.unit_models.nanofiltration_0D` diff --git a/setup.py b/setup.py index 9c30ceb1d5..a1489aef86 100644 --- a/setup.py +++ b/setup.py @@ -67,7 +67,7 @@ python_requires=">=3.9", install_requires=[ # primary requirements for unit and property models - "idaes-pse >=2.7.0,<2.8.0rc0", + "idaes-pse @ https://github.com/watertap-org/idaes-pse/archive/2.7.0.dev.watertap.2025.01.21.zip", "pyomo>=6.6.1", "watertap-solvers", "pyyaml", # watertap.core.wt_database diff --git a/watertap/property_models/multicomp_aq_sol_prop_pack.py b/watertap/property_models/multicomp_aq_sol_prop_pack.py index 63d6926efa..c8f76199e8 100644 --- a/watertap/property_models/multicomp_aq_sol_prop_pack.py +++ b/watertap/property_models/multicomp_aq_sol_prop_pack.py @@ -24,8 +24,6 @@ # -add viscosity as func of temp and concentration # Import Python libraries -import idaes.logger as idaeslog - from enum import Enum, auto # Import Pyomo libraries @@ -49,6 +47,7 @@ from pyomo.core.base.units_container import InconsistentUnitsError # Import IDAES cores +import idaes.logger as idaeslog from idaes.core import ( declare_process_block_class, MaterialFlowBasis, @@ -67,7 +66,6 @@ solve_indexed_blocks, ) from idaes.core.util.misc import add_object_reference -from watertap.core.solvers import get_solver from idaes.core.util.model_statistics import ( degrees_of_freedom, number_unfixed_variables, @@ -78,6 +76,10 @@ PropertyPackageError, ) import idaes.core.util.scaling as iscale +from idaes.core.scaling import CustomScalerBase, ConstraintScalingScheme +from idaes.core.initialization import BlockTriangularizationInitializer + +from watertap.core.solvers import get_solver from watertap.core.util.scaling import transform_property_constraints from watertap.core.util.chemistry import ( get_charge, @@ -120,6 +122,272 @@ class TransportNumberCalculation(Enum): ElectricalMobility = auto() +class MCASScaler(CustomScalerBase): + DEFAULT_SCALING_FACTORS = { + "act_coeff_phase_comp": 1, + "debye_huckel_constant": 100, + "dens_mass_phase": 1e-3, + "dens_mass_solvent": 1e-3, + "diffus_phase_comp": 1e10, + "enth_mass_phase": 1e-5, + "flow_mass_phase_comp": 1e3, + "flow_mol_phase_comp": 1e2, + "pressure_sat": 1e-3, + "trans_num_phase_comp": 10, + "visc_k_phase": 1e6, + } + + def variable_scaling_routine( + self, model, overwrite: bool = False, submodel_scalers: dict = None + ): + self.scale_variable_by_units(model.pressure, overwrite=overwrite) + self.scale_variable_by_units(model.temperature, overwrite=overwrite) + + if model.is_property_constructed("flow_mol_phase_comp"): + for j in model.component_list: + if j == "H2O": + self.set_variable_scaling_factor( + model.flow_mol_phase_comp["Liq", j], 0.1, overwrite=overwrite + ) + else: + self.scale_variable_by_default( + model.flow_mol_phase_comp["Liq", j], overwrite=overwrite + ) + + if model.is_property_constructed("flow_mass_phase_comp"): + for j in model.component_list: + if j == "H2O": + self.set_variable_scaling_factor( + model.flow_mass_phase_comp["Liq", j], 10, overwrite=overwrite + ) + else: + self.scale_variable_by_default( + model.flow_mass_phase_comp["Liq", j], overwrite=overwrite + ) + + if model.is_property_constructed("flow_equiv_phase_comp"): + for (p, j), v in model.flow_equiv_phase_comp.items(): + self.scale_variable_by_component( + target_variable=v, + scaling_component=model.flow_mass_phase_comp[p, j], + overwrite=overwrite, + ) + + for n in self.DEFAULT_SCALING_FACTORS.keys(): + if model.is_property_constructed(n): + c = model.find_component(n) + for v in c.values(): + self.scale_variable_by_default(v, overwrite=overwrite) + + if model.is_property_constructed("mole_frac_phase_comp"): + for (p, j), v in model.mole_frac_phase_comp.items(): + if j == "H2O": + self.set_variable_scaling_factor(v, 1, overwrite=overwrite) + else: + sf = self.get_scaling_factor( + model.flow_mol_phase_comp[p, j] + ) / self.get_scaling_factor(model.flow_mol_phase_comp["Liq", "H2O"]) + self.set_variable_scaling_factor(v, sf, overwrite=overwrite) + + if model.is_property_constructed("mass_frac_phase_comp"): + for (p, j), v in model.mass_frac_phase_comp.items(): + if j == "H2O": + self.set_variable_scaling_factor(v, 1, overwrite=overwrite) + else: + sf = self.get_scaling_factor( + model.flow_mass_phase_comp[p, j] + ) / self.get_scaling_factor( + model.flow_mass_phase_comp["Liq", "H2O"] + ) + self.set_variable_scaling_factor(v, sf, overwrite=overwrite) + + if model.is_property_constructed("mass_frac_phase_comp"): + for (p, j), v in model.mass_frac_phase_comp.items(): + if j == "H2O": + self.set_variable_scaling_factor(v, 1, overwrite=overwrite) + else: + sf = self.get_scaling_factor( + model.flow_mass_phase_comp[p, j] + ) / self.get_scaling_factor( + model.flow_mass_phase_comp["Liq", "H2O"] + ) + self.set_variable_scaling_factor(v, sf, overwrite=overwrite) + + if model.is_property_constructed("conc_mass_phase_comp"): + for (p, j), v in model.conc_mass_phase_comp.items(): + sf_dens = self.get_scaling_factor(model.dens_mass_phase[p]) + if j == "H2O": + # solvents typically have a mass fraction between 0.5-1 + self.set_variable_scaling_factor( + v, + sf_dens, + overwrite=overwrite, + ) + else: + sf_x = self.get_scaling_factor(model.mass_frac_phase_comp[p, j]) + self.set_variable_scaling_factor( + v, + sf_dens * sf_x, + overwrite=overwrite, + ) + + if model.is_property_constructed("conc_mol_phase_comp"): + for (p, j), v in model.conc_mol_phase_comp.items(): + if j == "H2O": + self.set_variable_scaling_factor( + v, + 1e-4, + overwrite=overwrite, + ) + else: + sf_x = self.get_scaling_factor(model.conc_mass_phase_comp[p, j]) + # Multiply by MW as scaling factor is inverse of value + self.set_variable_scaling_factor( + v, + value(sf_x * model.mw_comp[j]), + overwrite=overwrite, + ) + + if model.is_property_constructed("conc_equiv_phase_comp"): + for k, v in model.conc_equiv_phase_comp.items(): + self.scale_variable_by_component( + target_variable=v, + scaling_component=model.conc_mol_phase_comp[k], + overwrite=overwrite, + ) + + if model.is_property_constructed("pressure_osm_phase"): + sf_conc_mol = ( + sum( + self.get_scaling_factor(model.conc_mol_phase_comp["Liq", j]) ** -1 + for j in model.params.solute_set + ) + ) ** -1 + # R*T is generally of order 1e3 + self.set_variable_scaling_factor( + model.pressure_osm_phase["Liq"], sf_conc_mol * 1e-3, overwrite=overwrite + ) + + if model.is_property_constructed("elec_mobility_phase_comp"): + for ind, v in model.elec_mobility_phase_comp.items(): + if ( + model.params.config.elec_mobility_calculation + == ElectricalMobilityCalculation.EinsteinRelation + ): + sf = self.get_scaling_factor(model.diffus_phase_comp[ind]) / 40 + else: + sf = model.params.config.elec_mobility_data[ind] ** -1 + self.set_variable_scaling_factor(v, sf, overwrite=overwrite) + + if model.is_property_constructed("equiv_conductivity_phase"): + for v in model.equiv_conductivity_phase.values(): + if ( + model.params.config.equiv_conductivity_calculation + == EquivalentConductivityCalculation.ElectricalMobility + ): + sf = ( + 1e-5 + * sum( + self.get_scaling_factor( + model.elec_mobility_phase_comp["Liq", j] + ) + ** -1 + * self.get_scaling_factor( + model.conc_mol_phase_comp["Liq", j] + ) + ** -1 + for j in model.params.ion_set + ) + ** -1 + / sum( + self.get_scaling_factor(model.conc_mol_phase_comp["Liq", j]) + ** -1 + for j in model.params.cation_set + ) + ** -1 + ) + else: + sf = model.params.config.equiv_conductivity_phase_data[ind] ** -1 + self.set_variable_scaling_factor(v, sf, overwrite=overwrite) + + if model.is_property_constructed("elec_cond_phase"): + for ind, v in model.elec_cond_phase.items(): + sf_equiv_cond_phase = self.get_scaling_factor( + model.equiv_conductivity_phase[ind] + ) + sf_conc_mol_z = ( + sum( + self.get_scaling_factor(model.conc_mol_phase_comp["Liq", j]) + ** -1 + for j in model.params.cation_set + ) + ** -1 + ) + sf = sf_equiv_cond_phase * sf_conc_mol_z + self.set_variable_scaling_factor(v, sf, overwrite=overwrite) + + if model.is_property_constructed("flow_vol_phase"): + sf = ( + self.get_scaling_factor(model.flow_mol_phase_comp["Liq", "H2O"]) + / value(model.mw_comp["H2O"]) + / self.get_scaling_factor(model.dens_mass_phase["Liq"]) + ) + self.set_variable_scaling_factor( + model.flow_vol_phase["Liq"], sf, overwrite=overwrite + ) + + if model.is_property_constructed("molality_phase_comp"): + for (p, j), v in model.molality_phase_comp.items(): + sf = ( + self.get_scaling_factor(model.flow_mol_phase_comp[p, j]) + / self.get_scaling_factor(model.flow_mol_phase_comp[p, "H2O"]) + * value(model.mw_comp["H2O"]) + ) + self.set_variable_scaling_factor(v, sf, overwrite=overwrite) + + if model.is_property_constructed("ionic_strength_molal"): + sf = ( + sum( + 0.5 + * self.get_scaling_factor(model.molality_phase_comp["Liq", j]) ** -1 + * value(model.charge_comp[j]) ** 2 + for j in model.params.solute_set + ) + ** -1 + ) + self.set_variable_scaling_factor( + model.ionic_strength_molal, sf, overwrite=overwrite + ) + + if model.is_property_constructed("total_hardness"): + if value(model.total_hardness) == 0: + sf = 1 + else: + sf = 1 / value(model.total_hardness) + self.set_variable_scaling_factor( + model.total_hardness, sf, overwrite=overwrite + ) + + if model.is_property_constructed("total_dissolved_solids"): + if value(model.total_dissolved_solids) == 0: + sf = 1 + else: + sf = 1 / value(model.total_dissolved_solids) + self.set_variable_scaling_factor( + model.total_dissolved_solids, sf, overwrite=overwrite + ) + + def constraint_scaling_routine( + self, model, overwrite: bool = False, submodel_scalers: dict = None + ): + for c in model.component_data_objects(Constraint, descend_into=True): + self.scale_constraint_by_nominal_value( + c, + scheme=ConstraintScalingScheme.inverseMaximum, + overwrite=overwrite, + ) + + @declare_process_block_class("MCASParameterBlock") class MCASParameterData(PhysicalParameterBlock): CONFIG = PhysicalParameterBlock.CONFIG() @@ -689,6 +957,19 @@ class _MCASStateBlock(StateBlock): than individual elements of indexed Property Blocks. """ + default_initializer = BlockTriangularizationInitializer + default_scaler = MCASScaler + + def fix_initialization_states(self): + """ + Fixes state variables for state blocks. + + Returns: + None + """ + # Fix state variables + fix_state_vars(self) + def initialize( self, state_args=None, diff --git a/watertap/property_models/tests/test_multicomp_aq_sol_prop_pack.py b/watertap/property_models/tests/test_multicomp_aq_sol_prop_pack.py index 5ce08ad6c4..8262b31a37 100644 --- a/watertap/property_models/tests/test_multicomp_aq_sol_prop_pack.py +++ b/watertap/property_models/tests/test_multicomp_aq_sol_prop_pack.py @@ -12,6 +12,7 @@ import re import pytest + from pyomo.environ import ( ConcreteModel, assert_optimal_termination, @@ -23,6 +24,11 @@ Var, Constraint, ) +from pyomo.util.check_units import assert_units_consistent, assert_units_equivalent + +# Imports from idaes core +from idaes.core.base.components import Solvent, Solute, Cation, Anion +from idaes.core.base.phases import PhaseType as PT from idaes.core import ( FlowsheetBlock, MaterialFlowBasis, @@ -31,19 +37,6 @@ EnergyBalanceType, ) from idaes.core.util.scaling import calculate_scaling_factors, get_scaling_factor - -from pyomo.util.check_units import assert_units_consistent, assert_units_equivalent -from watertap.property_models.multicomp_aq_sol_prop_pack import ( - MCASParameterBlock, - MCASStateBlock, - ActivityCoefficientModel, - DensityCalculation, - DiffusivityCalculation, - ElectricalMobilityCalculation, - EquivalentConductivityCalculation, - TransportNumberCalculation, -) -from watertap.core.util.initialization import check_dof from idaes.core.util.model_statistics import ( number_variables, number_total_constraints, @@ -56,12 +49,11 @@ set_scaling_factor, ) from idaes.core.util.exceptions import ConfigurationError -from watertap.property_models.tests.property_test_harness import PropertyAttributeError -from watertap.core.solvers import get_solver - -# Imports from idaes core -from idaes.core.base.components import Solvent, Solute, Cation, Anion -from idaes.core.base.phases import PhaseType as PT +from idaes.core.scaling import set_scaling_factor +from idaes.core.initialization import ( + BlockTriangularizationInitializer, + InitializationStatus, +) # Imports from idaes generic models from idaes.models.properties.modular_properties.pure.ConstantProperties import Constant @@ -78,11 +70,41 @@ from idaes.models.unit_models.mixer import MixingType import idaes.logger as idaeslog +from watertap.property_models.multicomp_aq_sol_prop_pack import ( + MCASParameterBlock, + MCASStateBlock, + ActivityCoefficientModel, + DensityCalculation, + DiffusivityCalculation, + ElectricalMobilityCalculation, + EquivalentConductivityCalculation, + TransportNumberCalculation, + MCASScaler, +) +from watertap.core.util.initialization import check_dof +from watertap.property_models.tests.property_test_harness import PropertyAttributeError +from watertap.core.solvers import get_solver solver = get_solver() # ----------------------------------------------------------------------------- +@pytest.mark.unit +def test_default_initializer_scaler(): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["A"], + mw_data={"H2O": 0.018, "A": 0.01}, + charge={"A": 0}, + ) + + m.fs.state = m.fs.properties.build_state_block(m.fs.time) + + assert m.fs.state.default_initializer is BlockTriangularizationInitializer + assert m.fs.state.default_scaler is MCASScaler + + @pytest.mark.unit def test_config(): m = ConcreteModel() @@ -2485,3 +2507,625 @@ def test_automatic_charge_mw_population(): for comp, val in test_vals.items(): assert value(m.fs.stream2[0].mw_comp[comp]) == val[0] assert value(m.fs.stream2[0].charge_comp[comp]) == val[1] + + +class TestMCASInitializer: + @pytest.mark.component + def test_initialize(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Ca_2+", "SO4_2-", "Na_+", "Cl_-", "Mg_2+"], + material_flow_basis=MaterialFlowBasis.mass, + diffusivity_data={ + ("Liq", "Ca_2+"): 7.92e-10, + ("Liq", "SO4_2-"): 1.06e-09, + ("Liq", "Na_+"): 1.33e-09, + ("Liq", "Cl_-"): 2.03e-09, + ("Liq", "Mg_2+"): 7.06e-10, + }, + mw_data={ + "H2O": 0.018, + "Na_+": 0.023, + "Ca_2+": 0.04, + "Mg_2+": 0.024, + "Cl_-": 0.035, + "SO4_2-": 0.096, + }, + stokes_radius_data={ + "Na_+": 1.84e-10, + "Ca_2+": 3.09e-10, + "Mg_2+": 3.47e-10, + "Cl_-": 1.21e-10, + "SO4_2-": 2.3e-10, + }, + charge={"Na_+": 1, "Ca_2+": 2, "Mg_2+": 2, "Cl_-": -1, "SO4_2-": -2}, + elec_mobility_calculation=ElectricalMobilityCalculation.EinsteinRelation, + density_calculation=DensityCalculation.seawater, + activity_coefficient_model=ActivityCoefficientModel.davies, + ) + + m.fs.stream = stream = m.fs.properties.build_state_block( + [0], defined_state=True + ) + + mass_flow_in = 1 * pyunits.kg / pyunits.s + feed_mass_frac = { + "Na_+": 11122e-6, + "Ca_2+": 382e-6, + "Mg_2+": 1394e-6, + "SO4_2-": 2136e-6, + "Cl_-": 20300e-6, + } + for ion, x in feed_mass_frac.items(): + mass_comp_flow = x * pyunits.kg / pyunits.kg * mass_flow_in + + stream[0].flow_mass_phase_comp["Liq", ion].fix(mass_comp_flow) + + H2O_mass_frac = 1 - sum(x for x in feed_mass_frac.values()) + + stream[0].flow_mass_phase_comp["Liq", "H2O"].fix(H2O_mass_frac) + stream[0].temperature.fix(298.15) + stream[0].pressure.fix(101325) + + stream[0].assert_electroneutrality( + defined_state=True, tol=1e-2, adjust_by_ion="Cl_-" + ) + + metadata = m.fs.properties.get_metadata().properties + for v in metadata.list_supported_properties(): + getattr(stream[0], v.name) + assert stream[0].is_property_constructed(v.name) + + initializer = stream.default_initializer() + initializer.initialize(stream) + assert initializer.summary[stream]["status"] == InitializationStatus.Ok + + +class TestMCASScaler: + @pytest.mark.component + def test_scaler(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Ca_2+", "SO4_2-", "Na_+", "Cl_-", "Mg_2+"], + material_flow_basis=MaterialFlowBasis.mass, + diffusivity_data={ + ("Liq", "Ca_2+"): 7.92e-10, + ("Liq", "SO4_2-"): 1.06e-09, + ("Liq", "Na_+"): 1.33e-09, + ("Liq", "Cl_-"): 2.03e-09, + ("Liq", "Mg_2+"): 7.06e-10, + }, + mw_data={ + "H2O": 0.018, + "Na_+": 0.023, + "Ca_2+": 0.04, + "Mg_2+": 0.024, + "Cl_-": 0.035, + "SO4_2-": 0.096, + }, + stokes_radius_data={ + "Na_+": 1.84e-10, + "Ca_2+": 3.09e-10, + "Mg_2+": 3.47e-10, + "Cl_-": 1.21e-10, + "SO4_2-": 2.3e-10, + }, + charge={"Na_+": 1, "Ca_2+": 2, "Mg_2+": 2, "Cl_-": -1, "SO4_2-": -2}, + elec_mobility_calculation=ElectricalMobilityCalculation.EinsteinRelation, + density_calculation=DensityCalculation.seawater, + activity_coefficient_model=ActivityCoefficientModel.davies, + ) + + m.fs.stream = stream = m.fs.properties.build_state_block( + [0], defined_state=True + ) + + mass_flow_in = 1 * pyunits.kg / pyunits.s + feed_mass_frac = { + "Na_+": 11122e-6, + "Ca_2+": 382e-6, + "Mg_2+": 1394e-6, + "SO4_2-": 2136e-6, + "Cl_-": 20300e-6, + } + for ion, x in feed_mass_frac.items(): + mass_comp_flow = x * pyunits.kg / pyunits.kg * mass_flow_in + + stream[0].flow_mass_phase_comp["Liq", ion].fix(mass_comp_flow) + + H2O_mass_frac = 1 - sum(x for x in feed_mass_frac.values()) + + stream[0].flow_mass_phase_comp["Liq", "H2O"].fix(H2O_mass_frac) + stream[0].temperature.fix(298.15) + stream[0].pressure.fix(101325) + + stream[0].assert_electroneutrality( + defined_state=True, tol=1e-2, adjust_by_ion="Cl_-" + ) + + metadata = m.fs.properties.get_metadata().properties + for v in metadata.list_supported_properties(): + getattr(stream[0], v.name) + assert stream[0].is_property_constructed(v.name) + + initializer = stream.default_initializer() + initializer.initialize(stream) + assert initializer.summary[stream]["status"] == InitializationStatus.Ok + + scaler = stream.default_scaler() + + set_scaling_factor(stream[0].flow_mass_phase_comp["Liq", "H2O"], 1) + set_scaling_factor(stream[0].flow_mass_phase_comp["Liq", "Ca_2+"], 1e4) + set_scaling_factor(stream[0].flow_mass_phase_comp["Liq", "SO4_2-"], 1e3) + set_scaling_factor(stream[0].flow_mass_phase_comp["Liq", "Na_+"], 1e2) + set_scaling_factor(stream[0].flow_mass_phase_comp["Liq", "Cl_-"], 1e2) + set_scaling_factor(stream[0].flow_mass_phase_comp["Liq", "Mg_2+"], 1e3) + + scaler.scale_model(stream[0]) + + assert len(m.fs.stream[0].scaling_factor) == 154 + + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].temperature + ] == pytest.approx(0.01, rel=1e-3) + assert m.fs.stream[0].scaling_factor[m.fs.stream[0].pressure] == pytest.approx( + 0.00001, rel=1e-3 + ) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mass_phase_comp["Liq", "H2O"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mass_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(10000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mass_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mass_phase_comp["Liq", "Na_+"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mass_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mass_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mol_phase_comp["Liq", "H2O"] + ] == pytest.approx(0.0001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mol_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(0.4, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mol_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(0.096, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mol_phase_comp["Liq", "Na_+"] + ] == pytest.approx(0.0023, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mol_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(0.0035, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mol_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(0.024, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mass_phase_comp["Liq", "H2O"] + ] == pytest.approx(0.001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mass_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mass_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mass_phase_comp["Liq", "Na_+"] + ] == pytest.approx(0.1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mass_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(0.1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_mass_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].dens_mass_phase["Liq"] + ] == pytest.approx(0.001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mass_frac_phase_comp["Liq", "H2O"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mass_frac_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(10000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mass_frac_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mass_frac_phase_comp["Liq", "Na_+"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mass_frac_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mass_frac_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].dens_mass_solvent + ] == pytest.approx(0.001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].act_coeff_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].act_coeff_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].act_coeff_phase_comp["Liq", "Na_+"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].act_coeff_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].act_coeff_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].ionic_strength_molal + ] == pytest.approx(2.571, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].molality_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].molality_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].molality_phase_comp["Liq", "Na_+"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].molality_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].molality_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mol_phase_comp["Liq", "H2O"] + ] == pytest.approx(0.1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mol_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mol_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mol_phase_comp["Liq", "Na_+"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mol_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_mol_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].debye_huckel_constant + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].enth_mass_phase["Liq"] + ] == pytest.approx(0.00001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].total_dissolved_solids + ] == pytest.approx(0.0000278, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_vol_phase["Liq"] + ] == pytest.approx(5556, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mole_frac_phase_comp["Liq", "H2O"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mole_frac_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mole_frac_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mole_frac_phase_comp["Liq", "Na_+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mole_frac_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].mole_frac_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].pressure_osm_phase["Liq"] + ] == pytest.approx(0.00000129, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].pressure_sat + ] == pytest.approx(0.001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].visc_k_phase["Liq"] + ] == pytest.approx(1000000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_equiv_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(10000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_equiv_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_equiv_phase_comp["Liq", "Na_+"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_equiv_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].flow_equiv_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_equiv_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(0.4, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_equiv_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(0.096, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_equiv_phase_comp["Liq", "Na_+"] + ] == pytest.approx(0.0023, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_equiv_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(0.0035, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].conc_equiv_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(0.024, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].equiv_conductivity_phase["Liq"] + ] == pytest.approx(1545, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].elec_mobility_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(250000000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].elec_mobility_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(250000000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].elec_mobility_phase_comp["Liq", "Na_+"] + ] == pytest.approx(250000000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].elec_mobility_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(250000000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].elec_mobility_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(250000000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].elec_cond_phase["Liq"] + ] == pytest.approx(3.225, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].trans_num_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].trans_num_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].trans_num_phase_comp["Liq", "Na_+"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].trans_num_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].trans_num_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].total_hardness + ] == pytest.approx(0.0001443, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mol_phase_comp["Liq", "H2O"] + ] == pytest.approx(0.001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mol_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mol_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mol_phase_comp["Liq", "Na_+"] + ] == pytest.approx(0.1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mol_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(0.1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mol_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mass_phase_comp["Liq", "H2O"] + ] == pytest.approx(0.001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mass_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mass_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mass_phase_comp["Liq", "Na_+"] + ] == pytest.approx(0.1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mass_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(0.1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_mass_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_dens_mass_phase["Liq"] + ] == pytest.approx(0.001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mass_frac_phase_comp["Liq", "H2O"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mass_frac_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(2617, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mass_frac_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(468.1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mass_frac_phase_comp["Liq", "Na_+"] + ] == pytest.approx(89.89, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mass_frac_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(49.74, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mass_frac_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(717.2, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_dens_mass_solvent + ] == pytest.approx(0.001, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_act_coeff_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(65.09, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_act_coeff_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(65.09, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_act_coeff_phase_comp["Liq", "Na_+"] + ] == pytest.approx(260.4, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_act_coeff_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(260.4, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_act_coeff_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(65.09, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_ionic_strength_molal + ] == pytest.approx(2.571, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_molality_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_molality_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_molality_phase_comp["Liq", "Na_+"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_molality_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_molality_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(18, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_mol_phase_comp["Liq", "H2O"] + ] == pytest.approx(1.037, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_mol_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(2500, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_mol_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(468.2, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_mol_phase_comp["Liq", "Na_+"] + ] == pytest.approx(89.91, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_mol_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(49.75, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_mol_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(717.4, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_debye_huckel_constant + ] == pytest.approx(64.32, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_enth_mass_phase["Liq"] + ] == pytest.approx(0.0000007982, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_total_dissolved_solids + ] == pytest.approx(0.0000278, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_vol_phase["Liq"] + ] == pytest.approx(5489, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mole_frac_phase_comp["Liq", "H2O"] + ] == pytest.approx(1, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mole_frac_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mole_frac_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mole_frac_phase_comp["Liq", "Na_+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mole_frac_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_mole_frac_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_pressure_osm_phase["Liq"] + ] == pytest.approx(0.0000009278, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_pressure_sat + ] == pytest.approx(0.0003216, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_visc_k_phase["Liq"] + ] == pytest.approx(1000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_equiv_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(50, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_equiv_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(50, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_equiv_phase_comp["Liq", "Na_+"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_equiv_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(100, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_flow_equiv_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(50, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_equiv_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(0.2, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_equiv_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(0.048, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_equiv_phase_comp["Liq", "Na_+"] + ] == pytest.approx(0.0023, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_equiv_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(0.0035, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_conc_equiv_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(0.012, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_equiv_conductivity_phase["Liq"] + ] == pytest.approx(1545, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_elec_mobility_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(16220000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_elec_mobility_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(12120000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_elec_mobility_phase_comp["Liq", "Na_+"] + ] == pytest.approx(19320000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_elec_mobility_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(12660000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_elec_mobility_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(18200000, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_elec_cond_phase["Liq"] + ] == pytest.approx(3.225, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_trans_num_phase_comp["Liq", "Ca_2+"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_trans_num_phase_comp["Liq", "SO4_2-"] + ] == pytest.approx(10, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_trans_num_phase_comp["Liq", "Na_+"] + ] == pytest.approx(1.908, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_trans_num_phase_comp["Liq", "Cl_-"] + ] == pytest.approx(2.904, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_trans_num_phase_comp["Liq", "Mg_2+"] + ] == pytest.approx(9.956, rel=1e-3) + assert m.fs.stream[0].scaling_factor[ + m.fs.stream[0].eq_total_hardness + ] == pytest.approx(0.0001443, rel=1e-3) diff --git a/watertap/unit_models/__init__.py b/watertap/unit_models/__init__.py index 6b86ce0036..aa83d68d1d 100644 --- a/watertap/unit_models/__init__.py +++ b/watertap/unit_models/__init__.py @@ -16,6 +16,7 @@ from .osmotically_assisted_reverse_osmosis_0D import OsmoticallyAssistedReverseOsmosis0D from .nanofiltration_ZO import NanofiltrationZO from .nanofiltration_DSPMDE_0D import NanofiltrationDSPMDE0D +from .nanofiltration_0D import Nanofiltration0D from .pressure_exchanger import PressureExchanger from .pressure_changer import Pump, EnergyRecoveryDevice from .crystallizer import Crystallization diff --git a/watertap/unit_models/nanofiltration_0D.py b/watertap/unit_models/nanofiltration_0D.py new file mode 100644 index 0000000000..1603c244e1 --- /dev/null +++ b/watertap/unit_models/nanofiltration_0D.py @@ -0,0 +1,732 @@ +################################################################################# +# WaterTAP Copyright (c) 2020-2024, The Regents of the University of California, +# through Lawrence Berkeley National Laboratory, Oak Ridge National Laboratory, +# National Renewable Energy Laboratory, and National Energy Technology +# Laboratory (subject to receipt of any required approvals from the U.S. Dept. +# of Energy). All rights reserved. +# +# Please see the files COPYRIGHT.md and LICENSE.md for full copyright and license +# information, respectively. These files are also available online at the URL +# "https://github.com/watertap-org/watertap/" +################################################################################# +""" +Simple zero-dimensional Nanofiltration unit model. +""" + +from pyomo.environ import Block, Constraint, Set, value, Var +from pyomo.common.config import ConfigDict, ConfigValue, In + +from idaes.core import ( + declare_process_block_class, + EnergyBalanceType, + MomentumBalanceType, + useDefault, + UnitModelBlockData, +) +from idaes.core.initialization import ModularInitializerBase +from idaes.core.scaling import CustomScalerBase, ConstraintScalingScheme +from idaes.core.util.config import is_physical_parameter_block +from idaes.core.util.exceptions import ConfigurationError +import idaes.logger as idaeslog + + +_log = idaeslog.getLogger(__name__) + + +class Nanofiltration0DInitializer(ModularInitializerBase): + """ + Initializer for 0D Nanofiltration models. + + """ + + CONFIG = ModularInitializerBase.CONFIG() + + def initialize_main_model( + self, + model: Block, + ): + """ + Initialization routine for main Nanofiltration0D models. + + Args: + model: Pyomo Block to be initialized. + copy_inlet_state: bool (default=False). Whether to copy inlet state to other states or not. + Copying will generally be faster, but inlet states may not contain + all properties required elsewhere. + + Returns: + Pyomo solver results object. + + """ + # Get loggers + init_log = idaeslog.getInitLogger( + model.name, self.get_output_level(), tag="unit" + ) + solve_log = idaeslog.getSolveLogger( + model.name, self.get_output_level(), tag="unit" + ) + + # Create solver + solver = self._get_solver() + + # --------------------------------------------------------------------- + # Initialize inlet state + prop_init = self.get_submodel_initializer(model.properties_in) + + if prop_init is not None: + prop_init.initialize( + model=model.properties_in, + output_level=self.get_output_level(), + ) + else: + raise ValueError( + "No Initializer found for property package. Please set provide " + "sub-model initializers or assign a default Initializer." + ) + + init_log.info_high("Initialization Step 1 (Inlet State) Complete.") + + # --------------------------------------------------------------------- + # Estimate outlet states from inlet and initialize + for t, in_state in model.properties_in.items(): + state_vars = in_state.define_state_vars() + state_vars_ret = model.properties_retentate[t].define_state_vars() + state_vars_per = model.properties_permeate[t].define_state_vars() + + for n, sv in state_vars.items(): + sv_ret = state_vars_ret[n] + sv_per = state_vars_per[n] + + if sv.local_name == "pressure" and not sv_ret.fixed: + # For pressure, only retentate is linked to inlet + if hasattr(model, "deltaP"): + sv_ret.set_value(sv + model.deltaP[t]) + else: + sv_ret.set_value(sv) + elif "flow" in sv.local_name: + self._init_flow(model, t, in_state, sv, sv_ret, sv_per) + elif any( + sv.local_name.startswith(i) for i in ["mass_frac", "mole_frac"] + ): + self._init_frac(model, t, in_state, sv, sv_ret, sv_per) + elif sv.local_name.startswith("conc"): + self._init_conc(model, t, in_state, sv, sv_ret, sv_per) + else: + # For everything else, assume outlet similar to inlet + for k, svd in sv.items(): + if not sv_ret[k].fixed: + sv_ret[k].set_value(svd) + if not sv_per[k].fixed: + sv_per[k].set_value(svd) + + prop_init.initialize( + model=model.properties_retentate, + output_level=self.get_output_level(), + ) + prop_init.initialize( + model=model.properties_permeate, + output_level=self.get_output_level(), + ) + + init_log.info_high("Initialization Step 2 (Outlet States) Complete.") + + # --------------------------------------------------------------------- + # Solve full model + with idaeslog.solver_log(solve_log, idaeslog.DEBUG) as slc: + res = solver.solve(model, tee=slc.tee) + + init_log.info("Initialization Completed, {}".format(idaeslog.condition(res))) + + return res + + def _init_conc(self, model, t, in_state, sv, sv_ret, sv_per): + # Component indexed flow - need to apply split fractions + for k, svd in sv.items(): + if isinstance(k, str): + # Single indexing set, component is the index + j = k + else: + # Component is always the last index + j = k[-1] + + # Determine split based on type + comp = in_state.params.get_component(j) + + if comp.is_solvent(): + # Assume solvent concentrations unchanged + sv_ret[k].set_value(svd) + sv_per[k].set_value(svd) + else: + if j == model.config.electroneutrality_ion: + # Guess electroneutrality ion will be based on solvent recovery + split = model.recovery_solvent[t] + else: + # For initialization, assume density is roughly constant and rejection can be used as split fraction + split = 1 - model.rejection_comp[t, j] + + if not sv_ret[k].fixed: + sv_ret[k].set_value(svd * (1 - split)) + if not sv_per[k].fixed: + sv_per[k].set_value(svd * split) + + def _init_flow(self, model, t, in_state, sv, sv_ret, sv_per): + # Determine if it is a component flow or not + if sv.local_name.endswith("_comp"): + # Component indexed flow - need to apply split fractions + for k, svd in sv.items(): + if isinstance(k, str): + # Single indexing set, component is the index + j = k + else: + # Component is always the last index + j = k[-1] + + # Determine split based on type + comp = in_state.params.get_component(j) + + if comp.is_solvent() or j == model.config.electroneutrality_ion: + # Assume electroneutrality ion will have similar separation to solvents + split = model.recovery_solvent[t] + else: + # For initialization, assume density is roughly constant and rejection can be used as split fraction + split = 1 - model.rejection_comp[t, j] + + if not sv_ret[k].fixed: + sv_ret[k].set_value(svd * (1 - split)) + if not sv_per[k].fixed: + sv_per[k].set_value(svd * split) + else: + # Assume a total flow basis, and use solvent recovery + split = model.recovery_solvent[t] + for k, svd in sv.items(): + if not sv_ret[k].fixed: + sv_ret[k].set_value(svd * (1 - split)) + if not sv_per[k].fixed: + sv_per[k].set_value(svd * split) + + def _init_frac(self, model, t, in_state, sv, sv_ret, sv_per): + # First need to iterate over all indices to collect normalized flows + # Assume a basis of 1 mass or mole unit + # Also, only single phase property packages supported, so we only need to + # worry about the component index + nom_ret_comp_flow = {} + nom_per_comp_flow = {} + + for k, svd in sv.items(): + if isinstance(k, str): + # Single indexing set, component is the index + j = k + else: + # Component is always the last index + j = k[-1] + + # Determine split based on type + comp = in_state.params.get_component(j) + + if comp.is_solvent() or j == model.config.electroneutrality_ion: + # Assume electroneutrality ion will have similar separation to solvents + split = model.recovery_solvent[t] + else: + # For initialization, assume density is roughly constant and rejection can be used as split fraction + split = 1 - model.rejection_comp[t, j] + + nom_ret_comp_flow[j] = value(svd * (1 - split)) + nom_per_comp_flow[j] = value(svd * split) + + nom_ret_tot_flow = sum(f for f in nom_ret_comp_flow.values()) + nom_per_tot_flow = sum(f for f in nom_per_comp_flow.values()) + + # # Next set value based on fraction of nominal total flow + for k, svd in sv.items(): + if isinstance(k, str): + # Single indexing set, component is the index + j = k + else: + # Component is always the last index + j = k[-1] + + if not sv_ret[k].fixed: + sv_ret[k].set_value(nom_ret_comp_flow[j] / nom_ret_tot_flow) + if not sv_per[k].fixed: + sv_per[k].set_value(nom_per_comp_flow[j] / nom_per_tot_flow) + + +class Nanofiltration0DScaler(CustomScalerBase): + """ + Scaler class for Nanofiltration0D models. + """ + + DEFAULT_SCALING_FACTORS = { + "deltaP": 1e-4, + "rejection_comp": 1e2, + "recovery_solvent": 10, + } + + def variable_scaling_routine( + self, model, overwrite: bool = False, submodel_scalers: dict = None + ): + """ + Apply variable scaling to model. + + Args: + model: model to be scaled + overwrite: whether to overwrite existing scaling factors + submodel_scalers: ComponentMap of Scaler objects to be applied to sub-models + + Returns: + None + + """ + # Scale inlet state + self.call_submodel_scaler_method( + submodel=model.properties_in, + method="variable_scaling_routine", + submodel_scalers=submodel_scalers, + overwrite=overwrite, + ) + + # Propagate scaling from inlet state + self.propagate_state_scaling( + target_state=model.properties_retentate, + source_state=model.properties_in, + overwrite=overwrite, + ) + self.propagate_state_scaling( + target_state=model.properties_permeate, + source_state=model.properties_in, + overwrite=overwrite, + ) + + # Scale remaining states + self.call_submodel_scaler_method( + submodel=model.properties_retentate, + method="variable_scaling_routine", + submodel_scalers=submodel_scalers, + overwrite=overwrite, + ) + self.call_submodel_scaler_method( + submodel=model.properties_permeate, + method="variable_scaling_routine", + submodel_scalers=submodel_scalers, + overwrite=overwrite, + ) + + # Scale remaining variables + # Variables with default scaling factors + for n in self.DEFAULT_SCALING_FACTORS.keys(): + c = model.find_component(n) + for v in c.values(): + self.scale_variable_by_default(v, overwrite=overwrite) + + def constraint_scaling_routine( + self, model, overwrite: bool = False, submodel_scalers: dict = None + ): + """ + Apply constraint scaling to model. + + Args: + model: model to be scaled + overwrite: whether to overwrite existing scaling factors + submodel_scalers: ComponentMap of Scaler objects to be applied to sub-models + + Returns: + None + + """ + # Call scaling methods for sub-models + self.call_submodel_scaler_method( + submodel=model.properties_in, + method="constraint_scaling_routine", + submodel_scalers=submodel_scalers, + overwrite=overwrite, + ) + self.call_submodel_scaler_method( + submodel=model.properties_retentate, + method="constraint_scaling_routine", + submodel_scalers=submodel_scalers, + overwrite=overwrite, + ) + self.call_submodel_scaler_method( + submodel=model.properties_permeate, + method="constraint_scaling_routine", + submodel_scalers=submodel_scalers, + overwrite=overwrite, + ) + + # Scale remaining constraints + for c in model.component_data_objects(Constraint, descend_into=False): + self.scale_constraint_by_nominal_value( + c, + scheme=ConstraintScalingScheme.inverseMaximum, + overwrite=overwrite, + ) + + +@declare_process_block_class("Nanofiltration0D") +class Nanofiltration0DData(UnitModelBlockData): + """ + Nanofiltration model. + """ + + default_initializer = Nanofiltration0DInitializer + default_scaler = Nanofiltration0DScaler + + CONFIG = ConfigDict() + CONFIG.declare( + "dynamic", + ConfigValue( + domain=In([False]), + default=False, + description="Dynamic model flag - must be False", + doc="""Indicates whether this model will be dynamic or not, + **default** = False. Product blocks are always steady-state.""", + ), + ) + CONFIG.declare( + "has_holdup", + ConfigValue( + default=False, + domain=In([False]), + description="Holdup construction flag - must be False", + doc="""Product blocks do not contain holdup, thus this must be + False.""", + ), + ) + CONFIG.declare( + "property_package", + ConfigValue( + default=useDefault, + domain=is_physical_parameter_block, + description="Property package to use for mixer", + doc="""Property parameter object used to define property + calculations, + **default** - useDefault. + **Valid values:** { + **useDefault** - use default package from parent model or flowsheet, + **PropertyParameterObject** - a PropertyParameterBlock object.}""", + ), + ) + CONFIG.declare( + "property_package_args", + ConfigDict( + implicit=True, + description="Arguments to use for constructing property packages", + doc="""A ConfigDict with arguments to be passed to a property + block(s) and used when constructing these, + **default** - None. + **Valid values:** { + see property package for documentation.}""", + ), + ) + CONFIG.declare( + "momentum_balance_type", + ConfigValue( + default=MomentumBalanceType.pressureTotal, + domain=In(MomentumBalanceType), + description="Retentate momentum balance construction flag", + doc="""Indicates what type of momentum balance should be constructed + for the retentate side. Only MomentumBalanceType.none and + MomentumBalanceType.pressureTotal (default) are supported. + **Valid values:** { + **MomentumBalanceType.none** - exclude momentum balances, + **MomentumBalanceType.pressureTotal** - single pressure balance for material}""", + ), + ) + CONFIG.declare( + "has_pressure_change", + ConfigValue( + default=False, + domain=In([True, False]), + description="Retentate pressure change term construction flag", + doc="""Indicates whether terms for retentate pressure change should be + constructed, + **default** - False. + **Valid values:** { + **True** - include pressure change terms, + **False** - exclude pressure change terms.}""", + ), + ) + CONFIG.declare( + "energy_balance_type", + ConfigValue( + default=EnergyBalanceType.isothermal, + domain=In(EnergyBalanceType), + description="Indicated what type of energy balance should be constructed.", + doc="""Indicates what type of momentum balance should be constructed + for the retentate side. Only EnergyBalanceType.none and + EnergyBalanceType.isothermal (default) are supported.""", + ), + ) + # TODO: Support multiple electroneutrality ions + # This needs some rule for prioritizing each ion. + CONFIG.declare( + "electroneutrality_ion", + ConfigValue( + default="Cl_-", + domain=str, + description="Balancing ion to use to maintain electroneutrality", + doc="""Name of ion to use to maintain electroneutrality in + permeate stream. If an ion is provided, the split fraction + for this ion will be solved implicitly to ensure electroneutrality + is maintained in the permeate stream. If None, user must provide + spl;it fractions for all monovalent ions and electroneutrality is + not guaranteed.""", + ), + ) + CONFIG.declare( + "passing_species_list", + ConfigValue( + default=None, + domain=list, + description="List of solutes with low rejection", + doc="""List of solute species which pass through the nanofiltration membrane + with low rejection. If not provided (None), solutes will be classified based + on their charge, with monovalent and uncharged species considered passing.""", + ), + ) + CONFIG.declare( + "default_passing_rejection", + ConfigValue( + default=0.1, + domain=float, + description="Default value for rejection of passing species", + doc="""Default value to be assigned as the rejection fraction for species + identified as passing through the membrane with low rejection.""", + ), + ) + CONFIG.declare( + "default_excluded_rejection", + ConfigValue( + default=1 - 1e-10, + domain=float, + description="Default value for rejection of excluded species", + doc="""Default value to be assigned as the rejection fraction for species + identified as be excluded by membrane (high rejection).""", + ), + ) + + def build(self): + super().build() + + prop_params = self.config.property_package + + # Check that property package only supports single phase + # TODO: Extend the model to handle phase equilibrium + if len(prop_params.phase_list) > 1: + raise ConfigurationError( + "Nanofiltration0D model only supports single phase " + "property packages." + ) + + # Check that the electroneutrality ion is a valid passing species + # I.e. it must be in the passing species list or a monovalent ion + if self.config.electroneutrality_ion is not None: + try: + bal_ion = prop_params.get_component(self.config.electroneutrality_ion) + except AttributeError: + raise ConfigurationError( + f"electroneutrality_ion ({self.config.electroneutrality_ion}) " + "is not a valid component in property package." + ) + + if self.config.passing_species_list is not None: + if ( + self.config.electroneutrality_ion + not in self.config.passing_species_list + ): + raise ConfigurationError( + f"electroneutrality_ion ({self.config.electroneutrality_ion}) " + "must be a member of the passing species list." + ) + else: + # Check that balancing ion is monovalent + if ( + not hasattr(bal_ion.config, "charge") + or abs(bal_ion.config.charge) > 1 + ): + raise ConfigurationError( + f"electroneutrality_ion ({self.config.electroneutrality_ion}) " + "must be a monovalent ion." + ) + + # Check energy balance type is supported + if self.config.energy_balance_type not in ( + EnergyBalanceType.none, + EnergyBalanceType.isothermal, + ): + raise ConfigurationError( + f"Nanofiltration0D model only supports isothermal operation or no energy balance " + f"(assigned {self.config.energy_balance_type})" + ) + + # Check that pressure balance arguments are consistent + if self.config.momentum_balance_type not in ( + MomentumBalanceType.none, + MomentumBalanceType.pressureTotal, + ): + raise ConfigurationError( + f"Nanofiltration0D model only supports total pressure balance or no pressure balance " + f"(assigned {self.config.momentum_balance_type})" + ) + elif ( + self.config.has_pressure_change + and self.config.momentum_balance_type == MomentumBalanceType.none + ): + raise ConfigurationError( + "Inconsistent configuration arguments. has_pressure_change=True " + "requires that momentum_balance_type not equal MomentumBalanceType.none." + ) + + tmp_dict_in = dict(**self.config.property_package_args) + tmp_dict_in["has_phase_equilibrium"] = False + tmp_dict_in["defined_state"] = True + + self.properties_in = self.config.property_package.build_state_block( + self.flowsheet().time, doc="Material properties at inlet", **tmp_dict_in + ) + + tmp_dict_out = dict(**self.config.property_package_args) + tmp_dict_out["has_phase_equilibrium"] = False + tmp_dict_out["defined_state"] = False + + self.properties_retentate = self.config.property_package.build_state_block( + self.flowsheet().time, doc="Material properties at inlet", **tmp_dict_out + ) + + self.properties_permeate = self.config.property_package.build_state_block( + self.flowsheet().time, doc="Material properties at inlet", **tmp_dict_out + ) + + self.add_port("inlet", self.properties_in, doc="Inlet Port") + self.add_port("retentate", self.properties_retentate, doc="Retentate Port") + self.add_port("permeate", self.properties_permeate, doc="Permeate Port") + + # NF separation variables + self.recovery_solvent = Var(self.flowsheet().time, initialize=0.8) + + self._solute_set = Set( + initialize=[ + j + for j in prop_params.component_list + if not prop_params.get_component(j).is_solvent() + and j != self.config.electroneutrality_ion + ] + ) + + def _init_rejection(b, _, j): + if b.config.passing_species_list is not None: + if j in b.config.passing_species_list: + return b.config.default_passing_rejection + return b.config.default_excluded_rejection + else: + if j in prop_params.ion_set: + # Check charge + comp = prop_params.get_component(j) + if abs(comp.config.charge) == 1: + # Monovalent - include in set + return b.config.default_passing_rejection + else: + return b.config.default_excluded_rejection + else: + # Non-ionic solute - assume passing species + return b.config.default_passing_rejection + + self.rejection_comp = Var( + self.flowsheet().time, + self._solute_set, + initialize=_init_rejection, + bounds=(1e-10, 1), + doc="Solute rejection fractions", + ) + + units = self.config.property_package.get_metadata().derived_units + if self.config.has_pressure_change: + self.deltaP = Var( + self.flowsheet().time, + initialize=0, + units=units.PRESSURE, + doc="Retentate side pressure change", + ) + + # Material balance + @self.Constraint(self.flowsheet().time, prop_params.component_list) + def material_balances(b, t, j): + pset = b.properties_in[t].phase_list + pcset = b.properties_in[t].phase_component_set + return sum( + b.properties_in[t].get_material_flow_terms(p, j) + for p in pset + if (p, j) in pcset + ) == sum( + b.properties_retentate[t].get_material_flow_terms(p, j) + for p in pset + if (p, j) in pcset + ) + sum( + b.properties_permeate[t].get_material_flow_terms(p, j) + for p in pset + if (p, j) in pcset + ) + + @self.Constraint(self.flowsheet().time, self.properties_in.phase_component_set) + def separation_constraint(b, t, p, j): + comp = prop_params.get_component(j) + + if comp.is_solvent(): + # Permeate flows equal to recovery * inlet flow + return b.recovery_solvent[t] * b.properties_in[ + t + ].get_material_flow_terms(p, j) == b.properties_permeate[ + t + ].get_material_flow_terms( + p, j + ) + elif j == self.config.electroneutrality_ion: + # Rejection of electroneutrality ion will be solved using electroneutrality constraint + return Constraint.Skip + + # else: Rejection = 1 - C_permeate/C_feed + return ( + (1 - b.rejection_comp[t, j]) + * b.properties_in[t].conc_mol_phase_comp[p, j] + ) == b.properties_permeate[t].conc_mol_phase_comp[p, j] + + if self.config.electroneutrality_ion is not None: + # Add electroneutrality constraint for permeate stream + @self.Constraint(self.flowsheet().time) + def permeate_electronegativity(b, t): + perm_state = b.properties_permeate[t] + return 0 == sum( + perm_state.params.get_component(j).config.charge + * perm_state.flow_mol_phase_comp[p, j] + for p, j in perm_state.phase_component_set + if j in perm_state.params.ion_set + ) + + # Retentate pressure balance + if self.config.momentum_balance_type == MomentumBalanceType.pressureTotal: + + @self.Constraint(self.flowsheet().time, doc="Retentate pressure balance") + def retentate_pressure_balance(b, t): + expr = b.properties_retentate[t].pressure + if b.config.has_pressure_change: + expr += -b.deltaP[t] + return b.properties_in[t].pressure == expr + + # Temperature equalities + if self.config.energy_balance_type == EnergyBalanceType.isothermal: + + @self.Constraint( + self.flowsheet().time, doc="Retentate temperature equality" + ) + def retentate_temperature_equality(b, t): + return ( + b.properties_in[t].temperature + == b.properties_retentate[t].temperature + ) + + @self.Constraint(self.flowsheet().time, doc="Permeate temperature equality") + def permeate_temperature_equality(b, t): + return ( + b.properties_in[t].temperature + == b.properties_permeate[t].temperature + ) diff --git a/watertap/unit_models/nanofiltration_ZO.py b/watertap/unit_models/nanofiltration_ZO.py index 1e413df417..e3f78ca19f 100644 --- a/watertap/unit_models/nanofiltration_ZO.py +++ b/watertap/unit_models/nanofiltration_ZO.py @@ -21,6 +21,7 @@ units as pyunits, ) from pyomo.common.config import Bool, ConfigBlock, ConfigValue, In +from pyomo.common.deprecation import deprecated # Import IDAES cores from idaes.core import ( @@ -46,6 +47,12 @@ _log = idaeslog.getLogger(__name__) +@deprecated( + "The NanofiltrationZO model has been deprecated in favor of the Nanofiltration0D model. " + "The Nanofiltation0D model has most of the capabilities of the NanofiltrationZO model " + "with additional support for ensuring electroneutrality in outlets.", + version="1.3.0", +) @declare_process_block_class("NanofiltrationZO") class NanofiltrationData(InitializationMixin, UnitModelBlockData): """ diff --git a/watertap/unit_models/tests/test_nanofiltration_0D.py b/watertap/unit_models/tests/test_nanofiltration_0D.py new file mode 100644 index 0000000000..3393a86612 --- /dev/null +++ b/watertap/unit_models/tests/test_nanofiltration_0D.py @@ -0,0 +1,1254 @@ +################################################################################# +# WaterTAP Copyright (c) 2020-2024, The Regents of the University of California, +# through Lawrence Berkeley National Laboratory, Oak Ridge National Laboratory, +# National Renewable Energy Laboratory, and National Energy Technology +# Laboratory (subject to receipt of any required approvals from the U.S. Dept. +# of Energy). All rights reserved. +# +# Please see the files COPYRIGHT.md and LICENSE.md for full copyright and license +# information, respectively. These files are also available online at the URL +# "https://github.com/watertap-org/watertap/" +################################################################################# + +import re +import pytest + +from pyomo.environ import ( + assert_optimal_termination, + ConcreteModel, + Constraint, + Set, + Suffix, + TransformationFactory, + value, + Var, +) +from pyomo.network import Port + +from idaes.core import EnergyBalanceType, FlowsheetBlock, MomentumBalanceType +from idaes.core.initialization import InitializationStatus +from idaes.core.scaling import set_scaling_factor +from idaes.core.solvers import get_solver +from idaes.core.util.exceptions import ConfigurationError +from idaes.core.util.model_diagnostics import DiagnosticsToolbox +from idaes.core.util.testing import PhysicalParameterTestBlock + +from watertap.property_models.multicomp_aq_sol_prop_pack import ( + MCASParameterBlock, +) +from watertap.unit_models.nanofiltration_0D import ( + Nanofiltration0D, + Nanofiltration0DInitializer, + Nanofiltration0DScaler, +) + + +class TestBuild: + @pytest.mark.unit + def test_config(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Na_+", "Cl_-"], + charge={"Na_+": 1, "Cl_-": -1}, + ) + m.fs.unit = Nanofiltration0D(property_package=m.fs.properties) + + assert len(m.fs.unit.config) == 11 + assert not m.fs.unit.config.dynamic + assert not m.fs.unit.config.has_holdup + assert m.fs.unit.config.property_package is m.fs.properties + assert m.fs.unit.config.property_package_args == {} + assert ( + m.fs.unit.config.momentum_balance_type == MomentumBalanceType.pressureTotal + ) + assert not m.fs.unit.config.has_pressure_change + assert m.fs.unit.config.energy_balance_type == EnergyBalanceType.isothermal + assert m.fs.unit.config.electroneutrality_ion == "Cl_-" + assert m.fs.unit.config.passing_species_list is None + assert m.fs.unit.config.default_passing_rejection == pytest.approx( + 0.1, rel=1e-12 + ) + assert m.fs.unit.config.default_excluded_rejection == pytest.approx( + 1 - 1e-10, rel=1e-12 + ) + + @pytest.mark.unit + def test_default_initializer_and_scaler(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Na_+", "Cl_-"], + charge={"Na_+": 1, "Cl_-": -1}, + ) + m.fs.unit = Nanofiltration0D(property_package=m.fs.properties) + + assert m.fs.unit.default_initializer is Nanofiltration0DInitializer + assert m.fs.unit.default_scaler is Nanofiltration0DScaler + + @pytest.mark.unit + def test_multiphase(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = PhysicalParameterTestBlock() + + with pytest.raises( + ConfigurationError, + match="Nanofiltration0D model only supports single phase " + "property packages.", + ): + m.fs.unit = Nanofiltration0D(property_package=m.fs.properties) + + @pytest.mark.unit + def test_electroneutrality_invalid_component(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Na_+", "Cl_-"], + charge={"Na_+": 1, "Cl_-": -1}, + ) + + with pytest.raises( + ConfigurationError, + match=re.escape( + "electroneutrality_ion (foo) " + "is not a valid component in property package." + ), + ): + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + electroneutrality_ion="foo", + ) + + @pytest.mark.unit + def test_electroneutrality_not_in_passing_list(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Na_+", "Cl_-"], + charge={"Na_+": 1, "Cl_-": -1}, + ) + + with pytest.raises( + ConfigurationError, + match=re.escape( + "electroneutrality_ion (Na_+) " + "must be a member of the passing species list." + ), + ): + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + electroneutrality_ion="Na_+", + passing_species_list=["Cl_-"], + ) + + @pytest.mark.unit + def test_electroneutrality_divalent(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Na_+", "Cl_-"], + charge={"Na_+": 1, "Cl_-": -2}, + ) + + with pytest.raises( + ConfigurationError, + match=re.escape("electroneutrality_ion (Cl_-) must be a monovalent ion."), + ): + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + ) + + unsupported_balance_types = [ + MomentumBalanceType.pressurePhase, + MomentumBalanceType.momentumTotal, + MomentumBalanceType.momentumPhase, + ] + + @pytest.mark.unit + @pytest.mark.parametrize("balance_type", unsupported_balance_types) + def test_unsupported_momentum_balance(self, balance_type): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Na_+", "Cl_-"], + charge={"Na_+": 1, "Cl_-": -1}, + ) + + with pytest.raises( + ConfigurationError, + match=re.escape( + f"Nanofiltration0D model only supports total pressure balance or no pressure balance " + f"(assigned {balance_type})" + ), + ): + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + momentum_balance_type=balance_type, + ) + + @pytest.mark.unit + def test_invalid_pressure_change(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Na_+", "Cl_-"], + charge={"Na_+": 1, "Cl_-": -1}, + ) + + with pytest.raises( + ConfigurationError, + match="Inconsistent configuration arguments. has_pressure_change=True " + "requires that momentum_balance_type not equal MomentumBalanceType.none.", + ): + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + momentum_balance_type=MomentumBalanceType.none, + has_pressure_change=True, + ) + + unsupported_ebalance_types = [ + e + for e in EnergyBalanceType + if e not in (EnergyBalanceType.none, EnergyBalanceType.isothermal) + ] + + @pytest.mark.unit + @pytest.mark.parametrize("balance_type", unsupported_ebalance_types) + def test_unsupported_energy_balance(self, balance_type): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Na_+", "Cl_-"], + charge={"Na_+": 1, "Cl_-": -1}, + ) + + with pytest.raises( + ConfigurationError, + match=re.escape( + f"Nanofiltration0D model only supports isothermal operation or no energy balance " + f"(assigned {balance_type})" + ), + ): + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + energy_balance_type=balance_type, + ) + + @pytest.mark.unit + def test_default_build(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Ca_2+", "SO4_2-", "Mg_2+", "Na_+", "Cl_-"], + charge={"Ca_2+": 2, "Mg_2+": 2, "SO4_2-": -2, "Na_+": 1, "Cl_-": -1}, + ) + m.fs.unit = Nanofiltration0D(property_package=m.fs.properties) + + assert hasattr(m.fs.unit, "properties_in") + assert hasattr(m.fs.unit, "properties_retentate") + assert hasattr(m.fs.unit, "properties_permeate") + + assert isinstance(m.fs.unit.inlet, Port) + assert isinstance(m.fs.unit.retentate, Port) + assert isinstance(m.fs.unit.permeate, Port) + + assert isinstance(m.fs.unit._solute_set, Set) + assert len(m.fs.unit._solute_set) == 4 + for j in m.fs.unit._solute_set: + assert j in ["Na_+", "Ca_2+", "SO4_2-", "Mg_2+"] + + assert isinstance(m.fs.unit.recovery_solvent, Var) + assert value(m.fs.unit.recovery_solvent[0]) == pytest.approx(0.8, rel=1e-8) + + assert isinstance(m.fs.unit.rejection_comp, Var) + for (_, j), v in m.fs.unit.rejection_comp.items(): + if j in ["Ca_2+", "SO4_2-", "Mg_2+"]: + assert value(v) == pytest.approx( + m.fs.unit.config.default_excluded_rejection, abs=1e-12 + ) + else: + assert value(v) == pytest.approx( + m.fs.unit.config.default_passing_rejection, abs=1e-12 + ) + + assert not hasattr(m.fs.unit, "deltaP") + + assert isinstance(m.fs.unit.material_balances, Constraint) + assert len(m.fs.unit.material_balances) == 6 # 1 time point, 6 components + for (t, j), cd in m.fs.unit.material_balances.items(): + assert str(cd.expr) == str( + m.fs.unit.properties_in[t].flow_mol_phase_comp["Liq", j] + == m.fs.unit.properties_retentate[t].flow_mol_phase_comp["Liq", j] + + m.fs.unit.properties_permeate[t].flow_mol_phase_comp["Liq", j] + ) + + assert isinstance(m.fs.unit.separation_constraint, Constraint) + assert ( + len(m.fs.unit.separation_constraint) == 5 + ) # 1 time point, 5 components (Cl- is skipped) + for (t, p, j), cd in m.fs.unit.separation_constraint.items(): + assert j != "Cl_-" + + if j == "H2O": + assert str(cd.expr) == str( + m.fs.unit.recovery_solvent[t] + * m.fs.unit.properties_in[t].flow_mol_phase_comp[p, j] + == m.fs.unit.properties_permeate[t].flow_mol_phase_comp[p, j] + ) + else: + assert str(cd.expr) == str( + ( + (1 - m.fs.unit.rejection_comp[t, j]) + * m.fs.unit.properties_in[t].conc_mol_phase_comp[p, j] + ) + == m.fs.unit.properties_permeate[t].conc_mol_phase_comp[p, j] + ) + + assert isinstance(m.fs.unit.permeate_electronegativity, Constraint) + assert len(m.fs.unit.permeate_electronegativity) == 1 + assert str(m.fs.unit.permeate_electronegativity[0].expr) == str( + 0 + == sum( + m.fs.properties.get_component(j).config.charge + * m.fs.unit.properties_permeate[0].flow_mol_phase_comp["Liq", j] + for j in m.fs.properties.ion_set + ) + ) + + assert isinstance(m.fs.unit.retentate_pressure_balance, Constraint) + assert len(m.fs.unit.retentate_pressure_balance) == 1 + assert str(m.fs.unit.retentate_pressure_balance[0].expr) == str( + m.fs.unit.properties_in[0].pressure + == m.fs.unit.properties_retentate[0].pressure + ) + + assert isinstance(m.fs.unit.retentate_temperature_equality, Constraint) + assert len(m.fs.unit.retentate_temperature_equality) == 1 + assert str(m.fs.unit.retentate_temperature_equality[0].expr) == str( + m.fs.unit.properties_in[0].temperature + == m.fs.unit.properties_retentate[0].temperature + ) + + assert isinstance(m.fs.unit.permeate_temperature_equality, Constraint) + assert len(m.fs.unit.permeate_temperature_equality) == 1 + assert str(m.fs.unit.permeate_temperature_equality[0].expr) == str( + m.fs.unit.properties_in[0].temperature + == m.fs.unit.properties_permeate[0].temperature + ) + + @pytest.mark.unit + def test_w_passing_listd(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Ca_2+", "SO4_2-", "Mg_2+", "Na_+", "Cl_-"], + charge={"Ca_2+": 2, "Mg_2+": 2, "SO4_2-": -2, "Na_+": 1, "Cl_-": -1}, + ) + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + passing_species_list=["Mg_2+", "Na_+", "Cl_-"], + ) + + assert hasattr(m.fs.unit, "properties_in") + assert hasattr(m.fs.unit, "properties_retentate") + assert hasattr(m.fs.unit, "properties_permeate") + + assert isinstance(m.fs.unit.inlet, Port) + assert isinstance(m.fs.unit.retentate, Port) + assert isinstance(m.fs.unit.permeate, Port) + + assert isinstance(m.fs.unit._solute_set, Set) + assert len(m.fs.unit._solute_set) == 4 + for j in m.fs.unit._solute_set: + assert j in ["Na_+", "Ca_2+", "SO4_2-", "Mg_2+"] + + assert isinstance(m.fs.unit.recovery_solvent, Var) + assert value(m.fs.unit.recovery_solvent[0]) == pytest.approx(0.8, rel=1e-8) + + assert isinstance(m.fs.unit.rejection_comp, Var) + for (_, j), v in m.fs.unit.rejection_comp.items(): + if j in ["Ca_2+", "SO4_2-"]: + assert value(v) == pytest.approx( + m.fs.unit.config.default_excluded_rejection, abs=1e-12 + ) + else: + assert value(v) == pytest.approx( + m.fs.unit.config.default_passing_rejection, abs=1e-12 + ) + + assert not hasattr(m.fs.unit, "deltaP") + + assert isinstance(m.fs.unit.material_balances, Constraint) + assert len(m.fs.unit.material_balances) == 6 # 1 time point, 6 components + for (t, j), cd in m.fs.unit.material_balances.items(): + assert str(cd.expr) == str( + m.fs.unit.properties_in[t].flow_mol_phase_comp["Liq", j] + == m.fs.unit.properties_retentate[t].flow_mol_phase_comp["Liq", j] + + m.fs.unit.properties_permeate[t].flow_mol_phase_comp["Liq", j] + ) + + assert isinstance(m.fs.unit.separation_constraint, Constraint) + assert ( + len(m.fs.unit.separation_constraint) == 5 + ) # 1 time point, 5 components (Cl- is skipped) + for (t, p, j), cd in m.fs.unit.separation_constraint.items(): + assert j != "Cl_-" + + if j == "H2O": + assert str(cd.expr) == str( + m.fs.unit.recovery_solvent[t] + * m.fs.unit.properties_in[t].flow_mol_phase_comp[p, j] + == m.fs.unit.properties_permeate[t].flow_mol_phase_comp[p, j] + ) + else: + assert str(cd.expr) == str( + ( + (1 - m.fs.unit.rejection_comp[t, j]) + * m.fs.unit.properties_in[t].conc_mol_phase_comp[p, j] + ) + == m.fs.unit.properties_permeate[t].conc_mol_phase_comp[p, j] + ) + + assert isinstance(m.fs.unit.permeate_electronegativity, Constraint) + assert len(m.fs.unit.permeate_electronegativity) == 1 + assert str(m.fs.unit.permeate_electronegativity[0].expr) == str( + 0 + == sum( + m.fs.properties.get_component(j).config.charge + * m.fs.unit.properties_permeate[0].flow_mol_phase_comp["Liq", j] + for j in m.fs.properties.ion_set + ) + ) + + assert isinstance(m.fs.unit.retentate_pressure_balance, Constraint) + assert len(m.fs.unit.retentate_pressure_balance) == 1 + assert str(m.fs.unit.retentate_pressure_balance[0].expr) == str( + m.fs.unit.properties_in[0].pressure + == m.fs.unit.properties_retentate[0].pressure + ) + + assert isinstance(m.fs.unit.retentate_temperature_equality, Constraint) + assert len(m.fs.unit.retentate_temperature_equality) == 1 + assert str(m.fs.unit.retentate_temperature_equality[0].expr) == str( + m.fs.unit.properties_in[0].temperature + == m.fs.unit.properties_retentate[0].temperature + ) + + assert isinstance(m.fs.unit.permeate_temperature_equality, Constraint) + assert len(m.fs.unit.permeate_temperature_equality) == 1 + assert str(m.fs.unit.permeate_temperature_equality[0].expr) == str( + m.fs.unit.properties_in[0].temperature + == m.fs.unit.properties_permeate[0].temperature + ) + + @pytest.mark.unit + def test_build_w_pressure_drop(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Ca_2+", "SO4_2-", "Mg_2+", "Na_+", "Cl_-"], + charge={"Ca_2+": 2, "Mg_2+": 2, "SO4_2-": -2, "Na_+": 1, "Cl_-": -1}, + ) + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + has_pressure_change=True, + ) + + assert hasattr(m.fs.unit, "properties_in") + assert hasattr(m.fs.unit, "properties_retentate") + assert hasattr(m.fs.unit, "properties_permeate") + + assert isinstance(m.fs.unit.inlet, Port) + assert isinstance(m.fs.unit.retentate, Port) + assert isinstance(m.fs.unit.permeate, Port) + + assert isinstance(m.fs.unit._solute_set, Set) + assert len(m.fs.unit._solute_set) == 4 + for j in m.fs.unit._solute_set: + assert j in ["Na_+", "Ca_2+", "SO4_2-", "Mg_2+"] + + assert isinstance(m.fs.unit.recovery_solvent, Var) + assert value(m.fs.unit.recovery_solvent[0]) == pytest.approx(0.8, rel=1e-8) + + assert isinstance(m.fs.unit.rejection_comp, Var) + for (_, j), v in m.fs.unit.rejection_comp.items(): + if j in ["Ca_2+", "SO4_2-", "Mg_2+"]: + assert value(v) == pytest.approx( + m.fs.unit.config.default_excluded_rejection, abs=1e-12 + ) + else: + assert value(v) == pytest.approx( + m.fs.unit.config.default_passing_rejection, abs=1e-12 + ) + + assert isinstance(m.fs.unit.deltaP, Var) + + assert isinstance(m.fs.unit.material_balances, Constraint) + assert len(m.fs.unit.material_balances) == 6 # 1 time point, 6 components + for (t, j), cd in m.fs.unit.material_balances.items(): + assert str(cd.expr) == str( + m.fs.unit.properties_in[t].flow_mol_phase_comp["Liq", j] + == m.fs.unit.properties_retentate[t].flow_mol_phase_comp["Liq", j] + + m.fs.unit.properties_permeate[t].flow_mol_phase_comp["Liq", j] + ) + + assert isinstance(m.fs.unit.separation_constraint, Constraint) + assert ( + len(m.fs.unit.separation_constraint) == 5 + ) # 1 time point, 5 components (Cl- is skipped) + for (t, p, j), cd in m.fs.unit.separation_constraint.items(): + assert j != "Cl_-" + + if j == "H2O": + assert str(cd.expr) == str( + m.fs.unit.recovery_solvent[t] + * m.fs.unit.properties_in[t].flow_mol_phase_comp[p, j] + == m.fs.unit.properties_permeate[t].flow_mol_phase_comp[p, j] + ) + else: + assert str(cd.expr) == str( + ( + (1 - m.fs.unit.rejection_comp[t, j]) + * m.fs.unit.properties_in[t].conc_mol_phase_comp[p, j] + ) + == m.fs.unit.properties_permeate[t].conc_mol_phase_comp[p, j] + ) + + assert isinstance(m.fs.unit.permeate_electronegativity, Constraint) + assert len(m.fs.unit.permeate_electronegativity) == 1 + assert str(m.fs.unit.permeate_electronegativity[0].expr) == str( + 0 + == sum( + m.fs.properties.get_component(j).config.charge + * m.fs.unit.properties_permeate[0].flow_mol_phase_comp["Liq", j] + for j in m.fs.properties.ion_set + ) + ) + + assert isinstance(m.fs.unit.retentate_pressure_balance, Constraint) + assert len(m.fs.unit.retentate_pressure_balance) == 1 + assert str(m.fs.unit.retentate_pressure_balance[0].expr) == str( + m.fs.unit.properties_in[0].pressure + == m.fs.unit.properties_retentate[0].pressure - m.fs.unit.deltaP[0] + ) + + assert isinstance(m.fs.unit.retentate_temperature_equality, Constraint) + assert len(m.fs.unit.retentate_temperature_equality) == 1 + assert str(m.fs.unit.retentate_temperature_equality[0].expr) == str( + m.fs.unit.properties_in[0].temperature + == m.fs.unit.properties_retentate[0].temperature + ) + + assert isinstance(m.fs.unit.permeate_temperature_equality, Constraint) + assert len(m.fs.unit.permeate_temperature_equality) == 1 + assert str(m.fs.unit.permeate_temperature_equality[0].expr) == str( + m.fs.unit.properties_in[0].temperature + == m.fs.unit.properties_permeate[0].temperature + ) + + @pytest.mark.unit + def test_build_no_pt_electroneutrality(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Ca_2+", "SO4_2-", "Mg_2+", "Na_+", "Cl_-"], + charge={"Ca_2+": 2, "Mg_2+": 2, "SO4_2-": -2, "Na_+": 1, "Cl_-": -1}, + ) + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + electroneutrality_ion=None, + momentum_balance_type=MomentumBalanceType.none, + energy_balance_type=EnergyBalanceType.none, + ) + + assert hasattr(m.fs.unit, "properties_in") + assert hasattr(m.fs.unit, "properties_retentate") + assert hasattr(m.fs.unit, "properties_permeate") + + assert isinstance(m.fs.unit.inlet, Port) + assert isinstance(m.fs.unit.retentate, Port) + assert isinstance(m.fs.unit.permeate, Port) + + assert isinstance(m.fs.unit._solute_set, Set) + assert len(m.fs.unit._solute_set) == 5 + for j in m.fs.unit._solute_set: + assert j in ["Na_+", "Ca_2+", "SO4_2-", "Mg_2+", "Cl_-"] + + assert isinstance(m.fs.unit.recovery_solvent, Var) + assert value(m.fs.unit.recovery_solvent[0]) == pytest.approx(0.8, rel=1e-8) + + assert isinstance(m.fs.unit.rejection_comp, Var) + for (_, j), v in m.fs.unit.rejection_comp.items(): + if j in ["Ca_2+", "SO4_2-", "Mg_2+"]: + assert value(v) == pytest.approx( + m.fs.unit.config.default_excluded_rejection, abs=1e-12 + ) + else: + assert value(v) == pytest.approx( + m.fs.unit.config.default_passing_rejection, abs=1e-12 + ) + + assert not hasattr(m.fs.unit, "deltaP") + + assert isinstance(m.fs.unit.material_balances, Constraint) + assert len(m.fs.unit.material_balances) == 6 # 1 time point, 6 components + for (t, j), cd in m.fs.unit.material_balances.items(): + assert str(cd.expr) == str( + m.fs.unit.properties_in[t].flow_mol_phase_comp["Liq", j] + == m.fs.unit.properties_retentate[t].flow_mol_phase_comp["Liq", j] + + m.fs.unit.properties_permeate[t].flow_mol_phase_comp["Liq", j] + ) + + assert isinstance(m.fs.unit.separation_constraint, Constraint) + assert len(m.fs.unit.separation_constraint) == 6 # 1 time point, 6 components + for (t, p, j), cd in m.fs.unit.separation_constraint.items(): + if j == "H2O": + assert str(cd.expr) == str( + m.fs.unit.recovery_solvent[t] + * m.fs.unit.properties_in[t].flow_mol_phase_comp[p, j] + == m.fs.unit.properties_permeate[t].flow_mol_phase_comp[p, j] + ) + else: + assert str(cd.expr) == str( + ( + (1 - m.fs.unit.rejection_comp[t, j]) + * m.fs.unit.properties_in[t].conc_mol_phase_comp[p, j] + ) + == m.fs.unit.properties_permeate[t].conc_mol_phase_comp[p, j] + ) + + assert not hasattr(m.fs.unit, "permeate_electronegativity") + assert not hasattr(m.fs.unit, "retentate_pressure_balance") + assert not hasattr(m.fs.unit, "retentate_temperature_equality") + assert not hasattr(m.fs.unit, "permeate_temperature_equality") + + +@pytest.fixture +def mcas_model(): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Ca_2+", "SO4_2-", "Mg_2+", "Na_+", "Cl_-"], + diffusivity_data={ + ("Liq", "Ca_2+"): 9.2e-10, + ("Liq", "SO4_2-"): 1.06e-09, + ("Liq", "Mg_2+"): 7.06e-10, + ("Liq", "Na_+"): 1.33e-09, + ("Liq", "Cl_-"): 2.03e-09, + }, + mw_data={ + "H2O": 0.018, + "Ca_2+": 0.04, + "Mg_2+": 0.024, + "SO4_2-": 0.096, + "Na_+": 0.023, + "Cl_-": 0.035, + }, + stokes_radius_data={ + "Ca_2+": 3.09e-10, + "Mg_2+": 3.47e-10, + "SO4_2-": 2.3e-10, + "Cl_-": 1.21e-10, + "Na_+": 1.84e-10, + }, + charge={"Ca_2+": 2, "Mg_2+": 2, "SO4_2-": -2, "Na_+": 1, "Cl_-": -1}, + ) + + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + has_pressure_change=True, + ) + + return m + + +class TestNanofiltration0DInitializer: + @pytest.mark.component + def test_initialize(self, mcas_model): + # Fix DoF + mcas_model.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "H2O"].fix(53.6036) + mcas_model.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Ca_2+"].fix(0.00955) + mcas_model.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Mg_2+"].fix(0.5808) + mcas_model.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "SO4_2-"].fix(0.02225) + mcas_model.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Cl_-"].fix(1.61977) + mcas_model.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Na_+"].fix(0.48357) + mcas_model.fs.unit.inlet.temperature[0].fix(298.15) + mcas_model.fs.unit.inlet.pressure[0].fix(4e5) + + mcas_model.fs.unit.recovery_solvent.fix() + mcas_model.fs.unit.rejection_comp.fix() + mcas_model.fs.unit.permeate.pressure[0].fix(1e5) + mcas_model.fs.unit.deltaP.fix(-1e4) + + # Initialize + initializer = Nanofiltration0DInitializer() + + initializer.initialize(mcas_model.fs.unit) + + assert ( + initializer.summary[mcas_model.fs.unit]["status"] == InitializationStatus.Ok + ) + + for c in mcas_model.component_data_objects(Constraint, descend_into=True): + assert ( + abs(value(c.body - c.lower)) <= 1e-5 + and abs(value(c.upper - c.body)) <= 1e-5 + ) + + @pytest.mark.unit + def test_init_conc(self, mcas_model): + # Use MCAS model to dummy calls to initializer sub methods + initializer = Nanofiltration0DInitializer() + + mcas_model.in_var = Var(mcas_model.fs.properties.component_list, initialize=10) + mcas_model.ret_var = Var(mcas_model.fs.properties.component_list, initialize=0) + mcas_model.perm_var = Var(mcas_model.fs.properties.component_list, initialize=0) + + initializer._init_conc( + mcas_model.fs.unit, + 0, # time index + mcas_model.fs.unit.properties_in[0], + mcas_model.in_var, + mcas_model.ret_var, + mcas_model.perm_var, + ) + + assert value(mcas_model.ret_var["H2O"]) == pytest.approx( + 10, rel=1e-10 + ) # solvent concentration unchanged + assert value(mcas_model.ret_var["Cl_-"]) == pytest.approx( + 2, rel=1e-10 + ) # should use solvent recovery + assert value(mcas_model.ret_var["Na_+"]) == pytest.approx(1, rel=1e-10) + assert value(mcas_model.ret_var["Ca_2+"]) == pytest.approx( + 10 * (1 - 1e-10), rel=1e-10 + ) + assert value(mcas_model.ret_var["Mg_2+"]) == pytest.approx( + 10 * (1 - 1e-10), rel=1e-10 + ) + assert value(mcas_model.ret_var["SO4_2-"]) == pytest.approx( + 10 * (1 - 1e-10), rel=1e-10 + ) + + assert value(mcas_model.perm_var["H2O"]) == pytest.approx( + 10, rel=1e-10 + ) # solvent concentration unchanged + assert value(mcas_model.perm_var["Cl_-"]) == pytest.approx( + 8, rel=1e-10 + ) # should use solvent recovery + assert value(mcas_model.perm_var["Na_+"]) == pytest.approx(9, rel=1e-10) + assert value(mcas_model.perm_var["Ca_2+"]) == pytest.approx( + 10 * 1e-10, rel=1e-10 + ) + assert value(mcas_model.perm_var["Mg_2+"]) == pytest.approx( + 10 * 1e-10, rel=1e-10 + ) + assert value(mcas_model.perm_var["SO4_2-"]) == pytest.approx( + 10 * 1e-10, rel=1e-10 + ) + + @pytest.mark.unit + def test_init_flow_comp(self, mcas_model): + # Use MCA model to dummy calls to initializer sub methods + initializer = Nanofiltration0DInitializer() + + mcas_model.in_var_comp = Var( + mcas_model.fs.properties.component_list, initialize=10 + ) + mcas_model.ret_var = Var(mcas_model.fs.properties.component_list, initialize=0) + mcas_model.perm_var = Var(mcas_model.fs.properties.component_list, initialize=0) + + initializer._init_flow( + mcas_model.fs.unit, + 0, # time index + mcas_model.fs.unit.properties_in[0], + mcas_model.in_var_comp, + mcas_model.ret_var, + mcas_model.perm_var, + ) + + assert value(mcas_model.ret_var["H2O"]) == pytest.approx(2, rel=1e-10) + assert value(mcas_model.ret_var["Cl_-"]) == pytest.approx( + 2, rel=1e-10 + ) # should use solvent recovery + assert value(mcas_model.ret_var["Na_+"]) == pytest.approx(1, rel=1e-10) + assert value(mcas_model.ret_var["Ca_2+"]) == pytest.approx( + 10 * (1 - 1e-10), rel=1e-10 + ) + assert value(mcas_model.ret_var["Mg_2+"]) == pytest.approx( + 10 * (1 - 1e-10), rel=1e-10 + ) + assert value(mcas_model.ret_var["SO4_2-"]) == pytest.approx( + 10 * (1 - 1e-10), rel=1e-10 + ) + + assert value(mcas_model.perm_var["H2O"]) == pytest.approx(8, rel=1e-10) + assert value(mcas_model.perm_var["Cl_-"]) == pytest.approx( + 8, rel=1e-10 + ) # should use solvent recovery + assert value(mcas_model.perm_var["Na_+"]) == pytest.approx(9, rel=1e-10) + assert value(mcas_model.perm_var["Ca_2+"]) == pytest.approx( + 10 * 1e-10, rel=1e-10 + ) + assert value(mcas_model.perm_var["Mg_2+"]) == pytest.approx( + 10 * 1e-10, rel=1e-10 + ) + assert value(mcas_model.perm_var["SO4_2-"]) == pytest.approx( + 10 * 1e-10, rel=1e-10 + ) + + @pytest.mark.unit + def test_init_flow_vol(self, mcas_model): + # Use MCA model to dummy calls to initializer sub methods + initializer = Nanofiltration0DInitializer() + + mcas_model.in_var_vol = Var(initialize=10) + mcas_model.ret_var = Var(initialize=0) + mcas_model.perm_var = Var(initialize=0) + + initializer._init_flow( + mcas_model.fs.unit, + 0, # time index + mcas_model.fs.unit.properties_in[0], + mcas_model.in_var_vol, + mcas_model.ret_var, + mcas_model.perm_var, + ) + + assert value(mcas_model.ret_var) == pytest.approx(2, rel=1e-10) + assert value(mcas_model.perm_var) == pytest.approx(8, rel=1e-10) + + @pytest.mark.unit + def test_init_frac(self, mcas_model): + # Use MCA model to dummy calls to initializer sub methods + initializer = Nanofiltration0DInitializer() + + mcas_model.in_var = Var(mcas_model.fs.properties.component_list, initialize=10) + mcas_model.ret_var = Var(mcas_model.fs.properties.component_list, initialize=0) + mcas_model.perm_var = Var(mcas_model.fs.properties.component_list, initialize=0) + + initializer._init_frac( + mcas_model.fs.unit, + 0, # time index + mcas_model.fs.unit.properties_in[0], + mcas_model.in_var, + mcas_model.ret_var, + mcas_model.perm_var, + ) + + assert value(mcas_model.ret_var["H2O"]) == pytest.approx(2 / 35, rel=1e-6) + assert value(mcas_model.ret_var["Cl_-"]) == pytest.approx(2 / 35, rel=1e-6) + assert value(mcas_model.ret_var["Na_+"]) == pytest.approx(1 / 35, rel=1e-6) + assert value(mcas_model.ret_var["Ca_2+"]) == pytest.approx(10 / 35, rel=1e-6) + assert value(mcas_model.ret_var["Mg_2+"]) == pytest.approx(10 / 35, rel=1e-6) + assert value(mcas_model.ret_var["SO4_2-"]) == pytest.approx(10 / 35, rel=1e-6) + + assert value(mcas_model.perm_var["H2O"]) == pytest.approx(8 / 25, rel=1e-6) + assert value(mcas_model.perm_var["Cl_-"]) == pytest.approx(8 / 25, rel=1e-6) + assert value(mcas_model.perm_var["Na_+"]) == pytest.approx(9 / 25, rel=1e-6) + assert value(mcas_model.perm_var["Ca_2+"]) == pytest.approx(0, abs=1e-6) + assert value(mcas_model.perm_var["Mg_2+"]) == pytest.approx(0, abs=1e-6) + assert value(mcas_model.perm_var["SO4_2-"]) == pytest.approx(0, abs=1e-6) + + +class TestNanofiltration0DScaler: + @pytest.mark.unit + def test_variable_scaling_routine(self, mcas_model): + scaler = Nanofiltration0DScaler() + + scaler.variable_scaling_routine(mcas_model.fs.unit) + + # Check that submodels have scaling factors + assert isinstance(mcas_model.fs.unit.properties_in[0].scaling_factor, Suffix) + assert len(mcas_model.fs.unit.properties_in[0].scaling_factor) > 1 + + assert isinstance( + mcas_model.fs.unit.properties_retentate[0].scaling_factor, Suffix + ) + assert len(mcas_model.fs.unit.properties_retentate[0].scaling_factor) > 1 + + assert isinstance( + mcas_model.fs.unit.properties_permeate[0].scaling_factor, Suffix + ) + assert len(mcas_model.fs.unit.properties_permeate[0].scaling_factor) > 1 + + # Check scaling factors for unit variables + assert isinstance(mcas_model.fs.unit.scaling_factor, Suffix) + assert len(mcas_model.fs.unit.scaling_factor) == 6 + + assert mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.deltaP[0] + ] == pytest.approx(1e-4, rel=1e-8) + + for j in ["Ca_2+", "Mg_2+", "SO4_2-", "Na_+"]: + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.rejection_comp[0, j] + ] + == 1e2 + ) + assert ( + mcas_model.fs.unit.scaling_factor[mcas_model.fs.unit.recovery_solvent[0]] + == 10 + ) + + @pytest.mark.unit + def test_constraint_scaling_routine(self, mcas_model): + scaler = Nanofiltration0DScaler() + + scaler.constraint_scaling_routine(mcas_model.fs.unit) + + # Check that submodels have scaling factors + # No constraint scaling factors applied ot sub-models, so no check for length of Suffix + assert isinstance(mcas_model.fs.unit.properties_in[0].scaling_factor, Suffix) + assert isinstance( + mcas_model.fs.unit.properties_retentate[0].scaling_factor, Suffix + ) + assert isinstance( + mcas_model.fs.unit.properties_permeate[0].scaling_factor, Suffix + ) + + # Check scaling factors for unit variables + assert isinstance(mcas_model.fs.unit.scaling_factor, Suffix) + assert len(mcas_model.fs.unit.scaling_factor) == 15 + + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.material_balances[0.0, "Ca_2+"] + ] + == 10.0 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.material_balances[0.0, "Cl_-"] + ] + == 10.0 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.material_balances[0.0, "H2O"] + ] + == 10.0 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.material_balances[0.0, "Mg_2+"] + ] + == 10.0 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.material_balances[0.0, "Na_+"] + ] + == 10.0 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.material_balances[0.0, "SO4_2-"] + ] + == 10.0 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.permeate_electronegativity[0.0] + ] + == 5.0 + ) + assert mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.permeate_temperature_equality[0.0] + ] == pytest.approx(3.354016e-3, rel=1e-5) + + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.separation_constraint[0.0, "Liq", "Ca_2+"] + ] + == 0.002 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.separation_constraint[0.0, "Liq", "H2O"] + ] + == 10.0 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.separation_constraint[0.0, "Liq", "Mg_2+"] + ] + == 0.002 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.separation_constraint[0.0, "Liq", "Na_+"] + ] + == 0.002 + ) + assert ( + mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.separation_constraint[0.0, "Liq", "SO4_2-"] + ] + == 0.002 + ) + assert mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.retentate_pressure_balance[0.0] + ] == pytest.approx(9.869233e-06, rel=1e-5) + assert mcas_model.fs.unit.scaling_factor[ + mcas_model.fs.unit.retentate_temperature_equality[0.0] + ] == pytest.approx(3.354016e-3, rel=1e-5) + + +class TestMCAS: + @pytest.fixture(scope="class") + def mcas_case(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock( + solute_list=["Ca_2+", "SO4_2-", "Mg_2+", "Na_+", "Cl_-"], + diffusivity_data={ + ("Liq", "Ca_2+"): 9.2e-10, + ("Liq", "SO4_2-"): 1.06e-09, + ("Liq", "Mg_2+"): 7.06e-10, + ("Liq", "Na_+"): 1.33e-09, + ("Liq", "Cl_-"): 2.03e-09, + }, + mw_data={ + "H2O": 0.018, + "Ca_2+": 0.04, + "Mg_2+": 0.024, + "SO4_2-": 0.096, + "Na_+": 0.023, + "Cl_-": 0.035, + }, + stokes_radius_data={ + "Ca_2+": 3.09e-10, + "Mg_2+": 3.47e-10, + "SO4_2-": 2.3e-10, + "Cl_-": 1.21e-10, + "Na_+": 1.84e-10, + }, + charge={"Ca_2+": 2, "Mg_2+": 2, "SO4_2-": -2, "Na_+": 1, "Cl_-": -1}, + ) + + m.fs.unit = Nanofiltration0D( + property_package=m.fs.properties, + has_pressure_change=True, + ) + + # Fix other inlet state variables + m.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "H2O"].fix(53.6036) + m.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Ca_2+"].fix(0.00955) + m.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Mg_2+"].fix(0.5808) + m.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "SO4_2-"].fix(0.02225) + m.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Cl_-"].fix(1.61977) + m.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Na_+"].fix(0.48357) + m.fs.unit.inlet.temperature[0].fix(298.15) + m.fs.unit.inlet.pressure[0].fix(4e5) + + m.fs.unit.recovery_solvent.fix() + m.fs.unit.rejection_comp.fix() + m.fs.unit.deltaP.fix(-1e4) + m.fs.unit.permeate.pressure[0].fix(101325) + + return m + + @pytest.mark.component + def test_structural_issues(self, mcas_case): + dt = DiagnosticsToolbox(mcas_case) + + dt.assert_no_structural_warnings(ignore_evaluation_errors=True) + + @pytest.mark.component + def test_scaling(self, mcas_case): + # Scale model + set_scaling_factor( + mcas_case.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "H2O"], 1e-1 + ) + set_scaling_factor( + mcas_case.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Ca_2+"], 1e3 + ) + set_scaling_factor( + mcas_case.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Mg_2+"], 10 + ) + set_scaling_factor( + mcas_case.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "SO4_2-"], 1e2 + ) + set_scaling_factor( + mcas_case.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Cl_-"], 1 + ) + set_scaling_factor( + mcas_case.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", "Na_+"], 10 + ) + + scaler = Nanofiltration0DScaler() + scaler.scale_model(mcas_case.fs.unit) + + assert ( + mcas_case.fs.unit.scaling_factor[mcas_case.fs.unit.recovery_solvent[0]] + == 10 + ) + for j in ["Ca_2+", "Mg_2+", "SO4_2-", "Na_+"]: + assert ( + mcas_case.fs.unit.scaling_factor[mcas_case.fs.unit.rejection_comp[0, j]] + == 100 + ) + + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.deltaP[0.0] + ] == pytest.approx(0.0001, rel=1e-5) + + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.material_balances[0.0, "H2O"] + ] == pytest.approx(0.01866, rel=1e-3) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.material_balances[0.0, "Ca_2+"] + ] == pytest.approx(104.7, rel=1e-3) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.material_balances[0.0, "SO4_2-"] + ] == pytest.approx(44.94, rel=1e-3) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.material_balances[0.0, "Mg_2+"] + ] == pytest.approx(1.722, rel=1e-3) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.material_balances[0.0, "Na_+"] + ] == pytest.approx(2.068, rel=1e-3) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.material_balances[0.0, "Cl_-"] + ] == pytest.approx(0.6174, rel=1e-3) + + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.separation_constraint[0.0, "Liq", "H2O"] + ] == pytest.approx(0.0233193, rel=1e-5) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.separation_constraint[0.0, "Liq", "Ca_2+"] + ] == pytest.approx(0.004, rel=1e-5) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.separation_constraint[0.0, "Liq", "SO4_2-"] + ] == pytest.approx(0.00960, rel=1e-5) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.separation_constraint[0.0, "Liq", "Mg_2+"] + ] == pytest.approx(0.00240, rel=1e-5) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.separation_constraint[0.0, "Liq", "Na_+"] + ] == pytest.approx(0.0023, rel=1e-3) + + assert ( + mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.permeate_electronegativity[0.0] + ] + == 1 + ) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.retentate_pressure_balance[0.0] + ] == pytest.approx(0.0000025, rel=1e-3) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.retentate_temperature_equality[0.0] + ] == pytest.approx(0.003354, rel=1e-3) + assert mcas_case.fs.unit.scaling_factor[ + mcas_case.fs.unit.permeate_temperature_equality[0.0] + ] == pytest.approx(0.003354, rel=1e-3) + + @pytest.mark.component + def test_initialize(self, mcas_case): + initializer = Nanofiltration0DInitializer() + initializer.initialize(mcas_case.fs.unit) + + assert ( + initializer.summary[mcas_case.fs.unit]["status"] == InitializationStatus.Ok + ) + + @pytest.mark.component + def test_solve(self, mcas_case): + solver = get_solver( + "ipopt_v2", writer_config={"scale_model": True, "linear_presolve": True} + ) + results = solver.solve(mcas_case) + + assert_optimal_termination(results) + + @pytest.mark.component + def test_numerical_issues(self, mcas_case): + sm = TransformationFactory("core.scale_model").create_using( + mcas_case, rename=False + ) + + dt = DiagnosticsToolbox(sm.fs.unit) + dt.assert_no_numerical_warnings() + + @pytest.mark.component + def test_solution(self, mcas_case): + # Retentate stream + assert value( + mcas_case.fs.unit.retentate.flow_mol_phase_comp[0, "Liq", "H2O"] + ) == pytest.approx(0.2 * 53.6036, rel=1e-5) + assert value( + mcas_case.fs.unit.retentate.flow_mol_phase_comp[0, "Liq", "Ca_2+"] + ) == pytest.approx(0.00955, rel=1e-5) + assert value( + mcas_case.fs.unit.retentate.flow_mol_phase_comp[0, "Liq", "Mg_2+"] + ) == pytest.approx(0.5808, rel=1e-5) + assert value( + mcas_case.fs.unit.retentate.flow_mol_phase_comp[0, "Liq", "SO4_2-"] + ) == pytest.approx(0.02225, rel=1e-5) + assert value( + mcas_case.fs.unit.retentate.flow_mol_phase_comp[0, "Liq", "Cl_-"] + ) == pytest.approx(1.29167, rel=1e-5) + assert value( + mcas_case.fs.unit.retentate.flow_mol_phase_comp[0, "Liq", "Na_+"] + ) == pytest.approx(0.155472, rel=1e-5) + assert value(mcas_case.fs.unit.retentate.temperature[0]) == pytest.approx( + 298.15, rel=1e-5 + ) + assert value(mcas_case.fs.unit.retentate.pressure[0]) == pytest.approx( + 3.9e5, rel=1e-5 + ) + + # Permeate stream + assert value( + mcas_case.fs.unit.permeate.flow_mol_phase_comp[0, "Liq", "H2O"] + ) == pytest.approx(0.8 * 53.6036, rel=1e-5) + assert value( + mcas_case.fs.unit.permeate.flow_mol_phase_comp[0, "Liq", "Ca_2+"] + ) == pytest.approx(0, abs=1e-7) + assert value( + mcas_case.fs.unit.permeate.flow_mol_phase_comp[0, "Liq", "Mg_2+"] + ) == pytest.approx(0, abs=1e-7) + assert value( + mcas_case.fs.unit.permeate.flow_mol_phase_comp[0, "Liq", "SO4_2-"] + ) == pytest.approx(0, abs=1e-7) + assert value( + mcas_case.fs.unit.permeate.flow_mol_phase_comp[0, "Liq", "Cl_-"] + ) == pytest.approx(0.328098, rel=1e-5) + assert value( + mcas_case.fs.unit.permeate.flow_mol_phase_comp[0, "Liq", "Na_+"] + ) == pytest.approx(0.328098, rel=1e-5) + assert value(mcas_case.fs.unit.permeate.temperature[0]) == pytest.approx( + 298.15, rel=1e-5 + ) + assert value(mcas_case.fs.unit.permeate.pressure[0]) == pytest.approx( + 101325, rel=1e-5 + ) + + # Rejection + for j in ["Ca_2+", "Mg_2+", "SO4_2-"]: + assert value( + mcas_case.fs.unit.properties_permeate[0].conc_mol_phase_comp["Liq", j] + / mcas_case.fs.unit.properties_in[0].conc_mol_phase_comp["Liq", j] + ) == pytest.approx(0, abs=1e-8) + assert value( + mcas_case.fs.unit.properties_permeate[0].conc_mol_phase_comp["Liq", "Na_+"] + / mcas_case.fs.unit.properties_in[0].conc_mol_phase_comp["Liq", "Na_+"] + ) == pytest.approx(0.9, abs=1e-8) + + @pytest.mark.component + def test_conservation(self, mcas_case): + for j in mcas_case.fs.properties.component_list: + assert ( + abs( + value( + mcas_case.fs.unit.inlet.flow_mol_phase_comp[0, "Liq", j] + - mcas_case.fs.unit.retentate.flow_mol_phase_comp[0, "Liq", j] + - mcas_case.fs.unit.permeate.flow_mol_phase_comp[0, "Liq", j] + ) + ) + <= 1e-6 + )