An MILP Method for Design of Distributed Energy Resource System Considering Stochastic Energy Supply and Demand

A distributed energy resource (DER) system, which can be defined as a medium or small energy conversion and utilization system with various functions for meeting multiple targets, is directly oriented towards users and achieves on-site production and energy supply according to users’ demands. Optimization research on system construction has recently become an important issue. In this paper, simple stochastic mathematical equations were used to interpret the optimal design problem of a DER system, and based on this, a novel method for solving the optimization problem, which has multi-dimensional stochastic uncertainties (involving the price of input-energy and energy supply and demand), was put forward. A mixed-integer linear programming (MILP) model was established for the optimal design of the DER system by combining the ideas of mean value and variance, aiming to minimize the total costs, including facility costs, energy purchase costs, and loss caused by energy supply shortage, and considering the energy balance and facility performance constraints. In the end, a DER system design for an office building district in Xuzhou, China, was taken as an example to verify the model. The influences of uncertainty on the selection of system facilities and the economic evaluation were analyzed. The result indicated that uncertainty of energy demand played a significant role in optimal design, whereas energy price played a negligible role. With respect to economy, if uncertainties are not considered in system design, it will result in a short supply, and therefore the total cost will increase considerably. The calculation convergence was compared with previous work. The implementation results showed the practicality and efficiency of the proposed method.


Background
A distributed energy resource (DER) system, which is directly user-oriented and can produce energy in local areas, is a small or medium energy conversion and utilization system characterized by high versatility and the ability to meet multiple objectives [1].The advantage of a DER system is embodied in combined cooling, heating, and power (CCHP) [2].CCHP conforms to "cascade utilization" of the total energy system, and has good energy efficiency.Although large (heat) power plants can transmit energy over long distances, it is necessary to build power grids and substations, which leads to transmission losses.Furthermore, long-distance delivery is difficult for heating and cooling energy.In general, it is difficult to find enough appropriate users who need heating and cooling energy near plant locations by virtue of the requirements for the site selection of large plants, leading to ineffective CCHP.However, a DER system can overcome this issue by being set near the demand and cooperating with users as much as possible.In this way, long-distance heating and cooling delivery, as well as power-loss of grids, are avoided [3].A DER system allows the use of units that possess more capability for adjustment, control, and assurance, and promises sufficient supply of various secondary energies for the using units, which are applicable for power, heating and cooling energy supplies for cities [4], factories [5], commercial districts [6], large-scale public facilities [7,8], and transportation systems [9], and reduces environmental pressure.All in all, a DER system can provide the possibility for comprehensive cascade utilization of energy, opens up new directions for renewable energy, and makes a significant contribution to improving energy efficiency, security, and solving environmental pollution [10].
A DER system, however, is related to many kinds of resources and facilities, which have complicated relationships [11].Considering the variation of regional resource prices [12], demands [13], the environment [14], etc., the optimal design of a DER system is a complex task [15].Hence, the study of the optimal design of a DER system plays an important role in cutting down on facility costs, energy purchase costs, and losses caused by energy supply shortage, as well as enhancing energy efficiency.

Literature Review
Recently, many scholars have researched different optimization issues related to DER systems, such as heating networks [16], power networks [17,18], heating and cooling networks [19], heating and power networks [20,21], and complex networks [22].Previously, some scholars solved this kind of issue based on the linear programming (LP) [23,24] method.Later, many were inclined to use mixed integer linear programming (MILP) [25][26][27], mixed-integer nonlinear programming (MINLP) [28,29] and multi-objective programming (MOP) [30][31][32][33] methods for this issue.Pedrasa et al. [34] scheduled a smart-home energy services DER system using a decision-support tool that was able to optimize energy service provisions, and then schedule the available DER to maximize the net benefits.The corresponding optimization problem was solved using particle swarm optimization (PSO).Zhu and Tomsovic [35] proposed a three-stage optimal algorithm for optimal distribution power flow of a DER system.In that paper, a region of Japan was set as an example to verify the practicality of the algorithm.Calvillo et al. [36] proposed a hybrid energy system consisting of wind, photovoltaic (PV), and fuel cells, designed to supply continuous power to a load.A simple and economical control was used for maximum power point tracking and, hence, maximum power extraction from a wind turbine and PV array.Mehleri et al. [37] presented a new model for the optimal configuration of DER systems and the design of a heating pipeline network.The model was applied to a neighborhood of 10 buildings located in Central Greece, which demonstrated the practicality of the proposed model.Subsequently, Mehleri et al. [38] proposed an MILP model that considered the effect of PVs on system economy and the environment, and applied the model to a neighborhood of five buildings located in Athens (Greece).Yang et al. [39] put forward an MILP model to design centralized generation district-scale DER systems, which was applied to an urban area in Lianyungang City, China.Results demonstrated that, based on reduced total annual costs of energy, equipment and operation, in addition to the revenue from energy subsidies of 14.1%, the additional initial investment for the proposed DER systems could be paid back within three years.This research established a DER system from a different perspective, which was successfully put into practice.Buoro et al. [40] took both economic and environmental objectives into account, and presented a multi-objective model for the optimal design and operation of a district heating network.A case of a region in Italy was given to demonstrate the practicality of the proposed model.Ruan et al. [41] analyzed the energy consumption characteristics of different commercial buildings, and examined the variability of BCHP (building combined heat and power) in different commercial buildings.However, few have discussed the uncertainty of model parameters.
Load analysis is essential for DER system design.If the load analysis or evaluation for buildings is not detailed enough, it will exert a great deal of influence on operations [42,43].However, energy supply and demand are usually characterized by strong uncertainty [44].If the uncertainty is not taken into account, there will be a mismatch between supply and demand for a DER system and, thereby, the investment and operation and maintenance (O&M) costs will be increased.For example, in winter or summer, there may be insufficient energy supply.In transition seasons, there may be inadequate utilization of a DER system.In light of the above considerations, many scholars have conducted further research on DER system optimization with regard uncertainty.With regard to stochastic uncertain optimization, the Monte Carlo simulation is one of the most common methods to deal with models [45][46][47].The principle is that the uncertain parameters of the models are stochastically generated with several groups of values according to probability distributions.The groups of values are combined to divide the problem into a number of scenarios.In this way, an uncertain optimization problem is transformed into a certain optimization problem for multiple scenarios [48].
Marquez et al. [49] proposed a general approach to complex system reliability assessment, based on Monte Carlo simulations.Alarcon-Rodriguez et al. [50] presented a novel multi-objective planning framework for the integration of stochastic and controllable DER in a distribution grid.Calvillo et al. [51] put forward an optimization model based on stochastic price for the optimal planning and operation of an aggregated DER system.Mavrotas et al. [52] established a novel method for energy-planning optimization, considering uncertain economic parameters and coupled with MILP and Monte Carlo simulations.Similarly, Li et al. [53] proposed an optimization method coupling MINLP and Monte Carlo simulations, applied to building cooling, heating, and power (BCHP) systems.Yang et al. [15] and Zhou et al. [54] proposed an effective two-stage stochastic programming model to work out a DER system on the basis of Monte Carlo simulations, which was respectively verified through a real case.Ahmed et al. [55] developed a branch-and-bound algorithm to address two-stage stochastic integer programs.This method was well applicable since the explicit enumeration of all discontinuities of value functions could be avoided and the uncertainties in the cost parameters and constraint matrix could be allowed.Somma et al. [56] introduced a multi-objective linear programming model that could satisfy time-varying user demands and reduce both energy costs and environmental impacts.A hypothetical large hotel (16,000 m 2 ) located in Italy was set as an example to verify the practicality.Afterwards, Somma et al. [57] developed a multi-objective model, in which exergy was investigated in a DES (Distributed Energy System) design to attain a rational use of energy resources.A case study on a hypothetical cluster of 30 buildings, located in Turin (Italy), was presented to demonstrate that the maximum reduction in total annual costs and primary exergy inputs were equal to 33.5% and 35.8%, respectively, when using the proposed model.Additionally, the sensitivities of energy prices and demand scales were also examined to investigate their influence on the optimization model.Somma et al. [58] considered energy demand and supply uncertainty, established a stochastic programming model for a DER system with multiple energy devices, and finally provided good solutions for DER system operators, based on economic and environmental priorities.
The optimization model of a DER system coupled with Monte Carlo simulations has matured and obtain good results.However, the disadvantage is that the number of divided scenarios cannot be decided before solving models.If there are fewer scenarios, the solving time will be longer, and the strong uncertainty of model parameters requires many more scenarios to ensure the rationality of the model.If the scale of the study problem is too great, leading to strong uncertainty in model parameters, it will be difficult to obtain an optimal result in an ideal solving time.
This paper puts forward a novel approach to the optimization issue with multi-dimensional stochastic uncertainty.The probability distribution function of uncertain parameters was incorporated in the MILP model using appropriate logical processes.The uncertainty of the model was expressed by analytic equations, so as to transform the optimization model with multi-dimensional stochastic uncertainty into certainty.The model solution was carried out by the MILP algorithm; therefore, the reduction in computational efficiency, which is caused by Monte Carlo simulations introducing more scenarios into the model could be effectively avoided.Additionally, the effect of the uncertainties in energy demand and energy price on the optimal system design was investigated in detail.

•
A novel approach is put forward to solve optimization issues with multi-dimensional stochastic uncertainty.

•
Mean values and variances are introduced to build up a MILP model for optimal design of a DER system.

•
The practicality and solving efficiencies of the proposed model are verified through a real-world case study.

Paper Organization
This paper is divided into six sections.The adopted methodology and details of the mathematical model are given in Sections 2 and 3, respectively.Section 4 provides a detailed description of the example, Xuzhou, China, and provides four cases.Section 5 presents the implementation results.Conclusions are provided in Section 6.

DER System Configuration
This paper focuses on the optimal design of a DER system.The superstructure of the system is shown in Figure 1.There are two kinds of facilities in the system, used for energy generation and storage.The energy generation facilities consist of energy supply and conversion facilities, and the storage facilities usually include heating and cooling storage.Energy, such as purchased power, wind energy, solar energy, and natural gas, is input into the system as the basis for the normal operation of energy generation facilities, which is called input-energy for convenience.The system finally generates demanded energy, such as power, heating and cooling energy, which is called output-energy, while the heating and cooling energy need to be stored in storage facilities before being provided to the demand-side.In a DER system, various types of energy may be generated during the operation of energy-generation facilities.For example, some energy (such as exhaust gas energy) can be used by other facilities (such as absorption chiller/heater) to generate output-energy.This kind of energy is named "mid-energy" in the paper.

Problem Description
Input-energy is mainly composed of purchased power, wind energy, solar energy, natural gas, etc.In districts with a complete infrastructure, the supply of purchased power and natural gas in a DER system is relatively stable, which indicates a smaller variance in supply distribution function.However, for wind energy and solar energy, the supply is associated with climate, thereby leading to strong uncertainty, indicating a greater variance in the supply distribution function.It is the uncertainty of input-energy that causes output-energy in a DER system to be full of strong uncertainty.Meanwhile, climate, user numbers, units using quantity, and the randomness of accidents continue to strengthen the uncertainty of demand for output-energy, contributing to the complex issue of optimal design of a DER system.On the basis of the probabilistic method, this paper was able to obtain the distribution function of available energy by subtracting energy demand from energy supply, and then solving the probability density of energy sources, which is less than zero; thus, the probability of the energy supply shortage can be worked out.In this paper, the output-energy demand and supply are assumed to obey the normal distribution function, which represents the mean value and variance.The distribution function of the available energy sources can be obtained through subtraction between the two normal functions, i.e., as shown in Figure 2. To solve the probability density of this function, which is less than zero, the zero point of this function can be determined by the corresponding probability density of the standard normal function, which can be obtained by consulting the table of standard normal distribution.
Input-energy mainly consists of power, wind energy, natural gas, etc.; among which power and natural gas are usually provided by industry or government with a more stable supply, and therefore with lower uncertainty.For wind energy and solar energy, supply is mainly decided by local climate, such as in summer having stronger solar radiation supporting relatively abundant solar energy, while in winter there is weaker solar radiation, supporting relatively little solar energy.Therefore, significant uncertainty can result in uncertain output-energy and storage facilities.Meanwhile, changes in the price of input-energy, user numbers, unit usage quantity, and random factors such as accidents, further make the output-energy uncertainty stronger, contributing to the complex issue of optimal design for a DER system.
This study suggests a new way to solve the multi-dimensional uncertainty issue, i.e., by introducing variance to show various uncertainties in the system based on a stochastic method.
Assuming input-energy fulfills a normal distribution ( ) σ shows the variance of output-energy; then, the available system energy will follow a normal distribution ( ) , as shown in Figure 2. If the available system energy's average value is lower than 0, the system will have a short supply.Based on the probability density function of available system energy, the probability of system in short supply can be worked out.In this paper, the output-energy demand and supply are assumed to obey the normal distribution function, which represents the mean value and variance.The distribution function of the available energy sources can be obtained through subtraction between the two normal functions, i.e., as shown in Figure 2. To solve the probability density of this function, which is less than zero, the zero point of this function can be determined by the corresponding probability density of the standard normal function, which can be obtained by consulting the table of standard normal distribution.
Input-energy mainly consists of power, wind energy, natural gas, etc.; among which power and natural gas are usually provided by industry or government with a more stable supply, and therefore with lower uncertainty.For wind energy and solar energy, supply is mainly decided by local climate, such as in summer having stronger solar radiation supporting relatively abundant solar energy, while in winter there is weaker solar radiation, supporting relatively little solar energy.Therefore, significant uncertainty can result in uncertain output-energy and storage facilities.Meanwhile, changes in the price of input-energy, user numbers, unit usage quantity, and random factors such as accidents, further make the output-energy uncertainty stronger, contributing to the complex issue of optimal design for a DER system.
This study suggests a new way to solve the multi-dimensional uncertainty issue, i.e., by introducing variance to show various uncertainties in the system based on a stochastic method.Assuming input-energy fulfills a normal distribution (µ s1 , σ s1 ),in which µ s1 shows the average value of input-energy, and σ s1 shows the variance of input-energy; output-energy fulfills a normal distribution (µ s2 , σ s2 ), in which µ s2 shows the average value of output-energy, σ s2 shows the variance of output-energy; then, the available system energy will follow a normal distribution (µ s1 − µ s2 , σ s1 + σ s2 ), as shown in Figure 2. If the available system energy's average value is lower than 0, the system will have a short supply.Based on the probability density function of available system energy, the probability of system in short supply can be worked out.Based on such a theory, to reduce the probability of a short supply of energy, the average values of input-energy and output-energy (higher and as closer as possible to zero) are needed.Therefore, the problem in the optimal design of a DER system is to find the most economic energy distribution plan for users.

MILP Model
Building a MILP model is a typical way of solving the issue mentioned in 2.2.For such a question, the lowest total cost (i.e., facility cost, energy purchase cost, and loss caused by energy supply shortage) is the objective function.Meanwhile, constraint conditions, such as facilities, energy balance, etc., are also taken into consideration.Based on the in situ data and facility alternatives, the optimization of energy distributed system can be converted into a MILP model.In accordance with the branch-and-bound approach, the final result can be determined.
To describe the method proposed in this paper and facilitate understanding for in situ engineers, the model can be performed using the following compact form: , , 0 where d represents the decision variable for facilities; o represents the decision variable for system operation; p represents the model parameters; F represents the objective function; 1 f represents the facility purchase cost function, which is relative to d , the variable for system operation; 2 f represents the energy purchase cost, which is relative to the variable for system operation; 1 ϕ , 2 ϕ represents the inequality in energy balance between supply and demand; and 1 φ , 2 φ represent the facility performance inequality and equality constraints.Based on such a theory, to reduce the probability of a short supply of energy, the average values of input-energy and output-energy (higher and as closer as possible to zero) are needed.Therefore, the problem in the optimal design of a DER system is to find the most economic energy distribution plan for users.

MILP Model
Building a MILP model is a typical way of solving the issue mentioned in 2.2.For such a question, the lowest total cost (i.e., facility cost, energy purchase cost, and loss caused by energy supply shortage) is the objective function.Meanwhile, constraint conditions, such as facilities, energy balance, etc., are also taken into consideration.Based on the in situ data and facility alternatives, the optimization of energy distributed system can be converted into a MILP model.In accordance with the branch-and-bound approach, the final result can be determined.
To describe the method proposed in this paper and facilitate understanding for in situ engineers, the model can be performed using the following compact form: where d represents the decision variable for facilities; o represents the decision variable for system operation; p represents the model parameters; F represents the objective function; f 1 represents the facility purchase cost function, which is relative to d, the variable for system operation; f 2 represents the energy purchase cost, which is relative to the variable for system operation; ϕ 1 , ϕ 2 represents the inequality in energy balance between supply and demand; and φ 1 , φ 2 represent the facility performance inequality and equality constraints.
Considering energy input, demand and price with a definite uncertainty, combining the MILP modeling method with multiple scenarios is common.The structural framework of the model is shown as Equation (2).Uncertainty among models could be solved perfectly using such a method.However, for such a problem, energy input, demand and price strongly mix together, and if the model's uncertainty is strong, or there are multiple uncertain variables, more scenarios will be required to show such uncertainty factors.Despite having more scenarios, the larger a model's scale is, the slower the model converges.min where p DUs , p PUs represent the indefinite parameter of model energy supply and demand and price in a specific scenario; p C represents the definite parameter of the model, such as facility cost; o s represents the system operation related decision variable in a specific scenario; f 1s , f 2s represents the facility purchase cost function in a specific scenario; N S is the total number of scenarios.
The new method put forth in this paper is shown below: We discovered that the uncertainty parameter (p DUs ) was only shown in the model constraints, and was not shown in an objective function.Therefore, we were able to use various kinds of transfers to change the uncertainty parameter into a certainty constraint in order to reduce the coupling degree of uncertainty.Based on the analysis in the previous section, we introduce the mean value and variance into the constraint to show the uncertainty of input-energy and demand, which can be shown by the probability of energy supply shortage.Therefore, the uncertainty parameter (p DUs ) can be saved to reduce the coupling degree of the uncertainty parameter in the model; in other words, to reduce the number of scenarios in order to solve the model, and to improve the convergence and calculation rates.The specific model structure is: where f 3s is a penalty term of the objective function, introduced when a system is in a scenario of short supply; µ s , σ s represent mean value and variance of energy supply and demand.

Model Requirements
The model was formulated as a MILP, and the optimization was executed using MATLAB 2014.A detailed design and operation scheme of a DER system can be determined by solving the model. Given: • Numbers of seasons and corresponding days.

•
Costs of purchased facility and energy, as well as O&M and depreciation rate.

•
Input-energy supply and output-energy demand distributions at every time window of every season. Determine: • Design scheme of DER system: Facility specifications.

•
Operation scheme of DER system: Operation scheme of energy generation and storage facilities.

•
The probability of energy supply shortage of a DER system.

•
Mean O&M cost of the system per day in every season.
Objective: Minimize the total mean facility cost per day, energy purchase cost, and loss caused by energy supply shortage under various operational and technical constraints.
In order to build and solve the model effectively, the following assumptions were made: • This paper assumes that stochastic uncertainties of supply and demand conform to a normal distribution.

•
The characteristic equations of input-and output-energy are linear.

•
The influence of ambient air temperature on facility operation is not taken into consideration.

•
The facility capacity is optimal, based on investigations, which can maximize savings in facility investment costs, and the parameter representing the facility capacity is a fixed value in the model.

Objective Function
The objective function is min where f 1 represents the mean facility cost per day; f 2 represents the mean input-energy cost per day; and f 3 represents the mean cost of the mismatch between supply and demand per day.
The mean facility cost per day includes purchase depreciation cost of storage and generation facilities and O&M costs, where the daily depreciation rate is equal to the ratio of the sum of 1 and the residual rate (unit: %) to the product of facility service life (a) and 365; in general, the residual rate is 5%.
The mean input-energy cost per day is related to the used quantity and unit price of fundamental energy: The mean cost of the mismatch between supply and demand per day is related to the mean probability of demand exceeding supply per day.Generally, an energy system should avoid the condition where demand exceeds supply as much as possible; thus, C LERRORs is given a larger value in the model, with an order of magnitude of 10 6 :

Energy Balance Constraints
At the start of a time window, the mean value of available output-energy equals the mean value of available output-energy at the start of the last time window, minus the attrition rate, plus the mean value of output-energy generated by the operating facilities, and minus the mean value of energy demand and the excess energy released: The variance of available output-energy at the start of a time window is equal to the variance of available output-energy at the start of the last time window, multiplied by 1, minus attrition rate difference, plus the square of the variance of the output-energy generated by the operational facilities, and minus the variance of demand: Binary variables can be used to judge the ranges that the variance of available output-energy belong to: (B ABse,t,s,a − 1)Ma + A Bmina ≤ δ SRse,t,s , se ∈ SE, t ∈ T, s ∈ S, a ∈ AB (9) The variance of output-energy can only belong to one single range: Linearize the relationship between variance and standard deviations using the phased approach, and determine the variance ranges in order to select the linearized equation parameters: After obtaining the variance and standard deviation of available output-energy at the start of a time window, the probability of demand exceeding supply (i.e., the probability density of available output-energy, which is less than zero) can be calculated.In this paper, the probability distribution of the available output-energy is transformed into a standard normal distribution, and then the application of binary variables and a look-up table were employed to solve the probability density of available output-energy, which was less than zero: The mean value of the available output-energy at the start of a time window should be greater than zero: The paper defines one day as a study cycle.The mean value of available output-energy at the final time window, plus the mean generation value of output-energy at the final time window, and minus the demanded mean value, should equal the mean value of the available output-energy at the start of the first time window: µ SRsa,se,tm,s + µ SAsa,se,tm,s − µ SDsa,se,tm,s = µ SRsa,se,1,s , sa ∈ SA, se ∈ SE, s ∈ S (19)

Facility Performance Constraints
The mean value of available output-energy at the start of a time window should be lower than the upper limit of the storage facilities: The mean value of output-energy generated by generation facilities can be worked out using their characteristic equations, which are related to input-energy and facility performance parameters.
The variance of output-energy generated by generation facilities remains linear with the input-energy: If one generation facility is not selected, the mean value of output-energy generated by this facility (i.e., the equation) should be zero: In the system, there should be no more than one type of facility: In the system, the total usage of input-energy should equal the sum of input-energy that all the generation facilities need: The input-energy variance of generation facilities should equal the variance of input-energy in the system: The total generation of output-energy in the system should equal the sum of the energy generated by all the generation facilities: The total generation variance of every output-energy should equal the variance sum of the energy generated by all the generation facilities: If there is mid-energy applied in the system, the mean value of the output should equal the mean value of usage minus the energy released: If there is mid-energy applied in the system, the variance of usage should equal the variance of the output:

Example
To verify the model's practicality, an example-namely, a DER system design based on an office building district in Xuzhou, China-is presented in this paper.The required parameters, including climate, energy price and demand of the district, as well as the facility specifications of a DER system, are displayed below.The local climate, including ambient air temperature, solar radiation intensity, wind energy density, and the energy demand, hardly vary within a two-hour time window.Fewer time windows and a smaller model scale will result in a faster solving efficiency.Therefore, we determined 12 time windows every day, with each time window being 2 h.

Climate Data
The ambient air temperature, solar radiation intensity, and wind energy density of the city fluctuate throughout different seasons.The ambient air temperature and solar radiation intensity are higher in summer, but lower in winter than in spring and autumn; and the wind energy density is similar for spring and autumn.Therefore, the climate can be divided into three periods: Mid-season, summer, and winter, with durations of 92, 153, and 120 days, respectively.The ambient air temperature, solar radiation intensity, and wind energy density of the city, for different periods, are shown in Figure 3. Here, we assume that solar radiation intensity meets the normal distribution, [59] with the variance taken as the mean value multiplied by 0.5, and that wind energy density meets the Weibull distribution [60], which is as follows: where τ is the scale parameter, and represents the value of each point on the curves in Figure 3c; ν is the shape parameter, and equals 3.69.

Energy Demand
The heating and cooling energy, and the power demand of the district at different periods were obtained by site investigation, as shown in Figure 4.More energy is required during the daytime, while less is required during the night.In the mid-season, summer, and during winter days, the demand for energy differs.The demand for cooling energy in summer is far higher than in winter, while the demand for heating energy is far higher in winter than in summer.The demand for power experiences less fluctuation in the three periods.In this paper, we assume that the energy demand satisfies a normal distribution, and the variance is taken as the mean value multiplied by 0.5.

Energy Demand
The heating and cooling energy, and the power demand of the district at different periods were obtained by site investigation, as shown in Figure 4.More energy is required during the daytime, while less is required during the night.In the mid-season, summer, and during winter days, the demand for energy differs.The demand for cooling energy in summer is far higher than in winter, while the demand for heating energy is far higher in winter than in summer.The demand for power experiences less fluctuation in the three periods.In this paper, we assume that the energy demand satisfies a normal distribution, and the variance is taken as the mean value multiplied by 0.5.

Energy Demand
The heating and cooling energy, and the power demand of the district at different periods were obtained by site investigation, as shown in Figure 4.More energy is required during the daytime, while less is required during the night.In the mid-season, summer, and during winter days, the demand for energy differs.The demand for cooling energy in summer is far higher than in winter, while the demand for heating energy is far higher in winter than in summer.The demand for power experiences less fluctuation in the three periods.In this paper, we assume that the energy demand satisfies a normal distribution, and the variance is taken as the mean value multiplied by 0.5.

Energy Price
The input-energy consists of power, natural gas, wind energy, and solar energy, where the power price is 0.92 CNY/kWh, the natural gas price is 0.32 CNY/kWh, and the wind energy and solar energy price are 0. To study the price uncertainty of input-energy, many scholars have applied a uniform distribution, normal distribution, and triangular distribution to describe the energy prices.In this paper, uniform distribution was selected to deal with the uncertainty of power price.The upper and lower limits of power price are 1.03 CNY/kWh and 0.81 CNY/kWh, respectively.The triangular distribution is determined for the uncertainty of the natural gas price, with the upper and lower limits of 0.37 CNY/kWh and 0.27 CNY/kWh, and a mode of 0.35 CNY/kWh.

DER System and Facility Parameter
This paper denotes power, natural gas, wind energy, and solar energy as the input-energy, exhaust gas as the mid-energy, and power, heating and cooling energy as the output-energy.The energy supply facilities generate the mid-or output-energy through the input-energy, mainly including internal combustion engines, gas turbines, gas boilers, fuel cells, wind turbines, solar PV, and solar heat collectors, shown in Table 1.The energy conversion facilities generate the outputenergy by inputting the mid-energy, which includes absorption chillers/heaters, heat exchangers, electric heaters, water-source chillers/HPs and air-source chillers/HPs, and the detailed parameters are shown in Table 2. Detailed parameters of energy storage facilities are shown in Table 3.The operation principle of the energy flow of each facility is shown in Figure 5.

Energy Price
The input-energy consists of power, natural gas, wind energy, and solar energy, where the power price is 0.92 CNY/kWh, the natural gas price is 0.32 CNY/kWh, and the wind energy and solar energy price are 0. To study the price uncertainty of input-energy, many scholars have applied a uniform distribution, normal distribution, and triangular distribution to describe the energy prices.In this paper, uniform distribution was selected to deal with the uncertainty of power price.The upper and lower limits of power price are 1.03 CNY/kWh and 0.81 CNY/kWh, respectively.The triangular distribution is determined for the uncertainty of the natural gas price, with the upper and lower limits of 0.37 CNY/kWh and 0.27 CNY/kWh, and a mode of 0.35 CNY/kWh.

DER System and Facility Parameter
This paper denotes power, natural gas, wind energy, and solar energy as the input-energy, exhaust gas as the mid-energy, and power, heating and cooling energy as the output-energy.The energy supply facilities generate the mid-or output-energy through the input-energy, mainly including internal combustion engines, gas turbines, gas boilers, fuel cells, wind turbines, solar PV, and solar heat collectors, shown in Table 1.The energy conversion facilities generate the output-energy by inputting the mid-energy, which includes absorption chillers/heaters, heat exchangers, electric heaters, water-source chillers/HPs and air-source chillers/HPs, and the detailed parameters are shown in Table 2. Detailed parameters of energy storage facilities are shown in Table 3.The operation principle of the energy flow of each facility is shown in Figure 5.

Setting of Cases
In the model, there may be uncertainties, such as the price of input-energy, supply of wind and solar energy, and energy demand.This paper will discuss four cases, as follows: 1.All of the prices of input-energy, as well as the energy supply and demand, are fixed values.2. The price of input-energy is perceived as an uncertain factor, and the energy supply and demand are fixed values.3. The price of input-energy is a fixed value and the energy supply and demand are uncertain factors.

Setting of Cases
In the model, there may be uncertainties, such as the price of input-energy, supply of wind and solar energy, and energy demand.This paper will discuss four cases, as follows: 1.
All of the prices of input-energy, as well as the energy supply and demand, are fixed values.

2.
The price of input-energy is perceived as an uncertain factor, and the energy supply and demand are fixed values.

3.
The price of input-energy is a fixed value and the energy supply and demand are uncertain factors.4.
All of the prices of input-energy and the energy supply and demand are uncertain factors.

Convergence Analysis
During model solving, it is possible that the calculation results of different scenarios differ from one another.A proper scenario can be obtained by analyzing the tendency of the objective value as the scenario number varies, in order to avoid failing to search for the optimal solution due to fewer scenario numbers and larger amounts of time being spent.Since the uncertainties of energy supply and demand are expressed by variances, the convergence effect should be better in light of theoretical analysis.To verify the advantage of this method, the previous method (the uncertainties of energy supply and demand are expressed by scenarios) was employed for comparison.Table 4 shows the problem size and computing times for Cases 2, 3, and 4, with 10, 20, and 30 selected scenarios, respectively.In each case, all scenarios have an equal probability of occurrence.It is obvious that, in different scenarios, the problem size of the proposed method is larger than the previous method, causing longer solving times.
Figure 6 shows the mean values of the daily total costs of the two different methods in different scenarios.Figure 7 displays the computing times of two different methods in different scenarios.It can be seen that the mean value of the objective of the proposed method fluctuates dramatically, and when the scenario number exceeds 23, the mean value of the objective tends to be stable, converging on 160,733 CNY.The corresponding computing time is approximately two hours.Thus, the scenario number was determined and set at 25, to ensure an accurate mean value of the objective.With respect to the previous method, the scenario number should be 40, so as to converge to 160,517 CNY, and the corresponding computing time is far longer than four hours.In addition, the relative error of the mean value of the total daily cost is only 0.13%.In this way, the computing efficiency and effect of the proposed method are superior to the previous method, as fewer scenarios and shorter computing times when converging, as well as lower relative errors, can be obtained from the proposed method.

Uncertainty Impact on Selection of System Facility
Table 5 shows the selection of system facilities in different scenarios, and it was found that the selections of different cases were the same.Natural gas, wind energy, and power were regarded as input-energy; internal combustion engines, gas turbines, and wind turbines were regarded as energy supply facilities to generate power and heat; water-source chiller/HP and heat exchangers were regarded as energy conversion facilities to generate heating and cooling energy by using part of the power from supply energy facilities.Meanwhile, the system still determined heating and cooling storage facilities.When the price of input-energy and the energy supply and demand were uncertain factors, the selection program and energy flow of facilities are shown in Figure 8.

Uncertainty Impact on Selection of System Facility
Table 5 shows the selection of system facilities in different scenarios, and it was found that the selections of different cases were the same.Natural gas, wind energy, and power were regarded as input-energy; internal combustion engines, gas turbines, and wind turbines were regarded as energy supply facilities to generate power and heat; water-source chiller/HP and heat exchangers were regarded as energy conversion facilities to generate heating and cooling energy by using part of the power from supply energy facilities.Meanwhile, the system still determined heating and cooling storage facilities.When the price of input-energy and the energy supply and demand were uncertain factors, the selection program and energy flow of facilities are shown in Figure 8.

Uncertainty Impact on Selection of System Facility
Table 5 shows the selection of system facilities in different scenarios, and it was found that the selections of different cases were the same.Natural gas, wind energy, and power were regarded as input-energy; internal combustion engines, gas turbines, and wind turbines were regarded as energy supply facilities to generate power and heat; water-source chiller/HP and heat exchangers were regarded as energy conversion facilities to generate heating and cooling energy by using part of the power from supply energy facilities.Meanwhile, the system still determined heating and cooling storage facilities.When the price of input-energy and the energy supply and demand were uncertain factors, the selection program and energy flow of facilities are shown in Figure 8.When the selection of the system facilities was optimized, under the premise that a short supply should be avoided, the price of input-energy, conversion efficiency of generating power, heating and cooling energy, facility price, etc., were considered, and multiple energy supply modes and operation strategies were comprehensively compared to select the most economical operation mode to satisfy user demands.In Case 2, the selected types of system facilities did not change compared with Case 1, while the capacity of cooling storage facilities increased by 1900 kW, and it could be determined that the price uncertainty of input-energy played less of a role in the implementation results.In Case 3, the uncertainties of energy supply and demand had a great impact on the distribution of system facilities compared with Case 1, and the capacities of gas turbines, heating and cooling storage facilities increased by 800 kW, 2024 kWh, and 3709 kWh, respectively.It attributed (to the system) the requirement to compare energy supply and demand of multiple scenarios and determined the maximum for preventing supply from falling short of demand.In Case 4, the uncertainties in the price of input-energy and the energy supply and demand influenced the system more significantly compared to Case 1, and the capacities of gas turbines, water-source chillers/HPs, and heating and cooling storage facilities increased by 1200 kW, 800 kWh, 9540 kWh, and 6053 kWh, respectively.

Uncertainty Impact on Economic Evaluation
Table 6 shows the daily costs in different cases.It can be seen that, in Case 2, the total cost for purchased power and natural gas contributed the most in terms of total daily cost increases, while facility investment and O&M costs contributed less to the total daily cost increase.This was because the mean power purchase price and natural gas price in Case 2 were higher than in Case 1, which played a leading role in power and natural gas cost variations.For Case 3, power and natural gas purchase costs barely contributed to daily total cost increases, while facility investment and O&M costs contributed the most.This is because, considering the uncertainty of energy supply and demand, the optimal capacity of some facilities could increase.For Case 4, the total cost of purchased power and natural gas contributed almost the same as facility investment and O&M costs to daily total cost increases.Every kind of cost was higher than in Case 1 when considering uncertainties, indicating that it exerted the greatest impact on the total cost when simultaneously considering the uncertainties of the price of input-energy and the energy supply and demand, followed by the uncertainty of the price of input-energy, and, finally, the uncertainty of energy supply and demand.Furthermore, the two types of uncertainties had cumulative effects on the total daily costs.When the selection of the system facilities was optimized, under the premise that a short supply should be avoided, the price of input-energy, conversion efficiency of generating power, heating and cooling energy, facility price, etc., were considered, and multiple energy supply modes and operation strategies were comprehensively compared to select the most economical operation mode to satisfy user demands.In Case 2, the selected types of system facilities did not change compared with Case 1, while the capacity of cooling storage facilities increased by 1900 kW, and it could be determined that the price uncertainty of input-energy played less of a role in the implementation results.In Case 3, the uncertainties of energy supply and demand had a great impact on the distribution of system facilities compared with Case 1, and the capacities of gas turbines, heating and cooling storage facilities increased by 800 kW, 2024 kWh, and 3709 kWh, respectively.It attributed (to the system) the requirement to compare energy supply and demand of multiple scenarios and determined the maximum for preventing supply from falling short of demand.In Case 4, the uncertainties in the price of input-energy and the energy supply and demand influenced the system more significantly compared to Case 1, and the capacities of gas turbines, water-source chillers/HPs, and heating and cooling storage facilities increased by 1200 kW, 800 kWh, 9540 kWh, and 6053 kWh, respectively.

Uncertainty Impact on Economic Evaluation
Table 6 shows the daily costs in different cases.It can be seen that, in Case 2, the total cost for purchased power and natural gas contributed the most in terms of total daily cost increases, while facility investment and O&M costs contributed less to the total daily cost increase.This was because the mean power purchase price and natural gas price in Case 2 were higher than in Case 1, which played a leading role in power and natural gas cost variations.For Case 3, power and natural gas purchase costs barely contributed to daily total cost increases, while facility investment and O&M costs contributed the most.This is because, considering the uncertainty of energy supply and demand, the optimal capacity of some facilities could increase.For Case 4, the total cost of purchased power and natural gas contributed almost the same as facility investment and O&M costs to daily total cost increases.Every kind of cost was higher than in Case 1 when considering uncertainties, indicating that it exerted the greatest impact on the total cost when simultaneously considering the uncertainties of the price of input-energy and the energy supply and demand, followed by the uncertainty of the price of input-energy, and, finally, the uncertainty of energy supply and demand.Furthermore, the two types of uncertainties had cumulative effects on the total daily costs.

Sensitivity Analysis on Uncertainty
In order to verify the effects of uncertainties on the model solution, we conducted sensitivity analyses on the convergence rate for the demand and energy price uncertainties, respectively.If energy supply and demand are uncertain factors, when the previous method starts to converge, the scenario number will increase constantly as the uncertainties in energy supply and demand are strengthened.However, the proposed method embeds the demand uncertainty into the model in the form of variance; therefore, the increase in the demand uncertainty had little effect on the convergence rate of the model; the results are shown in Table 7. Supposing that the price of input-energy is an uncertain factor, when the two methods converge, both scenario numbers will increase constantly, as the uncertainty of the input energy price is strengthened.However, the proposed method's convergence rate was faster than that of the previous method; the results are shown in Table 8.When the scenario number exceeds 50, the traditional method still fails to converge, and since the server memory is insufficient, it is difficult to predict the scenario number when the convergence is finished.This indicated that the price uncertainty was more sensitive to the convergence rate.However, the convergence rate of the proposed method was superior to that of the previous method, which showed that the former was more applicable for dealing with systems with complex uncertainties.

Conclusions
This study applied simple stochastic mathematical equations to illustrate the design optimization problems of a DER system.Then, a novel method for solving the optimization problem with multi-dimensional stochastic uncertainty was proposed.A MILP model of optimal design of a DER system was established based on the introduction of the mean value and variance, taking the minimum total cost of facility cost (investment and O&M), energy purchase cost, and loss caused by energy supply shortage as the objective, and considering the energy balance and facility performance constraints.Finally, this study took a DER system design of an office building district in Xuzhou, China, as an example, and four cases were investigated to verify the model: (1) all of the prices of input-energy and the energy supply and demand are fixed values; (2) the price of input-energy is perceived as an uncertain factor and the energy supply and demand are fixed values; (3) the price of input-energy is a fixed value and the energy supply and demand are uncertain factors; and (4) all of the prices of input-energy and the energy supply and demand are uncertain factors.Compared with previous research, our method proved to have a better convergence.The implementation results reflected the practicality and efficiency of the proposed method.However, at present, this method is only applicable to conditions where supply and demand meet the normal distribution, and there are different derivation processes and modeling methods for other distributions.This problem will be the focus of future research.

Figure 1 .
Figure 1.Superstructure of a DER (distributed energy resource) system.

σμ
shows the variance of input-energy; output-energy fulfills a shows the average value of output-energy, 2 s

Figure 1 .
Figure 1.Superstructure of a DER (distributed energy resource) system.

Figure 4 .
Figure 4. Environmental parameters of energy demand.

Figure 4 .
Figure 4. Environmental parameters of energy demand.

Figure 5 .
Figure 5. Operation principle of energy flow of each facility.

Figure 5 .
Figure 5. Operation principle of energy flow of each facility.

Figure 6 .
Figure 6.Mean value of daily total cost in different scenarios.

Figure 7 .
Figure 7. Computing time in different mean value of daily total cost.

Figure 6 .
Figure 6.Mean value of daily total cost in different scenarios.

Figure 6 .
Figure 6.Mean value of daily total cost in different scenarios.

Figure 7 .
Figure 7. Computing time in different mean value of daily total cost.

Figure 7 .
Figure 7. Computing time in different mean value of daily total cost.

Figure 8 .
Figure 8. Selection program and energy flow of facilities of Case 4.

Figure 8 .
Figure 8. Selection program and energy flow of facilities of Case 4.

Table 3 .
Parameters of energy storage facility.

Table 3 .
Parameters of energy storage facility.

Table 4 .
Problem size and computing time of Cases 2-4 with different numbers of scenarios.

Table 5 .
Selection of system facilities in different cases.

Table 5 .
Selection of system facilities in different cases.

Table 5 .
Selection of system facilities in different cases.

Table 6 .
Daily cost of the system in different cases.

Table 7 .
Comparison of fluctuation result for energy supply and demand.

Table 8 .
Comparison of fluctuation result for input energy price.

F
OXo,m If the energy generation facility o operates with input-energy m, F OXo,m = 1.Otherwise, F OXo,m = 0. F OXMo,mo If the energy generation facility o operates with mid-energy mo, F OXMo,mo = 1 Otherwise, F OXMo,mo = 0. F OYo,n,s If the n th product of energy generation facility o is out-energy s, F OYo,n,s = 1.Otherwise, F OYo,n,s = 0. F OYMo,n,mo If the n th product of energy generation facility o is mid-energy mo, F OYMo,n,mo = 1.Otherwise, F OYMo,n,mo = 0.

B
Ss,isIf the storage of output-energy .. with capacity is is selected for energy distribution system, B Ss,is = 1.Otherwise, B Ss,is = 0.B Oo,ioIf the energy generation facility o with capacity io is selected for energy distribution system, B Oo,io = 1.Otherwise, B Oo,io = 0. B ABse,t,s,a If variance of available output-energy s in season se during time-window t is in region a, B ABse,t,s,a = 1.Otherwise, B ABse,t,s,a = 0. B APLsa,se,t,s,a If probability of supply shortage of output-energy s in season se during time-window t is in region a, B APLse,t,s,a = 1.Otherwise, B APLse,t,s,a = 0.