Elektrotehniški vestnik 86(4): 185-191, 2019 Original scientific paper Discrete event simulation of distributed energy supply systems Dušan Kragelj1, Tomaž Kramberger2, Simona Sinko2, Bojan Rupnik2 1 Celje General Hospital, Oblakova ulica 5, 3000 Celje, Slovenia 2 University of Maribor, Faculty of logistics, Mariborska cesta 7, 3000 Celje, Slovenia E-pošta: bojan.rupnik@um.si Abstract. The paper presents a Discrete Event Simulation (DES) model to analyse the Distributed Energy Supply Systems (DESS) in order to increase the overall system efficiency. DES is used for assessing the energy demand at each participating source and overall energy demand for the entire system. The energy-flow through the network of energy suppliers, consumers, heating and cooling storage facilities, and power transmission lines is considered as a flow of goods known from the material logistics. In the DES model different tariffs of purchasing and selling the energy in case of energy deficit or surplus are applied. The model considers different tariffs for the different time-frames during a particular day. To increase the DESS cost efficiency, the energy distribution is adapted to the energy demand, demand-meeting and cost depending on time-frames. The model is calibrated based on historical data and gives promising results for identifying the optimal time-frame for the various components operating in the energy supply network. Keywords: Co-generation, Discrete event simulation, Distributed energy supply system Simulacija diskretnih dogodkov sistemov razpršene oskrbe ž energijo V članku predstavljamo analizo sistemov razpršene oskrbe z energijo s pristopom simulacije diskretnih dogodkov. Simulacijo uporabimo za oceno energetskih obremenitev vsakega vira in skupne potrebe celotnega sistema. Na energetske tokove skozi omrezje virov, potrošnikov, hladilnih in ogrevalnih sistemov ter elektricšnega omrezšja gledamo iz vidika lo-gisticne oskrbe. Simulacijski model upošteva razlicne nabavne in prodajne tarife, ki pridejo v posštev tako pri presezških kot tudi pri energetskih primanjkljajih. Analiziran model upošteva razlicne tarife v casovnih okvirih tekom dneva. Za izboljšanje stroškovne ucinkovitosti energetskega sistema se porazdelitve obremenitev prilagodijo potrebam, strosškom in porabam posameznih virov glede na casovne okvirje. Simulacijski model smo kalibrirali na podlagi dolgotrajnih meritev, kar s cimer je omogocena identifikacija optimalnih casovnih okvirjev za upravljanje gradnikov omrezja energetske oskrbe. 1 Introduction Distributed generation can play an important role for future energy systems. It can provide economically effective reduction of the greenhouse gas emissions of an energy supply and possible diversification of the primary energy sources where alternative fuels, such as natural gas or even solar and wind energy, are utilised [1]. Co-generation of electricity, heat and cooling power represents an established option for development of highly efficient and cost-effective integrated energy sys- Received 7 June 2019 Accepted 14 August 2019 tems [2]. The multi-energy customer demand-meeting also promotes the multi-energy systems integration for energy supply and utilisation of renewable energy hybrid systems by the distributed energy sources [3]. One of the critical factors that affect the energy saving, economic and environmental performance of such a combined cooling, heating and power systems (CCHP) are operation strategies [4]. The operation strategies tend to reach a twofold objective: minimizing the total operational cost and finding the optimal size of a CCHP system. So, a rational method for determining the size and operation strategy of CCHP system has to be developed [5] in the near future. Since the typical co-generation users are commercial malls, hospitals, office buildings and other public facilities, it is not always possible to use same tools for optimizing energy utilizations as at the production facilities. Usually raw data from production planning and control systems can be used to perform the energy load prediction in production facilities. In those cases, the historical data of the energy demand-meeting are the only source of data. The paper presents a DES model for analysing the energy workload requirements for a particular time of the day and a particular exterior temperature. 2 Literature review In the recent years, distributed generation of electricity, heat and cooling power has become an interesting topic 186 KRAGELJ, KRAMBERGER, SINKO, RUPNIK for researchers. Mathematical programming as such has been the most widely applied tool in energy planning applications. Mavrotas et al. introduce an energy planning framework combining Mathematical Programming and Monte Carlo Simulation. It can be applied in buildings of the Services Sector like hospitals, hotels, sport centres, universities and elsewhere where the uncertainty in the cost parameters is expressed by probability distributions [6]. Simulations are widely employed to provide a means for optimization for the distributed energy resource management [7]. Wu et al. [8] present an approach to simulate the CHP operations and to verify simulation results through an experiment, while Farmani et al. [9] propose a conceptual CCHP simulation model. An optimized schedule for generators and other resources is investigated by Hawkes et.al [1]. They use a deterministic linear formulation of the problem and system design. A linear programming cost-minimization problem for the system design and corresponding unit commitment of generators and storage within a micro-grid is introduced. A two-stage decomposition-based solution strategy to solve the optimisation problem using a genetic algorithm performing the search in the first stage and a Monte Carlo simulation dealing with uncertainty in the second stage are introduced by Zhou and co-authors [10]. They claim that an optimal design under the uncertainty does not deviate from the deterministic one. The preferred deterministic model is due to its computational efficiency. On the other hand, a comprehensive input-output matrix approach for modelling of co-generation equipment, which takes into account the interactions between the plant components and external energy network, was presented by Chicco et al. [2]. The results obtained from a numerical example show the modelling effectiveness of the proposed formulation. The study performed by Li et al. [4] analyses five different operational strategies such as: following the electric load, the thermal load, a hybrid electric-thermal load, the seasonal operation strategy, and the electric-thermal load. The authors claim that the redundant energy is not the best criterion to evaluate the CCHP systems under the current configuration. Understanding the different energy demands during the day is analysed in [11] where the authors observe short-term and long-term loads of an airport CCHP and use linear optimisation to reduce the operation cost, with the demand-meeting monitored on an hourly basis and the data used for the demand-meeting forecasting. Among the other simulation techniques Kohl et al. [12] present a DES approach for a production-based energy demand-meeting. 3 Methodology The aim of our paper is to build a hybrid DESS model capable of predicting and showing the energy load of each component of the given DESS for any selected time window throughout the year, depending on the external temperature. The proposed hybrid model consists of a discrete event simulation (DES) model of the DESS system which models the system functioning and a statistical regression to calculate the dependence of the energy demand on the external temperature. 3.1 DES model When dealing with the energy load prediction in production facilities, the data used in production planning and control systems are employed. But in facilities like hospitals, schools and other public facilities, such data are not available. In those cases, historical energy demand-meeting data are the only source available. In our investigation this problem is addressed using a DES methodology. We propose an approach where the energy is presented by logistic packages, with every package presenting a unit of energy. The simulation data input presents the DESS-produced and consumed packages. Discrete events present any changes in the energy production output and demand-meeting at different timeframes. Figure 1. DESS configuration In the DES model (Figure 1) the inputs present the CHHP heat and electricity production and energy demand for heating, cooling, and electricity. The grid electricity supply or heating using the a gas-fired boiler is determined depending on the deficit or surplus of the CCHP energy production based on demand-meeting. While the system energy demand requires constant heating for technological services, the cooling demand-meeting mainly depends on the external temperature. Heating demand (DH) presents a combined demand for hot water for heating purposes and technological services. A high-efficiency gas-fired boiler (WHB) is used only when the CHP heating, (WH C) does not meet the heating demands, with the CHP running at a constant load to meet the electricity demand. Heating surplus (SH) is used for cooling purposes. The consumed heat- DISCRETE EVENT SIMULATION OF DISTRIBUTED ENERGY SUPPLY SYSTEMS 187 ing energy (WH) is then: Whc - S h; Whc > D h Whc + Whb ; Whc < D h Wh = (1) It is obvious that when WHC > DH, WHB = 0 and SH = Dh - WHC. But when WHC < DH, WHB = Dh - Whc and Sh = 0. The cooling demands (DC) can be met by either the absorption heat pump (WCAhp) or from additional cooling energy from the compression heat pump (WCChp). The cooling surplus (SC) may be generated by a sufficient thermal or electrical surplus directed to the compression heat pump. The consumed cooling energy (WC) is: Wc = WcAhp — Sc ; WcAhp + Wcchp; WcAhp > Dc WcAhp < Dc Wc c Ahp Se * k Ap (2) (3) We WEC — SE ; WEC > DE , Wcchp < Dc -Wec + WEGin; WEC < DE W, cAhp WE EOhp = Wc cchlp/kchp (4) (5) window is checked. Next, the linear and polynomial correlation are checked and the latter is used due to its greater r2. The total individual energy loads is assumed to equal the total individual energy demand. When the correlation is sufficient, the regression curve is calculated: Dh = f (T ),Dc = f (T ),De = f (T ) (6) When WcAhp > Dc, Sc = Dc - WcAhp. But when WcAhp < Dh, Wcchp = Dc - WcAhp and Sc = 0. The efficiency coefficient of the absorption cooler is kAhp = 0.72. The electricity demand (DE) depends on both the cooling as well as on the service demand. An insufficient thermal surplus to meet the cooling demand requires utilization of the compression heat pump. A sufficient electricity surplus (SE) from the CHP facility may cover a possible compression heat-pump energy demand-meeting, however, a deficit requires an additional electricity from the power grid. The electrical demand-meeting is then: At that point DES model inputs can be evaluated. Using a set of 24 temperatures (T) measured each hour of a particular day as an input of the regression model provides the demands of an individual energy demand distribution (DH,DC,DE) during the 24 hours of a particular day. 3.3 Hybrid model Both the DES and the regression model are combined in order to provide a workload simulation based on historical demand-meeting data (Figure 2). Input: External temperatures T in Celsius for each hour within a particular time window. Output: Energy load of each component of DESS for each hour within a particular time window. i. The regression model calculates the demands (Dh ,Dc ,De ); ii. The DES model calculates the energy loads of the DESS components according to the input demands; iii. Both calculations are repeated for each hour within a particular time window. Also, when Wec > De , Wcchp < Dc - Wcahp, Se = De - WEChp and when Wec < De, Wwoin = De - WEChp,SE = 0. The compressor heat-pump efficiency is kChp = 4. As seen in Figure 1, WEGCinpresents the required electricity for the compression heat pump and WEEUin required power for other electricity workload. The electricity surplus that is not used for the cooling purpose is returned to the grid (WEGout = SE - WEChp). 3.2 Regression model The inputs for the DES model presented above are the values of demand for an individual energy (Dh , DC, De) during 24 hours of a particular day. As the energy demand-meeting is assumed to be correlated with external temperature it is necessary such an assumption is checked. First the correlation between the series of measured external temperatures (T) and the measured total energy loads within the selected time Figure 2. Hybrid simulation model. 4 Results 4.1 Regression model results The hourly temperature readings taken in June 2015 are examined. This term is selected for being typical 188 KRAGELJ, KRAMBERGER, SINKO, RUPNIK for an early summer month. The input data are hourly temperatures for each single day. Based on these data mean temperature for each hour throughout the month is calculated. The monthly mean hourly temperatures are given in table 1. Table 1. The calculated mean monthly temperatures in °C Hour Mean temp Hour Mean temp. 0 17,3 12 23.5 1 16,8 13 24.1 2 16,3 14 24.5 3 15,8 15 24.9 4 15,3 16 25.1 5 15,0 17 24.7 6 15,2 18 23.9 7 16,3 19 22.7 8 20,1 20 21.5 9 23,9 21 19.8 10 22,8 22 18.7 11 23,3 23 17.9 A sample of the total individual energy demands (Dh ,Dc ,De) is taken. The Mean hourly energy demands are given in table 2. The energy units are given in kilowatt hours for each different energy type (cooling, heating, and electricity). After performing a correlation and regression analysis for each pair (T°C,De), (T°C,Dc), (T°C,Dh) the energy demand based on exterior temperatures is calculated as: (T°C, De) ^ y = 40, 8x + 42,1; R2 = 0,6804 (7) (T°C, DC) ^ y = 42,7x - 236,7; R2 = 0,7413 (8) (T°C, Dh) ^ y = -0,83x3 + 48, 73x2 - 920,48x + 6410 (9) Figure 3 shows the demands obtained with the comparison between the regression model and the measured demand-meeting. 4.2 DES model Values Dh ,Dc ,De , used as an input in the DES model and shown in Figure 3 present a typical thermal, cooling and electricity demand for an early summer working day. The result of applying the DES model for 24 hours are given in Figure 4. While the regression model uses the data for one month, the simulation model follows the modelled hourly temperature and the temperature-based demands. The CCHP unit workload Table 2. Calculated monthly energy demands in kWh Hour Mean DH Mean Dc Mean DH 0 679.0 513.7 740.1 1 652.4 481.0 717.8 2 639.3 450.7 690.1 3 626.5 440.2 646.2 4 622.1 407.5 704.0 5 642.8 381.3 725.8 6 743.4 397.7 787.0 7 934.9 385.5 851.1 8 1053.9 435.3 860.5 9 1113.7 507.2 842.5 10 1147.0 670.7 876.0 11 1148.6 722.7 860.5 12 1149.1 803.7 839.9 13 1138.5 818.2 800.1 14 1076.5 841.5 732.3 15 1032.1 862.0 707.7 16 969.8 836.2 723.5 17 917.1 845.0 726.7 18 949.3 930.7 735.3 19 919.8 901.7 745.0 20 860.8 853.0 741.9 21 770.6 735.8 728.7 22 742.4 678.7 742.5 23 707.0 585.0 748.6 1200 200 0123456789 10 1112 1314 1516 171819 20 2122 23 - De=f(T°C)-Dc-f(T°C) -Dh=f(T"C) De ----Dc ----Dh Figure 3. Comparison between the calculated demands and measured energy demand-meeting. is constant, because the installed CCHP unit has the best efficiency maximum power, which is a constant load of 1100 kWh. Using the data of the early summer month i.e. June, the heat demand (DH) (Figure 4) is not high enough to engage the gas-fired boiler. Therefore the produced energy (WHB) is constantly zero. Hence, there is the possibility to use it as additional energy to DISCRETE EVENT SIMULATION OF DISTRIBUTED ENERGY SUPPLY SYSTEMS 189 drive the absorber heat pump. Its workload (WCAhp) is between 100 and 300 kWh. Since the cooling demand is the highest between noon and 8 p.m., the workload of the compressor heat pump (WCchp) is also the highest. As at that time window the CCHP does not produce enough electricity to meet all the electricity demand (DE) and still produces enough surplus to drive compressor heat pump to produce WCChp, electricity has to be bought from the grid as shown in Figure 5. On the other hand, in the morning and through the night, the demand-meeting for electricity is so low, that the energy is constantly sold to the grid. Figure 4. Calculated energy workload of the DESS components. Assuming that the electricity surplus is sold to the grid between 21 p.m. and 8 a.m. of the next morning, to cool down a building more than usual, by the internal temperature staying on a still comfortable level a little longer, meaning that less cooling is needed between 8 a.m. and 11 a.m., when electricity is bought from the grid. In the afternoon, the desired temperature is simply raised for 1 or 2 degrees. The variation with the cooling demand is presented as a dotted line in Figure 6. Figure 6. Calculated and modified cooling demand-meeting. Performing the DESS simulation provides the energy values of the surplus heating and electricity energy, purchased electricity and used natural gas for heating by the gas-fired boiler in order to harmonize the energy demand-meeting and production output by selling and buying less energy to and from the grid. Balancing between the bought and sold electricity presented in Figure 7 shows that a small change on the demand side causes a change on the demand-meeting side. While less cheaper energy is sold during the night, so is less cheaper energy is purchased during the day, resulting in a lower total cost. Figure 5. Purchased (red) and sold (green) electricity in kWh. 4.3 Discussion Using the presented DES model provides a tool capable of modelling the energy flow through the network. It shows the difference between the electricity bought and sold from or to the grid. Since the balance between the bought and sold electricity depends entirely on energy demand-meeting within a certain time window during the day, the best option of changing the balance is changing the demand. The model shows the workload of individual DESS components and the total electricity balance when an individual energy needs are changed (Figure 5), and when it is only the cooling demands that is taken into account. The system cooling demand-meeting depends on the external temperature and on the desired internal temperature. So, the higher is the external temperature, and the lower is the inside temperatures, the higher is the energy demand-meeting. Figure 7. Expected electricity demand-meeting of the calculated and modified demand-meeting. The workloads of appliances given in Figure 8 shows that by modifying only the cooling demands, the workload of the absorption heat-pump that relies on the heating surplus remains the same when the cooling demand changes. As expected, the only change is seen in the workload of the compressor heat pump (WCChp), and the electricity surplus is entirely used to run the compressor heat pump. Comparing the cooling demand from Table 2 and the external temperature from Table 1 shows that they are linearly correlated. 190 KRAGELJ, KRAMBERGER, SINKO, RUPNIK D 1 2 3 4 5 6 7 3 9 ID 11 12 13 14 15 16 17 IS 13 20 21 22 23 Wc ahp -Wcchp -Whb -Whc ----Wc ahp modfied ----Wcchp modified Figure 8. Expected electricity demand-meeting of calculated and modified demands-meeting. A similar effect can be observed when shifting the heat demand. When a facility requires heating for its service as well as for the temperature maintenance, shifting the service demand-meeting into a later time i.e. postponing the service requiring heat into later hours (Figure 9) when the external temperature is lower, provides a heating surplus that can be used for the absorption heat pump, thus reducing the electricity load of the compressor heat pump. However, due to a lower efficiency of the absorption cooler, such shifting of the heat demand provides an insignificant cost difference for the electricity supply. D 1 2 3 4 5 6 7 S 3 10 11 12 13 14 15 16 17 IS 19 20 21 22 23 x Time[h] -Dh ----Dh modified Figure 9. Expected demand of the calculated and modified heating demand-meeting. Similarly, any service that requires electricity can be postponed into a lower tariff hour, thus providing lower electricity cost. As the facilities benefiting from DESS generally require services for their personnel and customers, demand-meeting should not be viewed only from the cost point of view but also in terms of personnel comfort and, moreover, postponing services may raise other facility operating costs. 4.4 Conclusion The presented DES model provides a means to analyse the required load of each system component and can be adapted to variations in a particular heating, cooling, or electricity demand-meeting or in an external temperature. Depending on the CHP production, time of the day, and external temperature, the model identifies electricity, cooling, and heating deficits or surpluses at any node. With the gas-fired boiler and CHP running on the natural gas, the operating costs can be directly expressed by the cost of either the natural gas or electricity, while the latter needs to be adapted to a particular time of the day. Cost optimization can be performed by observing both the energy surpluses and deficits and by comparing the cost of acquiring either the grid electricity or increasing the gas-fired boiler output. The model can be further modified to include the use of ice banks and batteries for storing the energy surpluses and for cost-benefit analysis by comparing the cost of any possible system variation. References [1] A. Hawkes and M. Leach, "Modelling high level system design and unit commitment for a microgrid," Applied Energy, vol. 86, pp. 1253-1265, jul 2009. [2] G. Chicco and P. Mancarella, "Matrix modelling of small-scale trigeneration systems and application to operational optimization," Energy, vol. 34, pp. 261-273, mar 2009. [3] L. Ju, Z. Tan, H. Li, Q. Tan, X. Yu, and X. Song, "Multi-objective operation optimization and evaluation model for CCHP and renewable energy based hybrid energy system driven by distributed energy resources in china," Energy, vol. 111, pp. 322340, sep 2016. [4] L. Li, H. Mu, N. Li, and M. Li, "Analysis of the integrated performance and redundant energy of CCHP systems under different operation strategies," Energy and Buildings, vol. 99, pp. 231-242, jul 2015. [5] A. Costa and A. Fichera, "A mixed-integer linear programming (MILP) model for the evaluation of CHP system in the context of hospital structures," Applied Thermal Engineering, vol. 71, pp. 921-929, oct 2014. [6] G. Mavrotas, K. Florios, and D. Vlachou, "Energy planning of a hospital using mathematical programming and monte carlo simulation for dealing with uncertainty in the economic parameters," Energy Conversion and Management, vol. 51, pp. 722-731, apr 2010. [7] Z. Jun, L. Junfeng, W. Jie, and H. Ngan, "A multi-agent solution to energy management in hybrid renewable energy generation system," Renewable Energy, vol. 36, pp. 1352-1363, may 2011. [8] J. Wu, J. Wang, S. Li, and R. Wang, "Experimental and simulative investigation of a micro-CCHP (micro combined cooling, heating and power) system with thermal management controller," Energy, vol. 68, pp. 444-453, apr 2014. [9] F. Farmani, M. Parvizimosaed, H. Monsef, and A. Rahimi-Kian, "A conceptual model of a smart energy management system for a residential building equipped with CCHP system," International Journal of Electrical Power & Energy Systems, vol. 95, pp. 523536, feb 2018. [10] Z. Zhou, J. Zhang, P. Liu, Z. Li, M. C. Georgiadis, and E. N. Pis-tikopoulos, "A two-stage stochastic programming model for the optimal design of distributed energy systems," Applied Energy, vol. 103, pp. 135-144, mar 2013. [11] E. Cardona, P. Sannino, A. Piacentino, and F. Cardona, "Energy saving in airports by trigeneration. part II: Short and long term planning for the malpensa 2000 CHCP plant," Applied Thermal Engineering, vol. 26, pp. 1437-1447, oct 2006. [12] J. Kohl, S. Spreng, and J. Franke, "Discrete event simulation of individual energy consumption for product-varieties," Procedia CIRP, vol. 17, pp. 517-522, 2014. Dusan Kragelj received his B.Sc and M.Sc degrees from the Faculty of Mechanical Engineering and Faculty of Economics and Business, University of Maribor in 1980 and 2004. He is currently a Ph.D. candidate at the Faculty of logistics, University of Maribor. He is employed with the General Hospital Celje as an assistant to the general manager for the field of supply, investments and maintenance. His research work covers logistics and power engineering management. DISCRETE EVENT SIMULATION OF DISTRIBUTED ENERGY SUPPLY SYSTEMS 191 Tomaz Kramberger received his a Ph.D. degree in computer science from the University of Maribor in 2010, where he is an associate professor and holds the chair for quantitative modelling in logistics at the Faculty of Logistics. His research areas are quantitative modelling in logistics, graph theory, arc and vehicle routing, and port choice and port throughput analysis. He has cooperated with researchers of universities at Antwerp, Belgium, as well as the George Mason University in the USA. He is a member of the Slovenian Logistics Association, German Logistics Association and International Association of Maritime Economists. In 2012, he was a senior visiting research fellow at the Centre for Maritime Studies at the National University of Singapore and since 2014 at George Mason University in Virginia. He is a joint author of Sustainable Logistics and Strategic Transportation Planning. Simona Sinko received her B.Sc and M.Sc degree from the Faculty of Logistics, University of Maribor in 2016 and 2018, where is a currently teaching assistant. Her main research area is quantitative modelling in logistics. Bojan Rupnik received his Ph.D. degree in Computer Science in 2016 from the Faculty of Electrical Engineering and Computer Science, University of Maribor. Currently, he is a teaching assistant at the Faculty of Logistics, University of Maribor. His research interests include computer algorithms, optimization methods, simulations, and GIS.