Heat Exchanger NTU Unit Model with Aqueous MEA System ===================================================== .. image:: heat_exchanger_4.svg **Problem Statement**: In this example, we will be transfering heat from hot to cold streams of MEA (monoethanolamine), simulating heat integration of an aqueous MEA system. **Hot Side Inlet** Flow Rate = 60.54879 mol/s Mole fraction (CO2) = 0.0158 Mole fraction (H2O) = 0.8747 Mole fraction (MEA) = 0.1095 Pressure = 202650 Pa Temperature = 392.23 K **Cold Side Inlet** Flow Rate = 63.01910 mol/s Mole fraction (CO2) = 0.0414 Mole fraction (H2O) = 0.8509 Mole fraction (MEA) = 0.1077 Pressure = 202650 Pa Temperature = 326.36 K This example will demonstrate the simulation of the NTU heat exchanger by fixing the following degrees of freedom: - heat transfer area - heat transfer coefficient - effectiveness - hot and cold side pressure changes IDAES documentation reference for heat exchanger NTU model: https://idaes-pse.readthedocs.io/en/latest/reference_guides/model_libraries/generic/unit_models/heat_exchanger_ntu.html This example utilizes the NTU heat exchanger for a shell and tube system. Importing necessary tools ------------------------- First, import the required IDAES and Pyomo modules: .. code:: ipython3 # Import pyomo package from pyomo.environ import ConcreteModel, Constraint, value, units # Import idaes logger to set output levels import idaes.logger as idaeslog # Import the main FlowsheetBlock from IDAES. The flowsheet block will contain the unit model from idaes.core import FlowsheetBlock # Import the MEA property package to create a property block for the flowsheet from idaes.models_extra.column_models.properties.MEA_solvent \ import configuration as aqueous_mea from idaes.models.properties.modular_properties.base.generic_property import ( GenericParameterBlock) # Import the degrees_of_freedom function from the idaes.core.util.model_statistics package from idaes.core.util.model_statistics import degrees_of_freedom # Import the default IPOPT solver from idaes.core.solvers import get_solver # Import a heat exchanger unit from idaes.models.unit_models.heat_exchanger_ntu import HeatExchangerNTU as HXNTU Setting up the flowsheet ------------------------ Then, we will create the ``ConcreteModel`` foundation, attach the steady state flowsheet, and declare the property parameter block that will used for the hot and cold sides. More information on this general workflow can be found here: https://idaes-pse.readthedocs.io/en/stable/how_to_guides/workflow/general.html .. code:: ipython3 m = ConcreteModel() m.fs = FlowsheetBlock(dynamic=False) m.fs.hotside_properties = GenericParameterBlock(**aqueous_mea) m.fs.coldside_properties = GenericParameterBlock(**aqueous_mea) Next, we specify hot and cold streams; note that heat flow from the hot to the cold stream is calculated using the transfer unit method for the ``HXNTU`` unit model. The unit is created below: .. code:: ipython3 # Create an instance of the heat exchanger unit, attaching it to the flowsheet # Specify that the property package to be used with the heat exchanger is the one we created earlier. m.fs.heat_exchanger = HXNTU( hot_side_name="shell", cold_side_name="tube", shell={"property_package": m.fs.hotside_properties, "has_pressure_change": True}, tube={"property_package": m.fs.coldside_properties, "has_pressure_change": True}) # Call the degrees_of_freedom function, get initial DOF DOF_initial = degrees_of_freedom(m) print("The initial DOF is {0}".format(DOF_initial)) .. parsed-literal:: The initial DOF is 15 .. code:: ipython3 assert DOF_initial == 15 Fixing input specifications --------------------------- For this problem, we will fix the inlet conditions, re-evaluate the degrees of freedom to ensure the problem is square (i.e. DOF=0), and run two different options for unit specifications: .. code:: ipython3 m.fs.heat_exchanger.hot_side_inlet.flow_mol[0].fix(60.54879) # mol/s m.fs.heat_exchanger.hot_side_inlet.temperature[0].fix(392.23) # K m.fs.heat_exchanger.hot_side_inlet.pressure[0].fix(202650) # Pa m.fs.heat_exchanger.hot_side_inlet.mole_frac_comp[0, "CO2"].fix(0.0158) # dimensionless m.fs.heat_exchanger.hot_side_inlet.mole_frac_comp[0, "H2O"].fix(0.8747) # dimensionless m.fs.heat_exchanger.hot_side_inlet.mole_frac_comp[0, "MEA"].fix(0.1095) # dimensionless m.fs.heat_exchanger.cold_side_inlet.flow_mol[0].fix(63.01910) # mol/s m.fs.heat_exchanger.cold_side_inlet.temperature[0].fix(326.36) # K m.fs.heat_exchanger.cold_side_inlet.pressure[0].fix(202650) # Pa m.fs.heat_exchanger.cold_side_inlet.mole_frac_comp[0, "CO2"].fix(0.0414) # dimensionless m.fs.heat_exchanger.cold_side_inlet.mole_frac_comp[0, "H2O"].fix(0.8509) # dimensionless m.fs.heat_exchanger.cold_side_inlet.mole_frac_comp[0, "MEA"].fix(0.1077) # dimensionless DOF_initial = degrees_of_freedom(m) print("The DOF is {0}".format(DOF_initial)) .. parsed-literal:: The DOF is 3 Option 1: Fix heat transfer coefficient (HTC) and dimensions of each domain ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Below, we fix the heat exchanger area, effectiveness factor and heat transfer coefficient, which yields a fully defined problem with zero degrees of freedom that may be initialized and solved: .. code:: ipython3 # Unit design variables m.fs.heat_exchanger.area.fix(100) # m2 m.fs.heat_exchanger.heat_transfer_coefficient.fix(200) # W/m2/K m.fs.heat_exchanger.effectiveness.fix(0.7) # dimensionless m.fs.heat_exchanger.hot_side.deltaP.fix(-2000) # Pa m.fs.heat_exchanger.cold_side.deltaP.fix(-2000) # Pa DOF_final = degrees_of_freedom(m) print("The DOF is {0}".format(DOF_final)) .. parsed-literal:: The DOF is 0 .. code:: ipython3 assert DOF_final == 0 Now that the problem is square (zero degress of freedom), we can initialize and solve the full model: .. code:: ipython3 # Initialize the flowsheet, and set the output at INFO m.fs.heat_exchanger.initialize(outlvl=idaeslog.INFO) # Solve the simulation using ipopt # Note: If the degrees of freedom = 0, we have a square problem opt = get_solver() solve_status = opt.solve(m, tee = True) # Display a readable report m.fs.heat_exchanger.report() .. parsed-literal:: 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.hot_side.properties_in: Starting initialization 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.hot_side.properties_in: State variable initialization completed. 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.hot_side.properties_in: Property initialization: optimal - Optimal Solution Found. 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.hot_side.properties_out: Starting initialization 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.hot_side.properties_out: State variable initialization completed. 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.hot_side.properties_out: Property initialization: optimal - Optimal Solution Found. 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.hot_side: Initialization Complete 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.cold_side.properties_in: Starting initialization 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.cold_side.properties_in: State variable initialization completed. 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.cold_side.properties_in: Property initialization: optimal - Optimal Solution Found. 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.cold_side.properties_out: Starting initialization 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.cold_side.properties_out: State variable initialization completed. 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.cold_side.properties_out: Property initialization: optimal - Optimal Solution Found. 2023-03-04 01:48:06 [INFO] idaes.init.fs.heat_exchanger.cold_side: Initialization Complete 2023-03-04 01:48:07 [INFO] idaes.init.fs.heat_exchanger: Initialization Completed, optimal - Optimal Solution Found Ipopt 3.13.2: nlp_scaling_method=gradient-based tol=1e-06 max_iter=200 ****************************************************************************** This program contains Ipopt, a library for large-scale nonlinear optimization. Ipopt is released as open source code under the Eclipse Public License (EPL). For more information visit http://projects.coin-or.org/Ipopt This version of Ipopt was compiled from source code available at https://github.com/IDAES/Ipopt as part of the Institute for the Design of Advanced Energy Systems Process Systems Engineering Framework (IDAES PSE Framework) Copyright (c) 2018-2019. See https://github.com/IDAES/idaes-pse. This version of Ipopt was compiled using HSL, a collection of Fortran codes for large-scale scientific computation. All technical papers, sales and publicity material resulting from use of the HSL codes within IPOPT must contain the following acknowledgement: HSL, a collection of Fortran codes for large-scale scientific computation. See http://www.hsl.rl.ac.uk. ****************************************************************************** This is Ipopt version 3.13.2, running with linear solver ma27. Number of nonzeros in equality constraint Jacobian...: 575 Number of nonzeros in inequality constraint Jacobian.: 0 Number of nonzeros in Lagrangian Hessian.............: 446 Total number of variables............................: 122 variables with only lower bounds: 28 variables with lower and upper bounds: 76 variables with only upper bounds: 0 Total number of equality constraints.................: 122 Total number of inequality constraints...............: 0 inequality constraints with only lower bounds: 0 inequality constraints with lower and upper bounds: 0 inequality constraints with only upper bounds: 0 iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls 0 0.0000000e+00 4.33e+02 1.00e+00 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0 1 0.0000000e+00 8.16e+00 2.63e+01 -1.0 1.01e-02 - 9.90e-01 9.81e-01h 1 2 0.0000000e+00 7.34e-02 3.03e+01 -1.0 4.29e-04 - 9.90e-01 9.91e-01h 1 3 0.0000000e+00 1.23e-09 4.64e+02 -1.0 3.96e-06 - 9.91e-01 1.00e+00h 1 Number of Iterations....: 3 (scaled) (unscaled) Objective...............: 0.0000000000000000e+00 0.0000000000000000e+00 Dual infeasibility......: 0.0000000000000000e+00 0.0000000000000000e+00 Constraint violation....: 2.8492324393151369e-12 1.2266809790162370e-09 Complementarity.........: 0.0000000000000000e+00 0.0000000000000000e+00 Overall NLP error.......: 2.8492324393151369e-12 1.2266809790162370e-09 Number of objective function evaluations = 4 Number of objective gradient evaluations = 4 Number of equality constraint evaluations = 4 Number of inequality constraint evaluations = 0 Number of equality constraint Jacobian evaluations = 4 Number of inequality constraint Jacobian evaluations = 0 Number of Lagrangian Hessian evaluations = 3 Total CPU secs in IPOPT (w/o function evaluations) = 0.002 Total CPU secs in NLP function evaluations = 0.000 EXIT: Optimal Solution Found. ==================================================================================== Unit : fs.heat_exchanger Time: 0.0 ------------------------------------------------------------------------------------ Stream Table Units Hot Inlet Hot Outlet Cold Inlet Cold Outlet Total Molar Flowrate mole / second 60.549 60.549 63.019 63.019 Total Mole Fraction H2O dimensionless 0.87470 0.87470 0.85090 0.85090 Total Mole Fraction MEA dimensionless 0.10950 0.10950 0.10770 0.10770 Total Mole Fraction CO2 dimensionless 0.015800 0.015800 0.041400 0.041400 Temperature kelvin 392.23 344.00 326.36 374.33 Pressure pascal 2.0265e+05 2.0065e+05 2.0265e+05 2.0065e+05 ==================================================================================== .. code:: ipython3 from pyomo.environ import assert_optimal_termination import pytest # Check if termination condition is optimal assert_optimal_termination(solve_status) assert value(m.fs.heat_exchanger.hot_side_outlet.temperature[0]) == pytest.approx(344.00, abs=1e-2) assert value(m.fs.heat_exchanger.cold_side_outlet.temperature[0]) == pytest.approx(374.33, abs=1e-2) Option 2: Unfix shell length and fix shell outlet temperatures ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ In the previous case, we fixed the heat exchanger area and overall heat transfer coefficient. However, given that the models in IDAES are equation oriented, we can fix the outlet variables. For example, we can fix the hot outlet temperature and heat exchanger length, and solve for the heat exchanger area and effectiveness factor that will satisfy that condition. .. code:: ipython3 m.fs.heat_exchanger.effectiveness.unfix() m.fs.heat_exchanger.area.unfix() m.fs.heat_exchanger.hot_side_outlet.temperature.fix(344.0) DOF_final = degrees_of_freedom(m) print("The DOF is {0}".format(DOF_final)) .. parsed-literal:: The DOF is 0 .. code:: ipython3 result = opt.solve(m) print(result) # Display a readable report m.fs.heat_exchanger.report() .. parsed-literal:: Problem: - Lower bound: -inf Upper bound: inf Number of objectives: 1 Number of constraints: 122 Number of variables: 122 Sense: unknown Solver: - Status: ok Message: Ipopt 3.13.2\x3a Optimal Solution Found Termination condition: optimal Id: 0 Error rc: 0 Time: 0.015180826187133789 Solution: - number of solutions: 0 number of solutions displayed: 0 ==================================================================================== Unit : fs.heat_exchanger Time: 0.0 ------------------------------------------------------------------------------------ Stream Table Units Hot Inlet Hot Outlet Cold Inlet Cold Outlet Total Molar Flowrate mole / second 60.549 60.549 63.019 63.019 Total Mole Fraction H2O dimensionless 0.87470 0.87470 0.85090 0.85090 Total Mole Fraction MEA dimensionless 0.10950 0.10950 0.10770 0.10770 Total Mole Fraction CO2 dimensionless 0.015800 0.015800 0.041400 0.041400 Temperature kelvin 392.23 344.00 326.36 374.33 Pressure pascal 2.0265e+05 2.0065e+05 2.0265e+05 2.0065e+05 ==================================================================================== .. code:: ipython3 # Check if termination condition is optimal assert_optimal_termination(result) assert value(m.fs.heat_exchanger.area) == pytest.approx(100, abs=1e-2) assert value(m.fs.heat_exchanger.effectiveness[0]) == pytest.approx(0.7, abs=1e-2) assert value(m.fs.heat_exchanger.cold_side_outlet.temperature[0]) == pytest.approx( 374.33, rel=1e-3)