A Matheuristic Approach for Resource Scheduling and Design of a

Multi-energy System

A. Bartolini

, G. Comodi

, F. Marinelli

, A. Pizzuti

and R. Rosetti

DIISM - Dipartimento di Ingegneria Industriale e Scienze Matematiche, Universit Politecnica delle Marche,

Via Brecce Bianche, I-60131 Ancona, Italia

DII - Dipartimento di Ingegneria dell’Informazione, Universit Politecnica delle Marche,

rosetti@dii.univpm.it

Keywords:

Multi-energy System Management, District Design, MILP, Matheuristic.

Abstract:

Modern energy system are evolving due to the opportunities and challenges that new technologies pose in the

energy sector. These changes create the requirements of decision tools able to effectively sustain the processes

of design and retroﬁt of energy systems. In this paper a multi-energy system management problem is taken into

account and a mixed integer linear programming (MILP) formulation is proposed to model both the design and

the resource scheduling of energy districts. However, since the size of the formulation restricts its applicability

to small cases far from the application of interest, a matheuristic based on constraint relaxations and variable

ﬁxing has been designed. Preliminary computational results show that the proposed solution strategy is able

to achieve good solutions (i.e., solutions with small optimality gaps) on restricted random instances, and to

solve in reasonable times instances derived from a real case study.

1 INTRODUCTION

The push towards a more rational and sustainable use

of primary energy resources is increasing worldwide

and many countries and organization take concrete ac-

tions and commitments in addressing such challenges,

e.g. (European Commission, 2012), (BNEF, 2017).

Among those actions, a common commitment is to

adopt higher shares of renewable energy generation

technologies into the energy mix to satisfy the users

demand for energy and/or other commodities. More-

over, several studies as (European Commision, 2017)

and (ISE, 2015) predict in the near future a strong

reduction of the purchase cost of such technologies

that would induces a dramatic change, with relative

challenges, in existing energy distribution infrastruc-

tures. Most of renewable energy generation technolo-

gies are in-fact characterized by an intermittent and

non-dispatchable electric generation capacity, making

the integration with current national and local distri-

bution grids challenging (California ISO, 2012).

A possible solution lies in regional large central-

ized energy storage systems, such as pumped hy-

dro, compressed air and large battery storage systems.

Such approach has been followed both in many stud-

ies and actual projects (J

ulch, 2016). Another emerg-

ing approach to address such challenging integration

resides in decentralization with the aid of microgrids

(Speer et al., 2015). While microgrids are already a

well established solution to accommodate the needs

of rural off-grid communities, their deployment in al-

ready urbanized contexts of developed countries is

also gaining attention from researchers (Enea, 2017),

(Center for Climate and Energy Solutions, 2017).

The basic idea behind microgrids is strictly re-

lated to electric energy management, but their poten-

tial to handle and optimize the employment of re-

newables and primary energy resources greatly im-

proves in multi-energy systems (MES), that is sys-

tems in which multiple energy vectors (i.e., carriers

of energy as electricity) are concurrently considered

(Mancarella, 2014). To capture the potential opportu-

nities for efﬁcient energy use (polygeneration, (Jana

et al., 2017)), MES models have to explicitly consider

the interactions across different energy vectors, in ac-

cordance with the real scenario under analysis, e.g. a

comprehensive assessment of users demanded com-

modities such as heat, cooling and electric energy.

In order to exploit all the potential energy savings

opportunities, such cross vector interactions have to

be taken into account since the design phase of both

new energy systems and the retroﬁt of already exist-

Bartolini, A., Comodi, G., Marinelli, F., Pizzuti, A. and Rosetti, R.

A Matheuristic Approach for Resource Scheduling and Design of a Multi-energy System.

DOI: 10.5220/0007574104510458

In Proceedings of the 8th International Conference on Operations Research and Enterprise Systems (ICORES 2019), pages 451-458

ISBN: 978-989-758-352-0

451

ing settings. Thus, MES models must preferably be

adapted to represent such opportunities and have a

correct assessment of what high integration of renew-

ables might imply.

In this paper a mixed integer linear programming

(MILP) formulation is proposed to model the design

of MES, where the optimal investment choices are

performed by evaluating also the required commodi-

ties from the exogenous energy networks. Speciﬁ-

cally, the formulation aim at minimizing the set-up,

maintenance and procurement costs by taking into ac-

count all the operational constraints, i.e., the schedul-

ing of the operation performed by the available tech-

nologies in order to convert energy vectors. Due to

the inner complexity of the addressed problem and the

size of the mathematical programming formulation,

a matheuristic scheme based on ﬁxing integer vari-

ables and relaxing constraints is designed to produce

good sub-optimal solutions in a reasonable amount of

time. To the best of our knowledge such scheme has

only been applied by (Triad

o-Aymerich et al., 2016)

for energy systems related problems, where the opti-

mal design for a rural electrical microgrid is evaluated

using different heuristic approaches, among which a

ﬁx-and-relax strategy. In the last few years the opti-

mal design and operation management of distributed

MES has been addressed by means of several methods

(Singh and Sharma, 2017). The problem has been for-

malized with both single or multi-objective functions,

and solved by techniques that can be mainly clas-

siﬁed into single monolithic and decomposition ap-

proaches. With respect to the former, (Mehleri et al.,

2012) and (Omu et al., 2013) exploit MILP formula-

tions to compute the optimal design of a small district

needing electricity and heat; the same is done in (Bis-

chi et al., 2014) where the demand of cooling power

is also considered. For the latter, decomposition ap-

proaches rely on keeping the planning phase, that is

the decision to purchase and deploy, separated by the

operational phase. In (Elsido et al., 2017) a mixed in-

teger non linear programming (MINLP) formulation

is presented and solved by a two phases approach,

where the ﬁrst is handled by means of two evolu-

tionary algorithms and a discrete variable relaxation

technique, while the second is solved by exploiting

a linearized MILP. In (Li et al., 2017) the design of

an hydrogen based microgrid is performed by using

a genetic algorithm to determine the size of the sys-

tems, while the operational variables are set by solv-

ing a MILP. A similar technique is described in (Sachs

and Sawodny, 2016), where the same framework is

applied for a multi-objective variant of the problem.

The rest of the paper is structured as follows: in

Section 2 the approach to the energy hub and tech-

nologies modeling is described, in Section 3 the prob-

lem is formalized and the mathematical model is pre-

sented. Section 4 describes the matheuristic used; in

Section 5 preliminary computational results are pre-

sented, based on both randomly generated instances

and a real case study. Finally in Section 6 conclusions

and future perspectives are presented.

2 ENERGY HUB AND

TECHNOLOGY MODELLING

In this section the concept of energy hub and the spe-

ciﬁc setting related to the case study in §5.2 are de-

scribed. Furthermore, the technologies that are con-

sidered within the energy hub and their general func-

tion are depicted.

2.1 The Energy Hub

Among the suitable tools for the management of

MES, a relevant option can be found in energy hub

(Mohammadi et al., 2017). An energy hub can be

thought as an integrated unit within which the con-

version and storage of different energy carriers is un-

dertaken in order to satisfy different commodities de-

manded by an user. The energy hub can get some of

such energy vectors from outside of its boundaries,

as for example the withdrawal of electric energy or

natural gas from the respective national energy distri-

bution infrastructures. Then, within the boundaries

of the energy hub such vectors can be transformed

in order to meet the demand of users, this provided

that the technology needed to operate the conversion

is installed. Such technologies can either be already

available or be deployed within the energy hub by sus-

taining the relative purchase and maintenance cost.

Thus, the choices regarding which kind of technolo-

gies to deploy within an energy hub have to be as-

sessed by considering both the monetary investment

needed to purchase such technologies and how such

systems will work, providing the needs of the com-

munity of users tied to the energy hub.

Energy hub is then a general scheme that can be

adapted to various settings and the case study pre-

sented in §5.2 reﬂects the structure presented in Fig-

ure 1. It serves a community consisting in buildings of

the tertiary sector, which will need to properly func-

tion electric energy, space heating for the winter and

cooling for the summer. It has access to both the nat-

ural gas network and the national electric energy grid;

also has access to a solar resource to be used to power

the district.

ICORES 2019 - 8th International Conference on Operations Research and Enterprise Systems

452

Figure 1: Proposed energy hub model.

2.2 Technologies Modeling

The technologies deployed in the energy hub convert

energy among the different vectors, or store it for a

later use. The technologies considered in this study

include the majority of the most mature ones available

on the market for small urban districts. These can be

listed as:

• Natural Gas Engine (CHP): operates in com-

bined heat and power (CHP) mode, thus produc-

ing electricity and simultaneously recovering the

waste heat, using natural gas as fuel;

• Natural Gas Boiler (GB): production of heat

from natural gas;

• Electric Chiller (EC): conversion of electric en-

ergy in cooling power;

• Absorption Chiller (AC): conversion of heat in

cooling energy;

• Heat Pumps (HP): allows the conversion of elec-

tric energy into (alternatively) cooling energy or

heat, based on the required output;

• Photovoltaic Panels (PV ): production of electric

energy from solar radiation;

• Batteries (EES): storage of electricity;

• Hot Thermal Storage (HT ES): storage of heat;

• Cold Thermal Storage (CT ES): storage of cool-

ing energy.

A small database of purchasable devices is con-

sidered for each of such technologies, distinguished

by size. Then for each particular model, parameters

representing costs and performance are speciﬁed. The

performance of conversion devices is described by the

energy conversion efﬁciency for each available size,

whereas the relevant parameters for the storage de-

vices are the efﬁciency of both charge and discharge

phases. Besides, conversion devices are enforced to

work above a ﬁxed threshold value for a certain times-

pan, this to represent a lower bound of partialisation

(i.e., the minimum output power) for the technology.

The only exception concerns the photovoltaic panels,

whose produced electric energy is directly dependent

by the total panel surface and the solar exposition.

Hence, in absence of solar radiation the output of the

photovoltaic can be zero.

Regarding the cost, each model of technology is

characterised by its own purchase cost related to its

size and performance. Such cost has to be sustained

in order to make use of the technology for a represen-

tative lifetime in years, after which the purchase cost

has to be re-sustained. Moreover, additional costs re-

lated to the maintenance operations are considered.

In this way, all of the costs related to purchasing and

operating phases of a certain technology are taken

into account. Finally, further costs sustained within

the energy hub are given by the possible withdrawal

of energy from the two national infrastructures. The

withdrawal is priced with a ﬁxed cost per kWh of pur-

chased energy, different for natural gas and electricity.

3 PROBLEM DEFINITION AND

MATHEMATICAL

FORMULATION

The problem can be summarized as follow: an en-

ergy district has to be designed and managed (by se-

lecting an equipment of technologies and by control-

ling the relative operations) in order to satisfy the de-

mand of different energy vectors (commodities). De-

cisions have to be performed such that the total costs

sustained by the energy district within a certain time

horizon are minimized, both in the planning and op-

erational phases. We assume that the time horizon is

cyclically repeated to describe the optimal dynamic of

a fully operational energy district, satisfying requests

that (hopefully) follow periodical patterns.

Let K be the set of commodities, with a ﬁxed unit

purchasing and storing price, that is respectively c

and h

. Then, let us deﬁne Q as the set of technol-

ogy models that can be deployed in the district, e.g.

heat pump 170kW or electric chiller 120kW. More-

over, let S ⊆ Q be the subset of all the storage systems

and S

⊆ S that of models that store the commodity

k. Parameters v

are deﬁned for each i ∈ Q, stating

the costs for buying, installing and maintaining the i-

th appliance type. This cost is normalized along the

considered period: indeed, the approximated proﬁle

of the demand of each commodity is reported for a

discrete time horizon T = {1,...,

T } with intervals of

one hour; namely d

is the demand of the commodity

k ∈ K at time t ∈ T .

A Matheuristic Approach for Resource Scheduling and Design of a Multi-energy System

453

Therefore, the matrix v is normalized in this way:

given the purchasing and installing cost (IC), main-

tenance cost per year (MCY) and estimated lifetime

(expressed in year and hours respectively as ELY and

ELH), each element of v is deﬁned by:

+ MCY

× ELY

ELH

Although ﬁxed costs should be fully sustained in the

same moment, the deﬁnition of v

includes only the

portions of hourly costs related to the considered time

horizon T . Hence, the ratio between ﬁxed and vari-

ables costs (e.g. c

) is contained, while the objective

function weights the total cost associated to the life

cycle of each device.

As previously mentioned, each device can require

and supply different types of commodity, according

to his features. For example a Gas Engine uses gas

and can supply heat or electricity, whereas an Elec-

tric Chiller needs electricity to work and supply cool.

A parameter ω

is then deﬁned only for conversion

technology i ∈ Q \ S, expressing the unitary conver-

sion multiplier to convert commodity h into commod-

ity k by means of appliance i. Some devices, such

as solar panels or wind farms, have a base production

that exploits renewable resources and are dependent

by weather conditions at a certain time. Thus, pa-

rameter b

states the base production for commodity

k ∈ K produced by the equipment i ∈ Q \ S at time t.

The amount of output k produced by the conversion

device i at each time t ranges within the maximum

rated power U

and the lower bound L

Regarding the storage technologies, each device i ∈ S

is able to collect commodity k up to a ﬁxed capacity

, with a common efﬁciency φ

We deﬁne non-negative variables r and p to model

the ﬂow of commodities converted by conversion de-

vices: at each time instant t, variable r

∈ R

ex-

presses the amount of commodity h converted into

commodity k by device i ∈ Q \ S; on the other hand,

variable p

∈ R

models the output commodity k re-

sulting by the conversion operation related to r

and

alternative base productions b

(e.g. solar energy).

Furthermore, for each k ∈ K and t ∈ T we deﬁne vari-

ables l

and f

, to describe for the time t respectively

the total amount of commodity k stored within the dis-

trict, and the quantity of k acquired by the exogenous

supplier. By assuming that each technology has an

estimated lifetime much larger than

T , we can use bi-

nary variables y

∀i ∈ Q to model the selection of de-

vices, i.e. y

= 1 if and only if the district is provided

with equipment i since the beginning of T . Finally, for

conversion technologies i ∈ Q \ S and for each t ∈ T ,

variables z

∈ {0, 1} are further employed to describe

the operation of device i at time t through z

= 1.

The proposed MILP formulation reads as follow:

min

∑

t∈T

∑

k∈K

+ h

) +

∑

i∈Q

(1)

+ l

∑

i∈Q\S

∑

h∈K:

h6=k

t−1

∑

i∈Q\S

+ f

∀k ∈ K,∀t ∈ T (2)

∑

h∈K:

h6=k

+ b

∀i ∈ Q \ S,∀k ∈ K, ∀t ∈ T

(3)

≤ p

≤ U

∀i ∈ Q \ S,∀k ∈ K, ∀t ∈ T

(4)

∑

t∈T

≤ |T |y

∀i ∈ Q \ S (5)

≤

∑

i∈S

∀k ∈ K,t ∈ T (6)

= l

∀k ∈ K (7)

, f

≥ 0 ∀k ∈ K,∀t ∈ T (8)

≥ 0 ∀i ∈ Q \ S,∀k ∈ K, ∀t ∈ T (9)

≥ 0 ∀i ∈ Q \ S,∀h, k ∈ K : h 6= k,∀t ∈ T

(10)

∈ {0,1} ∀i ∈ Q (11)

∈ {0,1} ∀i ∈ Q \ S,∀t ∈ T (12)

The objective function (1) aims to minimize the

total cost given by the cost related to deployed tech-

nologies, the cost of acquiring external commodi-

ties from public suppliers and the expense for storing

surpluses for further use. Constraints (2) are stock-

balancing constraints and state that the commodity k

demanded, stored and required by devices is equal to

the summation among the previously amount of com-

modity k stored at time t −1, the total one produced by

all deployed technologies and the quantity acquired

from external networks. The set of constraints (3)

deﬁnes the variables p, for each conversion technol-

ogy i ∈ Q \ S , as the summation between the required

commodity h necessary to produce commodity k and

a base production without consumption. This is due to

the feature of green technologies, such as photovoltaic

and wind farm, to produce commodities without con-

suming other resources. The bounds of variables p are

given by the set of constraints (4): if binary variable

is equal to 1, the corresponding variable p

is lim-

ited between values L

and U

. The link among vari-

ables z and y is coherently modelled by the set of con-

straints (5). Constraints (6) limit the stored quantity

of commodity k ∈ K at any time t ∈ T to the total sum

of capacities of the storage systems in S

installed. In

order to translate the cyclicality of the selected time

ICORES 2019 - 8th International Conference on Operations Research and Enterprise Systems

454

horizon, the quantities of commodities stored in the

ﬁrst timestep (t = 1) are balanced in (2) with the quan-

tities stored in the last timestep (t =

T ) by means of

equality (7). It is worth noting that whenever an ap-

pliance i is not able to convert a certain commodity h

into k, the corresponding parameter ω

is null. In all

such cases, ∀t ∈ T variables r

are ﬁxed to zero as a

preprocessing phase. Moreover, ∀t ∈ T variables f

are set equal to zero for any commodity k ∈ K that is

not purchasable by external suppliers.

4 MATHEURISTIC APPROACH

Computational experiments showed that the mathe-

matical formulation proposed in §3 cannot be ex-

ploited to solve real-size instances in a reasonable

amount of time by means of standard solvers. This

can be (partially) related to the size of the formula-

tion, as the number of constraints is O(|Q||K||T |) and

the number of binary variable is O(|Q||T |), whereas

continuous variables are O(|Q||K|

|T |). Therefore,

we design a matheuristic algorithm to approach in-

stances that are meaningful in real context.

The matheuristic used to solve the proposed MILP

is based on relaxing integrality constraints and ﬁx-

ing subsets of variables, similarly to the Fix & Re-

lax scheme formalized in (Escudero and Salmeron,

2005). Starting with the linear relaxation of the math-

ematical formulation, the idea is to restrict it by restor-

ing the integrality conditions only for a subset of the

original integer variables. Then, the restricted linear

relaxation is solved and the integer variables are ﬁxed

to their values in the resulting solution; hence, the

number of variables is reduced and the space of re-

search shrunk. This process is iterated until a feasible

integer solution for the original problem is found.

A common challenge for this kind of approach

is to ﬁnd an equilibrium between the algorithm efﬁ-

ciency and quality of the solution. Indeed, while the

former usually requires to proceed by small subset of

variables, the latter could be compromised by such

scheme as the ﬁrst steps can bring the current solution

far from the optimum. Moreover, a robust strategy is

one that can avoid the convergence to an unfeasible

solution, or is able to efﬁciently backtrack whenever

this occurs.

The matheuristic approach is actually a cus-

tomization of the Fix & Relax scheme for the energy

district design and management problem, obtained by

specifying the strategy of selection and ﬁxing of the

subset of variables. In particular, the subsets of bi-

nary variables are obtained by separating the corre-

sponding technologies for their capacity of convert-

ing or storing commodities. Indeed, each technology

rarely produce (store) more than one commodity, thus

creating a natural separation between the relative vari-

ables. Moreover, as detailed below the ﬁxing policy

is deﬁned to allow a certain grade of ﬂexibility to the

matheuristic and to reduce the possibility of reaching

unfeasible states.

The scheme can be summarized as follow. Let us

consider the mathematical formulation proposed in §3

as M . Let V

and V

be the sets of variables in M for

which the binary condition holds, related respectively

to y and z. In a similar way, let us deﬁne V

and V

the subsets of variables whose integrality is relaxed in

M . Finally, let K

be an auxiliary set initialized as K

and collecting all the commodities.

The algorithm starts by considering the linear re-

laxation of M , that is V

0 and V

0. Firstly, the

demands and costs for acquiring external resources

are used to deﬁne a priority among the commodities

k ∈ K. Speciﬁcally, at each iteration the commodity

s ∈ K

is selected accordingly to the following rule:

s = argmax

k∈K

,t∈T

} (13)

That is, the commodity chosen is that one with the

largest cost required to satisfy its own demand peak

through the external supply. Commodity s is then re-

moved from K

and the integrality is restored for the

variables y

∈ V

modeling appliances able to con-

vert (store) the commodity s. In other words, for all

the conversion technologies i ∈ Q \ S such that exists

> 0, with h ∈ K : h 6= s, the binary condition for

is restored, that is y

is removed from V

and added

to V

; the same is done for the storage systems i ∈ S

with positive C

. Formulation M is then solved and

for each variable y

just added to V

, its value is ﬁxed

to one if and only if y

= 1 in the resulting solution.

Furthermore, only for conversion devices i ∈ Q \ S,

the value of y

is set equal to zero if y

= 0 currently

holds. By excluding the variables related to the stor-

age devices, the ﬂexibility of the current solution is

partially preserved. In fact, if in a later iteration re-

lated to h 6= s the requests r

increase due to deploy-

ment of different appliances, then storage devices of

s are still valid options to avoid the increment of f

For each y

currently ﬁxed to 1, all the corresponding

in V

are ﬁxed to one if z

> 1 − ε (with ε reason-

ably small), and moved from V

to V

. This process is

iterated until K

is empty.

Afterwards, the integrality for all the variables left

in V

and V

is restored and all the subsets are coher-

ently updated. In order to reduce the number of binary

variables in V

, a further ﬁxing step is performed. For

each i ∈ Q and t ∈ T , let ¯z

be the values of variables

obtained in the last solution. Each ¯z

is compared

A Matheuristic Approach for Resource Scheduling and Design of a Multi-energy System

455

with a threshold θ ∈ [0, 1] and ﬁxed to one if ¯z

≥ θ,

then M is solved. If the current state is unfeasible,

θ is incremented by a small step α and all the z

for

which ¯z

< θ are unﬁxed; then M is solved again.

This is done until M is feasible or it results unfeasi-

ble for θ = 1, case in which the algorithm ends with

an unfeasible state. Although possible, during com-

putational experiments no unfeasibility state has been

reached and a feasible integer solution is achieved for

the starting value of θ. It is worth noting that for each

commodity k ∈ K that is not purchasable by an exter-

nal supplier, c

= ∞ is assumed to preserve the con-

sistency within (13). Moreover, selecting those com-

modities in later steps help to prevent unfeasibility, as

variables f

are enforced to be zero and cannot bal-

ance the possible increment of requests r

, due to the

selection of device i in further steps of the algorithm.

5 PRELIMINARY EXPERIMENTS

The code was implemented in AMPL v.20180308

(MS VC++ 10.0, 64-bit) and experiments were per-

formed on a Intel

Core i7-7500U 2.90 GHz with

16Gb RAM. All the MILP were solved by IBM

CPLEX

12.8.0.0.

The matheuristic performance were tested on two

instances derived from a real case study and ﬁve ran-

domly generated instances, all based on the energy

hub scheme depicted in §2: 4 energy vectors uni-

formly measured in kWh are taken into account, and

a database of 37 conversion appliances and 18 storage

devices are considered. For conversion technologies,

the lower bound of partialisation is set equal to 30%

of the maximum rated power. Parameters were set to

ε = 0.1, θ = 0.7 and α = 0.1.

5.1 Random Instances

The ﬁve random instances were generated by consid-

ering a time horizon T of 48 time steps, simulating the

requests of commodities within a period of two days.

The values of demands were randomly chosen within

[0, 250] kWh for each time unit t ∈ T , whereas solar

radiation was derived by a real proﬁle.

Table 1 shows the results obtained by the

matheuristic (H ) with respect to the optimal solution

provided by CPLEX after exactly solving the MILP

formulation (M ). In all the instances the percentage

optimality gap is smaller then 1%, with an average

value of 0.65%; then the CPU time required by H is

reasonable and exceeds two minutes only in one case,

while being widely smaller than directly solving M

on the mean.

Table 1: Solutions for random instances.

Instance gap% CPU

CPU

0.93 62.22 807.48

0.79 128.46 1592.91

0.36 83.45 1052.86

0.82 18.60 1481.12

0.37 65.42 2017.73

Average 0.65 71.63 1390.42

5.2 Case Study: Tertiary Urban District

In order to validate the proposed method, the data de-

scribing the needs of an existing urban district is used

to derive realistic instances. Simulated data based on

building energy consumption models is used to ob-

tain the demands of a district composed of tertiary use

buildings, speciﬁcally ten ofﬁces and a school. De-

mands are distinguished in three commodities: elec-

tricity, space heating and cooling. The user demands

are computed as if the buildings were situated in Italy,

speciﬁcally in Rome, for a timespan of one year with

a resolution of one hour. The building models are pro-

totypes developed by the U.S. Department of Energy

(U.S. Department of Energy, a) which are then sim-

ulated by means of the software Energy Plus (U.S.

Department of Energy, b). The solar radiation is also

given as an input by referring to an average year ob-

tained from historical data. While the withdrawal of

electric energy or natural gas from the respective na-

tional infrastructures comes with no limit in terms of

maximum capacity, this is not the same of the radia-

tion resource. The availability of this resource is then

represented by realistic data with the same hourly res-

olution as for the commodities demanded.

5.3 Realistic Instances

Given the scenario depicted as case study, two real-

istic instances are derived. The former instance r

counts 2160 hour as time horizon of simulated de-

mands, corresponding to three months of data that

cross the ﬁrst two quarter of the year and have re-

quirements linked to both the winter and summer pe-

riod. The detail of the simulated demands is shown in

kWh in Figure 2; the availability of the solar resource

for a surface with an inclination of 30

◦

in kW/m

shown in Figure 3. The latter instance r

is based on

the yearly demands aggregated within a daily resolu-

tion, so that the time horizon is translated into 365

days. This is done in order to keep manageable the

size of the MILP, while obtaining a reasonable coarse-

grained description of the urban district annual be-

haviour. To preserve the coherence, all technical pa-

rameters of systems are scaled within r

. While solv-

ICORES 2019 - 8th International Conference on Operations Research and Enterprise Systems

456

ing r

and r

, the relative optimality gap tolerance of

CPLEX was set to 2%.

(a) Electric demand

(b) Heating demand

Figure 2: Three-month energy demands of the district.

Table 2 summarizes the results obtained by the

matheuristic. For the two instances, the total cost Ω

and the CPU time (in seconds) are reported; more-

over, column gap% reports the percentage gap ob-

tained by comparing the heuristic integer solution

with the lower bound resulting by the solution of the

MILP in the ﬁrst step of H ; that is the linear relax-

ation of M restricted only by the integrality of the vari-

ables y

related to the technologies that convert (store)

the ﬁrst selected commodity. In Table 3 the list of ap-

pliances deployed in each solutions is reported, with

subscripts specifying the maximum rated power of the

models (or rated capacity, in case of storage devices).

The matheuristic was able to achieve an integer solu-

Figure 3: Three-month solar radiation.

Table 2: Solutions for realistic instances.

Instance Ω gap% CPU

22185 12.90 5804.09

74356 11.51 432.18

Table 3: Selected devices of realistic instance solutions.

Instance Devices

CHP

140

,EC

120

,EC

160

,HP

160

,HP

240

,PV

2000

HT ES

1000

,CT ES

500

,CT ES

1000

CHP

140

,EC

,AC

,HP

160

,HP

240

,PV

500

,PV

2000

HT ES

500

,CT ES

500

tion for the instance r

in a CPU time of 5804.09 sec-

onds. This solution has a total cost of 22185e, with

a gap% of 12.90%. For the instance r

the heuristic

solution has total cost of 74356efound in 432.18 sec-

onds, with a gap% of 11.51%. Given the dimension

of the instances and the size of the MILP on which the

matheuristic is based (e.g., only the binary variables

in instance r

are 79975), the time required to solve

the problem is reasonable and the values of gap% ob-

tained testify that the maximum distance from the true

optimal solution is not excessive.

6 CONCLUSION AND FUTURE

STEPS

In this paper a MILP formulation is described for the

optimal design of energy hubs, aimed at modeling the

energy systems of modern urban districts. The pro-

posed model wishes to represent the possible interac-

tions among different energy vectors, with the goal of

meeting the multiple commodity demands of users.

However, due to the size of the formulation, it is pos-

sible to directly exploit the MILP only to solve small

instances. In order to efﬁciently solve instances of

more relevant dimension, a matheuristic scheme re-

lying on the MILP has been designed. The strategy

behind the matheuristic is based on ﬁxing the values

A Matheuristic Approach for Resource Scheduling and Design of a Multi-energy System

457

of subsets of binary variables while other integrality

constraints are relaxed, iteratively moving up to a fea-

sible integer solution. Preliminary experiments on a

restricted group of instances show that the matheuris-

tic is able to solve instances with limited time horizon

in restricted time, providing a solution with strict op-

timality gap. Moreover, two instances derived from

real case scenario were solved in reasonable time with

a limited gap from the best lower bound.

The matheuristic presented in this paper can be

seen as a ﬁrst promising step in approaching a chal-

lenging problem as the management of MES, with a

strategy that can be customized in different ways. In-

deed, elements of a MILP can be divided into sub-

blocks following several criteria. Further work will

be dedicated to reﬁne the proposed framework based

on variable separation and ﬁxing, whereas alternative

methodology (e.g., the decomposition into linked sub-

blocks of constraints) will be explored. Concerning

the formulation, future efforts could aim to include

a more accurate modeling of the technology func-

tion, as for example switching from a performance

based on a single constant efﬁciency value to a piece-

wise linear efﬁciency-load curve. Finally, the prob-

lem could be further extended by considering the ex-

istence of multiple separated districts, thus evaluating

the deployment of technologies for connecting differ-

ent districts and modeling the energy exchange.

REFERENCES

Bischi, A., Taccari, L., Martelli, E., Amaldi, E., M. G.,

Silva, P., Campanari, S., and Macchi, E. (2014). A de-

tailed MILP optimization model for combined cool-

ing, heat and power system operation planning. En-

ergy, 74(C):12–26.

BNEF (2017). New Energy Outlook 2017. page 6.

California ISO (2012). What the duck curve tells us about

managing a green grid. Technical report.

Center for Climate and Energy Solutions (2017). Micro-

grids: What Every City Should Know. Technical re-

port.

Elsido, C., Bischi, A., Silva, P., and Martelli, E. (2017).

Two-stage MINLP algorithm for the optimal synthesis

and design of networks of CHP units. Energy, 121.

Enea (2017). Urban Microgrids. Technical Report January.

Escudero, L. and Salmeron, J. (2005). On a Fix-and-Relax

Framework for a Class of Project Scheduling Prob-

lems. J. Ann Oper Res, 140:163–188.

European Commision, J. S. H. (2017). Pv status report

2017. Technical report.

European Commission (2012). Roadmap 2050. Technical

Report April.

ISE, F. (2015). Current and future cost of photovoltaics,

long-term scenarios for market development, system

prices and lcoe of utility-scale pv systems. Technical

report.

Jana, K., Ray, A., Majoumerd, M. M., Assadi, M., and De,

S. (2017). Polygeneration as a future sustainable en-

ergy solution A comprehensive review. Applied En-

ergy, 202:88–111.

ulch, V. (2016). Comparison of electricity storage options

using levelized cost of storage ( LCOS ) method. Ap-

plied Energy, 183:1594–1606.

Li, B., Roche, R., Paire, D., and Miraoui, A. (2017). Sizing

of a stand-alone microgrid considering electric power,

cooling/heating, hydrogen loads and hydrogen storage

degradation. Applied Energy, 205(April):1244–1259.

Mancarella, P. (2014). MES (multi-energy systems): An

overview of concepts and evaluation models. Energy,

65:1–17.

Mehleri, E. D., Sarimveis, H., Markatosa, N. C., and Papa-

georgiou, L. G. (2012). A mathematical programming

approach for optimal design of distributed energy sys-

tems at the neighbourhood level. Energy, 44(1):96–

104.

Mohammadi, M., Noorollahi, Y., Mohammadi-ivatloo, B.,

and Youseﬁ, H. (2017). Energy hub: From a model

to a concept A review. Renewable and Sustainable

Energy Reviews, 80(December 2016):1512–1527.

Omu, A., Choudhary, R., and Boies, A. (2013). Distributed

energy resource system optimisation using mixed in-

teger linear programming. Energy Policy, 61:249–

266.

Sachs, J. and Sawodny, O. (2016). Multi-objective three

stage design optimization for island microgrids. Ap-

plied Energy, 165:789–800.

Singh, B. and Sharma, J. (2017). A review on distributed

generation planning. Renewable and Sustainable En-

ergy Reviews, 76(March):529–544.

Speer, B., Miller, M., Renewable, N., States, U., Schaffer,

W., Gmbh, S. N., Gueran, L., Reuter, A., and Jang,

B. (2015). The Role of Smart Grids in Integrating

Renewable Energy. Technical report.

Triad

o-Aymerich, J., Ferrer-Mart

ı, L., Garc

ıa-Villoria, A.,

and Pastor, R. (2016). MILP-based heuristics for

the design of rural community electriﬁcation projects.

Computers and Operations Research, 71:90–99.

U.S. Department of Energy. Commercial prototype building

models.

U.S. Department of Energy. Energy plus v8.9.0.

ICORES 2019 - 8th International Conference on Operations Research and Enterprise Systems

458