Optimal Scheduling of Combined Heat and Power Generation Units Using the Thermal Inertia of the Connected District Heating Grid as Energy Storage

A better integration across sectors is an essential element of 4th generation district heating and smart energy systems allowing to react to volatile renewable energy generation. This sector coupling enables to use more cost-efficient storage as storage prices differ for different forms of energy. Thermal energy for example can be stored in comparably cheap storage tanks. Besides such dedicated storage, the thermal inertia of a heating grid can be used as thermal storage as well. In this paper, a classic unit commitment optimization for scheduling of combined heat and power units not considering grid dynamics is extended to cover thermal dynamics of heating grids. First an outer approximation of the grid storage capabilities is developed. Second, a very efficient formulation for the storage dynamics of a heating grid is introduced and its capabilities are shown in a motivating case study. In this study additional savings of several thousand Euros per day are achieved using the thermal inertia of a heating grid as storage.


Introduction
One of the essential elements in 4th generation district heating is a transition to smart energy systems with a better integration across sectors allowing to react to fluctuating renewable energy sources [1].Besides work looking into providing flexibility for renewable electricity integration within the power sector, work on a cross-sectoral integration including all sectors and/or energy carriers has been published [2].Hence, tools and approaches integrating different energy systems are needed [3].For providing flexibility to the power system, energy storage systems are crucial.As thermal energy storage tanks are several orders of magnitude cheaper than pumped hydro storage, an integration of the heating sector with the power market offers economically interesting flexibility options [4].Indeed, the thermal inertia of a heating grid can be used as heat storage as well, if it is considered in heat generation planning.Not requiring any investment into dedicated storage tanks, this approach allows to leverage additional flexibility at minimum cost.
Optimal scheduling of heat and electricity production units is traditionally done using unit commitment models with one energy balance for heat and one energy balance for power, neglecting grid dynamics [5].Integration of heating grid dynamics into scheduling of combined heat and power plants (CHPs) has been studied in several past publications.Most start from detailed physical models of the heating grid trying to achieve a model suitable for optimization by simplifying the grid topology [6] or assuming simplifications like fixed mass flows or fixed time delays in an iterative optimization scheme [7] or introducing a multi-step approach combining different simulation and optimization steps like static hydraulic, static thermal and dynamic thermal model [8] or a unit commitment and a dispatch model [9].Li et al considers the building thermal inertia in addition to pipeline dynamics with fixed mass flows [10].Dai et al solve the same problem as Li considering pipeline dynamics with fixed mass flows and building thermal inertia in an iterative optimization scheme [11].As well considering thermal inertia of buildings and heating grid dynamics [12] proposes an optimization model with fixed transport delays calculated based on design flows which allows to use off-the-shelf solvers.A very simple model of the heating grid as one dynamic energy mass is proposed and compared with other methods in a simulation in [13].Groß starts from a detailed physical model as well, but in contrast to others trains a linear regression model using multiple simulations.This regression function represents a direct link between energy stored within the heating grid and supply temperatures at the heat producer allowing a fast and efficient optimization [14].Giraud uses a similar function linking supply temperature changes at the heat producer with changes of the heat energy stored in the grid for optimization of supply temperature reducing heat and pumping losses.Its parameters are iteratively updated by a detailed grid simulation [15].
All approaches offer feasible optimization solutions, however as first a detailed grid model is needed, engineering efforts to setup such a solution are high.In this paper we are introducing a new method to consider the thermal inertia of a heating grid as storage in CHP scheduling.Based on a classic unit commitment model for CHPs not considering grid dynamics presented at the beginning of Section 3, we will first introduce an outer approximation of the thermal inertia of a heating grid in Section 3.1.This formulation is inspired by models for dedicated energy storage (like batteries or heat storage tanks).In Section 3.2 the heat energy stored in the grid will be modeled based on past and present supply temperatures.In contrast to [14,15] this function is created using time delays to the consumers as well as their share of consumption, allowing a more fast modeling as the creation of a detailed physical grid model is avoided.The fundamental dynamics of heating grids are explained in Section 2. Section 4 presents a motivating case study using the formulations of Sections 3.2 and 5 discusses the results and shows an outlook on future developments and applications.

District Heating Grid Dynamics
In district heating grids heat energy is transported with water at different temperatures via a pipe network.Hence, hydraulic and temperature dynamics might be of interest when the thermal inertia of a heating grid should be used as energy storage.As hydraulic dynamics are very fast compared to the temperature dynamics in a heating grid, they can be neglected when looking at the energy dynamics of a heating grid [8].The temperature dynamics are however rather slow and have an important effect on the energy balances of a heating grid.
The water transporting the thermal energy is heated at production sites to a controllable supply temperature and transported to consumers via supply pipes.The consumers are controlling the mass flow of heated water depending on their heat demand and the supply temperature.This consumed water is cooled down and returned to the production sites via separate return pipes.The change of inner energy of water at producer and consumer in a time slot t can be described with: where .m t is the mass flow of water between supply and return in time step t, T supply t is the supply temperature, T return t is the return temperature and c p is the specific heat capacity of water (assumed to be constant).
As the hydraulic dynamics can be neglected, it is assumed that mass flow changes induced by the consumers instantly influence the mass flow and hence the heat production at the producers.Changes in supply temperature by the producers however, are not directly seen at consumer stations, as the changed temperature first needs to propagate to the consumers with the water in the supply pipes.
As the maximum velocity of water in a supply pipe is limited, there is a delay between an increase of supply temperature at the producer and the related increase of supply temperature at the consumer.This delay results in an inherent thermal energy storage, charging when the supply temperature is increased and discharging when the supply temperature is decreased.
This storage using the thermal inertia of a heating grid is explained with the simple example with one producer and one consumer having constant heat demand shown in Figure 1a.In t2 the supply temperature at the generation is increased and following Equation (1) the heat output of the generation is increasing proportionally.As production in the following time slots is bigger than consumption, the heating grid shows a charging behavior.The increased supply temperature takes until t4 to reach the consumer.Now the consumer having a constant demand reduces the mass flow to keep Equation (1) in balance.Due to the fast hydraulics this mass flow change directly propagates to the producer and hence production is reduced to the original level and charging of the grid ends.When the supply temperature at the consumer is reduced in t7 we can see a similar inversed behavior discharging the grid until t10.The surface between produced energy and consumed energy is the energy charged or discharged in the heating grid.
Energies 2019, 12, x FOR PEER REVIEW 3 of 9 consumer.This delay results in an inherent thermal energy storage, charging when the supply temperature is increased and discharging when the supply temperature is decreased.This storage using the thermal inertia of a heating grid is explained with the simple example with one producer and one consumer having constant heat demand shown in Figure 1a.In t2 the supply temperature at the generation is increased and following Equation ( 1) the heat output of the generation is increasing proportionally.As production in the following time slots is bigger than consumption, the heating grid shows a charging behavior.The increased supply temperature takes until t4 to reach the consumer.Now the consumer having a constant demand reduces the mass flow to keep Equation (1) in balance.Due to the fast hydraulics this mass flow change directly propagates to the producer and hence production is reduced to the original level and charging of the grid ends.When the supply temperature at the consumer is reduced in t7 we can see a similar inversed behavior discharging the grid until t10.The surface between produced energy and consumed energy is the energy charged or discharged in the heating grid.For two consecutive consumers we get a behavior as shown in Figure 1b.A change of supply temperature at the consumer in time t2 is having a similar effect of increased generation as before.When the supply temperature reaches the first consumer in t3 the production is decreased only by the share of the load/mass flow of this costumer (in our case 50%).When the second consumer is reached in t6 the remaining overproduction is reduced.A decrease in supply temperature at the producer in t8 shows as before an inverted storage profile discharging the heating grid.

Optimal Unit Commitment Models
Unit commitment problems are a very well-studied field for electric power generation [16].Mixed-integer linear programs have shown to be an efficient way to formulate and solve this problem with off-the-shelf open source or commercial solvers.Most unit commitment models are time discrete models.Hence, the time horizon is discretized into several time slots of equal size and for each time slot  the variables of the model are optimized.Key goal is to find optimal setpoints for the electric power output  , and running status  , (binary variable: 0:off/1:on) of each generator  in each time slot .Electric power output and running status are linked with: For two consecutive consumers we get a behavior as shown in Figure 1b.A change of supply temperature at the consumer in time t2 is having a similar effect of increased generation as before.When the supply temperature reaches the first consumer in t3 the production is decreased only by the share of the load/mass flow of this costumer (in our case 50%).When the second consumer is reached in t6 the remaining overproduction is reduced.A decrease in supply temperature at the producer in t8 shows as before an inverted storage profile discharging the heating grid.

Optimal Unit Commitment Models
Unit commitment problems are a very well-studied field for electric power generation [16].Mixed-integer linear programs have shown to be an efficient way to formulate and solve this problem with off-the-shelf open source or commercial solvers.Most unit commitment models are time discrete models.Hence, the time horizon is discretized into several time slots of equal size and for each time slot t the variables of the model are optimized.Key goal is to find optimal setpoints for the electric power output p i,t and running status n i,t (binary variable: 0:off/1:on) of each generator i in each time slot t.Electric power output and running status are linked with: For each generator i in each time slot t with P min i and P max i being the minimum and maximum possible electric output of the generator.It is assumed that energy demand p demand t for every time slot is known.If it is possible to sell and buy electricity from the electricity market, we introduce positive variables p buy t and p sell t being the electricity bought and sold to the energy markets.If no electric grid constraints are considered, we get the resulting electric energy balance: Classic goals of unit commitment are either minimization of cost or maximization of profit.As objective function minimizing the overall cost we get: With price el,buy t and price el,sell t being the time varying prices for selling or buying electricity at the energy markets and C el,const i being the cost for running generation i in a time slot and C el,var i being the variable cost of generation of generator i.
If heat generation of combined heat and power plants (CHPs) is considered, we get an additional variable q i,t representing the heat output of the CHP.Like Equation (2) we get: With Q max i being the maximum heat output of CHP i. Electric output and heat output of a generator are linked.For gas turbines the power-to-heat ratio is usually constant.For CHPs with extraction condensing turbines both are linked as shown in Figure 2 [17].
Energies 2019, 12, x FOR PEER REVIEW 4 of 9 For each generator  in each time slot  with  and  being the minimum and maximum possible electric output of the generator.It is assumed that energy demand  for every time slot is known.If it is possible to sell and buy electricity from the electricity market, we introduce positive variables  and  being the electricity bought and sold to the energy markets.If no electric grid constraints are considered, we get the resulting electric energy balance: Classic goals of unit commitment are either minimization of cost or maximization of profit.As objective function minimizing the overall cost we get: With  , and  , being the time varying prices for selling or buying electricity at the energy markets and  , being the cost for running generation  in a time slot and  , being the variable cost of generation of generator .
If heat generation of combined heat and power plants (CHPs) is considered, we get an additional variable  , representing the heat output of the CHP.Like Equation (2) we get: With  being the maximum heat output of CHP .Electric output and heat output of a generator are linked.For gas turbines the power-to-heat ratio is usually constant.For CHPs with extraction condensing turbines both are linked as shown in Figure 2   If heating grid dynamics are not considered, we can formulate a heat energy balance like the electric energy balance Equation (3): with  being the known overall heat demand and  being the thermal losses per time slot , which are neglected in the following optimization.With  , being the variable cost for heat generation the objective function becomes: If heating grid dynamics are not considered, we can formulate a heat energy balance like the electric energy balance Equation (3): with q demand t being the known overall heat demand and q loss t being the thermal losses per time slot t, which are neglected in the following optimization.With C heat,var i being the variable cost for heat generation the objective function becomes:

Outer Approximation of Heating Grid Dynamics Using State-of-Charge Formulation
The integration of dedicated thermal storage with storage tanks into a unit commitment problem can be achieved by introducing a state-of-charge variable SOC heat t at each time slot t as well as a variable for the charging and discharging heat power q store t per time slot t.The state-of-charge can implicitly be formulated as: where q store t is limited by the maximum charging and discharging power and SOC heat t is limited by the maximum and minimum energy storage capacity of the thermal storage tank.
Introducing the state-of-charge formulation the heat energy balance Equation ( 6) without consideration of heat losses becomes: With adapted limits of q store t and SOC heat t this formulation can approximate the heat storage capability of a heating grid as shown in Figure 1.For a given mass flow .m t the limits of the grid charging and discharging power can be calculated based on Equation (1): with T supply,max being the maximum possible supply temperature and T supply,min t being the minimum possible supply temperature in time slot t.The minimum possible supply temperature can be calculated upfront based on expected heat demand or outdoor temperature.
If the maximum time delay to the last consumer τ max is known, an upper limit of energy stored in the grid at time t can be formulated as a summation of the last upper bounds of charging power: m tt T supply,max − T supply,min tt (11) However, this formulation is only an upper bound of the maximum of the heat energy stored within the grid as only the last and not all consumers are feed with delay τ max (Figure 1b).Hence, the presented Equations ( 8)-( 11) leads to an outer approximation of the storage capabilities of a heating grid.It might be used for an approximation of economic potential but is not suitable for operations planning as results might be outside the feasible operation region.

Modeling of Heating Grid Dynamics with Estimation of Stored Energy
Looking at Figure 1 the heat stored within the heating grid seems to be linked to the supply temperature.Hence, a function directly linking T supply t to q stored t is possible for heating grids with one major generation site if mass flows vary little.Groß proposes to train such a model based on simulations using linear regression [14].Here we are deriving such a function based on transport delays to consumers and consumers' share in overall consumption.It is assumed that demand forecasts for all consumers are available.Additionally, it is assumed that initial supply temperatures as well as respective mass flows/transport delays are given.They could be calculated by a supply temperature optimization like [15].
As in the previous section we can formulate the heat energy balance not considering heat losses: However, to get a more exact representation of the thermal inertia of a heating grid, we need to calculate q store t differently.Assuming constant mass flow and hence constant transport delays, the supply temperature at the producer at time t2-t3 in Figure 1a is increased by: T supply t3 − T supply t2 (13) Leading to an instantaneous increase in CHP heat output by: As soon as the supply temperature increase reaches the consumer after delay τ = 2 in time t4-t5 the CHP output is decreased by: Not considering heat losses the energy stored in the heating grid follows: As the heat load is constant in Figure 1a, the changes in heat output of the CHP directly are equal to changes in the energy charged to or discharged from the grid.Hence, combining Equation (14) and Equation (15) we get: with T supply t being the supply temperature at the producer, τ being the time delay between producer and consumer and .m being a constant mass flow.If a time varying mass flow .m t is introduced as a parameter, Equation (17) becomes: To consider multiple consumers, as in the example in Figure 1b, we introduce a matrix M based on time delay to and share of consumption of a consumer.In this matrix row 1 has non-zero entries whenever a water volume starting at time one at the producer reaches a consumer; row 2 is having non-zero entries when a water volume starting at time two at the producer reaches a consumer et cetera.The non-zero values in M matrix correspond the share of consumption of the consumer being reached.Those shares are assumed to be constant.Hence every row sums up to a value of 1 if all consumers are reached within the time horizon.
Looking at the example given in Figure 1b an increase in supply temperature at the producer in t2 reaches the first consumer in t4 and the second consumer in t6.Having equal share in overall consumption non-zero entries in row t2 are 0.5 in column t4 and 0.5 in column t6 (Table 1).Having a row for each time instant allows to account for variable time delays if mass flows are not constant.Extending the calculation of q stored t Equation ( 18) with matrix M (Table 1) to allow multiple consumers we get: With Equation ( 12), Equation ( 19) and Table 1 we know have a model which allows to directly calculate the currently charged or discharged thermal energy in the heating grid based on a supply temperature profile, assuming demands of and transport delays to consumers as well as mass flows at the producer are known.

Motivating Case Study
As a motivating case study, we have a closer look at a heating grid with one CHP and two consumers (comparable to Figure 1b).The CHP is having a maximum electric output of 500 MW and a maximum heat output of 400 MW.Electricity is produced using an extraction condensing turbine and hence power-to-heat ratio is variable within the limits shown in Figure 2. Heat consumer 1 having a load share of about 45% and a peak load of 120 MW is connected to the CHP with a 20 km pipe with 0.7 m diameter.Heat consumer 2 having a load share of about 55% and a peak load of 145 MW is following consumer 1 with a 10 km pipe with 0.7 m diameter.Electricity prices are taken from EPEX Day-Ahead Spot Market for Germany.The optimization horizon is one day with 24 time slots (one per hour).Matrix M as well as mass flows are calculated assuming that water in the pipes is flowing at maximum velocity of 1.5 m/s.
If supply temperature can be increased for heat storage by a maximum of 30 K, the production cost can be reduced by 2.4% (several k€/day) compared to a scenario without heat storage.Corresponding results with supply temperature and heat power profiles are shown in Figure 3. temperature profile, assuming demands of and transport delays to consumers as well as mass flows at the producer are known.

Motivating Case Study
As a motivating case study, we have a closer look at a heating grid with one CHP and two consumers (comparable to Figure 1b).The CHP is having a maximum electric output of 500 MW and a maximum heat output of 400 MW.Electricity is produced using an extraction condensing turbine and hence power-to-heat ratio is variable within the limits shown in Figure 2. Heat consumer 1 having a load share of about 45% and a peak load of 120 MW is connected to the CHP with a 20 km pipe with 0.7 m diameter.Heat consumer 2 having a load share of about 55% and a peak load of 145 MW is following consumer 1 with a 10 km pipe with 0.7 m diameter.Electricity prices are taken from EPEX Day-Ahead Spot Market for Germany.The optimization horizon is one day with 24 time slots (one per hour).Matrix M as well as mass flows are calculated assuming that water in the pipes is flowing at maximum velocity of 1.5 m/s.Table 2 shows achievable savings for different maximum increases of supply temperature.As expected the savings increase, with increasing flexibility of supply temperature and hence increasing storage capacity of the heating grid.As shown in Table 3 variations of the demand share have a rather small influence on the results.Variations of velocity lead to a bit larger changes, however remain about one order of magnitude smaller than the achieved savings.In almost all cases the supply temperature profile does not change.Only for a reduction of velocity by 10% the supply temperature profile shows minor changes like a lower temperature increase of 22.6 K instead of 30 K in hour 6.The optimization problem is modeled in Julia/JuMP [18] and solved with Gurobi 8.0 [19] on an Intel®Xeon®v2 with 2.6 GHz (2 processors) and 4 GB RAM.Solution times for all problems shown in Tables 2 and 3 are around 16 ms never exceeding 20 ms, hence an application for online-optimization seems in reach even for larger heating grids and real-world cases.

Discussion/Comparison of Modeling Approaches
In this paper we introduced two different ways to model the thermal inertia of heating grids without consideration of thermal losses.First, an outer approximation was derived using a simple state-of-charge formulation.As this representation did not lead to accurate results, second, an estimation of the stored energy using past and present supply temperatures was introduced.This formulation is reducing the modeling effort compared to many other approaches, as only maximum supply temperature deviations, transport delays to consumers and their respective load shares need to be parametrized.Hence, an extension to real world district heating grids is possible without major efforts, if consumers are grouped into consumption areas having similar transport delay from the producer.This consumption areas and respective transport delays can be used to create the M matrix.
As some but not all approaches in literature, this formulation only allows us to model dynamics of a heating grid feed by one generation site.Another limitation is that constant transport times from heat source to consumers are assumed.This leads to an error as when an increased temperature reaches a consumer the corresponding mass flow and hence velocity and transport time are changed to keep Equation (1) in balance.If temperature increases are rather short and small, this is neglectable, but with long and large temperature increases this can decrease velocity importantly.One way to overcome this limitation could be to use iterative updates with a simulation model like in [7].
The results of the motivating case study show that the presented formulation is computationally very efficient and outperforms many existing algorithms.The model is quite robust against small parameter variations and it allows interesting cost savings of several thousand Euros per day as shown with calculations for different maximum supply temperature increases.As a high temperature variation leads to increased component fatigue (e.g., for pipes), in real world applications the maximum increase of supply temperature will most likely be limited to 10K or 20 K being intra-day temperature variations which can be seen in today's operations as well.
In the case study only one pipe with a large diameter is considered and hence heat losses can be neglected.However, if a distribution grid is included, increase of thermal losses is becoming more important.Hence, future work will try to overcome the inaccuracies of the fixed mass flows, address the integration of heat losses and apply the second approach to real-world heating grids.

Figure 1 .
Figure 1.Storage behavior of a heating grid over time for (a) one consumers; (b) two consumers with constant load.Upper graphs: supply temperature at producer (blue) and at consumer (orange, gray) in °C.Lower graphs: heat output of producer (black) and sum of load of consumers (orange) in MW.

Figure 1 .
Figure 1.Storage behavior of a heating grid over time for (a) one consumers; (b) two consumers with constant load.Upper graphs: supply temperature at producer (blue) and at consumer (orange, gray) in • C. Lower graphs: heat output of producer (black) and sum of load of consumers (orange) in MW. [17].

Figure 2 .
Figure 2. Feasible operating region of a CHP with extraction condensing turbine.

Figure 2 .
Figure 2. Feasible operating region of a CHP with extraction condensing turbine.

Figure 3 .
Figure 3. Optimization results for a maximum supply temperature increase of 30 K: (a) Day Ahead Electricity Price in €/MWh; (b) Supply Temperature increase at CHP in K; (c) Heat power profiles in MW.

Figure 3 .
Figure 3. Optimization results for a maximum supply temperature increase of 30 K: (a) Day Ahead Electricity Price in €/MWh; (b) Supply Temperature increase at CHP in K; (c) Heat power profiles in MW.

Table 1 .
Matrix M showing time delay and share of consumption.

Table 2 .
Savings for different maximum increases of supply temperature.

Table 2 .
Savings for different maximum increases of supply temperature.

Table 3 .
Variation of objective induced by parameter variations for a maximum supply temperature increase of 30 K.