BSM2 Anoxic CSTR Costing

watertap/costing/costing_base.py

@@ -13,12 +13,13 @@
 import pyomo.environ as pyo
 from idaes.core import declare_process_block_class
 from idaes.core.base.costing_base import FlowsheetCostingBlockData
-from idaes.models.unit_models import Mixer, HeatExchanger
+from idaes.models.unit_models import Mixer, HeatExchanger, CSTR
 
 from watertap.core.util.misc import is_constant_up_to_units
 
 from watertap.costing.unit_models.mixer import cost_mixer
 from watertap.costing.unit_models.heat_exchanger import cost_heat_exchanger
+from watertap.costing.unit_models.cstr import cost_cstr
 
 
 @declare_process_block_class("WaterTAPCostingBlock")
@@ -33,6 +34,7 @@ class WaterTAPCostingBlockData(FlowsheetCostingBlockData):
     unit_mapping = {
         Mixer: cost_mixer,
         HeatExchanger: cost_heat_exchanger,
+        CSTR: cost_cstr,
     }
 
     def build(self):

watertap/costing/unit_models/cstr.py

@@ -0,0 +1,69 @@
+#################################################################################
+# WaterTAP Copyright (c) 2020-2023, The Regents of the University of California,
+# through Lawrence Berkeley National Laboratory, Oak Ridge National Laboratory,
+# National Renewable Energy Laboratory, and National Energy Technology
+# Laboratory (subject to receipt of any required approvals from the U.S. Dept.
+# of Energy). All rights reserved.
+#
+# Please see the files COPYRIGHT.md and LICENSE.md for full copyright and license
+# information, respectively. These files are also available online at the URL
+# "https://github.com/watertap-org/watertap/"
+#################################################################################
+
+import pyomo.environ as pyo
+from ..util import (
+    register_costing_parameter_block,
+    make_capital_cost_var,
+)
+
+
+def build_cstr_cost_param_block(blk):
+    # Source: https://www.fwrj.com/articles/9812.pdf
+    blk.sizing_cost = pyo.Var(
+        initialize=0.34,
+        doc="Reactor sizing cost",
+        units=pyo.units.USD_1998 / pyo.units.m**3,
+    )
+
+
+@register_costing_parameter_block(
+    build_rule=build_cstr_cost_param_block,
+    parameter_block_name="cstr",
+)
+def cost_cstr(blk):
+    """
+    CSTR costing method
+    """
+    cost_cstr_capital(
+        blk,
+        blk.costing_package.cstr.sizing_cost,
+    )
+
+
+def cost_cstr_capital(blk, sizing_cost):
+    """
+    Generic function for costing an ElectroNP system.
+    """
+    make_capital_cost_var(blk)
+
+    blk.sizing_cost = pyo.Expression(expr=sizing_cost)
+
+    flow_in = pyo.units.convert(
+        blk.unit_model.control_volume.properties_in[0].flow_vol,
+        to_units=pyo.units.m**3 / pyo.units.hr,
+    )
+
+    HRT = pyo.units.convert(
+        blk.unit_model.volume[0]
+        / blk.unit_model.control_volume.properties_in[0].flow_vol,
+        to_units=pyo.units.hr,
+    )
+
+    print(f"base_currency: {blk.costing_package.base_currency}")
+    blk.capital_cost_constraint = pyo.Constraint(
+        expr=blk.capital_cost
+        == pyo.units.convert(
+            HRT * flow_in * blk.sizing_cost,
+            to_units=blk.costing_package.base_currency,
+        )
+    )

watertap/examples/flowsheets/case_studies/full_water_resource_recovery_facility/BSM2.py

@@ -26,6 +26,7 @@
 from watertap.unit_models.anaerobic_digestor import AD
 from watertap.unit_models.thickener import Thickener
 from watertap.unit_models.dewatering import DewateringUnit
+from watertap.unit_models.cstr import AnoxicCSTR
 
 from watertap.unit_models.translators.translator_asm1_adm1 import Translator_ASM1_ADM1
 from watertap.unit_models.translators.translator_adm1_asm1 import Translator_ADM1_ASM1
@@ -38,9 +39,6 @@
 from watertap.property_models.anaerobic_digestion.adm1_properties import (
     ADM1ParameterBlock,
 )
-from watertap.property_models.activated_sludge.asm1_properties import (
-    ASM1ParameterBlock,
-)
 from watertap.property_models.anaerobic_digestion.adm1_reactions import (
     ADM1ReactionParameterBlock,
 )
@@ -49,9 +47,8 @@
     ADM1_vaporParameterBlock,
 )
 
-from idaes.core import FlowsheetBlock
+from idaes.core import FlowsheetBlock, UnitModelCostingBlock
 from idaes.models.unit_models import (
-    CSTR,
     Feed,
     Mixer,
     Separator,
@@ -65,6 +62,7 @@
     ASM1ReactionParameterBlock,
 )
 from watertap.core.util.initialization import assert_degrees_of_freedom
+from watertap.costing import WaterTAPCosting
 from pyomo.util.check_units import assert_units_consistent
 
 
@@ -79,15 +77,14 @@ def main():
     results = solve(m)
 
     add_costing(m)
+    m.fs.costing.initialize()
     # Assert DOF = 0 after adding costing
     # assert_degrees_of_freedom(m, 0)
 
-    # TODO: initialize costing after adding to flowsheet
-    # m.fs.costing.initialize()
-
-    # results = solve(m)
+    results = solve(m)
 
     display_results(m)
+    display_costing(m)
 
     return m, results
 
@@ -113,11 +110,11 @@ def build():
         property_package=m.fs.props_ASM1, inlet_list=["feed_water", "recycle"]
     )
     # First reactor (anoxic) - standard CSTR
-    m.fs.R1 = CSTR(
+    m.fs.R1 = AnoxicCSTR(
         property_package=m.fs.props_ASM1, reaction_package=m.fs.ASM1_rxn_props
     )
     # Second reactor (anoxic) - standard CSTR
-    m.fs.R2 = CSTR(
+    m.fs.R2 = AnoxicCSTR(
         property_package=m.fs.props_ASM1, reaction_package=m.fs.ASM1_rxn_props
     )
     # Third reactor (aerobic) - CSTR with injection
@@ -442,7 +439,30 @@ def function(unit):
 
 def add_costing(m):
     # TODO: implement unit model and flowsheet level costing
-    pass
+    m.fs.costing = WaterTAPCosting()
+    m.fs.costing.base_currency = pyo.units.USD_2020
+
+    # Costing Blocks
+    m.fs.R1.costing = UnitModelCostingBlock(flowsheet_costing_block=m.fs.costing)
+    m.fs.R2.costing = UnitModelCostingBlock(flowsheet_costing_block=m.fs.costing)
+
+    m.fs.RADM.costing = UnitModelCostingBlock(flowsheet_costing_block=m.fs.costing)
+
+    # Leaving out mixer costs for now
+    # m.fs.MX1.costing = UnitModelCostingBlock(flowsheet_costing_block=m.fs.costing)
+    # m.fs.MX6.costing = UnitModelCostingBlock(flowsheet_costing_block=m.fs.costing)
+    # m.fs.MX2.costing = UnitModelCostingBlock(flowsheet_costing_block=m.fs.costing)
+    # m.fs.MX3.costing = UnitModelCostingBlock(flowsheet_costing_block=m.fs.costing)
+    # m.fs.MX4.costing = UnitModelCostingBlock(flowsheet_costing_block=m.fs.costing)
+
+    # process costing and add system level metrics
+    m.fs.costing.cost_process()
+    m.fs.costing.add_electricity_intensity(m.fs.FeedWater.properties[0].flow_vol)
+    m.fs.costing.add_annual_water_production(m.fs.Treated.properties[0].flow_vol)
+    m.fs.costing.add_LCOW(m.fs.Treated.properties[0].flow_vol)
+    m.fs.costing.add_specific_energy_consumption(m.fs.Treated.properties[0].flow_vol)
+
+    m.fs.objective = pyo.Objective(expr=m.fs.costing.LCOW)
 
 
 def solve(blk, solver=None):
@@ -484,5 +504,28 @@ def display_results(m):
         m.fs.component(u).report()
 
 
+def display_costing(m):
+    print(
+        "Energy consumption: %.1f kWh/m3"
+        % pyo.value(m.fs.costing.specific_energy_consumption)
+    )
+    print(
+        "Electricity intensity: %.1f kWh/m3"
+        % pyo.value(m.fs.costing.electricity_intensity)
+    )
+    print("Levelized cost of water: %.2f $/m3" % pyo.value(m.fs.costing.LCOW))
+
+    print(
+        "Total operating cost: %.2f $/yr" % pyo.value(m.fs.costing.total_operating_cost)
+    )
+
+    print("Total capital cost: %.2f $" % pyo.value(m.fs.costing.total_capital_cost))
+
+    print(
+        "Total annualized cost: %.2f $/yr"
+        % pyo.value(m.fs.costing.total_annualized_cost)
+    )
+
+
 if __name__ == "__main__":
     m, results = main()

watertap/unit_models/__init__.py

@@ -29,3 +29,4 @@
 from .thickener import Thickener
 from .dewatering import DewateringUnit
 from .electroNP_ZO import ElectroNPZO
+from .cstr import AnoxicCSTR

watertap/unit_models/cstr.py

@@ -0,0 +1,84 @@
+###############################################################################
+# WaterTAP Copyright (c) 2021, The Regents of the University of California,
+# through Lawrence Berkeley National Laboratory, Oak Ridge National
+# Laboratory, National Renewable Energy Laboratory, and National Energy
+# Technology Laboratory (subject to receipt of any required approvals from
+# the U.S. Dept. of Energy). All rights reserved.
+#
+# Please see the files COPYRIGHT.md and LICENSE.md for full copyright and license
+# information, respectively. These files are also available online at the URL
+# "https://github.com/watertap-org/watertap/"
+#
+###############################################################################
+"""
+Anoxic CSTR unit model for BSM2 and plant-wide wastewater treatment modeling.
+This unit inherits from the IDAES CSTR unit.
+"""
+
+# Import IDAES cores
+from idaes.core import (
+    declare_process_block_class,
+)
+from idaes.models.unit_models.cstr import CSTRData
+
+import idaes.logger as idaeslog
+
+from pyomo.environ import (
+    Constraint,
+    NonNegativeReals,
+    Var,
+    units as pyunits,
+)
+
+from watertap.costing.unit_models.cstr import cost_cstr
+
+__author__ = "Marcus Holly"
+
+
+# Set up logger
+_log = idaeslog.getLogger(__name__)
+
+
+@declare_process_block_class("AnoxicCSTR")
+class AnoxicCSTRData(CSTRData):
+    """
+    Anoxic CSTR unit block for BSM2
+    """
+
+    CONFIG = CSTRData.CONFIG()
+
+    def build(self):
+        """
+        Begin building model.
+        Args:
+            None
+        Returns:
+            None
+        """
+
+        # Call UnitModel.build to set up dynamics
+        super(AnoxicCSTRData, self).build()
+
+        self.HRT = Var(
+            self.flowsheet().time,
+            initialize=4,
+            domain=NonNegativeReals,
+            units=pyunits.s,
+            doc="Hydraulic retention time",
+        )
+
+        def CSTR_retention_time_rule(self, t):
+            return (
+                self.HRT[t]
+                == self.volume[t] / self.control_volume.properties_in[t].flow_vol
+            )
+
+        self.CSTR_retention_time = Constraint(
+            self.flowsheet().time,
+            rule=CSTR_retention_time_rule,
+            doc="Total CSTR retention time",
+        )
+
+    @property
+    def default_costing_method(self):
+        return cost_cstr