1D Membrane Model for CO2 Capture and Utilization

docs/reference_guides/model_libraries/models_extra/index.rst

@@ -6,3 +6,5 @@ Additional IDAES Model Libraries
 
     phe
     temperature_swing_adsorption/fixed_bed_tsa0d
+    membrane_model/1d_membrane
+

docs/reference_guides/model_libraries/models_extra/membrane_model/1d_membrane.rst

@@ -0,0 +1,14 @@
+One-dimensional membrane class for CO2 gas separation
+=============================================
+
+This is a one-dimensional model for gas separation in CO₂ capture applications.
+The model will be discretized in the flow direction, and it supports two flow patterns:
+counter-current flow and co-current flow. The model was customized for gas-phase separation
+in CO₂ capture with a single-layer design. If a multi-layer design is needed, multiple units
+can be connected for this application. The two sides of the membrane are called the feed side
+and sweep side. The sweep stream inlet is optional. The driving force across the membrane is the
+partial pressure difference in this gas separation application. Additionally, the energy balance
+assumes that temperature remains constant on each side of the membrane.
+
+
+

idaes/models_extra/co2_capture_and_utilization/unit_models/__init__.py

@@ -3,7 +3,7 @@
 # Framework (IDAES IP) was produced under the DOE Institute for the
 # Design of Advanced Energy Systems (IDAES).
 #
-# Copyright (c) 2018-2023 by the software owners: The Regents of the
+# Copyright (c) 2018-2024 by the software owners: The Regents of the
 # University of California, through Lawrence Berkeley National Laboratory,
 # National Technology & Engineering Solutions of Sandia, LLC, Carnegie Mellon
 # University, West Virginia University Research Corporation, et al.

idaes/models_extra/co2_capture_and_utilization/unit_models/membrane_1d.py

@@ -3,7 +3,7 @@
 # Framework (IDAES IP) was produced under the DOE Institute for the
 # Design of Advanced Energy Systems (IDAES).
 #
-# Copyright (c) 2018-2023 by the software owners: The Regents of the
+# Copyright (c) 2018-2024 by the software owners: The Regents of the
 # University of California, through Lawrence Berkeley National Laboratory,
 # National Technology & Engineering Solutions of Sandia, LLC, Carnegie Mellon
 # University, West Virginia University Research Corporation, et al.
@@ -15,16 +15,16 @@
 One-dimensional membrane class for CO2 gas separation
 """
 
-# pylint: disable=unused-import
+
 from enum import Enum
 from pyomo.common.config import Bool, ConfigDict, ConfigValue, In
 from pyomo.environ import (
     Constraint,
     Param,
     Var,
     units,
     Expression,
 )
 from pyomo.network import Port
 
 from idaes.core import (
@@ -37,6 +37,7 @@
 from idaes.core.util.config import is_physical_parameter_block
 from idaes.models.unit_models.mscontactor import MSContactor
 from idaes.core.util.exceptions import ConfigurationError
+from idaes.core.util.tables import create_stream_table_dataframe
 
 __author__ = "Maojian Wang"
 
@@ -49,7 +50,6 @@
 
     COUNTERCURRENT = 1
     COCURRENT = 2
-    CROSSFLOW = 3
 
 
 @declare_process_block_class("Membrane1D")
@@ -85,6 +85,7 @@
     see property package for documentation.}""",
         ),
     )
+
     Stream_Config.declare(
         "has_energy_balance",
         ConfigValue(
@@ -134,31 +135,6 @@
                                                     sweep side flows from 1 to 0  (default)""",
         ),
     )
-    CONFIG.declare(
-        "property_package",
-        ConfigValue(
-            default=None,
-            domain=is_physical_parameter_block,
-            description="Property package to use for control volume",
-            doc="""Property parameter object used to define property
-            calculations
-            (default = 'use_parent_value')
-            - 'use_parent_value' - get package from parent (default = None)
-            - a ParameterBlock object""",
-        ),
-    )
-    CONFIG.declare(
-        "property_package_args",
-        ConfigValue(
-            default={},
-            description="Arguments for constructing property package",
-            doc="""A dict of arguments to be passed to the PropertyBlockData
-            and used when constructing these
-            (default = 'use_parent_value')
-            - 'use_parent_value' - get package from parent (default = None)
-            - a dict (see property package for documentation)""",
-        ),
-    )
 
     for side_name in ["feed", "sweep"]:
         CONFIG.declare(
@@ -167,14 +143,19 @@
         )
 
     def build(self):
+        """
+        This is a one-dimensional model for gas separation in CO₂ capture applications.
+        The model will be discretized in the flow direction, and it supports two flow patterns:
+        counter-current flow and co-current flow. The model was customized for gas-phase separation
+        in CO₂ capture with a single-layer design. If a multi-layer design is needed, multiple units
+        can be connected for this application. The two sides of the membrane are called the feed side
+        and sweep side. The sweep stream inlet is optional. The driving force across the membrane is the
+        partial pressure difference in this gas separation application. Additionally, the energy balance
+        assumes that temperature remains constant on each side of the membrane.
+
+        """
         super().build()
 
-        if self.config.property_package is not None:
-            if self.config.feed_side.property_package == useDefault:
-                self.config.feed_side.property_package = self.config.property_package
-            if self.config.sweep_side.property_package == useDefault:
-                self.config.sweep_side.property_package = self.config.property_package
-
         feed_dict = dict(self.config.feed_side)
         sweep_dict = dict(self.config.sweep_side)
 
@@ -210,14 +191,14 @@
 
     def _make_geometry(self):
 
-        self.area = Var(initialize=100, units=units.cm**2, doc="The membrane area")
+        self.area = Var(
+            initialize=100, units=units.cm**2, doc="Area per cell (or finite element)"
+        )
 
         self.length = Var(initialize=100, units=units.cm, doc="The membrane length")
         self.cell_length = Expression(expr=self.length / self.config.finite_elements)
 
-        self.cell_area = Var(
-            initialize=100, units=units.cm**2, doc="The membrane area"
-        )
+        self.cell_area = Var(initialize=100, units=units.cm**2, doc="The membrane area")
 
         @self.Constraint()
         def area_per_cell(self):
@@ -227,50 +208,33 @@
         feed_side_units = (
             self.config.feed_side.property_package.get_metadata().derived_units
         )
+        crossover_component_list = list(
+            set(self.mscontactor.feed_side.component_list)
+            & set(self.mscontactor.sweep_side.component_list)
+        )
 
         self.permeance = Var(
             self.flowsheet().time,
             self.mscontactor.elements,
-            self.mscontactor.feed_side.component_list,
+            crossover_component_list,
             initialize=1,
-            doc="Values in Gas Permeance Unit(GPU)",
+            doc="Values in Gas Permeance Unit (GPU)",
             units=units.dimensionless,
         )
 
         self.gpu_factor = Param(
             default=10e-8 / 13333.2239,
             units=units.m / units.s / units.Pa,
             mutable=True,
+            # This is a coefficient that will convert the unit of permeability from GPU to SI units for further calculation"
         )
 
-        self.selectivity = Var(
-            self.flowsheet().time,
-            self.mscontactor.elements,
-            self.mscontactor.feed_side.component_list,
-            self.mscontactor.feed_side.component_list,
-            initialize=1,
-            units=units.dimensionless,
-        )
-
-        @self.Constraint(
-            self.flowsheet().time,
-            self.mscontactor.elements,
-            self.mscontactor.feed_side.component_list,
-            self.mscontactor.feed_side.component_list,
-            doc="permeance calculation",
-        )
-        def permeance_calculation(self, t, e, a, b):
-            return (
-                self.permeance[t, e, a] * self.selectivity[t, e, a, b]
-                == self.permeance[t, e, b]
-            )
-
         p_units = feed_side_units.PRESSURE
 
         @self.Constraint(
             self.flowsheet().time,
             self.mscontactor.elements,
-            self.mscontactor.feed_side.component_list,
+            crossover_component_list,
             doc="permeability calculation",
         )
         def permeability_calculation(self, t, s, m):
@@ -282,7 +246,9 @@
                 mb_units = feed_side_units.FLOW_MASS
                 rho = self.mscontactor.feed_side[t, s].dens_mass
             else:
-                raise TypeError("Undefined flow basis, please define the flow basis")
+                raise TypeError(
+                    "This model only supports MaterialFlowBasis equal to molar or mass"
+                )
 
             return self.mscontactor.material_transfer_term[
                 t, s, "feed_side", "sweep_side", m
@@ -307,10 +273,31 @@
         @self.Constraint(
             self.flowsheet().time,
             self.mscontactor.elements,
-            doc="Energy balance",
+            doc="isothermal constraint",
         )
         def energy_transfer(self, t, s):
             return (
                 self.mscontactor.feed_side[t, s].temperature
                 == self.mscontactor.sweep_side[t, s].temperature
             )
+
+    def _get_stream_table_contents(self, time_point=0):
+        if self.config.sweep_flow:
+            return create_stream_table_dataframe(
+                {
+                    "Feed Inlet": self.feed_side_inlet,
+                    "Feed Outlet": self.feed_side_outlet,
+                    "Permeate Inlet": self.sweep_side_inlet,
+                    "Permeate Outlet": self.sweep_side_outlet,
+                },
+                time_point=time_point,
+            )
+        else:
+            return create_stream_table_dataframe(
+                {
+                    "Feed Inlet": self.feed_side_inlet,
+                    "Feed Outlet": self.feed_side_outlet,
+                    "Permeate Outlet": self.sweep_side_outlet,
+                },
+                time_point=time_point,
+            )

idaes/models_extra/co2_capture_and_utilization/unit_models/tests/test_membrane_1d.py

@@ -3,7 +3,7 @@
 # Framework (IDAES IP) was produced under the DOE Institute for the
 # Design of Advanced Energy Systems (IDAES).
 #
-# Copyright (c) 2018-2023 by the software owners: The Regents of the
+# Copyright (c) 2018-2024 by the software owners: The Regents of the
 # University of California, through Lawrence Berkeley National Laboratory,
 # National Technology & Engineering Solutions of Sandia, LLC, Carnegie Mellon
 # University, West Virginia University Research Corporation, et al.
@@ -20,6 +20,7 @@
 
 from pyomo.environ import (
     check_optimal_termination,
+    assert_optimal_termination,
     ConcreteModel,
     value,
 )
@@ -50,9 +51,10 @@
 # Get default solver for testing
 solver = get_solver()
 
+
 # -----------------------------------------------------------------------------
 @pytest.mark.unit
-def test_config():
+def test_config_countercurrent():
     m = ConcreteModel()
     m.fs = FlowsheetBlock(dynamic=False)
 
@@ -66,17 +68,44 @@ def test_config():
         dynamic=False,
         sweep_flow=True,
         flow_type=MembraneFlowPattern.COUNTERCURRENT,
-        property_package=m.fs.properties,
+        feed_side={"property_package": m.fs.properties},
+        sweep_side={"property_package": m.fs.properties},
+    )
+
+    # Check unit config arguments
+    print("=====================================================")
+    print(len(m.fs.unit.config))
+    assert len(m.fs.unit.config) == 7
+    assert not m.fs.unit.config.dynamic
+    assert not m.fs.unit.config.has_holdup
+
+
+@pytest.mark.unit
+def test_congif_cocurrent():
+    m = ConcreteModel()
+    m.fs = FlowsheetBlock(dynamic=False)
+
+    m.fs.properties = GenericParameterBlock(
+        **get_prop(["CO2", "H2O", "N2"], ["Vap"], eos=EosType.IDEAL),
+        doc="Key flue gas property parameters",
+    )
+
+    m.fs.unit = Membrane1D(
+        finite_elements=3,
+        dynamic=False,
+        sweep_flow=True,
+        flow_type=MembraneFlowPattern.COCURRENT,
+        feed_side={"property_package": m.fs.properties},
+        sweep_side={"property_package": m.fs.properties},
     )
 
     # Check unit config arguments
-    assert len(m.fs.unit.config) == 9
+    assert len(m.fs.unit.config) == 7
     assert not m.fs.unit.config.dynamic
     assert not m.fs.unit.config.has_holdup
-    assert m.fs.unit.config.property_package is m.fs.properties
 
 
-class TestMembrane(object):
+class TestMembrane:
     @pytest.fixture(scope="class")
     def membrane(self):
         m = ConcreteModel()
@@ -91,7 +120,8 @@ def membrane(self):
             dynamic=False,
             sweep_flow=True,
             flow_type=MembraneFlowPattern.COUNTERCURRENT,
-            property_package=m.fs.properties,
+            feed_side={"property_package": m.fs.properties},
+            sweep_side={"property_package": m.fs.properties},
         )
 
         m.fs.unit.permeance[:, :, "CO2"].fix(1500)
@@ -151,8 +181,8 @@ def test_build(self, membrane):
         assert hasattr(membrane.fs.unit, "permeability_calculation")
         assert hasattr(membrane.fs.unit, "energy_transfer")
 
-        assert number_variables(membrane) == 184
-        assert number_total_constraints(membrane) == 116
+        assert number_variables(membrane) == 157
+        assert number_total_constraints(membrane) == 89
         assert number_unused_variables(membrane) == 28
 
     @pytest.mark.component
@@ -168,7 +198,7 @@ def test_solve(self, membrane):
         initializer.initialize(membrane.fs.unit)
         results = solver.solve(membrane)
         # Check for optimal solution
-        assert check_optimal_termination(results)
+        assert_optimal_termination(results)
 
     @pytest.mark.solver
     @pytest.mark.skipif(solver is None, reason="Solver not available")
@@ -255,3 +285,22 @@ def test_enthalpy_balance(self, membrane):
             )
             <= 1e-6
         )
+
+    @pytest.mark.solver
+    @pytest.mark.skipif(solver is None, reason="Solver not available")
+    @pytest.mark.component
+    def test_material_balance(self, membrane):
+
+        assert (
+            abs(
+                value(
+                    (
+                        membrane.fs.unit.feed_side_inlet.flow_mol[0]
+                        + membrane.fs.unit.sweep_side_inlet.flow_mol[0]
+                        - membrane.fs.unit.feed_side_outlet.flow_mol[0]
+                        - membrane.fs.unit.sweep_side_outlet.flow_mol[0]
+                    )
+                )
+            )
+            <= 1e-3
+        )