*Article* **Research on a Day-Ahead Grouping Coordinated Preheating Method for Large-Scale Electrified Heat Systems Based on a Demand Response Model**

**Guodong Guo \* and Yanfeng Gong**

State Key Laboratory of Alternate Electrical Power System with Renewable Energy Sources, North China Electric Power University, Changping District, Beijing 102206, China **\*** Correspondence: gbjdsf@163.com; Tel.: +86-183-22595221

**Abstract:** In recent years, the increasing winter load peak has brought great pressure on the operation of power grids. The demand response on the load side helps to alleviate the expansion of the power grid and promote the consumption of renewable energy. However, the response of large-scale electric heat loads to the same electricity price curve will lead to new load peaks and regulation failure. This paper proposes a grouping coordinated preheating framework based on a demand response model, which realizes the interaction of information between the central controller and each regulation group. The room thermal parameter model and the performance map of the inverter air conditioner/heat pump are integrated into the demand response model. In this framework, the coordination mechanism is adopted to avoid regulation failure, an edge computing structure is applied to consider the users' preferences and plans, the grouping and parallel computing structure is proposed to improve the computing efficiency. Users optimize their heat load curves based on a demand response model, which can consider travel planning and ensure user comfort. The central controller updates the marginal cost curve based on the predicted scenario set to coordinate the regulation groups and suppress the new peaks. The simulation results show that the proposed method can promote the consumption of renewable energy through coordinated preheating and reduce the system energy consumption cost and user bills. The parallel computing structure within the regulation group also ensures the computing efficiency under large-scale loads.

**Keywords:** demand response; coordinated preheating; inverter air conditioner; equivalent thermal parameter model; smart grid

## **1. Introduction**

In recent years, the peak power consumption in winter has become more and more obvious. The maximum load in winter exceeds that in summer in many southern provinces in many southern provinces of China [1], and Texas in the United States has also set a record for the peak power consumption in winter in 2021 [2]. Different from the central heating in the north, the household independent heating mode is more common in hot summer and cold winter regions, with the characteristics of intermittent heating on demand [3]. It is important to promote heating electrification to achieve carbon neutrality in winter heating systems. The heat pump, hot air blower/cooling and heating air conditioner follows the reverse Carnot cycle and transfers more heat with less energy consumption, which is highly energy efficient. Jiang [4] pointed out that the air source heat pump is the most important way of heating electrification. The European Commission has also set a heat pump development target, that is, 40% of residential buildings and 65% of commercial buildings are expected to achieve electric heating by 2030 [5]. With the increase in power demand in the winter peak period [6], efficient demand management technology can help to reduce system costs, promote the consumption of renewable energy, and achieve carbon neutrality.

**Citation:** Guo, G.; Gong, Y. Research on a Day-Ahead Grouping Coordinated Preheating Method for Large-Scale Electrified Heat Systems Based on a Demand Response Model. *Appl. Sci.* **2022**, *12*, 10758. https:// doi.org/10.3390/app122110758

Academic Editors: Chun Sing Lai, Kim-Fung Tsang and Yinhai Wang

Received: 24 September 2022 Accepted: 21 October 2022 Published: 24 October 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

Residential air conditioning plays an important role in demand response resources. A demand response model for single residential buildings was established in [7], in which heat pumps preheat at low electricity prices to reduce power demand at high electricity prices. The home energy management system with integrated intelligent heat load has been studied in [8,9] and the corresponding optimization problems were found, with the purpose of minimizing energy consumption cost and ensuring user comfort. Comfort is actually guaranteed by the precooling/preheating. Large scale day-ahead heat load regulation is an attractive solution for power system scheduling due to its economy and security, but there are many difficulties to be solved. The assessment of the maximum adjustable load of the air conditioning cluster can be used for real-time regulation [10,11], but it cannot be used for day ahead scheduling because of the neglect of load-time coupling. For example, the time point and adjustment amount of preheating/precooling cannot be determined by the above methods. Equivalent energy storage models and cluster models for air-conditioning clusters are studied for day-ahead scheduling [12,13], in which temperature setpoints are regarded as consistent and fixed. In fact, the temperature setpoints of users are different and time-variant, and are related to their respective travel schedules. Most importantly, the scheduling results are not practical due to the lack of instructions specific to each user. In view of the above problems, some studies have been carried out from the characteristics of intermittent heating/cooling [14–16] and the respective thermal parameters and thermal demands of households [17]. However, the above methods are only applicable to a single residential or commercial building and will encounter difficulties in solving large-scale problems. It is worth noting that if all participants adjust based on the same electricity price curve, it will fundamentally change the marginal cost curve of the system, which will lead to failure of the adjustment. A new coordinated preheating scheme based on game theory was proposed in [18] to ensure the effectiveness of large-scale family collective preheating, which effectively takes into account changes in marginal costs in the coordination process. However, households need to iterate one by one, and the computation time is linearly related to the household scale, which is unacceptable in large-scale problems. Intermittent heating, personalized thermal demand, and computational efficiency are the three major difficulties in large-scale day-ahead thermal load regulation.

Accurate house thermal models [19] and air conditioning models have a significant impact on the conditioning effect. In terms of air conditioning, constant frequency air conditioning has been gradually replaced by inverter air conditioning (IAC). The performance of IAC is related to the indoor/outdoor temperature and the compressor speed, and its steady-state model can be used to study the coupled dynamic characteristics of the room and IAC. Research and experiments on steady-state models are abundant [20–22], but their computation time is unacceptable for scheduling problems. In [16,23], the performance map based on the steady-state model was obtained for the direct control of IAC.

This paper proposes a day-ahead group coordinated preheating method based on demand response model for large-scale electric heating load, which is carried out in a framework composed of a central controller and several regulation groups. Under this framework, users can reduce electricity bills through the proposed demand response model, in which personalized settings such as travel schedules and user temperature demand curves can be fully considered. The central controller updates the marginal cost curve after each round of regulation and transmits it to the next group, and the interaction avoids new peaks. In each round of adjustment, households in the regulation group solve their respective optimization problems in parallel, ensuring computational efficiency. In addition, a room equivalent thermal parameter (ETP) model and an IAC model are established to form a single household demand response model. In order to quickly obtain the performance parameters of IAC under given conditions, this paper develops a performance map based on the steady-state model of IAC, which can be applied to the direct control of the compressor frequency to determine the power consumption. And the mapping can be easily transformed into piecewise linear constraints and added to the optimization problem.

This article is organized as follows. In Section 3, a detailed ETP model of residence and the performance map based on the IAC model are established. In Section 4, a single household demand response model is given, and the grouping coordinated preheating framework is proposed. In Section 5, the performance of the proposed methods are presented through numerical simulation. Finally, the main conclusions are discussed in Section 6.

#### **2. Methodology**

#### *2.1. Research Objectives*

The preheating model is an economical and efficient load side regulation method to realize economic savings of a single household. However, regulation failure will occur such as new peaks and increasing system costs when preheating without coordination, which is due to the lack of interaction between global interests and demand scheduling. The response of large-scale electric heat loads to the same electricity price curve will lead to new load peaks and regulation failure. And most of the existing studies take the thermal comfort into consideration by reducing the temperature deviation with the desired value. However, the user's temperature demand is time-variant, which is related to the user's travel planning. In other words, many existing studies cannot take into account users' personal preferences and travel plans. Finally, considering the practicability of the model, the solution time of large-scale preheating planning must meet the scheduling requirements, which is extremely challenging. The research objectives of this paper are summarized as follows:


#### *2.2. Research Method*

In order to solve the problem of regulation failure under large-scale preheating, a coordinated preheating mechanism is proposed which links the global interests with the demand side response. To formulate a customized heat consumption plan for each user according to their temperature preferences and travel planning, a kind of edge computing and central regulation framework is applied to the coordination mechanism. In addition, to meet the scheduling time requirements, a grouping and parallel computing structure is proposed. According to the above logic, the research architecture is shown in Figure 1.

To formulate the single household regulation model with IAC, a detailed ETP model of residence and the performance map of the IAC model need to be established, which can link the indoor temperature change with the user energy consumption scheduling. Then the grouping coordinated preheating mechanism based on edge computing and central regulation framework is proposed, which links the single household regulation model and the changes of marginal cost curve.

#### *2.3. Simulation Parameters*

The simulation is implemented on the python platform. CoolProp is called in the IAC system simulation to obtain physical properties, and Gurobi is called as the solver for the single household demand response optimization. The building parameters and thermal parameters are shown in Table 1. User status division and distribution of parameters can be seen in Section 5.

**Figure 1.** Research architecture.

**Table 1.** The parameters of the room thermal model.


In this paper, a certain type of apartment in southern China is used as the standard type. The geometric parameters are: length 13.6 m, width 8.6 m, and height 2.6 m. The residence has a north/south external wall, with one external wall on the east or west side, totaling three external walls. The window-to-wall ratio of each orientation is 0.25 in the north direction; 0.35 in the south direction; and 0.2 in other directions.

#### **3. Residential Electric Heating System Model**

#### *3.1. Room Thermal Model*

In order to link the indoor temperature change with the user energy consumption scheduling, a residential ETP model is established, which consists of the outdoor part, the exterior wall and the indoor part, as seen in Figure 2. There are three main modes of heat transfer in the model: heat conduction, heat convection and heat radiation. The outer wall is equivalent to a thermal resistance and two thermal capacitances, which are denoted as *R* and *C*, respectively. For the external surface of the wall, it obtains the heat gain *Q*s*olar*,*<sup>w</sup>* from solar radiation, obtains heat from the inside of the wall, and transfers heat to the outdoor air by convection. Equation (1) depicts the change rate of the exterior wall temperature *Twe*, which is related to heat transfer. Similarly, Equation (2) depicts the change rate of the internal wall surface temperature *Twi*, which is influenced by both indoor air convection and wall conduction heat transfer. Equation (3) gives the energy balance equation of indoor air. Indoor air has a certain heat storage capacity, and it exchanges heat with internal mass, the air outside walls and windows. In the model, indoor air also absorbs heat gain *Q*s*olar*,*<sup>a</sup>* from solar radiation, internal heat gain *Qgain*,*<sup>a</sup>* and heat *QAC* generated by IAC, but will leak heat *Qlk* due to the gap at the junction. Equation (4) gives the energy balance equation of internal mass. The internal mass such as room partitions and furniture have large thermal inertia, and their equivalent heat capacity *C<sup>m</sup>* can maintain the slow change of mass temperature *Tm*. In addition to the heat exchange with indoor air, the internal mass absorbs the solar radiation heat gain and internal heat gain, which are denoted as *Q*s*olar*,*<sup>m</sup>* and *Qgain*,*m*. Equations (5)–(7) are the expressions of solar heat gain absorbed by exterior walls, internal mass and indoor air, respectively. Equation (8) is the heat loss. Equations (9) and (10) represent the absorption of internal heat gain by mass and air respectively. Equation (11) represents three sources of internal heat gain: household appliances, lighting and human body heat radiation.

$$\mathcal{C}\_{w}\frac{dT\_{w\epsilon}}{dt} = \frac{T\_o - T\_{we}}{R\_{wo}} + \frac{T\_{wi} - T\_{we}}{R\_{w}} + Q\_{solar,w} \tag{1}$$

$$\mathcal{C}\_{w}\frac{dT\_{wi}}{dt} = \frac{T\_{we} - T\_{wi}}{R\_{w}} + \frac{T\_{a} - T\_{wi}}{R\_{wi}} \tag{2}$$

$$\mathcal{C}\_{a}\frac{dT\_{a}}{dt} = \frac{T\_{m} - T\_{a}}{R\_{am}} + \frac{T\_{wi} - T\_{a}}{R\_{wa}} + \frac{T\_{o} - T\_{a}}{R\_{win}} + \mathcal{Q}\_{solar,a} + \mathcal{Q}\_{gain,a} + \mathcal{Q}\_{AC} - \mathcal{Q}\_{lk} \tag{3}$$

$$\mathcal{C}\_m \frac{dT\_m}{dt} = \frac{T\_a - T\_m}{R\_{am}} + Q\_{solar,m} + Q\_{gain,m} \tag{4}$$

$$Q\_{solar,w} = k\_{solar} A\_w I\_{solar} \tag{5}$$

$$Q\_{solar,m} = f\_{solar,m} \times SHGC \times A\_{win} I\_{solar} \tag{6}$$

$$Q\_{solar,a} = f\_{solar,a} \times SHGC \times A\_{win} I\_{solar} \tag{7}$$

$$Q\_{leak} = \rho C\_p V\_{room} \times ACH \times \left(T\_a - T\_o\right) / 3600 \tag{8}$$

$$Q\_{\mathcal{Z}^{\text{ain,m}}} = f\_{\mathcal{Z}^{\text{ain,m}}} Q\_{\mathcal{Z}^{\text{ain}}} \tag{9}$$

$$Q\_{gain,a} = f\_{\mathcal{S}^{ain,a}} Q\_{gain} \tag{10}$$

$$Q\_{\text{gain}} = Q\_{\text{equip}} + Q\_{\text{lamp}} + Q\_{\text{occup}} \tag{11}$$

where *ksolar* denotes absorptance of surface for solar radiation; *A<sup>w</sup>* denotes geometric area of the exterior wall; *Isolar* denotes solar radiation; *SHGC* represents solar heat gain coefficient; *fsolar,m* and *fsolar,a* represent the radiative/convective split for the solar heat gain respectively; *A<sup>w</sup>* denotes geometric area of the window; *ρ* is the air density; *C<sup>p</sup>* denotes the specific heat of air at constant pressure. *Vroom* represents the room volume; *ACH* denotes the air exchange per hour; *fgain,m* and *fgain,a* represent the radiative/convective split for the internal heat gain respectively. The thermal dynamic model of the room can be added to the convex

optimization problem as constraint conditions, which can be easily solved with advanced optimization techniques. The discretization expression is shown in (12)–(15).

$$T\_{\rm PU}(t) = \left(1 - \frac{\Delta t}{\mathbb{C}\_{\rm u} R\_{\rm uv}} - \frac{\Delta t}{\mathbb{C}\_{\rm u} R\_{\rm uv}}\right) T\_{\rm uv}(t-1) + \frac{\Delta t}{\mathbb{C}\_{\rm u} R\_{\rm uv}} T\_0(t-1) + \frac{\Delta t}{\mathbb{C}\_{\rm u} R\_{\rm uv}} T\_{\rm uv}(t-1) + \frac{\Delta t}{\mathbb{C}\_{\rm u}} Q\_{\rm solr, w}(t-1) \tag{12}$$

$$T\_{\rm wi}(t) = \left(1 - \frac{\Delta t}{\mathcal{C}\_{\rm w} R\_{\rm w}} - \frac{\Delta t}{\mathcal{C}\_{\rm w} R\_{\rm wi}}\right) T\_{\rm wi}(t-1) + \frac{\Delta t}{\mathcal{C}\_{\rm w} R\_{\rm w}} T\_{\rm we}(t-1) + \frac{\Delta t}{\mathcal{C}\_{\rm w} R\_{\rm wi}} T\_{\rm d}(t-1) \tag{13}$$

$$\begin{array}{c} T\_{\mathfrak{a}}(t) = \left(1 - \frac{\Lambda t}{\mathbb{C}\_{\mathfrak{a}}\mathbb{R}\_{\text{un}}} - \frac{\Lambda t}{\mathbb{C}\_{\mathfrak{a}}\mathbb{R}\_{\text{un}}} - \frac{\Lambda t}{\mathbb{C}\_{\mathfrak{a}}\mathbb{R}\_{\text{un}}}\right) T\_{\mathfrak{a}}(t-1) + \frac{\Lambda t}{\mathbb{C}\_{\mathfrak{a}}\mathbb{R}\_{\text{un}}} T\_{\mathfrak{m}}(t-1) + \frac{\Lambda t}{\mathbb{C}\_{\mathfrak{a}}\mathbb{R}\_{\text{un}}} T\_{\text{un}}(t-1) \\ \qquad + \frac{\Lambda t}{\mathbb{C}\_{\mathfrak{a}}\mathbb{R}\_{\text{un}}} T\_{\mathfrak{o}}(t-1) + \frac{\Lambda t}{\mathbb{C}\_{\mathfrak{a}}} \left(Q\_{\text{solar},\mathfrak{a}} + Q\_{\text{gain},\mathfrak{a}} + Q\_{\text{AC}} - Q\_{\text{lk}}\right) \end{array} \tag{14}$$

$$T\_m(t) = \left(1 - \frac{\Delta t}{\mathbb{C}\_m \mathbb{R}\_{am}}\right) T\_m(t - 1) + \frac{\Delta t}{\mathbb{C}\_m \mathbb{R}\_{am}} T\_a(t - 1) + \frac{\Delta t}{\mathbb{C}\_m} \left(Q\_{solar, m} + Q\_{gain, m}\right) \tag{15}$$

**Figure 2.** ETP model of residences.

## *3.2. Performance Maps of IAC*

## 3.2.1. Performance Calculation of IAC

Steady-state IAC models can be used to study the cycling characteristics of IAC under different environment and control conditions. Specifically, the performance of IAC includes power consumption *PAC* and coefficient of performance *COP*, which are affected by compressor speed *Ncomp*, outdoor temperature *T<sup>o</sup>* and indoor temperature *Ta*. Figure 3 shows the main components of the air conditioning heating system and the cycle process, which meets the conservation of energy and mass. Starting from point *A*, the cycle is analyzed. The specific enthalpy *h<sup>A</sup>* of refrigerant steam is calculated from the evaporation temperature *Tevap* and the degree of superheat *SH*, and the steam changes into high-temperature and high-pressure steam through the compressor module (including the accumulator, suction pipe, compressor and exhaust pipe, etc.). Then the compressor power *Wcomp* and the mass flow rate *mcomp* can be obtained. The high-temperature and high-pressure steam releases heat to the indoor air through the condenser and cools to subcooled liquid. The degree of subcooling *SC,* is calculated by the condenser model. The refrigerant liquid with medium temperature and high pressure is depressurized by the thermal expansion valve, and the mass flow rate *mexp* is calculated according to the expansion valve model. The low-temperature and low-pressure liquid absorbs heat and turns into a gas through the evaporator, and the output specific enthalpy *hevap* is calculated according to the evaporator

model. There are three independent variables in the above cycle, namely [*Devap*, *Dcond*, *SH*], where *Devap* is the difference between the inlet temperature at the air side of the evaporator and the evaporation temperature, and *Dcond* is the difference between the condensation temperature and the inlet temperature at the air side of the condenser. Accordingly, the cycle should also meet three balance constraints, namely, specific enthalpy constraints at the start and end of the cycle, condenser subcooling constraints and mass flow balance constraints. The constraints are represented by residuals [∆1, ∆2, ∆3] respectively. When the steady state model converges, the residual values are all zero.

$$P\_{A\mathbb{C}} = \mathbb{W}\_{comp} + \mathbb{W}\_{cond,fan} + \mathbb{W}\_{evap,fan} \tag{16}$$

**Figure 3.** Schematic diagram of AC heating cycle.

Figure 4 shows the numerical calculation flow of the system, which is mainly divided into four parts: initial value calculation, independent variable update, pressure drop correction and performance parameter calculation. The operating condition parameters include the indoor air temperature *Ta*, the outdoor temperature *Tout* and the rotor speed *Ncomp* of the inverter compressor. The model first solves the independent variables in the minimal model as the initial value of the iteration. Secondly, the Broyden algorithm is used to solve the independent variable iteratively to make the residual close to zero, and the pressure drop at the high pressure side ∆*Ph+*, and low pressure side ∆*pl-* are calculated when converged. The pressure drop is introduced into the compressor model, and the pressure drop and independent variables are updated iteratively under this condition. The iteration stops when the pressure drop difference is less than the threshold. Finally, the performance parameters are calculated. As shown in Equations (16) and (17), the power consumption of IAC is mainly related to the compressor, condenser fan and evaporator fan. At the same time, the electric energy consumed by the condenser fan is eventually converted into the internal energy of the indoor air. Therefore, in addition to the heat release of the condenser, the numerator of *COP* also contains the internal energy.

#### 3.2.2. 3D Storage of the Mapping

The calculation of the steady state model is time-consuming, while the demand response problem needs to ensure the accuracy and speed of the calculation. To this end, this paper calculates the performance parameters of IAC at each operating point in advance, and uses the 3D map to store the mapping relationships.

**Figure 4.** Flow chart of air conditioning system simulation calculation.

As shown in Figure 5, the power consumption of IAC increases with the rotational speed, while *COP* increases first and then decreases with the speed. As shown in Figure 5b, the *COP* corresponding to 30 Hz is lower than that of 40 Hz. At the same frequency, the closer *T<sup>a</sup>* and *Tout* are, the higher the *COP*. According to Figure 5a,b, there is the following mapping: *QAC* = *f*(*PAC*, *Ta*, *Tout*). With *T<sup>a</sup>* and *Tout* determined, *QAC* can be expressed as a piecewise linear function of *PAC*, denoted as *fPQ*. It can be easily added to mixed integer linear problems as constraints. *Tout* is obtained from day-ahead weather forecasts, while *T<sup>a</sup>* is related to the temperature set by the household.

**Figure 5.** Performance maps of IAC.

#### **4. Day-Ahead Coordinated Preheating Control**

The preheating control utilizes the thermal inertia of the room to transfer heat loads to improve the economy of the heating system. Specifically, houses are preheated with cheap electricity during the high generation period of renewable energy or the low load period of the system, so as to reduce the energy consumption of marginal units during the high load period of the system. At present, most of the southern areas of China are in the mode of intermittent air conditioning heating, and the heat load is flexible. This paper firstly proposes the heat load scheduling of a single house to minimize the operating cost of IACs while meeting the user's temperature requirements. We then propose an efficient grouping coordinated control framework considering that the large-scale heat load regulation will change the marginal cost of the system.

#### *4.1. Single Household Demand Response Regulation*

The objective of single household demand response regulation (SDR) is to reduce user bills within the scheduling period while meeting the user's temperature demands. The optimization problem is formulated as Equations (16)–(18).

$$\min\_{P\_{A\mathbb{C}}} \sum\_{t \in \mathcal{T}} P\_{A\mathbb{C}}(t) \Delta t \cdot pr\_t + \rho\_{\mathcal{C}} e\_t \tag{18}$$

$$s.t.\text{ Eq.(12)} \sim \text{Eq.(15)}$$

$$Q\_{A\mathcal{C}}(t) = f\_{PQ}(P\_{A\mathcal{C}}(t), T\_a^\*(t), T\_{out}(t))\tag{19}$$

$$T\_{lb,t} + e\_t \le T\_{\mathfrak{a}}(t) \le T\_{\mathfrak{u}b,t} + e\_t \tag{20}$$

Equation (18) is the objective function, where *T* is the scheduling interval; ∆*t* is the scheduling time step; *pr<sup>t</sup>* is the electricity price at time *t*; *e<sup>t</sup>* ≥ 0 is a slack variable, used to avoid failure to solve optimization problems under tight constraints, and *ρ<sup>e</sup>* is the corresponding penalty factor. Equations (12)–(15) are the discretized equivalent heat balance constraints, and each time step corresponds to four heat balance equations. Equation (19) is the *COP* mapping relationship of IACs, which is discretized into a piecewise linear constraint of *PAC* → *QAC*. The addGenConstrPWL function in Gurobi makes it easy to add to constraints. Equation (20) represents the user's demand constraints at different times, that is, the room temperature should be within the upper and lower limits set by the user. The upper and lower temperature limits at time *t* are respectively denoted as *Tlb,t* and *Tub,t*. There are different settings for working and sleeping. In addition, considering the threshold of IAC in actual operation, *PAC* is set as a semi continuous variable, that is, it should meet *PAC*(*t*) = 0 or *Plb*,*<sup>t</sup>* ≤ *PAC*(*t*) ≤ *Pub*,*<sup>t</sup>* .

#### *4.2. Grouping Coordinated Preheating Framework*

Although the above optimization model can reduce user bills, if all users schedule based on the same price curve, the marginal cost of power generation will change, which cannot guarantee the reduction in total power generation cost. In order to solve this problem, this paper proposes a grouping coordinated preheating framework with edge computing and central coordination.

As shown in Figure 6, each household acts as an edge computing unit to solve the demand response problems (18)–(20) according to the received marginal electricity price. The central controller is the key to coordinate all groups, which contains three functional modules: user grouping, marginal cost calculation and specifying scheduling group. According to the central controller, the user adjusts the heat load planning based on the updated marginal cost. And in each round of adjustment, only a small number of users participate in scheduling so that the marginal cost will not be significantly changed.

**Figure 6.** Edge computing architecture for coordinated preheating.

The interactive process of coordinated preheating is shown in Figure 7, which consists of the original load upload and several scheduling rounds. First, users of all groups upload the original load curve shown in <sup>1</sup> . Next, Group # 1 ~ Group # G participates in each scheduling round in turn. Taking Group # 1 as an example, each round incorporates three steps: (1) The central controller sends the marginal generation cost curve to Group # 1 as process <sup>2</sup> . (2) All local controllers in Group # 1 receive marginal cost curve based on Equations (18)–(20) to update the heat load curve. (3) The local controller uploads the new load curve, namely process <sup>3</sup> . After the central controller updates the marginal cost curve, Group # 2 conducts the next round of adjustment. The iteration will continue until the groups balance or the set number of rounds are met. The process is summarized as Algorithm 1.

#### **Algorithm 1 Coordinated Preheating Control Algorithm**

**Grouping**: Group #1 ~ Group #G

**Local controller input**: Day-ahead travel planning and temperature range of each user **Central controller input**: Generation cost function

1: Each user sets the day-ahead travel plan and acceptable temperature range, then the local controller calculates and uploads the heat load curve based on the on-off control.

2: The central controller calculates the marginal generation cost curve λ**n** based on the aggregated load curve and passes it to Group #1.

3: **for** iteration = 1 to *N* **do**

4: After receiving the marginal cost curve, all users within selected group get their own optimized heat load curve based on Equations (18)–(20). The local controller uploads the curve. 5: The central controller updates the load curve, calculates the new marginal generation cost curve and passes it on to the next group.

6: **end for**

**Figure 7.** Interactive process of coordinated preheating.

#### **5. Simulation Results**

This paper will evaluate the effect of coordinated preheating from indicators such as energy consumption cost, user bills, renewable energy consumption, and iteration rounds.

The demo is a regional power system including 20,000 households, in which the installed capacity of wind power is 30 MW. All users are divided into ten groups to ensure that the adjustment of a single group has less impact on the system. *Tlb,t* and *Tub,t* in Equation (20) are related to the user status, which includes awake at home, out of home and sleeping. The user status is divided by the time points shown in Table 2, where *tup* denotes wake-up time, *tleave* denotes departure time, *treturn* denotes home time, *tdown* denotes bedtime. In addition, the relationship between these time points is as follows: *tup* = *tleave* − *dmor*, *tdown* = *treturn* + *deve*, in which *dmor* and *deve* are the awake time spent at home in the morning/evening respectively. This paper assumes that *tleave* and *treturn* follow truncated normal distribution, and *dmor* and *deve* follow uniform distribution, as shown in Table 3. The value of the time step is an integer between 0 and 96, and the time step is the scheduling interval, which is 15 min.

#### **Table 2.** User status division.


**Table 3.** Distribution of parameters.


The generation cost function required by the central controller can be in any form, including expressions, step function curve and even calculation programs. It is worth mentioning

that the central controller only needs to perform a limited number of marginal cost calculations based on the generation cost during iteration, which can be extended to any user scale and network topology. As in [24–26], this paper sets the cost function in the form of quadratic function, that is, *Ctotal* = *ax*<sup>2</sup> + *bx* + *c*. The coefficient of the quadratic term is 5 ¥/MW<sup>2</sup> h, the coefficient of the primary term is 200 ¥/MWh, and the constant term is 800 ¥.

#### *5.1. Coordinated Preheating Results and Comparison*

The simulation will be implemented within the time range of 04:00~24:00. In order to comprehensively analyze the performance of the coordinated preheating algorithm, this paper compares several cases, and the results are shown in Figure 8. The gray area in the figure represents the base load, and the yellow area represents the wind power output. The base load is not adjustable. The black solid line is the baseline heating schedule, in which all users set heating schedules by on-off mode according to the set temperature. The solid red line is the result of coordinated preheating. It can be seen that the curve is smoother than the baseline because part of the peak load is shifted to the valley. And the red solid line is almost all above the yellow area, so the wind power consumption is greatly promoted. Another two cases were analyzed to demonstrate the superiority of the coordinated preheating strategy. In case 1, all households schedule their heat load curves based on the baseline electricity price to reduce the energy bill on the premise of meeting the temperature demand. The baseline electricity price curve is the marginal generation cost curve in the baseline case. It can be seen that in case 1, a forward and higher peak is formed due to the neglect of marginal cost change, which leads to excessive preheating and heat loss instead. Therefore, concentrated scheduling in the absence of coordination will lead to preheating failure and increase the cost of power generation. To further illustrate the necessity of coordination, case 2 is used for comparison. In case 2, all households schedule their heat load curves based on the marginal cost curve of the coordinated preheating case, which can be regarded as the optimal electricity price curve. However, it also failed to reduce the peak, indicating that the scheduling based on a single price curve is not feasible and the coordination is necessary.

**Figure 8.** Total load under different conditions.

A total of 30 adjustment rounds were conducted in the simulation, namely Round # 1~Round # 30. In Round #1, Group #1 performs rescheduling based on the electricity price given by the central controller. Similarly, Group #10 performs rescheduling in Round #10. In Round # 11, it returns to Group # 1. The change of each evaluation index with the number of rounds, is shown in Figure 9. Figure 9a shows the change of total generation cost, which is reduced by 36.8% due to coordinated preheating. Figure 9b shows that the total electricity bill of users has decreased by 48.3%. In addition, it can be seen from Figure 9c that wind power curtailment is significantly reduced, and almost all wind power can be consumed. It is worth noting that the total user bill decreases more than the total

generation cost, which is attributed to the smoother load curve after preheating adjustment. The two values are quite close when reaching the group equilibrium, which indicates that the proposed algorithm achieves the optimal adjustment.

**Figure 9.** Changes of Evaluation Indexes with the Adjustment Rounds.

In addition, it can be found that all indicators have reached a balance in ten rounds. At this time, all groups performed a round of adjustment, and subsequent optimization produced only slight fluctuations. In the simulation, limited by the computing power of a single computer, users in the group need to schedule one by one, so it takes a long time. However, in practical applications, with the help of users' edge computing capabilities, users in the group solve optimization problems in parallel, and the adjustment time for each round is within 10 s. The whole coordinated preheating process can reach equilibrium within 2 min, which is suitable for day-ahead scheduling. In addition, the time-consuming is related to the number of groups, but not related to the user scale due to the group-bygroup coordination method, which is extensible for super-large-scale coordination. On the contrary, the one-by-one coordination method requires 20,000 rounds of adjustment in this example, whose calculation time is linearly related to the user size, so it is not scalable.

Table 4 presents the performance comparison under different models, and two conclusions can be drawn: (1) A coordination mechanism is necessary. Preheating without coordination will lead to the failure of adjustment. Its net peak and generation cost are even higher than those without preheating; (2) A grouping mechanism is necessary, which

**Generation Cost (10,000 ¥) User Bill (10,000 ¥) Net Peak (MW) Wind Power Curtailment (MW) Required Rounds** No preheating 15.90 21.92 73.33 88.34 — Preheating without coordination 22.42 33.16 92.89 50.32 <sup>1</sup> Coordinated preheating group by group 10.04 11.33 38.40 1.03 <sup>20</sup> Coordinated preheating one by one — — — — 20,000 (Unacceptable)

can achieve scale-independence, while coordination one by one is unacceptable in terms of computing time.

**Table 4.** Comparison of different preheating effects.

In order to analyze the impact of each round of adjustment on the system, the important process curves are given in Figure 10, and their meanings have been marked in the figure. It can be seen from Figure 10a that both the morning and evening peaks of the total load decrease with the adjustment. The peak decreases more obviously in the first few rounds, and gradually approaches the balance in the later rounds. It is worth noting that compared with the morning peak, the evening peak has dropped more significantly. This is because the preheating period before the evening peak is long and there is no specific heat demand (the user has not yet arrived home), so the adjustable range is large. The preheating period before the morning peak has strict temperature restrictions (to ensure a proper sleep temperature), and the adjustment range is limited. In addition, since the earlier the preheating is, the greater the heat loss will be, the preheating period will not be too long. The gray part in Figure 10b represents the wind power curtailment area, and it can be seen that the curve gradually leaves the area with adjustment. Figure 10c shows the change process of the marginal cost curve, which becomes smoother with adjustment.

**Figure 10.** Process display of each round.

#### *5.2. Analysis of Temperature and Energy Consumption in the Coordinated Preheating Process*

To further explore the details of the preheating process, Figure 11 shows the temperature and power curves of users participating in Round #1, including the adjusted and original curves. In the temperature chart, the blue dotted line represents the temperature limit when awaking at home; The red dotted line represents the temperature limit when sleeping at home; The green area represents the sleeping period set by the user; The grey area refers to the period of awaking at home. Therefore, the indoor temperature curve must be in the orange and blue areas to meet user demands. The original curve adopts an on–off adjustment mode to ensure that the temperature does not exceed the limit.

**Figure 11.** Comparison between optimization curves and original curves of a family participating in the Round #1.

Figure 11b is the original temperature curve, and Figure 11d is the corresponding IAC power curve. It can be seen that the black solid line representing the indoor air temperature is located at the bottom of the orange and blue areas, which both meet the user's temperature demands and save energy. From 04:00 to 07:15, the user was sleeping, so the IAC adjusts the power to ensure that the temperature does not exceed the lower limit. The power is increased around 07:00 to ensure a comfortable temperature when the user wakes up. After the user leaves home, the IAC stops running, so *T<sup>a</sup>* and *T<sup>m</sup>* starts to drop. *T<sup>m</sup>* decreases more slowly than *T<sup>a</sup>* because the heat capacity of internal mass is larger than that of air. At this time, the internal mass transfers heat to the air. As the temperature difference between *T<sup>a</sup>* and *Twi* gradually decreases, the heat loss becomes smaller, so the

temperature drops faster first and then slower. The IAC power starts to rise at about 18:15 to ensure a comfortable temperature when users get home. As the air heat capacity is small, *T<sup>a</sup>* rises rapidly. And due to the large temperature difference between the air and internal mass, the air begins to transfer a lot of heat to internal mass until temperature *T<sup>m</sup>* is close to temperature *Ta*. In addition, as the temperature demand decreases during sleep, although the IAC power drops to zero at 23:15, the temperature still meets the demands when the user sleeps (around 23:30).

Figure 11a,c are the temperature and power curves after preheating adjustment. The marginal cost curve received by the user is the baseline curve in Figure 10c. It can be seen that the marginal cost is lower in the period of 04:00~04:30 and the period out of home, so the two periods are selected for preheating. Compared with the original power curve, the power of the IAC is higher during the preheating period, and the room temperature also increases while keeping within the set range. It is worth noting that, in order to reduce heat loss, the preheating power increases gradually in the second half of the time out of home. When users get home, the IAC is shut down for a long time, while the indoor temperature is still in a comfortable range, thus reducing the energy consumption during the high marginal cost period.

#### **6. Conclusions**

In this paper, a grouping coordinated preheating method based on a demand response model is proposed for large-scale electric heating load, and the effectiveness of the proposed method is verified by simulation. The main conclusions are as follows:


**Author Contributions:** Conceptualization, G.G. and Y.G.; methodology, G.G. and Y.G.; investigation, G.G. and Y.G.; writing—original draft preparation, G.G. and Y.G.; writing—review and editing, G.G. and Y.G.; supervision, G.G. and Y.G. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** All data included in this study are available upon request by contact with the corresponding author.

**Conflicts of Interest:** The authors declare no conflict of interest.
