The control‐oriented model of coordinated control system based on stochastic differential equations

Considering the large‐scale development of renewable energy resources globally, control optimization of coal‐fired power plants is becoming increasingly crucial. The dynamics of the coordinated control system (CCS) must be studied prior to designing a controller. These dynamics of the subcritical units have been expressed through various modeling methods. However, in previous studies, ordinary differential equations (ODEs) have been primarily employed, which cannot reflect the uncertainties in the system. Therefore, a model with a greater accuracy should be observed to comprise uncorrelated residuals. In this study, the uncertainties in the calorific value of fired coal and combustion process were analyzed first. A normal distribution of disturbance was assumed in this process. The dynamics of the CCS were then described with stochastic differential equations (SDEs). Furthermore, a parameter estimation procedure was designed. The residual evaluation is employed to improve the model evaluation. After simulations, the SDE‐based model‐3 in this study elucidates the dynamics of the subcritical boiler‐turbine system better than the ODE‐based model and other SDE‐based models. The measured values can be regarded as a possible result of this model, rendering it a potential platform to employ stochastic model predictive control in CCS.


|
DU et al.
significant proportion of power supply in China, are required to change the role positions and face more challenges in daily operations. For system operations, coal-fired power plants are expected to have a higher adjustment rate, greater accuracy of adjustment load, and a wider adjustment range. These requirements can be achieved through automatic generation control (AGC) by system operators. These targets need hardware and software renovations and updates, among which the optimization of the CCS update is significant. 4,5 Prior to the control strategies, the dynamics of the coordinated control system need to be analyzed; subsequently, the model should be established. The modeling of CCS has received significant attentions in the previous decades. There are three types of common modeling methods: white-box, black-box, and gray-box. 6 White-box modeling, also referred to as mechanism modeling, relies on fundamental physical laws to master the dynamic behavior of thermodynamic variables. Herein, mass, energy, and momentum balances are primarily used. The parameters in the model can be calculated by design parameters according to their own physical meaning. Black-box modeling, also referred to as data-based modeling, requires the use of experimental or operating data sets, and always consists of either a form of a support vector machine, neural network, fuzzy structure, etc The effectiveness of this model is significantly influenced by the training data, which implies that the model may be invalid if another set of data are adopted as input. Gray-box modeling is a combination of white-box modeling and black-box modeling. The main structure with unknown parameters is built by the fundamental laws. Regarding the unknown parameters, optimization algorithms are employed to utilize the information in the operating data. [7][8][9] The drum boiler models 10,11 were derived from the first principles after several experiments; the unknown parameters were calculated by the design parameters. The models are apparently built through white-box modeling and are mainly composed of differential equations. The 4th order model has been widely selected as a simulation platform for advanced control algorithms. 12,13 In black-box modeling, the coordinated control systems were modeled through fuzzy-modeling or neural modeling 14,15 to prepare for optimized controllers. Developed by an advanced process simulation software, APROS, the dynamic simulation models 16 can simulate the dynamic behavior of the real power plant with a high accuracy, even at very low loads. 17 These simulated models have the advantages of higher accuracy but their complexity disables them from being used for CCS directly. Other simulation models of supercritical once-through units 18,19 have also been developed by other commercial software in previous studies. The aforementioned gray-box modeling has simultaneously made use of prior knowledge of the working mechanism and information in the historical operation data 20 ; therefore, gray-box modeling has received significant attentions. In previous studies, the control-oriented models of subcritical units [21][22][23][24] and once-through units 6,[25][26][27] can be found. These are all conducted through gray-box modeling. Considering the preliminary investigation above, it was found that the control-oriented models of the coordinated control systems in previous studies are mainly established based on ODEs, and a relatively small part has the form of partial differential equations (PDEs).
The models based on ODEs or PDEs indicate that the coordinated control system was treated as a deterministic system. 28 The same inputs would lead to unique and same outputs. However, besides measurement errors, there are still several undetermined factors working within the coordinated control system operation. 21,27 These factors would lead to the undetermined outputs with the same inputs. In previous studies, the variable residual variations always remain between the measured values and the ODE-based models. Sometimes, the residuals would become evident. If the residual variations can be further separated into diffusion and measurement noise, a more accurate model for the system can be obtained. In model evaluation, if the model is detailed enough to describe the dynamics of the system, the residuals should be uncorrelated. 29 The correlation tests of residuals would assist in evaluating the model accuracy. These types of model evaluations are common in other research areas; however, they are inadequate for the establishment of a model with coordinated control system. 21,22,[24][25][26][27][28][29][30] This gap in literature may mislead a more accurate model 30,31 and it should be employed for improvements.
The intuitive idea of SDEs is to add a noise with some probability distribution into the ODEs of a dynamic system. Each time, the output is influenced by the output at the previous time and the probability distribution at the present time. The use of the stochastic approach can help to distinguish the system dynamic information existing in the observations from measurement noise. SDEs have been used in the forecasting of power demand, 32 wind speed, 33 wind power, 34 vehicle flow, 35 etc They are also used for improving graybox models in energy systems 36,37 for obtaining a detailed stochastic gray-box model. These models can simultaneously capture the main dynamics and the randomness hidden in the systems. To the best of our knowledge, few studies have previously focused on SDE-based model for CCS. And the correlation tests of residuals in model evaluation should also be employed for improvements in CCS. Based on SDE-based models, a stochastic model predictive control can be employed [38][39][40] and the performance of control systems would be improved. Hence, the objective of this study is to assess the potential of an SDE-based model for a coordinated control system.
The remaining sections of this paper are organized as follows. Section 1 analyzes the uncertainties in the modeling procedure of CCS and the three SDE-based models of the subcritical unit have been derived. Section 2 introduces the | 131 DU et al.
procedure and algorithms regarding parameter estimation. Section 3 compares the effects between different types of models following model evaluations. Finally, the conclusions are provided in Section 4.

| MODEL DESCRIPTION
In this study, a 300 MW subcritical power plant unit in the North of China is selected as the object of study. This unit consists of typical components such as a single-reheat turbine, tangential firing furnace, MPS-type medium speed coal mills, three superheaters, and two reheaters. In the past decade, this unit has been operating for peak shaving and frequency modulation. With more integrated renewable energy, a greater number of coal-fired power plants would be positioned as this type of unit by system operators in the future. Therefore, the optimized control in this unit has a certain reference value to be carried out. The control-oriented nonlinear model of CCS in a subcritical power plant has been previously derived by Tian. 24 The model is verified through the history data in the time domain, and the controllers 24,41 have been implemented for real power plants based on this model. The simplified model is displayed as follows: where C 0 = inertia time of pulverizing system (s); r B = mass flow rate of coal into the furnace (kg/s); r m = coal flow rate command (kg/s); = the time delay of pulverizing system (s); C b = inertia time of boiler section (s); p b = steam pressure of drum (MPa); p = main steam pressure after superheaters (MPa); Q = the heat from the coal to working medium (kJ); K t = inertia time of turbine section (s); N e = the unit load (MW); u t = throttle valve opening fraction (0 ~ 1); k 3 , k 2 , b = the unknown parameters.
Due to the lack of direct measuring sensors, the throttle valve opening fraction is calculated as follows: where p m represents the pressure in the governing stage of the turbine.
The model structure is mainly derived from the laws of mass and energy balance; Equation (1) represents the dynamics of the coal mill. The dynamics of the boiler and turbine section are expressed in Equations (2) and (3), respectively. Furthermore, Equation (4) indicates the relationship between the drum pressure and the main steam pressure. This is fitted by data in several steady working conditions. The coal flow rate command and valve opening fraction are the inputs of model, while the unit load and main steam pressure are the outputs. The detail derivation process is presented in the previous study. 24 In this model, the coupled relationship between the unit load and the main steam pressure is revealed by simplified expressions. The model is apparently composed of three ODEs and a nonlinear algebraic equation. For example, the undetermined parameters are estimated through estimation process proposed in the previous study 24 and a model is established to represent CCS for this unit. The residuals between the measured values and the outputs of the model are shown in Figure 1. A considerable bias can be observed in the results. An autocorrelation function (ACF) is employed to evaluate the deviations and the results are shown in Figure 2. A significant correlation can be found in both, unit load and main steam pressure, which indicates that the system could have a more accurate model to capture the measured values. To model the unknown factors in the system and separate the residuals, new models would be derived based on the SDEs in the following.
(1) Considering Equation (2), the symbol Q represents the heat absorbed by the working medium from the coal flow rate, as fired. Note, the coal as fired would release the heat through the burning in the furnace. The heat is absorbed by the water wall. Then, the heat is transferred to the work medium from the metal. The major processes contained can be regarded as the combustion process and the heat transfer process. They are considerably simplified in the modeling. Due to several disturbances inside, the value of Q needs to be analyzed further.
Nowadays, because of the tension in the China coal market, the coal used in plants may be different in a period of time; this phenomenon would influence the calorific value of coal as fired. First, for a long time, there may be a significant difference in the mean calorific value of the coal as received. The length of time should be a week, at a minimum, in these situations. The air coal ratio and water coal ratio needs to be updated with a mean calorific value of coal. Second, the co-combustion of different types of coal would increase the randomness in the calorific value of coal as fired; the randomness of coal as fired naturally increases with an actuator operation error of coal mills. In a power plant, the calorific value of coal as fired would be tested every eight hours, as shown in Figure 3. The calorific value of coal as fired would fluctuate within a certain range. The maximum deviation is observed as approximately 10% of the calorific values form the figure. Hence, the uncertainty in calorific value of coal into the furnace should be considered. It is assumed that the coal quality does not change dramatically over a period of time. And the uncertainty in this study is assumed to be expressed by the following equation with a Gaussian disturbance: where 1,t ∼  (0, 1) and 21 represent unknown parameter. The combustion process in furnace can be concluded as following the stoichiometric chemical reactions: A proper air coal ratio is needed to provide a reasonable value of O 2 for combustion according to coal quantity. The ACFs of the residuals in ODE-based model air volume is adjusted by readjusting the angle position of the damper. The combustion status is fairly influenced by the oxygen content. In dynamic operating, the oxygen volume can be expressed by adding a Gaussian disturbance. 42 Due to the active participation of power plants in peak shaving and the current primary frequency regulation, the fluctuation of coal quality and air volume would become more volatile. Here, in order to represent the factors caused by the combustion status and damper control in dynamic operation, the heat produced in the furnace is derived as follows: where Q com represents the heat produced through combustion. 22 , 1 represents unknown parameters and 2,t ∼  (0, 1) .
The heat transfers from the furnace to the working medium, including the economizer, water wall, superheaters, and the reheaters, can be regarded as a proportion to Q com as follows: where 2 represents an unknown parameter. Combined with Equations (6)-(8), the following equation can be derived: where In summary, significant stochasticity exists in the boiler section as described above.
where w 1 and w 2 represent the Wiener processes.
However, considering the model of the turbine section, there are no significant features as the coal quality and combustion process. The inertia time of turbine section is significantly smaller than the boiler section, and the disturbances in the turbine section are not particularly apparent. Hence, the original equation, Equation (3), is adopted and the entire model is denoted as the SDE-based model-1. Considering the comprehensiveness of the problem, the model is supplemented with a diffusion term to represent the stochasticity in the turbine section, Equation (11), denoted as the SDE-based model-2.
Due to the significantly smaller section of the inertia time of turbine, this section can be modeled as a simple algebraic expression, which has been applied in once-through units in previous studies. 25,26 In this study, the coordinated control system with Equation (12) is denoted as the SDE-based model-3.
The observation equations can be expressed as follows: where e k1 and e k2 represent measurement error and This study presents 4 different models to represent the subcritical boiler-turbine system. The testing models are shown in Table 1 All the models in Table 1 can be expressed in the following form: where is the vector of parameters that should be estimated.
The Equation (14) represents a continuous state-space description of the coordinated control system with a stochastic term. The Equation (15) represents the discrete values of the outputs y of time t k . This type of model is a common SDEbased model and has been comprehensively studied, 31 but is rarely seen in a coordinated control system. In this study, the models in Table 1 will be applied to modeling the real plant, after which the statistical tests would be employed to evaluate the models and assist in choosing the more proper one.

| PARAMETERS ESTIMATION
After the modeling in the previous section, there are several unknown parameters to be determined. The procedure of parameter estimation is designed with the consideration of the model structure as follows: Step 1: Time delay , which reflects the time delay between the coal flow rate command and the coal as fired, can be regarded as empirical parameters to reduce the complexity of model. 21,27 Step 2: k 1 , k 2 , k 3 , and b are regarded as constant parameters and calculated through regression and statistical analysis of several steady working conditions, 6,26 due to E(d (t)) = 0.
Step 3: The rest parameters of the ODE-based model are estimated by a genetic algorithm (GA).
Step 4: This step is separated into two stages. In the first stage, the output is p , donated as L( ;p ). In the second stage, the output is N e , donated as L( ;N e ). The diagram can be seen in Figure 4.
The maximum likelihood estimation and extended Kalman filter 43,44 are employed to approximate the parameters here.
In this study, when the maximum deviations of the main steam pressure and real time load do not exceed 8% and 2% for 2 hours, respectively, this state is selected as the steady working condition. Based on the criteria, the measured values of the units in the steady working conditions from the historical database are shown in Table 2 in Step 2. The static parameters can then be estimated. Considering the rest parameters in the ODE-based model, the parameter identifications based on optimization are employed as Step 3. In addition, due to the disadvantages of local convergence, the optimizations would be conducted a few times to get a better result.
The objective function, for parameter estimation Step 3, is designed as follows: where p ⋀ and N ⋀ e are the outputs of the ODE-based model. p and N e represent the measured values.
In Step 4, the maximum likelihood estimation should be claimed first. A sequence of measurements y 0 , ..., y N sampled from DCS are notated as  N . The joint probability density of  N can be expressed by

ODE-based model
where p y 0 | represents the initial probability density. The Gaussian density of y k is assumed and the its mean and covariance can be notated by The joint probability density can be expressed by The form of the negative logarithm of the joint probability density is chosen and the unknown parameters are estimated by solving the following nonlinear optimization problem: The undetermined parameters, k and R k|k−1 , can be calculated recursively by the extended Kalman filter (EXF). In the Kalman filter, the essential part is the Kalman gain equation, combined with output prediction equations, state prediction equations, and updating equations; the EXF considers the nonlinearity in the model. The detailed information regarding EXF can be referred to in previous studies. 43 Thus, the estimation in Step 4 can be carried out.

DISCUSSION
The simulations and evaluations are conducted in MATlab and RStudio. Following parameter estimation, the unknown parameters in the models are shown in Table 3. The parameters of the boiler section are kept the same in each SDE-based model, but the parameters in the turbine section are different for each model.
The use of the likelihood function for parameter estimation is based on the assumption of the Gaussian distributed residuals, which should be uncorrected in model evaluation. Hence, the residuals should be tested for normality in order to verify the assumptions. In this part, the raw periodogram (RP), ACF, and cumulated periodogram (CP) are used to evaluate the residuals simultaneously.
RMSE is designed to be the numerical evaluation index of the models. RMSE is defined as follows: where n is the length of data. y The residuals evaluating the results of the main steam pressure of the SDE-based models are plotted in Figure 5. The spectrum is apparently the same during varying frequency. The CP plots indicate that the residuals are random, because the autocorrelations are within the 95% confidence bands. The values of ACF are within 95% confidence bands ( Figure 5C). Consequently, the residuals obtained can be regarded as white noise. The RMSE(p ) is equal to 0.1286, which is close to the results of the ODE-based model in Table 4. The results of the first stage are satisfied with the assumptions, and the existence of stochastic errors is verified inside the boiler section. Due to the same model of the boiler section, this result is adopted in all SDE-based models.
Regarding the unit load of the SDE-based model-1, the evaluation of residuals can be seen in Figure 6. The spectrum cannot be considered a constant. The CP plot and ACF plot show that the residuals are clearly not white noise because most of the values are outside of the 95% confidence bands. Hence, although the RMSE(N e ) is equal to 3.1614, the SDE-based model-1 cannot be considered useful to simulate the system. In the SDE-based model-2, RMSE(N e ) is 2.8533, which is a significant improvement compared to the ODE-based model. However, the evaluation of residuals shown in Figure 7 cannot meet the demand. Apparent trends can be found in RP. The majority of the CP and ACF results are not within the 95% confidence bands, which also cannot be considered because the statistical test results recommend discarding this model.
There is no inertial item adopted in the SDE-based model-3. However, RMSE(N e ) is 2.7559, which is the most considerable result of the unit load. In addition, the same intension of the spectrum can be found in the RP plot in Figure 8A. The values of CP and ACF are close to the 95% confidence bands, and most of them are within the 95% confidence band. These results are a significant improvement compared to the SDEbased model-1 and the SDE-based model-2. Consequently, the residuals obtained can be treated as white noise.
Based on the analysis above, SDE-based model-3 and ODE-based model are simulated with comparison to the measured values in Figures 9 and 10. The SDE-based model-3 results are obtained from the Monte Carlo simulations. Apparent errors can be found between the measured values (blue lines) and the ODE-based model (green lines), both, in the main steam pressure and unit load; this is a significant feature of ODE-based models. The measured values are well covered by the SDE-based model (red lines). The coordinated control system operation in this moment can be considered as a possible determination of the SDE-based model, due to a certain degree of stochasticity.
Besides, step signals are added to inputs of the SDEbased model-3 to verify open-loop responses. Measured values are deficient because of the lack of practical tests. But the correctness of the conclusion is not affected. Figure 11 shows the responses of the model with step disturbance on coal flow rate command. After a period of time delay, the heat absorbed by steam gradually increases due to the inertia effect of the boiler. Then, the steam pressure rise up to the corresponding level as well. Therefore, the unit load reaches a higher level compared with the original level. Figure 12 shows the responses of the model with step disturbance on the throttle valve opening. The main steam pressure declines immediately with the abrupt change of the opening. The reduction of the main steam pressure leads to the release of the boiler's heat storage. The load of the unit is increased first. Then the heat storage of boiler is consumed, and the unit load will recover to the original level when the amount of fuel and water flow is not changed. The open-loop responses of SDE-based model-3 are in accordance with the actual situations.

| CONCLUSION
To optimize the performance of CCS, the control-oriented model of the coordinated control system needs to be further studied. In this study, the comparison between the ODEbased model established in the previous study and the measured values is first analyzed through an example. To make up the residuals, stochastic factors are considered in modeling. The randomness values in the calorific value of coal as fired in the combustion process are analyzed; these are expressed by adding Gaussian white noise. Subsequently, the models are established with the Wiener process. In order to achieve a more appropriate model, the test contains three types of SDE-based models for comparison with the ODEbased model. All the SDE-based models have the form of a continuous-discrete state-space model with stochastic differential equations. Due to the distinct characteristic of the model structure, the procedure of parameter estimation is  For model evaluation, the raw periodogram, autocorrelation function, and cumulated periodogram are adopted to evaluate the residuals. Results show that the residuals of the SDE-based model-3 can be regarded as white noise, which satisfies the assumption. In addition, RMSE is chosen to be the numerical evaluation of the models. Results show that the SDE-based model-3 has the smallest bias, 0.1281 and 2.7559. Subsequently, the forecast effects of the ODE-based model and the SDE-based model-3 have been compared. In