An ice storage air conditioning system economy and responsiveness double target optimization model predictive control method

Abstract

The present application relates to a kind of ice storage air conditioning system economy and responsiveness dual-objective optimization model predictive control method, including S1, based on the building measured data and prediction algorithm after preprocessing, load prediction model is established;S2, establishment equipment energy consumption model based on least square method and equipment dynamic response model based on transfer function;S3, the model predictive control algorithm model of economy and responsiveness dual optimization target is established.The present application not only maximizes the building peak power consumption, relieves the power grid peak power consumption pressure, reduces the user power cost by responding to the electricity price policy when ice storage air conditioning system control strategy is formulated;And when strategy implementation control, after the optimal cold distribution demand response strategy is formulated, higher control accuracy and more ideal economic benefits can be realized by tracking optimization to control target.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of building energy systems, and particularly relates to a method for predicting control of an ice storage air conditioning system economic and responsiveness double target optimization model. BACKGROUND

[0002] As an efficient and economical, flexible and controllable active energy storage device, the ice storage air conditioning system can effectively alleviate the peak pressure of the power grid and the intermittent renewable energy consumption problem under the optimization control of building demand response. Compared with building thermal mass (passive energy storage), the energy storage system (active energy storage) in the ice storage air conditioning system has greater capacity and better controllability in reducing peak power demand in response to the requirements of the smart grid, and does not need to change the indoor temperature, thus having less negative impact on indoor thermal comfort. It reduces the peak power demand of the building refrigeration system on the power grid by producing and storing cold during off-peak hours and then using the stored cold during peak hours.

[0003] According to different purposes of model simulation, energy system modeling can be divided into steady-state model and dynamic model. The former is used for system optimization design and to obtain a feasible optimization control or regulation strategy to improve system energy efficiency, economy, and provide grid services; the latter is used to investigate the rapidity, stability and accuracy of local equipment in response to regulation. For the development of optimization strategies of energy systems involving ice storage, due to the diversity of objects participating in regulation and the different characteristics, most studies are based on steady-state device models, without considering the dynamic response performance during actual control. However, the ice storage system involves a complex phase change heat transfer process, and when responding to the upper regulation signal, the system undergoes a dynamic process from one steady state to another steady state. The steady-state mathematical model cannot reflect the dynamic response characteristics when the control signal changes, which may make the system unable to respond to dynamic regulation instructions in real time and accurately, thus leading to the actual effect of the optimization strategy deviating from the expected benefit.

[0004] The ultimate engineering significance of the optimization control strategy is to truly achieve how much reduction of the peak power demand of the power grid. When developing the optimization strategy, the system operation effect is usually analyzed based on ideal working conditions, however, the device layer control is the decomposition and implementation of the optimization strategy, and is the process of executing the regulation instructions. The actual benefit that can be achieved by the control process is greatly related to the performance of the controller, which affects the accuracy and rapid following ability of the system control response.

[0005] Traditional proportional-integral-derivative (PID) control has been widely used in HVAC applications for decades, which is simple in structure, widely used, and the parameters are usually adjusted by trial and error. One of the main shortcomings of traditional control technology is that it cannot compensate for thermal lag and changing set points, nor can it adapt to changing control signals. For complex processes with nonlinearity, time-varying, and multi-parameter coupling, it is difficult to adjust the parameters. Model predictive control technology has been applied to building HVAC systems to overcome some of the weaknesses of traditional controllers. The predictive control algorithm can provide thermal lag compensation, thereby greatly improving the stability of the closed loop, and its predictive characteristics can enable control actions to start earlier, thereby closely tracking changing set points.

[0006] However, in the current large number of studies on the development of optimal control strategies for ice storage systems, few studies have focused on and optimized the influence of the time lag of energy storage equipment on the dynamic response characteristics of the system during control, which is crucial for the system to accurately implement control according to the operation strategy and achieve the expected benefits of the optimization strategy. SUMMARY

[0007] The purpose of the present application is to overcome the shortcomings of the prior art and provide an ice storage air conditioning system economic and responsive dual-objective optimization model predictive control method.

[0008] The present application solves its technical problems by the following technical solutions:

[0009] An ice storage air conditioning system economic and responsive dual-objective optimization model predictive control method, characterized in that the steps of the method are:

[0010] S1, establishing a load prediction model based on preprocessed building measured data and prediction algorithms

[0011] By analyzing the correlation of the preprocessed building cooling load related data, determine reasonable and effective parameters as the input variables of the prediction model; then use the long short-term memory network model to train the measured data and perform accuracy evaluation to obtain a reliable building cooling load prediction model;

[0012] S2, establishing a device energy consumption model based on the least squares method and a device dynamic response model based on the transfer function

[0013] Based on device parameters and operation data, for the main energy-using equipment in the building ice storage air conditioning system, a mathematical model of its device energy consumption is established using the least squares method; for the energy storage equipment in the ice storage air conditioning system, a dynamic response model is established using the transfer function method;

[0014] S3, establishing an economic and responsive dual-optimization target model predictive control algorithm model

[0015] In the established cost function, the economic optimization objective is to minimize the total operation cost of the system in a day, which is related to the energy consumption of equipment and the electricity price in different time periods; the responsive optimization objective is to minimize the square deviation between the predicted value and the set value of the control output temperature, which ensures that the control temperature reaches the set value as quickly as possible, and then the rolling optimization and feedback correction of the controller model are performed based on the particle swarm algorithm to obtain better economic and control response effects.

[0016] Further, the step S1 is specifically:

[0017] (1) Correlation analysis

[0018] The Pearson correlation coefficient is used to analyze the correlation of the measured multiple variables, and the calculation formula is as follows:

[0019]

[0020] Wherein: n is the total number of sample data; i indicates the i th sample;

[0021] And are the sample average values of two types of parameters;

[0022] The value range of r is-1-1, a negative value indicates that the two variables are negatively correlated, a positive value indicates that the two variables are positively correlated, and when |r|<0.2, it is considered that the two variables are extremely weakly correlated or not correlated;

[0023] (2) Long short-term memory network model training

[0024] The LSTM model is used for load prediction, and the key processing layer of the LSTM model mainly includes three gates: forgetting gate, input gate and output gate, and a cell state is further included. At a certain time, the processing process of the LSTM algorithm is generally divided into three steps:

[0025] First, determine the information to be passed and discarded by the LSTM memory cell, and use the Sigmoid function in the forgetting gate to determine it, and the calculation process is as follows:

[0026] f t =σ(W f ·[h t-1 ,x t ]+b f ) (2)

[0027] Wherein: W f is the forgetting gate weight matrix;

[0028] b f is the forgetting gate bias vector;

[0029] h t-1 is the output at the previous moment;

[0030] x t is the current time input;

[0031] σ is the Sigmoid activation function, the value is between 0 and 1, which represents how much information can pass through, 1 means all information can pass through, 0 means all information cannot pass through;

[0032] Second, determine the information that needs to be added to the memory cell of LSTM, the calculation process is as follows:

[0033] i t = σ (W i · [h t-1 , x t ] + b i ) (3)

[0034]

[0035] Where: is the input information of the newly generated hidden layer of the tanh layer at the current time;

[0036] W i , W c are the input gate and memory unit weight matrix respectively;

[0037] b i , b c are the input gate and memory unit bias vector respectively;

[0038] The Sigmoid and tanh activation functions are respectively:

[0039]

[0040]

[0041] Third, determine the output information of the LSTM output gate. First, the Sigmoid layer gets the initial output, and then uses the tanh function to make C t ∈ [-1, 1], which is multiplied by the initial input to obtain the output of the model, the calculation process is as follows:

[0042] o t = σ (W o · [h t-1 , x t ] + b o ) (7)

[0043] h t = o t · tanh (C t ) (8)

[0044] where: W o ,b o are the weight matrix and bias vector of the output gate, respectively;

[0045] The RMSprop method is used as the optimization algorithm of the LSTM prediction model, and the algorithm design is as follows:

[0046]

[0047]

[0048] where: E[g 2 ] is the squared mean of the gradient at the tth iteration;

[0049] α is the momentum, set to 0.9;

[0050] W t is the model parameter at the tth iteration;

[0051] g y = ΔJ(W t ) is the gradient of the cost function of the tth iteration model with respect to the parameter W;

[0052] η0 is the global initial learning rate;

[0053] ε is a value that tends to 0 infinitely;

[0054] (3) Accuracy evaluation

[0055] To measure the accuracy of the model prediction results, the coefficient of variation of root mean square error (CV-RMSE) and the normalized mean bias error (NMBE) are selected as evaluation indexes, and their calculation formulas are as follows:

[0056]

[0057]

[0058] where: n is the number of prediction samples;

[0059] y i is the predicted value of the ith sample;

[0060] is the actual value of the ith sample,

[0061] According to ASHRAE14-2014 standard, when using hourly data to predict and verify, the results of CV-RMSE below 30% can be considered as accurate, and NMBE below 10% can be considered as accurate.

[0062] Moreover, the step S2 is specifically:

[0063] (1) Establishing the mathematical model of energy consumption of equipment

[0064] Chiller and water pump are the main energy-consuming equipment in ice storage air conditioning system, and the energy consumption of the system affects the amount of electricity taken from the power grid by the building and the operation cost of the user. The multiple polynomial model is the relationship between COP and three variables of chilled return liquid temperature and cooling water return water temperature, and the calculation accuracy is relatively high. The multiple polynomial is selected to establish the mathematical model of chiller, and the calculation formula is as follows:

[0065]

[0066] Wherein: b1~b 10 is the parameter to be identified;

[0067] Q e is the unit load, kW;

[0068] T ei is the chilled return liquid temperature, ℃;

[0069] T ci is the cooling water return water temperature, ℃;

[0070] Because the frequency data of the water pump is difficult to obtain in the actual test process, the water pump model is often simplified as a single variable parameter function about flow, therefore, the mathematical expression of energy consumption of the variable frequency water pump is as follows:

[0071] P pump =a1+a2·M+a3·M 2 (14)

[0072] (2) Establishing the dynamic response model of equipment

[0073] 1) Modeling of heat transfer differential equation

[0074] The reasons for the state change of the secondary refrigerant in the coil include heat transfer in the flow direction and the pipe diameter direction. The mass balance formula of the secondary refrigerant along the flow direction of the pipeline (axial direction) is as follows:

[0075] G hf,in =G hf,out =G hf (15)

[0076] Wherein: G hf,in is the flow of the secondary refrigerant at the inlet of the ice storage tank, kg / s;

[0077] G hf,out is the outlet flow rate of the ice storage tank, kg / s;

[0078] The energy balance equation of the secondary refrigerant side:

[0079]

[0080] wherein: c hf is the specific heat capacity of the secondary refrigerant, kJ / kg·℃;

[0081] p hf is the density of the secondary refrigerant, kg / m 3 ;

[0082] A p is the cross-sectional area of the coil, m 2 ;

[0083] L p is the length of the coil, m;

[0084] N is the number of nodes set at the end of the coil, i.e. the number of intermediate state variables set from the input to the output;

[0085] T hf,i and T hf,i+1 are the inlet and outlet temperatures of the i-th unit of the ice storage coil, respectively, ℃;

[0086] U is the heat transfer coefficient, W / (m 2 ·℃);

[0087] A l is the inner surface area of the coil tube wall, m 2 ;

[0088] T pcm is the temperature of the phase change material outside the coil, ℃;

[0089] The enthalpy method model represents the state change of the phase change medium in the form of enthalpy H in the energy balance equation, which is applicable to the entire phase change stage. According to the change of the phase change medium, the heat transfer process is divided into three stages: liquid phase stage, solid-liquid paste stage, and solid phase stage,

[0090]

[0091] wherein: H is the specific enthalpy, kJ / kg;

[0092] T m is the phase change temperature, 0℃;

[0093] ΔH m is the latent heat of phase change per unit mass of the substance, 334 kJ / kg;

[0094] c pcm Cp is the specific heat capacity of the phase change medium, kJ / (kg·℃);

[0095] The energy balance equation of the water outside the pipe is:

[0096]

[0097] Wherein: p pcm is the density of the phase change medium, kg / m 3 ;

[0098] V pcm is the volume of the phase change medium (ice-water paste area), m 3 ;

[0099] λ pcm is the thermal conductivity of the phase change medium, W / (m·k);

[0100] A pcm is the cross-sectional area of the phase change medium, m 2 ;

[0101] UA l is the product of the comprehensive heat transfer coefficient between the refrigerant inside the coil and the phase change medium and the heat transfer area, and its value is equal to the inverse of the total thermal resistance of the radial heat transfer of the pipe;

[0102]

[0103]

[0104]

[0105]

[0106] Wherein: R hf is the convective heat transfer thermal resistance between the refrigerant and the pipe wall when flowing inside the ice storage coil, ℃ / W;

[0107] R pipe is the thermal resistance of the pipe wall, ℃ / W;

[0108] R pcm is the thermal resistance of the phase change medium outside the coil, ℃ / W;

[0109] r i , r o , r pcm are the inner radius, outer radius of the coil and the radius of the phase change medium, respectively, m;

[0110] The convective heat transfer coefficient h hf between the refrigerant and the pipe wall is related to the Reynolds number Re, and is calculated by the following formula:

[0111]

[0112] where: p, v, d, μ are the density (kg / m 3 ), flow rate (m / s), diameter (m) and viscosity coefficient (m 2 ) of the coolant, respectively;

[0113] The Prandtl criterion number Pr is calculated by the formula:

[0114]

[0115] When Re > 2300, the coolant in the coil is in a turbulent state, and the Nusselt criterion number Nu is calculated according to the Dittus-Boelter formula:

[0116] Nu = 0.023Re 0.8 Pr 0.4 (25)

[0117] When Re ≤ 2300, the coolant in the pipe flows in a laminar state,

[0118]

[0119]

[0120] In the ice melting and cold releasing link, when the ice layer and the coil form a concentric water ring after a period of ice melting, in addition to the thermal resistance, the convective heat transfer resistance also needs to be considered, but since the calculation process is too complex, the concept of effective thermal conductivity λ ac related to the Rayleigh number R eff is established to simplify the calculation, and the empirical formula of the effective thermal conductivity λ

[0121]

[0122]

[0123] where: Ra is the Rayleigh number corresponding to the average temperature of the inlet and outlet of the coolant, with δ = r w -r o as the characteristic length;

[0124] 2) Transfer function conversion

[0125] For complex physical system, it is difficult to convert to transfer function directly, which may involve multiple input and output variables. State space equation is often used as a bridge between differential equation and transfer function. A set of first-order differential equations with state variables are used as mathematical model to describe the process of "input-state-output" by state space method in modern control theory. The state space equation can be easily converted to transfer function in s domain.

[0126] In order to get the expression of coefficient matrix of linear state space equation, the steady-state solution of corresponding steady-state differential equation is needed to be solved first. The steady-state equation is the algebraic equation corresponding to the case when the derivative term is 0. Each variable in the equation is expressed in the form of steady-state value and small increment, i.e.

[0127]

[0128] Where, symbol - represents the steady-state value of each variable, and symbol δ represents the small increment of each variable. The overall heat transfer coefficient is a dependent variable, which changes with the change of coolant flow rate. The overall heat transfer coefficient U related to G hf is linearized by using first-order Taylor expansion formula:

[0129]

[0130] The linearized derivation result of the heat transfer differential equation can be obtained by bringing the linear form of each variable into the heat transfer differential equation:

[0131] Liquid phase (water)

[0132] Energy balance equation of coolant side:

[0133]

[0134] Energy balance equation of water outside the tube:

[0135]

[0136] Solid-liquid paste phase (ice-water mixture)

[0137] Energy balance equation of coolant side:

[0138]

[0139] Energy balance equation of paste zone outside the tube:

[0140]

[0141] Solid phase (ice)

[0142] Energy balance equation of coolant side:

[0143]

[0144] The energy balance equation of the ice outside the pipe is:

[0145]

[0146] The state variable, input variable and output variable of the state space equation are selected, wherein the state variable can reflect the change of the internal state of the system, T hf,in , G hf , H in are input parameters, T hf,i and H i with differential terms are state variables in the heat transfer process, the outlet temperature T hf,out of the cold carrier in the ice storage tank coil is an output variable, and the above two heat transfer differential equations can be converted into the following linear state space model:

[0147]

[0148] T hf,out = CT (39)

[0149] The Laplace transform is performed on the state space equation of the system to obtain the input and output model of the ice storage coil as follows:

[0150] T hf,out = G1(s)G hf +G2(s)T hf,in +G3(s)H in (40)

[0151] Wherein: G1(s) is a transfer function of the cold carrier flow rate with respect to the outlet temperature of the cold carrier;

[0152] G2(s) is a transfer function of the inlet temperature of the cold carrier with respect to the outlet temperature of the cold carrier;

[0153] G3(s) is a transfer function of the initial enthalpy of the phase change medium with respect to the outlet temperature of the cold carrier;

[0154] 3) Time delay time calculation

[0155] In the process of participating in the energy system regulation and control of the ice storage system, the total time delay time can be defined as the time period between the time when the control parameter of the local energy storage device is adjusted according to the scheduling instruction and the time when the control output parameter at the end changes, which is the coupling effect of the two time delay factors of thermal inertia and flow delay, and the two time delay factors can be represented by the inertia link and the pure delay link in the transfer function.

[0156] (a) Thermal inertia time delay

[0157] The thermal inertia is mainly derived from the heat transfer medium, which is defined as the heat transfer capacity of the medium due to its own heat capacity. When the system exchanges heat with the outside world, the temperature of the medium rises or falls, and it takes a certain time to change, not immediately. This feature is called thermal inertia. In the transfer function, the inertia element (also known as viscous damping) is described as when the input quantity changes suddenly, the output quantity cannot change immediately, but only gradually changes according to the exponential law, which reflects the inertia of the controlled object. Any actual system with inertia element characteristics can store energy. To characterize thermal inertia, first establish the energy balance relationship of the heat transfer process, and then perform Laplace transform on the partial differential equation to obtain the transfer function, so that the inertia time constant τ of the medium in the transfer function can be calculated c :

[0158]

[0159]

[0160] (b) Pure time delay

[0161] The pure lag phenomenon of the system is mainly caused by the medium transmission. In the transfer function, the delay element is also called pure lag element, which indicates that in the steady state operation of the system, the input and output quantities are stable, and the steady state change of the input quantity needs to be delayed for a period of time before being transmitted to the output quantity, which can be expressed as:

[0162] y(t) = f(t-τ p ) (43)

[0163]

[0164] Where y and u are the output and input signals, respectively;

[0165] τ p is the pure time delay;

[0166] In the case of constant system flow rate, the pure delay time is almost fixed; while in the variable flow system, the pure delay time will change with the change of flow rate.

[0167] Moreover, the step S3 is specifically:

[0168] (1) Double target cost function establishment

[0169] The application focuses on the economy of the ice storage system according to the demand response of the time-of-use electricity price and the dynamic characteristics of the system when responding to the control signal, uses a prediction model and an optimization algorithm to obtain the economically optimal combination of the power contribution of the chiller and the ice storage tank in the prediction period, and on this basis, the cold quantity provided by the ice storage system is transmitted to the end user through the plate heat exchanger to meet the cooling load demand of the end user, and the purpose of on-demand cooling is achieved through the adjustment of the refrigerant flow; the controlled parameter is the secondary chilled water supply temperature, i.e. the set value to be tracked, in the cost function, the first term is the total operation cost of the system in a day, which is related to the equipment energy consumption and the electricity price in different time periods, which reflects the control target of saving the operation cost of the user and also reflects the flexibility of demand response; the second term is the square of the deviation between the predicted value of the control output temperature and the set value, which ensures that the control temperature reaches the set value as soon as possible, so that the system accuracy and followability are better, and thus the optimization objective function (also called cost function) can be established:

[0170]

[0171] Wherein: E ch and E pump are the energy consumptions of a single chiller and a water pump, respectively, kW;

[0172] N1 and N2 are the total number of chillers and water pumps, respectively;

[0173] is the set value of the chilled water supply temperature, ℃;

[0174] (2) Optimization algorithm

[0175] The particle swarm optimization algorithm (PSO) is used to solve the optimization problem in the above model predictive control, there are N particles in the solution space constantly moving and searching for the optimal solution, this process is realized through iteration, in the population, each particle is determined by the speed and position of the two parameters to determine the state of each movement, and the particle can record the search movement process, can record the corresponding optimal position of each particle in the search process, and then the position with the best comprehensive result in all particles can be obtained, in the iteration process, the mathematical expressions of the speed v ij and the position x ij of the i-th particle in the j-th dimension are as follows:

[0176] v ij (t+1)=ωv ij (t)+c1r 1j (t)(p ij -x ij (t))+c2r 2j (t)(p gjx ij (t)) (46)

[0177] x ij (t+1) = x ij (t) + v ij (t+1) (47)

[0178] Wherein: t is iteration times;

[0179] J is the number of variables;

[0180] ω is inertia weight;

[0181] c1, c2 are acceleration coefficients;

[0182] r 1j ,r 2j It is a random number uniformly distributed in [0, 1];

[0183] p i It is the optimal position of particle individual;p g It is the optimal position in all particles;

[0184] (3) Feedback correction

[0185] Based on the time delay time calculation method, a time delay compensation module is added, the delay time is calculated according to the control increment, that is, the state of the refrigerant temperature and flow at each control node, and the signal is transmitted to the control sequence signal, so that the time delay time can be eliminated in the next response, and the real-time performance and robustness of the control system are ensured through the rolling prediction and feedback correction link; as a kind of advanced control method, the model predictive controller can predict the control effect of the expected output signal of the input control in the future time period through the obtained system scheduling signal, improve the response speed of the system through advanced control and response optimization, and reduce the influence of time lag.

[0186] The advantages and beneficial effects of the application are:

[0187] The ice storage air conditioning system economic and responsive double target optimization model predictive control method of the application not only reduces the building peak power consumption to the maximum extent, relieves the power grid peak power consumption pressure and reduces the user power consumption cost through the response to the electricity price policy when the ice storage air conditioning system control strategy is formulated, and when the strategy is implemented and controlled, the optimal cold quantity distribution demand response strategy is formulated, higher control precision and more ideal economic benefits are realized by tracking and optimizing the control target. BRIEF DESCRIPTION OF DRAWINGS

[0188] Figure 1 It is a flow chart of the application;

[0189] Figure 2A graph of the hourly cooling load correlation analysis result of the present application;

[0190] Figure 3 A comparison graph of the predicted value of the load prediction model obtained by training of the present application and the verification set;

[0191] Figure 4 A graph of the matching situation of the hourly cooling supply and the building load demand under different controls of the present application;

[0192] Figure 5 A graph of the optimal cooling capacity distribution of the economic optimization of the combined cooling of the present application;

[0193] Figure 6 A comparison graph of the actual reduction of the demand response electric load under MPC and PID controls of the present application;

[0194] Figure 7 A graph of the hourly flexibility and cumulative flexibility of the electric load of the present application. DETAILED DESCRIPTION

[0195] The present application will be further described in detail through specific embodiments, and the following embodiments are only descriptive and not limiting, and cannot limit the protection scope of the present application.

[0196] As shown in Figure 1 , an ice storage air conditioning system economic and responsiveness dual target optimization model predictive control method, the innovation lies in that: the steps of the method are:

[0197] Step 1: Establish a load prediction model based on preprocessed building measured data and prediction algorithm

[0198] A certain comprehensive business district office building in Beijing is selected as a typical case, and the building uses ice storage air conditioning system cooling machine and ice storage tank combined cooling mode. The case test work is mainly carried out in the cooling season, and the test parameters include indoor dry bulb temperature, outdoor dry bulb temperature, wind speed, wind direction, air conditioning load side flow, air conditioning load side supply and return water temperature, and indoor equipment power consumption, etc. The test time interval is 1h. The data acquisition instrument equipment used in the test process of the case building and the precision are shown in Table 1.

[0199] Table 1 Test instrument equipment and precision of the case building

[0200]

[0201]

[0202] After initially collecting the measured data, data preprocessing is necessary to improve its quality. This mainly includes outlier removal and missing value imputation. Box plots are used to identify outliers; those significantly deviating from the normal data range are corrected or discarded as appropriate. To ensure the integrity of the model input data, missing data in the original data are imputed using a weighted average arithmetic method. Simultaneously, since the case study building is an office building, data from non-working days needs to be filtered out to prevent interference with the workload patterns on normal working days.

[0203] Step 2: Establish a load forecasting model based on preprocessed building measurement data and forecasting algorithms.

[0204] From the acquired cooling load data, Pearson correlation analysis was performed on the hourly cooling load with eight variables: outdoor dry-bulb temperature, water vapor pressure (hPa), wind speed, previous load value, indoor temperature, relative humidity, indoor lighting power, and horizontal visibility. The results are as follows: Figure 2 As shown, the correlations with hourly cooling load, ranked from highest to lowest, are: previous load value > lighting electricity consumption > outdoor dry-bulb temperature > water vapor pressure (hPa) > indoor dry-bulb temperature > 0.2 > outdoor relative humidity > horizontal visibility > maximum wind speed. Therefore, this study determines that the previous load value, lighting electricity consumption, outdoor dry-bulb temperature, water vapor pressure (hPa), and indoor dry-bulb temperature will be used as the five parameters for the cooling load prediction model.

[0205] Load prediction was performed using a long short-term memory network model. 980 valid measured data samples were selected, of which 75% were used as the training set and 25% were used as the validation set. Figure 3 The comparison between the measured and predicted values ​​on the validation set after model training is shown. The curves of the predicted and actual values ​​are in good agreement. The calculated error evaluation metrics, NMBE and CV-RMSE, are 7.3% and 11.9%, respectively, both meeting the accuracy evaluation requirements for hourly load values ​​in ASHRAE Guideline 14-2014, as shown in Table 2. Therefore, the model's prediction results are relatively accurate and can be used for further research in this project.

[0206] Table 2 Results of accuracy evaluation indicators for cooling load prediction models

[0207]

[0208] Step 3: Establish a device energy consumption model based on the least squares method and a device dynamic response model based on the transfer function.

[0209] With the relevant measured data and the established equipment model, the mathematical model of the main equipment of the building energy system is trained by using the least square method program to determine the to-be-fitted parameters in the corresponding expression. Among them, for the double-working-condition water chiller, the equipment characteristics of the ice storage working condition and the air conditioning working condition are different, and need to be fitted respectively. The to-be-fitted coefficients of the double-working-condition water chiller obtained by using the least square method multiple regression fitting are as shown in Table 3.

[0210] Table 3 to-be-fitted parameter values of the mathematical model of the double-working-condition water chiller

[0211]

[0212] There are three types of water pumps in the building of the case, which are the primary side ice melting refrigeration circulating pump, the primary side cooling water pump and the secondary side refrigerant pump. The least square method is used to train the established mathematical model of the water pump, and the mathematical expressions of various types of water pumps are as follows:

[0213] (1) Primary side ice melting refrigeration circulating pump

[0214] P pump,melting = 38.6317 + 0.1406·M - 0.00006·M 2 (48)

[0215] (2) Primary side cooling water pump

[0216] P pump,cooling = 45.0537 + 0.02617·M + 0.00002·M 2 (49)

[0217] (3) Secondary side refrigerant pump

[0218] P pump,chiller = 68.3123 - 0.0203·M + 0.00003·M 2 (50)

[0219] The ice storage tank in the building of the case is an internal ice melting coil type ice storage tank. According to its equipment physical parameters and operation data, the mathematical expression of the ice storage tank obtained by using the phase change heat transfer modeling method based on transfer function is as follows:

[0220]

[0221] Comparative analysis of the effect of the double-objective optimization model predictive control

[0222] To verify the effect of the proposed method on the demand response control of building air conditioning system, the PID control condition in the traditional strategy benefit single objective optimization MPC is taken as the comparison condition, which is referred to as the traditional PID control condition. The full cooling machine condition refers to the traditional non-optimized control strategy of building air conditioning system without participating in the demand response of power grid, which is taken as the benchmark condition to quantify the demand response flexibility of MPC condition after economic optimization combined with time-of-use electricity price, and then the actual control effect and demand response benefit of different control conditions are compared and analyzed.

[0223] a. Matching of cooling capacity and load demand

[0224] In order to more intuitively observe the load matching situation, Figure 4 The curve comparison of average hourly cooling capacity and hourly predicted cooling load demand is shown. The matching degree of cooling capacity and actual cooling load demand of the proposed double objective optimization MPC control condition is better, and the shadow overlap degree is better. The relative error of cooling capacity and actual cooling load demand in a day is only 2.1%, while the relative error of traditional PID control condition is 4.3%. The cumulative cooling load deviation of the proposed double objective optimization MPC control condition is 629.13kW·h, which is 43.4% lower than 1113.3kWh of the traditional PID control condition.

[0225] b. Comparison and analysis of demand response energy consumption and cost

[0226] The cooling capacity allocation results of MPC for economic optimization of ice storage tank and ice maker are shown in Figure 5 The optimal strategy of model predictive control is to use double-working-condition water chiller to store ice in ice storage tank during 0:00-4:00 period of valley electricity price, and all the cooling capacity of ice storage tank is allocated to release in high electricity price stage during the day, which makes the end user benefit and improves the grid friendliness of building energy system. After getting the hourly cooling capacity control instruction of the optimal strategy, the actual demand response control effect that can be achieved by different controllers is not the same.

[0227] The electricity load usage of the proposed double objective optimization MPC control condition and the traditional PID control condition under the same demand response optimization control strategy is shown in Figure 6The full-chiller cooling condition is taken as the benchmark condition to compare the effect of the electric load reduction under the chiller ice storage tank combined cooling condition. Compared with the full-chiller cooling condition, the traditional PID control condition and the proposed dual-objective optimization MPC control condition both effectively reduce the electric load in the high electricity price period within the two red boxes, but the reduction amount is significantly different. For the demand response optimization control strategy, the proposed dual-objective optimization MPC control condition has less electricity consumption than the traditional PID control condition. The largest difference between the two conditions occurs at 12 o'clock, at which time the traditional PID control condition reduces the electricity load by 126.91 kWh, and the proposed dual-objective optimization MPC control condition reduces the electricity load by 156.88 kWh, achieving a peak reduction effect of nearly 30 kWh more than the traditional PID control condition.

[0228] The electricity flexibility curve, cumulative energy consumption, and cost reduction effect under the two control conditions are shown in FIG. 6. Figure 7 Compared with the full-chiller cooling condition, the total energy consumption of the traditional PID control condition and the proposed dual-objective optimization MPC control condition under the demand response scenario increases by 20.0% and 19.2%, respectively, and the total electricity consumption of the two conditions is not significantly different. The proposed dual-objective optimization MPC control condition reduces the peak electricity consumption by 24.3%, reduces the total cost by 11.7% throughout the day, and reduces the cost in the peak period by 24.3%. The traditional PID control condition reduces the peak electricity consumption by 22.7%, reduces the total cost by 10.7% throughout the day, and reduces the cost in the peak period by 22.7%. Under the time-of-use electricity price mode, the demand response economic optimization strategy of the proposed dual-objective optimization MPC control condition increases the total energy consumption, but both controllers effectively achieve electric load transfer and reduction during the peak electricity consumption period. In addition, compared with the traditional PID control condition, the proposed dual-objective optimization MPC control condition saves 8.5% of the total electricity cost and shifts 6.5% of the peak load.

[0229] Although the embodiments of the present application and the drawings are disclosed for illustrative purposes, those skilled in the art can understand that various alternatives, changes and modifications are possible without departing from the spirit and scope of the present application and the appended claims, and therefore the scope of the present application is not limited to the disclosed embodiments and drawings.

Claims

1. An ice storage air conditioning system economy and responsiveness dual-objective optimization model predictive control method, characterized in that: The steps of the method are: S1, establishing a load prediction model based on pre-processed building measured data and a prediction algorithm Through correlation analysis on the pre-processed building cold load related data, reasonable and effective parameters are determined as input variables of the prediction model; then a long short-term memory network model is used to train the measured data, and accuracy evaluation is performed to obtain a reliable building cold load prediction model; S2, establishing a device energy consumption model based on the least square method and a device dynamic response model based on the transfer function Based on device parameters and operation data, for the main energy-using devices in the building ice storage air conditioning system, a mathematical model of device energy consumption is established by using the least square method; for the energy storage device in the ice storage air conditioning system, a dynamic response model is established by using the transfer function method; S3, establishing a model prediction control algorithm model with economic and responsive optimization objectives In the established cost function, the economic optimization objective refers to the minimum total operation cost of the system in a day, which is related to device energy consumption and electricity prices in different time periods; the responsive optimization objective refers to the minimum square deviation between the predicted value of the control output temperature and the set value, which ensures that the control temperature reaches the set value as soon as possible; and then based on the particle swarm algorithm, rolling optimization and feedback correction of the controller model are performed to obtain better economic and control response effects.

2. The ice storage air conditioning system economy and responsiveness dual objective optimization model predictive control method of claim 1, wherein: The step S1 is specifically: (1) Correlation analysis Pearson correlation coefficient is used to perform correlation analysis on the measured multiple variables, and the calculation formula is as follows: Wherein: n is the total number of sample data; i refers to the i th sample; and are the sample means of the two types of parameters; The value range of r is -1~1, negative value indicates negative correlation between two variables, positive value indicates positive correlation between two variables, and when |r|<0.2, it is considered that two variables are extremely weakly correlated or not correlated; (2) Long short-term memory network model training LSTM model is used for load prediction, and the key processing layer of LSTM mainly includes three gates: forgetting gate, input gate and output gate, and a cell state is further included. At a certain time, the processing process information of the LSTM algorithm is generally divided into three steps: First step, determine the information passed and discarded by the LSTM memory cell, use the Sigmoid function in the forgetting gate to determine, and the calculation process is as follows: f t = σ(W f · [h t-1 , x t ]+ b f ) (2) wherein: W f is a forget gate weight matrix; b f is a bias vector for the forget gate; h t-1 is the output of the previous time instant; x t is the current time input; σ is the Sigmoid activation function, the value is between 0 and 1, indicating how much information is passed, 1 indicates that all information can be passed, and 0 indicates that all information cannot be passed; Second step, determine the information that needs to be added to the memory cell of LSTM, and the calculation process is as follows: i t = σ(W i · [h t-1 , x t ]+b i ) (3) wherein: is the input information of the newly generated hidden layer of the tanh layer at the current time point; W i ,W c are input gates, memory cell weight matrices, respectively; b i ,b c are input gate, memory cell bias vector, respectively; The Sigmoid and tanh activation functions are respectively: Third step, determine the output information of LSTM output gate, first, Sigmoid layer gets initial output, then use tanh function to make C t ∈[-1, 1], multiply it with initial input to get the output of the model, and the calculation process is as follows: o t = σ(W o · [h t-1 , x t ]+b o ) (7) h t = o t tanh(C t ) (8) wherein: W o b o are the weight matrix and bias vector of the output gate, respectively; The RMSprop method is used as the optimization algorithm of the LSTM prediction model, and the algorithm design is as follows: E[g 2 ] is the squared mean of the gradient at the tth iteration; α is the power, set to 0.9; W t is the model parameter at the tth iteration; g t = ΔJ(W t ) is the gradient of the cost function of the t-th iteration model with respect to the parameters W. η0 is the global initial learning rate; ε is a value that tends to 0 infinitely; (3) Accuracy evaluation To measure the accuracy of the model prediction results, the coefficient of variation of root mean square error (CV-RMSE) and the normalized mean bias error (NMBE) are selected as evaluation indexes, and their calculation formulas are as follows: Wherein: n is the number of prediction samples; y i predi is the prediction value for the i-th sample; actual value for the i-th sample, According to the ASHRAE14-2014 standard, when using hourly data for prediction and verification, the CV-RMSE result is less than 30%, and the NMBE is less than 10%, which can be considered as accurate.

3. The ice storage air conditioning system economy and responsiveness dual objective optimization model predictive control method of claim 1, wherein: The step S2 is specifically: (1) Establishing a device energy consumption mathematical model The chiller and the water pump are the main energy-consuming devices in the ice storage air conditioning system, and the energy consumption of the system affects the amount of electricity taken from the power grid by the building and the operation cost of the user. The multiple polynomial model is the relationship between the COP and the three variables of the chilled liquid return temperature and the cooling water return temperature. The calculation accuracy is relatively high. The multiple polynomial is selected to establish the mathematical model of the chiller, and the calculation formula is as follows: wherein: b1~b 10 are parameters to be identified; Q e is the cooling load of the unit, kW; T ei is the refrigerant return temperature, °C; T ci is the cooling water return water temperature, °C; Due to the difficulty in obtaining the frequency data of the water pump during actual testing, the water pump model is often simplified as a single variable parameter function about the flow rate. Therefore, the energy consumption mathematical expression of the variable frequency water pump is as follows: P pump = a1+ a2- M + a3- M 2 (14) (2) Establishing a device dynamic response model 1) Modeling of heat transfer differential equation The reasons for the state change of the secondary refrigerant in the coil include heat transfer in the flow direction and the pipe diameter direction. The mass balance formula of the secondary refrigerant in the flow direction (axial direction) is as follows: G hf,in = G hf,out = G hf (15) where: G hf,in is the flow rate of the cold carrier at the inlet of the ice storage tank, kg / s; G hf,out Qc = 0.5 * Qs for the coolant at the outlet of the ice storage tank, kg / s; The energy balance equation of the secondary refrigerant side is as follows: wherein: c hf specific heat capacity of the secondary refrigerant, kJ / kg °C; p hf p 3 ; A p A is the cross-sectional area of the coil, m 2 ; L p L is the length of the coil, m; N is the number of nodes set at the end of the coil, that is, the number of state variables set from the input to the output; T hf,i and T hf,i+1 Ti and T0 are the inlet and outlet temperatures of the i-th ice storage coil unit, respectively, in °C. U is the heat transfer coefficient, W / (m 2 • °C); A l is the inner surface area of the coil tube wall, m 2 ; T pcm T for the temperature of the phase change material outside the coil, °C; The enthalpy method model represents the state change of the phase change medium in the form of enthalpy H in the energy balance equation, which is suitable for the entire phase change stage. According to the change of the phase change medium, the heat transfer process is divided into three stages: liquid phase stage, solid-liquid paste stage, and solid phase stage, Wherein: H is the specific enthalpy, kJ / kg; T m Phase transition temperature, 0 °C; ΔH m latent heat of phase change for the substance unit mass, 334 kJ / kg; c pcm Cp specific heat capacity of the phase change medium, kJ / (kg °C); The energy balance equation of the water outside the pipe is as follows: where: p pcm is the density of the phase change medium, kg / m 3 ; V pcm Volume of the phase change medium (ice-water paste region), m 3 ; λ pcm λ is the thermal conductivity of the phase change medium, W / (m-k); A pcm A is the cross-sectional area of the phase change medium, m 2 ; UA l The UA is the product of the overall heat transfer coefficient between the phase change medium and the coolant within the coil and the heat transfer area, and is equal to the inverse of the total thermal resistance of the radial heat transfer through the pipe. wherein: R hf is the convective heat transfer resistance between the flow of the refrigerant inside the ice storage coil and the coil wall, °C / W; R pipe Rthw is the thermal resistance of the tube wall, °C / W; R pcm Rth, °C / W, for the phase change medium outside the coil r i , r o , r pcm are the inner radius, the outer radius of the coil and the radius of the phase change medium, respectively, m; The convective heat transfer coefficient h between the coolant and the tube wall of the coil hf Related to the Reynolds number Re, calculated by the following formula: wherein: p, v, d, m are respectively the density (kg / m 3 ), the flow rate (m / s), the diameter (m) and the viscosity coefficient (m 2 ) of the coolant. The calculation formula of the Prandtl criterion number Pr is as follows: When Re>2300, the secondary refrigerant is in a turbulent state in the coil. The Nusselt criterion number Nu is calculated according to the Dittus-Boelter formula: Nu = 0.023Re 0.8 Pr 0.4 (25) When Re≤2300, the secondary refrigerant flows in a laminar state in the pipe, In the ice-melting release cold link, when the ice layer and the coil between the ice-melting time after the formation of concentric water ring, in addition to the thermal resistance, but also to consider the convective heat transfer resistance, but because of its calculation process is too complex, to establish and Rayleigh number R ac The concept of the relevant effective thermal conductivity to simplify the calculation, the effective thermal conductivity λ eff Empirical formula is: where: Ra is the Rayleigh number with δ = r w -r o as qualitative length, the average temperature of the inlet and outlet of the coolant as the qualitative temperature corresponding Rayleigh number; 2) Conversion of transfer function For a complex physical system, it may involve multiple input and output variables, and it is difficult to directly convert into a transfer function. At this time, the state space equation is often used as a bridge between the differential equation and the transfer function. A set of first-order differential equations with state variables are used as mathematical models to describe the "input-state-output" process. Based on the established multi-input multi-output state space equation, the transfer function in the s domain is conveniently converted. To get the expression of the coefficient matrix of the linear state space equation, first need to solve the steady-state solution corresponding to the steady-state differential equation, the steady-state equation is the algebraic equation corresponding to the guide number item is 0, the form of each variable in the equation is expressed as a steady-state value and a small increment, that is, where the symbol - denotes the steady state value of each variable, the symbol δ denotes a small increment of each variable; the overall heat transfer coefficient U is a dependent variable which varies with the change of the coolant flow rate G hf The overall heat transfer coefficient U is linearized using the first order Taylor expansion formula: The linear form of the above variables is brought into the heat transfer differential equation to obtain the linearization derivation result: liquid phase (water) Energy balance equation of the secondary side of the load: Energy balance equation of the water outside the tube: Solid-liquid paste phase (ice-water mixture) Energy balance equation of the secondary side of the load: Energy balance equation of the paste zone outside the tube: Solid phase (ice) Energy balance equation of the secondary side of the load: Energy balance equation of the ice outside the tube: Selecting state variables, input variables and output variables of the state space equation, wherein the state variables can reflect the change of the internal state of the system, T hf,in , G hf , H in are input parameters, T hf,i and H i with differential terms are state variables in the heat transfer process, the outlet temperature T hf,out of the cold carrier in the ice storage tank coil is an output variable, and the above two heat transfer differential equations can be converted into the following linear state space model: T hf,out = CT (39) The Laplace transform is performed on the state space equation of the system to obtain the input-output model of the ice storage coil as follows: T hf,out = G1(s)G hf + G2(s)T hf,in + G3(s)H in (40) Wherein: G1(s) is the transfer function of the secondary side of the load; G2(s) is the transfer function of the secondary side of the load; G3(s) is the transfer function of the initial enthalpy of the phase change medium with respect to the outlet temperature of the secondary side of the load; 3) Time delay time calculation In the process of participating in the energy system regulation and control, the total time delay time is defined as the time period between the time when the control parameter of the local energy storage device is adjusted according to the scheduling instruction and the time when the control output parameter changes, which is the coupling effect of the two time delay factors of thermal inertia and flow delay, which are represented by the inertia link and the pure delay link in the transfer function. (a) Thermal inertia time delay The thermal inertia is mainly derived from the heat transfer medium, which is defined as the heat transfer capacity of the medium due to its own heat capacity. When the system exchanges heat with the outside world, the temperature of the medium rises or falls, which needs a certain time, not immediately complete change. This characteristic is called thermal inertia. In the transfer function, the inertia element (also known as viscous damping) is described as when the input quantity changes suddenly, the output quantity cannot change suddenly, but can only change gradually according to the exponential law, which reflects that the control object has inertia. Any actual system with inertia element characteristics can store energy. In order to characterize the thermal inertia, first establish the energy balance relationship of the heat transfer process, and then perform Laplace transform on the partial differential equation to obtain the transfer function, so as to calculate the inertia time constant τ of the medium in the transfer function c : (b) Pure time delay The pure delay phenomenon of the system is mainly caused by medium transmission, and in the transfer function, the delay link is also called pure delay link, which represents that in the steady-state operation of the system, the state of the input and output is stable, and the steady-state change of the input needs to be delayed for a period of time before being transmitted to the output, which is represented as: y(t) = f(t - τ p ) (43) Wherein: y and u are output and input signals respectively; τ p is the pure delay time; Under the condition that the flow rate of the system is constant, the time of pure delay is almost fixed; while in the variable flow system, the time of pure delay will change with the change of the flow rate.

4. The ice storage air conditioning system economy and responsiveness dual objective optimization model predictive control method of claim 1, wherein: The step S3 is specifically: (1) Double target cost function establishment The present application simultaneously pays attention to the economy of the ice storage system according to the demand response of time-of-use electricity price and the dynamic characteristics of the system responding to the regulation and control signal, uses a prediction model and an optimization algorithm to obtain the economic best combination of the power contribution of the chiller and the ice storage tank in the prediction period, and on this basis, the cold quantity provided by the ice storage system is transmitted to the end user through the plate heat exchanger to meet the cold load demand of the end user, and the purpose of on-demand cooling is realized through the adjustment of the refrigerant flow rate; The controlled parameter is the secondary side chilled water supply temperature, that is, the set value to be tracked, in the cost function, the first term is the total operation cost of the system in a day, which is related to the equipment energy consumption and the electricity price in different time periods, which embodies the control objective of saving user operation cost and the flexibility of demand response; the second term is the square of the deviation between the predicted value of the control output temperature and the set value, which ensures that the control temperature reaches the set value as soon as possible, so that the system accuracy and followability are better, and thus the optimization objective function (also called cost function) is established: where: E ch and E pump are the energy consumptions of the individual chiller and water pump, respectively, in kW. N1 and N2 are the total number of water chillers and water pumps, respectively; Chilled water supply temperature setpoint, °C; (2) optimization algorithm Particle Swarm Optimization Algorithm (PSO) is used to solve the optimization problem in the model predictive control. There are N particles in the solution space constantly moving and searching for the optimal solution. This process is achieved through iteration. In the population, each particle is determined by the speed and position of the two parameters to move the state, and the particle can record the search process of each search movement. The corresponding optimal position of each particle in the search process is recorded, and then the best position of the comprehensive result of all particles is obtained. In the iteration process, the mathematical expressions of the speed v ij and position x ij of the i-th particle in the j-th dimension are as follows: v ij (t+1) = ωv ij (t) + c1r 1j (t) (p ij -x ij (t) ) + c2r 2j (t) (p gj -x ij (t) ) (46) x ij (t+1) = x ij (t) + v ij (t+1) (47) Where: t is the number of iterations; j is the number of variables; ω is the inertia weight; c1, c2 are acceleration coefficients; r 1j ,r 2j is a random number uniformly distributed in [0, 1]; p i is the optimal position for the particle individual; p g is the optimal position in the population of particles; (3) feedback correction Based on the time delay time calculation method, a time delay compensation module is added. At each control node, the delay time is calculated according to the control increment, i.e. the state of the temperature and flow of the secondary refrigerant, and the signal is transmitted to the control sequence signal, so as to eliminate the time delay time in the next response. Through the rolling prediction and feedback correction link, the real-time performance and robustness of the control system are ensured. As a kind of advanced control method, the model predictive controller predicts the control effect of the expected output signal of the input control in the future time period in advance through the obtained system scheduling signal, so as to improve the response speed of the system and reduce the influence of time lag.

Images