Combined multi-objective optimization and robustness analysis framework for building integrated energy system under uncertainty

Energy Conversion and Management 208 (2020) 112589

Contents lists available at ScienceDirect

Energy Conversion and Management journal homepage: www.elsevier.com/locate/enconman

Combined multi-objective optimization and robustness analysis framework for building integrated energy system under uncertainty

T

⁎

Meng Wanga, Hang Yua, , Rui Jingb, He Liuc, Pengda Chena, Chaoen Lia a

School of Mechanical Engineering, Tongji University, Shanghai, China Department of Chemical Engineering, Imperial College London, London, UK c State Grid Tianjin Information & Telecommunication Company, Tianjin, China b

A R T I C LE I N FO

A B S T R A C T

Keywords: Building integrated energy system Uncertainty Stochastic programming Multi-objective optimization Robustness analysis

The optimal design of building integrated energy system is sensitive to the variation of uncertain parameters. For addressing the tradeoff of uncertainty and optimality-robustness, this study proposes a combined multi-objective optimization and robustness analysis framework for optimal design of building integrated energy system. The proposed framework includes two parts. In the optimization part, on the basis of scenario generation for capturing the uncertainties of renewable energy sources and energy demands, two-stage multi-objective stochastic mixed-integer nonlinear programming is conducted to optimize the system‘s economic and environmental objectives. Two decision-making methods are introduced to identify the final optimum solution from the obtained Pareto frontier. In the robustness-analysis part, a combined Monte Carlo simulation and optimization method is implemented to verify the robustness of the optimal solutions. The two parts of the framework are integrated to investigate the case of a hotel in Beijing, China. The results indicate that compared with the stochastic model, the deterministic model underestimates the annual total cost. Achieving economic and environmental optimum is conflicting and needs a trade-off through decision making. Moreover, in the robustness analysis, an acceptable robustness value is identified, considering both the selected objectives and the operation constraints’ probability of failure. The Shannon-entropy-based final optimum solution exhibits the best comprehensive performance, with an annual total cost of $695 × 103/year, an annual carbon emissions of 2100 tons/year, and an 8.81% probability of failure.

1. Introduction The integrated energy system (IES) is composed of several key elements, including production, conversion, transmission, and consumption, in which multiple types of energy are coupled together and converted among each energy carrier. Over the past decade, the IES has attracted increasing attention from researchers seeking to address the energy supply-demand issue in urban areas, as it has several advantages, including renewable-energy access, low carbon emissions, good flexibility, low costs (with a feed-in tariff), and the development of distributed generation. The building integrated energy system (BIES) is commonly employed to fully utilize local renewable energy resources and efficiently recover waste heat, as well as for power generation. Studies are increasingly focusing on optimizing the design and dispatch of the BIES. Research efforts have been directed toward modeling the IES for individual buildings [1], particularly for public buildings such as

⁎

hospitals and commercial buildings with a relatively large energy demand. Razmi et al. [2] investigated the operation and optimization of combined cooling, heating, and power (CCHP) systems for buildings and concluded that it is ecofriendly and energy-efficient. Pagliarini et al. [3] examined the feasibility of CCHP systems to be used in a hospital through economic and exergy analysis. Jing et al. [4] selected a hotel in Beijing, China to model the BIES via mixed-integer nonlinear programming (MINLP). The results indicated that the fuel-cell-based CCHP system has good performance for buildings. In addition to the determination of a suitable objective and feasibility analysis, the decision making of multi-objective optimization is a significant issue. Currently, single-objective optimization cannot sufficiently characterize the comprehensive performance of IES systems. Researchers have introduced two or more indicators to evaluate the system characteristics. The economic indicator is the most extensively utilized optimization objective. Soheyli et al. [5] developed a renewable assisted CCHP system model in which multi-objective particle swarm

Corresponding author at: Kaiwu Building, Tongji University, 4800 Cao’an Road, Jiading District, Shanghai, China. E-mail addresses: [email protected], [email protected] (H. Yu).

https://doi.org/10.1016/j.enconman.2020.112589 Received 6 November 2019; Received in revised form 4 February 2020; Accepted 7 February 2020 0196-8904/ © 2020 Elsevier Ltd. All rights reserved.

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

Nomenclature

Greek symbols

Abbreviations

α η θ π Ψ ω

ATC ACE BIES CCHP CRF CV DP FOS MINLP SP SRI STD TOPSIS

annual total cost annual carbon emissions building integrated energy system combined cooling, heating, and power capital recovery factor coefficient of variation deterministic programming final optimum solution mixed-integer nonlinear programming stochastic programming solar-radiation intensity standard deviation technique for order preference by similarity to an ideal solution

Subscripts/superscript ac b cap c-in c-out cool ec ED ele ex fs h heat hp im max maint mgt min NG norm obj pv SE s st st-in st-out t wt

Symbols A AM CAP CAPEX d E EDi+ EDifij l OPEX OBJ PL prob Q r SECE T UC v

coefficient efficiency on/off status emission factor start limit variable weight

area air mass install capacity (kW) capital cost deviation index electrical power (kW) Euclidean distance to ideal point Euclidean distance to nadir point jth objective-function value for ith point project life (year) operation cost objective ratio of partial load probability thermal energy (kW) interest rate scenario‘s carbon emissions temperature (K) unit cost wind speed (m/s)

absorption chiller boiler capacity cut-in speed cut-out speed cooling energy electrical chiller Euclidean distance based electrical energy electricity export feasible solution time step (hour) heating energy heat pump electricity import maximum maintenance cost micro gas turbine minimum natural gas normalization objective photovoltaic panel Shannon entropy based stochastic scenario heat storage heat storage charge heat storage discharge each energy technology wind turbine

based CHP system. Their case study indicated that the ADMM algorithm has far better convergence than the particle swarm optimization algorithm for the operation optimization problem. There are three main types of optimization approaches with different computing frameworks: classical strategies, metaheuristic strategies, and distributed computing. The classical strategies employ a rigorous theory to ensure a global optimal and unique solution. The metaheuristic strategies usually have excellent generalization ability; thus, almost all optimization problems can be solved using heuristic algorithms. However, owing to the lack of a mathematical theory, the “premature” problem (i.e., being trapped at a local optimum) cannot always be avoided, and the computational cost is unpredictable. Additionally, the uncertainty and reliability issues for IES optimal planning have attracted considerable attention. Razmi et al. [12] performed an exergoeconomic analysis of cogeneration considering the reliability and availability via the Markov method and then provided a realistic and valid payback time. Mavromatidis et al. [13] used Monte Carlo simulation to trace the source of uncertainty and performed a global sensitivity analysis to evaluate the robustness of the system. On this basis, Mavromatidis et al. [14] proposed a single-objective twostage stochastic programming (SP) model to optimize the design and operation of a multi-building energy system. The results indicated that

optimization was employed to perform multi-objective optimization. Wei et al. [6] adopted non-dominated sorting genetic algorithm-2 (NSGA-2) to optimize the operation of a CCHP system, considering the energy conservation and total cost as two objectives for bi-objective optimization. Somma et al. [7] developed a mixed-integer linear programming (MILP) model to optimize the design and dispatch of a gasfired internal combustion engine-based CCHP system. In their model, the daily costs of the energy carriers and the carbon emissions were combined into a comprehensive objective function. Ju et al. [8] employed the goal programming method for multi-objective optimization of an MILP-based CCHP model. By restructuring the single objective, the model can be optimized with the energy utilization and operation cost, as well as the carbon emissions. Jing et al. [9] used the epsilonconstraints method to convert multi-objective optimization into multiple times single-objective optimization for obtaining a series of nondominated solutions, by changing the objective function into several constraints. Zeng et al. [10] used the weighted-sum method to address multi-objective optimization problems. By assigning the economic, environmental, and technical indicators the same weighted value, the original problem was converted into a single-objective optimization problem. Tran et al. [11] used the alternating direction method of multipliers (ADMM) to optimize the whole-day operation of a fuel-cell2

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

the deterministic model may underestimate the total cost. Lv et al. [15] investigated an energy-water nexus system by using the interval fuzzy chance-constrained method. The results of a case study for Hebei province in China indicated that the proposed programming method could address the shortage of the energy supply with the reduced installed capacity of coal-fired power stations. Akbari et al. [16] conducted robust optimization for the design of distributed energy systems under uncertainty. The concept of the conservatism level was introduced to express the decision makers’ potential choices, which were strongly correlated with the sizes of the installed equipment and the heating pipeline. Gabrielli et al. [17] performed the robust and optimal design of a multi-energy system through robust scenario generation with limited information. The results indicated that the optimization improved the robustness and reduced the overall costs. Generally, sampling-based scenario generation is widely used to express the uncertainty, which may cause the scenario to strongly depend on the probability distribution. Moreover, because the number of dimensions for a stochastic model is significantly larger than that for a deterministic model, the computational cost of a nonlinear stochastic model may be unacceptable for real-world project designers, particularly when they face multi-objective optimization problems. Although many studies have been performed on approaches for optimizing the design and dispatch of the BIES with various economic and environmental objectives, a few challenges and issues concerning multi-objective optimization of the BIES with multiple uncertainties still exist, according to the aforementioned literature:

Fig. 1. Proposed BIES design optimization method involving SP.

(1) In previous studies, researchers tended to figure out effective but simple approaches to capture the uncertainty of parameters, i.e., generating various scenarios through sampling or developing different frameworks of the programming model under uncertainty. However, few studies have focused on multi-objective optimization under uncertainty. (2) Compared with the IES, which is optimized at a large scale, the BIES always encounters variability and uncertainty in the time-series demands, because of the lack of complementarity among various categories of buildings. Hence, precise formulations are needed to accurately model the off-design or part-load performance of energy technologies, which makes it difficult to obtain a Pareto frontier with both good diversity and a uniform distribution. (3) In previous studies, researchers typically converted multi-objective optimization problems into single-objective problems through the weighted-sum method or by adding constraints. This may result in an insufficient number of alternatives for decision makers. Thus, an efficient robustness-analysis tool for the multi-objective optimization of the BIES design is needed.

(AUGMECON2). The final optimum solutions (FOSs) are obtained via two decision-making methods. (3) A robustness-analysis approach with low computational cost is proposed to verify the robustness of the results posteriorly by combining Monte Carlo simulation and operation optimization. Additionally, the robustness of the capacity configuration is evaluated for various optimal solutions under the stochastic variations at the operation stage. The remainder of this paper is organized as follows. Section 2 introduces the framework and methods used in the present study. Section 3 presents the details of the BIES model. Section 4 describes the case study performed to evaluate the proposed approach. Section 5 presents discussions and analysis of the results. Conclusions are drawn in Section 6. 2. Framework and methods Fig. 1 illustrates the proposed framework, which can be divided into five major steps: (1) input parameters, (2) uncertainty quantification, (3) model development, (4) multi-objective optimization and decision making, and (5) robustness analysis. In step (1), the parameters are classified into two groups: deterministic and uncertainty parameters. The deterministic parameters can be directly introduced into the BIES SP model, which is established in step (3). For the uncertainty parameters, the scenario generation in step (2) is performed to quantify the uncertainty of the hourly meteorological parameters and the energy demands by constructing a scenario tree, which can capture their stochastic characteristics via the clustering method. Additionally, distribution fitting is performed to simulate the time-series variation for different distribution functions. In step (3), a two-stage multi-objective SP model is established and formulated as an MINLP problem. In step (4), the approach for multi-objective tradeoff optimization and decision making is employed to obtain the Pareto frontier and identify the FOS, whose robustness is verified via the robustness-analysis approach in step (5). The aforementioned steps are explained in detail in the following subsections.

Although the robust optimization could optimize the uncertainty model considering optimality and robustness, it focuses on optimizing the worst-case objective function, which makes the optimal design overconservative. In the present study, to balance the optimality and robustness, the present study established an efficient and nonconservative optimization framework for the BIES under uncertainty, while ensuring that it is solvable by personal computers, making the proposed framework easy to implement by system designers in practice. The major contributions of the present study are as follows: (1) A dimensionality reduction and clustering approach is developed for capturing the characteristics of each uncertainty parameter (i.e., building energy demand, wind speed, and solar radiation) by generating a probabilistic stochastic scenario-based binary tree. (2) A two-stage multi-objective SP-based MINLP model is established for optimizing the design and dispatch of the BIES system. The Pareto frontier with economic and environmental objectives is obtained by combining the state-of-the-art MINLP solver and the improved version of the augmented ε-constraint method 3

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

stochastic scenarios with associated probabilities. The probability of each final scenario can be determined by multiplying the probabilities of the clusters in the layers (their sum must be equal to 1).

2.1. Uncertainty quantification The present SP model involves two types of uncertainty parameters: hourly meteorological data and energy demands. To capture their uncertainties, the scenario-generation approach is employed. The distribution-fitting method is utilized for quantitative description of their variations.

2.1.2. Probability distribution fitting On the basis of the assumption that the uncertainty parameters (including the meteorological data and energy demand) are mutually independent, probability distribution functions are introduced to describe the variation of the uncertainty parameter. The Weibull distribution, which is widely used to fit and forecast wind data [20], is selected for the wind-speed parameter. The SRI is assumed to be normally distributed in the time-series data for a day. The uncertainty of the electrical demand is assumed to follow a triangular distribution [21], and the heating and cooling demands both follow a normal distribution [22]. In the present study, the distribution fitting is conducted by calling the built-in module of dfittool (Distribution Fitter app) in MATLAB.

2.1.1. Scenario generation Owing to the large amount of hourly data, the computational cost and model dimension may be excessive if all the data are considered. Therefore, a scenario-generation method is necessary to convert the data into a series of stochastic scenarios with associated probabilities for characterizing their temporal uncertainty and reducing the computational cost. The structure of the scenario tree is illustrated in Fig. 2. For each season, the binary-tree is constructed layer-by-layer. Each node in the scenario tree is described that every time-series data of the uncertainty parameter has two scenario clusters: high-value and lowvalue. The primary structure of the scenario tree is established by different combinations of one classification layer (the season) and three clustering layers (the solar-radiation intensity (SRI), wind speed, and energy demand). Moreover, for minimizing the gap between the historical data and the clustering results and ensuring that none of the clusters is empty, the k-medoids clustering method is employed. The detailed procedure used to generate stochastic scenarios can be found in Ref. [18,19]. In this study, different types of energy demand (electrical, heating, and cooling) are considered, and each clustering vector of the energy demand had 72 dimensions, which could result in a high computational cost and the non-convergence problem. Therefore, principal component analysis is used to reduce the number of dimensions of the vectors. Thus far, 24 representative stochastic scenarios have been obtained to describe uncertainty and fluctuation. However, underestimation of the installed capacity is inevitable, as the scenarios may not cover the peak energy demand. Hence, three additional scenarios (characterizing the peak load) are added to the scenario tree. Finally, there are 27

2.2. Stochastic programming approach The two-stage SP approach, which is widely used in supply chain modeling and transportation planning to handle the uncertainty of traffic flows or weather variations, is employed to address the uncertainty problem. A common formulation is as follows [23]:

min cT x + Eξ (min{qT y| W y = h − Tx , y ⩾ 0}) s.t. Ax ⩾ b

(1a)

Eq. (1a) can be simplified into Eq. (1b), which is more applicable to optimal design and dispatch in the proposed model.

(1b) As indicated by Eq. (1b), the first-stage (design) decision to select and size the components of the energy technologies is composed of

Fig. 2. Scenario-tree structure for SP. 4

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

normalization, the technique for order preference by similarity to ideal solution (TOPSIS) is conducted. “Ideal” and “nadir” correspond to the absolute best (but nonexistent) point and the most undesirable point. respectively. On this basis, indices can be calculated using the Euclidean (geometric) distance, as follows [25,26]:

design variables and deterministic parameters, while the operation variables and uncertainty parameters constitute the second-stage (operation) decision to optimize the dispatch. The design variables include the continuous variable denoting the installed capacity, while the energy input and output of the installed components in the time series are the operation variables in which the uncertainty is embodied. The values of the design variables are fixed over the time horizon, whereas the values of the operation variables change at each time interval. Because these two groups of decision variables must be optimized simultaneously, the entire two-stage SP model is established as an MINLP model.

n

EDi + (i−) =

∑ (fij

− f jideal (nadir) ) (4)

j=1

The solution with the maximum relative distance is desired and is labeled as the FOS.

2.3. Multi-objective optimization and decision making

EDi − ED ⎞ iFOS = arg max ⎛ ⎝ EDi + + EDi − ⎠ ⎜

In this section, a multi-objective optimization approach is introduced to optimize two conflicting objectives. Then, two absolutely objective decision-making methods are described.

(

S

S

Sp rp

SEj = −

1 ln(m)

m

⎛

∑⎜ i=1

⎝

f ijnorm m ∑i = 1 f ijnorm

norm

× ln

⎛ f ij ⎞⎞ ⎜ ∑m f norm ⎟ ⎟ ⎝ i = 1 ij ⎠⎠

ωj =

1 − SEj n

∑ j = 1 (1 − SEj )

n

SE iFOS = arg min

⎛ ⎞ f × ωj ⎜ ∑ ij ⎟ j = 1 ⎝ ⎠

(8)

(2) 2.4. Robustness-analysis approach In this part, a simulation-optimization based Monte Carlo method is employed to investigate the robustness of the optimal solutions, which

2.3.2. Decision making methods Although the Pareto frontier always provides superfluous alternatives with various values of the given objectives, at least one suitable method is needed to address the multiple-criteria decision-making problem. In this study, two absolutely objective decision-making methods, i.e., “Shannon-entropy-based” and “Euclidean-distancebased” are employed to identify the FOS. (1) Euclidean-distance-based method In view of the various scales and dimensions for objective values, normalizations are inevitable in the first step. Eq. (3) provides the most common approach for converting all objectives into several dimensionless scalars.

fij − min(fij ) max(fij ) − min(fij )

(7)

Finally, the desired solution is obtained by minimizing the sum of the products of the objective values and their corresponding entropy weight values:

) ⎞⎠

Fig. 3 illustrates the generic form of the Pareto frontier in the search space with bi-objective simultaneous minimization. No solution is observed at the bottom-left of the Pareto frontier, indicating that no solution has two less objective values than the non-dominated solutions in the Pareto frontier. Moreover, the “best of the best” solution must be selected from the Pareto frontier, which necessitates a tradeoff.

f ijnorm =

(6)

Second, the entropy weight value of each objective is calculated:

f3 (x ) − S3 = e3 ... fp (x ) − Sp = ep x ∈ S, Si ∈ R+

(5)

(2) Shannon-entropy-based method The Shannon-entropy-based method exhibits good performance for quantifying the uncertainty of information sources. Each objective is assigned a weight value indicating its importance, according to the distribution of all the solutions in the Pareto frontier. Generally, an indicator with lower information entropy has greater variation and higher importance, and its weight should be larger. The procedure of the Shannon-entropy-based decision-making method is as follows [27]. First, the Shannon entropy of each objective is calculated:

2.3.1. Multi-objective optimization approach Typically, in multi-objective optimization, two or more objective functions must be optimized simultaneously. In this study, Pareto-optimality is introduced to deal with bi-objective tradeoff optimization, which uses the Pareto frontier composed of a series of non-dominated solutions to characterize the economic-environmental performance. To obtain the Pareto frontier, the AUGMECON2 is introduced, which can convert the multi-objective optimization into multiple times and constraints-added single-objective optimization, as indicated by Eq. (2). AUGMECON2 enhances the conventional ε-constraint method by addressing weaknesses such as the guarantee of Pareto-optimality and the increased solution time. Details and applications of AUGMECON2 can be found in the model libraries on the GAMS official website and in Ref. [24].

max ⎛f1 (x ) + eps × r 2 + 10−1 × r 3 +⋯+10−(p − 2) × 2 3 ⎝ s. t . f2 (x ) − S2 = e2

⎟

(3)

where subscript i denotes the number of the optimal solution, and the subscript j indicates the jth objective to be optimized. After the

Fig. 3. Diagrammatic sketch of Pareto-optimality and decision making. 5

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

can be regarded as a posteriori analysis method performed after the decision-making process. 2.4.1. Monte Carlo simulation Owing to the various scenarios and their associated probabilities, the optimal solutions for single-objective optimization, i.e., the two endpoints of the Pareto frontier, may be sensitive to the operation parameters, necessitating robustness analysis and validation. Monte Carlo simulation via random sampling is a typical method for evaluating robustness. However, when this method is introduced directly into the BIES model, the changing parameters (particularly the energy demand) may cause the optimization result for the dispatch strategy to not satisfy the hourly energy demand; that is, most simulations exhibit “failure of energy balance constraints” without any result. Additionally, in practice, the dispatch strategy cannot be identical to that for the assumed scenarios in the planning model. Accordingly, a Monte Carlo simulation is performed to simulate the variation and fluctuation of each uncertainty parameter at the operation stage, which can provide many results to indicate the sensitivity of system. The random sampling and probability distribution fitting methods are described in detail in Section 2.1.2.

Fig. 5. Superstructure of the proposed renewable assisted CCHP-based BIES.

constraints are analyzed via a statistical method. 3. Model development In accordance with the proposed technical roadmap, the overall superstructure of a renewable assisted CCHP-based BIES system is presented in Fig. 5. A micro-gas turbine (MGT) fueled by natural gas is selected as the prime mover, and the generated electricity is combined with that of the photovoltaic (PV) panel and wind turbine (WT), as well as the utility grid, to satisfy the electrical demand. Moreover, the localized heating network utilizes waste heat produced by the prime mover, gas boiler, air-source heat pump (HP), and heat storage (ST), and the localized cooling network comprises an electrical chiller (EC) and an absorption chiller (AC).

2.4.2. Robustness analysis This paper presents an improved robustness-analysis method that combines Monte Carlo simulation and single-objective optimization to verify the robustness of the optimal solutions under operation-parameter variations, as shown in Fig. 4. First, time-series random sampling of the uncertainty parameters in different scenarios is performed according to their probability distribution functions. Meanwhile, the design variables of the optimal solutions are fixed and converted into the given parameters, and the operation variables are cleared and re-optimized. The sampling results and fixed design variables are introduced into the proposed SP model to optimize the operation variables. Moreover, the priori sampling and posteriori optimization are repeated until the final iteration is complete. Finally, the distribution of the objective-function values and the probability of failure of the operation

3.1. Objective functions In the present study, the annual total cost (ATC) and annual carbon emissions (ACE) are employed to characterize the comprehensive performance of the proposed BIES system. These two objectives are optimized simultaneously.

Fig. 4. Procedure of the simulation-optimization-based robustness-analysis approach. 6

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

3.1.1. Economic objective The ATC is an economic index commonly used in the planning and evaluation of energy systems, and it can be divided into two parts: the capital cost and the sum of all the scenarios’ operation costs. The CAPEX includes the cost of all the devices and their auxiliary facilities, and the OPEX comprises the fuel cost, maintenance cost, electricity purchasing cost, and profit from the feed-in tariff [28].

ATC = CAPEX +

∑ OPEXs × probs

Edemand, s, h = Epv, s, h + E wt, s, h + Emgt, s, h + Eim, s, h − Eec, s, h − Ehp, s, h − Eex, s, h

heat heat heat heat heat heat heat Qdemand, s, h ⩽ Qmgt, s, h + Q b, s, h + Qst - out, s, h + Q hp, s, h − Qst - in, s, h − Qac, s, h

CAPEX = ∑ CAPt × UCtcap × CRF (amortization of capital cost)

cool cool cool Qdemand, s, h ⩽ Qac, s, h + Qec, s, h

OPEXs = ∑ (FCs, h + MCs, h + EPCs, h − FIPs, h) (operation cost)

(9c)

h

∀ s, h

(11c)

the subscripts h and s indicate the hour and scenario, respectively. In Eq. (12), the binary variables θim and θex represent the statuses for electricity purchasing and feed-in to the grid, and their sum is ≤1 [4].

(9b)

t

(11b)

∀ s, h

(9a)

s

(11a)

∀ s, h

max 0 ⩽ Eex, s, h ⩽ θex, s, h × Eex ∀ s, h

(12a)

max 0 ⩽ Eim, s, h ⩽ θim, s, h × Eim

(12b)

θex, s, h + θim, s, h ⩽ 1

FC s, h

∀ s, h

∀ s, h

(12c)

= (NGmgt, s, h + NG b, s, h) × UCNG (fuel cost)

3.3. Prime mover constraints

(9d)

The MGT is the core of the BIES model, and it satisfies a large amount of the electrical and heating demand simultaneously. Hence, several different types of operation constraints are introduced to describe the off-design performance, including general constraints and specific constraints.

MCs, h heat cool maint = (Epv/mgt/wt, s, h + Q hp/b/st, s, h + Qec/ac, s, h ) × UCt

(9e)

(maintencace cost)

EPCs, h = Eim, s, h × UCim

3.3.1. General operation constraints The general constraints provide the relationships (i.e., electrical and heating efficiencies) between the hourly MGT output and the consumption rate of natural gas. To avoid low-efficiency operation, the prime mover must be started when the electric load is > 30% of the capacity [28].

(electricity purchasing (9f)

cost)

FIPs, h = Eex, s, h × UCex

=

CRF

Emgt, s, h ⩽ CAPmgt

(9g)

(feed - in tariff)

Emgt, s, h =

r × (1 + r )l

where CRF represents the capital recovery factor.

∑ SECEs×probs ∑ Eim,s,h + πNG × ∑ (NGmgt,s,h + NG b,s,h) h

(13b)

heat heat Qmgt, s, h = ηmgt, s, h × NGmgt, s, h

∀ s, h

(13c)

h

∀ s, h

(13d)

Additionally, the hourly electrical and heating efficiencies can both be calculated using the part load ratio (PL), which can be fitted as a quadratic nonlinear equation with coefficients (α) that are set as −0.2, 0.4, and 0.1 for electricity generation and −0.101, 0.202, and 0.458 for heating [29,30]. ele/heat ele/heat ηmgt, × (PL mgt, s, h )2 + a2ele/heat × (PL mgt, s, h) + a3ele/heat s, h = a1

PL mgt, s, h = (10b)

where πgrid and πNG represent the carbon emission factors of grid electricity and natural gas, respectively.

Emgt, s, h CAPmgt

∀ (14a)

s, h

(10a)

s

SECEs = πgrid ×

∀ s, h

Emgt, s, h ⩾ 30% × CAPmgt × θmgt, s, h

3.1.2. Environmental objective The ACE index is the other objective, corresponding to the environmental performance. In contrast to the economic objective, the carbon emissions at the operation stage are considered; i.e., the ACE is equal to the sum of all the scenarios’ carbon emissions multiplied by the associated probabilities. Hence, the ACE can be expressed as follows [14]:

ACE =

(13a)

× NGmgt, s, h

(9h)

(1 + r )l − 1

ele ηmgt, s, h

∀ s, h

∀ s, h

(14b)

3.3.2. Specific operation constraints The ramping constraints (Eq. (15)) are introduced to avoid drastic variations of the operation [31]. The production variation among adjacent times must be limited to 50% of the installed capacity.

3.2. Energy balance

|Emgt, s, h + 1 − Emgt, s, h| ⩽ 50% × CAPmgt

According to the schematic of the BIES model in Fig. 5, the energybalance formulations include the electrical balance, heating balance, and cooling balance. The energy supply must be greater than or equal to the energy demand. For the electrical balance in each scenario, the electrical demand is covered by the MGT, the wind turbine, the PV panel, and the electricity imported from the grid, subtracted by the electrical consumption of the electrical chiller and the air-source heat pump. Additionally, the surplus electricity exported to the grid must be deducted.

∀ s, h

(15)

Moreover, start-stop constraints (Eq. (16)) are proposed to ensure that the operation is started at most once per day [9], which ensures that the prime mover can operate for the entire project.

∑ ψmgt,s,h ⩽ 1

∀ s, h (16a)

h

θmgt, s, h − θmgt, s, h − 1 ⩽ ψmgt, s, h ⩽ 1 − θmgt, s, h − 1 7

∀ s, h

(16b)

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

ψmgt, s, h ⩽ θmgt, s, h

∀ s, h

cool/heat Qec/ac/hp/b, s, h ⩽ θec/ac/hp/b, s, h × CAPec/ac/hp/b ∀ s , h

(16c)

3.4. Renewable energy module constraints In the BIES system, each component should be compact and spacesaving; thus, a vertical-axis wind turbine with a small rotation radius and low cut-in wind speed is introduced and installed on the roof. The PV panel and wind turbine must share the limited available roof space, as indicated by Eq. (17).

Awt + Apv ⩽ Aroof

0

otherwise

=

⎜

a2pv

+

⎟

⎡ pv ⎛ Ts, h ⎞ ⎢1 + a4 × T ⎝ 0 ⎠ ⎣ ⎜

a2pv

⎟

a3pv

+ a5pv ×

∀ s, h

(23c)

(24)

A case study involving a 24 h-operational hotel in Beijing is conducted to implement and verify the proposed framework for optimizing the BIES system under uncertainty.

AMs, h ⎤ ∀ s, h AM0 ⎥ ⎦

4.1. Overview (19a)

Epv, s, h = Apv × ηpv, s, h × SRIs, h ∀ s, h

Beijing is located in north China. It is in the “Cold” climate zone of China, which has distinctive seasons, i.e., a cold winter, a hot summer, and a transition season. Meteorological data, including the SRI and wind speed, are acquired from the China Integrated Meteorological Information Sharing System database. To determine the building‘s energy demand (i.e., electricity, cooling, and heating) with an hourly resolution, the building performance simulation software Designer’s Simulation Toolkit (DeST) [36] was employed to model the hotel. According to the structural drawings of the building, DeST simulated the energy demand by considering the building thermal processes and human behaviors. Details regarding the building are presented in Table 1.

(19b)

where ηpv represents the efficiency of the PV panel, and AM and T represent the air mass and ambient temperature, respectively. In accordance with [4,32], the coefficients (α) are set as 0.2820, 0.3967, −0.4473, −0.093, and 0.1601, and the values of SRI0, AM0, and T0 are set as 1000 W/m2, 25 °C, and 1.5, respectively. 3.5. Energy-storage constraints Owing to the economic attractiveness of thermal-storage technologies, a storage tank is introduced. Eq. (20) shows that the stored heat at time h is equal to that at time h-1, plus the heat charged at time h, minus the heat discharged at time h. To avoid simultaneous charging and discharging, the binary variables θst-in and θst-out are utilized to express the charging and discharging statuses, respectively [29,33]. heat max Qst, ∀ s, h s, h ⩽ CAPst

Qstheat - out, s, h ηst - out

max 0 ⩽ Qstheat ∀ s, h - in/st - out, s, h ⩽ θst - in/st - out, s, h × CAPst

θst - in, s, h + θst - out, s, h ⩽ 1 ∀ s, h

CAPac

(23b)

4. Case study

⎟

heat heat heat Qst, s, h = ηst × Qst, s, h − 1 + ηst - in × Qst - in, s, h −

cool Qac, s, h

∀ s, h

where the subscript t indicates the energy-conversion technology.

SRIs, h ⎞ ⎤ ×⎛ ⎥ ⎝ SRI0 ⎠ ⎦ ⎜

(23a)

CAPt ⩽ CAPtmax

(18)

ηpv, s, h ⎡ SRIs, h ⎞ × ⎢⎛ SRI0 ⎠ ⎣⎝

cool heat Qac, s, h = θac, s, h × ηac,s, h × Qac, s, h ∀ s , h

where NG represents the energy content of natural gas. Finally, Eq. (24) suggests that each energy-conversion technology has an upper bound for the available capacity.

where ρair represents the air density. The PV module can be expressed as follows:

a1pv

(22b)

PLac, s, h =

v c - in ⩽ vs, h ⩽

∀ s, h

Q b,heat s, h = θ b, s, h × η b × NG b, s, h ∀ s , h

× PLac, s, h + 0.8114

v c - out

⎨ ⎩

(22a)

= 0.0001 × (PLac, s, h )3 − 0.0061 × (PLac, s, h )2 + 0.0578

The module of the wind turbine is described by Eq. (18). The wind turbine only operates when the real-time wind speed is in the range of the cut-in speed to the cut-out speed [5].

E wt, s, h =

cool/heat Qec/hp, s, h = θec/hp, s, h × ηec/hp × Eec/hp, s, h ∀ s , h

ηac,s, h

(17)

⎧ ηwt × 0.5 × ρair × Ablade × min(v rated, vs, h)3

(21)

∀ s, h

4.2. Baseline condition

(20a)

The technical, economic, and environmental parameters are strongly correlated to the optimal results; thus, it is necessary to set these parameters carefully in accordance with state-of-the-art

(20b)

Table 1 Parameters for the building and air-conditioning system. Items

(20c) (20d)

where ηst, ηst-in, and ηst-out represent the storage efficiency of the intertemporal loss, the heat charging efficiency, and the discharging efficiency. 3.6. Other associated constraints The binary variable θ is introduced to indicate the on/off status for other energy-conversion technologies. The electrical chiller, the heat pump, and the gas boiler have constant efficiencies, as indicated by Eqs. (21) and (22). Moreover, the part-load constraints can be used for the absorption chiller [34,35], as indicated by Eq. (23). 8

Unit

Value

2

Total floor area Storeys Heat transfer coefficients of roof Heat transfer coefficients of external wall Heat transfer coefficients of window Window-wall ratio Winter Summer Transition season

m – W/(m2·K) W/(m2·K)

7300 5 0.525 0.593

W/(m2·K) – date date date

Operating time of air-conditioning system Shading coefficient of windows Indoor design temperature

hours

2.4 0.40 23th Oct.-24th Apr. 9th Jun.-3th Sep. 25th Apr.-8th Jun. & 4th Sep.-22th Oct. 24

– °C

0.7 22–25

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

5.1.2. Influence of uncertainty parameters on system performance A comparison of the Pareto frontiers of the DP and SP models can elucidate the influence of the uncertainty parameters on the selected objectives, as the variation of the uncertainty parameters is generally unpredictable and uncontrollable. As shown in Fig. 7, the non-dominated solutions of the DP model are generally better than those of the SP model, and the ranges of the two objective functions are wider for the DP model than for the SP model. DL11 represents the gap between the optimal ATCs of the DP and SP models in the optimization of the economic objective, and DL12 represents the gap between their corresponding ACEs. For economic optimization, the optimal ATC for the DP model is approximately 5% lower than that for the SP model. This is because the ACE depends on the operation stage rather than the design stage and reducing the number of constraints of the DP model makes the operation scheme more flexible with regard to reducing the operation cost at the expense of increased carbon emissions. The gap DL21 exists and is significantly smaller than DL11, owing to the underestimation of the devices‘ installed capacities. Similarly, the gap DL22 is attributed to the more feasible solutions with fewer constraints. Hence, the SP model has a smaller search space than the DP model, whereas its computing cost increases rapidly. For the DP model, the ranges of the non-dominated solutions for the ATC and ACE are $651000/year to $729000/year and 1674 to 2575 tons/year, respectively. For the SP model, the ranges of the non-dominated solutions for the ATC and ACE are $680000/year to $737000/year and 1740 to 2419 ton/year, respectively. To identify the tradeoff between the ATC and ACE of the BIES model, the fitted curve equations derived from the Pareto frontier are obtained via polynomial fitting, as shown in Table 3.

information. Owing to the mature technologies of the electrical chiller and gas boiler, the efficiency or COP is considered as a constant value, which is not be affected by the part-load ratio. The project lifetime is conservatively set as 15 years, and the grid electricity emission factor is set as 1.04 kgCO2/kWh, in accordance with the power source of the North China grid [31,35]. The other key parameters are presented in Table 2. According to the description of the model development in Section 3, the programming model can be simplified by significantly reducing the number of stochastic scenarios to three typical days representing the winter, summer, and transition season, whereby the deterministic programming (DP) model can be established as a baseline condition. The data (including the energy demand, SRI, wind speed, and time-ofuse price of energy carriers) are presented in Fig. 6. 4.3. Modeling environment The modeling, optimization, and simulation are implemented using an ordinary computing platform. The scenario generation and Monte Carlo simulation are conducted in MATLAB R2018b, and the SP model is established in the Lindo optimizer-based GAMS, which has a MATLAB-GAMS interface called GDXMRW. 5. Results and discussions As described in the previous sections, the BIES system is established as an MINLP-based programming model in GAMS. In this section, the results and discussions are analyzed and divided into five parts: optimization, decision making, capacity configuration, energy consumption shares breakdown, and robustness analysis.

5.2. Decision making 5.1. Multi-objective optimization In the proposed approach, two different types of decision making based on Euclidean distance and Shannon entropy are employed to identify the FOS from the Pareto frontier, addressing the tradeoff between the two conflicting objectives. Fig. 8 illustrates the similarities and differences between these two decision-making methods. For the Euclidean-distance-based method, the maximum relative distance to the nadir point can be defined as the optimum principle, which is calculated using the normalized distances to the ideal point and nadir point. In contrast, the Shannon-entropy-based method uses the concept of information entropy to describe the degree of the chaos of nondominated solutions; then, the solution with the smallest value of the index is considered as the FOS. Although these two approaches are absolutely objective, normalization is necessary to avoid potential issues caused by differences in the scale or units. Because these two Pareto frontiers have asymmetric shapes with variable curvature, the two different solutions are identified and denoted as B2 and B3, respectively. Because there are two different but similar FOSs, a new assessment parameter called the deviation index [41] is introduced to choose between B2 and B3, which is defined as follows:

To compare the DP model with the SP model, both models are optimized to obtain their Pareto frontiers via the approach described in Section 2.3.1. 5.1.1. Pareto frontiers Fig. 7 presents the Pareto frontiers of the two models, where a series of non-dominated solutions represent the optimal solutions. In the search space, a non-achievable solution exists with the minimum ATC and ACE, indicating that there is a tradeoff between these two objectives [40]. Additionally, the set of non-dominated solutions composes a curved line convex to the ideal point, indicating that improving one objective function degrades the other. Moreover, according to the foregoing decision-making methods, the FOS is identified and marked by the following roles: one is that the letters A and B represent the DP and SP models, respectively; the other is that the numbers 1 and 4 indicate the minimum ATC solution and minimum ACE solution, respectively, while the numbers 2 and 3 indicate the FOSs for the Shannon-entropy-based method and the Euclidean-distance-based method, respectively. Table 2a Technical parameters for the case study [4,9,29,35,37]. Items

Parameters

Value

Items

Parameters

Value

Gas boiler

Efficiency Maximum capacity Cut-in speed Cut-out speed Generator efficiency Maximum capacity Utilization ratio Maximum capacity

0.85 1600 kWh 1.3 m/s 45 m/s 0.53 1200 kWe 0.9 1200 kWc

Electrical chiller

COP Maximum capacity COP Maximum capacity Storage efficiency Heat charge efficiency Heat discharge efficiency Maximum capacity

4 1800 kWc 3.5 800 kWh 0.95 0.95 0.95 1000 kWh

Wind turbine

Micro gas turbine Waste-heat utilization Absorption chiller

Heat pump Heat storage

9

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

Table 2b Economic parameters for the case study [38,39]. Items

Capital cost $/kW

Maintenance cost $/kWh

Items

Parameter

Value

Micro gas turbine Gas boiler Absorption chiller Electrical chiller Heat storage Heat pump PV Wind turbine

550 50 230 150 25 200 1,650 1,150

0.005 0.0003 0.002 0.002 0.0003 0.002 0.002 0.003

System

Lifetime

15 years

Interest rate

6%

Natural gas emission factor

0.18 kg/kWh

Grid electricity emission factor

1.04 kg/kWh

Fig. 6. (a) Energy demands, (b) SRI and wind speed, and (c) time-of-use energy prices for the case study.

Fig. 7. Pareto frontiers and FOSs for the DP and SP models. Table 3 Fitted curves of the Pareto frontier for the DP and SP models. Category

Fitting equation

Range of ATC

R2

DP

ACE = 0.16269ATC2 − 234.52ATC + 86220

0.990

SP

ACE = 0.14778ATC2 − 220.44ATC + 83949

651–729 (103 $/year) 680–737 (103 $/year)

d+ (−) =

d=

(ATC − ATC ideal (nadir))2 + (ACE − ACE ideal (nadir))2

d+ d+ + d−

Fig. 8. Decision-making details of the two methods for the SP model.

5.3. Optimal design and capacity configuration

0.994

The combination of multi-objective optimization and decision making provides several optimal alternatives with different perspectives, for which the economic and environmental performance is evaluated. Table 4 presents detailed information for the points marked in Fig. 7. For the DP model, the installed capacity of the boiler and the electrical chiller are underestimated, owing to the limited number of typical days, which fails to cover all the characteristics for each uncertainty parameter. Even if the method of selecting typical days is improved or the day with the maximum load is added manually, the proportion of each typical day cannot be allocated accurately, leading

(25a)

(25b)

where d represents the deviation index, which indicates the deviation of each solution from the ideal solution. 10

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

Table 4 Optimal design and performance of each point. Category

Point

Optimizing principle

Installed capacity 2

DP

SP

A1 A2 A3 A4 B1 B2 B3 B4

Min ATC FOS (SE) FOS (ED) Min ACE Min ATC FOS (SE) FOS (ED) Min ACE

2

System performance 3

PV (m )

WT (m )

HP (kW)

MGT (kW)

Boiler (kW)

AC (kW)

EC (kW)

ST (kW)

ATC (10 $/ year)

ACE (tons/ year)

1200 1200 1200 1200 1119 627 672 1200

0 0 0 0 81 573 528 0

108 0 0 0 123 0 0 0

116 359 353 505 128 309 352 505

831 665 649 649 1103 928 907 850

128 315 319 585 141 337 381 753

655 467 564 601 1627 1430 1387 1215

213 700 713 212 315 516 583 1000

651 681 679 727 680 695 701 736

2575 1956 1979 1674 2418 2100 2037 1740

Deviation index (d)

– – – – 0.012 0.274 0.370 0.987

configuration. Additionally, as the location for the case study is in the “Cold” climate zone, where the temperature can hardly satisfy the requirement of the air-source heat pump in winter, the capacity of the heat pump is small at points A1 and B1. The installed capacities of the absorption chiller and the MGT tend to increase with the improvement of the ACE owing to the environmental friendliness and high efficiency. The deviation index (d) for the Shannon-entropy-based decision making is lower than that for the Euclidean-distance-based decision making. Therefore, the FOS from the Shannon-entropy-based decision making is selected as the most preferred optimum solution and is analyzed in the following sections. Fig. 9 presents the distributions of each design variable

to an oversimplified model. Additionally, the DP model tends to neglect the complementarity among renewable energy technologies; thus, the optimal capacity configuration frequently includes only one of them. As shown in Table 4, overall optimal solutions take full advantage of the roof space to install the PV and the wind turbine, indicating the significance of renewable-energy technology in the BIES system. Points B1 and B4 (associated with the single-objective optimization) preferred installation of PV panels rather than a wind turbine, because solar energy is more abundant and stable than wind power in Beijing. Moreover, points B2 and B3 correspond to the FOSs, achieving a tradeoff between the two conflicting objectives, where the capacities of the alternative technologies are evenly allocated in the optimal

Fig. 9. Scattered distributions of (a) the installed area of the PV panels and wind turbine, (b) the installed capacities of the heat pump and MGT, (c) the installed capacities of the gas boiler and heat storage, and (d) the installed capacities of the absorption and electrical chillers in the Pareto frontier. 11

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

total number of simulations [42].

corresponding to the Pareto frontier. The upper and lower bounds are marked by red and blue dashed lines, respectively. As shown, for the SP model, the optimal values of the PV panel-installed area are located in the upper half of the value range and are generally greater than those of the wind turbine, while the optimal values of these two design variables have a dispersive distribution in their ranges. Moreover, most of the optimal heat-pump capacities suggest the inapplicability of the heat pump in this case. Fig. 9(c) indicates that the optimal heat-storage capacities for the SP model are generally higher than those for the DP model, and some of them even reach the upper bound, indicating that the heat-storage technology can help to address the uncertainty of the BIES.

0 Ii = ⎧ ⎨ ⎩1

∃ Xfs failture

(26)

In this study, a statistical analysis of points B1, B2, and B4 is conducted to evaluate the robustness performance of the optimal capacity configuration. 5000 Monte Carlo simulations and dispatch optimizations for the SP model are performed via sampling according to each uncertainty parameter’s distribution. For B1 and B4, only one optimization in each simulation is conducted, according to the economic and environmental objectives, respectively. For the optimal solution B2, one objective function can be optimized by converting the other objective function into a series of constraints; thus, two optimizations with two different objectives are needed in each simulation. Fig. 12 presents the statistical analysis results for B1, B2, and B4 with regard to the ATC and ACE. The violin plots indicate the robustness performance of the design variables of these three optimal solutions under the uncertainty of the input parameters. The actual operation stage is usually unpredictable and cannot completely follow the dispatch schedule optimized in the planning stage, necessitating reoptimization of the operation variables in each simulation. In Fig. 12, the violin plots show the probability distributions of each index. The box indicates the boundaries of the 25th percentile and 75th percentile, and the white dot represents the mean value. The whisker indicates the range of the 5th percentile to the 95th percentile, and some outliers are marked on the line. Additionally, the results are presented in Table 5. In Table 5, the coefficient of variation (CV) is defined as the ratio of the standard deviation (STD) to the mean value, characterizing the dispersion degree of the probability distribution. A comparison of the performance indices for the three optimal solutions indicated that the mean values deviated from the original optimal values. The degrees of deviation are small and acceptable, particularly for the ATC indices of Min(ATC) and FOS(SE). This is because the fixed design variables determined the system’s capital cost, whereas carbon emissions are only generated in the operation stage. Consequently, the mean value of the economic index is closer to the original optimal value than that of the environmental index. Furthermore, the probability of failure is introduced to indicate whether the optimal capacity configuration can satisfy the constraints generated by random sampling. In Table 5, the probability of failure of Min(ACE) is 0, indicating that all the simulations have a feasible solution, and no constraints are violated. The probabilities of failure for FOS (SE) and Min(ATC) are 8.81% and 40.66%, respectively, indicating that there is no feasible solution in some simulations. Owing to the various capacity configurations among the three optimal solutions, the probabilities of failure different and depend significantly on the heat-storage

5.4. Energy consumption shares breakdown According to the energy interaction with the external environment, the energy consumption shares of the proposed system involve five aspects: the natural gas shares of the MGT and the gas boiler, the grid shares of the electricity feed-in and purchase, and the renewable share. Fig. 10 presents the energy shares of selected solutions for the DP and SP models. The electricity feed-in shares at all the points are close to 0, as they are impacted by both the low feed-in tariff and inadequate renewable energy resources, whereas the energy shares for electricity purchase generally locate in the range of 0%–25%. The results indicate that the bulk grid is of vital importance to the on-grid IES and is also conducive to cost reduction. Additionally, the maximum energy share of all the technologies is observed in the DP model, while the energy shares of the SP model tend to exhibit a more balanced allocation, as the scenarios provide more constraints to avoid overdependence on a certain technology. Additional insight regarding the operation details for the SP model can be obtained from the ternary plots of each point, as shown in Fig. 11, which presents the distributions of the energy carrier consumption shares for each stochastic scenario (including natural gas, grid, and renewable energy). The Min(ACE) with the best environmental performance tends to push the operating region toward lower natural gas shares and higher renewable shares for many cases in Europe [14]. However, the grid-emission factor for North China is significantly higher than that for Europe, because coal-fired generation accounts for more than half of the power supply of the North China grid. Owing to the strict limit for carbon emissions, the variations of the energy consumption shares are concentrated in the bottom-right region to avoid the purchase of electricity from the grid. In comparison, coalfired generation may contribute to a lower purchase price of electricity; thus, the region of energy carrier shares for the Min(ATC) solution becomes significantly wider and moves toward higher grid shares, as shown in Fig. 11(a). As shown in Fig. 11(b), FOS(SE) have a moderate distribution region compared with the two single-objective optimization solutions. According to the analysis of the energy shares, natural gas is the main input energy carrier, particularly for a low shape coefficient building with high energy consumption. 5.5. Robustness analysis To verify the robustness of the proposed model, a Monte Carlo based robustness analysis combined with optimization and simulation is performed using MATLAB and GAMS. The analysis procedure is described in Section 2.4, and the distribution of each uncertainty parameter is presented in Section 2.1. Owing to the random sampling and fixed capacity, some of the time-series uncertainty parameters may cause a failure to satisfy the operation constraints, resulting in the absence of feasible solutions and the inability to identify the unsatisfied constraints. Hence, a binary indicator is introduced to express the success or failure status for the whole system. Its value is 1 in the case of a failure and 0 otherwise. The probability of failure is defined as the number of simulations whose binary-indicator value is 1 divided by the

Fig. 10. Radar figure of energy consumption shares for each optimal solution. 12

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

Fig. 11. Distribution of the SP model’s energy carrier shares for each stochastic scenario: (a) Min(ATC), (b) FOS(SE), and (c) Min(ACE).

Fig. 12. Robustness performance of each optimal point with regard to the ATC and ACE.

economic and environmental objectives, indicating that it is the best choice for system designers and stakeholders.

and MGT capacities, indicating the importance of the energy-storage technology for the BIES system. Increasing the capacity of the energystorage equipment may reduce the probability of failure. Therefore, the optimal solution of Min(ATC) has a relatively high risk of operationconstraint violation, so that the expectation of planning will not occur in one or more years of project life. In contrast, the optimal solution of Min(ACE) completely passes the robustness validation, with zero probability of failure, at the expense of a high ATC. FOS(SE) exhibits a lower and acceptable probability of failure, with a tradeoff between the

6. Conclusions An MINLP-based two-stage multi-objective SP framework combining optimization and robustness analysis is developed for optimal design and dispatch of a BIES system under uncertainty. The decision variables are divided into design variables and operation variables,

Table 5 Mean, STD, and CV of ATC and ACE for the optimal solutions. Item Min ATC FOS (SE) Min ACE

Performance index ATC ACE ATC ACE ATC ACE

Mean

STD 3

3

683.81 × 10 $/year 2,537.3 tons/year 693.89 × 103$/year 2,060.5 tons/year 760.11 × 103$/year 1.696.9 tons/year

4.01 × 10 $/year 15.86 tons/year 4.65 × 103$/year 61.16 tons/year 3.66 × 103$/year 12.63 tons/year

13

CV

Probability of failure

0.59% 0.65% 0.67% 2.97% 0.48% 0.74%

40.66% 8.81% 0.00%

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

Declaration of Competing Interest

considering the multiple uncertainties caused by the SRI, wind speed, electrical demand, heating demand, and cooling demand. To capture the characteristics of these uncertainty parameters, a series of stochastic scenarios with corresponding probabilities are generated as a binarytree structure using k-medoids clustering method. The proposed model allows multi-objective optimization, and the Pareto frontier is obtained for the economic and environmental performance via the AUGMECON2 approach. Furthermore, decision-making methods, i.e., Shannon-entropy-based and Euclidean-distance-based methods, are employed to identify the FOS from the Pareto frontier. Finally, a simulation-optimization-based Monte Carlo method is proposed and implemented to verify the robustness performance for the economically optimal solution, the environmentally optimal solution, and the FOS. The following conclusions are drawn.

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Acknowledgements This work was supported by the National Key R&D Program of China [grant numbers 2018YFC0704602 and 2016YFC0700305], as well as the Global Environment Facility [grant number 1-A-CS-014]. References [1] Zhang Y, Campana PE, Lundblad A, Zheng W, Yan J. Planning and operation of an integrated energy system in a Swedish building. Energy Convers Manag 2019;199:111920https://doi.org/10.1016/j.enconman.2019.111920. [2] Razmi A, Soltani M, Torabi M. Investigation of an efficient and environmentallyfriendly CCHP system based on CAES, ORC and compression-absorption refrigeration cycle: Energy and exergy analysis. Energy Convers Manag 2019;195:1199–211. https://doi.org/10.1016/j.enconman.2019.05.065. [3] Pagliarini G, Corradi C, Rainieri S. Hospital CHCP system optimization assisted by TRNSYS building energy simulation tool. Appl Therm Eng 2012;44:150–8. https:// doi.org/10.1016/j.applthermaleng.2012.04.001. [4] Jing R, Wang M, Zhang Z, Liu J, Liang H, Meng C, et al. Comparative study of posteriori decision-making methods when designing building integrated energy systems with multi-objectives. Energy Build 2019;194:123–39. https://doi.org/10. 1016/j.enbuild.2019.04.023. [5] Soheyli S, Shafiei Mayam MH, Mehrjoo M. Modeling a novel CCHP system including solar and wind renewable energy resources and sizing by a CC-MOPSO algorithm. Appl Energy 2016;184:375–95. https://doi.org/10.1016/j.apenergy.2016.09.110. [6] Wei D, Chen A, Sun B, Zhang C. Multi-objective optimal operation and energy coupling analysis of combined cooling and heating system. Energy 2016;98:296–307. https://doi.org/10.1016/j.energy.2016.01.027. [7] Di Somma M, Yan B, Bianco N, Luh PB, Graditi G, Mongibello L, et al. Multi-objective operation optimization of a Distributed Energy System for a large-scale utility customer. Appl Therm Eng 2016;101:752–61. https://doi.org/10.1016/j. applthermaleng.2016.02.027. [8] Ju L, Tan Z, Li H, Tan Q, Yu X, Song X. Multi-objective operation optimization and evaluation model for CCHP and renewable energy based hybrid energy system driven by distributed energy resources in China. Energy 2016;111:322–40. https:// doi.org/10.1016/j.energy.2016.05.085. [9] Jing R, Wang M, Wang W, Brandon N, Li N, Chen J, et al. Economic and environmental multi-optimal design and dispatch of solid oxide fuel cell based CCHP system. Energy Convers Manag 2017;154:365–79. https://doi.org/10.1016/j. enconman.2017.11.035. [10] Zeng R, Li H, Jiang R, Liu L, Zhang G. A novel multi-objective optimization method for CCHP-GSHP coupling systems. Energy Build 2016;112:149–58. https://doi.org/ 10.1016/j.enbuild.2015.11.072. [11] Tran HN, Narikiyo T, Kawanishi M, Kikuchi S, Takaba S. Whole-day optimal operation of multiple combined heat and power systems by alternating direction method of multipliers and consensus theory. Energy Convers Manag 2018;174:475–88. https://doi.org/10.1016/j.enconman.2018.08.046. [12] Razmi AR, Janbaz M. Exergoeconomic assessment with reliability consideration of a green cogeneration system based on compressed air energy storage (CAES). Energy Convers Manag 2019;112320. https://doi.org/10.1016/j.enconman.2019.112320. [13] Mavromatidis G, Orehounig K, Carmeliet J. Uncertainty and global sensitivity analysis for the optimal design of distributed energy systems. Appl Energy 2018;214:219–38. https://doi.org/10.1016/j.apenergy.2018.01.062. [14] Mavromatidis G, Orehounig K, Carmeliet J. Design of distributed energy systems under uncertainty: A two-stage stochastic programming approach. Appl Energy 2018;222:932–50. https://doi.org/10.1016/j.apenergy.2018.04.019. [15] Lv J, Li YP, Shan BG, Jin SW, Suo C. Planning energy-water nexus system under multiple uncertainties – A case study of Hebei province. Appl Energy 2018;229:389–403. https://doi.org/10.1016/j.apenergy.2018.08.010. [16] Akbari K, Jolai F, Ghaderi SF. Optimal design of distributed energy system in a neighborhood under uncertainty. Energy 2016;116:567–82. https://doi.org/10. 1016/j.energy.2016.09.083. [17] Gabrielli P, Fürer F, Mavromatidis G, Mazzotti M. Robust and optimal design of multi-energy systems with seasonal storage through uncertainty analysis. Appl Energy 2019;238:1192–210. https://doi.org/10.1016/j.apenergy.2019.01.064. [18] Sharifzadeh M, Lubiano-Walochik H, Shah N. Integrated renewable electricity generation considering uncertainties: The UK roadmap to 50% power generation from wind and solar energies. Renew Sustain Energy Rev 2017;72:385–98. https:// doi.org/10.1016/j.rser.2017.01.069. [19] Gianniou P, Liu X, Heller A, Nielsen PS, Rode C. Clustering-based analysis for residential district heating data. Energy Convers Manag 2018;165:840–50. https:// doi.org/10.1016/j.enconman.2018.03.015. [20] Soulouknga MH, Doka SY, Revanna N, Djongyang N, Kofane TC. Analysis of wind speed data and wind energy potential in Faya-Largeau, Chad, using Weibull distribution. Renew Energy 2018;121:1–8. https://doi.org/10.1016/j.renene.2018.01.

(1) The proposed two-stage multi-objective SP model allows the design and operation optimization of BIES considering uncertainty, whereby the Pareto frontier with two competing objectives (ATC and ACE) can be obtained. The ATC values of the optimal solutions range from $680 × 103/year to $737 × 103/year, and the ACE values are in the range of 1740 to 2419 tons/year. By employing two decision-making methods and the deviation index, the FOS is identified as the desired solution, with an ATC of $695 × 103/year and an ACE of 2100 tons/year. (2) With regard to the optimal capacity configurations, the SP model tends to deploy both PV panels and wind turbines, whereas the DP model deploys only PV panels, failing to fully utilize the local renewable energy sources. Additionally, the DP model exhibits a lower capacity of the gas boiler and the electrical chiller in the optimal capacity configuration; thus, the ATC is typically underestimated. (3) According to the energy consumption shares breakdown for the operation stage, the environmentally optimal configuration tends to utilize natural gas rather than the power grid. In contrast, the economically optimal configuration tends to comprehensively utilize each energy carrier. (4) A Monte Carlo simulation-based robustness analysis verifies the robustness performance of the operation stage for the FOS, the economically optimal solution, and the environmentally optimal solution. The objective values of the simulation-optimization are concentrated around the original optimal objective value, and the CV is < 1%. According to the probability of failure, the economically optimal capacity configuration has a relatively high risk (40.66%) of operation-constraint violation, while the FOS (with a value of 8.81%) not only represents the optimal tradeoff between the economic and environmental objectives but also exhibits good robustness performance. (5) The proposed framework fully considers the main uncertainties and provides a series of optimal alternatives for stakeholders to deal with the optimality-robustness tradeoff. Additionally, according to the specific needs and the platforms’ computing performance in practical applications, designers can flexibly select the uncertain factors to be addressed and then modify the scenario structure for modeling. Furthermore, the robustness of optimal designs can be analyzed within the proposed framework.

CRediT authorship contribution statement Meng Wang: Conceptualization, Methodology, Software, Writing original draft. Hang Yu: Supervision, Project administration. Rui Jing: Methodology, Investigation, Writing - review & editing. He Liu: Resources, Investigation. Pengda Chen: Visualization, Software. Chaoen Li: Resources.

14

Energy Conversion and Management 208 (2020) 112589

M. Wang, et al.

[32] Wang J, Lu Y, Yang Y, Mao T. Thermodynamic performance analysis and optimization of a solar-assisted combined cooling, heating and power system. Energy 2016;115:49–59. https://doi.org/10.1016/j.energy.2016.08.102. [33] Bischi A, Taccari L, Martelli E, Amaldi E, Manzolini G, Silva P, et al. A detailed MILP optimization model for combined cooling, heat and power system operation planning. Energy 2014;74:12–26. https://doi.org/10.1016/j.energy.2014.02.042. [34] Chen Y, Zhang T, Yang H, Peng J. Study on energy and economic benefits of converting a combined heating and power system to a tri-generation system for sewage treatment plants in subtropical area. Appl Therm Eng 2016;94:24–39. https://doi. org/10.1016/j.applthermaleng.2015.10.078. [35] Jing R, Zhu X, Zhu Z, Wang W, Meng C, Shah N, et al. A multi-objective optimization and multi-criteria evaluation integrated framework for distributed energy system optimal planning. Energy Convers Manag 2018;166:445–62. https://doi. org/10.1016/j.enconman.2018.04.054. [36] You T, Wang B, Wu W, Shi W, Li X. Performance analysis of hybrid ground-coupled heat pump system with multi-functions. Energy Convers Manag 2015;92:47–59. https://doi.org/10.1016/j.enconman.2014.12.036. [37] Li M, Mu H, Li N, Ma B. Optimal design and operation strategy for integrated evaluation of CCHP (combined cooling heating and power) system. Energy 2016;99:202–20. https://doi.org/10.1016/j.energy.2016.01.060. [38] Wu Q, Ren H, Gao W, Ren J. Multi-criteria assessment of building combined heat and power systems located in different climate zones: Japan-China comparison. Energy 2016;103:502–12. https://doi.org/10.1016/j.energy.2016.02.156. [39] Wang X, Yang C, Huang M, Ma X. Multi-objective optimization of a gas turbinebased CCHP combined with solar and compressed air energy storage system. Energy Convers Manag 2018;164:93–101. https://doi.org/10.1016/j.enconman.2018.02. 081. [40] Wang M, Jing R, Zhang H, Meng C, Li N, Zhao Y. An innovative Organic Rankine Cycle (ORC) based Ocean Thermal Energy Conversion (OTEC) system with performance simulation and multi-objective optimization. Appl Therm Eng 2018;145:743–54. https://doi.org/10.1016/j.applthermaleng.2018.09.075. [41] Feng Y, Hung TC, Zhang Y, Li B, Yang J, Shi Y. Performance comparison of lowgrade ORCs (organic Rankine cycles) using R245fa, pentane and their mixtures based on the thermoeconomic multi-objective optimization and decision makings. Energy 2015;93:2018–29. https://doi.org/10.1016/j.energy.2015.10.065. [42] Gholaminezhad I, Jafarpur K, Paydar MH, Karimi G. Multi-scale multi-objective optimization and uncertainty analysis of methane-fed solid oxide fuel cells using Monte Carlo simulations. Energy Convers Manag 2017;153:175–87. https://doi. org/10.1016/j.enconman.2017.10.011.

002. [21] Kang J, Wang S. Robust optimal design of distributed energy systems based on lifecycle performance analysis using a probabilistic approach considering uncertainties of design inputs and equipment degradations. Appl Energy 2018;231:615–27. https://doi.org/10.1016/j.apenergy.2018.09.144. [22] Gang W, Wang S, Shan K, Gao D. Impacts of cooling load calculation uncertainties on the design optimization of building cooling systems. Energy Build 2015;94:1–9. https://doi.org/10.1016/j.enbuild.2015.02.032. [23] Zhou Z, Zhang J, Liu P, Li Z, Georgiadis MC, Pistikopoulos EN. A two-stage stochastic programming model for the optimal design of distributed energy systems. Appl Energy 2013;103:135–44. https://doi.org/10.1016/j.apenergy.2012.09.019. [24] Mavrotas G, Florios K. An improved version of the augmented s-constraint method (AUGMECON2) for finding the exact pareto set in multi-objective integer programming problems. Appl Math Comput 2013;219:9652–69. https://doi.org/10. 1016/j.amc.2013.03.002. [25] Feng Y, Hung TC, Greg K, Zhang Y, Li B, Yang J. Thermoeconomic comparison between pure and mixture working fluids of organic Rankine cycles (ORCs) for low temperature waste heat recovery. Energy Convers Manag 2015;106:859–72. https://doi.org/10.1016/j.enconman.2015.09.042. [26] Luo Z, Yang S, Xie N, Xie W, Liu J, Souley Agbodjan Y, et al. Multi-objective capacity optimization of a distributed energy system considering economy, environment and energy. Energy Convers Manag 2019;200:112081https://doi.org/10. 1016/j.enconman.2019.112081. [27] Lotfi, Hosseinzadeh F, Fallahnejad, Reza. Entropy, Vol. 12, Pages 53-62: Imprecise Shannon’s Entropy and Multi Attribute Decision Making 2010. [28] Jing R, Wang M, Zhang Z, Wang X, Li N, Shah N, et al. Distributed or centralized? Designing district-level urban energy systems by a hierarchical approach considering demand uncertainties. Appl Energy 2019;252:113424https://doi.org/10. 1016/j.apenergy.2019.113424. [29] Zhu X, Zhan X, Liang H, Zheng X, Qiu Y, Lin J, et al. The optimal design and operation strategy of renewable energy-CCHP coupled system applied in five building objects. Renew Energy 2020;146:2700–15. https://doi.org/10.1016/j. renene.2019.07.011. [30] Zeng R, Zhang X, Deng Y, Li H, Zhang G. An off-design model to optimize CCHPGSHP system considering carbon tax. Energy Convers Manag 2019;189:105–17. https://doi.org/10.1016/j.enconman.2019.03.062. [31] Jing R, Wang M, Liang H, Wang X, Li N, Shah N, et al. Multi-objective optimization of a neighborhood-level urban energy network: considering Game-theory inspired multi-benefit allocation constraints. Appl Energy 2018;231:534–48. https://doi. org/10.1016/j.apenergy.2018.09.151.

15