From ab9870b9263238ba6d8e457612e0335a0a8709a3 Mon Sep 17 00:00:00 2001 From: luohezhiming <98901358+luohezhiming@users.noreply.github.com> Date: Mon, 25 Mar 2024 20:12:37 -0400 Subject: [PATCH] Add more detailed pressure exchanger (#1264) * add * delete redundant files * Add more detailed pressure exchanger * code linting * minor changes * add test for multi_comp property package * revise documentation * change efficiency to universal variable * temporary version * Revise doc * code pylint * Resolve test problem * New version of PX * revise bounds for LPD * revise leakage equations * update doc * code linting * revise typo * revise typos * Update docs/technical_reference/unit_models/pressure_exchanger.rst Co-authored-by: Adam Atia * Update docs/technical_reference/unit_models/pressure_exchanger.rst Co-authored-by: Adam Atia * Update docs/technical_reference/unit_models/pressure_exchanger.rst Co-authored-by: Adam Atia * Update docs/technical_reference/unit_models/pressure_exchanger.rst Co-authored-by: Adam Atia * Update docs/technical_reference/unit_models/pressure_exchanger.rst Co-authored-by: Adam Atia * Update watertap/unit_models/pressure_exchanger.py Co-authored-by: Adam Atia * Update watertap/unit_models/pressure_exchanger.py Co-authored-by: Adam Atia * add comments * resolve issue with inden * rename LPS to feed_side and HPS to brine_side * resolve tests problems * revise test file * efficiency equation change to ener_out/ener_in * has_leakage for initialization * Adjust dye desal with RO test solutions --------- Co-authored-by: Adam Atia Co-authored-by: MarcusHolly Co-authored-by: MarcusHolly <96305519+MarcusHolly@users.noreply.github.com> --- .../unit_models/pressure_exchanger.rst | 66 +- .../costing/unit_models/pressure_exchanger.py | 2 +- .../RO_with_energy_recovery.py | 24 +- .../default_configuration.yaml | 52 +- ...test_RO_with_energy_recovery_simulation.py | 32 +- .../seawater_RO_desalination.py | 16 +- .../dye_desalination_withRO.py | 16 +- .../tests/test_dye_desalination_withRO.py | 4 +- watertap/unit_models/pressure_exchanger.py | 391 +++++++---- .../tests/test_pressure_exchanger.py | 639 ++++++++++++------ 10 files changed, 826 insertions(+), 416 deletions(-) diff --git a/docs/technical_reference/unit_models/pressure_exchanger.rst b/docs/technical_reference/unit_models/pressure_exchanger.rst index 9df0d91bc0..6493205645 100644 --- a/docs/technical_reference/unit_models/pressure_exchanger.rst +++ b/docs/technical_reference/unit_models/pressure_exchanger.rst @@ -5,8 +5,7 @@ This pressure exchanger unit model: * is isothermal * supports a single liquid phase only * supports steady-state only - * assumes no mixing or leakage between the low and high pressure side - * assumes equal flowrates on both sides + * supports leakage and mixing effect .. index:: pair: watertap.unit_models.pressure_exchanger;pressure_exchanger @@ -32,9 +31,16 @@ Where the system is also subject to following constraints: Figure 1. Schematic representation of an energy recovery system using a pressure exchanger. +When setting the ``has_mixing`` configuration option to ``True``, there is 1 additional variable ``mixing_vol`` that must be fixed. + +When setting the ``has_leakage`` configuration option to ``True``, there is 1 additional variable ``leakage_vol`` that must be fixed. + +When setting the ``pressure_exchange_calculation`` configuration option to ``PressureExchangeType.high_pressure_difference``, +there are 2 additional variables ``high_pressure_difference`` and ``low_pressure_difference`` that must be fixed. Instead, ``efficiency_pressure_exchanger`` is unfixed. + Model Structure ------------------ -The pressure exchanger model consists of 2 `ControlVolume0DBlocks`: one for the low-pressure side and high-pressure side. +The pressure exchanger model consists of 2 `ControlVolume0DBlocks`: one for the feed side and brine side. Sets ---- @@ -56,8 +62,14 @@ The pressure exchanger unit model includes the following variables: :header: "Description", "Symbol", "Variable Name", "Index", "Units", "Pyomo Type" "Efficiency", ":math:`\eta`", "efficiency_pressure_exchanger", "[t]", ":math:`\text{dimensionless}`", "Var" + "Volumetric leakage fraction", ":math:`\delta`", "leakage_vol", "[t]", ":math:`\text{dimensionless}`", "Var" + "Volumetric mixing fraction", ":math:`\chi`", "mixing_vol", "[t]", ":math:`\text{dimensionless}`", "Var" + "High pressure difference*", ":math:`HPD`", "high_pressure_difference", "[t]", ":math:`\text{Pa}`", "Var" + "Low pressure difference*", ":math:`LPD`", "low_pressure_difference", "[t]", ":math:`\text{Pa}`", "Var" + +\*High pressure difference and low pressure difference are non-negative values -Each control volume (i.e. `low_presssure_side`, and `high_pressure_side`) has the following variables of interest: +Each control volume (i.e. `feed_side`, and `brine_side`) has the following variables of interest: .. csv-table:: :header: "Description", "Symbol", "Variable Name", "Index", "Units", "Pyomo Type" @@ -72,7 +84,9 @@ Each property block on both control volumes (i.e. `properties_in` and `propertie .. csv-table:: :header: "Description", "Symbol", "Variable Name", "Index", "Units", "Pyomo Type" + "Mass flowrate", ":math:`M`", "flow_mass_phase_comp", "[p, j]", "\*", "Var" "Volumetric flowrate", ":math:`Q`", "flow_vol", "None", "\*", "Var" + "Concentration", ":math:`C`", "conc_mass_phase_comp", "[p, j]", "\*", "Var" "Temperature", ":math:`T`", "temperature", "[t]", "\*", "Var" "Pressure", ":math:`P`", "pressure", "[t]", "\*", "Var" @@ -81,6 +95,7 @@ Each property block on both control volumes (i.e. `properties_in` and `propertie Equations ----------- +if ``has_leakage`` and ``has_mixing`` are set to default (``False``): .. csv-table:: :header: "Description", "Equation" @@ -88,10 +103,47 @@ Equations "Mass balance for each side", ":math:`M_{out, j} = M_{in, j}`" "Momentum balance for each side", ":math:`P_{out} = P_{in} + ΔP`" "Isothermal assumption for each side", ":math:`T_{out} = T_{in}`" - "Equal volumetric flowrate*", ":math:`Q_{LPS} = Q_{HPS}`" - "Pressure transfer*", ":math:`ΔP_{LPS} = - \eta ΔP_{HPS}`" + "Equal volumetric flowrate*", ":math:`Q_{out, F} = Q_{in, F}`" + "Equal pressure*", ":math:`P_{out, B} = P_{in, F}`" + "Pressure transfer*", ":math:`ΔP_{F} = - \eta ΔP_{B}`" + +\* F stands for feed side, B stands for brine side + +if ``has_leakage`` is set to ``True``, then the equal volumetric flowrate equation is replaced by: + +.. csv-table:: + :header: "Description", "Equation" + + "Equal volumetric flowrate", ":math:`Q_{out, F} = (1 - \delta) Q_{in, B}`" + +if ``has_mixing`` is set to ``True``, the mass balance equations for each side become: + +.. csv-table:: + :header: "Description", "Equation" + + "Mass balance for each side*", ":math:`M_{out, j} = M_{in, j} + MTT_{j}`" + +\* MTT is mass transfer term into the control volume + +and there are 3 additional constraints: + +.. csv-table:: + :header: "Description", "Equation" + + "Mixing effect of solute*", ":math:`C_{out, F} = C_{in, F} (1-\chi) + C_{in, B} \chi`" + "Linking mass transfer terms", ":math:`MTT_{j, F} = -MTT_{j, B}`" + "Equal feed side volumetric flowrate", ":math:`Q_{out, F} = Q_{in, F}`" + +\* C represents solute concentration + +if ``pressure_change_calculation`` is set to ``PressureExchangeType.high_pressure_difference``, +then there is 1 additional constraint and the equal pressure equation is replaced: + +.. csv-table:: + :header: "Description", "Equation" -\* LPS stands for low pressure side, HPS stands for high pressure side + "Pressure drop across B,in and F, out", ":math:`P_{out, F} + HPD = P_{in, B}`" + "Pressure drop across B, out and F, in", ":math:`P_{out, B} = P_{in, F} + LPD`" Class Documentation ------------------- diff --git a/watertap/costing/unit_models/pressure_exchanger.py b/watertap/costing/unit_models/pressure_exchanger.py index 371b8caeac..a2dffc3dcc 100644 --- a/watertap/costing/unit_models/pressure_exchanger.py +++ b/watertap/costing/unit_models/pressure_exchanger.py @@ -37,7 +37,7 @@ def cost_pressure_exchanger(blk): blk, blk.costing_package.pressure_exchanger.cost, pyo.units.convert( - blk.unit_model.low_pressure_side.properties_in[0].flow_vol, + blk.unit_model.feed_side.properties_in[0].flow_vol, (pyo.units.meter**3 / pyo.units.hours), ), ) diff --git a/watertap/examples/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py b/watertap/examples/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py index d90bd46146..79c2723acc 100644 --- a/watertap/examples/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py +++ b/watertap/examples/flowsheets/RO_with_energy_recovery/RO_with_energy_recovery.py @@ -152,14 +152,10 @@ def build(erd_type=ERDtype.pressure_exchanger): m.fs.s03 = Arc(source=m.fs.P1.outlet, destination=m.fs.M1.P1) m.fs.s04 = Arc(source=m.fs.M1.outlet, destination=m.fs.RO.inlet) m.fs.s05 = Arc(source=m.fs.RO.permeate, destination=m.fs.product.inlet) - m.fs.s06 = Arc( - source=m.fs.RO.retentate, destination=m.fs.PXR.high_pressure_inlet - ) - m.fs.s07 = Arc( - source=m.fs.PXR.high_pressure_outlet, destination=m.fs.disposal.inlet - ) - m.fs.s08 = Arc(source=m.fs.S1.PXR, destination=m.fs.PXR.low_pressure_inlet) - m.fs.s09 = Arc(source=m.fs.PXR.low_pressure_outlet, destination=m.fs.P2.inlet) + m.fs.s06 = Arc(source=m.fs.RO.retentate, destination=m.fs.PXR.brine_inlet) + m.fs.s07 = Arc(source=m.fs.PXR.brine_outlet, destination=m.fs.disposal.inlet) + m.fs.s08 = Arc(source=m.fs.S1.PXR, destination=m.fs.PXR.feed_inlet) + m.fs.s09 = Arc(source=m.fs.PXR.feed_outlet, destination=m.fs.P2.inlet) m.fs.s10 = Arc(source=m.fs.P2.outlet, destination=m.fs.M1.P2) elif erd_type == ERDtype.pump_as_turbine: m.fs.s01 = Arc(source=m.fs.feed.outlet, destination=m.fs.P1.inlet) @@ -184,8 +180,8 @@ def build(erd_type=ERDtype.pressure_exchanger): if erd_type == ERDtype.pressure_exchanger: iscale.set_scaling_factor(m.fs.P2.control_volume.work, 1e-3) - iscale.set_scaling_factor(m.fs.PXR.low_pressure_side.work, 1e-3) - iscale.set_scaling_factor(m.fs.PXR.high_pressure_side.work, 1e-3) + iscale.set_scaling_factor(m.fs.PXR.feed_side.work, 1e-3) + iscale.set_scaling_factor(m.fs.PXR.brine_side.work, 1e-3) # touch properties used in specifying and initializing the model m.fs.S1.mixed_state[0].mass_frac_phase_comp m.fs.S1.PXR_state[0].flow_vol_phase["Liq"] @@ -407,7 +403,7 @@ def initialize_pressure_exchanger(m, optarg): # ---initialize splitter and pressure exchanger--- # pressure exchanger high pressure inlet propagate_state(m.fs.s06) # propagate to PXR high pressure inlet from RO retentate - m.fs.PXR.high_pressure_side.properties_in.initialize(optarg=optarg) + m.fs.PXR.brine_side.properties_in.initialize(optarg=optarg) # splitter inlet propagate_state(m.fs.s01) # propagate to splitter inlet from feed @@ -421,7 +417,7 @@ def initialize_pressure_exchanger(m, optarg): "flow_vol_phase", "Liq", ): value( # same volumetric flow rate as PXR high pressure inlet - m.fs.PXR.high_pressure_side.properties_in[0].flow_vol_phase["Liq"] + m.fs.PXR.brine_side.properties_in[0].flow_vol_phase["Liq"] ), ("mass_frac_phase_comp", ("Liq", "NaCl")): value( m.fs.S1.mixed_state[0].mass_frac_phase_comp["Liq", "NaCl"] @@ -627,12 +623,12 @@ def print_state(s, b): print_state("Split 1 ", m.fs.S1.P1) print_state("P1 out ", m.fs.P1.outlet) print_state("Split 2 ", m.fs.S1.PXR) - print_state("PXR LP out", m.fs.PXR.low_pressure_outlet) + print_state("PXR feed out", m.fs.PXR.feed_outlet) print_state("P2 out ", m.fs.P2.outlet) print_state("Mix out ", m.fs.M1.outlet) print_state("RO perm ", m.fs.RO.permeate) print_state("RO reten ", m.fs.RO.retentate) - print_state("PXR HP out", m.fs.PXR.high_pressure_outlet) + print_state("PXR brine out", m.fs.PXR.brine_outlet) if __name__ == "__main__": diff --git a/watertap/examples/flowsheets/RO_with_energy_recovery/default_configuration.yaml b/watertap/examples/flowsheets/RO_with_energy_recovery/default_configuration.yaml index 7605398c03..e396c32005 100644 --- a/watertap/examples/flowsheets/RO_with_energy_recovery/default_configuration.yaml +++ b/watertap/examples/flowsheets/RO_with_energy_recovery/default_configuration.yaml @@ -61,32 +61,32 @@ fs.P2.costing.capital_cost: 610.4424522403153 fs.P2.ratioP[0.0]: 1.0725115194945298 fs.P2.work_fluid[0.0]: 255.95071372759554 fs.PXR.costing.capital_cost: 973.1020678190223 -fs.PXR.high_pressure_side.deltaP[0.0]: -7247342.029398318 -fs.PXR.high_pressure_side.properties_in[0.0].dens_mass_phase[Liq]: 1045.9125190112295 -fs.PXR.high_pressure_side.properties_in[0.0].flow_mass_phase_comp[Liq,H2O]: 0.49285444999999994 -fs.PXR.high_pressure_side.properties_in[0.0].flow_mass_phase_comp[Liq,NaCl]: 0.035587728093265734 -fs.PXR.high_pressure_side.properties_in[0.0].flow_vol_phase[Liq]: 0.0005052451027097726 -fs.PXR.high_pressure_side.properties_in[0.0].mass_frac_phase_comp[Liq,H2O]: 0.9326553981332943 -fs.PXR.high_pressure_side.properties_in[0.0].mass_frac_phase_comp[Liq,NaCl]: 0.06734460186670563 -fs.PXR.high_pressure_side.properties_in[0.0].pressure: 7348667.029398319 -fs.PXR.high_pressure_side.properties_in[0.0].temperature: 298.1500072638073 -fs.PXR.high_pressure_side.properties_out[0.0].flow_mass_phase_comp[Liq,H2O]: 0.49285444999999994 -fs.PXR.high_pressure_side.properties_out[0.0].flow_mass_phase_comp[Liq,NaCl]: 0.035587728093265734 -fs.PXR.high_pressure_side.properties_out[0.0].pressure: 101325.0 -fs.PXR.high_pressure_side.properties_out[0.0].temperature: 298.1500072638073 -fs.PXR.low_pressure_side.deltaP[0.0]: 6884974.927928402 -fs.PXR.low_pressure_side.properties_in[0.0].dens_mass_phase[Liq]: 1021.46 -fs.PXR.low_pressure_side.properties_in[0.0].flow_mass_phase_comp[Liq,H2O]: 0.498024594422437 -fs.PXR.low_pressure_side.properties_in[0.0].flow_mass_phase_comp[Liq,NaCl]: 0.018063068191487355 -fs.PXR.low_pressure_side.properties_in[0.0].flow_vol_phase[Liq]: 0.0005052451027097726 -fs.PXR.low_pressure_side.properties_in[0.0].mass_frac_phase_comp[Liq,H2O]: 0.965 -fs.PXR.low_pressure_side.properties_in[0.0].mass_frac_phase_comp[Liq,NaCl]: 0.03500000000000001 -fs.PXR.low_pressure_side.properties_in[0.0].pressure: 101325.0 -fs.PXR.low_pressure_side.properties_in[0.0].temperature: 298.15 -fs.PXR.low_pressure_side.properties_out[0.0].flow_mass_phase_comp[Liq,H2O]: 0.498024594422437 -fs.PXR.low_pressure_side.properties_out[0.0].flow_mass_phase_comp[Liq,NaCl]: 0.018063068191487355 -fs.PXR.low_pressure_side.properties_out[0.0].pressure: 6986299.927928403 -fs.PXR.low_pressure_side.properties_out[0.0].temperature: 298.15 +fs.PXR.brine_side.deltaP[0.0]: -7247342.029398318 +fs.PXR.brine_side.properties_in[0.0].dens_mass_phase[Liq]: 1045.9125190112295 +fs.PXR.brine_side.properties_in[0.0].flow_mass_phase_comp[Liq,H2O]: 0.49285444999999994 +fs.PXR.brine_side.properties_in[0.0].flow_mass_phase_comp[Liq,NaCl]: 0.035587728093265734 +fs.PXR.brine_side.properties_in[0.0].flow_vol_phase[Liq]: 0.0005052451027097726 +fs.PXR.brine_side.properties_in[0.0].mass_frac_phase_comp[Liq,H2O]: 0.9326553981332943 +fs.PXR.brine_side.properties_in[0.0].mass_frac_phase_comp[Liq,NaCl]: 0.06734460186670563 +fs.PXR.brine_side.properties_in[0.0].pressure: 7348667.029398319 +fs.PXR.brine_side.properties_in[0.0].temperature: 298.1500072638073 +fs.PXR.brine_side.properties_out[0.0].flow_mass_phase_comp[Liq,H2O]: 0.49285444999999994 +fs.PXR.brine_side.properties_out[0.0].flow_mass_phase_comp[Liq,NaCl]: 0.035587728093265734 +fs.PXR.brine_side.properties_out[0.0].pressure: 101325.0 +fs.PXR.brine_side.properties_out[0.0].temperature: 298.1500072638073 +fs.PXR.feed_side.deltaP[0.0]: 6884974.927928402 +fs.PXR.feed_side.properties_in[0.0].dens_mass_phase[Liq]: 1021.46 +fs.PXR.feed_side.properties_in[0.0].flow_mass_phase_comp[Liq,H2O]: 0.498024594422437 +fs.PXR.feed_side.properties_in[0.0].flow_mass_phase_comp[Liq,NaCl]: 0.018063068191487355 +fs.PXR.feed_side.properties_in[0.0].flow_vol_phase[Liq]: 0.0005052451027097726 +fs.PXR.feed_side.properties_in[0.0].mass_frac_phase_comp[Liq,H2O]: 0.965 +fs.PXR.feed_side.properties_in[0.0].mass_frac_phase_comp[Liq,NaCl]: 0.03500000000000001 +fs.PXR.feed_side.properties_in[0.0].pressure: 101325.0 +fs.PXR.feed_side.properties_in[0.0].temperature: 298.15 +fs.PXR.feed_side.properties_out[0.0].flow_mass_phase_comp[Liq,H2O]: 0.498024594422437 +fs.PXR.feed_side.properties_out[0.0].flow_mass_phase_comp[Liq,NaCl]: 0.018063068191487355 +fs.PXR.feed_side.properties_out[0.0].pressure: 6986299.927928403 +fs.PXR.feed_side.properties_out[0.0].temperature: 298.15 fs.RO.feed_side.K[0.0,0.0,NaCl]: 3.4512803512388215e-05 fs.RO.feed_side.K[0.0,1.0,NaCl]: 2.7060824531310153e-05 fs.RO.feed_side.N_Re[0.0,0.0]: 339.4894225951861 diff --git a/watertap/examples/flowsheets/RO_with_energy_recovery/tests/test_RO_with_energy_recovery_simulation.py b/watertap/examples/flowsheets/RO_with_energy_recovery/tests/test_RO_with_energy_recovery_simulation.py index c6a3bc4593..dc3948de27 100644 --- a/watertap/examples/flowsheets/RO_with_energy_recovery/tests/test_RO_with_energy_recovery_simulation.py +++ b/watertap/examples/flowsheets/RO_with_energy_recovery/tests/test_RO_with_energy_recovery_simulation.py @@ -123,10 +123,10 @@ def test_build(self, system_frame): fs.s03: (fs.P1.outlet, fs.M1.P1), fs.s04: (fs.M1.outlet, fs.RO.inlet), fs.s05: (fs.RO.permeate, fs.product.inlet), - fs.s06: (fs.RO.retentate, fs.PXR.high_pressure_inlet), - fs.s07: (fs.PXR.high_pressure_outlet, fs.disposal.inlet), - fs.s08: (fs.S1.PXR, fs.PXR.low_pressure_inlet), - fs.s09: (fs.PXR.low_pressure_outlet, fs.P2.inlet), + fs.s06: (fs.RO.retentate, fs.PXR.brine_inlet), + fs.s07: (fs.PXR.brine_outlet, fs.disposal.inlet), + fs.s08: (fs.S1.PXR, fs.PXR.feed_inlet), + fs.s09: (fs.PXR.feed_outlet, fs.P2.inlet), fs.s10: (fs.P2.outlet, fs.M1.P2), } for arc, port_tpl in arc_dict.items(): @@ -197,32 +197,30 @@ def test_initialize_system(self, system_frame): # check results across pressure exchanger, proxy for both upstream and downstream of RO # high pressure inlet assert value( - m.fs.PXR.high_pressure_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + m.fs.PXR.brine_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] ) == pytest.approx(0.4928, rel=1e-3) assert value( - m.fs.PXR.high_pressure_inlet.flow_mass_phase_comp[0, "Liq", "NaCl"] + m.fs.PXR.brine_inlet.flow_mass_phase_comp[0, "Liq", "NaCl"] ) == pytest.approx(3.561e-2, rel=1e-3) - assert value(m.fs.PXR.high_pressure_inlet.temperature[0]) == pytest.approx( + assert value(m.fs.PXR.brine_inlet.temperature[0]) == pytest.approx( 298.15, rel=1e-3 ) - assert value(m.fs.PXR.high_pressure_inlet.pressure[0]) == pytest.approx( + assert value(m.fs.PXR.brine_inlet.pressure[0]) == pytest.approx( 7.242e6, rel=1e-3 ) # low pressure inlet assert value( - m.fs.PXR.low_pressure_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + m.fs.PXR.feed_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] ) == pytest.approx(0.4980, rel=1e-3) assert value( - m.fs.PXR.low_pressure_inlet.flow_mass_phase_comp[0, "Liq", "NaCl"] + m.fs.PXR.feed_inlet.flow_mass_phase_comp[0, "Liq", "NaCl"] ) == pytest.approx(1.806e-2, rel=1e-3) - assert value(m.fs.PXR.low_pressure_inlet.temperature[0]) == pytest.approx( + assert value(m.fs.PXR.feed_inlet.temperature[0]) == pytest.approx( 298.15, rel=1e-3 ) - assert value(m.fs.PXR.low_pressure_inlet.pressure[0]) == pytest.approx( - 101325, rel=1e-3 - ) + assert value(m.fs.PXR.feed_inlet.pressure[0]) == pytest.approx(101325, rel=1e-3) # low pressure outlet - assert value(m.fs.PXR.low_pressure_outlet.pressure[0]) == pytest.approx( + assert value(m.fs.PXR.feed_outlet.pressure[0]) == pytest.approx( 6.885e6, rel=1e-3 ) @@ -310,12 +308,12 @@ def test_display_state(self, system_frame, capsys): Split 1 : 0.505 kg/s, 35000 ppm, 1.0 bar P1 out : 0.505 kg/s, 35000 ppm, 74.9 bar Split 2 : 0.516 kg/s, 35000 ppm, 1.0 bar -PXR LP out: 0.516 kg/s, 35000 ppm, 68.9 bar +PXR feed out: 0.516 kg/s, 35000 ppm, 68.9 bar P2 out : 0.516 kg/s, 35000 ppm, 74.9 bar Mix out : 1.021 kg/s, 35000 ppm, 74.9 bar RO perm : 0.493 kg/s, 240 ppm, 1.0 bar RO reten : 0.528 kg/s, 67424 ppm, 72.4 bar -PXR HP out: 0.528 kg/s, 67424 ppm, 1.0 bar +PXR brine out: 0.528 kg/s, 67424 ppm, 1.0 bar """ ) diff --git a/watertap/examples/flowsheets/case_studies/seawater_RO_desalination/seawater_RO_desalination.py b/watertap/examples/flowsheets/case_studies/seawater_RO_desalination/seawater_RO_desalination.py index 09f64bddfe..17837f3295 100644 --- a/watertap/examples/flowsheets/case_studies/seawater_RO_desalination/seawater_RO_desalination.py +++ b/watertap/examples/flowsheets/case_studies/seawater_RO_desalination/seawater_RO_desalination.py @@ -288,16 +288,12 @@ def eq_flow_mass_comp(blk, j): # pylint: disable=function-redefined desal.s01 = Arc(source=desal.S1.P1, destination=desal.P1.inlet) desal.s02 = Arc(source=desal.P1.outlet, destination=desal.M1.P1) desal.s03 = Arc(source=desal.M1.outlet, destination=desal.RO.inlet) - desal.s04 = Arc( - source=desal.RO.retentate, destination=desal.PXR.high_pressure_inlet - ) - desal.s05 = Arc(source=desal.S1.PXR, destination=desal.PXR.low_pressure_inlet) - desal.s06 = Arc( - source=desal.PXR.low_pressure_outlet, destination=desal.P2.inlet - ) + desal.s04 = Arc(source=desal.RO.retentate, destination=desal.PXR.brine_inlet) + desal.s05 = Arc(source=desal.S1.PXR, destination=desal.PXR.feed_inlet) + desal.s06 = Arc(source=desal.PXR.feed_outlet, destination=desal.P2.inlet) desal.s07 = Arc(source=desal.P2.outlet, destination=desal.M1.P2) m.fs.s_disposal = Arc( - source=desal.PXR.high_pressure_outlet, destination=m.fs.disposal.inlet + source=desal.PXR.brine_outlet, destination=m.fs.disposal.inlet ) elif erd_type == "pump_as_turbine": m.fs.s_tb_desal = Arc( @@ -343,8 +339,8 @@ def eq_flow_mass_comp(blk, j): # pylint: disable=function-redefined iscale.set_scaling_factor(desal.RO.area, 1e-4) if erd_type == "pressure_exchanger": iscale.set_scaling_factor(desal.P2.control_volume.work, 1e-5) - iscale.set_scaling_factor(desal.PXR.low_pressure_side.work, 1e-5) - iscale.set_scaling_factor(desal.PXR.high_pressure_side.work, 1e-5) + iscale.set_scaling_factor(desal.PXR.feed_side.work, 1e-5) + iscale.set_scaling_factor(desal.PXR.brine_side.work, 1e-5) elif erd_type == "pump_as_turbine": iscale.set_scaling_factor(desal.ERD.control_volume.work, 1e-5) # calculate and propagate scaling factors diff --git a/watertap/examples/flowsheets/case_studies/wastewater_resource_recovery/dye_desalination/dye_desalination_withRO.py b/watertap/examples/flowsheets/case_studies/wastewater_resource_recovery/dye_desalination/dye_desalination_withRO.py index dd6567c299..722ebaa326 100644 --- a/watertap/examples/flowsheets/case_studies/wastewater_resource_recovery/dye_desalination/dye_desalination_withRO.py +++ b/watertap/examples/flowsheets/case_studies/wastewater_resource_recovery/dye_desalination/dye_desalination_withRO.py @@ -314,15 +314,11 @@ def eq_flow_mass_comp(blk, j): desal.s01 = Arc(source=desal.S1.P2, destination=desal.P2.inlet) desal.s02 = Arc(source=desal.P2.outlet, destination=desal.M1.P2) desal.s03 = Arc(source=desal.M1.outlet, destination=desal.RO.inlet) - desal.s04 = Arc( - source=desal.RO.retentate, destination=desal.PXR.high_pressure_inlet - ) - desal.s05 = Arc(source=desal.S1.PXR, destination=desal.PXR.low_pressure_inlet) - desal.s06 = Arc(source=desal.PXR.low_pressure_outlet, destination=desal.P3.inlet) + desal.s04 = Arc(source=desal.RO.retentate, destination=desal.PXR.brine_inlet) + desal.s05 = Arc(source=desal.S1.PXR, destination=desal.PXR.feed_inlet) + desal.s06 = Arc(source=desal.PXR.feed_outlet, destination=desal.P3.inlet) desal.s07 = Arc(source=desal.P3.outlet, destination=desal.M1.P3) - m.fs.s_disposal = Arc( - source=desal.PXR.high_pressure_outlet, destination=m.fs.brine.inlet - ) + m.fs.s_disposal = Arc(source=desal.PXR.brine_outlet, destination=m.fs.brine.inlet) m.fs.s_permeate = Arc(source=desal.RO.permeate, destination=m.fs.permeate.inlet) @@ -341,8 +337,8 @@ def eq_flow_mass_comp(blk, j): ) iscale.set_scaling_factor(desal.RO.area, 1e-4) iscale.set_scaling_factor(desal.P3.control_volume.work, 1e-5) - iscale.set_scaling_factor(desal.PXR.low_pressure_side.work, 1e-5) - iscale.set_scaling_factor(desal.PXR.high_pressure_side.work, 1e-5) + iscale.set_scaling_factor(desal.PXR.feed_side.work, 1e-5) + iscale.set_scaling_factor(desal.PXR.brine_side.work, 1e-5) # calculate and propagate scaling factors iscale.calculate_scaling_factors(m) diff --git a/watertap/examples/flowsheets/case_studies/wastewater_resource_recovery/dye_desalination/tests/test_dye_desalination_withRO.py b/watertap/examples/flowsheets/case_studies/wastewater_resource_recovery/dye_desalination/tests/test_dye_desalination_withRO.py index 7b6384d606..841c0f04c1 100644 --- a/watertap/examples/flowsheets/case_studies/wastewater_resource_recovery/dye_desalination/tests/test_dye_desalination_withRO.py +++ b/watertap/examples/flowsheets/case_studies/wastewater_resource_recovery/dye_desalination/tests/test_dye_desalination_withRO.py @@ -380,7 +380,7 @@ def test_solve(self, system_frame): m.fs.treated.flow_mass_phase_comp[0, "Liq", "H2O"] ) - assert pytest.approx(0.000504058, rel=1e-3) == value( + assert pytest.approx(0.000485234, rel=1e-3) == value( m.fs.treated.flow_mass_phase_comp[0, "Liq", "dye"] ) @@ -392,7 +392,7 @@ def test_solve(self, system_frame): m.fs.adsorbed_dye.flow_mass_phase_comp[0, "Liq", "H2O"] ) - assert pytest.approx(0.0149224, rel=1e-3) == value( + assert pytest.approx(0.0149412, rel=1e-3) == value( m.fs.adsorbed_dye.flow_mass_phase_comp[0, "Liq", "dye"] ) diff --git a/watertap/unit_models/pressure_exchanger.py b/watertap/unit_models/pressure_exchanger.py index 4c2795fe07..98449ece12 100644 --- a/watertap/unit_models/pressure_exchanger.py +++ b/watertap/unit_models/pressure_exchanger.py @@ -10,6 +10,8 @@ # "https://github.com/watertap-org/watertap/" ################################################################################# +from enum import Enum, auto + # Import Pyomo libraries from pyomo.common.config import Bool, ConfigBlock, ConfigValue, In from pyomo.environ import ( @@ -42,6 +44,11 @@ _log = idaeslog.getLogger(__name__) +# --------------------------------------------------------------------- +class PressureExchangeType(Enum): + efficiency = auto() + high_pressure_difference = auto() + @declare_process_block_class("PressureExchanger") class PressureExchangerData(InitializationMixin, UnitModelBlockData): @@ -163,13 +170,37 @@ class PressureExchangerData(InitializationMixin, UnitModelBlockData): ) CONFIG.declare( - "has_mass_transfer", + "has_leakage", ConfigValue( default=False, domain=In([True, False]), - description="Defines if there is mass transport between high- and low-pressure sides", - doc="""Indicates whether pressure exchanger solution mass transfer terms should be constructed or not. - **default** - False.""", + description="Defines if there is leakage", + doc="""Indicates whether pressure exchanger has leakage. + **default** - False.""", + ), + ) + + CONFIG.declare( + "has_mixing", + ConfigValue( + default=False, + domain=In([True, False]), + description="Defines if there is mixing between high- and low-pressure side", + doc="""Indicates whether pressure exchanger has mixing. + **default** - False.""", + ), + ) + CONFIG.declare( + "pressure_exchange_calculation", + ConfigValue( + default=PressureExchangeType.efficiency, + domain=In(PressureExchangeType), + description="Pressure exchanger calculation method", + doc="""Indicates what type of pressure exchange calculation method should be used. + **default** - PressureExchangeType.efficiency. + **Valid values:** { + **PressureExchangeType.efficiency** - calculates momentum transfer by pressure exchanger efficiency, + **PressureExchangeType.high_pressure_difference** - calculates momentum transfer by high pressure difference}""", ), ) @@ -210,134 +241,203 @@ def build(self): doc="Pressure exchanger efficiency", ) - if self.config.has_mass_transfer: - self.mass_transfer_fraction_comp = Var( + if ( + self.config.pressure_exchange_calculation + is PressureExchangeType.high_pressure_difference + ): + self.high_pressure_difference = Var( self.flowsheet().config.time, - self.config.property_package.component_list, - initialize=0.05, + initialize=0.8e5, + bounds=(0, None), + domain=NonNegativeReals, + units=pyunits.Pa, + doc="High pressure difference", + ) + self.low_pressure_difference = Var( + self.flowsheet().config.time, + initialize=0.4e5, + bounds=(0, None), + domain=NonNegativeReals, + units=pyunits.Pa, + doc="Low pressure difference", + ) + + if self.config.has_leakage: + self.leakage_vol = Var( + self.flowsheet().config.time, + initialize=0.01, bounds=(1e-6, 1), domain=NonNegativeReals, units=pyunits.dimensionless, - doc="The fraction of solution transfering from high to low pressure side", + doc="The volumetric leakage fraction of the brine side to the feed side", ) - # Build control volume for high pressure side - self.high_pressure_side = ControlVolume0DBlock( + if self.config.has_mixing: + self.mixing_vol = Var( + self.flowsheet().config.time, + initialize=0.035, + bounds=(1e-6, 1), + domain=NonNegativeReals, + units=pyunits.dimensionless, + doc="The volumetric mixing fraction of the brine side and feed side", + ) + + # Build control volume for brine side + self.brine_side = ControlVolume0DBlock( dynamic=False, has_holdup=False, property_package=self.config.property_package, property_package_args=self.config.property_package_args, ) - self.high_pressure_side.add_state_blocks(has_phase_equilibrium=False) + self.brine_side.add_state_blocks(has_phase_equilibrium=False) - self.high_pressure_side.add_material_balances( + self.brine_side.add_material_balances( balance_type=self.config.material_balance_type, - has_mass_transfer=self.config.has_mass_transfer, + has_mass_transfer=self.config.has_mixing, ) - self.high_pressure_side.add_energy_balances( + self.brine_side.add_energy_balances( balance_type=self.config.energy_balance_type, has_enthalpy_transfer=False, ) if self.config.is_isothermal: - self.high_pressure_side.add_isothermal_assumption() + self.brine_side.add_isothermal_assumption() - self.high_pressure_side.add_momentum_balances( + self.brine_side.add_momentum_balances( balance_type=self.config.momentum_balance_type, has_pressure_change=True ) - @self.high_pressure_side.Expression( + @self.brine_side.Expression( self.flowsheet().config.time, - doc="Work transferred to high pressure side fluid (should be negative)", + doc="Work transferred to brine side fluid (should be negative)", ) def work(b, t): return b.properties_in[t].flow_vol * b.deltaP[t] - # Build control volume for low pressure side - self.low_pressure_side = ControlVolume0DBlock( + # Build control volume for feed side + self.feed_side = ControlVolume0DBlock( dynamic=False, has_holdup=False, property_package=self.config.property_package, property_package_args=self.config.property_package_args, ) - self.low_pressure_side.add_state_blocks(has_phase_equilibrium=False) + self.feed_side.add_state_blocks(has_phase_equilibrium=False) - self.low_pressure_side.add_material_balances( + self.feed_side.add_material_balances( balance_type=self.config.material_balance_type, - has_mass_transfer=self.config.has_mass_transfer, + has_mass_transfer=self.config.has_mixing, ) - self.low_pressure_side.add_energy_balances( + self.feed_side.add_energy_balances( balance_type=self.config.energy_balance_type, has_enthalpy_transfer=False, ) if self.config.is_isothermal: - self.low_pressure_side.add_isothermal_assumption() + self.feed_side.add_isothermal_assumption() - self.low_pressure_side.add_momentum_balances( + self.feed_side.add_momentum_balances( balance_type=self.config.momentum_balance_type, has_pressure_change=True ) - @self.low_pressure_side.Expression( + @self.feed_side.Expression( self.flowsheet().config.time, - doc="Work transferred to low pressure side fluid", + doc="Work transferred to feed side fluid", ) def work(b, t): # pylint: disable=function-redefined return b.properties_in[t].flow_vol * b.deltaP[t] # Add Ports - self.add_inlet_port(name="high_pressure_inlet", block=self.high_pressure_side) - self.add_outlet_port(name="high_pressure_outlet", block=self.high_pressure_side) - self.add_inlet_port(name="low_pressure_inlet", block=self.low_pressure_side) - self.add_outlet_port(name="low_pressure_outlet", block=self.low_pressure_side) + self.add_inlet_port(name="brine_inlet", block=self.brine_side) + self.add_outlet_port(name="brine_outlet", block=self.brine_side) + self.add_inlet_port(name="feed_inlet", block=self.feed_side) + self.add_outlet_port(name="feed_outlet", block=self.feed_side) # Performance equations @self.Constraint(self.flowsheet().config.time, doc="Pressure transfer") def eq_pressure_transfer(b, t): return ( - b.low_pressure_side.deltaP[t] - == b.efficiency_pressure_exchanger[t] * -b.high_pressure_side.deltaP[t] + b.feed_side.deltaP[t] + == b.efficiency_pressure_exchanger[t] * -b.brine_side.deltaP[t] ) - @self.Constraint(self.flowsheet().config.time, doc="Equal volumetric flow rate") - def eq_equal_flow_vol(b, t): + if ( + self.config.pressure_exchange_calculation + is PressureExchangeType.high_pressure_difference + ): - return ( - b.high_pressure_side.properties_out[t].flow_vol - == b.low_pressure_side.properties_in[t].flow_vol + @self.Constraint( + self.flowsheet().config.time, + doc="Pressure transfer by high pressure difference", ) + def eq_pressure_difference(b, t): + return b.brine_side.properties_in[ + t + ].pressure == b.feed_side.properties_out[t].pressure + pyunits.convert( + b.high_pressure_difference[t], to_units=pyunits.Pa + ) - @self.Constraint( - self.flowsheet().config.time, doc="Equal low pressure on both sides" - ) - def eq_equal_low_pressure(b, t): - return ( - b.high_pressure_side.properties_out[t].pressure - == b.low_pressure_side.properties_in[t].pressure + @self.Constraint( + self.flowsheet().config.time, doc="Equal low pressure on both sides" + ) + def eq_equal_low_pressure(b, t): + return b.brine_side.properties_out[ + t + ].pressure == b.feed_side.properties_in[t].pressure + pyunits.convert( + b.low_pressure_difference[t], to_units=pyunits.Pa + ) + + else: + + @self.Constraint( + self.flowsheet().config.time, doc="Equal low pressure on both sides" + ) + def eq_equal_low_pressure(b, t): + return ( + b.brine_side.properties_out[t].pressure + == b.feed_side.properties_in[t].pressure + ) + + if self.config.has_leakage: + + @self.Constraint( + self.flowsheet().config.time, doc="Equal volumetric flow rate" + ) + def eq_equal_flow_vol(b, t): + return ( + b.feed_side.properties_out[t].flow_vol + == (1 - b.leakage_vol[t]) * b.brine_side.properties_in[t].flow_vol + ) + + else: + + @self.Constraint( + self.flowsheet().config.time, doc="Equal volumetric flow rate" ) + def eq_equal_flow_vol(b, t): + return ( + b.feed_side.properties_out[t].flow_vol + == b.brine_side.properties_in[t].flow_vol + ) - if self.config.has_mass_transfer: + if self.config.has_mixing: @self.Constraint( self.flowsheet().config.time, self.config.property_package.phase_list, - self.config.property_package.component_list, - doc="Mass transfer from high pressure side", + self.config.property_package.solute_set, + doc="Mixing effect of the unit", ) - def eq_mass_transfer_from_high_pressure_side(b, t, p, j): - comp = self.config.property_package.get_component(j) - return b.high_pressure_side.mass_transfer_term[ - t, p, j - ] == -b.mass_transfer_fraction_comp[ - t, j - ] * b.high_pressure_side.properties_in[ - t - ].get_material_flow_terms( - p, j + def eq_mixing(b, t, p, j): + return ( + b.feed_side.properties_out[t].conc_mass_phase_comp[p, j] + == b.feed_side.properties_in[t].conc_mass_phase_comp[p, j] + * (1 - b.mixing_vol[t]) + + b.brine_side.properties_in[t].conc_mass_phase_comp[p, j] + * b.mixing_vol[t] ) @self.Constraint( @@ -348,8 +448,19 @@ def eq_mass_transfer_from_high_pressure_side(b, t, p, j): ) def eq_mass_transfer_term(b, t, p, j): return ( - b.high_pressure_side.mass_transfer_term[t, p, j] - == -b.low_pressure_side.mass_transfer_term[t, p, j] + b.brine_side.mass_transfer_term[t, p, j] + == -b.feed_side.mass_transfer_term[t, p, j] + ) + + @self.Constraint( + self.flowsheet().config.time, + doc="Equal volumetric flow rate on low-pressure side", + ) + def eq_equal_LPS_flow_vol(b, t): + + return ( + b.feed_side.properties_out[t].flow_vol + == b.feed_side.properties_in[t].flow_vol ) def initialize_build( @@ -385,14 +496,14 @@ def initialize_build( # Set solver and options opt = get_solver(solver, optarg) # initialize inlets - flags_low_in = self.low_pressure_side.properties_in.initialize( + flags_low_in = self.feed_side.properties_in.initialize( outlvl=outlvl, optarg=optarg, solver=solver, state_args=state_args, hold_state=True, ) - flags_high_in = self.high_pressure_side.properties_in.initialize( + flags_high_in = self.brine_side.properties_in.initialize( outlvl=outlvl, optarg=optarg, solver=solver, @@ -403,31 +514,32 @@ def initialize_build( init_log.info_high("Initialize inlets complete") # check that inlets are feasible - if value(self.low_pressure_side.properties_in[0].pressure) > value( - self.high_pressure_side.properties_in[0].pressure + if value(self.feed_side.properties_in[0].pressure) > value( + self.brine_side.properties_in[0].pressure ): raise ConfigurationError( "Initializing pressure exchanger failed because " - "the low pressure side inlet has a higher pressure " - "than the high pressure side inlet" + "the feed side inlet has a higher pressure " + "than the brine side inlet" ) # only needed when there is no mass trnasfer if ( abs( - value(self.low_pressure_side.properties_in[0].flow_vol) - - value(self.high_pressure_side.properties_in[0].flow_vol) + value(self.feed_side.properties_in[0].flow_vol) + - value(self.brine_side.properties_in[0].flow_vol) ) - / value(self.high_pressure_side.properties_in[0].flow_vol) + / value(self.brine_side.properties_in[0].flow_vol) > 1e-4 - and not self.config.has_mass_transfer + and not self.config.has_mixing + and not self.config.has_leakage ): # flow_vol values are not within 0.1% raise ConfigurationError( "Initializing pressure exchanger failed because " "the volumetric flow rates are not equal for both inlets " - + str(value(self.high_pressure_side.properties_out[0].flow_vol)) + + str(value(self.brine_side.properties_out[0].flow_vol)) + "," - + str(value(self.low_pressure_side.properties_in[0].flow_vol)) + + str(value(self.feed_side.properties_in[0].flow_vol)) ) else: # volumetric flow is equal, deactivate flow constraint for the solve self.eq_equal_flow_vol.deactivate() @@ -443,29 +555,35 @@ def propogate_state(sb1, sb2): else: state_dict_2[k].value = state_dict_1[k].value - # low pressure side + # feed side propogate_state( - self.low_pressure_side.properties_in[0], - self.low_pressure_side.properties_out[0], + self.feed_side.properties_in[0], + self.feed_side.properties_out[0], ) - self.low_pressure_side.properties_out[ - 0 - ].pressure = self.low_pressure_side.properties_in[ - 0 - ].pressure.value + self.efficiency_pressure_exchanger[ - 0 - ].value * ( - self.high_pressure_side.properties_in[0].pressure.value - - self.low_pressure_side.properties_in[0].pressure.value - ) - # high pressure side + if self.config.pressure_exchange_calculation is PressureExchangeType.efficiency: + self.feed_side.properties_out[0].pressure = self.feed_side.properties_in[ + 0 + ].pressure.value + self.efficiency_pressure_exchanger[0].value * ( + self.brine_side.properties_in[0].pressure.value + - self.feed_side.properties_in[0].pressure.value + ) + elif ( + self.config.pressure_exchange_calculation + is PressureExchangeType.high_pressure_difference + ): + self.feed_side.properties_out[0].pressure = ( + self.brine_side.properties_in[0].pressure.value + - self.high_pressure_difference[0].value + ) + + # brine side propogate_state( - self.high_pressure_side.properties_in[0], - self.high_pressure_side.properties_out[0], + self.brine_side.properties_in[0], + self.brine_side.properties_out[0], ) - self.high_pressure_side.properties_out[ + self.brine_side.properties_out[0].pressure.value = self.feed_side.properties_in[ 0 - ].pressure.value = self.low_pressure_side.properties_in[0].pressure.value + ].pressure.value init_log.info_high("Initialize outlets complete") # Solve unit @@ -474,8 +592,8 @@ def propogate_state(sb1, sb2): init_log.info("Initialization complete: {}".format(idaeslog.condition(res))) # release state of fixed variables - self.low_pressure_side.properties_in.release_state(flags_low_in) - self.high_pressure_side.properties_in.release_state(flags_high_in) + self.feed_side.properties_in.release_state(flags_low_in) + self.brine_side.properties_in.release_state(flags_high_in) if not check_optimal_termination(res): raise InitializationError(f"Unit model {self.name} failed to initialize") @@ -490,69 +608,80 @@ def calculate_scaling_factors(self): if iscale.get_scaling_factor(self.efficiency_pressure_exchanger) is None: # efficiency should always be between 0.1-1 iscale.set_scaling_factor(self.efficiency_pressure_exchanger, 1) - if hasattr(self, "mass_transfer_fraction_comp"): - if iscale.get_scaling_factor(self.mass_transfer_fraction_comp) is None: - iscale.set_scaling_factor(self.mass_transfer_fraction_comp, 1) + if hasattr(self, "high_pressure_difference"): + if iscale.get_scaling_factor(self.high_pressure_difference) is None: + iscale.set_scaling_factor(self.high_pressure_difference, 1e-5) + if hasattr(self, "low_pressure_difference"): + if iscale.get_scaling_factor(self.low_pressure_difference) is None: + iscale.set_scaling_factor(self.low_pressure_difference, 1e-5) + if hasattr(self, "leakage_vol"): + if iscale.get_scaling_factor(self.leakage_vol) is None: + iscale.set_scaling_factor(self.leakage_vol, 1) + if hasattr(self, "mixing_vol"): + if iscale.get_scaling_factor(self.mixing_vol) is None: + iscale.set_scaling_factor(self.mixing_vol, 1) # scale expressions - if iscale.get_scaling_factor(self.low_pressure_side.work) is None: - sf = iscale.get_scaling_factor( - self.low_pressure_side.properties_in[0].flow_vol - ) - sf = sf * iscale.get_scaling_factor(self.low_pressure_side.deltaP[0]) - iscale.set_scaling_factor(self.low_pressure_side.work, sf) + if iscale.get_scaling_factor(self.feed_side.work) is None: + sf = iscale.get_scaling_factor(self.feed_side.properties_in[0].flow_vol) + sf = sf * iscale.get_scaling_factor(self.feed_side.deltaP[0]) + iscale.set_scaling_factor(self.feed_side.work, sf) - if iscale.get_scaling_factor(self.high_pressure_side.work) is None: - sf = iscale.get_scaling_factor( - self.high_pressure_side.properties_in[0].flow_vol - ) - sf = sf * iscale.get_scaling_factor(self.high_pressure_side.deltaP[0]) - iscale.set_scaling_factor(self.high_pressure_side.work, sf) + if iscale.get_scaling_factor(self.brine_side.work) is None: + sf = iscale.get_scaling_factor(self.brine_side.properties_in[0].flow_vol) + sf = sf * iscale.get_scaling_factor(self.brine_side.deltaP[0]) + iscale.set_scaling_factor(self.brine_side.work, sf) # transform constraints for t, c in self.eq_pressure_transfer.items(): - sf = iscale.get_scaling_factor(self.low_pressure_side.deltaP[t]) + sf = iscale.get_scaling_factor(self.efficiency_pressure_exchanger[t]) iscale.constraint_scaling_transform(c, sf) for t, c in self.eq_equal_flow_vol.items(): - sf = iscale.get_scaling_factor( - self.low_pressure_side.properties_in[t].flow_vol - ) + sf = iscale.get_scaling_factor(self.brine_side.properties_in[t].flow_vol) iscale.constraint_scaling_transform(c, sf) for t, c in self.eq_equal_low_pressure.items(): - sf = iscale.get_scaling_factor( - self.low_pressure_side.properties_in[t].pressure - ) + sf = iscale.get_scaling_factor(self.feed_side.properties_in[t].pressure) iscale.constraint_scaling_transform(c, sf) - if hasattr(self, "eq_mass_transfer_from_high_pressure_side"): - for (t, p, j), c in self.eq_mass_transfer_from_high_pressure_side.items(): + if hasattr(self, "eq_equal_LPS_flow_vol"): + for t, c in self.eq_equal_LPS_flow_vol.items(): + sf = iscale.get_scaling_factor(self.feed_side.properties_in[t].flow_vol) + iscale.constraint_scaling_transform(c, sf) + + if hasattr(self, "eq_mixing"): + for (t, p, j), c in self.eq_mixing.items(): sf = iscale.get_scaling_factor( - self.high_pressure_side.properties_in[t].get_material_flow_terms( - p, j - ) + self.feed_side.properties_in[t].conc_mass_phase_comp[p, j] ) iscale.constraint_scaling_transform(c, sf) if hasattr(self, "eq_mass_transfer_term"): for (t, p, j), c in self.eq_mass_transfer_term.items(): sf = iscale.get_scaling_factor( - self.high_pressure_side.mass_transfer_term[t, p, j] + self.brine_side.mass_transfer_term[t, p, j] ) iscale.constraint_scaling_transform(c, sf) sf = iscale.get_scaling_factor( - self.low_pressure_side.mass_transfer_term[t, p, j] + self.feed_side.mass_transfer_term[t, p, j] + ) + iscale.constraint_scaling_transform(c, sf) + + if hasattr(self, "eq_pressure_difference"): + for t, c in self.eq_pressure_difference.items(): + sf = iscale.get_scaling_factor( + self.brine_side.properties_in[t].pressure ) iscale.constraint_scaling_transform(c, sf) def _get_stream_table_contents(self, time_point=0): return create_stream_table_dataframe( { - "HP Side In": self.high_pressure_inlet, - "HP Side Out": self.high_pressure_outlet, - "LP Side In": self.low_pressure_inlet, - "LP Side Out": self.low_pressure_outlet, + "Brine Side In": self.brine_inlet, + "Brine Side Out": self.brine_outlet, + "Feed Side In": self.feed_inlet, + "Feed Side Out": self.feed_outlet, }, time_point=time_point, ) @@ -562,12 +691,12 @@ def _get_performance_contents(self, time_point=0): return { "vars": { "Efficiency": self.efficiency_pressure_exchanger[t], - "HP Side Pressure Change": self.high_pressure_side.deltaP[t], - "LP Side Pressure Change": self.low_pressure_side.deltaP[t], + "Brine Side Pressure Change": self.brine_side.deltaP[t], + "Feed Side Pressure Change": self.feed_side.deltaP[t], }, "exprs": { - "HP Side Mechanical Work": self.high_pressure_side.work[t], - "LP Side Mechanical Work": self.low_pressure_side.work[t], + "Brine Side Mechanical Work": self.brine_side.work[t], + "Feed Side Mechanical Work": self.feed_side.work[t], }, "params": {}, } diff --git a/watertap/unit_models/tests/test_pressure_exchanger.py b/watertap/unit_models/tests/test_pressure_exchanger.py index 2802ee0283..237fc2f623 100644 --- a/watertap/unit_models/tests/test_pressure_exchanger.py +++ b/watertap/unit_models/tests/test_pressure_exchanger.py @@ -18,6 +18,7 @@ TerminationCondition, SolverStatus, value, + units as pyunits, ) from pyomo.network import Port from idaes.core import ( @@ -26,10 +27,12 @@ EnergyBalanceType, MomentumBalanceType, ) -from watertap.unit_models.pressure_exchanger import PressureExchanger +from watertap.unit_models.pressure_exchanger import ( + PressureExchanger, + PressureExchangeType, +) import watertap.property_models.seawater_prop_pack as props -import watertap.property_models.seawater_ion_prop_pack as property_seawater_ions - +import watertap.property_models.multicomp_aq_sol_prop_pack as props_multi from idaes.core.util.model_statistics import ( degrees_of_freedom, number_variables, @@ -61,7 +64,7 @@ def test_config_no_mass_transfer(): m.fs.unit = PressureExchanger(property_package=m.fs.properties) # check unit config arguments - assert len(m.fs.unit.config) == 9 + assert len(m.fs.unit.config) == 11 assert not m.fs.unit.config.dynamic assert not m.fs.unit.config.has_holdup @@ -71,7 +74,7 @@ def test_config_no_mass_transfer(): assert m.fs.unit.config.property_package is m.fs.properties # verify no mass transfer - assert not m.fs.unit.config.has_mass_transfer + assert not m.fs.unit.config.has_mixing @pytest.mark.unit @@ -79,28 +82,28 @@ def test_config_mass_transfer(): m = ConcreteModel() m.fs = FlowsheetBlock(dynamic=False) m.fs.properties = props.SeawaterParameterBlock() - m.fs.unit = PressureExchanger( - property_package=m.fs.properties, has_mass_transfer=True - ) + m.fs.unit = PressureExchanger(property_package=m.fs.properties, has_mixing=True) # check mass_transfer is added - assert m.fs.unit.config.has_mass_transfer + assert m.fs.unit.config.has_mixing @pytest.mark.unit -def test_build(has_mass_transfer=False, extra_variables=0, extra_constraint=0): +def test_build( + has_leakage=False, has_mixing=False, extra_variables=0, extra_constraint=0 +): m = ConcreteModel() m.fs = FlowsheetBlock(dynamic=False) m.fs.properties = props.SeawaterParameterBlock() m.fs.unit = PressureExchanger( - property_package=m.fs.properties, has_mass_transfer=has_mass_transfer + property_package=m.fs.properties, has_leakage=has_leakage, has_mixing=has_mixing ) # test ports and state variables port_lst = [ - "low_pressure_inlet", - "low_pressure_outlet", - "high_pressure_inlet", - "high_pressure_outlet", + "feed_inlet", + "feed_outlet", + "brine_inlet", + "brine_outlet", ] port_vars_lst = ["flow_mass_phase_comp", "pressure", "temperature"] @@ -117,10 +120,15 @@ def test_build(has_mass_transfer=False, extra_variables=0, extra_constraint=0): # test unit variables assert hasattr(m.fs.unit, "efficiency_pressure_exchanger") assert isinstance(m.fs.unit.efficiency_pressure_exchanger, Var) - if not has_mass_transfer: - assert not hasattr(m.fs.unit, "mass_transfer_fraction_comp") + if not has_leakage: + assert not hasattr(m.fs.unit, "leakage_vol") else: - assert isinstance(m.fs.unit.mass_transfer_fraction_comp, Var) + assert isinstance(m.fs.unit.leakage_vol, Var) + + if not has_mixing: + assert not hasattr(m.fs.unit, "mixing_vol") + else: + assert isinstance(m.fs.unit.mixing_vol, Var) # test unit constraints @@ -130,8 +138,10 @@ def test_build(has_mass_transfer=False, extra_variables=0, extra_constraint=0): "eq_equal_low_pressure", ] - if has_mass_transfer: - unit_cons_lst.append("eq_mass_transfer_from_high_pressure_side") + if has_leakage or has_mixing: + # unit_cons_lst.append("eq_mass_transfer_from_brine_side") + unit_cons_lst.append("eq_leakage") + unit_cons_lst.append("eq_mixing") unit_cons_lst.append("eq_mass_transfer_term") for c in unit_cons_lst: @@ -141,9 +151,9 @@ def test_build(has_mass_transfer=False, extra_variables=0, extra_constraint=0): # test control volumes, only terms directly used by pressure exchanger - cv_list = ["low_pressure_side", "high_pressure_side"] + cv_list = ["feed_side", "brine_side"] cv_var_lst = ["deltaP"] - if has_mass_transfer: + if has_mixing: cv_var_lst.append("mass_transfer_term") cv_exp_lst = ["work"] @@ -182,22 +192,17 @@ def test_build(has_mass_transfer=False, extra_variables=0, extra_constraint=0): @pytest.mark.unit -def test_build_with_mass_transfer(): - # 2 variables for mass transfer fractions (solute,solvent) - # 4 extra vars for mass_transfer_term on LP and HP side, and solute, solvent - # 4 extra constraints for mass transsfer from HP side (soute, solvent) - test_build(has_mass_transfer=True, extra_variables=6, extra_constraint=4) +def test_build_without_mass_transfer(): + test_build(has_leakage=False, has_mixing=False) -class TestPressureExchanger: +class TestPressureExchanger_without_mass_transfer: @pytest.fixture(scope="class") def unit_frame(self): m = ConcreteModel() m.fs = FlowsheetBlock(dynamic=False) m.fs.properties = props.SeawaterParameterBlock() - m.fs.unit = PressureExchanger( - property_package=m.fs.properties, has_mass_transfer=True - ) + m.fs.unit = PressureExchanger(property_package=m.fs.properties) # Specify inlet conditions temperature = 25 + 273.15 @@ -207,52 +212,35 @@ def unit_frame(self): highP_mass_frac_TDS = 0.07 highP_pressure = 50e5 - solute_transfer = 0.05 - solvent_transfer = 0.1 - m.fs.unit.mass_transfer_fraction_comp[0, "H2O"].fix(solvent_transfer) - m.fs.unit.mass_transfer_fraction_comp[0, "TDS"].fix(solute_transfer) - - m.fs.unit.low_pressure_side.properties_in[0].flow_vol_phase["Liq"].fix(flow_vol) - m.fs.unit.low_pressure_side.properties_in[0].mass_frac_phase_comp[ - "Liq", "TDS" - ].fix(lowP_mass_frac_TDS) + m.fs.unit.feed_side.properties_in[0].flow_vol_phase["Liq"].fix(flow_vol) + m.fs.unit.feed_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].fix( + lowP_mass_frac_TDS + ) - m.fs.unit.low_pressure_side.properties_in[0].pressure.fix(lowP_pressure) - m.fs.unit.low_pressure_side.properties_in[0].temperature.fix(temperature) + m.fs.unit.feed_side.properties_in[0].pressure.fix(lowP_pressure) + m.fs.unit.feed_side.properties_in[0].temperature.fix(temperature) - m.fs.unit.high_pressure_side.properties_in[0].flow_vol_phase["Liq"].fix( - flow_vol + m.fs.unit.brine_side.properties_in[0].flow_vol_phase["Liq"].fix(flow_vol) + m.fs.unit.brine_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].fix( + highP_mass_frac_TDS ) - m.fs.unit.high_pressure_side.properties_in[0].mass_frac_phase_comp[ - "Liq", "TDS" - ].fix(highP_mass_frac_TDS) - m.fs.unit.high_pressure_side.properties_in[0].pressure.fix(highP_pressure) - m.fs.unit.high_pressure_side.properties_in[0].temperature.fix(temperature) + m.fs.unit.brine_side.properties_in[0].pressure.fix(highP_pressure) + m.fs.unit.brine_side.properties_in[0].temperature.fix(temperature) # solve inlet conditions and only fix state variables (i.e. unfix flow_vol and mass_frac_phase) - results = solver.solve(m.fs.unit.high_pressure_side.properties_in[0]) + results = solver.solve(m.fs.unit.brine_side.properties_in[0]) assert results.solver.termination_condition == TerminationCondition.optimal - m.fs.unit.high_pressure_side.properties_in[0].flow_mass_phase_comp[ - "Liq", "TDS" - ].fix() - m.fs.unit.high_pressure_side.properties_in[0].flow_vol_phase["Liq"].unfix() - m.fs.unit.high_pressure_side.properties_in[0].mass_frac_phase_comp[ - "Liq", "TDS" - ].unfix() - - results = solver.solve(m.fs.unit.low_pressure_side.properties_in[0]) + m.fs.unit.brine_side.properties_in[0].flow_mass_phase_comp["Liq", "TDS"].fix() + m.fs.unit.brine_side.properties_in[0].flow_vol_phase["Liq"].unfix() + m.fs.unit.brine_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].unfix() + + results = solver.solve(m.fs.unit.feed_side.properties_in[0]) assert results.solver.termination_condition == TerminationCondition.optimal - m.fs.unit.low_pressure_side.properties_in[0].flow_mass_phase_comp[ - "Liq", "H2O" - ].fix() - m.fs.unit.low_pressure_side.properties_in[0].flow_mass_phase_comp[ - "Liq", "TDS" - ].fix() - m.fs.unit.low_pressure_side.properties_in[0].flow_vol_phase["Liq"].unfix() - m.fs.unit.low_pressure_side.properties_in[0].mass_frac_phase_comp[ - "Liq", "TDS" - ].unfix() + m.fs.unit.feed_side.properties_in[0].flow_mass_phase_comp["Liq", "H2O"].fix() + m.fs.unit.feed_side.properties_in[0].flow_mass_phase_comp["Liq", "TDS"].fix() + m.fs.unit.feed_side.properties_in[0].flow_vol_phase["Liq"].unfix() + m.fs.unit.feed_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].unfix() # Specify unit efficiency = 0.95 @@ -277,6 +265,7 @@ def test_calculate_scaling(self, unit_frame): unscaled_constraint_list = list(unscaled_constraints_generator(m)) assert len(unscaled_constraint_list) == 0 + @pytest.mark.requires_idaes_solver @pytest.mark.component def test_initialize(self, unit_frame): m = unit_frame @@ -288,6 +277,7 @@ def test_var_scaling(self, unit_frame): badly_scaled_var_lst = list(badly_scaled_var_generator(m)) assert badly_scaled_var_lst == [] + @pytest.mark.requires_idaes_solver @pytest.mark.component def test_solve(self, unit_frame): m = unit_frame @@ -297,56 +287,333 @@ def test_solve(self, unit_frame): assert results.solver.termination_condition == TerminationCondition.optimal assert results.solver.status == SolverStatus.ok + @pytest.mark.requires_idaes_solver @pytest.mark.component def test_solution(self, unit_frame): m = unit_frame assert pytest.approx(0.9877, rel=1e-3) == value( - m.fs.unit.low_pressure_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + m.fs.unit.feed_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] ) assert pytest.approx(3.582e-2, rel=1e-3) == value( - m.fs.unit.low_pressure_inlet.flow_mass_phase_comp[0, "Liq", "TDS"] + m.fs.unit.feed_inlet.flow_mass_phase_comp[0, "Liq", "TDS"] ) assert pytest.approx(4.755e6, rel=1e-3) == value( - m.fs.unit.low_pressure_outlet.pressure[0] + m.fs.unit.feed_outlet.pressure[0] ) - assert pytest.approx(0.9877 * 1.1, rel=1e-3) == value( - m.fs.unit.high_pressure_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + assert pytest.approx(0.9767, rel=1e-3) == value( + m.fs.unit.brine_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] ) assert pytest.approx(7.352e-2, rel=1e-3) == value( - m.fs.unit.high_pressure_inlet.flow_mass_phase_comp[0, "Liq", "TDS"] + m.fs.unit.brine_inlet.flow_mass_phase_comp[0, "Liq", "TDS"] + ) + assert pytest.approx(4.654e6, rel=1e-3) == value(m.fs.unit.feed_side.deltaP[0]) + assert pytest.approx(4.654e3, rel=1e-3) == value(m.fs.unit.feed_side.work[0]) + assert pytest.approx(-4.899e6, rel=1e-3) == value( + m.fs.unit.brine_side.deltaP[0] + ) + assert pytest.approx(-4.899e3, rel=1e-3) == value(m.fs.unit.brine_side.work[0]) + + # testing solvent transfer + assert pytest.approx(0.9767, rel=1e-3) == value( + m.fs.unit.brine_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] + ) + assert pytest.approx(0.9877, rel=1e-3) == value( + m.fs.unit.feed_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] + ) + # testing solute transfer + assert pytest.approx(7.352e-2, rel=1e-3) == value( + m.fs.unit.brine_outlet.flow_mass_phase_comp[0, "Liq", "TDS"] + ) + assert pytest.approx(3.582e-2, rel=1e-3) == value( + m.fs.unit.feed_outlet.flow_mass_phase_comp[0, "Liq", "TDS"] + ) + + @pytest.mark.unit + def test_report(self, unit_frame): + unit_frame.fs.unit.report() + + +class TestPressureExchanger_with_high_pressure_difference: + @pytest.fixture(scope="class") + def unit_frame(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = props.SeawaterParameterBlock() + m.fs.unit = PressureExchanger( + property_package=m.fs.properties, + pressure_exchange_calculation=PressureExchangeType.high_pressure_difference, ) - assert pytest.approx(4.654e6, rel=1e-3) == value( - m.fs.unit.low_pressure_side.deltaP[0] + + # Specify inlet conditions + temperature = 25 + 273.15 + flow_vol = 1e-3 + lowP_mass_frac_TDS = 0.035 + lowP_pressure = 101325 + highP_mass_frac_TDS = 0.07 + highP_pressure = 50e5 + + m.fs.unit.feed_side.properties_in[0].flow_vol_phase["Liq"].fix(flow_vol) + m.fs.unit.feed_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].fix( + lowP_mass_frac_TDS ) - assert pytest.approx(4.654e3, rel=1e-3) == value( - m.fs.unit.low_pressure_side.work[0] + + m.fs.unit.feed_side.properties_in[0].pressure.fix(lowP_pressure) + m.fs.unit.feed_side.properties_in[0].temperature.fix(temperature) + + m.fs.unit.brine_side.properties_in[0].flow_vol_phase["Liq"].fix(flow_vol) + m.fs.unit.brine_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].fix( + highP_mass_frac_TDS ) + m.fs.unit.brine_side.properties_in[0].pressure.fix(highP_pressure) + m.fs.unit.brine_side.properties_in[0].temperature.fix(temperature) + + # solve inlet conditions and only fix state variables (i.e. unfix flow_vol and mass_frac_phase) + + results = solver.solve(m.fs.unit.brine_side.properties_in[0]) + assert results.solver.termination_condition == TerminationCondition.optimal + m.fs.unit.brine_side.properties_in[0].flow_mass_phase_comp["Liq", "TDS"].fix() + m.fs.unit.brine_side.properties_in[0].flow_vol_phase["Liq"].unfix() + m.fs.unit.brine_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].unfix() + + results = solver.solve(m.fs.unit.feed_side.properties_in[0]) + assert results.solver.termination_condition == TerminationCondition.optimal + m.fs.unit.feed_side.properties_in[0].flow_mass_phase_comp["Liq", "H2O"].fix() + m.fs.unit.feed_side.properties_in[0].flow_mass_phase_comp["Liq", "TDS"].fix() + m.fs.unit.feed_side.properties_in[0].flow_vol_phase["Liq"].unfix() + m.fs.unit.feed_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].unfix() + + # Specify unit + high_pressure_difference = 0.8 * pyunits.bar + low_pressure_difference = 0 * pyunits.bar + m.fs.unit.high_pressure_difference.fix(high_pressure_difference) + m.fs.unit.low_pressure_difference.fix(low_pressure_difference) + return m + + @pytest.mark.unit + def test_dof(self, unit_frame): + assert degrees_of_freedom(unit_frame) == 0 + + @pytest.mark.unit + def test_calculate_scaling(self, unit_frame): + m = unit_frame + calculate_scaling_factors(m) + + # check that all variables have scaling factors + unscaled_var_list = list( + unscaled_variables_generator(m.fs.unit, include_fixed=True) + ) + assert len(unscaled_var_list) == 0 + # check that all constraints have been scaled + unscaled_constraint_list = list(unscaled_constraints_generator(m)) + assert len(unscaled_constraint_list) == 0 + + @pytest.mark.component + def test_initialize(self, unit_frame): + m = unit_frame + initialization_tester(unit_frame) + + @pytest.mark.component + def test_var_scaling(self, unit_frame): + m = unit_frame + badly_scaled_var_lst = list(badly_scaled_var_generator(m)) + assert badly_scaled_var_lst == [] + + @pytest.mark.component + def test_solve(self, unit_frame): + m = unit_frame + results = solver.solve(m) + + # Check for optimal solution + assert results.solver.termination_condition == TerminationCondition.optimal + assert results.solver.status == SolverStatus.ok + + @pytest.mark.component + def test_solution(self, unit_frame): + m = unit_frame + assert pytest.approx(0.9877, rel=1e-3) == value( + m.fs.unit.feed_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + ) + assert pytest.approx(3.582e-2, rel=1e-3) == value( + m.fs.unit.feed_inlet.flow_mass_phase_comp[0, "Liq", "TDS"] + ) + assert pytest.approx(4.92e6, rel=1e-3) == value( + m.fs.unit.feed_outlet.pressure[0] + ) + assert pytest.approx(0.9767, rel=1e-3) == value( + m.fs.unit.brine_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + ) + assert pytest.approx(7.352e-2, rel=1e-3) == value( + m.fs.unit.brine_inlet.flow_mass_phase_comp[0, "Liq", "TDS"] + ) + assert pytest.approx(4.819e6, rel=1e-3) == value(m.fs.unit.feed_side.deltaP[0]) + assert pytest.approx(4.819e3, rel=1e-3) == value(m.fs.unit.feed_side.work[0]) assert pytest.approx(-4.899e6, rel=1e-3) == value( - m.fs.unit.high_pressure_side.deltaP[0] + m.fs.unit.brine_side.deltaP[0] + ) + assert pytest.approx(-4.899e3, rel=1e-3) == value(m.fs.unit.brine_side.work[0]) + + # testing solvent transfer + assert pytest.approx(0.9767, rel=1e-3) == value( + m.fs.unit.brine_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] + ) + assert pytest.approx(0.9877, rel=1e-3) == value( + m.fs.unit.feed_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] + ) + # testing solute transfer + assert pytest.approx(7.352e-2, rel=1e-3) == value( + m.fs.unit.brine_outlet.flow_mass_phase_comp[0, "Liq", "TDS"] ) - assert pytest.approx(-5.436e3, rel=1e-3) == value( - m.fs.unit.high_pressure_side.work[0] + assert pytest.approx(3.582e-2, rel=1e-3) == value( + m.fs.unit.feed_outlet.flow_mass_phase_comp[0, "Liq", "TDS"] ) + @pytest.mark.unit + def test_report(self, unit_frame): + unit_frame.fs.unit.report() + + +class TestPressureExchanger_with_mass_transfer: + @pytest.fixture(scope="class") + def unit_frame(self): + m = ConcreteModel() + m.fs = FlowsheetBlock(dynamic=False) + m.fs.properties = props.SeawaterParameterBlock() + m.fs.unit = PressureExchanger( + property_package=m.fs.properties, + has_leakage=True, + has_mixing=True, + pressure_exchange_calculation=PressureExchangeType.high_pressure_difference, + ) + + # Specify inlet conditions + temperature = 25 + 273.15 + flow_vol = 1e-3 + lowP_mass_frac_TDS = 0.035 + lowP_pressure = 101325 + highP_mass_frac_TDS = 0.07 + highP_pressure = 50e5 + leakage_vol = 0.01 + mixing_vol = 0.035 + + m.fs.unit.feed_side.properties_in[0].flow_vol_phase["Liq"].fix(flow_vol) + m.fs.unit.feed_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].fix( + lowP_mass_frac_TDS + ) + + m.fs.unit.feed_side.properties_in[0].pressure.fix(lowP_pressure) + m.fs.unit.feed_side.properties_in[0].temperature.fix(temperature) + + m.fs.unit.brine_side.properties_in[0].flow_vol_phase["Liq"].fix(flow_vol) + m.fs.unit.brine_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].fix( + highP_mass_frac_TDS + ) + m.fs.unit.brine_side.properties_in[0].pressure.fix(highP_pressure) + m.fs.unit.brine_side.properties_in[0].temperature.fix(temperature) + + m.fs.unit.leakage_vol[0].fix(leakage_vol) + m.fs.unit.mixing_vol[0].fix(mixing_vol) + + # solve inlet conditions and only fix state variables (i.e. unfix flow_vol and mass_frac_phase) + + results = solver.solve(m.fs.unit.brine_side.properties_in[0]) + assert results.solver.termination_condition == TerminationCondition.optimal + m.fs.unit.brine_side.properties_in[0].flow_mass_phase_comp["Liq", "TDS"].fix() + m.fs.unit.brine_side.properties_in[0].flow_vol_phase["Liq"].unfix() + m.fs.unit.brine_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].unfix() + + results = solver.solve(m.fs.unit.feed_side.properties_in[0]) + assert results.solver.termination_condition == TerminationCondition.optimal + m.fs.unit.feed_side.properties_in[0].flow_mass_phase_comp["Liq", "H2O"].fix() + m.fs.unit.feed_side.properties_in[0].flow_mass_phase_comp["Liq", "TDS"].fix() + m.fs.unit.feed_side.properties_in[0].flow_vol_phase["Liq"].unfix() + m.fs.unit.feed_side.properties_in[0].mass_frac_phase_comp["Liq", "TDS"].unfix() + + # Specify unit + high_pressure_difference = 0.8 * pyunits.bar + low_pressure_difference = 0 * pyunits.bar + m.fs.unit.high_pressure_difference.fix(high_pressure_difference) + m.fs.unit.low_pressure_difference.fix(low_pressure_difference) + return m + + @pytest.mark.unit + def test_dof(self, unit_frame): + assert degrees_of_freedom(unit_frame) == 0 + + @pytest.mark.unit + def test_calculate_scaling(self, unit_frame): + m = unit_frame + calculate_scaling_factors(m) + + # check that all variables have scaling factors + unscaled_var_list = list( + unscaled_variables_generator(m.fs.unit, include_fixed=True) + ) + assert len(unscaled_var_list) == 0 + # check that all constraints have been scaled + unscaled_constraint_list = list(unscaled_constraints_generator(m)) + assert len(unscaled_constraint_list) == 0 + + @pytest.mark.requires_idaes_solver + @pytest.mark.component + def test_initialize(self, unit_frame): + m = unit_frame + initialization_tester(unit_frame) + + @pytest.mark.requires_idaes_solver + @pytest.mark.component + def test_solve(self, unit_frame): + m = unit_frame + results = solver.solve(m) + + # Check for optimal solution + assert results.solver.termination_condition == TerminationCondition.optimal + assert results.solver.status == SolverStatus.ok + + @pytest.mark.requires_idaes_solver + @pytest.mark.component + def test_solution(self, unit_frame): + m = unit_frame + assert pytest.approx(0.9877, rel=1e-3) == value( + m.fs.unit.feed_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + ) + assert pytest.approx(3.582e-2, rel=1e-3) == value( + m.fs.unit.feed_inlet.flow_mass_phase_comp[0, "Liq", "TDS"] + ) + assert pytest.approx(4.92e6, rel=1e-3) == value( + m.fs.unit.feed_outlet.pressure[0] + ) + assert pytest.approx(0.9868, rel=1e-3) == value( + m.fs.unit.brine_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + ) + assert pytest.approx(7.352e-2, rel=1e-3) == value( + m.fs.unit.brine_inlet.flow_mass_phase_comp[0, "Liq", "TDS"] + ) + assert pytest.approx(4.819e6, rel=1e-3) == value(m.fs.unit.feed_side.deltaP[0]) + assert pytest.approx(4.819e3, rel=1e-3) == value(m.fs.unit.feed_side.work[0]) + assert pytest.approx(-4.899e6, rel=1e-3) == value( + m.fs.unit.brine_side.deltaP[0] + ) + assert pytest.approx(-4.948e3, rel=1e-3) == value(m.fs.unit.brine_side.work[0]) + # testing solvent transfer - assert pytest.approx(-0.9877 * 1.1 * 0.1, rel=1e-3) == value( - m.fs.unit.high_pressure_side.mass_transfer_term[0, "Liq", "H2O"] + assert pytest.approx(3.559e-4, rel=1e-3) == value( + m.fs.unit.brine_side.mass_transfer_term[0, "Liq", "H2O"] ) - assert pytest.approx(0.9778, rel=1e-3) == value( - m.fs.unit.high_pressure_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] + assert pytest.approx(0.9871, rel=1e-3) == value( + m.fs.unit.brine_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] ) - assert pytest.approx(0.9877 + 0.9877 * 1.1 * 0.1, rel=1e-3) == value( - m.fs.unit.low_pressure_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] + assert pytest.approx(0.9874, rel=1e-3) == value( + m.fs.unit.feed_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] ) # testing solute transfer - assert pytest.approx(-7.352e-2 * 0.05, rel=1e-3) == value( - m.fs.unit.high_pressure_side.mass_transfer_term[0, "Liq", "TDS"] + assert pytest.approx(-1.293e-3, rel=1e-3) == value( + m.fs.unit.brine_side.mass_transfer_term[0, "Liq", "TDS"] ) - assert pytest.approx(7.352e-2 - 7.352e-2 * 0.05, rel=1e-3) == value( - m.fs.unit.high_pressure_outlet.flow_mass_phase_comp[0, "Liq", "TDS"] + assert pytest.approx(7.222e-2, rel=1e-3) == value( + m.fs.unit.brine_outlet.flow_mass_phase_comp[0, "Liq", "TDS"] ) - assert pytest.approx(3.582e-2 + 7.352e-2 * 0.05, rel=1e-3) == value( - m.fs.unit.low_pressure_outlet.flow_mass_phase_comp[0, "Liq", "TDS"] + assert pytest.approx(3.712e-2, rel=1e-3) == value( + m.fs.unit.feed_outlet.flow_mass_phase_comp[0, "Liq", "TDS"] ) @pytest.mark.unit @@ -359,22 +626,31 @@ class TestPressureExchanger_with_ion_prop_pack: def unit_frame(self): m = ConcreteModel() m.fs = FlowsheetBlock(dynamic=False) - m.fs.properties = property_seawater_ions.PropParameterBlock() + m.fs.properties = props_multi.MCASParameterBlock( + solute_list=["Na", "Ca", "Mg", "SO4", "Cl"], + mw_data={ + "H2O": 0.018, + "Na": 0.023, + "Ca": 0.04, + "Mg": 0.024, + "SO4": 0.096, + "Cl": 0.035, + }, + ignore_neutral_charge=True, + ) m.fs.unit = PressureExchanger( - property_package=m.fs.properties, has_mass_transfer=True + property_package=m.fs.properties, + has_leakage=True, + has_mixing=True, ) # Specify inlet conditions temperature = 25 + 273.15 lowP_pressure = 101325 - highP_mass_frac_TDS = 0.07 highP_pressure = 50e5 - solvent_transfer = 0.05 - solute_transfer = {"Na": 0.1, "Ca": 0.11, "Mg": 0.12, "SO4": 0.13, "Cl": 0.14} - - feed_flow_mass = 1 + feed_flow_mass = 1 * pyunits.kg / pyunits.s lowP_mass_frac = { "Na": 11122e-6, "Ca": 382e-6, @@ -391,23 +667,31 @@ def unit_frame(self): "Cl": 2 * 20300e-6, } + leakage_vol = 0.01 + mixing_vol = 0.035 + # Specify unit efficiency = 0.95 m.fs.unit.efficiency_pressure_exchanger.fix(efficiency) - m.fs.unit.mass_transfer_fraction_comp[0, "H2O"].fix(solvent_transfer) - for s in solute_transfer: - m.fs.unit.mass_transfer_fraction_comp[0, s].fix(solute_transfer[s]) + m.fs.unit.leakage_vol[0].fix(leakage_vol) + m.fs.unit.mixing_vol[0].fix(mixing_vol) + + LPS_in = m.fs.unit.feed_side.properties_in[0] + for ion, x in lowP_mass_frac.items(): + mol_comp_flow = ( + x * pyunits.kg / pyunits.kg * feed_flow_mass / LPS_in.mw_comp[ion] + ) + + LPS_in.flow_mol_phase_comp["Liq", ion].fix(mol_comp_flow) - m.fs.unit.low_pressure_side.properties_in[0].flow_mass_phase_comp[ - "Liq", "H2O" - ].fix(feed_flow_mass * (1 - sum(x for x in lowP_mass_frac.values()))) - for s in lowP_mass_frac: - m.fs.unit.low_pressure_side.properties_in[0].flow_mass_phase_comp[ - "Liq", s - ].fix(feed_flow_mass * lowP_mass_frac[s]) - m.fs.unit.low_pressure_side.properties_in[0].pressure.fix(lowP_pressure) - m.fs.unit.low_pressure_side.properties_in[0].temperature.fix(temperature) + LPS_in.flow_mol_phase_comp["Liq", "H2O"].fix( + (feed_flow_mass * (1 - sum(x for x in lowP_mass_frac.values()))) + / LPS_in.mw_comp["H2O"] + ) + + m.fs.unit.feed_side.properties_in[0].pressure.fix(lowP_pressure) + m.fs.unit.feed_side.properties_in[0].temperature.fix(temperature) ### NOTE THE PRESSURE EXHCANGER EXPECTS EQUAL VOLUMETRIC FLOW RATES # so one of mass fraction (here we unfix H2O mass fraction) @@ -418,17 +702,22 @@ def unit_frame(self): # a -1 DOF error and should not occur in practice where typicly at least one of # the inlet mass flow rates is colaculated. - m.fs.unit.high_pressure_side.properties_in[0].flow_mass_phase_comp[ - "Liq", "H2O" - ] = feed_flow_mass * (1 - sum(x for x in highP_mass_frac.values())) + HPS_in = m.fs.unit.brine_side.properties_in[0] + for ion, x in highP_mass_frac.items(): + mol_comp_flow = ( + x * pyunits.kg / pyunits.kg * feed_flow_mass / HPS_in.mw_comp[ion] + ) - for s in highP_mass_frac: - m.fs.unit.high_pressure_side.properties_in[0].flow_mass_phase_comp[ - "Liq", s - ].fix(feed_flow_mass * highP_mass_frac[s]) + HPS_in.flow_mol_phase_comp["Liq", ion].fix(mol_comp_flow) - m.fs.unit.high_pressure_side.properties_in[0].pressure.fix(highP_pressure) - m.fs.unit.high_pressure_side.properties_in[0].temperature.fix(temperature) + HPS_in.flow_mol_phase_comp["Liq", "H2O"] = ( + feed_flow_mass + * (1 - sum(x for x in highP_mass_frac.values())) + / HPS_in.mw_comp["H2O"] + ) + + m.fs.unit.brine_side.properties_in[0].pressure.fix(highP_pressure) + m.fs.unit.brine_side.properties_in[0].temperature.fix(temperature) return m @@ -450,17 +739,13 @@ def test_calculate_scaling(self, unit_frame): unscaled_constraint_list = list(unscaled_constraints_generator(m)) assert len(unscaled_constraint_list) == 0 + @pytest.mark.requires_idaes_solver @pytest.mark.component def test_initialize(self, unit_frame): m = unit_frame initialization_tester(unit_frame) - @pytest.mark.component - def test_var_scaling(self, unit_frame): - m = unit_frame - badly_scaled_var_lst = list(badly_scaled_var_generator(m)) - assert badly_scaled_var_lst == [] - + @pytest.mark.requires_idaes_solver @pytest.mark.component def test_solve(self, unit_frame): m = unit_frame @@ -470,91 +755,49 @@ def test_solve(self, unit_frame): assert results.solver.termination_condition == TerminationCondition.optimal assert results.solver.status == SolverStatus.ok + @pytest.mark.requires_idaes_solver @pytest.mark.component def test_solution(self, unit_frame): m = unit_frame - solvent_transfer = 0.05 - solute_transfer = {"Na": 0.1, "Ca": 0.11, "Mg": 0.12, "SO4": 0.13, "Cl": 0.14} - - feed_flow_mass = 1 - lowP_mass_frac = { - "Na": 11122e-6, - "Ca": 382e-6, - "Mg": 1394e-6, - "SO4": 2136e-6, - "Cl": 20300e-6, - } - highP_mass_frac = { - "Na": 2 * 11122e-6, - "Ca": 2 * 382e-6, - "Mg": 2 * 1394e-6, - "SO4": 2 * 2136e-6, - "Cl": 2 * 20300e-6, - } - - assert pytest.approx( - feed_flow_mass * (1 - sum(x for x in lowP_mass_frac.values())), rel=1e-3 - ) == value(m.fs.unit.low_pressure_inlet.flow_mass_phase_comp[0, "Liq", "H2O"]) - - for s in lowP_mass_frac: - assert pytest.approx(feed_flow_mass * lowP_mass_frac[s], rel=1e-3) == value( - m.fs.unit.low_pressure_inlet.flow_mass_phase_comp[0, "Liq", s] - ) - assert pytest.approx(4.755e6, rel=1e-3) == value( - m.fs.unit.low_pressure_outlet.pressure[0] + m.fs.unit.feed_outlet.pressure[0] ) - assert pytest.approx(0.9876, rel=1e-3) == value( - m.fs.unit.high_pressure_inlet.flow_mass_phase_comp[0, "Liq", "H2O"] + assert pytest.approx(4.654e6, rel=1e-3) == value(m.fs.unit.feed_side.deltaP[0]) + assert pytest.approx(4.654e3, rel=1e-3) == value(m.fs.unit.feed_side.work[0]) + assert pytest.approx(-4.899e6, rel=1e-3) == value( + m.fs.unit.brine_side.deltaP[0] ) + assert pytest.approx(-4.948e3, rel=1e-3) == value(m.fs.unit.brine_side.work[0]) - for s in highP_mass_frac: - assert pytest.approx( - feed_flow_mass * highP_mass_frac[s], rel=1e-3 - ) == value(m.fs.unit.high_pressure_inlet.flow_mass_phase_comp[0, "Liq", s]) - - assert pytest.approx(4.654e6, rel=1e-3) == value( - m.fs.unit.low_pressure_side.deltaP[0] + # testing solvent transfer + assert pytest.approx(0.06733, rel=1e-3) == value( + m.fs.unit.brine_side.mass_transfer_term[0, "Liq", "H2O"] ) - assert pytest.approx(4.654e3, rel=1e-3) == value( - m.fs.unit.low_pressure_side.work[0] + assert pytest.approx(0.5, rel=1e-3) == value( + m.fs.unit.brine_side.properties_out[0].flow_mass_phase_comp["Liq", "H2O"] ) - assert pytest.approx(-4.899e6, rel=1e-3) == value( - m.fs.unit.high_pressure_side.deltaP[0] - ) - assert pytest.approx(-5.184e3, rel=1e-3) == value( - m.fs.unit.high_pressure_side.work[0] + assert pytest.approx(0.9634, rel=1e-3) == value( + m.fs.unit.feed_side.properties_out[0].flow_mass_phase_comp["Liq", "H2O"] ) - # testing solvent transfer - assert pytest.approx(-0.9876 * solvent_transfer, rel=1e-3) == value( - m.fs.unit.high_pressure_side.mass_transfer_term[0, "Liq", "H2O"] + # testing solute transfer + assert pytest.approx(-0.01659, rel=1e-3) == value( + m.fs.unit.brine_side.mass_transfer_term[0, "Liq", "Na"] ) - assert pytest.approx(0.9876 - 0.9876 * solvent_transfer, rel=1e-3) == value( - m.fs.unit.high_pressure_outlet.flow_mass_phase_comp[0, "Liq", "H2O"] + assert pytest.approx(-3.276e-4, rel=1e-3) == value( + m.fs.unit.brine_side.mass_transfer_term[0, "Liq", "Ca"] + ) + assert pytest.approx(-1.992e-3, rel=1e-3) == value( + m.fs.unit.brine_side.mass_transfer_term[0, "Liq", "Mg"] + ) + assert pytest.approx(-7.632e-4, rel=1e-3) == value( + m.fs.unit.brine_side.mass_transfer_term[0, "Liq", "SO4"] + ) + assert pytest.approx(-1.989e-2, rel=1e-3) == value( + m.fs.unit.brine_side.mass_transfer_term[0, "Liq", "Cl"] ) - assert pytest.approx( - feed_flow_mass * (1 - sum(x for x in lowP_mass_frac.values())) - + 0.9876 * solvent_transfer, - rel=1e-3, - ) == value(m.fs.unit.low_pressure_outlet.flow_mass_phase_comp[0, "Liq", "H2O"]) - # testing solute transfer - for s in highP_mass_frac: - assert pytest.approx( - -feed_flow_mass * highP_mass_frac[s] * solute_transfer[s], rel=1e-3 - ) == value(m.fs.unit.high_pressure_side.mass_transfer_term[0, "Liq", s]) - assert pytest.approx( - feed_flow_mass * highP_mass_frac[s] - - feed_flow_mass * highP_mass_frac[s] * solute_transfer[s], - rel=1e-3, - ) == value(m.fs.unit.high_pressure_outlet.flow_mass_phase_comp[0, "Liq", s]) - assert pytest.approx( - feed_flow_mass * lowP_mass_frac[s] - + feed_flow_mass * highP_mass_frac[s] * solute_transfer[s], - rel=1e-3, - ) == value(m.fs.unit.low_pressure_outlet.flow_mass_phase_comp[0, "Liq", s]) @pytest.mark.unit def test_report(self, unit_frame):