##############################################################################
# Institute for the Design of Advanced Energy Systems Process Systems
# Engineering Framework (IDAES PSE Framework) Copyright (c) 2018, 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. All rights reserved.
#
# Please see the files COPYRIGHT.txt and LICENSE.txt for full copyright and
# license information, respectively. Both files are also available online
# at the URL "https://github.com/IDAES/idaes".
##############################################################################
"""
Basic IDAES 1D Heat Exchanger Model.
Single pass shell and tube HX model with 0D/1D wall conduction model
"""
from __future__ import division
# Import Python libraries
import logging
import math
# Import Pyomo libraries
from pyomo.environ import SolverFactory, Var, Param, Constraint, Expression,\
value, TransformationFactory
from pyomo.opt import TerminationCondition
from pyomo.common.config import ConfigValue, In
from pyomo.dae import DerivativeVar, ContinuousSet
# Import IDAES cores
from idaes.core import UnitBlockData, declare_process_block_class, \
Holdup1D, CONFIG_Base_1D
from idaes.core.util.config import is_parameter_block, list_of_strings
from idaes.core.util.misc import add_object_ref
__author__ = "Jaffer Ghouse"
# Set up logger
logger = logging.getLogger('idaes.unit_model')
[docs]@declare_process_block_class("HeatExchanger1D")
class HeatExchanger1dData(UnitBlockData):
"""Standard Heat Exchanger 1D Unit Model Class."""
CONFIG = CONFIG_Base_1D()
# Set default values of inherited attributes
CONFIG.get("has_heat_transfer")._default = True
# Removing redundant phase equilibrium attribute
# Workaround for now - del attribute being added to Pyutilib by Pyomo
CONFIG._declared.remove("has_phase_equilibrium")
CONFIG._decl_order.remove("has_phase_equilibrium")
del CONFIG._data["has_phase_equilibrium"]
# Add unit model attributes
CONFIG.declare("shell_has_phase_equilibrium", ConfigValue(
default=True,
domain=In([True, False]),
description="Phase equilibrium for shell side",
doc="""Argument to enable phase equilibrium on the shell side.
- True - include phase equilibrium term
- False - do not include phase equilinrium term"""))
CONFIG.declare("tube_has_phase_equilibrium", ConfigValue(
default=True,
domain=In([True, False]),
description="Phase equilibrium for tube side",
doc="""Argument to enable phase equilibrium on the tube side.
- True - include phase equilibrium term
- False - do not include phase equilinrium term"""))
CONFIG.declare("shell_has_mass_diffusion", ConfigValue(
default=False,
domain=In([True, False]),
description="Mass diffusion flag",
doc="""Flag indicating whether mass diffusion/dispersion should be
included in material balance equations (default=False)"""))
CONFIG.declare("shell_has_energy_diffusion", ConfigValue(
default=False,
domain=In([True, False]),
description="Energy diffusion flag",
doc="""Flag indicating whether energy diffusion/dispersion should be
included in energy balance equations (default=False)"""))
CONFIG.declare("tube_has_mass_diffusion", ConfigValue(
default=False,
domain=In([True, False]),
description="Mass diffusion flag",
doc="""Flag indicating whether mass diffusion/dispersion should be
included in material balance equations (default=False)"""))
CONFIG.declare("tube_has_energy_diffusion", ConfigValue(
default=False,
domain=In([True, False]),
description="Energy diffusion flag",
doc="""Flag indicating whether energy diffusion/dispersion should be
included in energy balance equations (default=False)"""))
CONFIG.declare("shell_property_package", ConfigValue(
default=None,
domain=is_parameter_block,
description="Property package to use for shell holdup",
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("shell_property_package_args", ConfigValue(
default={},
description="Arguments for constructing shell 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)"""))
CONFIG.declare("tube_property_package", ConfigValue(
default=None,
domain=is_parameter_block,
description="Property package to use for tube holdup",
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("tube_property_package_args", ConfigValue(
default={},
description="Arguments for constructing tube 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)"""))
CONFIG.declare("shell_inlet_list", ConfigValue(
domain=list_of_strings,
description="List of inlet names for shell side of heat exchanger",
doc="""A list containing names of inlets for shell of heat exchanger
(default = None)
- None - default single inlet
- list - a list of names for inlets"""))
CONFIG.declare("shell_num_inlets", ConfigValue(
domain=int,
description="Number of inlets to shell side",
doc="""Argument indication number (int) of inlets to construct
(default = None). Not used if inlet_list arg is provided.
- None - use inlet_list arg instead
- int - Inlets will be named with sequential numbers from 1
to num_inlets"""))
CONFIG.declare("shell_outlet_list", ConfigValue(
domain=list_of_strings,
description="List of outlet names for shell side of heat exchanger",
doc="""A list containing names of outlets for shell of heat exchanger
(default = None)
- None - default single inlet
- list - a list of names for inlets"""))
CONFIG.declare("shell_num_outlets", ConfigValue(
domain=int,
description="Number of outlets to shell side",
doc="""Argument indication number (int) of outlets to construct
(default = None). Not used if outlet_list arg is provided.
- None - use outlet_list arg instead
- int - Outlets will be named with sequential numbers from 1
to num_outlets"""))
CONFIG.declare("tube_inlet_list", ConfigValue(
domain=list_of_strings,
description="List of inlet names for tube side of heat exchanger",
doc="""A list containing names of inlets for tube of heat exchanger
(default = None)
- None - default single inlet
- list - a list of names for inlets"""))
CONFIG.declare("tube_num_inlets", ConfigValue(
domain=int,
description="Number of inlets to tube side",
doc="""Argument indication number (int) of inlets to construct
(default = None). Not used if inlet_list arg is provided.
- None - use inlet_list arg instead
- int - Inlets will be named with sequential numbers from 1
to num_inlets"""))
CONFIG.declare("tube_outlet_list", ConfigValue(
domain=list_of_strings,
description="List of outlet names for tube side of heat exchanger",
doc="""A list containing names of outlets for tube of heat exchanger
(default = None)
- None - default single inlet
- list - a list of names for inlets"""))
CONFIG.declare("tube_num_outlets", ConfigValue(
domain=int,
description="Number of outlets to tube side",
doc="""Argument indication number (int) of outlets to construct
(default = None). Not used if outlet_list arg is provided.
- None - use outlet_list arg instead
- int - Outlets will be named with sequential numbers from 1
to num_outlets"""))
# TODO : We should probably think about adding a consistency check for the
# TODO : discretisation methdos as well.
CONFIG.declare("shell_discretization_method", ConfigValue(
default="OCLR",
domain=In(['OCLR', 'OCLL', 'BFD', 'FFD']),
description="Discretization method to apply for length domain",
doc="""Method to be used by DAE transformation when discretizing length
domain on shell side (default = `OCLR`).
- 'OCLR' - orthogonal collocation (Radau roots)
- 'OCLL' - orthogonal collocation (Legendre roots)
- 'BFD' - backwards finite difference (1st order)
- 'FFD' - forwards finite difference (1st order)"""))
CONFIG.declare("tube_discretization_method", ConfigValue(
default="OCLR",
domain=In(['OCLR', 'OCLL', 'BFD', 'FFD']),
description="Discretization method to apply for length domain",
doc="""Method to be used by DAE transformation when discretizing
length
domain on tube side (default = `OCLR`).
- 'OCLR' - orthogonal collocation (Radau roots)
- 'OCLL' - orthogonal collocation (Legendre roots)
- 'BFD' - backwards finite difference (1st order)
- 'FFD' - forwards finite difference (1st order)"""))
CONFIG.declare("finite_elements", ConfigValue(
default=20,
domain=int,
description="Number of finite elements length domain",
doc="""Number of finitie elements to use when discretizing length
domain (default=20)"""))
CONFIG.declare("collocation_points", ConfigValue(
default=3,
domain=int,
description="Number of collocation points per finite element",
doc="""Number of collocation points to use per finite element when
discretizing length domain (default=3)"""))
CONFIG.declare("flow_type", ConfigValue(
default="co_current",
domain=In(['co_current', 'counter_current']),
description="Flow configuration of heat exchanger",
doc="""Flow configuration of heat exchanger
- co_current: shell and tube flows from 0 to 1
- counter_current: shell side flows from 0 to 1
tube side flows from 1 to 0"""))
CONFIG.declare("has_wall_conduction", ConfigValue(
default="none",
domain=In(["none", "1D", "2D"]),
description="Conduction model for tube wall",
doc="""Argument to enable type of wall heat conduction model.
- none - 0D wall model
- 1D - 1D wall model along the thickness of the tube
- 2D - 2D wall model along the lenghth and thickness
of the tube"""))
[docs] def build(self):
"""
Begin building model (pre-DAE transformation).
Args:
None
Returns:
None
"""
# Call UnitModel.build to setup dynamics
super(HeatExchanger1dData, self).build()
# Set flow directions for holdup blocks
if self.config.flow_type is "co_current":
set_direction_shell = "forward"
set_direction_tube = "forward"
else:
set_direction_shell = "forward"
set_direction_tube = "backward"
# Holdup 1D for shell
self.shell = Holdup1D(
has_phase_equilibrium=self.config.shell_has_phase_equilibrium,
has_mass_diffusion=self.config.shell_has_mass_diffusion,
has_energy_diffusion=self.config.shell_has_energy_diffusion,
property_package=self.config.shell_property_package,
property_package_args=self.config.shell_property_package_args,
flow_direction=set_direction_shell,
discretization_method=self.config.shell_discretization_method,
finite_elements=self.config.finite_elements,
collocation_points=self.config.collocation_points)
# Holdup 1D for tube
self.tube = Holdup1D(
has_phase_equilibrium=self.config.tube_has_phase_equilibrium,
has_mass_diffusion=self.config.tube_has_mass_diffusion,
has_energy_diffusion=self.config.tube_has_energy_diffusion,
property_package=self.config.tube_property_package,
property_package_args=self.config.tube_property_package_args,
flow_direction=set_direction_tube,
discretization_method=self.config.tube_discretization_method,
finite_elements=self.config.finite_elements,
collocation_points=self.config.collocation_points)
# Set Unit Geometry and holdup Volume
self._set_geometry()
# Construct the unit model
self._make_performance()
[docs] def post_transform_build(self):
"""
Continue model construction after DAE transformation.
Args:
None
Returns:
None
"""
# Construct Inlets
self.build_inlets(holdup='shell',
inlets=self.config.shell_inlet_list,
num_inlets=self.config.shell_num_inlets)
self.build_inlets(holdup='tube',
inlets=self.config.tube_inlet_list,
num_inlets=self.config.tube_num_inlets)
# Build Outlets
self.build_outlets(holdup='shell',
outlets=self.config.shell_outlet_list,
num_outlets=self.config.shell_num_outlets)
self.build_outlets(holdup='tube',
outlets=self.config.tube_outlet_list,
num_outlets=self.config.tube_num_outlets)
def _set_geometry(self):
"""
Define the geometry of the unit and link to holdup volume.
Args:
None
Returns:
None
"""
# TODO : For consistency with arguments, shell_area might be better.
# Create references to all geometry variables
# Shell side
add_object_ref(self, "shell_area", self.shell.area)
add_object_ref(self, "shell_length", self.shell.length)
# Tube side
add_object_ref(self, "tube_area", self.tube.area)
add_object_ref(self, "tube_length", self.tube.length)
def _make_performance(self):
"""
Constraints for unit model.
Args:
None
Returns:
None
"""
# Unit model parameters
self.pi = Param(initialize=math.pi, doc="pi")
# Unit model variables
# HX dimensions
self.d_shell = Var(initialize=1, doc="Diameter of shell")
self.d_tube_outer = Var(initialize=0.011, doc="Outer diameter of tube")
self.d_tube_inner = Var(initialize=0.010, doc="Inner diameter of tube")
self.N_tubes = Var(initialize=1, doc="Number of tubes")
# Note: In addition to the above variables, "shell_length" and
# "tube_length" need to be fixed at the flowsheet level
# Performance variables
self.shell_heat_transfer_coefficient = Var(self.time,
self.shell.ldomain,
initialize=50,
doc="Heat transfer "
"coefficient")
self.tube_heat_transfer_coefficient = Var(self.time,
self.tube.ldomain,
initialize=50,
doc="Heat transfer "
"coefficient")
# Wall 0D model (Q_in = Q_out)
if self.config.has_wall_conduction is "none":
self.temperature_wall = Var(self.time, self.shell.ldomain,
initialize=298.15)
# Performance equations
# Energy transfer between shell and tube wall
@self.Constraint(self.time, self.tube.ldomain,
doc="Heat transfer between shell and tube")
def shell_heat_transfer_eq(self, t, x):
return self.shell.heat[t, x] == - self.N_tubes *\
(self.shell_heat_transfer_coefficient[t, x] *
self.pi * self.d_tube_inner *
(self.shell.properties[t, x].temperature -
self.temperature_wall[t, x]))
# Energy transfer between tube wall and tube
@self.Constraint(self.time, self.shell.ldomain,
doc="Convective heat transfer")
def tube_heat_transfer_eq(self, t, x):
return self.tube.heat[t, x] == \
self.tube_heat_transfer_coefficient[t, x] *\
self.pi * self.d_tube_inner * \
(self.temperature_wall[t, x] -
self.tube.properties[t, x].temperature)
# Wall 0D model
if self.config.has_wall_conduction is "none":
@self.Constraint(self.time, self.shell.ldomain,
doc="wall 0D model")
def wall_0D_model(self, t, x):
return self.tube.heat[t, x] == -(self.shell.heat[t, x] /
self.N_tubes)
# wall 1D model (1D across thickness of wall)
if self.config.has_wall_conduction is "1D":
# Parameters specific for 1D wall model
self.xdomain = ContinuousSet(bounds=(0, 1),
doc="domain for wall thickness")
# Variables specific for 1D wall model
self.therm_cond_wall = Var(initialize=28.5,
doc="thermal conductivity of wall")
self.dens_wall = Var(initialize=7880,
doc="density of wall")
self.cp_wall = Var(initialize=741,
doc="specific heat capacity of wall")
self.temperature_wall = Var(self.time, self.shell.ldomain,
self.xdomain, initialize=298.15,
doc="wall temperature")
self.dTdx = DerivativeVar(self.temperature_wall,
wrt=self.xdomain,
doc="""first derivative of wall temperature
with respect to xdomain""")
self.dTdx2 = DerivativeVar(self.temperature_wall,
wrt=(self.xdomain, self.xdomain),
doc="""second derivative of wall
temperature with respect
to xdomain""")
# Discretisation for the wall along the thickness (xdomain)
self.wall_discretizer = \
TransformationFactory('dae.finite_difference')
self.nfe = Param(initialize=5,
doc="number of finite elements")
self.wall_discretizer.apply_to(self, nfe=self.nfe.value,
wrt=self.xdomain,
scheme='CENTRAL')
# Expression
self.thickness_wall = Expression(expr=0.5*(self.d_tube_outer -
self.d_tube_inner))
# Constraints
@self.Constraint(self.time, self.shell.ldomain, self.xdomain,
doc="wall 1D model across the thickness")
def wall_1D_model(self, t, z, x):
if x == 0 or x == 1:
return Constraint.Skip
return 0 == (self.therm_cond_wall/(self.thickness_wall**2 *
self.dens_wall*self.cp_wall)) *\
self.dTdx2[t, z, x]
# Boundary condition
# At x=0
@self.Constraint(self.time, self.shell.ldomain, self.xdomain,
doc="wall 1D model across the thickness")
def wall_1D_outer_BC(self, t, z, x):
if x == 0:
return -(self.therm_cond_wall/self.thickness_wall) * \
(self.temperature_wall[t, z, self.xdomain.next(0)] -
self.temperature_wall[t, z, 0]) == (1/self.nfe) *\
self.shell_heat_transfer_coefficient[t, z] *\
(self.shell.properties[t, z].temperature -
self.temperature_wall[t, z, 0])
return Constraint.Skip
# At x=L
@self.Constraint(self.time, self.tube.ldomain, self.xdomain,
doc="wall 1D model across the thickness")
def wall_1D_inner_BC(self, t, z, x):
if x == 1:
return -(self.therm_cond_wall/self.thickness_wall) * \
(self.temperature_wall[t, z, 1] -
self.temperature_wall[t, z,
self.xdomain.prev(1)]) == \
(1/self.nfe) *\
self.tube_heat_transfer_coefficient[t, z] * \
(self.temperature_wall[t, z, 1] -
self.tube.properties[t, z].temperature)
return Constraint.Skip
# Energy transfer between shell and tube wall
@self.Constraint(self.time, self.tube.ldomain,
doc="Heat transfer between shell and tube")
def shell_heat_transfer_eq(self, t, z):
return self.shell.heat[t, z] == - self.N_tubes *\
(self.shell_heat_transfer_coefficient[t, z] * self.pi *
self.d_tube_outer*(self.shell.properties[t, z].
temperature -
self.temperature_wall[t, z, 0]))
# Energy transfer between tube wall and tube
@self.Constraint(self.time, self.shell.ldomain,
doc="Convective heat transfer")
def tube_heat_transfer_eq(self, t, z):
return self.tube.heat[t, z] == \
self.tube_heat_transfer_coefficient[t, z] *\
self.pi * self.d_tube_inner * \
(self.temperature_wall[t, z, 1] -
self.tube.properties[t, z].temperature)
# Define tube area in terms of tube diameter
self.area_calc_tube = Constraint(expr=4*self.tube_area == self.pi *
self.d_tube_inner**2)
# Define shell area in terms of shell and tube diameter
self.area_calc_shell = Constraint(expr=4*self.shell_area ==
self.pi * (self.d_shell**2 -
self.N_tubes *
self.d_tube_outer**2))
[docs] def model_check(blk):
"""
Model checks for unit (call model check on holdups)
Args:
None
Returns:
None
"""
# Run holdup block model checks
blk.shell.model_check()
blk.tube.model_check()
[docs] def initialize(blk, shell_state_args={}, tube_state_args={}, outlvl=0,
solver='ipopt', optarg={'tol': 1e-6}):
"""
Initialisation routine for isothermal unit (default solver ipopt).
Keyword Arguments:
state_args : a dict of arguments to be passed to the property
package(s) to provide an initial state for
initialization (see documentation of the specific
property package) (default = {}).
outlvl : sets output level of initialisation routine
* 0 = no output (default)
* 1 = return solver state for each step in routine
* 2 = return solver state for each step in subroutines
* 3 = include solver output infomation (tee=True)
optarg : solver options dictionary object (default={'tol': 1e-6})
solver : str indicating whcih solver to use during
initialization (default = 'ipopt')
Returns:
None
"""
# Set solver options
if outlvl > 3:
stee = True
else:
stee = False
opt = SolverFactory(solver)
opt.options = optarg
# ---------------------------------------------------------------------
# Initialize shell block
flags = blk.shell.initialize(outlvl=outlvl-1,
optarg=optarg,
solver=solver,
state_args=shell_state_args)
flags = blk.tube.initialize(outlvl=outlvl-1,
optarg=optarg,
solver=solver,
state_args=tube_state_args)
if outlvl > 0:
logger.info('{} Initialisation Step 1 Complete.'.format(blk.name))
# ---------------------------------------------------------------------
# Solve unit
# Wall 0D
if blk.config.has_wall_conduction is "none":
for t in blk.time:
for z in blk.shell.ldomain:
blk.temperature_wall[t, z].fix(value(
0.5*(blk.shell.properties[t, 0].temperature +
blk.tube.properties[t, 0].temperature)))
blk.tube.deactivate()
blk.tube_heat_transfer_eq.deactivate()
blk.wall_0D_model.deactivate()
results = opt.solve(blk, tee=stee)
if outlvl > 0:
if results.solver.termination_condition \
== TerminationCondition.optimal:
logger.info('{} Initialisation Step 2 Complete.'
.format(blk.name))
else:
logger.warning('{} Initialisation Step 2 Failed.'
.format(blk.name))
blk.tube.activate()
blk.tube_heat_transfer_eq.activate()
results = opt.solve(blk, tee=stee)
if outlvl > 0:
if results.solver.termination_condition \
== TerminationCondition.optimal:
logger.info('{} Initialisation Step 3 Complete.'
.format(blk.name))
else:
logger.warning('{} Initialisation Step 3 Failed.'
.format(blk.name))
blk.wall_0D_model.activate()
blk.temperature_wall.unfix()
results = opt.solve(blk, tee=stee)
if outlvl > 0:
if results.solver.termination_condition \
== TerminationCondition.optimal:
logger.info('{} Initialisation Step 4 Complete.'
.format(blk.name))
else:
logger.warning('{} Initialisation Step 4 Failed.'
.format(blk.name))
# Wall 0D
if blk.config.has_wall_conduction is "1D":
for t in blk.time:
blk.temperature_wall[t, :, :].fix(value(
0.5*(blk.shell.properties[t, 0].temperature +
blk.tube.properties[t, 0].temperature)))
blk.tube.deactivate()
blk.tube_heat_transfer_eq.deactivate()
blk.wall_1D_model.deactivate()
blk.wall_1D_outer_BC.deactivate()
blk.wall_1D_inner_BC.deactivate()
results = opt.solve(blk, tee=stee)
if outlvl > 0:
if results.solver.termination_condition \
== TerminationCondition.optimal:
logger.info('{} Initialisation Step 2 Complete.'
.format(blk.name))
else:
logger.warning('{} Initialisation Step 2 Failed.'
.format(blk.name))
blk.tube.activate()
blk.tube_heat_transfer_eq.activate()
results = opt.solve(blk, tee=stee)
if outlvl > 0:
if results.solver.termination_condition \
== TerminationCondition.optimal:
logger.info('{} Initialisation Step 3 Complete.'
.format(blk.name))
else:
logger.warning('{} Initialisation Step 3 Failed.'
.format(blk.name))
blk.wall_1D_model.activate()
blk.wall_1D_outer_BC.activate()
blk.wall_1D_inner_BC.activate()
blk.temperature_wall[:, :, :].unfix()
results = opt.solve(blk, tee=stee)
if outlvl > 0:
if results.solver.termination_condition \
== TerminationCondition.optimal:
logger.info('{} Initialisation Step 4 Complete.'
.format(blk.name))
else:
logger.warning('{} Initialisation Step 4 Failed.'
.format(blk.name))
# ---------------------------------------------------------------------
# Release Inlet state
blk.shell.release_state(flags, outlvl-1)
blk.tube.release_state(flags, outlvl-1)
if outlvl > 0:
logger.info('{} Initialisation Complete.'.format(blk.name))