##############################################################################
# Institute for the Design of Advanced Energy Systems Process Systems
# Engineering Framework (IDAES PSE Framework) Copyright (c) 2018-2019, 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-pse".
##############################################################################
"""This is an example supercritical pulverized coal (SCPC) power plant steam cycle
model. This model doesn't represent any specific power plant, but it represents
what could be considered a typical SCPC plant, producing around 620 MW gross.
This model is for demonstration and tutorial purposes only. Before looking at the
model, it may be useful to look at the process flow diagram (PFD).
"""
__author__ = "John Eslick"
# Import Python libraries
from collections import OrderedDict
import argparse
import logging
# Import Pyomo libraries
import pyomo.environ as pyo
from pyomo.network import Arc, Port
# IDAES Imports
from idaes.core import FlowsheetBlock # Flowsheet class
from idaes.core.util import model_serializer as ms # load/save model state
from idaes.core.util.misc import svg_tag # place numbers/text in an SVG
from idaes.property_models import iapws95 # steam properties
from idaes.unit_models.power_generation import ( # power generation unit models
TurbineMultistage,
FWH0D,
)
from idaes.unit_models import ( # basic IDAES unit models, and enum
Mixer,
HeatExchanger,
PressureChanger,
MomentumMixingType, # Enum type for mixer pressure calculation selection
)
from idaes.core.util import copy_port_values as _set_port # for model intialization
from idaes.core.util.model_statistics import degrees_of_freedom
from idaes.core.util.tables import create_stream_table_dataframe # as Pandas DataFrame
# Callback used to construct heat exchangers with the Underwood approx. for LMTD
from idaes.unit_models.heat_exchanger import delta_temperature_underwood_callback
# Pressure changer type (e.g. adiabatic, pump, isentropic...)
from idaes.unit_models.pressure_changer import ThermodynamicAssumption
from idaes.logger import getModelLogger, getInitLogger, init_tee, condition
_log = getModelLogger(__name__, logging.INFO)
[docs]def create_model():
"""Create the flowsheet and add unit models. Fixing model inputs is done
in a separate function to try to keep this fairly clean and easy to follow.
Args:
None
Returns:
(ConcreteModel) Steam cycle model
"""
############################################################################
# Flowsheet and Properties #
############################################################################
m = pyo.ConcreteModel(name="Steam Cycle Model")
m.fs = FlowsheetBlock(default={"dynamic": False}) # Add steady state flowsheet
# A physical property parameter block for IAPWS-95 with pressure and enthalpy
# (PH) state variables. Usually pressure and enthalpy state variables are
# more robust especially when the phases are unknown.
m.fs.prop_water = iapws95.Iapws95ParameterBlock(
default={"phase_presentation": iapws95.PhaseType.LG}
)
# A physical property parameter block with temperature, pressure and vapor
# fraction (TPx) state variables. There are a few instances where the vapor
# fraction is known and the temperature and pressure state variables are
# preferable.
m.fs.prop_water_tpx = iapws95.Iapws95ParameterBlock(
default={
"phase_presentation": iapws95.PhaseType.LG,
"state_vars": iapws95.StateVars.TPX,
}
)
############################################################################
# Turbine with fill-in reheat constraints #
############################################################################
# The TurbineMultistage class allows creation of the full turbine model by
# providing several configuration options, including: throttle valves;
# high-, intermediate-, and low-pressure sections; steam extractions; and
# pressure driven flow. See the IDAES documentation for details.
m.fs.turb = TurbineMultistage(
default={
"property_package": m.fs.prop_water,
"num_parallel_inlet_stages": 4, # number of admission arcs
"num_hp": 7, # number of high-pressure stages
"num_ip": 10, # number of intermediate-pressure stages
"num_lp": 11, # number of low-pressure stages
"hp_split_locations": [4, 7], # hp steam extraction locations
"ip_split_locations": [5, 10], # ip steam extraction locations
"lp_split_locations": [4, 8, 10, 11], # lp steam extraction locations
"hp_disconnect": [7], # disconnect hp from ip to insert reheater
"ip_split_num_outlets": {10: 3}, # number of split streams (default is 2)
}
)
# This model is only the steam cycle, and the reheater is part of the boiler.
# To fill in the reheater gap, a few constraints for the flow, pressure drop,
# and outlet temperature are added. A detailed boiler model can be coupled later.
#
# hp_split[7] is the splitter directly after the last HP stage. The splitter
# outlet "outlet_1" is always taken to be the main steam flow through the turbine.
# When the turbine model was instantiated the stream from the HP section to the IP
# section was omitted, so the reheater could be inserted.
# The flow constraint sets flow from outlet_1 of the splitter equal to
# flow into the IP turbine.
@m.fs.turb.Constraint(m.fs.time)
def constraint_reheat_flow(b, t):
return b.ip_stages[1].inlet.flow_mol[t] == b.hp_split[7].outlet_1.flow_mol[t]
# Create a variable for pressure change in the reheater (assuming
# reheat_delta_p should be negative).
m.fs.turb.reheat_delta_p = pyo.Var(m.fs.time, initialize=0)
# Add a constraint to calculate the IP section inlet pressure based on the
# pressure drop in the reheater and the outlet pressure of the HP section.
@m.fs.turb.Constraint(m.fs.time)
def constraint_reheat_press(b, t):
return (
b.ip_stages[1].inlet.pressure[t]
== b.hp_split[7].outlet_1.pressure[t] + b.reheat_delta_p[t]
)
# Create a variable for reheat temperature and fix it to the desired reheater
# outlet temperature
m.fs.turb.reheat_out_T = pyo.Var(m.fs.time, initialize=866)
# Create a constraint for the IP section inlet temperature.
@m.fs.turb.Constraint(m.fs.time)
def constraint_reheat_temp(b, t):
return (
b.ip_stages[1].control_volume.properties_in[t].temperature
== b.reheat_out_T[t]
)
############################################################################
# Add Condenser/hotwell/condensate pump #
############################################################################
# Add a mixer for all the streams coming into the condenser. In this case the
# main steam, and the boiler feed pump turbine outlet go to the condenser
m.fs.condenser_mix = Mixer(
default={
"momentum_mixing_type": MomentumMixingType.none,
"inlet_list": ["main", "bfpt"],
"property_package": m.fs.prop_water,
}
)
# The pressure in the mixer comes from the connection to the condenser. All
# the streams coming in and going out of the mixer are equal, but we created
# the mixer with no calculation for the unit pressure. Here a constraint that
# specifies that the mixer pressure is equal to the main steam pressure is
# added. There is also a constraint that specifies the that BFP turbine outlet
# pressure is the same as the condenser pressure. Combined with the stream
# connections between units, these constraints effectively specify that the
# mixer inlet and outlet streams all have the same pressure.
@m.fs.condenser_mix.Constraint(m.fs.time)
def mixer_pressure_constraint(b, t):
return b.main_state[t].pressure == b.mixed_state[t].pressure
# The condenser model uses the physical property model with TPx state
# variables, while the rest of the model uses PH state variables. To
# translate between the two property calculations, an extra port is added to
# the mixer which contains temperature, pressure, and vapor fraction
# quantities. The first step is to add references to the temperature and
# vapor fraction expressions in the IAPWS-95 property block. The references
# are used to handle time indexing in the ports by using the property blocks
# time index to create references that appear to be time indexed variables.
# These references mirror the references created by the framework automatically
# for the existing ports.
m.fs.condenser_mix._outlet_temperature_ref = pyo.Reference(
m.fs.condenser_mix.mixed_state[:].temperature
)
m.fs.condenser_mix._outlet_vapor_fraction_ref = pyo.Reference(
m.fs.condenser_mix.mixed_state[:].vapor_frac
)
# Add the new port with the state information that needs to go to the
# condenser
m.fs.condenser_mix.outlet_tpx = Port(
initialize={
"flow_mol": m.fs.condenser_mix._outlet_flow_mol_ref,
"temperature": m.fs.condenser_mix._outlet_temperature_ref,
"pressure": m.fs.condenser_mix._outlet_pressure_ref,
"vapor_frac": m.fs.condenser_mix._outlet_vapor_fraction_ref,
}
)
# Add the heat exchanger model for the condenser.
m.fs.condenser = HeatExchanger(
default={
"delta_temperature_callback": delta_temperature_underwood_callback,
"shell": {"property_package": m.fs.prop_water_tpx},
"tube": {"property_package": m.fs.prop_water},
}
)
m.fs.condenser.delta_temperature_out.fix(5)
# Everything condenses so the saturation pressure determines the condenser
# pressure. Deactivate the constraint that is used in the TPx version vapor
# fraction constraint and fix vapor fraction to 0.
m.fs.condenser.shell.properties_out[:].eq_complementarity.deactivate()
m.fs.condenser.shell.properties_out[:].vapor_frac.fix(0)
# There is some subcooling in the condenser, so we assume the condenser
# pressure is actually going to be slightly higher than the saturation
# pressure.
m.fs.condenser.pressure_over_sat = pyo.Var(
m.fs.time,
initialize=500,
doc="Pressure added to Psat in the condeser. This is to account for"
"some subcooling. (Pa)",
)
# Add a constraint for condenser pressure
@m.fs.condenser.Constraint(m.fs.time)
def eq_pressure(b, t):
return (
b.shell.properties_out[t].pressure
== b.shell.properties_out[t].pressure_sat + b.pressure_over_sat[t]
)
# Extra port on condenser to hook back up to pressure-enthalpy properties
m.fs.condenser._outlet_1_enth_mol_ref = pyo.Reference(
m.fs.condenser.shell.properties_out[:].enth_mol
)
m.fs.condenser.outlet_1_ph = Port(
initialize={
"flow_mol": m.fs.condenser._outlet_1_flow_mol_ref,
"pressure": m.fs.condenser._outlet_1_pressure_ref,
"enth_mol": m.fs.condenser._outlet_1_enth_mol_ref,
}
)
# Add the condenser hotwell. In steady state a mixer will work. This is
# where makeup water is added if needed.
m.fs.hotwell = Mixer(
default={
"momentum_mixing_type": MomentumMixingType.none,
"inlet_list": ["condensate", "makeup"],
"property_package": m.fs.prop_water,
}
)
# The hotwell is assumed to be at the same pressure as the condenser.
@m.fs.hotwell.Constraint(m.fs.time)
def mixer_pressure_constraint(b, t):
return b.condensate_state[t].pressure == b.mixed_state[t].pressure
# Condensate pump
m.fs.cond_pump = PressureChanger(
default={
"property_package": m.fs.prop_water,
"thermodynamic_assumption": ThermodynamicAssumption.pump,
}
)
############################################################################
# Add low pressure feedwater heaters #
############################################################################
# All the feedwater heater sections will be set to use the Underwood
# approximation for LMTD, so create the fwh_config dict to make the config
# slightly cleaner
fwh_config = {"delta_temperature_callback": delta_temperature_underwood_callback}
# The feedwater heater model allows feedwater heaters with a desuperheat,
# condensing, and subcooling section to be added an a reasonably simple way.
# See the IDAES documentation for more information of configuring feedwater
# heaters
m.fs.fwh1 = FWH0D(
default={
"has_desuperheat": False,
"has_drain_cooling": False,
"has_drain_mixer": True,
"property_package": m.fs.prop_water,
"condense": fwh_config,
}
)
# pump for fwh1 condensate, to pump it ahead and mix with feedwater
m.fs.fwh1_pump = PressureChanger(
default={
"property_package": m.fs.prop_water,
"thermodynamic_assumption": ThermodynamicAssumption.pump,
}
)
# Mix the FWH1 drain back into the feedwater
m.fs.fwh1_return = Mixer(
default={
"momentum_mixing_type": MomentumMixingType.none,
"inlet_list": ["feedwater", "fwh1_drain"],
"property_package": m.fs.prop_water,
}
)
# Set the mixer pressure to the feedwater pressure
@m.fs.fwh1_return.Constraint(m.fs.time)
def mixer_pressure_constraint(b, t):
return b.feedwater_state[t].pressure == b.mixed_state[t].pressure
# Add the rest of the low pressure feedwater heaters
m.fs.fwh2 = FWH0D(
default={
"has_desuperheat": True,
"has_drain_cooling": True,
"has_drain_mixer": True,
"property_package": m.fs.prop_water,
"desuperheat": fwh_config,
"cooling": fwh_config,
"condense": fwh_config,
}
)
m.fs.fwh3 = FWH0D(
default={
"has_desuperheat": True,
"has_drain_cooling": True,
"has_drain_mixer": True,
"property_package": m.fs.prop_water,
"desuperheat": fwh_config,
"cooling": fwh_config,
"condense": fwh_config,
}
)
m.fs.fwh4 = FWH0D(
default={
"has_desuperheat": True,
"has_drain_cooling": True,
"has_drain_mixer": False,
"property_package": m.fs.prop_water,
"desuperheat": fwh_config,
"cooling": fwh_config,
"condense": fwh_config,
}
)
############################################################################
# Add deaerator and boiler feed pump (BFP) #
############################################################################
# The deaerator is basically an open tank with multiple inlets. For steady-
# state, a mixer model is sufficient.
m.fs.fwh5_da = Mixer(
default={
"momentum_mixing_type": MomentumMixingType.none,
"inlet_list": ["steam", "drain", "feedwater"],
"property_package": m.fs.prop_water,
}
)
@m.fs.fwh5_da.Constraint(m.fs.time)
def mixer_pressure_constraint(b, t):
# Not sure about deaerator pressure, so assume same as feedwater inlet
return b.feedwater_state[t].pressure == b.mixed_state[t].pressure
# Add the boiler feed pump and boiler feed pump turbine
m.fs.bfp = PressureChanger(
default={
"property_package": m.fs.prop_water,
"thermodynamic_assumption": ThermodynamicAssumption.pump,
}
)
m.fs.bfpt = PressureChanger(
default={
"property_package": m.fs.prop_water,
"compressor": False,
"thermodynamic_assumption": ThermodynamicAssumption.isentropic,
}
)
# The boiler feed pump outlet pressure is the same as the condenser
@m.fs.Constraint(m.fs.time)
def constraint_out_pressure(b, t):
return (
b.bfpt.control_volume.properties_out[t].pressure
== b.condenser.shell.properties_out[t].pressure
)
# Instead of specifying a fixed efficiency, specify that the steam is just
# starting to condense at the outlet of the boiler feed pump turbine. This
# ensures approximately the right behavior in the turbine. With a fixed
# efficiency, depending on the conditions you can get odd things like steam
# fully condensing in the turbine.
@m.fs.Constraint(m.fs.time)
def constraint_out_enthalpy(b, t):
return (
b.bfpt.control_volume.properties_out[t].enth_mol
== b.bfpt.control_volume.properties_out[t].enth_mol_sat_phase["Vap"] - 200
)
# The boiler feed pump power is the same as the power generated by the
# boiler feed pump turbine. This constraint determines the steam flow to the
# BFP turbine. The turbine work is negative for power out, while pump work
# is positive for power in.
@m.fs.Constraint(m.fs.time)
def constraint_bfp_power(b, t):
return 0 == b.bfp.control_volume.work[t] + b.bfpt.control_volume.work[t]
############################################################################
# Add high pressure feedwater heaters #
############################################################################
m.fs.fwh6 = FWH0D(
default={
"has_desuperheat": True,
"has_drain_cooling": True,
"has_drain_mixer": True,
"property_package": m.fs.prop_water,
"desuperheat": fwh_config,
"cooling": fwh_config,
"condense": fwh_config,
}
)
m.fs.fwh7 = FWH0D(
default={
"has_desuperheat": True,
"has_drain_cooling": True,
"has_drain_mixer": True,
"property_package": m.fs.prop_water,
"desuperheat": fwh_config,
"cooling": fwh_config,
"condense": fwh_config,
}
)
m.fs.fwh8 = FWH0D(
default={
"has_desuperheat": True,
"has_drain_cooling": True,
"has_drain_mixer": False,
"property_package": m.fs.prop_water,
"desuperheat": fwh_config,
"cooling": fwh_config,
"condense": fwh_config,
}
)
############################################################################
# Additional Constraints/Expressions #
############################################################################
# Add a few constraints to allow a for complete plant results despite the
# lack of a detailed boiler model.
# Boiler pressure drop
m.fs.boiler_pressure_drop_fraction = pyo.Var(
m.fs.time,
initialize=0.01,
doc="Fraction of pressure lost from boiler feed pump and turbine inlet",
)
@m.fs.Constraint(m.fs.time)
def boiler_pressure_drop(b, t):
return (
m.fs.bfp.control_volume.properties_out[t].pressure
* (1 - b.boiler_pressure_drop_fraction[t])
== m.fs.turb.inlet_split.mixed_state[t].pressure
)
# Again, since the boiler is missing, set the flow of steam into the turbine
# equal to the flow of feedwater out of the last feedwater heater.
@m.fs.Constraint(m.fs.time)
def close_flow(b, t):
return (
m.fs.bfp.control_volume.properties_out[t].flow_mol
== m.fs.turb.inlet_split.mixed_state[t].flow_mol
)
# Calculate the amount of heat that is added in the boiler, including the
# reheater.
@m.fs.Expression(m.fs.time)
def boiler_heat(b, t):
return (
b.turb.inlet_split.mixed_state[t].enth_mol
* b.turb.inlet_split.mixed_state[t].flow_mol
- b.fwh8.desuperheat.tube.properties_out[t].enth_mol
* b.fwh8.desuperheat.tube.properties_out[t].flow_mol
+ b.turb.ip_stages[1].control_volume.properties_in[t].enth_mol
* b.turb.ip_stages[1].control_volume.properties_in[t].flow_mol
- b.turb.hp_split[7].outlet_1.enth_mol[t]
* b.turb.hp_split[7].outlet_1.flow_mol[t]
)
# Calculate the efficiency of the steam cycle. This doesn't account for
# heat loss in the boiler, so actual plant efficiency would be lower.
@m.fs.Expression(m.fs.time)
def steam_cycle_eff(b, t):
return -100 * b.turb.power[t] / b.boiler_heat[t]
############################################################################
## Create the stream Arcs ##
############################################################################
############################################################################
# Connect turbine and condenser units #
############################################################################
m.fs.EXHST_MAIN = Arc(
source=m.fs.turb.outlet_stage.outlet, destination=m.fs.condenser_mix.main
)
m.fs.condenser_mix_to_condenser = Arc(
source=m.fs.condenser_mix.outlet_tpx, destination=m.fs.condenser.inlet_1
)
m.fs.COND_01 = Arc(
source=m.fs.condenser.outlet_1_ph, destination=m.fs.hotwell.condensate
)
m.fs.COND_02 = Arc(source=m.fs.hotwell.outlet, destination=m.fs.cond_pump.inlet)
############################################################################
# Low pressure FWHs #
############################################################################
m.fs.EXTR_LP11 = Arc(
source=m.fs.turb.lp_split[11].outlet_2, destination=m.fs.fwh1.drain_mix.steam
)
m.fs.COND_03 = Arc(
source=m.fs.cond_pump.outlet, destination=m.fs.fwh1.condense.inlet_2
)
m.fs.FWH1_DRN1 = Arc(
source=m.fs.fwh1.condense.outlet_1, destination=m.fs.fwh1_pump.inlet
)
m.fs.FWH1_DRN2 = Arc(
source=m.fs.fwh1_pump.outlet, destination=m.fs.fwh1_return.fwh1_drain
)
m.fs.FW01A = Arc(
source=m.fs.fwh1.condense.outlet_2, destination=m.fs.fwh1_return.feedwater
)
# fwh2
m.fs.FW01B = Arc(
source=m.fs.fwh1_return.outlet, destination=m.fs.fwh2.cooling.inlet_2
)
m.fs.FWH2_DRN = Arc(
source=m.fs.fwh2.cooling.outlet_1, destination=m.fs.fwh1.drain_mix.drain
)
m.fs.EXTR_LP10 = Arc(
source=m.fs.turb.lp_split[10].outlet_2,
destination=m.fs.fwh2.desuperheat.inlet_1,
)
# fwh3
m.fs.FW02 = Arc(
source=m.fs.fwh2.desuperheat.outlet_2, destination=m.fs.fwh3.cooling.inlet_2
)
m.fs.FWH3_DRN = Arc(
source=m.fs.fwh3.cooling.outlet_1, destination=m.fs.fwh2.drain_mix.drain
)
m.fs.EXTR_LP8 = Arc(
source=m.fs.turb.lp_split[8].outlet_2, destination=m.fs.fwh3.desuperheat.inlet_1
)
# fwh4
m.fs.FW03 = Arc(
source=m.fs.fwh3.desuperheat.outlet_2, destination=m.fs.fwh4.cooling.inlet_2
)
m.fs.FWH4_DRN = Arc(
source=m.fs.fwh4.cooling.outlet_1, destination=m.fs.fwh3.drain_mix.drain
)
m.fs.EXTR_LP4 = Arc(
source=m.fs.turb.lp_split[4].outlet_2, destination=m.fs.fwh4.desuperheat.inlet_1
)
############################################################################
# FWH5 (Deaerator) and boiler feed pump (BFP) #
############################################################################
m.fs.FW04 = Arc(
source=m.fs.fwh4.desuperheat.outlet_2, destination=m.fs.fwh5_da.feedwater
)
m.fs.EXTR_IP10 = Arc(
source=m.fs.turb.ip_split[10].outlet_2, destination=m.fs.fwh5_da.steam
)
m.fs.FW05A = Arc(source=m.fs.fwh5_da.outlet, destination=m.fs.bfp.inlet)
m.fs.EXTR_BFPT_A = Arc(
source=m.fs.turb.ip_split[10].outlet_3, destination=m.fs.bfpt.inlet
)
m.fs.EXHST_BFPT = Arc(source=m.fs.bfpt.outlet, destination=m.fs.condenser_mix.bfpt)
############################################################################
# High-pressure feedwater heaters #
############################################################################
# fwh6
m.fs.FW05B = Arc(source=m.fs.bfp.outlet, destination=m.fs.fwh6.cooling.inlet_2)
m.fs.FWH6_DRN = Arc(
source=m.fs.fwh6.cooling.outlet_1, destination=m.fs.fwh5_da.drain
)
m.fs.EXTR_IP5 = Arc(
source=m.fs.turb.ip_split[5].outlet_2, destination=m.fs.fwh6.desuperheat.inlet_1
)
# fwh7
m.fs.FW06 = Arc(
source=m.fs.fwh6.desuperheat.outlet_2, destination=m.fs.fwh7.cooling.inlet_2
)
m.fs.FWH7_DRN = Arc(
source=m.fs.fwh7.cooling.outlet_1, destination=m.fs.fwh6.drain_mix.drain
)
m.fs.EXTR_HP7 = Arc(
source=m.fs.turb.hp_split[7].outlet_2, destination=m.fs.fwh7.desuperheat.inlet_1
)
# fwh8
m.fs.FW07 = Arc(
source=m.fs.fwh7.desuperheat.outlet_2, destination=m.fs.fwh8.cooling.inlet_2
)
m.fs.FWH8_DRN = Arc(
source=m.fs.fwh8.cooling.outlet_1, destination=m.fs.fwh7.drain_mix.drain
)
m.fs.EXTR_HP4 = Arc(
source=m.fs.turb.hp_split[4].outlet_2, destination=m.fs.fwh8.desuperheat.inlet_1
)
############################################################################
# Turn the Arcs into constraints and return the model #
############################################################################
pyo.TransformationFactory("network.expand_arcs").apply_to(m.fs)
return m
def _stream_dict(m):
"""Adds _streams to m, which contains a dictionary of streams for display
Args:
m (ConcreteModel): A Pyomo model from create_model()
Returns:
None
"""
# Put all the top-level streams in the stream dictionary
m._streams = dict(
[(c.getname(), c) for c in m.fs.component_objects(Arc, descend_into=False)]
)
# There are some additional streams we are interested in that are either
# inlets or outlets, where there is no arc or arcs buried in unit models.
# Next those are added specifically.
m._streams.update(
{
"STEAM_MAIN": m.fs.turb.inlet_split.mixed_state,
"THRTL1": m.fs.turb.inlet_stage[1].control_volume.properties_in,
"THRTL2": m.fs.turb.inlet_stage[2].control_volume.properties_in,
"THRTL3": m.fs.turb.inlet_stage[3].control_volume.properties_in,
"THRTL4": m.fs.turb.inlet_stage[4].control_volume.properties_in,
"RHT_COLD": m.fs.turb.hp_split[7].outlet_1_state,
"RHT_HOT": m.fs.turb.ip_stages[1].control_volume.properties_in,
"STEAM_LP": m.fs.turb.ip_split[10].outlet_1_state,
"CW01": m.fs.condenser.tube.properties_in,
"CW02": m.fs.condenser.tube.properties_out,
"MAKEUP_01": m.fs.hotwell.makeup_state,
"FW08": m.fs.fwh8.desuperheat.tube.properties_out,
}
)
# Alphabetize to find streams easier in the tabular stream table view.
m._streams = OrderedDict(sorted(m._streams.items()))
# Now all the model input has been specified.
[docs]def initialize(m, fileinput=None, outlvl=3):
""" Initialize a mode from create_model(), set model inputs before
initializing.
Args:
m (ConcreteModel): A Pyomo model from create_model()
fileinput (str|None): File to load initialized model state from. If a
file is supplied skip initialization routine. If None, initialize.
Returns:
solver: A Pyomo solver object, that can be used to solve the model.
"""
init_log = getInitLogger(m.name, outlvl)
solver = pyo.SolverFactory("ipopt")
solver.options = {
"tol": 1e-7,
"linear_solver": "ma27",
"max_iter": 40,
}
if fileinput is not None:
init_log.log(5, "Loading initial values from file: {}".format(fileinput))
ms.from_json(m, fname=fileinput)
return solver
init_log.log(5, "Starting initialization")
############################################################################
# Initialize turbine #
############################################################################
# Extraction rates are calculated from the feedwater heater models, so to
# initialize the turbine fix some initial guesses. They get unfixed after
# solving the turbine
m.fs.turb.outlet_stage.control_volume.properties_out[:].pressure.fix(3500)
m.fs.turb.lp_split[11].split_fraction[:, "outlet_2"].fix(0.04403)
m.fs.turb.lp_split[10].split_fraction[:, "outlet_2"].fix(0.04025)
m.fs.turb.lp_split[8].split_fraction[:, "outlet_2"].fix(0.04362)
m.fs.turb.lp_split[4].split_fraction[:, "outlet_2"].fix(0.08025)
m.fs.turb.ip_split[10].split_fraction[:, "outlet_2"].fix(0.045)
m.fs.turb.ip_split[10].split_fraction[:, "outlet_3"].fix(0.04)
m.fs.turb.ip_split[5].split_fraction[:, "outlet_2"].fix(0.05557)
m.fs.turb.hp_split[7].split_fraction[:, "outlet_2"].fix(0.09741)
m.fs.turb.hp_split[4].split_fraction[:, "outlet_2"].fix(0.0740)
# Put in a rough initial guess for the IP section inlet, since it is
# disconnected from the HP section for the reheater.
ip1_pin = 5.35e6
ip1_hin = iapws95.htpx(T=866, P=ip1_pin)
ip1_fin = pyo.value(m.fs.turb.inlet_split.inlet.flow_mol[0])
m.fs.turb.ip_stages[1].inlet.enth_mol[:].value = ip1_hin
m.fs.turb.ip_stages[1].inlet.flow_mol[:].value = ip1_fin
m.fs.turb.ip_stages[1].inlet.pressure[:].value = ip1_pin
# initialize turbine
assert degrees_of_freedom(m.fs.turb) == 0
m.fs.turb.initialize(outlvl=5, optarg=solver.options)
# The turbine outlet pressure is determined by the condenser once hooked up
m.fs.turb.outlet_stage.control_volume.properties_out[:].pressure.unfix()
# Extraction rates are calculated once the feedwater heater models are
# added so unfix the splits.
m.fs.turb.lp_split[11].split_fraction[:, "outlet_2"].unfix()
m.fs.turb.lp_split[10].split_fraction[:, "outlet_2"].unfix()
m.fs.turb.lp_split[8].split_fraction[:, "outlet_2"].unfix()
m.fs.turb.lp_split[4].split_fraction[:, "outlet_2"].unfix()
m.fs.turb.ip_split[10].split_fraction[:, "outlet_3"].unfix()
m.fs.turb.ip_split[5].split_fraction[:, "outlet_2"].unfix()
m.fs.turb.hp_split[7].split_fraction[:, "outlet_2"].unfix()
m.fs.turb.hp_split[4].split_fraction[:, "outlet_2"].unfix()
# Initialize the boiler feed pump turbine.
_set_port(m.fs.bfpt.inlet, m.fs.turb.ip_split[10].outlet_3)
m.fs.bfpt.initialize(outlvl=5, optarg=solver.options)
############################################################################
# Condenser section #
############################################################################
# initialize condenser mixer
_set_port(m.fs.condenser_mix.main, m.fs.turb.outlet_stage.outlet)
_set_port(m.fs.condenser_mix.bfpt, m.fs.bfpt.outlet)
m.fs.condenser_mix.initialize(outlvl=5, optarg=solver.options)
# initialize condenser hx
_set_port(m.fs.condenser.inlet_1, m.fs.condenser_mix.outlet_tpx)
_set_port(m.fs.condenser.outlet_2, m.fs.condenser.inlet_2)
# This still has the outlet vapor fraction fixed at 0, so if that isn't
# true for the initial pressure guess this could fail to initialize
m.fs.condenser.inlet_1.fix()
m.fs.condenser.inlet_1.pressure.unfix()
res = solver.solve(m.fs.condenser, tee=init_tee(init_log))
init_log.log(4, "Condenser initialized: {}".format(condition(res)))
init_log.log(3, "Init condenser pressure: {}".format(
m.fs.condenser.shell.properties_in[0].pressure.value
)
)
m.fs.condenser.inlet_1.unfix()
# initialize hotwell
_set_port(m.fs.hotwell.condensate, m.fs.condenser.outlet_1_ph)
m.fs.hotwell.initialize(outlvl=5, optarg=solver.options)
m.fs.hotwell.condensate.unfix()
# initialize condensate pump
_set_port(m.fs.cond_pump.inlet, m.fs.hotwell.outlet)
m.fs.cond_pump.initialize(outlvl=5, optarg=solver.options)
############################################################################
# Low-pressure FWH section #
############################################################################
# fwh1
m.fs.fwh1.drain_mix.drain.flow_mol[:] = 1000
m.fs.fwh1.drain_mix.drain.pressure[:] = 1e5
m.fs.fwh1.drain_mix.drain.enth_mol[:] = 6117
_set_port(m.fs.fwh1.condense.inlet_2, m.fs.cond_pump.outlet)
_set_port(m.fs.fwh1.drain_mix.steam, m.fs.turb.lp_split[11].outlet_2)
m.fs.fwh1.initialize(outlvl=5, optarg=solver.options)
# initialize fwh1 drain pump
_set_port(m.fs.fwh1_pump.inlet, m.fs.fwh1.condense.outlet_1)
m.fs.fwh1_pump.initialize(outlvl=5, optarg=solver.options)
# initialize mixer to add fwh1 drain to feedwater
_set_port(m.fs.fwh1_return.feedwater, m.fs.fwh1.condense.outlet_2)
_set_port(m.fs.fwh1_return.fwh1_drain, m.fs.fwh1.condense.outlet_1)
m.fs.fwh1_return.initialize(outlvl=5, optarg=solver.options)
m.fs.fwh1_return.feedwater.unfix()
m.fs.fwh1_return.fwh1_drain.unfix()
# fwh2
m.fs.fwh2.drain_mix.drain.flow_mol[:] = 100
m.fs.fwh2.drain_mix.drain.pressure[:] = 1.5e5
m.fs.fwh2.drain_mix.drain.enth_mol[:] = 7000
_set_port(m.fs.fwh2.cooling.inlet_2, m.fs.fwh1_return.outlet)
_set_port(m.fs.fwh2.desuperheat.inlet_1, m.fs.turb.lp_split[10].outlet_2)
m.fs.fwh2.initialize(outlvl=5, optarg=solver.options)
# fwh3
m.fs.fwh3.drain_mix.drain.flow_mol[:] = 100
m.fs.fwh3.drain_mix.drain.pressure[:] = 2.5e5
m.fs.fwh3.drain_mix.drain.enth_mol[:] = 8000
_set_port(m.fs.fwh3.cooling.inlet_2, m.fs.fwh2.desuperheat.outlet_2)
_set_port(m.fs.fwh3.desuperheat.inlet_1, m.fs.turb.lp_split[8].outlet_2)
m.fs.fwh3.initialize(outlvl=5, optarg=solver.options)
# fwh4
_set_port(m.fs.fwh4.cooling.inlet_2, m.fs.fwh3.desuperheat.outlet_2)
_set_port(m.fs.fwh4.desuperheat.inlet_1, m.fs.turb.lp_split[4].outlet_2)
m.fs.fwh4.initialize(outlvl=5, optarg=solver.options)
############################################################################
# boiler feed pump and deaerator #
############################################################################
_set_port(m.fs.fwh5_da.feedwater, m.fs.fwh4.desuperheat.outlet_2)
_set_port(m.fs.fwh5_da.steam, m.fs.turb.ip_split[10].outlet_2)
m.fs.fwh5_da.drain.flow_mol[:] = 2000
m.fs.fwh5_da.drain.pressure[:] = 3e6
m.fs.fwh5_da.drain.enth_mol[:] = 9000
m.fs.fwh5_da.initialize(outlvl=5, optarg=solver.options)
_set_port(m.fs.bfp.inlet, m.fs.fwh5_da.outlet)
m.fs.bfp.control_volume.properties_out[:].pressure.fix()
m.fs.bfp.initialize(outlvl=5, optarg=solver.options)
m.fs.bfp.control_volume.properties_out[:].pressure.unfix()
############################################################################
# High-pressure feedwater heaters #
############################################################################
# fwh6
m.fs.fwh6.drain_mix.drain.flow_mol[:] = 1000
m.fs.fwh6.drain_mix.drain.pressure[:] = 1e7
m.fs.fwh6.drain_mix.drain.enth_mol[:] = 9500
_set_port(m.fs.fwh6.cooling.inlet_2, m.fs.bfp.outlet)
_set_port(m.fs.fwh6.desuperheat.inlet_1, m.fs.turb.ip_split[5].outlet_2)
m.fs.fwh6.initialize(outlvl=5, optarg=solver.options)
# fwh7
m.fs.fwh7.drain_mix.drain.flow_mol[:] = 2000
m.fs.fwh7.drain_mix.drain.pressure[:] = 1e7
m.fs.fwh7.drain_mix.drain.enth_mol[:] = 9500
_set_port(m.fs.fwh7.cooling.inlet_2, m.fs.fwh6.desuperheat.outlet_2)
_set_port(m.fs.fwh7.desuperheat.inlet_1, m.fs.turb.hp_split[7].outlet_2)
m.fs.fwh7.initialize(outlvl=5, optarg=solver.options)
# fwh8
_set_port(m.fs.fwh8.cooling.inlet_2, m.fs.fwh7.desuperheat.outlet_2)
_set_port(m.fs.fwh8.desuperheat.inlet_1, m.fs.turb.hp_split[4].outlet_2)
m.fs.fwh8.initialize(outlvl=5, optarg=solver.options)
res = solver.solve(m, tee=init_tee(init_log))
init_log.log(5, "Initialization Complete: {}".format(condition(res)))
return solver
[docs]def pfd_result(m, df, svg):
"""Insert model results into the PFD and return a new SVG string, which
can be displayed, further edited, or saved to a file.
Args:
m (ConcreteModel): A steam cycle model
df (Pandas DataFrame): Stream table
svg (FILE*, str, bytes): Origianl svg svg as either a file-like object,
a string, or a byte array.
Returns:
(str): SVG content.
"""
tags = {} # dict of tags and data to insert into SVG
for i in df.index: # Create entires for streams
tags[i + "_F"] = df.loc[i, "Molar Flow (mol/s)"]
tags[i + "_T"] = df.loc[i, "T (K)"]
tags[i + "_P"] = df.loc[i, "P (Pa)"]
tags[i + "_X"] = df.loc[i, "Vapor Fraction"]
# Add some addtional quntities from the model to report
tags["gross_power"] = -pyo.value(m.fs.turb.power[0])
tags["gross_power_mw"] = -pyo.value(m.fs.turb.power[0]) * 1e-6
tags["steam_mass_flow"] = df.loc["STEAM_MAIN", "Mass Flow (kg/s)"]
tags["sc_eff"] = pyo.value(m.fs.steam_cycle_eff[0])
tags["boiler_heat"] = pyo.value(m.fs.boiler_heat[0]) * 1e-6
tags["steam_pressure"] = df.loc["STEAM_MAIN", "P (Pa)"] / 1000.0
tags["cond_pressure"] = df.loc["EXHST_MAIN", "P (Pa)"] / 1000.0
tags["bfp_power"] = pyo.value(m.fs.bfp.work_mechanical[0])
tags["bfp_eff"] = pyo.value(m.fs.bfp.efficiency_pump[0]) * 100
tags["bfpt_power"] = pyo.value(m.fs.bfpt.work_mechanical[0])
tags["bfpt_eff"] = pyo.value(m.fs.bfpt.efficiency_isentropic[0]) * 100
return svg_tag(tags, svg=svg)
[docs]def main(initialize_from_file=None, store_initialization=None):
""" Create and initalize a model and solver
Args:
None
Returns:
A tuple of a model and solver
"""
m = create_model()
_stream_dict(m)
set_model_input(m)
if initialize_from_file is None:
solver = initialize(m)
else:
solver = initialize(m, fileinput=initialize_from_file)
solver.solve(m, tee=True)
if store_initialization is not None:
ms.to_json(m, fname=store_initialization)
return m, solver
if __name__ == "__main__":
parser = argparse.ArgumentParser()
parser.add_argument(
"--initialize_from_file",
help="File from which to load initialized values. If specified, the "
"initialization proceedure will be skipped.",
default=None,
)
parser.add_argument(
"--store_initialization",
help="If specified, store" " the initialized model values, to reload later.",
default=None,
)
args = parser.parse_args()
m, solver = main(
initialize_from_file=args.initialize_from_file,
store_initialization=args.store_initialization,
)