Market Analysis with Integrated Ultra-supercritical Power Plant Model: Pricetaker Assumption

Author: Naresh Susarla (naresh.susarla@netl.doe.gov)

This notebook presents a market analysis considering the pricetaker assumption for the integrated system of ultra-supercritical power plant and Solar salt as the storage fluid in thermal energy storage system. The electricity prices or LMP (locational marginal prices) are given and remain constant through the assumed time horizon. The prices used in this study are obtained from a synthetic database, i.e. RTS-GMLC.

To start this analysis, first import all packages and libraries required.

[1]:
from pyomo.environ import Param, Objective, Expression, SolverFactory, value, Constraint
import numpy as np
from dispatches.case_studies.fossil_case.ultra_supercritical_plant.storage.\
    multiperiod_integrated_storage_usc import create_multiperiod_usc_model

# For plots
from matplotlib import pyplot as plt
import matplotlib
matplotlib.rc('font', size=24)
plt.rc('axes', titlesize=24)
WARNING: DEPRECATED: The 'pyomo.contrib.incidence_analysis.util' module has
    been moved to 'pyomo.contrib.incidence_analysis.scc_solver'. However, we
    recommend importing this functionality (e.g.
    solve_strongly_connected_components) directly from
    'pyomo.contrib.incidence_analysis'.  (deprecated in 6.5.0) (called from
    c:\gmlc_fe\source_code\idaes-
    pse\idaes\core\initialization\block_triangularization.py:18)

Get the LMP data from a database or define a list of electricity prices. In this analysis, we define a list and populate it with custom electricity prices. This list is named price.

[2]:
# Select LMP source data and scaling factor
use_rts_data = False
use_mod_rts_data = True
if use_rts_data:
    print('>>>>>> Using RTS lmp data')
    with open('rts_results_all_prices_base_case.npy', 'rb') as f:
        price = np.load(f)
elif use_mod_rts_data:
    price = [22.9684, 21.1168, 20.4, 20.419,
             20.419, 21.2877, 23.07, 25,
             18.4634, 0, 0, 0,
             0, 0, 0, 0,
             19.0342, 23.07, 200, 200,
             200, 200, 200, 200]
else:
    print('>>>>>> Using NREL lmp data')
    price = np.load("nrel_scenario_average_hourly.npy")

if use_rts_data:
    lmp = price[0:number_hours].tolist()
elif use_mod_rts_data:
    lmp = price

Define the time horizon and period for the market analysis using the pricetaker assumption and add the initial status of the Solar salt storage tank.

[3]:
# Add number of days and hours per week
ndays = 1
nweeks = 1
number_hours = 24 * ndays
n_time_points = nweeks * number_hours

# Add status for storage tanks in kg
tank_status = "hot_empty"
tank_min = 1
tank_max = 6739292

Construct a time-indexed multiperiod model using the multiperiod method with the integrated ultra-supercritical power plant. During this step, the model is initialized for each time period.

[4]:
# Create the multiperiod model object
m = create_multiperiod_usc_model(
    n_time_points=n_time_points, pmin=None, pmax=None
)

# Retrieve pyomo model and active time blocks
blks = [m.period[t] for t in m.set_period]

# Add periodic linking constraint
m.periodic_variable_pair_costraint = Constraint(
    expr=(m.period[24].fs.salt_inventory_hot ==
         m.period[1].fs.previous_salt_inventory_hot))

[+   0.00] Beginning the formulation of the multiperiod optimization problem.
2023-03-28 13:53:20 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[1]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[1]
2023-03-28 13:53:29 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[2]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[2]
2023-03-28 13:53:38 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[3]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[3]
2023-03-28 13:53:48 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[4]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[4]
2023-03-28 13:53:55 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[5]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[5]
2023-03-28 13:54:03 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[6]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[6]
2023-03-28 13:54:11 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[7]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[7]
2023-03-28 13:54:18 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[8]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[8]
2023-03-28 13:54:25 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[9]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[9]
2023-03-28 13:54:32 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[10]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[10]
2023-03-28 13:54:39 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[11]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[11]
2023-03-28 13:54:45 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[12]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[12]
2023-03-28 13:54:53 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[13]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[13]
2023-03-28 13:54:59 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[14]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[14]
2023-03-28 13:55:06 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[15]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[15]
2023-03-28 13:55:14 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[16]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[16]
2023-03-28 13:55:20 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[17]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[17]
2023-03-28 13:55:27 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[18]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[18]
2023-03-28 13:55:36 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[19]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[19]
2023-03-28 13:55:43 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[20]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[20]
2023-03-28 13:55:50 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[21]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[21]
2023-03-28 13:55:57 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[22]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[22]
2023-03-28 13:56:05 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[23]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[23]
2023-03-28 13:56:13 [INFO] idaes.apps.grid_integration.multiperiod.multiperiod: ...Constructing the flowsheet model for period[24]
INFO:idaes.apps.grid_integration.multiperiod.multiperiod:...Constructing the flowsheet model for period[24]
[+ 180.71] Completed the formulation of the multiperiod optimization problem.
2023-03-28 13:56:34 [INFO] idaes.init.fs.boiler.control_volume: Initialization Complete
2023-03-28 13:56:35 [INFO] idaes.init.fs.boiler: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:35 [INFO] idaes.init.fs.turbine_splitter[1]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:35 [INFO] idaes.init.fs.turbine_splitter[2]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:35 [INFO] idaes.init.fs.reheater[1].control_volume: Initialization Complete
2023-03-28 13:56:35 [INFO] idaes.init.fs.reheater[1]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:35 [INFO] idaes.init.fs.turbine_splitter[3]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:35 [INFO] idaes.init.fs.turbine_splitter[4]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:35 [INFO] idaes.init.fs.reheater[2].control_volume: Initialization Complete
2023-03-28 13:56:35 [INFO] idaes.init.fs.reheater[2]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:35 [INFO] idaes.init.fs.turbine_splitter[5]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:36 [INFO] idaes.init.fs.turbine_splitter[6]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:36 [INFO] idaes.init.fs.turbine_splitter[7]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:36 [INFO] idaes.init.fs.turbine_splitter[8]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:36 [INFO] idaes.init.fs.turbine_splitter[9]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:36 [INFO] idaes.init.fs.turbine_splitter[10]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:37 [INFO] idaes.init.fs.condenser_mix: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:37 [INFO] idaes.init.fs.condenser.control_volume: Initialization Complete
2023-03-28 13:56:37 [INFO] idaes.init.fs.condenser: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh_mixer[1]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh[1].hot_side: Initialization Complete
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh[1].cold_side: Initialization Complete
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh[1]: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh_mixer[2]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh[2].hot_side: Initialization Complete
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh[2].cold_side: Initialization Complete
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh[2]: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh_mixer[3]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh[3].hot_side: Initialization Complete
2023-03-28 13:56:37 [INFO] idaes.init.fs.fwh[3].cold_side: Initialization Complete
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[3]: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh_mixer[4]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[4].hot_side: Initialization Complete
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[4].cold_side: Initialization Complete
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[4]: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[5].hot_side: Initialization Complete
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[5].cold_side: Initialization Complete
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[5]: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:38 [INFO] idaes.init.fs.deaerator: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh_mixer[6]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[6].hot_side: Initialization Complete
2023-03-28 13:56:38 [INFO] idaes.init.fs.fwh[6].cold_side: Initialization Complete
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[6]: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh_mixer[7]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[7].hot_side: Initialization Complete
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[7].cold_side: Initialization Complete
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[7]: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh_mixer[8]: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[8].hot_side: Initialization Complete
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[8].cold_side: Initialization Complete
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[8]: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[9].hot_side: Initialization Complete
2023-03-28 13:56:39 [INFO] idaes.init.fs.fwh[9].cold_side: Initialization Complete
2023-03-28 13:56:40 [INFO] idaes.init.fs.fwh[9]: Initialization Completed, optimal - Optimal Solution Found
Model Initialization =  optimal
*******************  USC Model Initialized   ********************
2023-03-28 13:56:43 [INFO] idaes.init.fs.ess_hp_split: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:43 [INFO] idaes.init.fs.hxc.hot_side: Initialization Complete
2023-03-28 13:56:43 [INFO] idaes.init.dispatches.properties.solarsalt_properties: fs.hxc.cold_side.properties_in Initialisation Step 1 Complete.
2023-03-28 13:56:43 [INFO] idaes.init.dispatches.properties.solarsalt_properties: Initialization Step 1 Complete.
2023-03-28 13:56:43 [INFO] idaes.init.dispatches.properties.solarsalt_properties: fs.hxc.cold_side.properties_out Initialisation Step 1 Complete.
2023-03-28 13:56:43 [INFO] idaes.init.dispatches.properties.solarsalt_properties: Initialization Step 1 Complete.
2023-03-28 13:56:43 [INFO] idaes.init.dispatches.properties.solarsalt_properties: State Released.
2023-03-28 13:56:43 [INFO] idaes.init.fs.hxc.cold_side: Initialization Complete
2023-03-28 13:56:43 [INFO] idaes.init.dispatches.properties.solarsalt_properties: State Released.
2023-03-28 13:56:43 [INFO] idaes.init.fs.hxc: Initialization Completed, optimal - Optimal Solution Found
2023-03-28 13:56:43 [INFO] idaes.init.fs.cooler.control_volume: Initialization Complete
2023-03-28 13:56:43 [INFO] idaes.init.fs.cooler: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:43 [INFO] idaes.init.fs.hx_pump.control_volume: Initialization Complete
2023-03-28 13:56:44 [INFO] idaes.init.fs.hx_pump: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:44 [INFO] idaes.init.fs.ess_bfp_split: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:44 [INFO] idaes.init.fs.recycle_mixer: Initialization Complete: optimal - Optimal Solution Found
2023-03-28 13:56:44 [INFO] idaes.init.dispatches.properties.solarsalt_properties: fs.hxd.hot_side.properties_in Initialisation Step 1 Complete.
2023-03-28 13:56:44 [INFO] idaes.init.dispatches.properties.solarsalt_properties: Initialization Step 1 Complete.
2023-03-28 13:56:44 [INFO] idaes.init.dispatches.properties.solarsalt_properties: fs.hxd.hot_side.properties_out Initialisation Step 1 Complete.
2023-03-28 13:56:44 [INFO] idaes.init.dispatches.properties.solarsalt_properties: Initialization Step 1 Complete.
2023-03-28 13:56:44 [INFO] idaes.init.dispatches.properties.solarsalt_properties: State Released.
2023-03-28 13:56:44 [INFO] idaes.init.fs.hxd.hot_side: Initialization Complete
2023-03-28 13:56:44 [INFO] idaes.init.fs.hxd.cold_side: Initialization Complete
2023-03-28 13:56:44 [INFO] idaes.init.dispatches.properties.solarsalt_properties: State Released.
2023-03-28 13:56:44 [INFO] idaes.init.fs.hxd: Initialization Completed, optimal - Optimal Solution Found
Integrated Model Initialization =  optimal
***************   Integrated Model Initialized   ***************
Cost Initialization =  optimal
******************** Costing Initialized *************************


[+  28.16] Created an instance of the flowsheet and initialized it.
[+   5.37] Initialized the entire multiperiod optimization model.
[+   0.02] Unfixed the degrees of freedom from each period model.

Use the LMP data from the price list above to assign the electricity prices to each time period. Also, include revenue and operating cost functions to calculate a total cost.

[5]:
count = 0
for blk in blks:
    blk.lmp_signal = Param(default=0, mutable=True)
    blk.revenue = lmp[count]*blk.fs.net_power
    blk.operating_cost = Expression(
        expr=(
            (blk.fs.operating_cost
             + blk.fs.plant_fixed_operating_cost
             + blk.fs.plant_variable_operating_cost) / (365 * 24)
        )
    )
    blk.cost = Expression(expr=-(blk.revenue - blk.operating_cost))
    count += 1

Add the total cost expression as the objective function.

[6]:
m.obj = Objective(expr=sum([blk.cost for blk in blks]))

Add initial state for the two linking variables in the model: the Solar salt tank for different scenarios and the power.

[7]:
# Initial state for Solar salt tank for differrent tank scenarios:"hot_empty","hot_full","hot_half_full"
if tank_status == "hot_empty":
    blks[0].fs.previous_salt_inventory_hot.fix(1103053.48)
    blks[0].fs.previous_salt_inventory_cold.fix(tank_max-1103053.48)
elif tank_status == "half_full":
    blks[0].fs.previous_salt_inventory_hot.fix(tank_max/2)
    blks[0].fs.previous_salt_inventory_cold.fix(tank_max/2)
elif tank_status == "hot_full":
    blks[0].fs.previous_salt_inventory_hot.fix(tank_max-tank_min)
    blks[0].fs.previous_salt_inventory_cold.fix(tank_min)
else:
    print("Unrecognized scenario! Try hot_empty, hot_full, or half_full")

# Initial state for power linking variable
blks[0].fs.previous_power.fix(447.66)

Finally, solve the entire multi-period model and save the results in lists to be used later to plot the results.

[8]:
# Declare NLP solver
opt = SolverFactory('ipopt')

# Declare lists to save results
hot_tank_level = []
cold_tank_level = []
net_power = []
hxc_duty = []
hxd_duty = []
for week in range(nweeks):
    print()
    print(">>>>>> Solving for week {}: {} hours of operation in {} day(s) "
          .format(week + 1, number_hours, ndays))

    # Solve the multi-period model
    opt.solve(m, tee=True)

    hot_tank_level.append([(value(blks[i].fs.salt_inventory_hot)) * 1e-3
                           for i in range(n_time_points)])
    cold_tank_level.append([(value(blks[i].fs.salt_inventory_cold)) * 1e-3
                            for i in range(n_time_points)])
    net_power.append([value(blks[i].fs.net_power)
                      for i in range(n_time_points)])
    hxc_duty.append([value(blks[i].fs.hxc.heat_duty[0]) * 1e-6
                     for i in range(n_time_points)])
    hxd_duty.append([value(blks[i].fs.hxd.heat_duty[0]) * 1e-6
                     for i in range(n_time_points)])

>>>>>> Solving for week 1: 24 hours of operation in 1 day(s)
Ipopt 3.13.2:

******************************************************************************
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...:    38348
Number of nonzeros in inequality constraint Jacobian.:      284
Number of nonzeros in Lagrangian Hessian.............:    14520

Total number of variables............................:    14685
                     variables with only lower bounds:      144
                variables with lower and upper bounds:    10053
                     variables with only upper bounds:        0
Total number of equality constraints.................:    14567
Total number of inequality constraints...............:      168
        inequality constraints with only lower bounds:       24
   inequality constraints with lower and upper bounds:        0
        inequality constraints with only upper bounds:      144

iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
   0 -3.7718386e+05 6.50e+07 1.00e+02  -1.0 0.00e+00    -  0.00e+00 0.00e+00   0
   1 -3.7219484e+05 4.05e+07 3.06e+03  -1.0 1.93e+07    -  2.70e-01 4.26e-01h  1
   2 -3.6953074e+05 2.36e+07 1.67e+03  -1.0 1.04e+07    -  5.80e-01 4.48e-01h  1
   3 -3.6941476e+05 2.19e+07 1.54e+03  -1.0 6.16e+06    -  7.68e-01 7.18e-02h  1
   4 -3.7476158e+05 1.39e+06 2.64e+03  -1.0 7.33e+06    -  6.99e-01 1.00e+00f  1
   5 -3.9019830e+05 3.95e+05 6.37e+02  -1.0 2.29e+07    -  5.52e-01 7.21e-01f  1
   6 -3.9081574e+05 3.88e+05 6.27e+02  -1.0 4.31e+07    -  2.24e-02 1.68e-02f  1
   7 -3.9163604e+05 3.78e+05 6.09e+02  -1.0 3.71e+07    -  8.28e-01 2.75e-02f  1
   8 -4.2481146e+05 6.44e+05 4.10e+02  -1.0 2.68e+08    -  1.82e-01 1.76e-01f  1
   9 -4.2517311e+05 6.43e+05 4.09e+02  -1.0 3.27e+08    -  1.70e-02 1.67e-03f  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  10 -4.2528820e+05 6.41e+05 4.07e+02  -1.0 6.85e+07    -  4.45e-01 3.80e-03f  1
  11 -4.2753226e+05 5.42e+05 3.45e+02  -1.0 5.76e+07    -  2.47e-02 1.54e-01f  1
  12 -4.2799776e+05 5.16e+05 3.28e+02  -1.0 5.55e+07    -  6.64e-01 4.80e-02f  1
  13 -4.3006960e+05 4.31e+05 2.74e+02  -1.0 1.07e+08    -  1.42e-01 1.65e-01f  1
  14 -4.3130567e+05 3.89e+05 2.47e+02  -1.0 1.20e+08    -  8.88e-01 9.74e-02f  1
  15 -4.3186028e+05 4.50e+05 2.46e+02  -1.0 7.55e+09    -  1.62e-02 5.28e-03f  1
  16 -4.3268106e+05 4.45e+05 2.43e+02  -1.0 1.72e+09    -  2.21e-02 1.16e-02f  1
  17 -4.3375779e+05 4.38e+05 2.39e+02  -1.0 1.48e+09    -  7.74e-02 1.58e-02f  1
  18 -4.3544275e+05 4.25e+05 2.32e+02  -1.0 8.01e+08    -  8.55e-02 2.91e-02f  1
  19 -4.3740807e+05 4.09e+05 2.24e+02  -1.0 5.57e+08    -  1.65e-01 3.65e-02f  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  20 -4.4205866e+05 9.97e+05 2.05e+02  -1.0 5.07e+08    -  1.68e-01 8.65e-02f  1
  21 -4.4231008e+05 9.95e+05 2.04e+02  -1.0 4.83e+08    -  2.93e-01 4.79e-03f  1
  22 -4.4719532e+05 3.10e+06 1.86e+02  -1.0 4.73e+08    -  8.29e-02 8.55e-02f  1
  23 -4.4758765e+05 3.10e+06 1.84e+02  -1.0 1.13e+08    -  1.22e-03 1.21e-02f  1
  24 -4.4845636e+05 3.23e+06 1.79e+02  -1.0 4.36e+08    -  5.13e-02 2.73e-02f  1
  25 -4.5168962e+05 5.43e+06 1.62e+02  -1.0 4.21e+08    -  2.71e-03 9.36e-02f  1
  26 -4.5185267e+05 5.41e+06 1.61e+02  -1.0 4.13e+08    -  1.66e-01 5.23e-03f  1
  27 -4.5364183e+05 6.01e+06 1.52e+02  -1.0 2.70e+08    -  3.78e-01 5.76e-02f  1
  28 -4.5770806e+05 1.06e+07 1.31e+02  -1.0 1.12e+08    -  4.86e-01 1.40e-01f  1
  29 -4.7540023e+05 1.01e+08 4.36e+01  -1.0 8.52e+07    -  6.12e-01 7.75e-01f  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  30 -4.7693908e+05 1.64e+07 1.60e+01  -1.0 2.40e+07    -  8.20e-01 8.51e-01f  1
  31 -4.7830655e+05 2.14e+07 6.29e+00  -1.0 3.58e+07    -  4.58e-01 1.00e+00f  1
  32 -4.7842684e+05 1.07e+07 1.19e+00  -1.0 1.39e+07    -  9.28e-01 1.00e+00h  1
  33 -4.7847765e+05 8.86e+06 1.13e-01  -1.0 1.45e+07    -  1.00e+00 1.00e+00h  1
  34 -4.7848098e+05 3.21e+06 1.12e-02  -1.0 5.53e+06    -  1.00e+00 1.00e+00h  1
  35 -4.7848125e+05 3.31e+05 2.02e-03  -1.0 1.21e+06    -  1.00e+00 1.00e+00h  1
  36 -4.7848127e+05 3.67e+03 1.85e-03  -1.0 1.15e+05    -  1.00e+00 1.00e+00h  1
  37 -4.7849227e+05 1.23e+05 1.11e+01  -2.5 1.70e+06    -  1.00e+00 7.52e-01f  1
  38 -4.7849558e+05 4.71e+02 1.98e-03  -2.5 1.56e+05    -  1.00e+00 1.00e+00h  1
  39 -4.7849557e+05 7.10e-02 8.03e-03  -2.5 1.02e+03    -  1.00e+00 1.00e+00h  1
iter    objective    inf_pr   inf_du lg(mu)  ||d||  lg(rg) alpha_du alpha_pr  ls
  40 -4.7849596e+05 4.94e+01 7.96e-03  -3.8 2.44e+04    -  1.00e+00 1.00e+00f  1
  41 -4.7849596e+05 3.19e-03 1.51e-09  -3.8 2.21e+02    -  1.00e+00 1.00e+00h  1
  42 -4.7849598e+05 1.47e-01 4.75e-08  -5.7 1.33e+03    -  1.00e+00 1.00e+00f  1
  43 -4.7849598e+05 2.09e-05 1.79e-11  -8.6 1.58e+01    -  1.00e+00 1.00e+00h  1

Number of Iterations....: 43

                                   (scaled)                 (unscaled)
Objective...............:  -2.3924799034252550e+05   -4.7849598068505101e+05
Dual infeasibility......:   1.7907030466656758e-11    3.5814060933313515e-11
Constraint violation....:   4.4237822294235229e-09    2.0943582057952881e-05
Complementarity.........:   2.5088487928275600e-09    5.0176975856551201e-09
Overall NLP error.......:   4.4237822294235229e-09    2.0943582057952881e-05


Number of objective function evaluations             = 44
Number of objective gradient evaluations             = 44
Number of equality constraint evaluations            = 44
Number of inequality constraint evaluations          = 44
Number of equality constraint Jacobian evaluations   = 44
Number of inequality constraint Jacobian evaluations = 44
Number of Lagrangian Hessian evaluations             = 43
Total CPU secs in IPOPT (w/o function evaluations)   =      2.468
Total CPU secs in NLP function evaluations           =    117.869

EXIT: Optimal Solution Found.

Add maximum and minimum values for the power plant power production and storage heat duty. This is needed when plotting the results.

[9]:
max_power = 436
max_power_storage = 30
max_power_total = max_power + max_power_storage
min_storage_heat_duty = 10
max_storage_heat_duty = 200

In Figure 1, plot the Solar salt levels in the storage tank for the entire time horizon.

[10]:
hours = np.arange(n_time_points * nweeks)
lmp_array = np.asarray(lmp[0:n_time_points])
hot_tank_array = np.asarray(hot_tank_level[0:nweeks]).flatten()
cold_tank_array = np.asarray(cold_tank_level[0:nweeks]).flatten()

# Convert array to list to include hot tank level at time zero
hot_tank_array0 = value(blks[0].fs.previous_salt_inventory_hot) * 1e-3
cold_tank_array0 = value(blks[0].fs.previous_salt_inventory_cold) * 1e-3
hours_list = hours.tolist() + [n_time_points]
hot_tank_list = [hot_tank_array0] + hot_tank_array.tolist()
cold_tank_list = [cold_tank_array0] + cold_tank_array.tolist()

# Declare settings for plot
font = {'size': 16}
plt.rc('font', **font)
fig1, ax1 = plt.subplots(figsize=(12, 8))
color = ['r', 'b', 'tab:green', 'k', 'tab:orange']

ax1.set_xlabel('Time Period (hr)')
ax1.set_ylabel('Salt Tank Level (metric ton)',
               color=color[3])
ax1.spines["top"].set_visible(False)
ax1.spines["right"].set_visible(False)
ax1.grid(linestyle=':', which='both',
         color='gray', alpha=0.30)
plt.axhline(tank_max * 1e-3,
            ls=':', lw=1.75,
            color=color[4])
plt.text(n_time_points / 2 - 1.5,
         tank_max * 1e-3 + 100,
         'max salt',
         color=color[4])
ax1.step(hours_list, hot_tank_list,
         marker='^', ms=4,
         lw=1, color=color[0],
         label='Hot Salt')
ax1.step(hours_list, cold_tank_list,
         marker='v', ms=4,
         lw=1, color=color[1],
         label='Cold Salt')
ax1.legend(loc="center right", frameon=False)
ax1.tick_params(axis='y')
ax1.set_xticks(np.arange(0, n_time_points * nweeks + 1, step=2))

ax2 = ax1.twinx()
ax2.set_ylabel('LMP ($/MWh)',
               color=color[2])
ax2.step([x + 1 for x in hours], lmp_array,
         marker='o', ms=3, alpha=0.5,
         ls='-', lw=1,
         color=color[2])
ax2.tick_params(axis='y',
                labelcolor=color[2])
../_images/examples_pricetaker_integrated_ultrasupercritical_power_plant_22_0.png

In Figure 2, plot the operating profile of the power plant in terms of total power for the entire time horizon.

[11]:
font = {'size': 18}
plt.rc('font', **font)

power_array = np.asarray(net_power[0:nweeks]).flatten()

# Convert array to list to include net power at time zero
power_array0 = value(blks[0].fs.previous_power)
power_list = [power_array0] + power_array.tolist()

fig2, ax3 = plt.subplots(figsize=(12, 8))
ax3.set_xlabel('Time Period (hr)')
ax3.set_ylabel('Net Power Output (MW)',
               color=color[1])
ax3.spines["top"].set_visible(False)
ax3.spines["right"].set_visible(False)
ax3.grid(linestyle=':', which='both',
         color='gray', alpha=0.30)
plt.text(n_time_points / 2 - 3, max_power - 5.5,
         'max plant power',
         color=color[4])
plt.axhline(max_power,
            ls='-.', lw=1.75,
            color=color[4])
ax3.step(hours_list, power_list,
         marker='o', ms=4,
         lw=1, color=color[1])
ax3.tick_params(axis='y',
                labelcolor=color[1])
ax3.set_xticks(np.arange(0, n_time_points * nweeks + 1, step=2))

ax4 = ax3.twinx()
ax4.set_ylabel('LMP ($/MWh)',
               color=color[2])
ax4.step([x + 1 for x in hours], lmp_array,
         marker='o', ms=3, alpha=0.5,
         ls='-', lw=1,
         color=color[2])
ax4.tick_params(axis='y',
                labelcolor=color[2])
../_images/examples_pricetaker_integrated_ultrasupercritical_power_plant_24_0.png

In Figure 3, plot the storage heat exchanger operating profiles in terms of heat duties for the entire time horizon.

[12]:
zero_point = True
hxc_array = np.asarray(hxc_duty[0:nweeks]).flatten()
hxd_array = np.asarray(hxd_duty[0:nweeks]).flatten()
hxc_duty0 = 0
hxc_duty_list = [hxc_duty0] + hxc_array.tolist()
hxd_duty0 = 0
hxd_duty_list = [hxd_duty0] + hxd_array.tolist()

fig3, ax5 = plt.subplots(figsize=(12, 8))
ax5.set_xlabel('Time Period (hr)')
ax5.set_ylabel('Storage Heat Duty (MW)',
               color=color[3])
ax5.spines["top"].set_visible(False)
ax5.spines["right"].set_visible(False)
ax5.grid(linestyle=':', which='both',
         color='gray', alpha=0.30)
plt.text(n_time_points / 2 - 2.2,
         max_storage_heat_duty + 1,
         'max storage',
         color=color[4])
plt.text(n_time_points / 2 - 2,
         min_storage_heat_duty - 6.5,
         'min storage',
         color=color[4])
plt.axhline(max_storage_heat_duty,
            ls=':', lw=1.75,
            color=color[4])
plt.axhline(min_storage_heat_duty,
            ls=':', lw=1.75,
            color=color[4])
if zero_point:
    ax5.step(hours_list, hxc_duty_list,
             marker='^', ms=4,
             color=color[0], label='Charge',)
    ax5.step(hours_list, hxd_duty_list,
             marker='v', ms=4,
             color=color[1], label='Discharge')
else:
    ax5.step([x + 1 for x in hours], hxc_array,
             marker='^', ms=4, lw=1,
             color=color[0], label='Charge')
    ax5.step([x + 1 for x in hours], hxd_array,
             marker='v', ms=4, lw=1,
             color=color[1], label='Discharge')
ax5.legend(loc="center right", frameon=False)
ax5.tick_params(axis='y',
                labelcolor=color[3])
ax5.set_xticks(np.arange(0, n_time_points * nweeks + 1, step=2))

ax6 = ax5.twinx()
ax6.set_ylabel('LMP ($/MWh)',
               color=color[2])
ax6.step([x + 1 for x in hours], lmp_array,
         marker='o', ms=3,
         alpha=0.5, ls='-',
         color=color[2])
ax6.tick_params(axis='y',
                labelcolor=color[2])
../_images/examples_pricetaker_integrated_ultrasupercritical_power_plant_26_0.png