Flexibility management with virtual batteries of thermostatically controlled loads: real-time control system and potential in Spain

Virtual batteries composed of aggregated thermostatically controlled loads are able to provide real-time frequency regulation to electrical grids. Load flexibility management can be helpful in solving the problem of balancing generation and demand, which is becoming more complex due to the variability of renewable energies. A real-time virtual battery control system is presented in this paper. As an example, a virtual battery of 1000 thermostatically controlled loads is operated. In order to quantify the potential of virtual batteries, a study focused on residential thermostatically controlled loads in Spain is reported.


Introduction
A sustainable and environmentally friendly electric system requires the increase in generation coming from renewable energy sources, such as solar and wind energy. The increase in renewable energy sources in the electric system reduces emissions of greenhouse gases and pollutants. However, this type of energy is uncertain and intermittent, and therefore not completely predictable. Consequently, the task of balancing electricity demand and generation is becoming more complicated. As the penetration of generation based on renewable energy sources grows, more advanced and efficient regulation strategies are required to balance the power of the system. In this scenario, flexibility management (FM) could become a suitable costeffective alternative to overcome these technical difficulties.
Nowadays, an increasing number of houses have several electric appliances that work using their own thermal inertia or that of the home. These appliances are called "thermostatically controlled loads" (TCLs), and they are important FM resources. Their operation can be controlled in order to provide frequency regulation to the electrical grid they are connected to, while respecting at every moment their own physical constraints and the required user comfort. TCLs are also available in many industrial processes. Some examples are cold rooms and certain chemical reactors. In order to promote their operation as balancing resources for distribution system operators, TCLs can be aggregated into virtual batteries (VBs). As with conventional batteries, they are defined in terms of global quantities such as capacity, state of charge (SOC) and power, but with some differences, as is explained later in this paper. Several studies characterize VBs involving homogeneous or heterogeneous TCLs [1,2,3]. The last one is more realistic, has a greater significance and a broader application possibilities.
A great deal of research on TCLs and VBs has been carried out in recent years. Flexible loads, such as electrical vehicles (EVs) or TCLs, have been extensively studied as elements that could be aggregated to act as a virtual battery [4,5,6]. A theoretical characterization of the aggregate power and energy capacities for a collection of TCLs is reported in [3], as well as a priority-stack-based control strategy for frequency regulation using TCLs. A generalized battery model that can be directly controlled by an aggregator is introduced in [7]. The two abovementioned papers establish a TCL behaviour model and a VB control system. Another paper that experimentally models VBs is [8], where a first-order model is adjusted using binary search algorithms. The operating reserve hold by an aggregation of heterogeneous TCLs can also be modelled probabilistically, as explained in [9,10,11,12,13], using Markov chains, as in [14,15], or applying model predictive control [16,17]. Several of these papers do not take into account the real short-term instantaneous status of TCLs, and some of them are focused on changing the set-point temperature of the devices, which could imply reducing users comfort. In [18], the author proposes a stochastic battery model and a control system which performs TCLs monitoring by taking data at each instant, but it is not able to accurately fit the grid requirements. In [19], the technical viability of frequency regulation by TCLs is analysed, and it is proved that the load contribution can be significant. These thermal devices can also be the basis for a voltage control system [20]. TCLs are useful in microgrids, as shown in [21], and in virtual energy storage systems (VESS) [22], where they are able to provide ancillary services. Another use of TCLs is the management of intraday wholesale energy market prices [23].
Furthermore, the capabilities of TCLs have been studied in several regions such as Denmark, [24], Switzerland [25], Germany, [26], California [27] and Sardinia [28]. All these works show that TCLs can provide huge potential to mitigate the negative effects of the increase in renewable energy sources in the electrical system.
Two main contributions are reported in this paper. On the one hand, an accurate real-time VB control system is developed and its performance is tested by simulation. In the new control system, the grid operator determines the power required at every time instant based on the grid status. The resulting tracking order of error is of an individual TCL rated power, which is extremely low in comparison with other previously proposed methods. On the other hand, a study based on the residential TCLs population in Spain is carried out. The territory of this country can be classified into three different climate areas [29], and has a great renewable energy potential (high solar irradiance, strong wind in the central plateau and on the coasts, etc.). For these reasons, Spain exemplifies how VBs would help countries and regions to continue developing the implementation of this type of energy.
The remainder of this paper is organised as follows. In Section 2, models of an individual TCL and a VB obtained by aggregating a large number of TCLs are introduced. These models are used to develop a new VB control system that improves the accuracy in power delivery as compared to previously designed control systems. In Section 3, the VB and its control systems are used to study the capabilities of VBs as a cost-effective approach to balance power generation and consumption in electric systems with a large penetration of renewable energy sources. The study is carried out for residential TCLs in Spain. Section 4 shows the results obtained operating the VB control system and the results of the study in Spain, followed by a discussion. Finally, the conclusion is given in Section 5.

The Virtual Battery
A VB is defined as an aggregated collection of TCLs operated by a suitable control system to provide power and frequency regulation to the electrical grid. VBs are combined with adequate regulation policies to increase penetration of renewable energy generation sources in the grid by balancing power demand and generation in a cost-effective manner.
A new real-time VB control system, which accurately follows the operator signal, is presented in this section.  The control strategy complies with users and device constraints, such as TCL temperature bounds or short-cycling prevention. The goal is achieved by anticipating variations imposed by the abovementioned constraints. Our control strategy is an improvement of that reported in [30]. The models of TCL and VB used in our control system are explained below. They quantify the stored energy and the maximum power the VB may provide. The control methodology is also expounded and implemented in discrete-time for a certain sampling time denoted by h.

The TCL model
A model allows the controller to predict the behavior of each TCL and to estimate its stored energy. The parameters describing a TCL are collected in Table 1. The nameplate power P i is a positive quantity in cooling devices and negative in heating devices. The typical values of these parameters for residential TCLs are shown in Table 2. Reversible heat pumps are classified into two different types: 'reversible heat pumps (heat)' and 'reversible heat pumps (cold)'. The motivation for this classification is that some reversible heat pumps are only used for heating or cooling, not for both tasks. Besides, they have a different operation mode in the VB control algorithm.
Several discrete-time models of a TCL have been proposed in the literature [3,18]. Here, we use the model proposed in [18] as it is more realistic with a small cost of increasing complexity. In this model, the temperature θ i of each TCL i is calculated for every time instant k using the following equation where TCL Reversible heat pump (heat) 1.5 -2.  Each TCL has a set-point temperature (θ si ), which is set by the user, and a temperature dead-band (∆ i ) determined by the design of the TCL or by the user. Both parameters define the comfort band, whose bounds are θ si ± ∆ i . The temperature θ i must always be between these two bounds.
Consequently, the control system aims to maintain the temperature θ i in the comfort band. The average power to get the objective (P 0i ) is given by In order to maintain the temperature in the comfort band, each load can be switched on and off by the control system depending on the temperature and its own status.
Each device is characterized by four binary-valued variables: device type (φ i ), status (u i ), availability (δ i ) and total availability (γ i ). The device type φ i is 0 for a cooling TCL and 1 for a heating TCL. The device status u i is 1 when it is working, i.e., the motor is on with fixed power consumption P i , and 0 when off, with zero power consumption. The device availability δ i is 1 when the TCL is available for the control system and 0 otherwise. A TCL is available when θ i lies within the comfort band, and the time elapsed after its last change of status is large enough to avoid short-cycling. The full availability of the device γ i is 0 when the ambient forecast temperatureθ ai is lower (in cooling devices) or higher (in heating devices) than any TCL temperature bound. This last parameter indicates whether the TCL is ready to operate or not.
Short-cycling has harmful effects for electrical TCL systems. It causes damage in the electromechanical components, decreases their expected lifetime and can also significantly drive up the energy consumption. In order to avoid short-cycling, a minimum cycle elapsed time κ i is defined. Let ζ i be the elapsed time after the last status change, then the TCL device is available if The time to reach the temperature bound (λ i ) is the time that a TCL can be off without reaching that bound. It is calculated by simulating the evolution of θ i over time using the TCL dynamic model given by (1). If the TCL is not available, i.e. if δ i = 0, then λ i also equals 0.

Symbol Meaning
N Set of TCLs with cardinality N C c /C d Charging/Discharging capacity (kWh) SOC c /SOC d Charging/Discharging state of charge (kWh) n + /n − Maximum charging/discharging power (kW) n + /n − Maximum available charging/discharging power (kW)

The VB model
The VB model is an abstraction that represents a large number of TCLs with a reduced number of parameters. It allows the grid operator to make decisions in order to efficiently manage the grid. The set of parameters that describe a VB is given in Table 3.
The capacity of the VB is represented by two different parameters, C c and C d , to model the lack of symmetry in the dynamics of the charging and discharging processes of a TCL. In one case, the temperature change is forced by mechanical or electrical components (increasing θ in heating devices or reducing it in cooling ones), and in the other case it is not. Likewise, the state of charge (SOC c and SOC d ), the maximum power (n + and n − ), and the maximum available power (n + and n − ) are also duplicated.
The charging capacity (C c ) is the maximum energy that can be used from every TCL i ∈ N whose status u i is 1. If the VB is charging, then the charging state of charge (SOC c ) is the current energy reserve of the VB for charging. Likewise, the discharging capacity (C d ) and the discharging state of charge (SOC d ) represent the maximum energy that can be used for discharging and the current energy reserve for discharging. The capacities C c and C d are calculated by using the models to obtain the time that total available TCLs (i.e., {i ∈ N : γ i = 1}) take to evolve from one temperature bound to the other. The states of charge SOC c and SOC d also use the TLC models, but operating from the actual temperature θ i to the corresponding temperature bound. The algorithms developed for computing the capacities and states of charge require long-term temperature forecastsθ ai to produce accurate results. Accurate forecasts are assumed to be available in this paper.
The maximum charging power (n + ) is the maximum power a VB can provide for charging, assuming that all TCLs are available (δ i is 1): The maximum available charging power (n + ) is the maximum power that a VB can provide in the following instant, considering only the available TCLs. It is calculated by excluding the TCLs that are unavailable for charging: (n + ) k+1 = n k+1 where P + is the sum of power of the TCLs which are unavailable for charging. TCLs are unavailable for charging when their status cannot be on in the next time instant. The maximum available discharging power (n − ) is the same as n + , but for the discharging process: The maximum available discharging power (n − ) excludes the power of the TCLs which are unavailable where P − is the sum of power of TCLs that are unavailable for discharging. TCLs are unavailable for discharging when their status cannot be off in the following instant. The maximum powers n + , n − , n + and n − are calculated considering only TCLs that are totally available (γ i is 1).

The VB controller
The real-time VB control algorithm is composed of three blocks: Check of TCLs, Aggregation, and Priority Control. The first block is executed individually by each TCL, while the others are managed by a centralised controller. The communication between TCLs and the aggregator is done by using the method explained in [32], which requires an internet connection. An accurate long-term weather forecast is required to obtain reliable results. Figure 1 shows the complete controller structure where the interconnection of the three blocks is displayed. The set of variables that are used in the control algorithm are listed in Table 4.

Checking of TCL
The objective of this block is to compute the variables of a TCL at the next time instant and to avoid its temperature from leaving the comfort band. The process that a TCL i executes is detailed in Algorithm 1, see the Appendix. In this algorithm, a TCL with variables that are close to violating their constraints are identified and marked, which means that δ i and/or γ i become 1. The powers P + , P − and the extra power (P extra ) are calculated as well. P extra is the VB total power that will be System operator signal (kW) P agg Aggregated power (kW) Power switched (kW) switched on or off at the next time instant in order to avoid violating the TCLs constraints.
The computed variables of each TCL should be automatically communicated to the centralised controller every time instant h. If the length of the sampling time interval h is too small, then the communication system would require a high bandwidth and very expensive hardware due to the high volume of data. Notwithstanding, the variables must be sent as often as possible, so the controller may quickly correct any power deviation. In our studies, we have used h = 10 seconds.

Aggregation
In this controller block, the data obtained from each device are aggregated and compared to the system operator signal (r). This value is the amount of power that the grid operator wants to be provided. It can be either positive or negative, and has been previously determined by the grid operator using the current grid status and information about the VB capacity, state of charge and maximum power. If r is positive, then the VB is ordered to absorb energy. It could be caused by an excess in the renewable energy production or a strong power demand fall. However, if r is negative, then the VB has to stop consuming energy. This case corresponds to a higher power demand or a lower electric generation than expected.
System deviation (ψ) is the difference between the current power consumed by the VB, and its base power. The system deviation is obtained as where As a result of the aggregation process, the regulation signal ( ) is given by This variable will be used at the Priority Control block, and it determines the power that finally has to be provided. The aggregation block is implemented in Algorithm 2. See the Appendix.

Priority Control
The last block of the control system decides which TCLs change their status variables u i in such a way that the output VB power tracks the system operator signal r. The selected TCLs must be available and its temperature be far enough from the temperature bounds. The decision making process has been implemented in Algorithm 3, see the Appendix. The variable p tracks the total power obtained during the execution of the algorithm. The process stops when the absolute value of p exceeds the absolute value of or when every TCL has been assessed. When a TCL changes its status u i , then the cycle elapsed time ζ i is reset. The algorithm only works properly when the requirements of the grid operator are between the previously calculated capacities and the maximum available power. In that case, the maximum error between ψ and r is the greatest P i in the aggregation of TCLs.

Residential Virtual Battery Potential in Spain
The model of a VB and its control system developed in Section 2 are used here to perform studies about the implementation of residential VBs and their capabilities. The study focus on Spain for several reasons: its great potential for renewable energy penetration in the electric grid, the presence of different climate areas, and the availability of data. The main source of information about residential TCLs in Spain is the SECH-SPAHOUSEC project [29]. It was developed by IDAE (Instituto para la Diversificación y Ahorro de la Energía), and it estimates the penetration per house of heat pumps, cold pumps, electric water heaters and refrigerators in the different climate areas of the country in 2010. A summary of that information is presented in Table 5. According to [29], Spain is geographically classified in three different climate areas: North Atlantic, Continental and Mediterranean, see  Table 6, the variation in the number of houses between 2010 and 2019 in each climate area is considered. This information is obtained from the ECH survey (Encuesta Continua de Hogares) [33], developed by the Spanish INE (Instituto Nacional de Estadística).
In the North Atlantic climate area, the deployment of cold pumps and reversible heat pumps for air conditioning North Atlantic Continental Mediterranean Figure 2: Territorial distribution of climatic zones in Spain [29] is very reduced. The reason is that the temperature is not so hot in summer as in other climate areas.
Unfortunately, to the best of our knowledge, there does not exist publicly available databases iconcerning thermal and electrical characteristics of TCLs in Spain. Thus, the study has been performed by sampling random data from Table 2.

Results and discussion
In this section, an example of the VB controller operation and the residential VB charging capacity, discharging capacity, maximum charging power, and maximum discharging power of each Spain climate area are reported.

VB controller operation
The control system developed in Section 2 has been implemented in MATLAB [34] and its performance is studied for a VB of 1000 TCLs classified as follows: 125 reversible heat pumps (cold), 125 cold pumps, 125 reversible heat pumps (heat), 125 non-reversible heat pumps, 250 refrigerators and 250 electric water heaters. The parameters of each TCL have been randomly selected by using a Gaussian probability distribution centered on the mean values given in Table 2 and standard deviation 0.1. The simulation covers 200 time instants with a sampling interval h of 10 seconds, which amounts to a total of 2000 seconds. The minimum cycle elapsed time, κ, is given a value of 60 seconds. No disturbance has been considered, i.e. ω = 0. The initial status of each TCL has also been conveniently randomized using a binary probability distribution. The estimated ambient temperatureθ a is fixed at 20 • C in case of refrigerators and electric water heaters. For the remaining TCLs,θ a is variable. Figures 3-8 show the obtained results.
The system deviation ψ accurately tracks the system operator power signal r most of the time. The exception is the time interval between 510 and 750 seconds, where the maximum ψ stays at the value given by n + , the maximum       available charging power. The reason is that the grid operator established a value of r that the VB is not capable of producing as it is outside the range previously defined by n + and n − . This situation should be avoided by taking into account the capability, SOC and power information of the VB. When the signal r changes, several TCLs change their status, which means that they become unavailable for 1 minute. This explains the peaks and valleys in n + and n − when r is modified. The absolute error is lower than 6 kW in absolute value whenever r is between n + and n − . Otherwise, the VB cannot follow the system operator signal and the absolute error is either |r − n + | or |n − − r|. The charging/discharging capacities C c /C d and the maximum charging/discharging powers n + /n − appear to be constant in Figures 3, 6 and 7. In fact, they evolve with time but the changes are very small because the simulation time is not large enough to notice great changes in the ambient temperatures. The greatest values of these variables are given in Table 7. The charging capacity C d is about 6.5 times greater than the discharging capacity C c at every time instant. The reason is the lack of symmetry in the charging and discharging process, as the TCL takes more time to change its temperature θ when it is not forced to by electromechanical components, i.e., when it is switched off (u i = 0). Regarding the maximum charging/discharging power, n + is about twice as large as n − , since the sum of the nameplate powers P i doubles the sum of the average powers P 0i for every TCL i ∈ N . Figure 8 demonstrates the improvement achieved by the new controller proposed in this paper. In this figure, ψ is compared to ψ * between the seconds 270 and 510, which is the system deviation signal when the anticipation to the variations imposed by TCLs constraints is not included in the controller. As mentioned above, ψ follows r with an absolute error not greater than an individual TCL nameplate power; while ψ * is several times greater, as it depends on how many TCLs have to change their status forced by any of their inner constraints, which is almost random. Specifically, in this case, the maximum absolute error ψ * becomes around 15 times as higher as the maximum absolute error ψ, in absolute value.
Besides performing short-duration regulation, the VB is able to maintain a system deviation of -300 kW for more than 15 minutes, as can be seen from time 1010 to 2000 seconds. This means that this VB can provide ancillary services to the grid similar to the primary and secondary regulation, according to the current definition given by the Spanish normative [35], which is followed by the grid operator of the country, Red Eléctrica de España. Using European Network of Transmission System Operators (ENTSO-E) terminology [36], this VB could supply regulation similar to Frequency Containment Reserves (FCR) and Frequency Restoration Reserves (FRR) services. The main difference between the regulation provided by the VB and the primary and secondary regulation is that the first acts on the demand side, whereas the others act on the generation side. As a general rule, the VB will fail to follow r earlier if the amount of power requested is larger.

VB potential in Spain
The results of the study of the VB potential in Spain are presented here. Every existing TCL in each of the three climate areas of this country are assumed be part of the corresponding VB. The charging capacity, discharging capacity, maximum charging power, and maximum discharging power of each climate area of Spain throughout a natural year are shown in Figures 10-13. The legend and colour codes for these plots are given in Figure 9. The ratios of greatest capacity and maximum power per home are shown in Table 8. Information about the contribution of each residential type of TCL to the VB potential of each climate region is collected in Tables 9-12. The average values from Table 2 for the parameters of each type of TCL have been used in the study, along with the ambient temperature of each climate area, obtained from the hourly temperature model [37] of the most populated cities in each climate area during 2019: Bilbao (North Atlantic), Madrid (Continental) and Barcelona (Mediterranean). In the case of refrigerators and electric water heaters, the ambient temperature is assumed to be constant and equal to 20 • C.
Some interesting results should be highlighted. Mediterranean is the climate area with the highest VB potential, while North Atlantic has the lowest. Electric water heaters and refrigerators have the highest contribution to the charging and discharging capacity in the three climate areas. Moreover, this contribution is constant during the year, because the ambient temperature inside home is considered constant and equal to 20 • C. Refrigerators have more charging capacity in every climate area, while water heaters have more discharging capacity. These constant capacities can be used for demand management at any time of the year. The remaining devices have a variable capacity contribution, as can be easily seen by the numerous peaks and valleys appearing in Figures 10-13.        This is a consequence of the variability of the temperature along different hours, days and seasons. The capacity and maximum power of heat and cold pumps depend on the weather. In general, heat pumps can only be used on cold days and most of them are in winter. In contrast, cold pumps are mostly used on hot days, typically during summer. Reversible heat pumps and cold pumps have a greater contribution in Continental and Mediterranean areas than in the North Atlantic area. The reason is the greater number of them in these areas (see Table 6). Consequently, the variability in capacity and maximum power potential is more important in those climatic areas. Some summer days, the contribution of pumps to maximum charging and discharging power exceeds 50%. During these days, most of the cold pumps are working. Finally, the lack of symmetry in the charging and discharging process, explained in Section 4.1, clearly manifests here as a greater discharging capacity than charging capacity, as well as a greater maximum charging power than maximum discharging power, both for every climate area.

Discussion
A huge amount of energy and power flexibility of TCLs can be efficiently managed at real-time in a cost-effective way using a virtual battery. The control method proposed in this article to perform power regulation, control, and communication between TCLs and the aggregator improves the accuracy in tracking the system operator command power signal and does not require great investments in hardware or electronics because most of them are already installed, including the TCLs, thermal sensors and Internet connection. This is an efficient and cost-effective alternative to conventional batteries or fossil fuel solutions that require complete new installations. Furthermore, VB potential is expected to increase in the next few years due to the progressive electrification of heating devices.
Another important benefit that VBs provide is the decentralization and empowerment of consumers in the goal of balancing the grid. They become potential participants in balancing markets. Although these markets are not currently completely developed in Europe, the European Union has launched Regulation 2019/943 [38] and Directive 2019/944 [39], where balancing markets are considered. As mentioned in [38], these markets can either be individually accessed or be part of an aggregation, ensuring non-discriminatory access to all participants and respecting the need to adapt to the increasing share of variable generation and increased demand responsiveness. VBs fit properly with these requirements.
In Spain, net balancing energy amounted to approximately 1209 GWh in 2019 [40]. This means an average power of 138 MW. However, positive and negative hourly peaks of more than 4000 MW are produced [41]. The study performed in this paper shows that VBs clearly help to achieve power regulation goals in Spain, especially in extreme situations. The remuneration for participating in balancing markets must be important to encourage TCL owners to become contributors. Although the Mediterranean climate area has the highest residential power and capacity potential in Spain, VBs management should be profitable in any climate region if the economic profit is sufficient.

Conclusion
A VB real-time control system has been presented in this paper. Its operation shows low levels of error when VB constraints are fulfilled. The results prove that VBs can be very useful in providing regulation support to electrical grids. If the communication between TCLs and aggregators is fast, VBs can become an important short-term energy storage technology which, in combination with other existing ones (SMES, flywheels, etc.), provides renewable energy penetration. Furthermore, the capacity and maximum power of the aggregation of TCLs in the different climate areas of Spain has been computed. The obtained results show the potential of the different climate areas of the country. Virtual batteries are a powerful instrument for implementing demand side management programs. Their utility is closely related to the speed and reliability of communication networks, the quality of TCLs models and the ambient temperature prediction. Improving these aspects will be the main motivation for our future work. if φ i = 0 then 22: