From 2819acdedb5432f377ffdd1eb2a26f4dca4c934d Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Sun, 23 Nov 2025 15:26:46 -0700 Subject: [PATCH 01/15] create unit_models dir --- watertap/tools/unit_models/__init__.py | 0 1 file changed, 0 insertions(+), 0 deletions(-) create mode 100644 watertap/tools/unit_models/__init__.py diff --git a/watertap/tools/unit_models/__init__.py b/watertap/tools/unit_models/__init__.py new file mode 100644 index 0000000000..e69de29bb2 From 6d9e6d403999927e1d056127f57781ee7f7477f7 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Sun, 23 Nov 2025 15:32:06 -0700 Subject: [PATCH 02/15] add calculate_operating_pressure helper function --- watertap/tools/unit_models/__init__.py | 13 ++++ watertap/tools/unit_models/reverse_osmosis.py | 72 +++++++++++++++++++ 2 files changed, 85 insertions(+) create mode 100644 watertap/tools/unit_models/reverse_osmosis.py diff --git a/watertap/tools/unit_models/__init__.py b/watertap/tools/unit_models/__init__.py index e69de29bb2..649a45a002 100644 --- a/watertap/tools/unit_models/__init__.py +++ b/watertap/tools/unit_models/__init__.py @@ -0,0 +1,13 @@ +################################################################################# +# 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/" +################################################################################# + +from .reverse_osmosis import * diff --git a/watertap/tools/unit_models/reverse_osmosis.py b/watertap/tools/unit_models/reverse_osmosis.py new file mode 100644 index 0000000000..b113ba8748 --- /dev/null +++ b/watertap/tools/unit_models/reverse_osmosis.py @@ -0,0 +1,72 @@ +################################################################################# +# 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/" +################################################################################# +""" +This module contains functions to be used with WaterTAP +ReverseOsmosis0D or ReverseOsmosis1D unit models. +""" + +from pyomo.environ import ( + ConcreteModel, + assert_optimal_termination, + value, +) +from idaes.core.util.initialization import solve_indexed_blocks + + +__all__ = ["calculate_operating_pressure"] + + +def calculate_operating_pressure( + feed_state_block=None, + over_pressure=0.15, + water_recovery_mass=0.5, + NaCl_passage=0.01, + solver=None, +): + """ + Estimate operating pressure for RO unit model given the following arguments: + + Arguments: + feed_state_block: the state block of the RO feed that has the non-pressure state + variables initialized to their values (default=None) + over_pressure: the amount of operating pressure above the brine osmotic pressure + represented as a fraction (default=0.15) + water_recovery_mass: the mass-based fraction of inlet H2O that becomes permeate + (default=0.5) + NaCl_passage: the mass-based fraction of inlet NaCl that becomes permeate + (default=0.01) + solver: solver object to be used (default=None) + """ + t = ConcreteModel() # create temporary model + prop = feed_state_block.config.parameters + t.brine = prop.build_state_block([0]) + + # specify state block + t.brine[0].flow_mass_phase_comp["Liq", "H2O"].fix( + value(feed_state_block.flow_mass_phase_comp["Liq", "H2O"]) + * (1 - water_recovery_mass) + ) + t.brine[0].flow_mass_phase_comp["Liq", "NaCl"].fix( + value(feed_state_block.flow_mass_phase_comp["Liq", "NaCl"]) * (1 - NaCl_passage) + ) + t.brine[0].pressure.fix( + 101325 + ) # valid when osmotic pressure is independent of hydraulic pressure + t.brine[0].temperature.fix(value(feed_state_block.temperature)) + # calculate osmotic pressure + # since properties are created on demand, we must touch the property to create it + t.brine[0].pressure_osm_phase + # solve state block + results = solve_indexed_blocks(solver, [t.brine]) + assert_optimal_termination(results) + + return value(t.brine[0].pressure_osm_phase["Liq"]) * (1 + over_pressure) From 969a60bb3cbf5548037f1709225eb947daa47e08 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Sun, 23 Nov 2025 15:32:21 -0700 Subject: [PATCH 03/15] import calculate_operating_pressure --- .../RO_with_energy_recovery.py | 50 +------------------ 1 file changed, 2 insertions(+), 48 deletions(-) diff --git a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py index 3366297b8f..85453332d0 100644 --- a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py +++ b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py @@ -24,7 +24,6 @@ from idaes.core import FlowsheetBlock from watertap.core.solvers import get_solver from idaes.core.util.model_statistics import degrees_of_freedom -from idaes.core.util.initialization import solve_indexed_blocks, propagate_state from idaes.models.unit_models import Mixer, Separator, Product, Feed from idaes.models.unit_models.mixer import MomentumMixingType from idaes.core import UnitModelCostingBlock @@ -42,6 +41,7 @@ from watertap.unit_models.pressure_exchanger import PressureExchanger from watertap.unit_models.pressure_changer import Pump, EnergyRecoveryDevice from watertap.core.util.initialization import assert_degrees_of_freedom +from watertap.tools.unit_models import calculate_operating_pressure from watertap.costing import WaterTAPCosting @@ -255,7 +255,7 @@ def set_operating_conditions( operating_pressure = calculate_operating_pressure( feed_state_block=m.fs.feed.properties[0], over_pressure=over_pressure, - water_recovery=water_recovery, + water_recovery_mass=water_recovery, NaCl_passage=0.01, solver=solver, ) @@ -312,52 +312,6 @@ def set_operating_conditions( ) -def calculate_operating_pressure( - feed_state_block=None, - over_pressure=0.15, - water_recovery=0.5, - NaCl_passage=0.01, - solver=None, -): - """ - Estimate operating pressure for RO unit model given the following arguments: - - Arguments: - feed_state_block: the state block of the RO feed that has the non-pressure state - variables initialized to their values (default=None) - over_pressure: the amount of operating pressure above the brine osmotic pressure - represented as a fraction (default=0.15) - water_recovery: the mass-based fraction of inlet H2O that becomes permeate - (default=0.5) - NaCl_passage: the mass-based fraction of inlet NaCl that becomes permeate - (default=0.01) - solver: solver object to be used (default=None) - """ - t = ConcreteModel() # create temporary model - prop = feed_state_block.config.parameters - t.brine = prop.build_state_block([0]) - - # specify state block - t.brine[0].flow_mass_phase_comp["Liq", "H2O"].fix( - value(feed_state_block.flow_mass_phase_comp["Liq", "H2O"]) - * (1 - water_recovery) - ) - t.brine[0].flow_mass_phase_comp["Liq", "NaCl"].fix( - value(feed_state_block.flow_mass_phase_comp["Liq", "NaCl"]) * (1 - NaCl_passage) - ) - t.brine[0].pressure.fix( - 101325 - ) # valid when osmotic pressure is independent of hydraulic pressure - t.brine[0].temperature.fix(value(feed_state_block.temperature)) - # calculate osmotic pressure - # since properties are created on demand, we must touch the property to create it - t.brine[0].pressure_osm_phase - # solve state block - results = solve_indexed_blocks(solver, [t.brine]) - assert_optimal_termination(results) - - return value(t.brine[0].pressure_osm_phase["Liq"]) * (1 + over_pressure) - def solve(blk, solver=None, tee=False, check_termination=True): if solver is None: From 6794e0a3251041b434e13493c8cfdef1e43162dc Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Sun, 23 Nov 2025 15:33:43 -0700 Subject: [PATCH 04/15] add back propagate_state import --- .../RO_with_energy_recovery/RO_with_energy_recovery.py | 1 + 1 file changed, 1 insertion(+) diff --git a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py index 85453332d0..772c878ea5 100644 --- a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py +++ b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py @@ -23,6 +23,7 @@ from pyomo.network import Arc from idaes.core import FlowsheetBlock from watertap.core.solvers import get_solver +from idaes.core.util.initialization import propagate_state from idaes.core.util.model_statistics import degrees_of_freedom from idaes.models.unit_models import Mixer, Separator, Product, Feed from idaes.models.unit_models.mixer import MomentumMixingType From 57e66f81e0711273d8f02131cca9b59e5f69a538 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Sun, 23 Nov 2025 15:51:18 -0700 Subject: [PATCH 05/15] initial test --- watertap/tools/unit_models/reverse_osmosis.py | 45 +++++++++++++------ 1 file changed, 31 insertions(+), 14 deletions(-) diff --git a/watertap/tools/unit_models/reverse_osmosis.py b/watertap/tools/unit_models/reverse_osmosis.py index b113ba8748..cc0dd6d5bf 100644 --- a/watertap/tools/unit_models/reverse_osmosis.py +++ b/watertap/tools/unit_models/reverse_osmosis.py @@ -21,6 +21,9 @@ ) from idaes.core.util.initialization import solve_indexed_blocks +from watertap.property_models.seawater_prop_pack import SeawaterStateBlockData +from watertap.property_models.NaCl_prop_pack import NaClStateBlockData + __all__ = ["calculate_operating_pressure"] @@ -29,7 +32,7 @@ def calculate_operating_pressure( feed_state_block=None, over_pressure=0.15, water_recovery_mass=0.5, - NaCl_passage=0.01, + salt_passage=0.01, solver=None, ): """ @@ -42,31 +45,45 @@ def calculate_operating_pressure( represented as a fraction (default=0.15) water_recovery_mass: the mass-based fraction of inlet H2O that becomes permeate (default=0.5) - NaCl_passage: the mass-based fraction of inlet NaCl that becomes permeate + salt_passage: the mass-based fraction of inlet salt that becomes permeate (default=0.01) solver: solver object to be used (default=None) """ - t = ConcreteModel() # create temporary model + if not any( + isinstance(feed_state_block, cls) + for cls in [SeawaterStateBlockData, NaClStateBlockData] + ): + raise TypeError( + "feed_state_block argument must be a SeawaterStateBlockData or NaClStateBlockData object" + ) + + if isinstance(feed_state_block, NaClStateBlockData): + comp = "NaCl" + if isinstance(feed_state_block, SeawaterStateBlockData): + comp = "TDS" + + tmp = ConcreteModel() # create temporary model prop = feed_state_block.config.parameters - t.brine = prop.build_state_block([0]) + tmp.brine = prop.build_state_block([0]) + tmp.brine[0].pressure_osm_phase # specify state block - t.brine[0].flow_mass_phase_comp["Liq", "H2O"].fix( + tmp.brine[0].flow_mass_phase_comp["Liq", "H2O"].fix( value(feed_state_block.flow_mass_phase_comp["Liq", "H2O"]) * (1 - water_recovery_mass) ) - t.brine[0].flow_mass_phase_comp["Liq", "NaCl"].fix( - value(feed_state_block.flow_mass_phase_comp["Liq", "NaCl"]) * (1 - NaCl_passage) + tmp.brine[0].flow_mass_phase_comp["Liq", comp].fix( + value(feed_state_block.flow_mass_phase_comp["Liq", comp]) * (1 - salt_passage) ) - t.brine[0].pressure.fix( + tmp.brine[0].pressure.fix( 101325 ) # valid when osmotic pressure is independent of hydraulic pressure - t.brine[0].temperature.fix(value(feed_state_block.temperature)) - # calculate osmotic pressure - # since properties are created on demand, we must touch the property to create it - t.brine[0].pressure_osm_phase + tmp.brine[0].temperature.fix(value(feed_state_block.temperature)) + # solve state block - results = solve_indexed_blocks(solver, [t.brine]) + results = solve_indexed_blocks(solver, [tmp.brine]) assert_optimal_termination(results) - return value(t.brine[0].pressure_osm_phase["Liq"]) * (1 + over_pressure) + op_pressure = value(tmp.brine[0].pressure_osm_phase["Liq"]) * (1 + over_pressure) + + return op_pressure From d6a172ec53d9dadf06d8e7ddc653584df00ffaf4 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Sun, 23 Nov 2025 15:51:40 -0700 Subject: [PATCH 06/15] generic salt_passage as kwarg --- .../RO_with_energy_recovery/RO_with_energy_recovery.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py index 772c878ea5..7b7a8c1fdb 100644 --- a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py +++ b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py @@ -257,7 +257,7 @@ def set_operating_conditions( feed_state_block=m.fs.feed.properties[0], over_pressure=over_pressure, water_recovery_mass=water_recovery, - NaCl_passage=0.01, + salt_passage=0.01, solver=solver, ) m.fs.P1.control_volume.properties_out[0].pressure.fix(operating_pressure) From e97d19686a9a85496505fd123c1c75063469ab05 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Sun, 23 Nov 2025 15:51:48 -0700 Subject: [PATCH 07/15] initial test --- .../tests/test_reverse_osmosis_tools.py | 60 +++++++++++++++++++ 1 file changed, 60 insertions(+) create mode 100644 watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py diff --git a/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py b/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py new file mode 100644 index 0000000000..7ffcf06f46 --- /dev/null +++ b/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py @@ -0,0 +1,60 @@ +################################################################################# +# 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 pytest + +from pyomo.environ import ConcreteModel +from idaes.core import FlowsheetBlock +from idaes.core.util.model_statistics import degrees_of_freedom +from idaes.models.unit_models import Feed + +from watertap.property_models.seawater_prop_pack import ( + SeawaterParameterBlock, + SeawaterStateBlockData, +) +from watertap.tools.unit_models import calculate_operating_pressure +from watertap.unit_models import ReverseOsmosis0D +from watertap.core.solvers import get_solver + +solver = get_solver() + + +def build_seawater_prop_model(): + """ + Create feed model using seawater property package + """ + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = SeawaterParameterBlock() + m.fs.feed = Feed(property_package=m.fs.properties) + m.fs.feed.properties[0].pressure_osm_phase + m.fs.feed.properties[0].conc_mass_phase_comp + m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "H2O"].fix(0.965) + m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "TDS"].fix(0.035) + m.fs.feed.properties[0].temperature.fix(273 + 25) + m.fs.feed.properties[0].pressure.fix(101325) + + # Initialize and solve for the initial conditions + m.fs.feed.initialize() + results = solver.solve(m) + print(f"Solve termination {results.solver.termination_condition}") + # Unfix TDS mass flow + m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "TDS"].unfix() + + return m + + +m = build_seawater_prop_model() +osm = calculate_operating_pressure( + m.fs.feed.properties[0], water_recovery_mass=0.5, solver=solver +) +print(f"Calculated osmotic pressure: {osm:.2f} Pa") From 24bc795505242dedfbb371b57eea9a9d1f52e8ef Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Mon, 24 Nov 2025 10:22:00 -0700 Subject: [PATCH 08/15] all passage of MembraneChannel blocks also; error handling for other input bounds --- watertap/tools/unit_models/reverse_osmosis.py | 50 ++++++++++++++----- 1 file changed, 37 insertions(+), 13 deletions(-) diff --git a/watertap/tools/unit_models/reverse_osmosis.py b/watertap/tools/unit_models/reverse_osmosis.py index cc0dd6d5bf..44a162378d 100644 --- a/watertap/tools/unit_models/reverse_osmosis.py +++ b/watertap/tools/unit_models/reverse_osmosis.py @@ -23,6 +23,11 @@ from watertap.property_models.seawater_prop_pack import SeawaterStateBlockData from watertap.property_models.NaCl_prop_pack import NaClStateBlockData +from watertap.property_models.NaCl_T_dep_prop_pack import ( + NaClStateBlockData as NaClTDepStateBlockData, +) +from watertap.core import MembraneChannel0DBlock, MembraneChannel1DBlock +from watertap.core.solvers import get_solver __all__ = ["calculate_operating_pressure"] @@ -49,41 +54,60 @@ def calculate_operating_pressure( (default=0.01) solver: solver object to be used (default=None) """ + + if any( + isinstance(feed_state_block, cls) + for cls in [MembraneChannel0DBlock, MembraneChannel1DBlock] + ): + feed_state_block = feed_state_block.properties[0, 0] + if not any( isinstance(feed_state_block, cls) - for cls in [SeawaterStateBlockData, NaClStateBlockData] + for cls in [SeawaterStateBlockData, NaClStateBlockData, NaClTDepStateBlockData] ): raise TypeError( - "feed_state_block argument must be a SeawaterStateBlockData or NaClStateBlockData object" + "feed_state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock" ) - if isinstance(feed_state_block, NaClStateBlockData): + if not 1e-3 < salt_passage < 0.999: + raise ValueError("salt_passage argument must be between 0.001 and 0.999") + + if not 1e-3 < water_recovery_mass < 0.999: + raise ValueError("water_recovery_mass argument must be between 0.001 and 0.999") + + if isinstance(feed_state_block, NaClStateBlockData) or isinstance( + feed_state_block, NaClTDepStateBlockData + ): comp = "NaCl" if isinstance(feed_state_block, SeawaterStateBlockData): comp = "TDS" - + + if solver is None: + solver = get_solver() + tmp = ConcreteModel() # create temporary model prop = feed_state_block.config.parameters - tmp.brine = prop.build_state_block([0]) - tmp.brine[0].pressure_osm_phase + + tmp.feed = prop.build_state_block([0]) + tmp.feed[0].pressure_osm_phase # specify state block - tmp.brine[0].flow_mass_phase_comp["Liq", "H2O"].fix( + tmp.feed[0].flow_mass_phase_comp["Liq", "H2O"].fix( value(feed_state_block.flow_mass_phase_comp["Liq", "H2O"]) * (1 - water_recovery_mass) ) - tmp.brine[0].flow_mass_phase_comp["Liq", comp].fix( + tmp.feed[0].flow_mass_phase_comp["Liq", comp].fix( value(feed_state_block.flow_mass_phase_comp["Liq", comp]) * (1 - salt_passage) ) - tmp.brine[0].pressure.fix( + tmp.feed[0].pressure.fix( 101325 ) # valid when osmotic pressure is independent of hydraulic pressure - tmp.brine[0].temperature.fix(value(feed_state_block.temperature)) - + tmp.feed[0].temperature.fix(value(feed_state_block.temperature)) + # solve state block - results = solve_indexed_blocks(solver, [tmp.brine]) + results = solve_indexed_blocks(solver, [tmp.feed]) assert_optimal_termination(results) - op_pressure = value(tmp.brine[0].pressure_osm_phase["Liq"]) * (1 + over_pressure) + op_pressure = value(tmp.feed[0].pressure_osm_phase["Liq"]) * (1 + over_pressure) return op_pressure From 85e04e5e5f3b5a464fe4a3c86fefeb5522e930c4 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Mon, 24 Nov 2025 10:22:34 -0700 Subject: [PATCH 09/15] test errors and other valid prop models --- .../tests/test_reverse_osmosis_tools.py | 185 ++++++++++++++++-- 1 file changed, 168 insertions(+), 17 deletions(-) diff --git a/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py b/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py index 7ffcf06f46..ef22c192b9 100644 --- a/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py +++ b/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py @@ -14,15 +14,17 @@ from pyomo.environ import ConcreteModel from idaes.core import FlowsheetBlock -from idaes.core.util.model_statistics import degrees_of_freedom from idaes.models.unit_models import Feed -from watertap.property_models.seawater_prop_pack import ( - SeawaterParameterBlock, - SeawaterStateBlockData, +from watertap.property_models.seawater_prop_pack import SeawaterParameterBlock +from watertap.property_models.NaCl_prop_pack import NaClParameterBlock +from watertap.property_models.NaCl_T_dep_prop_pack import ( + NaClParameterBlock as NaClTDepParameterBlock, ) +from watertap.property_models.multicomp_aq_sol_prop_pack import MCASParameterBlock + from watertap.tools.unit_models import calculate_operating_pressure -from watertap.unit_models import ReverseOsmosis0D +from watertap.unit_models import ReverseOsmosis0D, ReverseOsmosis1D from watertap.core.solvers import get_solver solver = get_solver() @@ -36,25 +38,174 @@ def build_seawater_prop_model(): m.fs = FlowsheetBlock(dynamic=False) m.fs.properties = SeawaterParameterBlock() m.fs.feed = Feed(property_package=m.fs.properties) - m.fs.feed.properties[0].pressure_osm_phase - m.fs.feed.properties[0].conc_mass_phase_comp m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "H2O"].fix(0.965) m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "TDS"].fix(0.035) m.fs.feed.properties[0].temperature.fix(273 + 25) m.fs.feed.properties[0].pressure.fix(101325) - # Initialize and solve for the initial conditions m.fs.feed.initialize() - results = solver.solve(m) - print(f"Solve termination {results.solver.termination_condition}") - # Unfix TDS mass flow - m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "TDS"].unfix() return m -m = build_seawater_prop_model() -osm = calculate_operating_pressure( - m.fs.feed.properties[0], water_recovery_mass=0.5, solver=solver -) -print(f"Calculated osmotic pressure: {osm:.2f} Pa") +def build_nacl_prop_model(): + """ + Create feed model using NaCl property package + """ + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = NaClParameterBlock() + m.fs.feed = Feed(property_package=m.fs.properties) + m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "H2O"].fix(0.965) + m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "NaCl"].fix(0.035) + m.fs.feed.properties[0].temperature.fix(273 + 25) + m.fs.feed.properties[0].pressure.fix(101325) + + m.fs.feed.initialize() + + return m + + +def build_nacl_tdep_prop_model(): + """ + Create feed model using NaCl temp dependence property package + """ + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = NaClTDepParameterBlock() + m.fs.feed = Feed(property_package=m.fs.properties) + m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "H2O"].fix(0.965) + m.fs.feed.properties[0].flow_mass_phase_comp["Liq", "NaCl"].fix(0.035) + m.fs.feed.properties[0].temperature.fix(273 + 25) + m.fs.feed.properties[0].pressure.fix(101325) + + m.fs.feed.initialize() + + return m + + +def build_ro0d_model(): + """ + Create RO0D model for testing + """ + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = SeawaterParameterBlock() + + m.fs.RO = ReverseOsmosis0D( + property_package=m.fs.properties, + concentration_polarization_type="none", + mass_transfer_coefficient="none", + has_pressure_change=False, + ) + + m.fs.RO.inlet.flow_mass_phase_comp[0, "Liq", "H2O"].fix(0.965) + m.fs.RO.inlet.flow_mass_phase_comp[0, "Liq", "TDS"].fix(0.035) + m.fs.RO.inlet.temperature.fix(273 + 25) + m.fs.RO.inlet.pressure.fix(101325) + + return m + + +def build_ro1d_model(): + """ + Create RO1D model for testing + """ + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = NaClTDepParameterBlock() + + m.fs.RO = ReverseOsmosis1D( + property_package=m.fs.properties, + concentration_polarization_type="none", + mass_transfer_coefficient="none", + has_pressure_change=False, + ) + + m.fs.RO.inlet.flow_mass_phase_comp[0, "Liq", "H2O"].fix(0.965) + m.fs.RO.inlet.flow_mass_phase_comp[0, "Liq", "NaCl"].fix(0.035) + m.fs.RO.inlet.temperature.fix(273 + 25) + m.fs.RO.inlet.pressure.fix(101325) + + return m + + +@pytest.mark.component +def test_calculate_operating_pressure_sw(): + + m = build_seawater_prop_model() + osm = calculate_operating_pressure(m.fs.feed.properties[0]) + assert pytest.approx(osm, rel=1e-3) == 6034067.12 + + +@pytest.mark.component +def test_calculate_operating_pressure_nacl(): + + m = build_nacl_prop_model() + osm = calculate_operating_pressure(m.fs.feed.properties[0]) + assert pytest.approx(osm, rel=1e-3) == 6624988.52 + + +@pytest.mark.component +def test_calculate_operating_pressure_nacl_tdep(): + + m = build_nacl_tdep_prop_model() + osm = calculate_operating_pressure(m.fs.feed.properties[0]) + assert pytest.approx(osm, rel=1e-3) == 6607568.15 + + +@pytest.mark.component +def test_calculate_operating_pressure_ro0d(): + + m = build_ro0d_model() + osm1 = calculate_operating_pressure(m.fs.RO.feed_side.properties[0, 0]) + osm2 = calculate_operating_pressure(m.fs.RO.feed_side) + assert osm1 == osm2 + assert pytest.approx(osm1, rel=1e-3) == 6034067.12 + + +@pytest.mark.component +def test_calculate_operating_pressure_ro1d(): + + m = build_ro1d_model() + osm1 = calculate_operating_pressure(m.fs.RO.feed_side.properties[0, 0]) + osm2 = calculate_operating_pressure(m.fs.RO.feed_side) + assert osm1 == osm2 + assert pytest.approx(osm1, rel=1e-3) == 6607568.14 + + +@pytest.mark.unit +def test_calculate_operating_pressure_errors(): + + with pytest.raises( + TypeError, + match="feed_state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock", + ): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = MCASParameterBlock(solute_list=["Na_+"]) + m.fs.feed = Feed(property_package=m.fs.properties) + calculate_operating_pressure(feed_state_block=m.fs.feed.properties[0]) + + with pytest.raises( + TypeError, + match="feed_state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock", + ): + m = build_ro0d_model() + calculate_operating_pressure(feed_state_block=m.fs.RO.inlet) + + with pytest.raises( + ValueError, match="salt_passage argument must be between 0.001 and 0.999" + ): + m = build_seawater_prop_model() + calculate_operating_pressure( + feed_state_block=m.fs.feed.properties[0], salt_passage=1.1 + ) + + with pytest.raises( + ValueError, match="water_recovery_mass argument must be between 0.001 and 0.999" + ): + m = build_nacl_prop_model() + calculate_operating_pressure( + feed_state_block=m.fs.feed.properties[0], water_recovery_mass=2.5 + ) From 8ccdb287f61b32b3f10f04638b96774bb7c47366 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Mon, 24 Nov 2025 10:35:52 -0700 Subject: [PATCH 10/15] black --- .../RO_with_energy_recovery/RO_with_energy_recovery.py | 1 - 1 file changed, 1 deletion(-) diff --git a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py index 7b7a8c1fdb..b1fadafd7a 100644 --- a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py +++ b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py @@ -313,7 +313,6 @@ def set_operating_conditions( ) - def solve(blk, solver=None, tee=False, check_termination=True): if solver is None: solver = get_solver() From 91b60b7cd751e75958cd2a41916dc82cf2346146 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Mon, 24 Nov 2025 11:11:28 -0700 Subject: [PATCH 11/15] comments: kwargs to state_block, over_pressure_factor; find comp dynamically; check over_pressure_factor >= 1; raise if not optimal termination --- watertap/tools/unit_models/reverse_osmosis.py | 57 +++++++++++-------- .../tests/test_reverse_osmosis_tools.py | 21 +++++-- 2 files changed, 48 insertions(+), 30 deletions(-) diff --git a/watertap/tools/unit_models/reverse_osmosis.py b/watertap/tools/unit_models/reverse_osmosis.py index 44a162378d..ab117914f0 100644 --- a/watertap/tools/unit_models/reverse_osmosis.py +++ b/watertap/tools/unit_models/reverse_osmosis.py @@ -16,7 +16,7 @@ from pyomo.environ import ( ConcreteModel, - assert_optimal_termination, + check_optimal_termination, value, ) from idaes.core.util.initialization import solve_indexed_blocks @@ -34,8 +34,8 @@ def calculate_operating_pressure( - feed_state_block=None, - over_pressure=0.15, + state_block=None, + over_pressure_factor=1.15, water_recovery_mass=0.5, salt_passage=0.01, solver=None, @@ -44,9 +44,9 @@ def calculate_operating_pressure( Estimate operating pressure for RO unit model given the following arguments: Arguments: - feed_state_block: the state block of the RO feed that has the non-pressure state + state_block: the state block of the RO feed that has the non-pressure state variables initialized to their values (default=None) - over_pressure: the amount of operating pressure above the brine osmotic pressure + over_pressure_factor: the amount of operating pressure above the brine osmotic pressure represented as a fraction (default=0.15) water_recovery_mass: the mass-based fraction of inlet H2O that becomes permeate (default=0.5) @@ -56,17 +56,17 @@ def calculate_operating_pressure( """ if any( - isinstance(feed_state_block, cls) + isinstance(state_block, cls) for cls in [MembraneChannel0DBlock, MembraneChannel1DBlock] ): - feed_state_block = feed_state_block.properties[0, 0] + state_block = state_block.properties[0, 0] if not any( - isinstance(feed_state_block, cls) + isinstance(state_block, cls) for cls in [SeawaterStateBlockData, NaClStateBlockData, NaClTDepStateBlockData] ): raise TypeError( - "feed_state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock" + "state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock" ) if not 1e-3 < salt_passage < 0.999: @@ -75,39 +75,48 @@ def calculate_operating_pressure( if not 1e-3 < water_recovery_mass < 0.999: raise ValueError("water_recovery_mass argument must be between 0.001 and 0.999") - if isinstance(feed_state_block, NaClStateBlockData) or isinstance( - feed_state_block, NaClTDepStateBlockData - ): - comp = "NaCl" - if isinstance(feed_state_block, SeawaterStateBlockData): - comp = "TDS" + if not over_pressure_factor >= 1.0: + raise ValueError( + "over_pressure_factor argument must be greater than or equal to 1.0" + ) + + print(over_pressure_factor, water_recovery_mass, salt_passage) + + comp = [c for c in state_block.component_list if c != "H2O"][0] + + if comp not in ["NaCl", "TDS"]: + raise ValueError( + f"salt_passage calculation only supported for NaCl or TDS components but found {comp}" + ) if solver is None: solver = get_solver() tmp = ConcreteModel() # create temporary model - prop = feed_state_block.config.parameters + prop = state_block.config.parameters tmp.feed = prop.build_state_block([0]) tmp.feed[0].pressure_osm_phase # specify state block tmp.feed[0].flow_mass_phase_comp["Liq", "H2O"].fix( - value(feed_state_block.flow_mass_phase_comp["Liq", "H2O"]) + value(state_block.flow_mass_phase_comp["Liq", "H2O"]) * (1 - water_recovery_mass) ) tmp.feed[0].flow_mass_phase_comp["Liq", comp].fix( - value(feed_state_block.flow_mass_phase_comp["Liq", comp]) * (1 - salt_passage) + value(state_block.flow_mass_phase_comp["Liq", comp]) * (1 - salt_passage) ) - tmp.feed[0].pressure.fix( - 101325 - ) # valid when osmotic pressure is independent of hydraulic pressure - tmp.feed[0].temperature.fix(value(feed_state_block.temperature)) + tmp.feed[0].temperature.fix(value(state_block.temperature)) + tmp.feed[0].pressure.fix(101325) # solve state block results = solve_indexed_blocks(solver, [tmp.feed]) - assert_optimal_termination(results) - op_pressure = value(tmp.feed[0].pressure_osm_phase["Liq"]) * (1 + over_pressure) + if not check_optimal_termination(results): + raise RuntimeError( + "Failed to solve temporary state block for operating pressure" + ) + + op_pressure = value(tmp.feed[0].pressure_osm_phase["Liq"]) * over_pressure_factor return op_pressure diff --git a/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py b/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py index ef22c192b9..dc8a5a60a2 100644 --- a/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py +++ b/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py @@ -179,27 +179,27 @@ def test_calculate_operating_pressure_errors(): with pytest.raises( TypeError, - match="feed_state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock", + match="state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock", ): m = ConcreteModel() m.fs = FlowsheetBlock(dynamic=False) m.fs.properties = MCASParameterBlock(solute_list=["Na_+"]) m.fs.feed = Feed(property_package=m.fs.properties) - calculate_operating_pressure(feed_state_block=m.fs.feed.properties[0]) + calculate_operating_pressure(state_block=m.fs.feed.properties[0]) with pytest.raises( TypeError, - match="feed_state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock", + match="state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock", ): m = build_ro0d_model() - calculate_operating_pressure(feed_state_block=m.fs.RO.inlet) + calculate_operating_pressure(state_block=m.fs.RO.inlet) with pytest.raises( ValueError, match="salt_passage argument must be between 0.001 and 0.999" ): m = build_seawater_prop_model() calculate_operating_pressure( - feed_state_block=m.fs.feed.properties[0], salt_passage=1.1 + state_block=m.fs.feed.properties[0], salt_passage=1.1 ) with pytest.raises( @@ -207,5 +207,14 @@ def test_calculate_operating_pressure_errors(): ): m = build_nacl_prop_model() calculate_operating_pressure( - feed_state_block=m.fs.feed.properties[0], water_recovery_mass=2.5 + state_block=m.fs.feed.properties[0], water_recovery_mass=2.5 + ) + + with pytest.raises( + ValueError, + match="over_pressure_factor argument must be greater than or equal to 1.0", + ): + m = build_nacl_prop_model() + calculate_operating_pressure( + state_block=m.fs.feed.properties[0], over_pressure_factor=0.9 ) From a78f2f0f735d912773c86affeab3d94c98d115ab Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Mon, 24 Nov 2025 11:17:19 -0700 Subject: [PATCH 12/15] default salt_passage=0; update test --- watertap/tools/unit_models/reverse_osmosis.py | 6 +++--- .../unit_models/tests/test_reverse_osmosis_tools.py | 12 ++++++------ 2 files changed, 9 insertions(+), 9 deletions(-) diff --git a/watertap/tools/unit_models/reverse_osmosis.py b/watertap/tools/unit_models/reverse_osmosis.py index ab117914f0..9368a3a897 100644 --- a/watertap/tools/unit_models/reverse_osmosis.py +++ b/watertap/tools/unit_models/reverse_osmosis.py @@ -37,7 +37,7 @@ def calculate_operating_pressure( state_block=None, over_pressure_factor=1.15, water_recovery_mass=0.5, - salt_passage=0.01, + salt_passage=0, solver=None, ): """ @@ -69,8 +69,8 @@ def calculate_operating_pressure( "state_block must be created with SeawaterParameterBlock, NaClParameterBlock, or NaClTDepParameterBlock" ) - if not 1e-3 < salt_passage < 0.999: - raise ValueError("salt_passage argument must be between 0.001 and 0.999") + if not 0 <= salt_passage < 0.999: + raise ValueError("salt_passage argument must be between 0 and 0.999") if not 1e-3 < water_recovery_mass < 0.999: raise ValueError("water_recovery_mass argument must be between 0.001 and 0.999") diff --git a/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py b/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py index dc8a5a60a2..40dc94d711 100644 --- a/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py +++ b/watertap/tools/unit_models/tests/test_reverse_osmosis_tools.py @@ -135,7 +135,7 @@ def test_calculate_operating_pressure_sw(): m = build_seawater_prop_model() osm = calculate_operating_pressure(m.fs.feed.properties[0]) - assert pytest.approx(osm, rel=1e-3) == 6034067.12 + assert pytest.approx(osm, rel=1e-3) == 6098990.3 @pytest.mark.component @@ -143,7 +143,7 @@ def test_calculate_operating_pressure_nacl(): m = build_nacl_prop_model() osm = calculate_operating_pressure(m.fs.feed.properties[0]) - assert pytest.approx(osm, rel=1e-3) == 6624988.52 + assert pytest.approx(osm, rel=1e-3) == 6695261.0 @pytest.mark.component @@ -151,7 +151,7 @@ def test_calculate_operating_pressure_nacl_tdep(): m = build_nacl_tdep_prop_model() osm = calculate_operating_pressure(m.fs.feed.properties[0]) - assert pytest.approx(osm, rel=1e-3) == 6607568.15 + assert pytest.approx(osm, rel=1e-3) == 6677498.9 @pytest.mark.component @@ -161,7 +161,7 @@ def test_calculate_operating_pressure_ro0d(): osm1 = calculate_operating_pressure(m.fs.RO.feed_side.properties[0, 0]) osm2 = calculate_operating_pressure(m.fs.RO.feed_side) assert osm1 == osm2 - assert pytest.approx(osm1, rel=1e-3) == 6034067.12 + assert pytest.approx(osm1, rel=1e-3) == 6098990.3 @pytest.mark.component @@ -171,7 +171,7 @@ def test_calculate_operating_pressure_ro1d(): osm1 = calculate_operating_pressure(m.fs.RO.feed_side.properties[0, 0]) osm2 = calculate_operating_pressure(m.fs.RO.feed_side) assert osm1 == osm2 - assert pytest.approx(osm1, rel=1e-3) == 6607568.14 + assert pytest.approx(osm1, rel=1e-3) == 6677498.9 @pytest.mark.unit @@ -195,7 +195,7 @@ def test_calculate_operating_pressure_errors(): calculate_operating_pressure(state_block=m.fs.RO.inlet) with pytest.raises( - ValueError, match="salt_passage argument must be between 0.001 and 0.999" + ValueError, match="salt_passage argument must be between 0 and 0.999" ): m = build_seawater_prop_model() calculate_operating_pressure( From d2dd303febf45cb6231378fbef8fabbf0b8d3fad Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Mon, 24 Nov 2025 11:29:03 -0700 Subject: [PATCH 13/15] apply changes to RO w ERD fs and test --- .../RO_with_energy_recovery/RO_with_energy_recovery.py | 6 +++--- .../tests/test_RO_with_energy_recovery_simulation.py | 4 ++-- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py index b1fadafd7a..266ee2a770 100644 --- a/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py +++ b/watertap/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py @@ -202,7 +202,7 @@ def build(erd_type=ERDtype.pressure_exchanger): def set_operating_conditions( m, water_recovery=0.5, - over_pressure=0.3, + over_pressure_factor=1.3, flow_vol=1e-3, salt_mass_conc=35e-3, solver=None, @@ -254,8 +254,8 @@ def set_operating_conditions( # pump 1, high pressure pump, 2 degrees of freedom (efficiency and outlet pressure) m.fs.P1.efficiency_pump.fix(0.80) # pump efficiency [-] operating_pressure = calculate_operating_pressure( - feed_state_block=m.fs.feed.properties[0], - over_pressure=over_pressure, + state_block=m.fs.feed.properties[0], + over_pressure_factor=over_pressure_factor, water_recovery_mass=water_recovery, salt_passage=0.01, solver=solver, diff --git a/watertap/flowsheets/RO_with_energy_recovery/tests/test_RO_with_energy_recovery_simulation.py b/watertap/flowsheets/RO_with_energy_recovery/tests/test_RO_with_energy_recovery_simulation.py index 936aa9078c..4320ff9e3b 100644 --- a/watertap/flowsheets/RO_with_energy_recovery/tests/test_RO_with_energy_recovery_simulation.py +++ b/watertap/flowsheets/RO_with_energy_recovery/tests/test_RO_with_energy_recovery_simulation.py @@ -142,7 +142,7 @@ def test_set_operating_conditions(self, system_frame): m = system_frame set_operating_conditions( - m, water_recovery=0.5, over_pressure=0.3, solver=solver + m, water_recovery=0.5, over_pressure_factor=1.3, solver=solver ) # check fixed variables @@ -513,7 +513,7 @@ def test_set_operating_conditions(self, system_frame): m = system_frame set_operating_conditions( - m, water_recovery=0.5, over_pressure=0.3, solver=solver + m, water_recovery=0.5, over_pressure_factor=1.3, solver=solver ) # check fixed variables From ea50b84e22cff56d3fd6aa5971a0fce764b78b80 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Mon, 24 Nov 2025 11:35:46 -0700 Subject: [PATCH 14/15] better way to find comp --- watertap/tools/unit_models/reverse_osmosis.py | 18 ++++++------------ 1 file changed, 6 insertions(+), 12 deletions(-) diff --git a/watertap/tools/unit_models/reverse_osmosis.py b/watertap/tools/unit_models/reverse_osmosis.py index 9368a3a897..298ec85b58 100644 --- a/watertap/tools/unit_models/reverse_osmosis.py +++ b/watertap/tools/unit_models/reverse_osmosis.py @@ -44,15 +44,11 @@ def calculate_operating_pressure( Estimate operating pressure for RO unit model given the following arguments: Arguments: - state_block: the state block of the RO feed that has the non-pressure state - variables initialized to their values (default=None) - over_pressure_factor: the amount of operating pressure above the brine osmotic pressure - represented as a fraction (default=0.15) - water_recovery_mass: the mass-based fraction of inlet H2O that becomes permeate - (default=0.5) - salt_passage: the mass-based fraction of inlet salt that becomes permeate - (default=0.01) - solver: solver object to be used (default=None) + state_block: the state block of the RO feed that has the non-pressure state variables set to desired values (default=None) + over_pressure_factor: the amount of operating pressure above the brine osmotic pressure represented as a fraction (default=1.15) + water_recovery_mass: the mass-based fraction of inlet H2O that becomes permeate (default=0.5) + salt_passage: the mass-based fraction of inlet salt that becomes permeate (default=0) + solver: solver object to be used (default=None) """ if any( @@ -80,9 +76,7 @@ def calculate_operating_pressure( "over_pressure_factor argument must be greater than or equal to 1.0" ) - print(over_pressure_factor, water_recovery_mass, salt_passage) - - comp = [c for c in state_block.component_list if c != "H2O"][0] + comp = state_block.params.solute_set.first() if comp not in ["NaCl", "TDS"]: raise ValueError( From 02e94fba796916fd9de100a1608ff1b6db4f4e43 Mon Sep 17 00:00:00 2001 From: kurbansitterley Date: Mon, 24 Nov 2025 11:50:05 -0700 Subject: [PATCH 15/15] update kwarg --- .../RO_with_energy_recovery/monte_carlo_sampling_RO_ERD.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/watertap/flowsheets/RO_with_energy_recovery/monte_carlo_sampling_RO_ERD.py b/watertap/flowsheets/RO_with_energy_recovery/monte_carlo_sampling_RO_ERD.py index 87a3051d0d..9f4b9a8f9c 100644 --- a/watertap/flowsheets/RO_with_energy_recovery/monte_carlo_sampling_RO_ERD.py +++ b/watertap/flowsheets/RO_with_energy_recovery/monte_carlo_sampling_RO_ERD.py @@ -69,7 +69,9 @@ def build_model( # Build, set, and initialize the system (these steps will change depending on the underlying model) m = build() - set_operating_conditions(m, water_recovery=0.5, over_pressure=0.3, solver=solver) + set_operating_conditions( + m, water_recovery=0.5, over_pressure_factor=1.3, solver=solver + ) initialize_system(m, solver=solver) # Check if we need to read in the default model values from a file