Flexibility Assessment of Multi-Energy Residential and Commercial Buildings

: With the growing concern about decreasing CO 2 emissions, renewable energy sources are being vastly integrated in the energy systems worldwide. This will bring new challenges to the network operators, which will need to ﬁnd sources of ﬂexibility to cope with the variable-output nature of these technologies. Demand response and multi-energy systems are being widely studied and considered as a promising solution to mitigate possible problems that may occur in the energy systems due to the large-scale integration of renewables. In this work, an optimal model to manage the resources and loads within residential and commercial buildings was developed, considering consumers preferences, electrical network restrictions and CO 2 emissions. The ﬂexibility that these buildings can provide was analyzed and quantiﬁed. Additionally, it was shown how this model can be used to solve technical problems in electrical networks, comparing the performance of two scenarios of ﬂexibility provision: ﬂexibility obtained only from electrical loads vs. ﬂexibility obtained from multi-energy loads. It was proved that multi-energy systems bring more options of ﬂexibility, as they can rely on non-electrical resources to supply the same energy needs and thus relieve the electrical network. It was also found that commercial buildings can offer more ﬂexibility during the day, while residential buildings can offer more during the morning and evening. Nonetheless, Multi-Energy System (MES) buildings end up having higher CO 2 emissions due to a higher consumption of natural gas.


Motivation and Aim
The growing worldwide concern with global warming is leading to the implementation of crosscutting policies to reduce greenhouse gas emissions, in particular in the transport and electricity and heat generation sectors, which produced two-thirds of the world's CO 2 emissions in 2016 [1].
The increasing pressure to reduce CO 2 emissions in these sectors is leading policy-makers to incentivize the adoption of emissions-free technologies, such as electrical vehicles (EVs) and renewable energy sources [2]. Despite the evident environmental benefits that these technologies yield, they also pose new challenges to network operators, which are facing higher uncertainties due to the volatile behavior of renewable energy technologies.
As a consequence, the network operation paradigm has been changing in the last decades, relying increasingly in flexibility provided by generation and demand side assets. As the integration of renewables is expected to boost in the coming years, it is of utmost importance to find new ways of increasing energy systems flexibility.

Contributions and Advantages of the Proposed Model
This paper proposes an innovative approach to optimize the flexibility provided by multi-energy residential and commercial buildings with the aim of solving technical problems in electricity networks operation. The model generally contributes to increase networks robustness and may be used by system operators to increase the quality of service indexes. It aims at minimizing costs, considering the prices of energy, the costs of using DR programs and the costs of CO 2 emissions, while keeping the network operating within the specified technical limits.
The proposed model considers different sources of flexibility from residential and commercial buildings. This flexibility comes from changing the flexible load patterns or by managing the resources within the buildings. The flexible loads include thermal loads, like space heating and cooling, water heating or refrigeration systems, and also other loads, like EVs or lighting. Space heating and cooling is assumed to be controlled through temperature set-points defined for centralized HVAC systems and thermostats. Although only some type of resources are used in this work, such as electric boilers (EB), gas boilers (GB), CHPs, air conditioners (AC) and heat pumps (HP), the proposed model is prepared to include other resources like P2G or combined cooling, heat and power (CCHP). The studies performed look forward to a future where these technologies are widely used, since although these resources and their combinations are not extensively used nowadays, the situation may change in the coming years as the costs of these technologies have been showing a clear decreasing trend.
The optimization problem used assumes that consumers provide the parameters necessary for its resolution, like temperature ranges, lighting levels or EVs charging needs. The optimization problem differs from [7][8][9][10]14,[16][17][18][19][20][21][22][23] by allowing the participation of different sources of load flexibility (e.g., PVs, EVs, HVAC systems, thermal loads, lighting) in DR programs and by considering the constraints imposed by the consumers needs and by the technical characteristics of the electrical network.
In the literature, the models developed for HVAC systems considering thermal models are usually non-linear models [24,25] or quadratic [26]. In order to linearize these models, some works use a sequential linearization approach [17] or assume a linear model with parameters calculated using regression methods [16]. A new linear model for HVAC systems was developed in this work. This linear model differs from the other linearized models by discretizing the air flow rate supply variable. This way, the approach for the linearization of the models is simpler as it does not require to be done individually for each building.
With this model, the authors studied and compared the flexibility that can be provided by buildings to electrical network operators considering two scenarios: the buildings only have electrical loads (only-electrical buildings) and the buildings have loads that consume different types of energy (MES buildings). In the literature there are several studies regarding the advantages of MES, but this study specifically compares the same buildings with two different structures: one purely electrified and other with MES. Electrification and MES are two different perspectives of looking into future scenarios of the power systems development. They have aroused great interest in the scientific community and this study offers a new level of comparison of both of them.
Furthermore, the flexibility offered by commercial and residential buildings is studied and compared. There are several works in the literature that study the flexibility that commercial [16,17] and residential [18,19] buildings are able to offer to network operators or aggregators but these works only focus on each type of buildings individually. Commercial and residential buildings have different patterns of consumption which may reflect the availability of flexibility in different times of the day. In this work, it was considered the operation of both type of buildings and the flexibility that each building is able to provide was compared.
The proposed approach was applied to a case study that includes an electrical network with technical problems. The results achieved show how the model manages the flexibility from the resources available to overcome the identified technical problems. In addition to the global flexibility quantification, the contribution of different flexible resources was also analyzed. The importance of a detailed model to explore every source of flexibility within buildings in order to mitigate any problem present in power networks is thus highlighted.
Lastly, an annual analysis of the different cases was also performed with emphasis on the costs, energy consumed and CO 2 emissions.

Paper Organization
The rest of the paper is organized as follows: Section 2 explains the problem and the model developed in this study; in Section 3, the mathematical model of the proposed approach is presented; Section 4 provides the case studies and numerical results; finally, Section 5 presents the main findings of the work developed.

Problem and Model Description
This work focus on the flexibility that can be provided by residential and commercial buildings to a system operator. This flexibility is obtained by changing the resources' point of operation. The buildings' flexibility will be quantified and analyzed for two scenarios: when they have resources that only use electricity and when they have multi-energy resources. The model developed allows the system operator to optimize the flexibility of the buildings' resources through a DR-based program in order to solve network technical problems. Therefore, the model assumes that the system operator is able to control the resources and loads of the EHs by communicating with their EH management system.

System Operator
The system operator is responsible to satisfy the energy needs of the consumers, while maintaining the network operating within the predefined technical limits. When the network is operating in normal conditions, systems operators may use the clients' flexibility to minimize operating costs (OPEX) or to maximize RES integration. When technical problems occur, flexibility may be used to solve them. Of course that consumers should be paid for the flexibility provided. The cost of this DR-based service is assumed to be established in a contract set between the system operator and the consumer.
The consumers present in the network are commercial and residential buildings. They are represented by EHs and it is assumed that they have an EHMS integrated. The EHMS is responsible for the control of all the assets within the EH, to fulfill the consumers' load requirements and optimize the EH. When needed, the system operator communicates with the EHMSs to request the delivery of flexibility.
EHs are connected to the electrical network in a certain point and they consume or inject electricity in that same point. This way, the system operator is able to manage the consumption or injection of electricity in that point by by using the flexibility of the EHs. When needed, the system operator orders the resources to decrease or increase the consumption of one energy vector or change the load patterns through the DR program, always taking into account the limits predefined by the client.

Energy Hub
EHs represent commercial and residential buildings and they are formed by input and output energy points, conversion technologies, storage systems and different types of loads ( Figure 1).
The input energy points considered are the electricity consumed from the electrical network and from local PVs and the gas consumed from the gas network. The output point is the electricity injected into the electrical network.
Conversion technologies can transform one type of energy to another and storage systems can store energy to be used in another time.
The loads to be satisfied include electricity, gas, heating and cooling loads. They can be inflexible (e.g., electrical equipment and cookers) or flexible (lighting, space heating and cooling, water heating, refrigeration and EVs). Inflexible loads have to be completely satisfied while flexible loads can be modified, according to the flexible limits imposed by the consumers, in order to optimize the EH.

Energy Hub Management System
As previously explained, an EHMS is responsible for controlling the assets of the EH in order to fulfill the requirements of the consumers. They are also responsible to communicate with the system operator to inform their energy needs and their availability to participate in the DR program. These requirements are initially set in the EHMS and they involve defining different restrictions and flexibility levels for each type of load. The definition of requirements for each type of load, resources and storage system is explained in the following paragraphs. A scheme of the EHMS is presented in Figure 2. This concept was adapted from [27] and extended for the application in an EH.
The space heating and cooling temperature requirements for each space are initially communicated by consumers to the EHMS and the actual temperature of each room is calculated by the EHMS using a thermal model. The HVAC systems control the heating and cooling power of commercial buildings and thermostats control the heating and cooling power of residential buildings. These systems are responsible for activating the necessary resources to provide the cooling and heating power according to the temperature requirements of the consumers and the actual room temperature calculated by the EHMS.
The consumers communicate the temperature for water heating and refrigeration and the range of acceptable temperatures (i.e., their temperature flexibility) to the EHMS. Water heating and refrigeration loads are also defined by thermodynamic models and calculated by the EHMS. Electric equipment, lighting and gas cooking loads are defined according to the profiles sent by the consumers. Lighting in commercial buildings has a flexibility range that is also defined by the clients.
The restrictions of EVs are defined by the consumers and they involve defining the hours when the EVs will connect and disconnect from the EH and the state-of-charge (SOC) desired at the moment of disconnection.
Storage systems, PVs and other conversion resources have no requirements and can be managed as better fits the EHMS. The only limits imposed are related with their technical characteristics. In summary, the inputs that the consumers need to provide to configure the EHMS are the following: • Spaces temperature set points and flexibility bands; These inputs need to be provided only once. After that, the consumer may provide new inputs if he wants to reconfigure the EHMS.
With all this information, the EHMS will calculate the energy required by the EH and manage the resources, storage systems and flexible loads in order to satisfy all the requirements imposed by the consumers.

Model Constraints
The model developed in this work and presented in this section has the goal of optimizing the costs of the resources present in the network in order to fulfill all technical requirements from the network and from the resources.

Objective Function
The objective function in Equation (1) considers the costs of buying electricity as seen in Equation (2) and gas in Equation (3), the costs of DR in Equation (4), the profit of selling electricity in Equation (5), the costs of CO 2 emissions (6) and penalization costs for voltage violations in Equation (7). The costs of DR consider the modification of the profile of lighting and temperature of rooms, water and refrigeration systems. ).π DR t , ∀j ∈ {EB} (4) The objective function is subjected to operational network constraints which are presented in the following subsections.

Network Constraints
In this section, the operation of an electrical network is presented. Some of the variables presented are linked with the variables of the EH modeling in Section 3.3, as it will be discussed.

Power Flow Constraints
The considered network is a radial distribution network. A branch power flow model (DistFlow) [28,29] was considered and modeled through Equations (8)- (13).
The mathematical formulation presented above is non-linear and non-convex. This formulation was simplified using a linear version of the DistFlow, the LinDistFlow [30,31] in Equations (14)- (17). In this linear approximation, the branch loss terms (r 2 i,j + x 2 i,j ).I 2 i,j were neglected. This is possible if it is assumed that the line losses are much smaller than the branch power terms P i,j and Q i,j . Voltages are also nearly balanced and (V i − V 0 ) 2 ≈ 0. These assumptions will lead to the expression The substitution of this term in Equation (10) leads to Equation (16). Differently from the DC-power flow, this model considers lines resistance and models active and reactive power flow. Variable V viol j represents a voltage violation and it is used to allow the problem to converge even when a solution with all the voltages inside the allowable limits (lower limit-V j ; upper limit-V j ) does not exist. This value has a high penalization in the objective function.
The power flow was calculated for each instant t. It should be mentioned that this model can only be applied in radial networks.
Equation (18) defines that the power P j in each bus is related with the power of all EHs in that bus. The calculation of the power I w n,t for each EH is presented in Section 3.3.

Generator Units Constraints
Equations (19) and (20) define the limits of PV and wind generators.

Energy Hub Modelling Constraints
In this section, the operation for an EH is presented. Some of the variables are linked with the variables of the network modelling in Section 3.2 and of the resources modelling in Section 3.4, as it will be explained.

Load
Equations (21)-(24) represent the consumption of each type of load. The load of each energy vector is equal to the sum of the energy coming from the energy network, storage systems or the converters i generating that energy vector.
The electrical load L w t is equal to the sum of the electricity needed for EVs, HVAC systems, refrigeration systems, electrical equipment and lighting, as seen in Equation (25). The heating load L h t is equal to the sum of the heat needed to supply the heating for spaces and water as seen in Equation (26). The cooling load L c t is equal to the sum of the cooling needs for space cooling as seen in Equation (27). The gas load L g t is equal to the sum of gas needed for cooking as seen in Equation (28).

Energy Hub Inputs
Equation (29)-(30) represents the energy inputs of each EH. Each energy input is equal to the sum of the energy that goes from the energy network to the storage systems, loads and converters i that consume that energy vector.

Energy Hub Outputs
Equations (31) represents the electrical output of each EH, i.e., the electricity sold to the network. The energy output is equal to the sum of the energy that goes from the EVs, storage systems or converters i generating electricity to the electrical network.

Energy Converters and Storage Systems Input
Equations (32)-(35) represent the energy inputs of each type of storage system or converter inside the EH. The energy input of each system j is equal to the sum of the energy that goes from the network, storage systems or converter i generating that energy vector to that system.
Equation (36) represent the limits of the energy inputs of each storage system or converter j.

Energy Resources Output
Equations (37) to (39) represent the output of the storage systems or resources present in the EH. The output of each system is the sum of energy that goes from that system i to the network, storage system, load or resource j that consumes that energy vector. It should be noted that the electricity sold to the electrical network is considered in Equation (37).

Load, EVs, Storage and PV Contraints
In this section the constraints of the different loads, EVs, storage systems and PVs are presented. As mentioned, some of the variables are linked with the variables from the energy huh modelling in Section 3.3. In Table A1 of Appendix A, it is synthesized the values of the main parameters used in the models presented in this section.

Thermal Model
The temperature in a room is calculated by Equation (40) [19]. Equation (41) set the temperature at which the room has to be, taking into account the temperature set by the consumer Θ room i,t and the temperature changes due to the participation in the DR program. Equations (42) and (43)  of the temperature in a room for the participation in the DR program.
Equation (44) defines the infiltration rate present in the rooms [32].
In order to calculate the energy necessary for the cooling and heating needs, it is necessary to define first the thermal parameters for each room. Equation (45) defines the thermal resistance [33,34].
The thermal reactance is calculated by using the heat capacity from the air in the room and from the walls. For the sake of simplicity, the outside walls and windows parameters were not considered in this study. The heat capacity contribution of the materials inside the room were not considered as well. Equation (46) defines the thermal reactance [35].

HVAC System
Commercial buildings have an HVAC system installed. Their typical configuration is represented in Figure 3. The cooling coil provides the cooling power and the heating coils in each room provide the heating power. There is a supply fan which moves the cooling or heating air throughout the entire system. For each room there is also an individual variable air volume (VAV) box connected to the HVAC system, which regulates the temperature inside a room. Equations (47)-(50) define the cooling and heating power provided by the respective coils [16].
The coils are air-water heat exchangers and Equations (51) and (52) define the power needed to power the cooling and heating coils [26], respectively, from the heating and cooling resources of the EH.
The relationship between the fan power and the airflow is between a quadratic and a cubic function, depending on the fan system [33]. This behavior is explained by the fan needing more power to compensate the leakage in the system, when an increase in airflow is necessary. In this work, it was assumed that the fan system model used is the same as the one presented in [16], as it is a model for commercial buildings. This model is defined in Equation (53) and it is applied to the commercial buildings considered in this work. The values of b 1 , b 2 , b 3 and b 4 are taken from the same reference. Nonetheless, if the required data was available, they could be identified using several techniques, including the ordinary least squares technique [25].
The total airflow supplied by the fan system is equal to the summation of the individual airflow of each room as seen in Equation (54).
The value of η h and η c is 0.98 and the value of COP h and COP c is 5.92. The discharge air cooling temperature Θ C,c t is set to 12.5 • C and the discharge air heating temperature Θ C,h t is set to 50 • C. The values to calculate the fan power are the following: b 1 = 2.57 × 10 −12 , b 2 = −4.45 × 10 −9 , b 3 = 1.46 × 10 −4 and b 4 = 4.71 × 10 −3 . Each fan has a limit of 20,000 cfm. Equation (49) present a non-linearity between two continuous variables m i,t and Θ room i,t . Variable m i,t was discretized by assuming binary variables b f an j,i,t that are multiplied by a small value of 100 cfm. This way, the problem can optimize the supply air flow rate with a rate of 100 cfm. By assuming a maximum of 20,000 cfm for a fan, there would exist 200 new binary variables. This would increase the complexity of the problem. Thus, instead of increasing the values by 100 each time, the supply air flow rate can be increase by 100, 200, 300, 400, 2000, 3000, 4000 and 10,000 cfm. This way, the problem can still look for solutions between 0 and 20000 cfm while needing only 8 new binary variables. This linearity is present in Equation (55).
Nonetheless, with the discretization of m i,t there is still a non-linearity due to the multiplication of the binary variables b n i,t and Θ room i,t . This non-linearity can be solved by replacing that multiplication by a new continuous variable y i,t , as represented in Equation (56), and by adding the restrictions represented by Equations (57)-(60).
With the discretization of m i,t , Equation (53) also presents several non-linearities due to the multiplication of several binary variables. It was considered an additional binary variable y that replaces the multiplication of the two binary variables. In Equations (61)-(63) it is presented a general approach for the linearization of all these non-linearities.

Thermostats
Contrary to the HVAC systems, thermostats used in residential buildings control directly the heating and cooling resources according to the thermal model defined in Equation (40). This way, P c it and P h i,t are the power needs that are requested from the heating and cooling converters as seen in Equations (64) and (65). P c,resources

Water Heating
The power necessary to heat the water to the required temperature is presented in Equation (66) [20]. The required temperature of the water is calculated according to the temperature defined by the consumer Θ water i,t and the changes from the DR program as seen in Equation (67).
The consumers also define the temperature limits of the DR program Θ DR,water,− i,t+1 and Θ DR,water,+ i,t+1 , as present in Equations (68) and (69).
The value of c w is 4.18 J/g· • C (considering water at 18 C • ) and initial temperature Θ water,i j is set to 15 • C in Winter. The mass of water m i is defined according to the water consumption profile.

Refrigeration
The model used in this work to represent a refrigerator is represented in Equations (70) and (71) [22]. The required temperature of the refrigeration system is calculated according to the temperature defined by the consumer Θ  (73) and (74).
The coefficient of performance COP re f r was assumed to be 3 and the insulation a was assumed to be 15 kJ/ • C for residential systems and 850 kJ/ • C for supermarket systems.

Other Loads
The power necessary to supply energy to cooking and other electrical equipment, according to their load profiles, is defined in Equations (75) and (76).
The power necessary to supply energy for lighting is defined in Equation (77), which already takes into account the changes due to the participation in the DR program. This participation is constrained by the limits imposed by the consumers P   This network has several loads, including residential buildings, commercial buildings, an hospital, warehouses, medium and small offices and a supermarket. The location and number of buildings in the network is presented in Table 1. This table also presents the type of loads of these buildings, besides space heating and cooling, lighting and electrical equipment. The supermarket has 10 refrigeration systems, while the residential buildings only have 1 system each. Refrigeration systems have a maximum power of 5 kW for supermarkets and 1 kW for residential buildings.

Case Study Description
These buildings are represented by the models shown in Sections 2 and 3. They are supplied by the electrical network, the gas grid and by local PVs. The technical restrictions of the gas grid were not considered in this study. The resources assumed for each case, i.e., only electrical buildings and MES buildings, are presented in Table 2. All the buildings also have heat storage systems installed. The efficiency of EB is 0.95, GB is 0.9 and the efficiency of CHP for electricity is 0.35 and for heat is 0.45. The efficiency of heating storage systems is 0.88. The COP of both HP and AC systems is 3.45.

Buildings
Electrical In every building, 10% of the occupants have an EV. EVs have a charging and discharging power limit of 10 kW and an efficiency of 0.9. The battery capacity of these vehicles varies between 35 kWh and 100 kWh.

MES EB HP AC GB CHP HP AC
PVs are installed in every building: 10 kW in small offices, 50 kW in medium offices, 20 kW in warehouses, 20 kW in the supermarket and 200 kW in the hospital. The systems' efficiency is 0.227.
The parameters and profiles of each building are taken from [36]. The assumptions for the parameters of the models were already presented in Section 3.
The wind and PV generation profiles are presented in Figure 5 and they were taken from [38], for the day 7 January 2019. The outside temperature is presented in Figure 6 and was taken from [38] for the same day. The average monthly wind and PV generation profiles are presented in Figures A1 and A2 in the Appendix B and they were calculated with the information taken from [38] for the years 2017 and 2018. The average monthly outside temperature is presented in Figure A3    In Figure 7, the prices of electricity and gas are shown. They were taken from ENTSO-E [39] and from the MIBGAS platform [40], respectively, for the day January 7th of 2019. The price of the electricity sold to the grid was assumed to be 0.5 of the electricity price in order to enforce self-consumption and avoid injection into the grid. The average monthly electricity and gas prices are presented in Figures A4 and A5 in the Appendix B and they were calculated with the information from [39] and [40], respectively, for the years 2017 and 2018. The price of DR services related with temperatures flexibility was assumed to be 0.5 e/ • C. This figure was assumed for the purpose of this work, as no supporting reference was found in the literature. Although it is outside the scope of this work, it is a very interesting topic that should be further studied in a dedicated work.
The price of each voltage violations was assumed to be very high to force the optimisation model to give priority to the resolution of these technical problems. This value was set to 10,000,000 e.
The cost of the CO 2 emissions from the natural gas is calculated by considering the average market price of CO 2 emissions and a factor representing the total emissions per MWh of natural gas consumed. The CO 2 emissions market price considered was 25 e/tonCO 2 for every month and it was taken from [41], whereas the emission factor was set to 0.2 tonCO 2 /MWh and it was taken from [42]. This means that the CO 2 emissions cost for natural gas was set to 5 e/MWh for the entire year.
The cost of the CO 2 emissions from electricity is calculated as follows: • First, the percentage of coal and natural gas in the electricity generation mix is determined, defining a daily profile for each of these resources; • Then, the daily profiles are multiplied by the CO 2 emissions market price for that day and by a factor representing the total emissions of each resource per MWh of generated electricity.
The daily profiles of coal and natural gas in the electricity generation mix were calculated with the information from [43] for the years 2016, 2017 and 2018 and are presented in Figures A6 and A7 of Appendix B. The emission factors per MWh of electricity generated for coal and natural gas were set to 1 tonCO 2 /MWh and 0.41 tonCO 2 /MWh and they were taken from [44]. The average monthly cost of the CO 2 emissions from the electricity generated is presented in Figure A8 of Appendix B.
The cost of the CO 2 emissions, for both electricity and natural gas, will not have a significant impact in the results as they are much lower than electricity and gas prices.
The optimization problem described in the previous sections is a mixed-integer linear programming problem (MILP) and was solved with pycharm [45], using the CPLEX optimizer in a core i7, 3.5 GHz PC.

Results
In this section it is studied how the flexibility provided by the different buildings can be used to solve voltage problems that may appear in the network. To this end, the following cases were created:

•
Only electrical buildings without optimization (no DR); • Only electrical buildings with optimization (with DR); • MES buildings without optimization (no DR); • MES buildings with optimization (with DR).
With these cases it will be possible to study and compare the flexibility provided by electrical and MES buildings and by residential and commercial buildings, as well as their annual costs, energy consumption and CO 2 emissions.
The organization of this section is as follows:

Grid Technical Problems
In this section it is presented the technical results of the network on the 7 January 2019. The voltages in all buses at peak hour and the voltages at bus 20 during all day are presented and analyzed. In Figure 8 it is presented the voltages in every bus of the network at 6h. It is possible to observe that buses 18, 19, 20 and 21 present undervoltage problems in the scenario with "electrical buildings without optimization", while in the scenario with "electrical buildings with optimization" only bus 20 presents a slight undervoltage problem. The other scenarios do not present this type of problems. These problems occur due to the fact that these buses have high consumption of electricity at this hour, as it will be analyzed in Section 4.2.2, and also because they are the farthest from the feeding point, as it is possible to observe in Figure 4.  Figure 9 shows the daily profile of the voltages in bus 20. It is possible to observe that severe voltage problems occur in the "electrical buildings without optimization" scenario. Buses 18, 19 and 21 also reveal similar undervoltage problems. The importance of DR programs is easily perceptible in the "electrical buildings with optimization" scenario, as almost all undervoltage problems cease to exist. Nevertheless, there are still some minor problems in buses 20 and 21. In both MES scenarios, with and without DR, no voltage problems occur. Even without using the DR program, the system was never with technical problems as part of the MES buildings consumption relies on other energy vectors besides electricity.
In total, in the scenario with electrical buildings without optimization, there were 54 cases of undervoltage affecting 4 buses, while in the scenario with electrical buildings with optimization there were 14 cases affecting 2 buses. This means that flexibility activation through DR programs reduce the number of undervoltage problems by 74%.

Daily Load Diagrams
In this section it is presented the electricity and gas consumed in the network during the day of 7 January 2019. The changes in the load due to the participation in DR programs are analyzed.
In Figure 10 it is presented the electrical consumption in the network. It is possible to observe that without optimization, the case with "only electrical" buildings consumes more electricity than the case with MES buildings. The peak consumption at 6h in the scenario with electrical buildings without optimization is due to the fact that electricity prices are cheaper at this hour and clients start heating the spaces, especially in commercial buildings. This higher consumption of electricity causes several voltage problems in the network, as it was seen in the previous section. The utilization of flexibility to solve the referred problems impose changes to the daily load profiles. It is possible to observe that the main changes occur in the scenario with only electrical buildings with optimization as this case presents several voltage problems during the entire day. This way, the electrical load was decreased in this case.
In the cases without optimization, the peak load is 85% lower in the scenario MES buildings in comparison with the scenario only electrical buildings; comparing the same cases with optimization, it is observable that the peak load is reduced by 67%. Figure 11 shows the overall gas consumption in the network. The peak consumption is at 7 h as most of the buildings with gas also start heating their spaces at this hour.

Flexibility of Only Electrical and MES Buildings
In this section, the flexibility that electrical and MES buildings are able to offer is analyzed and quantified. The minimum and maximum flexibility represents the minimum and maximum consumption that each building is able to consume at each hour. Figure 12 presents the minimum and maximum flexibility (red lines) that the buildings are capable of providing in the only electrical and MES scenarios. Eventually, this flexibility may be activated if proper DR-based programs are implemented. The most important conclusion to retain is that buildings have the capability to provide a significant band of flexibility that may be very important to system operators, as they can use it to solve low voltage and/or branches' overload problems in distribution networks, as it will be seen in Section 4.2.1.
Although the electrical buildings scenario have a larger flexibility band, their minimum limit is slightly higher as all the loads rely entirely on electricity. In the MES scenario, the flexibility band is generally straighter since there are less loads that use electricity. However, in this case, the overall buildings' consumption is lower than in the only electrical scenario, what may be an interesting situation for system operators since networks utilization will be lower.
It is possible to observe that only electrical buildings have major modifications in their profiles during the first 7 h of the day. This happens to avoid higher peaks of consumption at 7 h, when buildings start heating their spaces. Buildings accumulate heat in the previous hours decreasing their energy needs at 7 h. This is also shown in Section 4.2.5.
In average, the total minimum consumption of MES buildings is 40% lower than electrical buildings while the maximum consumption is 37% lower.

Flexibility of Residential and Commercial Buildings
In this section, the flexibility that residential and commercial buildings are able to offer is analyzed and quantified.
In Figure 13 it is presented the range of flexibility of residential and commercial buildings. It is possible to observe that commercial buildings have a higher range of flexibility than residential buildings. They have a higher flexibility during the day as they are intensively used during these periods, with a higher number of appliances and loads connected and capable of increasing or decreasing their energy consumption. Air heating/cooling and water heating are good examples of loads that are a great source of flexibility in these periods.
Residential buildings have higher flexibility after 16h and during the night. This happens since residential buildings do not usually have EVs connected between 8 h-16 h, which are a significant source of flexibility. Additionally, the utilization of other residential equipment (e.g., heating and cooling) is also reduced during this period, contributing to decrease the flexibility margin. It is possible to observe that the main modifications in the load profiles with the activation of flexibility occurred in the case with "only electrical" buildings, to both residential and commercial buildings. Nonetheless, the flexibility utilization has a higher impact in residential buildings during the first 7 h of the day.

Impact on the Daily Load Profiles of Flexible Loads
The flexibility that buildings are able to offer comes essentially from using resources that consume gas instead of electricity and from DR programs. DR programs are capable of changing the load patterns of lighting, refrigeration, water heating and space heating. In this section, it is analyzed the flexibility that each source of flexibility is able to provide to the system operator. This highlights the importance of considering every possible source of flexibility within a building. Figure 14 shows the lighting consumption in a supermarket. The lighting systems can be decreased by 10% in relation to their initial value.
It is possible to observe a decrease in the lighting consumption in the scenario electrical buildings with optimization. This is the only scenario where the DR program was activated for lighting. In the MES scenarios, changes in the lighting consumption patterns did not happen (i.e., DR programs were not activated) since all the network problems detected were solved by replacing part of the electricity consumption by gas. This is the reason why the series of the scenarios electrical buildings without optimization, MES buildings without optimization and MES buildings with optimization overlap. In Figure 15 it is presented the temperature and the energy consumed for space heating in a room of the supermarket. Each type of building and each room has its own temperature requirements and different schedules of the HVAC system. In order to simplify the analysis of the results, only the requirements for the temperature of a room in a supermarket are presented, although the requirements for the other buildings are similar: • HVAC system is always "on"; • Between 6 h and 23 h: The temperature can vary between 15 • C and 25 • C.

• Between 0 h and 6 h and after 23 h:
-There is no set-point; - The temperature can vary between 15 • C and 25 • C.
As Figure 15 shows, the temperature requirements were met in all the scenarios. It is possible to observe that the case electrical buildings with optimization was the only case where flexibility of the space temperature was used. In this case, the DR program was activated and the temperature in the room decreased below the set-point but remained inside the allowed interval. This happened since the DR program was activated and the energy consumed for heating this room was decreased to zero. In the other scenarios, the temperature strictly followed the set-point defined. It is also possible to observe that between 0 h and 6 h and after 23 h, in all scenarios, the temperature was different from 21 • C due to the absence of set-point. Figure 16 shows the temperature profile of the refrigeration systems of the supermarket as well as the respective power profile (the series of the scenarios without optimization are overlapped). The temperature requirements of all refrigeration systems present in this network are the following: • The temperature is set to −5 • C for the entire day; • The temperature can vary between −6 • C and −3 • C by using the DR program.
More flexibility was used in the scenario electrical buildings with optimization, as the temperature increases due to the reduction of the power consumption to provide flexibility to the grid. In the scenario MES buildings with optimization, the system operator did not request any flexibility. However, some flexibility was used to decrease energy costs at the expense of a slightly higher temperature in the refrigeration system. An example of the flexibility utilization for this purpose is shown in Figure 16, where more power was consumed between 5-7 h to decrease temperature, as the energy cost in this period was lower.
In the other scenarios the temperature is the same as the set-point as the DR program is not available and their series are overlapped. Figure 17 presents the water temperature in the hospital at bus 13 and the respective power consumption. The water temperature requirements for all buildings present in the network are the following: • The temperature is set to 49 • C for the entire day; • The temperature can vary between 44 • C and 54 • C by using the DR program.
As it is possible to observe, the scenario electrical buildings with optimization was the only case where there was a decrease in the water temperature, with the correspondent power decrease to provide flexibility to the grid.
It is important to notice that the power necessary to heat the water depends on the set-point of the water temperature, which is a flexible value, and the water volume consumed, which is not flexible and it is set by the consumer. This way, the flexibility provided by heating the water is provided only by changing the water temperature and not the used volume. In Figures 18 and 19 it is presented the charging and discharging of EVs in a residence and a medium size office. It is possible to observe that the charging and discharging patterns in the scenario with electrical buildings with optimization differ from both MES buildings scenarios and the scenario with electrical buildings without optimization. This reflects the utilization of EVs flexibility to mitigate the undervoltage problems in the network.   Figure 20 presents the electricity and heat production profiles of the CHP unit in the supermarket. As it is possible to observe, there is a peak at 7 h. Comparing these figures with Figure 10, it is possible to observe that the supermarket consumes much less electricity from the network at this hour and relies more on the gas network and the CHP unit for the provision of electricity and heat. This is because gas prices are lower than electricity prices and by consuming less electricity the power network is not so overload and its limits are more relaxed.

Annual Results
In this section, it is first presented the annual consumption of electricity and gas as well as the annual emissions of CO 2 . After, the costs of buying electricity and gas, the costs of the emitted CO 2 and the costs of using DR programs are presented. The results obtained in this section allow the generalization of the results presented in the previous sections.
From Figure 21 and Table 3 it is possible to observe that the cases with only electrical buildings consume naturally more electricity than the cases with MES buildings. In total, the electricity consumed in the case with MES buildings without optimization is 33% lower than the case with only electrical buildings, while in the cases with optimization it is 30% lower. The utilization of flexibility (in the scenario with optimization) decreases the consumption of electricity by 6% in the only-electrical buildings cases, while in the MES buildings it only decreases 1%. This happens since in the cases with only electrical buildings, DR program are activated more often due to the technical problems in the network.
The gas consumption is very low in only electrical buildings cases, as the gas is only used for cooking in residential buildings. In the cases with MES buildings, more gas is used during the winter as it is the main source of heating.  In Figure 22 the annual emissions of CO 2 for each case are presented. It is possible to observe that MES buildings have naturally higher emissions than only electrical buildings. Although they consume less electricity, the gas consumption is much higher which makes CO 2 emissions higher. It is also possible to observe that the cases with optimization have lower emissions as they also consume less energy due to the activation of downward flexibility. These results show that DR programs can have an impact on the CO 2 emissions from buildings. In total, only electrical buildings have 40% and 44% lower emissions of CO 2 than MES buildings in the cases without optimization and with optimization, respectively. DR programs (i.e., optimization scenarios) lowered the emissions of CO 2 by 6% in the case with only electrical buildings and 4% in the case with MES buildings.
The costs of electricity and gas consumption, emissions of CO 2 and DR programs are presented in In Figure 23. As expected from the previous results, the costs of electricity are higher in the cases with only electrical buildings as more electricity is consumed than in the MES buildings cases. The costs of gas and CO 2 emissions are higher in the MES buildings scenarios. The case with only electrical buildings with optimization has the higher costs related with flexibility activation (DR programs), as this is the case that presents more technical problems in the electrical network.
It is also possible to observe that the total costs at the end-of the year are lower in the case with MES buildings with optimization, as they consume less electricity that have higher prices when compared to natural gas.
In total, without optimization, the annual costs of electricity in the cases with only electrical buildings are 32% higher than the cases with MES buildings, the CO 2 costs are 40% lower and the total annual costs are 10% higher. With optimization, the annual electricity costs are 29% higher, the costs of CO 2 emissions are 44% lower and the total annual costs are 29% higher.  Figure 23. Annual costs of electricity, gas, CO 2 emissions and DR programs.

Conclusions
This work presented a model to optimize energy consumption in buildings, aiming at minimizing costs while satisfying the technical constraints of the power network. The model developed is capable of controlling a wide variety of loads (thermal loads, lighting, EVs, PVs, etc.) taking into account the flexibility that their owners are willing to provide. This flexibility can be used for technical matters, i.e., to improve networks operation, or for economic purposes, i.e., decrease overall energy cost, including CO 2 costs. A novel linearized model for the HVAC system was also developed in this work, decreasing the complexity of the problem. The model was used with a test network to quantify and compare the capability of only electrical buildings and MES buildings, both commercial and residential, to offer flexibility. Annual costs, energy consumption and CO 2 emissions were also studied.
Analyzing the grid technical problems for the day 7 January 2019, in the scenarios with only electrical buildings, there are a number of undervoltage problems detected in the network and they cannot be entirely solved even when flexibility is activated through DR programs. These stressful operating conditions are not verified in the MES scenarios, where no technical problems were detected. In the only electrical buildings scenarios, the voltage problems detected were greatly reduced after the activation of the DR program-the number of buses with voltages below 0.9 p.u. was reduced from 54 to 14, representing a 74% reduction. The flexibility provided in the DR scenario was obtained from electrical and thermal loads, which can offer flexibility through thermal storage, lightings, EVs or CHP units and thus contribute to reduce the power absorbed from the grid. It is thus possible to conclude that DR programs are an important tool for system operators as they can enable important changes in the load profiles to avoid voltage problems or overloaded branches.
As expected, it is in the MES scenarios that electricity consumption is lower. One of the outcomes of this study is the quantification of the network load reduction that can be achieved in the MES scenarios in comparison with the only electrical buildings scenarios. In average, in the scenarios without DR, the load is 85% lower during the peak hour while in the scenarios with DR it is 67% lower. In the only electrical scenarios, the DR program was capable of reducing the peak load by 55%, although it was not enough to mitigate all voltage problems.
With lower electricity consumption, it is not surprising that buildings can reach lower absolute values of consumed power in MES buildings scenarios. In the comparison between commercial and residential buildings, the former revealed overall wider ranges of flexibility. Commercial buildings can generally provide more flexibility during the day, especially in the MES scenarios, while residential buildings are more flexible during the night and late afternoon.
The last outcome from this work is related with the annual costs, energy consumption and CO 2 emissions. It was concluded that MES consume ca. 33% less electricity, although they consume more gas. This way, annual energy costs are lower when considering MES buildings, presenting a reduction of ca. 32% in electricity costs and ca. 10% in total costs. In the cases with only electrical buildings, system operators need to solve technical problems through the activation of flexibility. This will increase total energy costs by ca. 25%, although the number of technical problems will be decreased by 74%. It was also concluded that only electrical buildings have ca. 44% and 40% lower emissions of CO 2 , with and without optimization, respectively. This is due to the fact that MES buildings consume more natural gas, while only electrical buildings depend strictly on electricity. If the mix of electricity generation is highly based on renewable sources, the emissions of CO 2 by consuming electricity will be lower than if consuming natural gas. The annual results follow the trend of the results obtained for the day 7 January 2019.
Future work should address the integration of the restrictions of the gas network and the models of shiftable loads, such as washing machines and dryers. Other energy vectors could also be integrated, like hydrogen or district heating. The utilization of the flexibility provided by MES buildings in energy markets, trough multi-energy aggregators, should also be studied. Other points to be taken in consideration and that should be further studied are the fixed costs in the consumer 's tariffs and the costs of DR programs. For planning purposes, investment costs should also be a matter of study.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.