oemof / tespy

Thermal Engineering Systems in Python (TESPy). This package provides a powerful simulation toolkit for thermal engineering plants such as power plants, district heating systems or heat pumps.
https://tespy.readthedocs.io
MIT License
272 stars 85 forks source link

Request for help with offdesign calculations #256

Closed govind-menon110 closed 3 years ago

govind-menon110 commented 3 years ago

I have an issue with turning a design run into an offdesign one. The jacobian matrix is always ending up forming a singularity. I feel I am doing something wrong with what parameters I omit during offdesign. Any help is highly appreciated.

What I am trying to do is the follows:

  1. If I change the mass flow rate, I want the turbine efficiencies to change and correspondingly its output powers.
  2. If it is possible, I also want to be able to set a pressure change at crucial junctures and see its effect on the power output

Attached herewith is the design code. Do let me know how do I turn it into an offdesign case.

# %% Starting the Network
fluid_list = ['Ar', 'N2', 'H2O', 'CO2', 'O2']
nw = Network(fluids=fluid_list, T_unit='C', p_unit='bar', h_unit='kJ / kg')

# Parameters for fluid flow
flueflow = 604.8

hptsink = 45.3
waterflow = 425.8

iptsploffset = 0 # +-30% 190.25 | max value is 0.3*(waterflow-hptsink)/2 | min is -0.3*(waterflow-hptsink)/2

ipt1psink = 10.5 
ipt2psink = 10.5 

iptlptsink = 42.6

lptsploffset = 0 # Can be maximum of +-0.3*(waterflow-hptsink-ipt1psink-ipt2psink-iptlptsink)/2

lpt1spl = 6.1
lpt2spl = 6.1

lpt1pspl = 11.65
lpt2pspl = 11.65

lpt1dpspl = 5.1
lpt2dpspl = 5.1

# %% Instantiating Components

#src sink
wtr_src = Source('Water Source')
flue = Source('Flue Gas Source')
stack = Sink('Flue Gas Stack')

hpt_sink = Sink('HPT Sink')
ipt_sink_1 = Sink('IPT 1 & 2 Sink')
ipt_sink_2 = Sink('IPT 1p & 2p Sink')
lpt_sink_1 = Sink('LPT 1 & 2 Sink')
lpt_sink_2 = Sink('LPT 1p & 2p Sink')
lpt_sink_3 = Sink('LPT 1dp & 2dp Sink')
lpt_sink_4 = Sink('LPT 1tp & 2tp Sink')

#Heat Exhangers

evap_sup = SteamGenerator('Evaporator & Superheater')
reh = SteamGenerator('Reheater')

evap_sup.set_attr(pr1=0.98, pr2=0.888, design=['pr2'], offdesign=['zeta2','kA_char2'])
reh.set_attr(pr1=0.99, pr2=0.91, design=['pr2'], offdesign=['zeta2', 'kA_char2'])

#Splitter & Merger
hpt_spl = Splitter('HPT Splitter')

ipt_spl = Splitter('IPT Splitter')
ipt1_spl = Splitter('IPT1 Splitter')
ipt2_spl = Splitter('IPT2 Splitter')
ipt12_mrg = Merge('IPT1 & 2 Merger')
ipt12p_mrg = Merge('IPT1p & 2p Merger')

ipt_lpt_spl = Splitter('IPT to LPT Splitter')

lpt_spl = Splitter('LPT Splitter')
lpt1_spl = Splitter('LPT1 Splitter')
lpt1p_spl = Splitter('LPT1p Splitter')
lpt1dp_spl = Splitter('LPT1dp Splitter')
lpt2_spl = Splitter('LPT2 Splitter')
lpt2p_spl = Splitter('LPT2p Splitter')
lpt2dp_spl = Splitter('LPT2dp Splitter')

lpt12_mrg = Merge('LPT1 & 2 Merger')
lpt12p_mrg = Merge('LPT1p & 2p Merger')
lpt12dp_mrg = Merge('LPT1dp & 2dp Merger')
lpt12tp_mrg = Merge('LPT1tp & 2tp Merger')

#Turbines
hpt = Turbine('HPT', eta_s=0.89, design=['eta_s'], offdesign=['eta_s_char', 'cone'])

ipt1 = Turbine('IPT 1', eta_s=0.903, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
ipt2 = Turbine('IPT 2', eta_s=0.903, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
ipt1p = Turbine('IPT 1p', eta_s=0.903, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
ipt2p = Turbine('IPT 2p', eta_s=0.903, design=['eta_s'], offdesign=['eta_s_char', 'cone'])

lpt1 = Turbine('LPT 1', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt2 = Turbine('LPT 2', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt1p = Turbine('LPT 1p', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt2p = Turbine('LPT 2p', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt1dp = Turbine('LPT 1dp', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt2dp = Turbine('LPT 2dp', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt1tp = Turbine('LPT 1tp', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt2tp = Turbine('LPT 2tp', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])

# %% Setting the Connections
# Air/Flue

flue_evap = Connection(flue, 'out1', evap_sup, 'in1')
evap_reh = Connection(evap_sup, 'out1', reh, 'in1')
reh_stack = Connection(reh, 'out1', stack, 'in1')

# Water
src_evap = Connection(wtr_src, 'out1', evap_sup, 'in2')
evap_hpt = Connection(evap_sup, 'out2', hpt, 'in1')

hpt_hptspl = Connection(hpt, 'out1', hpt_spl, 'in1')
hptspl_sink = Connection(hpt_spl, 'out1', hpt_sink, 'in1')
hptspl_reh = Connection(hpt_spl, 'out2', reh, 'in2')

reh_iptspl = Connection(reh, 'out2', ipt_spl, 'in1')
iptspl_1 = Connection(ipt_spl, 'out1', ipt1, 'in1')
iptspl_2 = Connection(ipt_spl, 'out2', ipt2, 'in1')

ipt1_1spl = Connection(ipt1, 'out1', ipt1_spl, 'in1')
ipt1spl_1p = Connection(ipt1_spl, 'out1', ipt1p, 'in1')
ipt1spl_mrg = Connection(ipt1_spl, 'out2', ipt12_mrg, 'in1')

ipt2_2spl = Connection(ipt2, 'out1', ipt2_spl, 'in1')
ipt2spl_2p = Connection(ipt2_spl, 'out1', ipt2p, 'in1')
ipt2spl_mrg = Connection(ipt2_spl, 'out2', ipt12_mrg, 'in2')

ipt12mrg_s = Connection(ipt12_mrg, 'out1', ipt_sink_1, 'in1')

ipt1p_mrg = Connection(ipt1p, 'out1', ipt12p_mrg, 'in1')
ipt2p_mrg = Connection(ipt2p, 'out1', ipt12p_mrg, 'in2')
ipt12pmrg_spl = Connection(ipt12p_mrg, 'out1', ipt_lpt_spl, 'in1')

iptlptspl_spl = Connection(ipt_lpt_spl, 'out1', lpt_spl, 'in1')
iptlptspl_s = Connection(ipt_lpt_spl, 'out2', ipt_sink_2, 'in1')

lptspl_1 = Connection(lpt_spl, 'out1', lpt1, 'in1')
lptspl_2 = Connection(lpt_spl, 'out2', lpt2, 'in1')

lpt1_1spl = Connection(lpt1, 'out1', lpt1_spl, 'in1')
lpt1spl_1p = Connection(lpt1_spl, 'out1', lpt1p, 'in1')
lpt1spl_mrg = Connection(lpt1_spl, 'out2', lpt12_mrg, 'in1')

lpt2_2spl = Connection(lpt2, 'out1', lpt2_spl, 'in1')
lpt2spl_2p = Connection(lpt2_spl, 'out1', lpt2p, 'in1')
lpt2spl_mrg = Connection(lpt2_spl, 'out2', lpt12_mrg, 'in2')

lpt12mrg_sink = Connection(lpt12_mrg, 'out1', lpt_sink_1, 'in1')

lpt1p_1spl = Connection(lpt1p, 'out1', lpt1p_spl, 'in1')
lpt1pspl_1dp = Connection(lpt1p_spl, 'out1', lpt1dp, 'in1')
lpt1pspl_mrg = Connection(lpt1p_spl, 'out2', lpt12p_mrg, 'in1')

lpt2p_2spl = Connection(lpt2p, 'out1', lpt2p_spl, 'in1')
lpt2pspl_2dp = Connection(lpt2p_spl, 'out1', lpt2dp, 'in1')
lpt2pspl_mrg = Connection(lpt2p_spl, 'out2', lpt12p_mrg, 'in2')

lpt12pmrg_sink = Connection(lpt12p_mrg, 'out1', lpt_sink_2, 'in1')

lpt1dp_1spl = Connection(lpt1dp, 'out1', lpt1dp_spl, 'in1')
lpt1dpspl_1tp = Connection(lpt1dp_spl, 'out1', lpt1tp, 'in1')
lpt1dpspl_mrg = Connection(lpt1dp_spl, 'out2', lpt12dp_mrg, 'in1')

lpt2dp_2spl = Connection(lpt2dp, 'out1', lpt2dp_spl, 'in1')
lpt2dpspl_2tp = Connection(lpt2dp_spl, 'out1', lpt2tp, 'in1')
lpt2dpspl_mrg = Connection(lpt2dp_spl, 'out2', lpt12dp_mrg, 'in2')

lpt12dpmrg_sink = Connection(lpt12dp_mrg, 'out1', lpt_sink_3, 'in1')

lpt1tp_mrg = Connection(lpt1tp, 'out1', lpt12tp_mrg, 'in1')
lpt2tp_mrg = Connection(lpt2tp, 'out1', lpt12tp_mrg, 'in2')
lpt12tpmrg_sink = Connection(lpt12tp_mrg, 'out1', lpt_sink_4, 'in1')

# %% Setting Connection Parameters
# Combustion
flue_evap.set_attr(p=1.01, T=1782.7, m=flueflow, fluid={'CO2': 0.25, 'N2': 0.709, 'Ar': 0.0102, 
                                            'O2': 0.0348, 'H2O': 0.028})

evap_reh.set_attr(T=802.3)
reh_stack.set_attr(T=529.8)

# Water
src_evap.set_attr(p=193.7, T=324, m=waterflow, fluid={'Ar': 0, 'N2': 0, 'H2O': 1,
                                                    'CO2': 0, 'O2': 0})

hpt_hptspl.set_attr(p=44.1)

hptspl_sink.set_attr(m=hptsink)

iptspl_1.set_attr(m=(waterflow-hptsink)/2+iptsploffset)

ipt1_1spl.set_attr(p=17.4)

ipt1spl_mrg.set_attr(m=ipt1psink)
ipt2spl_mrg.set_attr(m=ipt2psink)

ipt1p_mrg.set_attr(p=7.93)

iptlptspl_s.set_attr(m=iptlptsink)

lptspl_1.set_attr(m=(waterflow-hptsink-ipt1psink-ipt2psink-iptlptsink)/2+lptsploffset)

lpt1_1spl.set_attr(p=2.6)

lpt1spl_mrg.set_attr(m=lpt1spl)
lpt2spl_mrg.set_attr(m=lpt2spl)

lpt1p_1spl.set_attr(p=1.3)

lpt1pspl_mrg.set_attr(m=lpt1pspl)
lpt2pspl_mrg.set_attr(m=lpt2pspl)

lpt1dp_1spl.set_attr(p=0.26)

lpt1dpspl_mrg.set_attr(m=lpt1dpspl)
lpt2dpspl_mrg.set_attr(m=lpt2dpspl)

lpt1tp_mrg.set_attr(p=0.103)

nw.add_conns(flue_evap, evap_reh, reh_stack)
nw.add_conns(src_evap, evap_hpt, hpt_hptspl, hptspl_sink, hptspl_reh, reh_iptspl,
             iptspl_1, iptspl_2, ipt1_1spl, ipt1spl_1p, ipt1spl_mrg, ipt2_2spl, ipt2spl_2p, ipt2spl_mrg, ipt12mrg_s,
             ipt1p_mrg, ipt2p_mrg, ipt12pmrg_spl, iptlptspl_s, iptlptspl_spl,
             lptspl_1, lptspl_2, lpt1_1spl, lpt1spl_1p, lpt1spl_mrg, lpt2_2spl, lpt2spl_2p, lpt2spl_mrg, lpt12mrg_sink,
             lpt1p_1spl, lpt1pspl_1dp, lpt1pspl_mrg, lpt2p_2spl, lpt2pspl_2dp, lpt2pspl_mrg, lpt12pmrg_sink,
             lpt1dp_1spl, lpt1dpspl_1tp, lpt1dpspl_mrg, lpt2dp_2spl, lpt2dpspl_2tp, lpt2dpspl_mrg, lpt12dpmrg_sink,
             lpt1tp_mrg, lpt2tp_mrg, lpt12tpmrg_sink)

# %% Solving Plant for design
nw.solve('design')
nw.print_results()
nw.save(path='./my_char')

If it helps I have also attached a diagram of the plant image

Note: The combustor was removed because of inadequate data I had

fwitte commented 3 years ago

Hi, the imports are missing in your script, can you provide them? Is the SteamGenerator a custom built component? Without those I cannot run the design case either.

Generally, I recommend starting with a minimal example and build from there step by step testing offdesign calculation as well asking what parameters should be fixed for offdesign and what parameters should be free instead (e.g. calculate efficiency from characteristic function instead of specified value or calculate turbine inlet pressure from cone law instead of directly specified value). Also, a way to test, if your specifications for design/offdesign switching are correct is to run the offdesign case with no changes to the design case. The solver should easily find the solution and the result should be identical to the design case.

Best regards

govind-menon110 commented 3 years ago

Apologies for not mentioning. Actually SteamGenerator is HeatExchanger itself. Here is the code according to tespy components. Running the offdesign case without any changes also doesn't happen with the jacobian matrix failing.

# %% Starting the Network
fluid_list = ['Ar', 'N2', 'H2O', 'CO2', 'O2']
nw = Network(fluids=fluid_list, T_unit='C', p_unit='bar', h_unit='kJ / kg')

flueflow = 604.8

hptsink = 45.3
waterflow = 425.8

iptsploffset = 0 # +-30% 190.25 | max value is 0.3*(waterflow-hptsink)/2 | min is -0.3*(waterflow-hptsink)/2

ipt1psink = 10.5 
ipt2psink = 10.5 

iptlptsink = 42.6

lptsploffset = 0 # Can be maximum of +-0.3*(waterflow-hptsink-ipt1psink-ipt2psink-iptlptsink)/2

lpt1spl = 6.1
lpt2spl = 6.1

lpt1pspl = 11.65
lpt2pspl = 11.65

lpt1dpspl = 5.1
lpt2dpspl = 5.1

# %% Instantiating Components

#src sink
wtr_src = Source('Water Source')
flue = Source('Flue Gas Source')
stack = Sink('Flue Gas Stack')

hpt_sink = Sink('HPT Sink')
ipt_sink_1 = Sink('IPT 1 & 2 Sink')
ipt_sink_2 = Sink('IPT 1p & 2p Sink')
lpt_sink_1 = Sink('LPT 1 & 2 Sink')
lpt_sink_2 = Sink('LPT 1p & 2p Sink')
lpt_sink_3 = Sink('LPT 1dp & 2dp Sink')
lpt_sink_4 = Sink('LPT 1tp & 2tp Sink')

#Heat Exhangers

evap_sup = HeatExchanger('Evaporator & Superheater')
reh = HeatExchanger('Reheater')

evap_sup.set_attr(pr1=0.98, pr2=0.888, design=['pr2'], offdesign=['zeta2','kA_char2'])
reh.set_attr(pr1=0.99, pr2=0.91, design=['pr2'], offdesign=['zeta2', 'kA_char2'])

#Splitter & Merger
hpt_spl = Splitter('HPT Splitter')

ipt_spl = Splitter('IPT Splitter')
ipt1_spl = Splitter('IPT1 Splitter')
ipt2_spl = Splitter('IPT2 Splitter')
ipt12_mrg = Merge('IPT1 & 2 Merger')
ipt12p_mrg = Merge('IPT1p & 2p Merger')

ipt_lpt_spl = Splitter('IPT to LPT Splitter')

lpt_spl = Splitter('LPT Splitter')
lpt1_spl = Splitter('LPT1 Splitter')
lpt1p_spl = Splitter('LPT1p Splitter')
lpt1dp_spl = Splitter('LPT1dp Splitter')
lpt2_spl = Splitter('LPT2 Splitter')
lpt2p_spl = Splitter('LPT2p Splitter')
lpt2dp_spl = Splitter('LPT2dp Splitter')

lpt12_mrg = Merge('LPT1 & 2 Merger')
lpt12p_mrg = Merge('LPT1p & 2p Merger')
lpt12dp_mrg = Merge('LPT1dp & 2dp Merger')
lpt12tp_mrg = Merge('LPT1tp & 2tp Merger')

#Turbines
hpt = Turbine('HPT', eta_s=0.89, design=['eta_s'], offdesign=['eta_s_char', 'cone'])

ipt1 = Turbine('IPT 1', eta_s=0.903, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
ipt2 = Turbine('IPT 2', eta_s=0.903, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
ipt1p = Turbine('IPT 1p', eta_s=0.903, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
ipt2p = Turbine('IPT 2p', eta_s=0.903, design=['eta_s'], offdesign=['eta_s_char', 'cone'])

lpt1 = Turbine('LPT 1', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt2 = Turbine('LPT 2', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt1p = Turbine('LPT 1p', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt2p = Turbine('LPT 2p', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt1dp = Turbine('LPT 1dp', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt2dp = Turbine('LPT 2dp', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt1tp = Turbine('LPT 1tp', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])
lpt2tp = Turbine('LPT 2tp', eta_s=0.851, design=['eta_s'], offdesign=['eta_s_char', 'cone'])

# %% Setting the Connections
# Air/Flue

flue_evap = Connection(flue, 'out1', evap_sup, 'in1')
evap_reh = Connection(evap_sup, 'out1', reh, 'in1')
reh_stack = Connection(reh, 'out1', stack, 'in1')

# Water
src_evap = Connection(wtr_src, 'out1', evap_sup, 'in2')
evap_hpt = Connection(evap_sup, 'out2', hpt, 'in1')

hpt_hptspl = Connection(hpt, 'out1', hpt_spl, 'in1')
hptspl_sink = Connection(hpt_spl, 'out1', hpt_sink, 'in1')
hptspl_reh = Connection(hpt_spl, 'out2', reh, 'in2')

reh_iptspl = Connection(reh, 'out2', ipt_spl, 'in1')
iptspl_1 = Connection(ipt_spl, 'out1', ipt1, 'in1')
iptspl_2 = Connection(ipt_spl, 'out2', ipt2, 'in1')

ipt1_1spl = Connection(ipt1, 'out1', ipt1_spl, 'in1')
ipt1spl_1p = Connection(ipt1_spl, 'out1', ipt1p, 'in1')
ipt1spl_mrg = Connection(ipt1_spl, 'out2', ipt12_mrg, 'in1')

ipt2_2spl = Connection(ipt2, 'out1', ipt2_spl, 'in1')
ipt2spl_2p = Connection(ipt2_spl, 'out1', ipt2p, 'in1')
ipt2spl_mrg = Connection(ipt2_spl, 'out2', ipt12_mrg, 'in2')

ipt12mrg_s = Connection(ipt12_mrg, 'out1', ipt_sink_1, 'in1')

ipt1p_mrg = Connection(ipt1p, 'out1', ipt12p_mrg, 'in1')
ipt2p_mrg = Connection(ipt2p, 'out1', ipt12p_mrg, 'in2')
ipt12pmrg_spl = Connection(ipt12p_mrg, 'out1', ipt_lpt_spl, 'in1')

iptlptspl_spl = Connection(ipt_lpt_spl, 'out1', lpt_spl, 'in1')
iptlptspl_s = Connection(ipt_lpt_spl, 'out2', ipt_sink_2, 'in1')

lptspl_1 = Connection(lpt_spl, 'out1', lpt1, 'in1')
lptspl_2 = Connection(lpt_spl, 'out2', lpt2, 'in1')

lpt1_1spl = Connection(lpt1, 'out1', lpt1_spl, 'in1')
lpt1spl_1p = Connection(lpt1_spl, 'out1', lpt1p, 'in1')
lpt1spl_mrg = Connection(lpt1_spl, 'out2', lpt12_mrg, 'in1')

lpt2_2spl = Connection(lpt2, 'out1', lpt2_spl, 'in1')
lpt2spl_2p = Connection(lpt2_spl, 'out1', lpt2p, 'in1')
lpt2spl_mrg = Connection(lpt2_spl, 'out2', lpt12_mrg, 'in2')

lpt12mrg_sink = Connection(lpt12_mrg, 'out1', lpt_sink_1, 'in1')

lpt1p_1spl = Connection(lpt1p, 'out1', lpt1p_spl, 'in1')
lpt1pspl_1dp = Connection(lpt1p_spl, 'out1', lpt1dp, 'in1')
lpt1pspl_mrg = Connection(lpt1p_spl, 'out2', lpt12p_mrg, 'in1')

lpt2p_2spl = Connection(lpt2p, 'out1', lpt2p_spl, 'in1')
lpt2pspl_2dp = Connection(lpt2p_spl, 'out1', lpt2dp, 'in1')
lpt2pspl_mrg = Connection(lpt2p_spl, 'out2', lpt12p_mrg, 'in2')

lpt12pmrg_sink = Connection(lpt12p_mrg, 'out1', lpt_sink_2, 'in1')

lpt1dp_1spl = Connection(lpt1dp, 'out1', lpt1dp_spl, 'in1')
lpt1dpspl_1tp = Connection(lpt1dp_spl, 'out1', lpt1tp, 'in1')
lpt1dpspl_mrg = Connection(lpt1dp_spl, 'out2', lpt12dp_mrg, 'in1')

lpt2dp_2spl = Connection(lpt2dp, 'out1', lpt2dp_spl, 'in1')
lpt2dpspl_2tp = Connection(lpt2dp_spl, 'out1', lpt2tp, 'in1')
lpt2dpspl_mrg = Connection(lpt2dp_spl, 'out2', lpt12dp_mrg, 'in2')

lpt12dpmrg_sink = Connection(lpt12dp_mrg, 'out1', lpt_sink_3, 'in1')

lpt1tp_mrg = Connection(lpt1tp, 'out1', lpt12tp_mrg, 'in1')
lpt2tp_mrg = Connection(lpt2tp, 'out1', lpt12tp_mrg, 'in2')
lpt12tpmrg_sink = Connection(lpt12tp_mrg, 'out1', lpt_sink_4, 'in1')

# %% Setting Connection Parameters
# Combustion
flue_evap.set_attr(p=1.01, T=1782.7, m=flueflow, fluid={'CO2': 0.25, 'N2': 0.709, 'Ar': 0.0102, 
                                            'O2': 0.0348, 'H2O': 0.028})

evap_reh.set_attr(T=802.3)
reh_stack.set_attr(T=529.8)

# Water
src_evap.set_attr(p=193.7, T=324, m=waterflow, fluid={'Ar': 0, 'N2': 0, 'H2O': 1,
                                                    'CO2': 0, 'O2': 0})

hpt_hptspl.set_attr(p=44.1)

hptspl_sink.set_attr(m=hptsink)

iptspl_1.set_attr(m=(waterflow-hptsink)/2+iptsploffset)

ipt1_1spl.set_attr(p=17.4)

ipt1spl_mrg.set_attr(m=ipt1psink)
ipt2spl_mrg.set_attr(m=ipt2psink)

ipt1p_mrg.set_attr(p=7.93)

iptlptspl_s.set_attr(m=iptlptsink)

lptspl_1.set_attr(m=(waterflow-hptsink-ipt1psink-ipt2psink-iptlptsink)/2+lptsploffset)

lpt1_1spl.set_attr(p=2.6)

lpt1spl_mrg.set_attr(m=lpt1spl)
lpt2spl_mrg.set_attr(m=lpt2spl)

lpt1p_1spl.set_attr(p=1.3)

lpt1pspl_mrg.set_attr(m=lpt1pspl)
lpt2pspl_mrg.set_attr(m=lpt2pspl)

lpt1dp_1spl.set_attr(p=0.26)

lpt1dpspl_mrg.set_attr(m=lpt1dpspl)
lpt2dpspl_mrg.set_attr(m=lpt2dpspl)

lpt1tp_mrg.set_attr(p=0.103)

nw.add_conns(flue_evap, evap_reh, reh_stack)
nw.add_conns(src_evap, evap_hpt, hpt_hptspl, hptspl_sink, hptspl_reh, reh_iptspl,
             iptspl_1, iptspl_2, ipt1_1spl, ipt1spl_1p, ipt1spl_mrg, ipt2_2spl, ipt2spl_2p, ipt2spl_mrg, ipt12mrg_s,
             ipt1p_mrg, ipt2p_mrg, ipt12pmrg_spl, iptlptspl_s, iptlptspl_spl,
             lptspl_1, lptspl_2, lpt1_1spl, lpt1spl_1p, lpt1spl_mrg, lpt2_2spl, lpt2spl_2p, lpt2spl_mrg, lpt12mrg_sink,
             lpt1p_1spl, lpt1pspl_1dp, lpt1pspl_mrg, lpt2p_2spl, lpt2pspl_2dp, lpt2pspl_mrg, lpt12pmrg_sink,
             lpt1dp_1spl, lpt1dpspl_1tp, lpt1dpspl_mrg, lpt2dp_2spl, lpt2dpspl_2tp, lpt2dpspl_mrg, lpt12dpmrg_sink,
             lpt1tp_mrg, lpt2tp_mrg, lpt12tpmrg_sink)

# %% Solving Plant for design
nw.solve('design')
nw.print_results()
nw.save(path='./my_char')

Thank you so much for your patience @fwitte

fwitte commented 3 years ago

When I run the offdesign case right after, there are missing a couple of parameters. How I said before, I strongly recommend building the offdesign case step by step. For example consider this simple setup:

# %% Importing modules
from tespy.networks import Network
from tespy.components import Sink, Source
from tespy.components import HeatExchangerSimple, Turbine
from tespy.connections import Connection

# %% Starting the Network
fluid_list = ['water']
nw = Network(fluids=fluid_list, T_unit='C', p_unit='bar', h_unit='kJ / kg', p_range=[0.01, 250])

# %% Instantiating Components

#src sink
so = Source('Source')
si = Sink('Sink')

#Heat Exhangers

b = HeatExchangerSimple('Boiler')
t = Turbine('Turbine')

b.set_attr(pr=0.98, design=['pr'], offdesign=['zeta'])
t.set_attr(eta_s=0.9, design=['eta_s'], offdesign=['cone', 'eta_s_char'])

c1 = Connection(so, 'out1', b, 'in1')
c2 = Connection(b, 'out1', t, 'in1')
c3 = Connection(t, 'out1', si, 'in1')

c1.set_attr(T=150, m=425.8, fluid={'water': 1}, design=['T'])
c2.set_attr(T=650, p=190)
c3.set_attr(p=10)

nw.add_conns(c1, c2, c3)

# %% Solving Plant
nw.solve('design')
nw.save('tmp')
nw.solve('offdesign', design_path='tmp')

The offdesign case will not be solved as the isentropic efficiency characteristic of the turbine and the cone low impose two equations on the change of state over the turbine while the inlet state is fully fixed by p, T specification and the outlet state is missing one parameter (p specified only, thus requiring one equation). The specification of the boiler inlet temperature as offdesign parameter simply does not make sense in this example.

Therefore, the outlet pressure or the turbine inlet temperature would be kind of reasonable parameters to unset for offdesign calculation (depending on how the plant is actually controlled).

# %% Importing modules
from tespy.networks import Network
from tespy.components import Sink, Source
from tespy.components import HeatExchangerSimple, Turbine
from tespy.connections import Connection

# %% Starting the Network
fluid_list = ['water']
nw = Network(fluids=fluid_list, T_unit='C', p_unit='bar', h_unit='kJ / kg', p_range=[0.01, 250])

# %% Instantiating Components

#src sink
so = Source('Source')
si = Sink('Sink')

#Heat Exhangers

b = HeatExchangerSimple('Boiler')
t = Turbine('Turbine')

b.set_attr(pr=0.98, design=['pr'], offdesign=['zeta'])
t.set_attr(eta_s=0.9, design=['eta_s'], offdesign=['cone', 'eta_s_char'])

c1 = Connection(so, 'out1', b, 'in1')
c2 = Connection(b, 'out1', t, 'in1')
c3 = Connection(t, 'out1', si, 'in1')

c1.set_attr(T=150, m=425.8, fluid={'water': 1})
c2.set_attr(T=650, p=190)
c3.set_attr(p=10, design=['p'])

nw.add_conns(c1, c2, c3)

# %% Solving Plant
nw.solve('design')
nw.save('tmp')
nw.solve('offdesign', design_path='tmp')

The user needs to make sure the input parameters for the offdesign case match. Best regards

govind-menon110 commented 3 years ago

Will try and update. Thanks @fwitte

fwitte commented 3 years ago

@govind-menon110, you might be interested in joining us in the upcoming oemof user meeting: https://github.com/oemof/oemof/issues/93.

govind-menon110 commented 3 years ago

@fwitte thank you so much, how do I join?

fwitte commented 3 years ago

For now, just comment that you like to participate. We will publish a schedule and more detailed information in a separate wiki https://github.com/oemof/oemof/wiki.

govind-menon110 commented 3 years ago

Hi @fwitte, Actually I have tried multiple times to form the off design case for the power plant but failed. It works on a smaller scale but as soon as I want to add Splitter, Merge and a few more components, it isn't able to converge. Is there a talk that has taken place in the past on how to form such a case? Unable to get a hang of it. Also the above simple example that you provided works for the offdesign only if all the parameters are exactly the same. If I change for example 'm' a little, it pushes the fluid out of range and fails.

fwitte commented 3 years ago

For completion, solution is to fix the lowest turbine stage outlet pressure for offdesign:

# %% Importing modules

from tespy.networks import Network
from tespy.connections import Connection
from tespy.components import Sink, Source
from tespy.components import HeatExchangerSimple, Turbine
from tespy.components import Merge, Splitter

import numpy as np

# %% Network

fluid_list = ['BICUBIC::H2O']
nw = Network(fluids=fluid_list, T_unit='C', p_unit='bar', h_unit='kJ / kg')

# %% Starting the components

wtr = Source("Water")
sink_1 = Sink("HPT Sink")
sink_2 = Sink("IPT 1 Sink")
sink_3 = Sink("LPT 1 Sink")
sink_4 = Sink("LPT 2 Sink")
sink_5 = Sink("LPT 3 Sink")
sink_6 = Sink("LPT 4 Sink")

heat_ex_1 = HeatExchangerSimple('Evaporator and Superheater', pr=0.98, design=['pr'], offdesign=["zeta"])
heat_ex_2 = HeatExchangerSimple("Reheater", pr=0.98, design=['pr'], offdesign=["zeta"])

turbine_1 = Turbine('HPT', eta_s=0.89, design=["eta_s"], offdesign=["eta_s_char", "cone"])
turbine_2 = Turbine('IPT 1', eta_s=0.903, design=["eta_s"], offdesign=["eta_s_char", "cone"])
turbine_3 = Turbine('IPT 2', eta_s=0.903, design=["eta_s"], offdesign=["eta_s_char", "cone"])
turbine_4 = Turbine('LPT 1', eta_s=0.85,design=["eta_s"], offdesign=["eta_s_char", "cone"])
turbine_5 = Turbine('LPT 2', eta_s=0.85, design=["eta_s"], offdesign=["eta_s_char", "cone"])
turbine_6 = Turbine('LPT 3', eta_s=0.85, design=["eta_s"], offdesign=["eta_s_char", "cone"])
turbine_7 = Turbine('LPT 4', eta_s=0.85, design=["eta_s"], offdesign=["eta_s_char", "cone"])

spl_1 = Splitter('HPT Splitter')
spl_2 = Splitter('IPT Splitter')
spl_3 = Splitter('LPT 1 Splitter')
spl_4 = Splitter('LPT 2 Splitter')
spl_5 = Splitter('LPT 3 Splitter')

# %% Building Connections
water = Connection(wtr, 'out1', heat_ex_1, 'in1')
evap_hpt = Connection(heat_ex_1, 'out1', turbine_1, 'in1')
hpt_spl = Connection(turbine_1, 'out1', spl_1, 'in1')
spl1_water = Connection(spl_1, 'out1', sink_1, 'in1')
spl1_reh = Connection(spl_1, 'out2', heat_ex_2, 'in1')
reh_ipt = Connection(heat_ex_2, 'out1', turbine_2, 'in1')
ipt1_spl = Connection(turbine_2, 'out1', spl_2, 'in1')
spl2_water = Connection(spl_2, 'out1', sink_2, 'in1')
spl2_ipt2 = Connection(spl_2, 'out2', turbine_3, 'in1')
ipt2_lpt1 = Connection(turbine_3, 'out1', turbine_4, 'in1')
lpt1_spl = Connection(turbine_4, 'out1', spl_3, 'in1')
spl3_water = Connection(spl_3, 'out1', sink_3, 'in1')
spl3_lpt2 = Connection(spl_3, 'out2', turbine_5, 'in1')
lpt2_spl = Connection(turbine_5, 'out1', spl_4, 'in1')
spl4_water = Connection(spl_4, 'out1', sink_4, 'in1')
spl4_lpt3 = Connection(spl_4, 'out2', turbine_6, 'in1')
lpt3_spl = Connection(turbine_6, 'out1', spl_5, 'in1')
spl5_water = Connection(spl_5, 'out1', sink_5, 'in1')
spl5_lpt4 = Connection(spl_5, 'out2', turbine_7, 'in1')
lpt4_water = Connection(turbine_7, 'out1', sink_6, 'in1')

# %% Setting Connection Params
water.set_attr(T=170, m=425.8, fluid={'H2O': 1})
evap_hpt.set_attr(p=187, T=537, design=["p"])
hpt_spl.set_attr(p=44.1, design=["p"])
spl1_water.set_attr(m=45)
reh_ipt.set_attr(T=537)
ipt1_spl.set_attr(p=17.4, design=["p"])
spl2_water.set_attr(m=40)
ipt2_lpt1.set_attr(p=7.93, design=["p"])
lpt1_spl.set_attr(p=2.6, design=["p"])
spl3_water.set_attr(m=24)
lpt2_spl.set_attr(p=1.3, design=["p"])
lpt3_spl.set_attr(p=0.26, design=["p"])
spl4_water.set_attr(m=15)
spl5_water.set_attr(m=10)
lpt4_water.set_attr(p=0.103)

# %% Adding Connections to Network

nw.add_conns(water, evap_hpt,hpt_spl, spl1_water, spl1_reh, reh_ipt,
            ipt1_spl, spl2_water, spl2_ipt2, ipt2_lpt1, lpt1_spl,
            spl3_water, spl3_lpt2, lpt2_spl, spl4_water, spl4_lpt3,
            lpt3_spl, spl5_water, spl5_lpt4, lpt4_water)

# %% Solving the plant

nw.solve('design')
nw.save('./srinath')
nw.print_results()
# %% Turbine Power
print(turbine_1.P.val + turbine_2.P.val + turbine_3.P.val + turbine_4.P.val
        + turbine_5.P.val + turbine_6.P.val + turbine_7.P.val)

# %% Offdesign
for m in np.linspace(24, 34, 11):
    spl3_water.set_attr(m=m)
    nw.solve(mode='offdesign', design_path='./srinath', max_iter=75)
    print(lpt4_water.h.val)
govind-menon110 commented 3 years ago

Hi @fwitte , Thank you so much. Works now. I'll put all my doubts and learnings in a doc and share when I finish my project. Probably that'll help people that want to start off with tespy.

Regards, Govind

fwitte commented 3 years ago

Actually, I sent you an email with the above solution on the 18th of April with a short explanation. Did you receive it?

fwitte commented 3 years ago

Hi Govind,

usually, the pressure at the lowest turbine stage is not a result of the cone law, but rather a result of the cooling water specs at the condenser. The cone law does take the outlet pressure into account, but it has a very small influence (see https://tespy.readthedocs.io/en/master/api/tespy.components.html#tespy.components.turbomachinery.turbine.Turbine.cone_func), and thus it is very difficult for the solver to find a solution (if even possible) in this case. Usually, the inlet pressure is subject to change by the cone law in power plant design, therefore I suggest changing the script in the way I did in the script appended. But I am not an expert on this particular subject, I would need to consult someone myself to be fully sure.

BTW: If your simulations will be run-time critical you might want to check, if the results you get with BICUBIC::water instead of HEOS::water (default CoolProp back-end) are consistent (they should be) because the solver is significantly faster. I did so in the script.

Best regards

Francesco

fwitte commented 3 years ago

Thank you so much. Works now. I'll put all my doubts and learnings in a doc and share when I finish my project. Probably that'll help people that want to start off with tespy.

Some kind of FAQ for beginner issues could be nice. Also is related to https://github.com/oemof/oemof/issues/94.

govind-menon110 commented 3 years ago

Your email did not reach my inbox @fwitte. Thank you for providing the email body here :) Will read and understand. It'll also be helpful for the others I feel.

Regards, Govind