1 Introduction

With the development of smart grid, demand side management (DSM) has been considered as an effective way to improve energy efficiency and optimize grid operation. DSM allows customers to adjust their energy consumption based on certain price signals or incentive schemes [1], so as to accomplish specific grid-level objectives, such as peak load reduction and frequency regulation. In summer, residential air-conditioning loads (ACLs) contribute to an increasing proportion of residential demand. To deal with the short term but sharp peak demand caused by the ACLs, enormous investments have to be made on upgrading the power infrastructure. In recent years, many DSM techniques have been proposed to re-shape the energy consumption profiles of residential appliance loads. These techniques can be categorized into two classes: the price-based DSM schemes (also known as indirect load control) [2–4] and incentive-based DSM schemes (or direct load control) [5, 6].

In addition to the DSM schemes, distributed renewable energy sources have also been widely deployed in distribution networks [7]. As one of the most promising renewable energy sources, wind power has gained widespread application in recent decades [8]. However, wind power generation is variable and stochastic in nature, which has negative impacts on the grid (e.g. power quality and protection settings). To satisfy the system operational constraints, in some cases wind power needs to be curtailed [9], which reduces wind power utilization to some extent.

In the literature, many efforts have been made to improve wind power penetration level while accomplishing peak load reduction in distribution networks. In [10], the authors presented an approach to minimize the expected operational cost of the microgrid and power losses while accommodating the intermittent nature of renewable energy resources. In [11], energy storage was utilized to increase wind power penetration level and curtail peak load with revenue optimization taken into consideration. In [12], the authors proposed a direct load control scheme for the ACLs based on the distributed imperialist competitive algorithm to minimize total operation cost over the whole dispatch horizon. In [13], the authors exploited the thermal characteristics of buildings to compensate for the renewable power fluctuations and reduce the renewable energy curtailment. In [14], the authors proposed a method to minimize the interconnection point power flow fluctuation by controlling the distributed controllable loads and renewable resources. In [15], the authors coordinated demand response programs to reduce conventional energy storage systems for the large scale integration of renewable energy resources. In [16], the authors developed model predictive control application for load management under high wind penetration, which effectively realized load shifting and real-time load balance. In [17], the authors reported a novel methodology of sizing energy storage systems to accommodate high penetration of variable energy resources to maintain power system stability.

By reviewing the literature, we find that many of the existing works focus on improving wind power penetration levels combined with optimal allocation of energy storage systems. Few works have attempted to establish an accurate ACL model to accommodate the intermittent renewable resources, especially small wind turbine resources at community group level, while taking the forecasting errors into account. Therefore, the major contribution of this paper is to propose a coordinated control scheme for wind generation units, battery energy storage systems (BESSs), and controllable ACLs to minimize system operation cost while reducing the system peak load. The proposed control scheme is modelled as a mix integer linear programming (MILP) problem, and rolling horizon optimization (RHO) is applied to alleviate initial prediction uncertainties. Meanwhile, in order to minimize customers’ thermal discomfort, an advanced two-parameter thermal inertia model is applied to more precisely model the thermal transition process of the indoor environment.

This paper is organized as follows. In Section 2, the modelling of thermal appliances is presented. In Section 3, MILP-based control model is proposed. The performance of proposed scheme is verified by simulation on a five node radial distribution network. In Section 4, the simulation studies are performed. Finally, Section 5 concludes this paper and future work is described.

2 Thermal modelling

2.1 Two-parameter thermal model

A full understanding of the dynamic thermal behavior of an air conditioned household is required if a DSM scheme is to maintain residents’ thermal comfort. Different thermal models have been developed and applied to gain information about the thermal process of a household and characterize its thermal energy consumption [18]. A one-parameter thermal model, shown in Fig. 1a, has been widely used in many previous smart home management studies [19, 20]. In this one-parameter model, the authors assume that the heat gain/loss of a house is only related to internal and ambient temperature difference, and the building envelope equivalent thermal resistance. It neglects any thermal capacitance effect associated with the walls. After turning on a room heating system (e.g. central heating), a thermal transient phenomenon occurs on the wall–room system, until it reaches a final equilibrium state [21].

Fig. 1
figure 1

Schematic of thermal models

Full size image

A more complex two-parameter thermal model considering the thermal capacitance of walls is presented in Fig. 1b [22]. The thermal process of a house is composed of two components. One component is the thermal mass inside the house and the other is the thermal mass of the walls with a notably different thermal capacity. According to the research in [23], the house indoor air temperature change can be significantly different when taking thermal capacity of the walls into account. Since the thermal models’ complexity can pose significant impact on the accuracy of calculating cooling energy, the more complex and accurate model in [22] is chosen here.

Compared with the thermal models in [19, 22], we innovatively improve the model performance by calculating air conditioner coefficient of performance (COP) shown in (5), which defines the ratio of cooling or heating energy supplied by the air conditioner over its consumed electrical energy. The efficiency of this novel two-parameter thermal model has been shown in our previous published works [12, 18, 24]. The dynamic two-parameter thermal model can be explained below:

$$\frac{{{\text{d}}T_{r} \left( t \right)}}{{{\text{d}}t}} = \frac{1}{{M_{a} \times Cp_{a} }}\left[ {\frac{{{\text{d}}Q_{gain\_a} \left( t \right)}}{{{\text{d}}t}} - \frac{{{\text{d}}Q_{ex\_w\_r} \left( t \right)}}{{{\text{d}}t}} - \frac{{{\text{d}}Q_{ac} \left( t \right)}}{{{\text{d}}t}}} \right]$$
(1)
$$\frac{{{\text{d}}T_{w} \left( t \right)}}{{{\text{d}}t}} = \frac{1}{{M_{w} \times Cp_{w} }}\left[ {\frac{{{\text{d}}Q_{gain\_w} \left( t \right)}}{{{\text{d}}t}} + \frac{{{\text{d}}Q_{ex\_w\_r} \left( t \right)}}{{{\text{d}}t}}} \right]$$
(2)
$$\frac{{{\text{d}}Q_{gain\_a} \left( t \right)}}{{{\text{d}}t}} = \frac{{T_{amb} - T_{r} }}{{R_{eq} }}$$
(3)
$$\frac{{{\text{d}}Q_{ex\_w\_r} \left( t \right)}}{{{\text{d}}t}} = \frac{{T_{w} - T_{r} }}{{R_{wr} }}$$
(4)
$$\frac{{{\text{d}}Q_{ac} \left( t \right)}}{{{\text{d}}t}} = COP \times P_{ac}$$
(5)
$$\frac{{{\text{d}}Q_{gain\_w} \left( t \right)}}{{{\text{d}}t}} = \frac{{T_{amb} - T_{w} }}{{R_{wa} }}$$
(6)

where T r (t) is indoor air temperature of the house at time t (°C); T w (t) is wall temperature of the house at time t (°C); Cp a and Cp w are heat capacity of the air and the wall (J/kg·°C); Q gain_a is heat gain which the ambient brings to indoor air (J); Q gain_w is heat gain which the ambient brings to the wall (J); Q ex_w_r is heat exchange between indoor air and indoor walls (J); Q ac is air conditioner cooling energy (J); M a and M w are the mass of indoor house air and walls(kg); T amb (t) is ambient temperature at time t (°C); R eq is house envelope equivalent thermal resistance; R wr is equivalent thermal resistance between house indoor air and wall inner surface; COP is coefficient of performance of air-conditioner; P ac is rated power of air conditioner (kW); R wa is equivalent thermal resistance between the ambient and wall outer surface.

Equations (1) and (2) refer to indoor air temperature and wall temperature change rate. Equations (3)–(6) represent the heat gain change rate that the ambient brings to the indoor air, heat gain change rate between indoor air and walls, air conditioners cooling energy change rate, and heat gain change rate that the ambient brings to the wall respectively.

According to [21], these parameters listed above can be easily estimated from the physical data of a building. Alternatively, it is also possible to use a parameter identification technique to obtain their values.

2.2 Thermal model linearization

To calculate indoor air and wall temperature variations more conveniently, the thermal dynamic model can be linearized. The whole operation time period t is split into N time steps, and t refers to 24 hour in this paper. Assuming that N is large enough, the ambient temperature, the wall temperature, and indoor air temperature are assumed constant within any time step [19]. Elaborate description about this linearization process can be found in [24]. The linearization is expressed in (7)–(10):

$$\begin{aligned} T_{r} \left( t \right) & = \left( {1 - \frac{1}{{M_{a} \times Cp_{a} \times R_{eq} }}} \right) \times T_{r\_init} + \frac{1}{{M_{a} \times Cp_{a} \times R_{eq} }} \times T_{amb\_init} \\ & \quad + \frac{{T_{w\_init} - T_{r\_init} }}{{M_{a} \times Cp_{a} \times R_{wr} }} - S_{ac\_init} \times \frac{{Q_{ac} }}{{M_{a} \times Cp_{a} }},\;t = 1 \\ \end{aligned}$$
(7)
$$\begin{aligned} T_{r} \left( t \right) & = \left( {1 - \frac{1}{{M_{a} \times Cp_{a} \times R_{eq} }}} \right) \times T_{r} \left( {t - 1} \right) + \frac{{T_{amb} \left( {t - 1} \right)}}{{M_{a} \times Cp_{a} \times R_{eq} }} \\ & \quad + \frac{{T_{w} \left( {t - 1} \right) - T_{r} \left( {t - 1} \right)}}{{M_{a} \times Cp_{a} \times R_{wr} }} \\&- S_{ac} (t) \times \frac{{Q_{ac} \left( {t - 1} \right)}}{{M_{a} \times Cp_{a} }},\;\forall t \in [2,N] \\ \end{aligned}$$
(8)
$$\begin{aligned}T_{w} \left( t \right) = T_{w\_init} + \frac{{T_{amb\_init} - T_{w\_init} }}{{M_{w} \times Cp_{w} \times R_{wa} }} \\+\frac{{T_{r\_init} - T_{w\_init} }}{{M_{w} \times Cp_{w} \times R_{wr} }},t = 1\end{aligned}$$
(9)
$$\begin{aligned} T_{w} \left( t \right) & = T_{w} \left( {t - 1} \right) + \frac{{T_{amb} \left( {t - 1} \right) - T_{w} \left( {t - 1} \right)}}{{M_{w} \times Cp_{w} \times R_{wa} }} \\ & \quad + \frac{{T_{r} \left( {t - 1} \right) - T_{w} \left( {t - 1} \right)}}{{M_{w} \times Cp_{w} \times R_{wr} }},\forall t \in [2,N] \\ \end{aligned}$$
(10)

where T r_init is indoor air temperature of the house at initial time (°C); T amb_init is ambient temperature at initial time (°C); \(T_{{w\_\text{in} it}}\) is wall temperature of the house at initial time (°C); S ac_init and S ac (t) are operation status of air conditioners at initial time and at time t.

Equations (7) and (9) are included for accurately calculating the indoor and wall temperature at the initial time step within the operating time horizon.

3 MILP-based rolling horizon optimization model

The participation of end users in demand response programs can effectively accommodate renewable energy integration in distribution systems. In general, loads within residential buildings can be classified in two categories: uncontrollable loads (e.g. lighting, refrigerator), and controllable loads (e.g., water heater, air conditioner). In this paper, we primarily consider air-conditioning loads in cooling mode as the controllable load. We study the control of aggregated ACLs through direct load control to achieve a certain grid-level objective, while ensuring that the room temperatures remain within the residents’ thermal comfort range. In this section, a MILP-based RHO model is proposed to control the ACLs to improve the wind power penetration level and reduce system peak load. Participating customers’ electricity cost is minimized on the premise of satisfying relevant distribution network device constraints and customers’ thermal comfort constraints.

Figure 2 depicts a schematic of the network configuration studied in this paper. In this network, the net feed-in-tariff for wind resources is lower than the electricity retail rate. Therefore, it is economically favorable for the system to consume local renewable energy as much as possible. Considering the variable nature of wind power, a BESS is directly charged by the wind turbines to smooth the power output available to end users. When the BESS is unable to fully supply the demand, customers need to buy electricity from the external grid to meet their electricity requirement. We assume that there is no direct power flow from the external grid to charge the battery. Conversely, the BESS could sell surplus power to the grid when customers are at a low demand level. The wind turbine and BESS are regulated by a micro source controller (MC). The ACL group is regulated by a controllable load controller (LC). MC and LC are further regulated by a microgrid control center (MGCC). For better depicting the power flow directions of system devices, the MC and LC are omitted in this figure but will be mentioned in the Case Study Section. Currently, we utilize an AC microgrid to organize wind turbine and BESS and we will design a DC microgrid where small wind turbines and BESSs have parallel connection [25] in future work to improve energy efficiency.

Fig. 2
figure 2

Simplified system network configuration

Full size image

3.1 Objective function

The objective function is represented as (11), aiming to minimize the system operation costs. The net cost in this distribution network is decided by the real-time power exchange between external grid and demand side, and real-time electricity retail price with net feed in tariff. In this case, the objective is formulated to minimize the total operation cost over all nodes and entire time steps.

$$\hbox{min} {\kern 1pt} {\kern 1pt} \sum\limits_{t = 1}^{N} {\sum\limits_{i = 1}^{N} {\left[ {P_{g\_buyi} \left( t \right) \times C_{buyi} \left( t \right) - P_{g\_selli} \left( t \right) \times C_{selli} \left( t \right)} \right] \cdot \tau } }$$
(11)

where P g_buyi (t) and P g_selli (t) represent the amount of power bought from external grid and sold to external grid at node i at time t; C buyi (t) and C selli (t) represent electricity purchase and selling price at node i at time t; τ is defined as the time step during which the parameters are assumed as constants.

It should be noted that the decision variables in (11) are P g_buyi (t) and P g_selli (t). These variables are time varying and determine the net cost in the distribution network over the whole time period.

3.2 Constraints

Objective (11) is subjected to following constraints.

3.2.1 Load balance constraints

The power flow between the distribution network and the external grid must match the local load balance at any time. This means the sum of the uncontrollable load and air-conditioning load is equal to BESS discharge power and power bought from external grid. Moreover, the amount of power bought from grid should not exceed grid maximum power capacity. Power sold to the grid should be less than the battery maximum discharge capacity. These constraints can be mathematically expressed as follows (12)–(14):

$$P_{uncontro\_load} \left( t \right) + P_{ac} \times S_{ac} \left( t \right) = P_{bat} \left( t \right) + P_{g\_buy} \left( t \right)$$
(12)
$$0 \le P_{g\_sell} \left( t \right) \le P_{bat}^{\hbox{max} \_discharge}$$
(13)
$$0 \le P_{g\_buy} \left( t \right) \le P_{g\hbox{max} }$$
(14)

where P uncontro_load (t) is uncontrollable load at time t (kW); P bat (t) is power flow from battery to consumers at time t(kW); P g_buy (t) is power bought from grid at time t (kW); P g_sell (t) is power sold to the grid at time t (kW); P max_discharge bat is maximum discharge power from battery (kW); P gmax is maximum grid power capacity (kW).

3.2.2 Wind power constraints

The output of the wind turbine largely depends on wind speed. In general, generated power increases as the cube of wind speed until it reaches rated speed. Beyond rated speed the wind turbine will generate a fixed amount of power while ever wind speed remains below the cut-off speed. At the other end of the scale, if wind speed is lower than cut-in speed or higher than cut-off speed, there would be no wind power generation. The wind power generation model is given in (15):

$$P_{w} \left( t \right) = \left\{ {\begin{array}{*{20}c} {\frac{1}{2}\rho \pi R^{2} v_{w}^{3} \left( t \right)c_{p} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} v_{c} \le v_{w} \left( t \right) \le v_{r} {\kern 1pt} {\kern 1pt} } \\ {P_{rated} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} v_{r} \le v_{w} \left( t \right) \le v_{f} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} } \\ {0{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} v_{w} \left( t \right) \le v_{c} \cup v_{w} \left( t \right) > v_{f} } \\ \end{array} } \right.$$
(15)

where P w is forecast wind turbine power (kW); ρ is outdoor air density (kg/m3); R is blade radius of wind turbine (m); c p is efficiency coefficient of wind turbine; v w , v c , v r and v f are wind speed, cut-in wind speed, rated wind speed and cut-off wind speed (m/s) respectively; P rated is rated power of wind turbine (kW).

3.2.3 BESS constraints

The state of charge (SoC) limits and rated power constraints of the BESS are considered below,

$$\begin{aligned} SoC\left( t \right) &= \frac{{\left[ {P_{w} \left( {t - 1} \right) - P_{bat} \left( {t - 1} \right) - P_{g\_sell} \left( {t - 1} \right)} \right] \cdot n_{effi} \cdot \tau }}{{E_{BESS} }} \hfill \\ &\quad+ SoC\left( {t - 1} \right) \hfill \\ \end{aligned}$$
(16)
$$0 \le SoC\left( t \right) \le 1$$
(17)
$$0 \le P_{bat} \left( t \right) \le P_{bat}^{\hbox{max} \_discharge}$$
(18)

where SoC(t) is battery state of charge level at time t (%); n effi is battery use efficiency during charge and discharge scheme; τ is time interval between time t and time t + 1 (hour); E BESS is battery capacity (kWh).

For simplicity, battery use efficiency n effi is assumed to be constant at 0.95 over the whole charging/discharging cycle. In a real scenario, n effi is closely related to battery SoC level. Due to the limited page length requirement, this discussion is beyond this paper’s scope.

3.2.4 ACL group thermal comfort constraints

The thermal comfort constraint for occupants is described in (19), where the lower and upper thermal comfort levels are presented in terms of temperatures. It is assumed that air-conditioners are operating at rated power once switched on. Air conditioners’ status is represented by the binary integer: when air conditioners are ON, S ac  = 1, otherwise, S ac  = 0. The lower and upper bounds for the walls temperature are formulated in (20):

$$T_{r}^{\hbox{min} } \le T_{r} \left( t \right) \le T_{r}^{\hbox{max} }$$
(19)
$$T_{w}^{\hbox{min} } \le T_{w} \left( t \right) \le T_{w}^{\hbox{max} }$$
(20)
$$S_{ac} \left( t \right) = \left\{ \begin{aligned} 1{\text{ ON}} \hfill \\ 0{\text{ OFF}} \hfill \\ \end{aligned} \right.$$
(21)

where T min r and T max r are lower and upper limit of indoor air temperature; T min w and T max w are lower and upper limit of wall temperature.

3.3 Rolling horizon optimization

Considering the uncertainties associated with electricity price, ambient temperature and wind speed prediction, a RHO strategy has been employed. By performing RHO strategy, model inputs are updated at each time step. Therefore, deviations from the initial forecast can be readily accounted for and the effects of prediction errors can be mitigated. The proposed MILP control model is implemented for a 24-hour time period, with time steps of 15 minutes. The main RHO procedures are given below.

  1. 1)

    At the 1st time step, the parameters including electricity price, ambient temperature, and wind power outputs are determined based on day-ahead forecast data. The MILP model calculates and generates a set of parameters (e.g., T r and T w ) based on the minimum operation cost objective.

  2. 2)

    At the next time step, cost function is optimized utilizing the MILP model based on updated input parameters, including real-time price (RTP), newly forecast ambient temperature, updated future wind power output, previous step generated data T r , etc., and generates a new set of parameters for next control window

  3. 3)

    At each time step, move forward the control windows, repeat the above procedure until the last time step of the planning horizon is finished.

Detailed procedures of the RHO strategy are shown in Fig. 3. It is assumed in this paper that new forecasts are available every 15 minutes in advance. The flow chart in Fig. 4 gives an overview of MILP-based RHO control scheme.

Fig. 3
figure 3

Representation of rolling horizon optimization process

Full size image
Fig. 4
figure 4

Flow chart of MILP-based rolling horizon optimization control scheme

Full size image

4 Case studies

The suburb of Randwick in Sydney has a 4.8 m/s average wind speed on site and is considered to be a good location to install small wind turbines for community group use. In 2010, the government installed a Skystream 2.4 kW wind turbine to educate the community about renewable energy alternatives and reducing carbon emissions [26]. It is expected more small wind turbines would be installed in this suburb to serve the community. Under this circumstance, this paper aims to provide pioneer research, for any interested stakeholder, on coordinating residential side controllable load and micro sources to improve renewable energy penetration levels and curtail peak load demand. The main objective of this section is to investigate the performance of the proposed control scheme under day-ahead control scenario and rolling horizon control scenario. Case studies are based on a five node radial distribution network with small wind turbines and Li-Ion batteries installed in the community group. The proposed MILP control scheme is conducted using MATLAB software combined with the Mosek tool box [27]. The simulation program is executed on a 4 core, 64-bit DELL Desktop with Intel Core i5-2400 CPU and RAM 4 Giga-byte.

4.1 Day ahead control

4.1.1 Experimental setup

  1. 1)

    The system network configuration is shown in Fig. 5. 2.4 kW Skystream wind turbines and 5 kWh Li-Ion batteries are installed near the community group households. As mentioned in Section 3, the capacity of small wind turbines and Li-Ion BESSs are aggregated to be controlled and regulated by micro source controllers (MC) at each community group. Air-conditioning controllable load in each group are gathered and managed in controllable load controller (LC). Local MC and LC are further regulated by decentralized micro grid control center (MGCC) for exchanging information with upper layer distribution management system (DMS). Total capacity of small wind turbines, BESSs, and ACL groups are given in Table 1.

    Fig. 5
    figure 5

    Configuration of radial distribution network with five nodes

    Full size image
    Table 1 Aggregated devices’ parameters on each node
    Full size table
  2. 2)

    Wind speed measurement is collected from day-ahead wind data [9]. Based on the wind speed, wind power output can be calculated by (15). Cut in wind speed, rated wind speed and cut off wind speed for Skystream wind turbine respectively are 3.5 m/s, 9 m/s and 25 m/s.

  3. 3)

    The five ACL groups are assumed to have distinct building environments and different thermal comfort requirements. Table 2 lists the information of the five ACL groups. In this research, the aggregation of thermal loads is simplified to treat households in each group with same building environment and thermal comfort preferences. Two methods can be used in future work to aggregate thermal loads. One is historical data based modelling [28]. It requires the specified type of household model to run for a long duration (e.g. a typical summer) so that seasonal nature of heating and cooling by householders can be captured. Another method is Monte Carlo simulation [29]. By assuming uncertainties in different households’ parameters (e.g. width, length and building materials) to closely simulate real-life conditions, Monte Carlo simulation is performed by repeated sampling of stochastic parameters to generate heterogeneous operating scenarios.

    Table 2 ACL groups’ information
    Full size table
  4. 4)

    The ambient temperature, electricity purchase price, and selling price for a typical summer day in Sydney [30, 31] are shown in Fig. 6, which is used for the forecasted day-ahead data in our research. It should be noted that only summer season data is used here for conveniently testing air conditioners operation status in a day. For real-time data, it is generated by forecasting methods. Here, we randomly generate three profiles of Gaussian noisy variables with small deviations and add them on the day-ahead temperature profile, electricity purchase price profile and wind speed profile respectively. The day-ahead data and real-time data used in the rolling horizon optimization are sent from DMS to MGCC, and from there relayed to local MC and LC. It is worth noting that selling price here refers to the net feed-in-tariff for wind resource. The electricity retail price is always higher than selling price. The peak electricity retail rate period is 17:00–21:00. Ambient temperature peak time in one day occurs from 13:00 to 16:00.

    Fig. 6
    figure 6

    Ambient temperature and day-ahead electricity price

    Full size image

4.1.2 Simulation results

Figures 7 and 8 indicate the scheduled air conditioner’s operation status and indoor air temperature in the five ACL groups during one day. Since ACL groups have different building characteristics and thermal comfort preferences, the air conditioners have different operation status in one day. As indicated in Fig. 8, the five groups’ room temperature is well controlled within the pre-set indoor temperature comfort zones. Peak load demand occurs during late afternoon when householders get back to home from work. As it is shown in Fig.7, air conditioners in five ACL groups are pre-opened to keep the house cool before peak load demand periods and shut down during peak load periods. By selectively turning off air conditioners at certain time, the peak load demand of this system is reduced to a lower level without compromising customers’ thermal comfort.

Fig. 7
figure 7

Air conditioners operation status in five groups

Full size image
Fig. 8
figure 8

Room temperature change in five groups

Full size image

Figure 9 outlines the power flow of the BESSs to consumers and corresponding SOC profiles. It is observed that the five BESSs discharge power at a relatively high rate during periods where electricity retail prices are high (6:00–10:00 and 17:00–22:00). BESSs are charged by wind turbines during low electricity retail price periods 1:00–6:00 and 11:00–16:00. It should be noted that wind power generation during these periods is completely stored by the BESSs without any curtailment. Therefore, the wind power utilization level is greatly improved by the proposed control method.

Fig. 9
figure 9

BESS operation conditions

Full size image

Figure 10 reports the power exchanges between the community group and external grid. According to the analysis in Section 3, the best operational mode for proposed system is to consume renewable generation directly from local micro sources. It can be clearly seen that the peak electricity purchase occurs in the period between 10:00–16:00, when the battery operates in charge mode. Due to high electricity demand between 17:00–20:00, local renewable resources are not enough to meet the requirements and are supplemented with the external grid. The simulation results imply that as long as the distributed renewable energy can be consumed locally through an appropriate control scheme, it won’t impact on system stability. On the contrary, it helps support grid.

Fig. 10
figure 10

Power bought from external grid in five groups

Full size image

The total calculation time for the benchmark system is 27.02 seconds, which is fast enough for 15 minutes time step control window. It should be noted that under the cost minimization objective, all the variables (e.g.S ac (t), T r (t) and T w (t)) are acquired by the Mosek tool box except the actual cost amount. In order to show the cost effectiveness of the proposed control scheme, we compare the system operation cost in our system with the situation that all the electricity is bought from external grid. By consuming the local power generation from wind turbines instead of buying electricity from external grid, the cost savings for one day in our system is $765.89, which accounts for 31.4% of original system cost $2436.6 (i.e., all the electricity is bought from external grid).

4.2 Rolling horizon optimization

Next, we validate the effectiveness of the RHO strategy on alleviating the day-ahead prediction errors. Based on the same network configuration in Fig. 5, continuous real-time information is updated and the control window is proceeded. The updates of the wind speed, ambient temperature, and electricity price on 15-minute basis are shown in Fig. 11. Two randomly selected groups are used to show the difference with day-ahead and real-time controls.

Fig. 11
figure 11

15-minute ahead data update

Full size image

As shown in Fig. 12a, air conditioners in group three and group five are switched on more frequently under the real-time scenario and operated for a longer time period. This is because the actual ambient temperature is slightly higher than the day-ahead forecast data, and longer operating time is therefore needed under this condition. It can be clearly seen that in Fig. 12b, the indoor temperatures are kept within the customers’ pre-set temperature limits under both scenarios, indicating that the customers’ thermal comforts are respected. The results also show that the system operation costs in the real-time control are higher than those in the day-ahead stage. This is because longer air conditioner operation time consumes more electricity from the external grid.

Fig. 12
figure 12

Air conditioners operation status and indoor temperature change

Full size image

By utilizing the forecasted day-ahead data and real-time data, we can compare the control outcomes under the two different scenarios. The result indicates the application of RHO scheme shows a more accurate control effect. With the wide deployment of bidirectional smart metering among end users, proposed methodology can effectively utilize real-time information of electricity pricing, forecasted weather condition, and residents’ thermal preferences to optimize the system operation.

To deal with the uncertainties brought by real-time data prediction, future work will compare another two existing techniques to validate the efficiency of the RHO scheme. One is a scenario-based approach [32], in which the possible uncertainties are realized by generating multiple scenarios. Scenario-generation methods include Monte Carlo sampling, moment matching principles and other methods motivated by stability analysis. The limitation of this method is the significant increase in the computation burden and scale. Another method is interval optimization [33], in which confidence intervals in terms of upper and lower bounds are used to represent the uncertainty spectrum. While this method does not need assumed probability distribution for uncertainties, it needs to carefully select the uncertainty intervals.

5 Conclusion

A MILP-based RHO control scheme is proposed in this paper. The objective of the proposed model is to improve wind power utilization and minimize the total operation costs. Interruptible ACLs are scheduled without significantly sacrificing the customers’ thermal comfort. In order to mitigate the negative effects brought by the uncertainties associated with some variables, the RHO strategy is employed.

The proposed control scheme is tested on a five node radial distribution network. Two cases in Section 4 validate the effectiveness of the proposed method. The simulation results prove that the wind power penetration level is guaranteed, and the system operation costs are also minimized while different operational constraints are satisfied. Meanwhile, the results show that the RHO strategy can effectively update the real-time information and re-schedule the resources based on the day-ahead decisions.

Future work would focus on the following aspects: 1) Consider the battery life influences caused by frequent charging and discharging of renewable power. It can be achieved by setting the battery operational constraint of depth of discharge (DOD) and locating in a suitable environment. 2) Design a DC microgrid for renewable energy and battery with parallel connection instead of current AC microgrid to improve energy use efficiency. 3) Employ historical data based modelling and Monte Carlo simulation to aggregate controllable thermal loads. 4) Compare the rolling horizon optimization scheme with scenario-based approach and interval optimization to validate its efficiency.