Next Article in Journal
Reflection-Like Maps in High-Dimensional Euclidean Space
Next Article in Special Issue
Multiscale Model Reduction of the Unsaturated Flow Problem in Heterogeneous Porous Media with Rough Surface Topography
Previous Article in Journal
Erratum: Al-Marzouki, S., et al. Type II Topp Leone Power Lomax Distribution with Applications. Mathematics 2020, 8, 4
Previous Article in Special Issue
Linear Model Predictive Control for a Coupled CSTR and Axial Dispersion Tubular Reactor with Recycle
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Dynamic Optimal Dispatch of Energy Systems with Intermittent Renewables and Damage Model

Department of Chemical and Biomedical Engineering, West Virginia University, P.O. Box 6102, Morgantown, WV 26506, USA
*
Author to whom correspondence should be addressed.
Submission received: 18 April 2020 / Revised: 18 May 2020 / Accepted: 21 May 2020 / Published: 28 May 2020

Abstract

:
With the increasing penetration of intermittent renewable energy sources into the grid, there is a growing need for process systems-based strategies that integrate dispatchable and variable energy systems for supplying the demand while maintaining grid reliability. The proposed framework corresponds to a dynamic mixed-integer linear programming optimization approach that integrates coal-fired and natural gas-fired power plants, NaS batteries for energy storage, and solar/wind energy to supply the demand. This optimization approach considers an economic goal and constraints to provide power balance while maintaining the overall damage of the natural gas combined cycle (NGCC) power plant drum under a maximum stress as well as avoiding the overheating of the NGCC superheater and reheater. Renewable curtailment levels are also retained at minimum levels. Case studies are analyzed considering different loads and renewable penetration levels. The results show that the demand was met for all cases. Grid flexibility was mostly provided by the NGCC, while the batteries were used sparingly. In addition, considering a CO2 equivalent analysis, the environmental performance was intrinsically connected to grid flexibility and the level of renewable penetration. Stress analysis results reinforced the necessity for an equipment health-related constraint.

Graphical Abstract

1. Introduction

Changes in global climate, societal and cultural pressures, and a rapid decrease in the cost of renewable technologies have stimulated the increased participation of intermittent renewable energy sources in the power generation grid [1]. Analysis and projections of electricity generation fuels show an increase of natural gas and renewables participation, along with a decrease in coal usage. Specifically, by 2050, it is expected for the natural gas participation to increase from 34% to 39%, while the renewables participation to increase from 18% to 31%, and the coal participation to decrease from 28% to 17% [2]. Within the renewable energy sources, the two technologies enabling the rapid deployment are wind and solar [2]. Intermittent or variable renewable energy sources, such as solar and wind, are known to be non-dispatchable, i.e., the energy is only available when the fuel source is available and cannot be produced as required by the demand. Furthermore, solar and wind are also known to have limited alignment with electricity demand [3], e.g., the period when the availability of the renewable energy sources peaks is not necessarily the same period when the electricity demand peaks. Thus, baseload thermal power plants, such as natural gas-fired and coal-fired power plants, are often allocated to supply the net load demand [4]. The behavior of ramping up/down the baseload power plants to accommodate variability, which could have various sources such as demand variability or renewables variability, will be considered from now on in the paper as cycling of the power plant. Some of the main concerns associated with an increase in cycling of baseload power plants are the power plant minimum operating loads, environmental performance during cycling, and the potential of prohibitive ramping rates. Therefore, at high penetration rates, the integration of variable renewable energy sources into the grid becomes even more challenging.
Regarding coal, although its usage is predicted to decrease, it still has a significant participation [2] in electricity generation. The competition between coal, natural gas, and renewable sources [5] has discouraged the launching of new coal-fired power plants. Thus, there is a high possibility that most of the coal-fired power plants being utilized are originally designed to be operated in a baseload manner. Furthermore, these coal-fired power plants are mostly operated at off-design conditions to provide flexible power generation [6]. Enforcing cycling operation on designed baseload power plants may introduce inefficiencies to the systems, as these plants are no longer operating at their optimal design points. Cycling may also increase the emissions from the power plant side and decrease the power plant lifetime due to unexpected wear and tear [4,6].
One way to attenuate the increase in wear and tear in power plants would be to consider the damage that the equipment could suffer when subjected to cycling operations. Researches regarding damage models show that fatigue and creep damages can be estimated considering start-up [7] as well as partial load operations [8]. Consequently, start-up procedure or the design optimization of specific components can be determined to improve the lifetime under cycling. However, currently, there is a gap in how the wear and tear knowledge may be applied during the optimal dispatch of power plants under cycling operation.
Additionally, without a sufficiently flexible grid, thermal power plants may not be able to reduce their power output, and wind and solar generation need to be restricted [3]. For instance, a key limiting factor for deploying photovoltaics (PV) is curtailment, i.e., the generated PV energy would need to be rejected by system operators to maintain the supply/demand balance of the system [9]. As curtailment increases, the solar/wind energy offsets less fossil-based generation. Demand response and storage are enabling technologies that could reduce curtailment and facilitate higher penetrations of intermittent renewables into the grid [3,10]. In addition, if natural gas prices were to rise in the future, a large-scale energy system based on such technologies may be used in utility and industrial scales at a lower level of CO2 emissions [11]. From the grid operation side, within the energy storage technologies, sodium sulfur batteries are a promising technology. At the grid level, sodium sulfur batteries have high potential for electric storage due to their high energy density, the low cost of the reactants, high open-circuit voltage, and they are also one of the most mature technologies [12,13]. Some challenges associated with these batteries are the exclusivity of the technology (as it is currently produced by only one company), the high installation and operating cost, as well as the thermal management challenge that arises during the operation of several sodium sulfur cells [13,14]. Due to the challenges in quantifying the point at which energy storage becomes the least-costly flexibility option, evaluating the role of storage interventions in power systems with a high penetration of intermittent renewables requires continued analysis, improved data, and the implementation of advanced process systems techniques.
Flexibility requirements of the grid that facilitate the integration of high penetration levels of renewables have been a subject of recent studies [9,15,16]. The impact of cycling on individual thermal generators under high flexibility conditions is also a topic of concern [6,17,18,19,20] regarding economic, environmental, and power plant health objectives. The integration of different energy systems in a microgrid considering an economic objective as well as demand response battery installations are of interest [10,16,21]. However, there is a scarcity of studies on methods that analyze grid-level dynamic power plant models and energy storage units with different levels of variable renewable integration considering the damage of cycling on the power plants.
In this work, a mixed-integer linear programming (MILP)-based dynamic optimization strategy is explored to contribute in closing the aforementioned gap by determining the optimal dispatch and the integration of fossil-fueled power plants, energy storage units, and renewable energy production considering the damage of cycling. Natural gas-fired [17] and coal-fired [22] power plants correspond to the fossil-fueled power plants. The energy storage unit is represented by the sodium sulfur batteries [12]. The MILP optimization objective is to minimize the overall economic cost. The curtailment of variable renewable energy is kept as minimal as possible. In addition, a damage model for the stressed components in the natural gas-fired power plant is considered in order to determine the optimal ramping rates and dispatch while maintaining acceptable component life. Ultimately, the generated optimal policy of different energy modules could be sent to advanced model-based controllers [23,24].
This paper is organized into four sections. In Section 2, the energy systems models, the main constraints, and the MILP-based optimization problem statement are introduced. Then, the results of the case studies and their respective implications are presented in Section 3. Finally, the paper conclusions are summarized in Section 4.

2. Modeling and Optimization Approach

The mathematical definition of optimization is to find the best solution for a set of systems within the process constraints in a systematic manner. The quality of the solution is evaluated by the objective function. Dynamic optimization considers the process dynamics to determine the optimality [19]. In particular, the optimal dispatch problem addressed here is formulated as a dynamic MILP optimization that includes integer and continuous decision variables, and the linear dynamic models and grid power balance as some of the constraints. The general dynamic MILP optimization formulation is expressed in Equations (1)–(4):
min u Φ ( x , y , u , k ) ,
subject to:
x ( k ) = x ( k 1 ) + f ( x , y , u , k ) × Δ t
g ( x , y , u , k ) = 0 ,
h ( x , y , u , k ) 0 , ( x , u , k ) R , y     Z
in which Φ is the objective function to be minimized, f is the dynamic model that describes the system, g are the equality constraints, h are the inequality constraints, x are the state variables, u are the input/decision variables, y are the integer variables, and k is the current time step of the discretization. The posed optimization problem was solved using branch and bound in the intlinprog subroutine available in MATLAB 2017a® considering every 6-min time step dispatch scenario for a 30-day time horizon. The optimality gap tolerances were set to the subroutine default values. A perfect demand forecast of at least 6 min ahead of time is assumed. The overall framework could be easily adapted to incorporate the most updated forcing function data accordingly. The energy systems are assumed to be always committed at this level by an independent scheduler; thus, shutdown/startup procedures are not considered in this study during the dispatch.
The objective is to minimize the overall cost of dispatch defined in Equation (5):
m i n   Φ = C o s t E S R e v e n u e E S ,
in which ES represents the energy systems in the grid. The “CostES” and “RevenueES” variables are further defined in Equations (6) and (7):
C o s t E S = V a r i a b l e   o p e r a t i o n   a n d   m a i n t e n a n c e E S + C o s t   o f   f u e l E S + P e n a l t y E S + E l e c t r i c i t y   C o n s u m p t i o n E S ,
R e v e n u e E S = W h o l e s a l e   o f   g e n e r a t e d   p o w e r E S ,
The “PenaltyES” variable in Equation (6) represents the penalty from a potential policy that penalizes the curtailment of renewables or the investment cost per cycle of the battery.
The dynamic models of the energy systems considered are linearized models of the high-fidelity first-principle models of the charge and discharge models of the sodium sulfur (NaS) batteries [12], natural gas-fired combined cycle (NGCC) power plant [17,18], and supercritical pulverized coal-fired (SCPC) power plant [22], as depicted in Figure 1, Figure 2 and Figure 3, respectively. The wind and solar power generation are considered as forcing functions based on historical data [25]. Thus, the energy system set is defined as ES = {NaS batteries, NGCC, SCPC}.
The data used to obtain the parameters for the linear models of the NaS batteries for charge and discharge cycles are shown in Figure 4 and Figure 5.
The NaS battery linear charging model parameters expressed in Equations (8) and (9) are calculated by linear regression of the charging data obtained by exciting the high-fidelity model in Aspen Custom Modeler [12]:
S O D = 0.055 t + 0.46910
C a p a c i t y = 162.13   S O D + 144.6
in which SOD is the state of discharge (fractional), t is the time in hour, and Capacity is the capacity in Wh.
The NaS battery linear discharging model parameters expressed in Equations (10) and (11) are calculated by considering the charging data and adding a deviation term ( ε ) to shift the model predictions for alignment with the discharging data. The discharging data were also obtained from the high-fidelity model in Aspen Custom Modeler [12]. A limitation in the linear representation of the discharging model is observed regarding the mismatch between the prediction and the data at lower levels of SOD. This issue arises due to the phase change characteristics of NaS batteries during charging/discharging that alter their dynamic behavior [26].
S O D = 0.0771 t + 0.9027
C a p a c i t y = 177.98   S O D 70.486
The calculated R2 considering all of the linear regressions performed for the NaS batteries is higher than 95%.
The reduced-order linear model for the NGCC, obtained by linearization of a nonlinear high-fidelity model in Aspen Dynamics [17], is a state-space model as shown in Equations (12) and (13). For conciseness, parameters of the NGCC linear model are provided in the Supplementary Materials.
x N G C C ( k + 1 ) = A N G C C x N G C C ( k ) + B N G C C u N G C C ( k )
z N G C C ( k ) = C N G C C x N G C C ( k ) + D N G C C u N G C C ( k )
In Equations (12) and (13), k is the time step, zNGCC is the output vector, and (ANGCC, BNGCC, CNGCC, DNGCC) are the NGCC system matrices. The NGCC input is the natural gas flowrate (kg/h) and the outputs are the gross power output (MWe), gas turbine power output (MWe), the main steam temperature (°C), the reheat steam temperature (°C), the drum body mechanical stress (MPa), and the thermal drum body stress (MPa). In particular, the formulation for the stress components is explained later.
The reduced-order linear model for the SCPC obtained by linear regression using the generated data from the high-fidelity model in Aspen Dynamics [22] is an autoregressive with exogenous inputs (ARX) model, as shown in Equation (14):
A S C P C ( z ) z S C P C ( k ) = B S C P C ( z ) u S C P C ( k ) + e ( k )
in which k is the time step, zSCPC is the output vector, and (ASCPC, BSCPC) are the SCPC system matrices. The SCPC regressors are composed of the past and current coal feedrates (kg/s) and the output vector contains the past and current gross power outputs (kWe). Further details regarding parameters and regressors of the SCPC linear model are also provided in the Supplementary Materials. The datasets used for the linear regression of the power plants were obtained by exciting the input and recording the output responses. Linearization of the high-fidelity nonlinear models inherently introduced inaccuracies in the representation of the energy systems. In particular, the mechanical stress is underestimated at full load, and the thermal stress is underestimated at peaks. However, the overall dynamic behavior of the integrated energy system is still captured by the reduced-order model.
The nominal power generation for the NGCC is 641 MWe, for the SCPC is 600 MWe, and for the NaS batteries is 99.4 MW/845 MWh of storage capacity with the considered capacity of a single sodium sulfur cell being 102 Wh. The maximum discharge and charge rates of a single NaS cell are considered to be 12.5 Wh and 8.9 Wh, respectively. Furthermore, the optimization problem statement also includes the constraints expressed in Equations (15)–(28).
The smart-grid power balance is expressed in Equation (15).
i E S ( E l e c t r i c i t y   p r o d u c e d i E l e c t r i c i t y   c o n s u m e d i ) + R e n e w a b l e s = D e m a n d
The “Electricity consumedi” variable refers to the case in which the energy system requires electricity from the grid that cannot be supplied from within the energy system (for instance, when the battery requires electricity to charge). The “Renewables” variable considers the wind and solar power generation as forcing function profiles built from data that were obtained from the PJM Regional Transmission Organization [25]. More specifically, the data used corresponds to 30 days of hourly aggregate demand loads and the electricity price in August 2019, which were interpolated with the chosen time step of 6 min. The absolute values of the considered demand and the renewable power generation were adjusted to be proportional to the available energy sources’ nominal maximum and minimum power generations, while maintaining the load demand slopes. These actions were taken to guarantee feasibility to the current forcing functions that use 2019 information. However, further modifications to the forcing functions (the modified 2019 load and the 2050 renewables projection) might render the problem infeasible. Thus, to partially manage this issue, renewables curtailment was allowed and flexibility requirements of the grid analyzed. The proportionality factors used for the scaling of renewable generation in 2050 were based on information in Table 1.
A logical constraint using disjunctive programming is also included, as the NaS batteries cannot be charged and discharged at the same time. This constraint is expressed in Equations (16) and (17).
y Z | y d i s c h a r g e + y c h a r g e 1 ,   y { 0 , 1 }
y d i s c h a r g e / c h a r g e = { 0 1 d i s c h a r g e / c h a r g e   m o d e l   n o t   s e l e c t e d d i s c h a r g e / c h a r g e   m o d e l   s e l e c t e d
The NaS batteries range of operation considering the regression data fit is shown in Equation (18).
40 % S t a t e   o f   D i s c h a r g e 89 %
Finally, constraints regarding the damage of the drum, the superheater, and reheater are also included. Two distinct mechanisms can influence and change the material lifetime, namely creep and fatigue. Creep is due to the prolonged exposure to high temperature, while fatigue is caused by the repeated and fluctuating stress that results from both temperature and pressure transients under load-following conditions [7].
In the optimization problem, the maximum stress that is allowed for the drum is considered as a constraint along with a constraint on the maximum temperature to avoid overheating of the NGCC superheater and reheater.
Stress analysis of the NGCC drum was based on the European norm [27], design rules for steam boilers [28], and available literature in this area [7]. Dimensions of the drum main body and drum branches as well as the material properties of the drum are shown in Table 2.
The tangential stress defined in Equations (19)–(21) is used to represent the limiting stress as it is the highest stress component. These equations were added in the high-fidelity model.
σ θ , b o d y = σ θ , b o d y P + σ θ , b o d y T
σ θ , b o d y P = ( p p 0 ) ( 2 r i n m s + s m s ) 2 s m s  
σ θ , b o d y T = α m E m 1 ν ( 1 α w Φ w s m s 2 δ T δ t     )
In Equation (19), σ θ , b o d y is the overall tangential stress on the body of the drum, σ θ , b o d y P is the tangential stress due to pressure further defined in Equation (20), and σ θ , b o d y T is the tangential stress due to temperature also further defined in Equation (21). In Equation (20), p is the operating pressure (MPa) and p0 is the outside pressure (MPa), while T t is the temperature transient (K·h−1) in Equation (21).
To define the threshold for the stress constraint on the body of the drum, the maximum allowable stress (MAS) is defined in Equation (22). Specifically, the material SA-515 grade 70 requires the tensile strength to be within the range of 485–620 MPa and the minimum yield strength of 260 MPa [29].
M A S = t e n s i l e   s t r e n g t h m i n 4
The calculated MAS value is about 120   M P a , which is more conservative than the threshold considered in a similar work [7]. The MAS value is used to determine the maximum tangential stress value during dispatch according to Equations (23) and (24):
P c = 2 M A S 1 + 2 r i n m s s m s E
T h r e s h o l d = ( P c p 0 ) ( 2 r i n m s + s m s ) 2 s m s   .
In Equation (23), E is the joint efficiency of the shell–branch connection, with E = 0.95.
Then, the constraint on the overall stress that the drum body endures is given by Equation (25).
σ θ , b o d y T h r e s h o l d
Finally, to avoid overheating of the superheater and reheater, temperature constraints given by Equations (26) and (28) are considered:
T s u p e r h e a t e r T m a x
T r e h e a t T m a x
T m a x = 635   .
Additionally, two post-optimization analyses are carried out: one considering the stress of the most stressed part of the drum and the other on the CO2 equivalent assessment. The post-optimization analysis of the most stressed part of the drum is performed by using Equations (29)–(36) along with Equation (19) after the optimal dispatch is determined. Equations (29) and (30) are developed by modifying Equations (21) and (22) using the stress concentration factors due to the shell–branch connections (e.g., drum–downcomer junction). A diagram of the steam drum is shown in Figure 6.
σ θ , m o s t   s t r e s s e d P = λ m ( p p 0 ) ( 2 r i n m s + s m s ) 2 s m s  
σ θ , m o s t   s t r e s s e d T = λ t α m E m 1 ν ( 1 α w Φ w s m s 2 δ T δ t )
In Equations (29) and (30) λ m is the stress concentration factor due to pressure and is calculated using Equations (31)–(34), while λ t is the stress concentration factor due to temperature and is calculated using Equations (35) and (36).
λ m = 2.2 + e A ζ B  
A = 1.14 ( s b s m s ) 2 0.89 ( s b s m s ) + 1.43
B = 0.326 ( s b s m s ) 2 0.59 ( s b s m s ) + 1.08
ζ = d m b d m s d m s 2 s m s
λ t = { [ 2 h + 2700 h + 1700 z + h h + 1700 ( exp 7 z 1 ) ] 2 + 0.81 z 2 } 1 / 2  
z = d m b d m s  
In Equations (34) and (36), dmb and dms are the branch and main body mean diameters.
Regarding the environmental CO2 equivalent analysis, complete combustion of the fuel (natural gas and coal) is assumed for the calculation of the CO2 emissions. The coal considered is a bituminous Illinois no. 6 with Table 3 showing its composition [30]. Table 4 shows the composition used for the natural gas [30].
Considering the statement above, the optimization problem is solved for different renewable power generations and power grid components.

3. Results and Discussion

Overall, five case studies are considered. Initially, four case studies are assessed considering the dispatch of NGCC, SCPC, and NaS batteries to accommodate different levels of renewable penetrations into the grid, while supplying the demand using the optimization problem statement in Section 2. Figure 7 and Figure 8 show the different power demands and renewable penetration levels considered in these four case studies. Then, for the fifth case study, a more challenging renewable penetration scenario is elaborated, and an in-depth stress analysis is performed.
A technology-invariant curtailment penalty of 110 US$/MWh based on available literature [31] was considered in all case studies. The main purpose of this penalty is to maintain the curtailment at minimum levels, i.e., the wind/solar supply is only curtailed when the lack of curtailment renders the dispatch infeasible. The battery lifecycle is assumed to be relatively shorter than the lifetime of the power plants, as it depends on the degradation and the number of cycles. In this study, a penalty term is considered in the objective to take into account each cycle completed by the NaS batteries with a cost of 300 US$/kWh and a maximum number of 2500 cycles. A complete charge/discharge of the battery is considered as a cycle.
Figure 9, Figure 10, Figure 11, Figure 12 and Figure 13 show the results as stacked plots for the four case studies, i.e., the plot is shown as cumulative in the order of appearance of each technology (SCPC, NGCC, NGCC 2, NaS batteries, and renewables). For instance, if at some time instance a second line superimposes the first line, the energy system of the second line is not providing energy at that instance. The amount of energy being provided by a specific energy system is the subtraction between the previous and current energy system. This type of plot is chosen to rapidly verify if the power of the grid is balanced. If the line that represents the last energy source to be added matches the contour of the gray shaded area (demand), that means the power is balanced.
Case 1 considers the year 2019 for electricity price, demand, and solar/wind penetration into the grid. The power grid in this case is composed of SCPC, NGCC, and NaS batteries. The main reason for simulating Case 1 is to check the optimizer’s ability to dispatch the energy systems considering the power demand without curtailment of renewables. Results depicted in Figure 9 show that for Case 1, the load demand is supplied by the SCPC and NGCC with no curtailment of the renewable energy sources while considering potential damage for the NGCC power plant. The NaS batteries are not utilized at any point. The results obtained are in agreement with the usual strategy employed nowadays, as there is no large-scale deployment of NaS batteries. Essentially, the results indicated that at the 2019 renewables penetration levels, the already installed SCPC and NGCC are sufficient to supply the energy within the acceptable planned damage maximum without requiring curtailment or introduction of NaS batteries.
Case 2 considers a minor modification to the 2019 demand by increasing the minimum and maximum loads, i.e., the valleys and peaks, which would require faster ramping abilities from the grid. The electricity price and solar/wind penetrations into the grid are considered to be the same as in Case 1. The rationale in Case 2 is to verify the grid’s capability to supply the demand under a more challenging demand peak/valley scenario. In this case, the results depicted in Figure 10 show that the load demand is supplied mainly by the SCPC and NGCC, while the NaS batteries are used only twice (as highlighted in Figure 10) to smooth the peak at the approximate time of 450 h and to smooth the valley at the approximate time of 580 h. Specifically, at approximately 580 h, the sum of energy provided by the SCPC, NGCC, and renewables (green dashed line) surpasses the total load demand (shaded area). The reason is that there is an excess of electricity in the grid, which is used to recharge the batteries (orange line matches with shaded area). Due to the battery operating costs, they are used sparingly to prevent infeasibility. The detected behavior indicates that in case an unexpected demand was to happen, that would lead to suboptimal dispatch, and the batteries could be used more often for frequency regulation to satisfy the demand. As a perfect demand forecast of 6 min is assumed by the dispatcher, it is considered that any sudden deviation that was to happen could be supplied by spinning reserves. The results obtained in Case 2 are in accordance with the literature, in which the batteries are often employed in the absence of transmission structure or isolated areas where energy sources are scarce [13].
In Case 3, the electricity price, grid components, and demand are held constant and equal to those in Case 1, while the solar/wind penetrations are switched from 2019 to 2050 values. The 2050 wind/solar penetrations are shown in Figure 8. The 2050 intermittent renewable energy penetration is based on the 2019 penetration and Table 1 proportionality factors. The main goal of simulating Case 3 is to verify the potential of the 2019 grid components to supply the demand considering the high penetration of wind/solar sources. The results for Case 3 depicted in Figure 11 show that the load demand is supplied by the SCPC and NGCC components. However, to avoid infeasibility during dispatch, a 45% maximum curtailment is allowed. Figure 12 shows the curves for the accepted and curtailed renewable levels in Case 3. Under these circumstances, the batteries are not enough to avoid curtailment. To lower levels of curtailment, the optimizer stopped mostly after utilizing the energy available in the battery. This result indicates that considering the 2050 renewables penetration projection, the current grid is not sufficiently flexible to prevent curtailment. The reasons may be because the power plants are not capable to ramp up/down fast enough, the minimum load of the power plant is too high, and the batteries do not have enough capacity available to manage the variability introduced by the renewables.
Finally, Case 4 considers the same forcing functions as in Case 3; however, the grid is revised by retiring the SCPC and including another 641 MWe NGCC. The goal of this case study is to assess how the change in grid components may impact renewable curtailment. The results depicted in Figure 13 show that the demand load is supplied by the two NGCC power plants: one operated as a baseload power plant and the other operated as a load-following power plant, while achieving 0% renewable curtailment levels for 2050. This result indicates that for the grid to absorb higher variable renewable integration, there is a need to substantially increase grid flexibility. In Case 4, the simulated scenario is to replace the SCPC by an NGCC, but other alternatives would include retrofitting the SCPC, increasing NaS batteries capacity—albeit not being necessarily cost effective, or introducing an alternative energy storage system with higher capacity than the NaS batteries for a lower cost, such as pumped hydro. Table 5 shows a summary of the four analyzed cases regarding the year of demand, electricity price, wind/solar production, grid components, and a brief description of findings.
Overall, it is observed that the SCPC is dispatched as a baseload power plant while the NGCC is dispatched as load-following power plant when both technologies are available. This behavior is expected as the SCPC is a less flexible technology than the NGCC. The average calculation time to determine the optimal dispatch for Case 1 to 4 was approximately 1200 s for each entire profile. The computational time considering the optimal dispatch problem is promising for extending the overall problem into a more complex unit commitment problem, which would include managing energy systems reserves for coping with uncertainties regarding demand and renewable levels.
As mentioned above, after the optimal dispatch is determined, CO2 equivalent and concentrated stress analyses are then performed for each case study. The calculated CO2 equivalent is for the 30-day timespan and is depicted in Figure 14.
In Figure 14, note that most of the CO2 emissions derive from the SCPC power plant, as the NGCC has a lower rate of CO2 emissions per generated unit of power. In particular, even though Case 3 has the 2050 renewable power generation, the amount of CO2 emissions is similar to that in Cases 1 and 2, as renewables curtailment was allowed. This result implies that increasing the renewable production should be done along with increasing grid flexibility; otherwise, the environmental performance may not improve as expected. A possibility to lower the emissions would be to combine the SCPC with a carbon capture technology. The carbon capture unit would have an energy penalty for the SCPC [20,22], but it could substantially decrease the CO2 emissions. In addition, a multi-objective optimization considering economic and environmental objectives could also be applied to target emission reduction while determining the optimal dispatch.
Figure 15 shows both stress results on the body and at the most stressed part of the NGCC drum for all case studies. Specifically, as Case 4 has two NGCCs in operation, stress analysis is performed to both power plants.
In Figure 15, there are two distinct regions, namely, “Body” (highlighted in pink in Figure 15) and “Concentrated” (highlighted in orange in Figure 15). The “Body” section refers to the overall stress that the body of the drum is subjected to, and the maximum value of this stress is considered as a constraint in the optimal dispatch algorithm, as mentioned in Section 2. Although close to the threshold value, the optimal dispatch was determined without the body stress constraint ever becoming active. The “Concentrated” section in Figure 15 refers to the post-optimization stress analysis in the most stressed part of the drum. The stress in this specific part remained higher than the body stress threshold as well as the yield strength (260 MPa), while it was lower than the tensile strength. Thus, local uniform plastic deformation is expected to occur during a number of initial cycles in the most stressed part of the NGCC drum. However, the stress does not exceed the elastic shakedown limit, which is defined as twice the yield strength (520 MPa), and fully elastic behavior is anticipated in subsequent cycles.
Based on Case 4 results, Case 5 was elaborated under an even more challenging scenario regarding renewable penetration that demands up to 30%/min ramping rates from the power plants and shorter time steps (36 sec). Therefore, a shorter overall time horizon load and renewable penetration curves are created and shown in Figure 16 and Figure 17. The grid components are the same as Case 4, i.e., the two NGCCs and the NaS batteries.
The goal of this new case study is to give an insight about the NGCC drum thickness influence on the optimization objective function as well as the stresses on the NGCC’s drums under high ramping rates scenario. Therefore, a sensitivity study is conducted initially considering a range of [0.1461 m, 0.1905 m] for both drums’ thicknesses (from NGCC and NGCC 2), while holding all other forcing functions constant. The results of the sensitivity study are shown in Figure 18. The results indicate that for cases in which the cycling leads to prohibitive stresses, the optimization objective showed some sensitivity to the drum thickness. It is expected for the objective to become more sensitive as the number of NGCC’s considered in the grid increases. It was also observed that the thicknesses are not the same for baseload and load-following power plants to achieve the strict minimum value in the sensitivity analysis. At drum thicknesses higher than 0.1905 m, the dispatch was not completed, as prohibitive stresses reduced the feasible space to null at certain time instances.
After the sensitivity analysis was complete, a point for the normalized objective equal to zero is chosen for the optimal dispatch simulation, as shown in Figure 16. The drum thickness for the NGCC is 0.17145 m and for the NGCC 2, it is equal to 0.1461 m. The results of the optimal dispatch and corresponding stresses are shown in Figure 19 and Figure 20. Figure 19 shows that under a highly variable renewables scenario that demands faster ramping rates, the NaS batteries are used more often, discharging at approximately t = 9.6 min and t = 24 min. These results indicate that under more aggressive ramp rates (> 30%/min), the use of batteries could have a larger participation in the grid. Curtailment is still avoided under this scenario, while NGCC acts as baseload power plant, ramping down when there is excess of renewable energy, and NGCC 2 acts as load following power plant, cycling its load more often than the other NGCC. Figure 20 shows that the drum body stress constraint was active at times for NGCC 2, while for the other NGCC, the stress value was close but below the threshold. In addition, as the threshold is dependent on drum thickness, the thresholds for NGCC and NGCC 2 are slightly different.
Further analysis is performed considering the individual components of the body stress depicted in Figure 21 and Figure 22. For the NGCC, the mechanical stress, which is tensile in nature, was mostly at its design point, with the thermal stress frequently being compressive and relieving the overall body stress. This behavior agrees with the NGCC dispatch, which is employed often as baseload with minimal cycling. For the NGCC 2, the mechanical stress was often lower than its design point, which allowed a larger margin for thermal stress to vary between tensile and compressive stresses. This behavior agrees with the NGCC 2 dispatch, which is employed often as load-following and with more frequent cycling.

4. Conclusions

In this work, a dynamic mixed-integer linear programming optimization approach was proposed and successfully implemented to integrate different dispatchable energy systems with variable renewable energies to supply the demand. The dynamic optimization with an economic objective was used to determine the optimal dispatch of SCPC, NGCC, and NaS batteries under different levels of renewable penetrations. Overall, five case studies were presented. Four case studies were performed with different grid components with the integration of 2019 and 2050 variable renewable penetrations, with an additional case study in which challenging demand and renewable penetrations were considered.
The NGCC was successfully scheduled considering the damage model, and grid flexibility was mostly provided by the NGCC while either SCPC (Cases 1, 2, and 3) or a second NGCC (Cases 4 and 5) performed as baseload power plant. The NaS batteries were only utilized when the ramping rates were extreme (Cases 2 and 5), which is often avoided due to their high installation and operating costs. In addition, it was observed that a change in the grid composition is required to avoid curtailment in high-variable renewable penetration scenarios. Two of the case studies (Cases 4 and 5) explored a potential structure in which the SCPC was retired and replaced by NGCC. The results indicated that the proposed change prevented renewables curtailment. Other alternatives would be to retrofit the SCPC to accommodate cycling and/or add an energy system with larger capacity at lower cost than the NaS batteries. Additional environmental and equipment health analyses were performed on the results of the case studies. The environmental analysis corroborates that an increase of the renewable generation with no regard to grid flexibility is not sufficient to achieve desired environmental performance. The stress analysis indicated that under challenging ramping rates, the equipment health constraint was often active, reinforcing its importance for cycling applications. Additional sensitivity studies regarding drum thickness showed that different designs may have to be considered for power plants with different purposes, e.g., baseload and load-following.
Overall, the proposed framework shows that the grid operators could be benefitted in terms of guiding the decision-making process of dispatching energy systems. However, one of the limitations of the framework is that synergy between fatigue and creep is not considered, as well as potential corrosion due to tension. In addition, due to the use of linear models for linear programming, inaccuracies when compared to the high-fidelity nonlinear models were introduced. For simplicity, only the tangential direction of stress was considered. Future work would include a more complex damage model as well as converting the problem statement to nonlinear programming. The inclusion of a fatigue-based damage model during a unit commitment problem employing dynamic models is also of interest, as it would consider the complete cycle damage during the scheduling of power plants. Finally, the framework could be extended to include more energy systems and account for start-up/shutdown costs and time constraints for scheduling purposes.

Supplementary Materials

The following are available online at https://0-www-mdpi-com.brum.beds.ac.uk/2227-7390/8/6/868/s1, “Power_plants_model_parameters” and “Power plants model parameters–Help Guide”.

Author Contributions

Conceptualization, R.K., Y.W., S.P.V., D.B., F.V.L., and R.T.; methodology, R.K., Y.W., and S.P.V.; formal analysis, R.K., Y.W., and S.P.V.; investigation, R.K., Y.W., and S.P.V.; resources, D.B., F.V.L., and R.T.; data curation, R.K., Y.W., and S.P.V.; writing—original draft preparation, R.K.; writing—review and editing, R.K., Y.W., S.P.V., D.B., F.V.L., and R.T.; visualization, R.K.; supervision, D.B., F.V.L., and R.T.; project administration, D.B., F.V.L., and R.T.; funding acquisition, D.B., F.V.L., and R.T. All authors have read and agreed to the published version of the manuscript.

Funding

The authors gratefully acknowledge the financial support from U.S.-China Clean Energy Research Center (CERC) program and U.S. Department of Energy under Cooperative Agreement DEPI0000017.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. REN21. Renewable Energy Policy Network for 21st Century. Available online: http://www.ren21.net/status-of-renewables/global-status-report/ (accessed on 16 April 2018).
  2. U.S. Energy Information Administration. Annual Energy Outlook 2019 with Projections to 2050. Available online: https://www.eia.gov/outlooks/aeo/pdf/aeo2019.pdf (accessed on 12 December 2019).
  3. Greening the Grid. Demand Response and Storage. Available online: http://www.greeningthegrid.org/integration-in-depth/demand-response-and-storage (accessed on 7 April 2019).
  4. Bentek Energy LLC. How Less Became More… Wind, Power and Unintended Consequences in the Colorado Energy Market. Prepared for Independent Petroleum Association of Mountain States. 16 April 2010. Available online: https://docs.wind-watch.org/BENTEK-How-Less-Became-More.pdf (accessed on 22 May 2020).
  5. U.S. Energy Information Administration. More U.S. Coal-Fired Power Plants Are Decommissioning as Retirements Continue. Available online: https://www.eia.gov/todayinenergy/detail.php?id=40212 (accessed on 12 December 2019).
  6. Tuttle, J.F.; Powell, K.M. Analysis of thermal generator’s participation in the Western Energy Imbalance Market and the resulting effects on overall performance and emissions. Electr. J. 2019, 32, 38–46. [Google Scholar] [CrossRef]
  7. Kim, T.S.; Lee, D.K.; Ro, S.T. Analysis of thermal stress evolution in the steam drum during start-up of a heat recovery steam generator. Appl. Therm. Eng. 2010, 20, 977–992. [Google Scholar] [CrossRef]
  8. Fontaine, I.P.; Galopin, I.J. HRSG optimisation for cycling duty based on Euro Norm EN 12952-3. OMMI 2008, 5, 1–14. [Google Scholar]
  9. Denholm, P.; Margolis, R. Energy Storage Requirements for Achieving 50% Solar Photovoltaic Energy Penetration in California; NREL Technical Report, NREL/TP-6A20-66595; 2016. Available online: https://www.nrel.gov/docs/fy16osti/66595.pdf (accessed on 22 May 2020).
  10. Sufyan, M.; Rahim, N.A.; Tan, C.; Muhammad, M.A.; Raihan, S.R.S. Optimal Sizing and Energy Scheduling of Isolated Microgrid Considering the Battery Lifetime Degradation. PLoS ONE 2019, 14. [Google Scholar] [CrossRef] [PubMed]
  11. Adams II, T.A.; Nease, J.; Tucker, D.; Barton, P.I. Energy Conversion with Solid Oxide Fuel Cell Systems: A Review of Concepts and Outlooks for the Short- and Long-Term. Ind. Eng. Chem. Res. 2013, 52, 3089–3111. [Google Scholar] [CrossRef]
  12. Schaefer, S.; Vudata, S.P.; Bhattacharyya, D.; Turton, R. Transient modeling and simulation of a nonisothermal sodium-sulfur cell. J. Power Sources 2020, 453. [Google Scholar] [CrossRef]
  13. International Renewable Energy Agency. Battery Storage for Renewables: Market Status and Technology Outlook. January 2015. Available online: https://www.irena.org/documentdownloads/publications/irena_battery_storage_report_2015.pdf (accessed on 12 December 2019).
  14. Vudata, S.P.; Bhattacharyya, D.; Turton, R. Optimal Thermal Management of a High-Temperature Sodium Sulphur Battery. In Proceedings of the 2018 AIChE Annual Meeting, Pittsburgh, PA, USA, 28 October 2018; p. 49e. [Google Scholar]
  15. Denholm, P.; Hand, M. Grid flexibility and storage required to achieve very high penetration of variable renewable electricity. Energy Policy 2011, 39, 1817–1830. [Google Scholar] [CrossRef]
  16. Trifkovic, M.; Marvin, W.A.; Daoutidis, P.; Sheikhzadeh, M. Dynamic Real-Time Optimization and Control of a Hybrid Energy System. AIChE J. 2014, 60, 2546–2556. [Google Scholar] [CrossRef]
  17. Wang, Y.; Bhattacharyya, D.; Turton, R. Dynamic Modeling and Control of a Natural Gas Combined Cycle Power Plant for Load-following Operation. Comput. Aided Chem. 2019, 47, 101–106. [Google Scholar]
  18. Wang, Y.; Bhattacharyya, D.; Turton, R. Evaluation of Novel Configurations of Natural Gas Combined Cycle (NGCC) Power Plants for Load-Following Operation using Dynamic Modeling and Optimization. Energy Fuels 2020, 34, 1053–1070. [Google Scholar] [CrossRef]
  19. Kim, R.; Lima, F.V. A Tchebycheff-based multi-objective combined with a PSO-SQP dynamic real-time optimization framework for cycling energy systems. Chem. Eng. Res. Des. 2020, 156, 180–194. [Google Scholar] [CrossRef]
  20. He, X.; Wang, Y.; Bhattacharyya, D.; Lima, F.V.; Turton, R. Dynamic modeling and advanced control of post-combustion CO2 capture plants. Chem. Eng. Res. Des. 2018, 131, 430–439. [Google Scholar] [CrossRef]
  21. Corengia, M.; Torres, A.I. Effect of Tariff Policy and Battery Degradation on Optimal Energy Storage. Processes 2018, 6, 204. [Google Scholar] [CrossRef] [Green Version]
  22. Zhang, Q.; Turton, R.; Bhattacharyya, D. Development of Model and Model-Predictive Control of an MEA-Based Postcombustion CO₂ Capture Process. Ind. Eng. Chem. Res. 2016, 55, 1292–1308. [Google Scholar] [CrossRef]
  23. He, X.; Lima, F.V. Development and implementation of advanced control strategies for power plant cycling with carbon capture. Comput. Chem. Eng. 2019, 121, 497–509. [Google Scholar] [CrossRef]
  24. He, X.; Lima, F.V. A modified SQP-based model predictive control algorithm: Application to supercritical coal-fired power plant cycling. Ind. Eng. Chem. Res 2020. submitted for publication. [Google Scholar]
  25. PJM. Markets & Operations. Available online: https://www.pjm.com/markets-and-operations.aspx (accessed on 16 December 2019).
  26. McKubre, M.C.H.; Tanzella, F.L.; Smedley, S.I. The Electromotive Force of the Na/S Cell. J. Electrochem. Soc. 1989, 136, 303–305. [Google Scholar] [CrossRef]
  27. European Committee for Standardization. EN 13345 Part 3, Unfired Pressure Vessels, Clause 17; Simplified Assessment of Fatigue Life, and Clause 18; Detailed Assessment of Fatigue Life; 2002; Available online: https://www.unm.fr/medias/files/107_FICHIER_0.pdf (accessed on 23 May 2020).
  28. Technical Rules for Steam Boilers (TRD 301, Annex 1, Design). Calculation for Cyclic Loading Due to Pulsating Internal Pressure or Combined Changes of Internal Pressure and Temperature. Available online: https://kupdf.net/download/trd-301-annex-1-design_58c5ef24dc0d609909339028_pdf (accessed on 23 May 2020).
  29. Henan Steel Guang International Trade Co, LTD. Available online: http://www.steelplatecuttingpart.com/products/asme-sa515-boiler-or-pressure-vessel-steel-plate.html (accessed on 4 January 2020).
  30. National Energy Technology Laboratory. Cost and Performance Baseline for Fossil Energy Plants Volume 1: Bituminous Coal and Natural Gas to Electricity, NETL-PUB-22638. 2019. Available online: https://netl.doe.gov/projects/files/CostAndPerformanceBaselineForFossilEnergyPlantsVol1BitumCoalAndNGtoElectBBRRev4-1_092419.pdf. (accessed on 8 April 2020).
  31. Schill, W.; Pahle, M.; Gambardella, C. Start-up costs of thermal power plants in markets with increasing shares of variable renewable generation. Supplementary Information. Nat. Energy 2017, 2. [Google Scholar] [CrossRef] [Green Version]
Figure 1. Schematic of a single sodium sulfur cell that composes the battery.
Figure 1. Schematic of a single sodium sulfur cell that composes the battery.
Mathematics 08 00868 g001
Figure 2. Schematic of the natural gas-fired combined cycle power plant (NGCC).
Figure 2. Schematic of the natural gas-fired combined cycle power plant (NGCC).
Mathematics 08 00868 g002
Figure 3. Schematic of the supercritical pulverized coal-fired power plant (SCPC).
Figure 3. Schematic of the supercritical pulverized coal-fired power plant (SCPC).
Mathematics 08 00868 g003
Figure 4. Plot of capacity (Wh) vs. state of discharge of a single NaS cell.
Figure 4. Plot of capacity (Wh) vs. state of discharge of a single NaS cell.
Mathematics 08 00868 g004
Figure 5. Plot of state of discharge vs. time for a single NaS cell.
Figure 5. Plot of state of discharge vs. time for a single NaS cell.
Mathematics 08 00868 g005
Figure 6. Diagram of a steam drum (including the drum–downcomer junction).
Figure 6. Diagram of a steam drum (including the drum–downcomer junction).
Mathematics 08 00868 g006
Figure 7. 2019 total and modified electricity demand load.
Figure 7. 2019 total and modified electricity demand load.
Mathematics 08 00868 g007
Figure 8. 2019 and 2050 variable renewable energy penetrations (solar/wind).
Figure 8. 2019 and 2050 variable renewable energy penetrations (solar/wind).
Mathematics 08 00868 g008
Figure 9. Power dispatch optimization result for Case 1.
Figure 9. Power dispatch optimization result for Case 1.
Mathematics 08 00868 g009
Figure 10. Power dispatch optimization result for Case 2.
Figure 10. Power dispatch optimization result for Case 2.
Mathematics 08 00868 g010
Figure 11. Power dispatch optimization result for Case 3.
Figure 11. Power dispatch optimization result for Case 3.
Mathematics 08 00868 g011
Figure 12. 2050 renewable levels accepted and curtailed for Case 3.
Figure 12. 2050 renewable levels accepted and curtailed for Case 3.
Mathematics 08 00868 g012
Figure 13. Power dispatch optimization result for Case 4.
Figure 13. Power dispatch optimization result for Case 4.
Mathematics 08 00868 g013
Figure 14. CO2 equivalent emissions by case study.
Figure 14. CO2 equivalent emissions by case study.
Mathematics 08 00868 g014
Figure 15. Drum body and concentrated tangential stresses during dispatch: Cases 1 to 4.
Figure 15. Drum body and concentrated tangential stresses during dispatch: Cases 1 to 4.
Mathematics 08 00868 g015
Figure 16. Modified electricity demand load for Case 5.
Figure 16. Modified electricity demand load for Case 5.
Mathematics 08 00868 g016
Figure 17. Modified variable renewable energy penetrations (solar/wind) for Case 5.
Figure 17. Modified variable renewable energy penetrations (solar/wind) for Case 5.
Mathematics 08 00868 g017
Figure 18. Sensitivity analysis of optimization objective due to NGCC and NGCC 2 drums’ thicknesses.
Figure 18. Sensitivity analysis of optimization objective due to NGCC and NGCC 2 drums’ thicknesses.
Mathematics 08 00868 g018
Figure 19. Power dispatch optimization result for Case 5.
Figure 19. Power dispatch optimization result for Case 5.
Mathematics 08 00868 g019
Figure 20. Drum body and concentrated tangential stresses during dispatch for Case 5.
Figure 20. Drum body and concentrated tangential stresses during dispatch for Case 5.
Mathematics 08 00868 g020
Figure 21. Mechanical and thermal stresses for NGCC: Case 5.
Figure 21. Mechanical and thermal stresses for NGCC: Case 5.
Mathematics 08 00868 g021
Figure 22. Mechanical and thermal stresses for NGCC 2: Case 5.
Figure 22. Mechanical and thermal stresses for NGCC 2: Case 5.
Mathematics 08 00868 g022
Table 1. Electricity generation from selected fuels [1].
Table 1. Electricity generation from selected fuels [1].
Fuel20182050
Natural gas34%39%
Coal28%17%
Renewables18%31%
        Wind        6.66%        7.77%
        Solar        2.34%        14.88%
          Other 1        9.00%        8.35%
Nuclear19%12%
Total99%99%
1 Other renewable sources, e.g., hydro and geothermal.
Table 2. High-pressure drum specification with MOC 1 being SA-515 grade 70.
Table 2. High-pressure drum specification with MOC 1 being SA-515 grade 70.
ParameterSymbolValueUnit
Branch radiusrinb0.20m
Main body radiusrinms0.8382m
Branch thicknesssb0.04m
Main body thicknesssms0.1524m
Density ρ 7753kg·m−3
Differential thermal expansion coefficient α m 1.4 × 10−5K−1
Modulus of elasticityEm207 × 103MPa
Poisson’s ratio ν 0.292-
Shape factor Φ w −0.358-
Thermal diffusivity α w 0.0342m2·h−1
Thermal transmittance for steamh1000Wm−2K−1
1 MOC: Material of Construction.
Table 3. Bituminous coal composition.
Table 3. Bituminous coal composition.
Rank
Seam
Bituminous
Illinois No. 6
ComponentWeight (%)
Moisture11.12
Carbon63.75
Hydrogen4.50
Nitrogen1.25
Chlorine0.15
Sulfur2.51
Ash9.70
Oxygen7.02
Table 4. Natural gas composition.
Table 4. Natural gas composition.
ComponentVolumetric (%)
Methane93.1
Ethane3.2
Propane0.7
n-Butane0.4
Carbon Dioxide1.0
Nitrogen1.6
Methanethiol(A)0.00000575
Table 5. Case studies according to their main differences and findings.
Table 5. Case studies according to their main differences and findings.
Case NumberDemandElectricity Price Wind/Solar Penetration Grid Components Findings
1201920192019SCPC, NGCC, NaS BatteriesDispatch the NGCC and SCPC. NaS batteries not dispatched. No renewables curtailment.
22019 120192019SCPC, NGCC, NaS BatteriesDispatch the NGCC and SCPC. NaS batteries dispatched to smooth peaks and valleys. No renewables curtailment.
3201920192050 2SCPC, NGCC, NaS BatteriesDispatch the NGCC and SCPC. NaS batteries not dispatched. Up until 45% renewables curtailment allowed.
42019201920502 NGCCs, NaS BatteriesDispatch the 2 NGCC’s. NaS batteries not dispatched. No renewables curtailment.
1 Modified; 2 Curtailed.

Share and Cite

MDPI and ACS Style

Kim, R.; Wang, Y.; Vudata, S.P.; Bhattacharyya, D.; Lima, F.V.; Turton, R. Dynamic Optimal Dispatch of Energy Systems with Intermittent Renewables and Damage Model. Mathematics 2020, 8, 868. https://0-doi-org.brum.beds.ac.uk/10.3390/math8060868

AMA Style

Kim R, Wang Y, Vudata SP, Bhattacharyya D, Lima FV, Turton R. Dynamic Optimal Dispatch of Energy Systems with Intermittent Renewables and Damage Model. Mathematics. 2020; 8(6):868. https://0-doi-org.brum.beds.ac.uk/10.3390/math8060868

Chicago/Turabian Style

Kim, Rebecca, Yifan Wang, Sai Pushpitha Vudata, Debangsu Bhattacharyya, Fernando V. Lima, and Richard Turton. 2020. "Dynamic Optimal Dispatch of Energy Systems with Intermittent Renewables and Damage Model" Mathematics 8, no. 6: 868. https://0-doi-org.brum.beds.ac.uk/10.3390/math8060868

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop