Operational Planning of Energy for Non ‐ Interconnected Zones: A Simulation ‐ Optimization Approach and a Case Study to Tackle Energy Poverty in Colombia

: There is a need to look for alternative sources of renewable energy, especially in zones where people continue to live under energy poverty conditions. Consequently, to enhance the per ‐ formance of energy systems, algorithms to support planning decisions are required. This article proposes a simulation ‐ optimization framework to solve the stochastic version of the integrated en ‐ ergy dispatch and unit commitment problem for a solar radiation system operating in non ‐ intercon ‐ nected zones. Our study was motivated by challenges faced by a rural school located in Cundina ‐ marca, Colombia. Particularly, a simulation with optimization ‐ based iterations approach is used, modeling solar radiation as a random variable. The optimization phase uses a heuristic procedure that enables good solutions to be found in short computational times. To test our method, compu ‐ tational experiments were conducted using a set of randomly generated cases. The results suggest that our approach is useful and able to handle the random nature of the process for the school “Vol ‐ canes”. Additionally, we were able to quantify the impact that using a deterministic approach has on service levels for such systems. The novelty of the article lies in the proposed method and its application to a rural school with a low ‐ budget system.


Introduction
Renewable sources have been proven to be an efficient solution to increase the available energy in isolated zones [1][2][3] and reduce energy poverty [4,5]. According to the Sustainable Development Goals, by 2030, we should be able to provide affordable and clean energy for all [6]. However, with only ten years left, energy poverty continuous to be a concern [7]. In Latin America, between 10% and 20% of the population does not have electricity, and in Colombia, close to 1 million families in the rural sector fall into this category [8]. Consequently, renewable energies provide an alternative solution while avoiding the impact that climate change has on health and economic growth. For example, the Word Health Organization estimates that approximately seven million people die from air pollution every year and fossil fuels are the main source of air pollutants [9][10][11]. In addition, the use of renewable energies can reduce carbon dioxide emissions and improve people's health and quality of life [12,13].
In this context, solar energy has great potential in Colombia considering that average radiation is 4.5 kWh/m 2 [14]. However, policies to encourage its use are recent [15], and different barriers for their implementation have been identified [16]. Table 1 shows the effective capacity for electricity production by type of generation in Colombia. Particularly, in the department of Cundinamarca (Colombia), between 5% and 20% of the rural population does not have access to electricity. Therefore, security and sustainable development dimensions such as equal access to education and health services, among others, are at risk. Thus, local governments have recognized that there is a need for alternative ways to generate energy through sustainable technologies and social innovation mechanisms in this department. This work is motivated by the implementation of a solar energy generation system in a rural school called "Volcanes". The school does not have electric power and it is located at the municipality of Caparrapí, Cundinamarca. Figure 1 is a graphical representation of the proposed system. The source of energy is given by solar panels which can supply the batteries or loads through direct current/alternating current (DC/AC) power converters. The converter is designed to step up the voltage specifically for this particular case. The usual on/off loads needed in a rural school are considered, e.g., bulbs, computers, etc. The idea is to turn on/off the loads when it is necessary, taking into account the current energy. We propose an algorithm to accomplish two tasks: (i) to decide if energy comes from the batteries or directly from the solar panels and when the batteries are going to be charged (dispatching phase) and (ii) to address the problem of maximizing the number of on/off loads (devices) that can be connected considering a demand profile. This must take into account the constraints of charge in the batteries (unit commitment). Several optimization models have been used to support dispatching decisions. In Ashok [18], the author aims at finding the optimal combination of different energy sources to ensure balance and availability. In Pezzini et al. [19], the authors deal with energy distribution techniques in order to improve energy efficiency, and in Cormio et al. [20], the authors propose a methodology for regional energy planning. Although exact approaches are effective, recent research efforts have been focused in the development of heuristic and metaheuristic algorithms to obtain satisfactory solutions in good computational times [21][22][23][24][25]. In Zhou et al. [26], for example, the authors solve a problem of energy generation from wind and batteries based on genetic algorithms. More recently, in Khare and Rangnekar [27], the authors discuss the use of particle swarm optimization algorithms to find economic solutions in urban and rural regions.
Although renewable energy is the key to the transition to a carbon-neutral economy, its intermittent generation behavior, produced by the variability in climatic conditions, leads to the need for storage [28][29][30]. Therefore, uncertainty in solar radiation is a critical issue that should be researched in more depth [22,[31][32][33]. Despite having the potential of reducing energy poverty, solar energy has the disadvantage of being variable in a way that is not controllable [34]. In this context, assuming solar radiation is a deterministic parameter could impact service levels [35] and prevent nationwide projects from being implemented [36]. As a consequence, more research in the stochastic version of the problem is needed in order to support the decision-making process both at the strategic and operational levels [31,37].
Recently, in Mallor et al. [31], the authors proposed a mathematical model to find the optimal performance of a system with renewable sources and energy storage. Additionally, in Jeon and Shin [35], the authors use Monte Carlo simulation to manage photovoltaic energy and to perform long-term evaluation. We argue that the use of hybrid simulationoptimization approaches to model uncertainty in solar radiation is a promising research area. This method is a very effective mechanism to tackle variability in combinatorial optimization problems and its applications have been widely documented [38]. In Zeng et al. [39], the authors use this framework to reduce the risk of decision-making under uncertainties in the context of a water-food-energy plan, and in Ferrara et al. [40], the authors propose a similar approach to verify energy demand and supply it to a multi-family building. Lastly, in Roberts et al. [41], the authors incorporate a probabilistic simulation optimization approach to finding the optimal design of an energy system, taking into account the uncertainties in the sources, the load demand, and the unavailability of the components subjected to failure.
In this work, we present a novel solution approach for the stochastic version of the integrated problem of dispatching and unit commitment of energy generated by a photovoltaic system when it is not connected to the grid. Moreover, the solution is applied to a rural school with a low budget. Here, demand for electrical energy and its production are determined, taking into account solar radiation estimation and the storage capacity of the system. To the authors' knowledge, this simulation-optimization method has not been previously used to solve energy dispatch problems. Moreover, our deterministic solution approach, embedded in the simulation-optimization framework, is able to compute good solutions in much less time than exact optimization methods. This is a very important feature when taking into account the variability in solar radiation and the computing constraints of an isolated setting.
This article is organized as follows: In Section 2, we start with the description of the modeling framework and the mathematical model. Then, Sections 3 and 4 present the solution approach and computational experiments, respectively. Lastly, while Section 5 discuss the case study at "Volcanes", Section 6 presents our concluding remarks.

Problem Context and Mathematical Model
In this section, we first present the problem context and describe a decision framework for the planning process. Then, we formulate a linear programing model to support the decision-making process. The main limitation of using linear programming as a modeling approach is the need to know (or, at least, be able to accurately estimate) the solar radiation in the planning horizon. As stated before, there are a number of problems with this assumption; therefore, we propose to deal with this limitation by designing a simulation-optimization framework to solve the stochastic version of the problem and apply it to a real case in the school "Volcanes". However, the mathematical model allows us to assess the quality of the heuristic procedure (in terms of the objective function) and to quantify the impact of including the variability in the solution design (in terms of the service level).

Problem Context
Operational planning can be divided in two phases: dispatching and unit commitment. Figure 2 shows a graphic representation of the different stages of the implemented system. At the dispatching stage, the decision of energy generation for each planning horizon is modeled. Therefore, the amount of energy that will be sent to the devices and how much will be stored in the batteries need to be decided. Depending on the context, the planning horizon can range from one to seven days [39,40]. As can be seen in Figure 2, this model uses two known parameters: the capacity of electric power generation and the storage capacity of the batteries. Additionally, two parameters need to be estimated. First, a forecast model for the energy demand is produced based on historical data. Second, different scenarios of solar radiation are generated considering medium and high variability. Lastly, the performance of the systems is assessed using the service level and the percentage of missing of fulfilment.
At the operational level, the unit commitment must be executed. This consists of setting when loads have to turn on and off during the day. Consumption profiles and nominal consumption for each load are required as input parameters. Likewise, the priority of each load (ranging from one to three) is assumed to be known. For example, on a school day, between 7 am and 3 pm, computers would have priority over light bulbs; and on a non-school day, the additional room socket would have priority over the exterior bulbs. Since the capacity of the system may not be enough to meet all the requirements, priority plays an important role in user satisfaction.
The model framework shown in Figure 2 can be represented in a mathematical optimization model. However, it should be noted that radiation on the solar panels is a priori unknown. Therefore, it is modeled as a random variable [41], and it is required that the deterministic formulation be extended to a stochastic formulation. In order to do so, we used a hybrid simulation-optimization approach based on a heuristic procedure and Monte Carlo simulation [40]. Finding the optimal solution of the mathematical model, for all the solar radiation scenarios, requires high computational time. This disadvantage becomes important in the case of schools not connected to the grid. Therefore, finding good solutions in less time would ensure that priority is given to the essential loads in the school.

Mathematical Model
The dispatching phase can be modeled as a network flow problem on a directed graph, in which the panels, batteries, and loads are represented by nodes. Four sets are defined: days in the planning horizon ( ), periods within each day which can be minutes ( ), batteries (B), and loads ( ). In each period, energy supply ( ) and energy demand ( , ) are known. With this information, a binary parameter , can be built which takes a value of one when the demand is greater than the supply and battery energy is needed. Finally, it is necessary to know the maximum state of charge of the batteries and, also, the minimum . In this context, it is necessary to decide the amount of energy that will be sent between the panels and the batteries ( ), between the panels and the loads ( ), and between the batteries and the loads ( ). Additionally, auxiliary variables are constructed to represent the state of charge of the batteries ( ) and the energy available for the loads ( ). Finally, the binary variable takes a value of 1 if the state of charge of a battery b is greater than . Equations (1)-(7) represent the dispatching phase. First, binary variable is computed using constraints (1) to (3). Particularly, constraint (3) guarantees energy will be sent between the batteries and the loads only in those periods in which the charge in the batteries is greater than . is a big arbitrary number which is set for optimization purposes. C The following constraints ensure that the batteries deliver energy only in those periods in which the demand for energy is greater than the one given by the panels. Constraint (4) is only active when the parameter is equal to 1. If the battery can be used ( 1), the amount of energy given by the batteries to the loads must be less than the energy available minus the minimum state of charge. When the battery is not charged above the required limit ( = 0), no energy can be supplied.
Constraint (5) computes the state of charge of each battery for each period. Constraint (6) calculates the energy that will be available for the loads. Finally, constraint (7) ensures that the battery charge does not exceed its maximum capacity.
For the unit commitment phase, the status of each load (on or off) in each period must be decided. In order to do that, a binary variable and a set Ψ of constraints (8) and (9) are used. They allow us to coordinate this decision with the dispatching phase. Two parameters associated with the loads are assumed to be known: (i) the demand ( ) in each period and (ii) the priority they have ( ). The priority of loads can change between periods depending on the criticality of the task it supports. Constraint (8) ensures that only loads whose aggregate demand is less than the energy available in the period are turned on. Constraint (9) ensures that loads that are not connected cannot be turned on.
Finally, the objective function (OF) represented in Equation (10) maximizes the weighted satisfaction of the demand. The prioritization parameter, used to weight the demand, is of high importance in restricted resource scenarios in which it is not possible to turn on all the loads that have positive demand.
Summarizing, in order to determine the energy that is going to be used taking into account the supply, there is need for a dispatching phase. This is given by the cost in Equation (10) with the constraints (1)- (7). In order to determine which loads are turned, constraints (8) and (9) are added. This is the unit commitment phase.
The previous notation will be used in the algorithms of the next sections; however, we used decision structures in the algorithms to facilitate understanding of the implemented process and how it is related to the preceding equations.

Solution Approach
This section starts with the simulation-optimization framework to solve the stochastic version of the problem. The deterministic version then follows, presenting the designed heuristic algorithm.

Solving the Stochastic Problem: A Simulation-Optimization Framework
We run an optimization procedure for each iteration of a simulation. Figure 3 shows how deterministic and stochastic phases interact, in each iteration of the algorithm, in order to solve the stochastic version of the problem. First, we assume that solar radiation can be modeled as a random variable. In fact, other studies have found that solar radiation follows a beta distribution [42,43]. Based on this information, we generate realization of such a variable and solve the optimization problem using the heuristic algorithm. Thus, at this stage, the stochastic component of the problem is neglected. This solution is stored and compared with the ones obtained in subsequent iterations. In this context, the stochastic solution of the problem will be the most frequent deterministic solution in a set of simulation iterations. Therefore, we are aiming at finding the deterministic solution with the best performance for most of the simulation iterations. This form of interaction is called simulation with optimization-based iterations [44]. For comparison purposes, we decided that two solutions are similar if 90% or more (similarity threshold) of the loads share the same state (on/off) for all the time slots of the planning horizon. We conducted a parametric analysis to quantify the impact that changing the similarity threshold has in the number of required iterations for the algorithm. Therefore, the value that guarantees solution stability and minimizes the number of required simulation iterations was selected.
In order to assess the performance of the stochastic solution, a new simulation is needed. Therefore, for each iteration of the simulation, the confidence interval for the objective function is computed. The total number of required iterations for this simulation depends on the amplitude of such an interval. Additionally, the same procedure is used to assess the performance of deterministic solutions. In order to reduce variability of the comparison process, we used common random numbers [43] to quantify the impact of using the simulation-optimization approach.

Solving the Deterministic Problem: A Heuristic Approach
The deterministic version of the optimization problem can be solved using exact methods. However, computational time is a critical issue in scenarios where the decision will be made using a low-performance central processing unit [21]. Moreover, when tackling the stochastic version of a combinatorial problem, the optimization algorithm must use as few computational resources as possible [38]. Therefore, a heuristic procedure is designed to solve the optimization problem.
In this heuristic, dispatching and unit commitment phases are solved sequentially. Algorithm 1 starts by computing the difference between the available energy (provided by solar radiation) and the demand of the loads for each time slot. If the former is greater, the state of charge in batteries is increased. Otherwise, some batteries will supply energy to the loads while ensuring the minimum level of charge . Therefore, at the end of this stage, it is possible to know if the system will have enough energy to meet the total demand. Likewise, periods of charge and use of each battery are decided during this phase. Using this information as input, Algorithm 2 determines which loads are turned on and which are accordingly turned off for each time slot. First, loads are sorted in descending order of consumption. Then, following constructive logic, the algorithm decides to turn on the first load of the list and updates the available level of energy for the time slot. This procedure is repeated until all the loads are turned on or there is not more energy in the system. Lastly, loads are sorted according with the priority parameter.
In this section, algorithms are presented in decision structures because it makes easier to replicate the method.

Computational Experiments
In this section, we present the results of computational experimentation designed to test two hypotheses: (i) the heuristic procedure is able to find satisfactory solutions with low computational times and (ii) the use of a simulation-optimization framework improves the service level of the system.

Heuristic Performance
The mathematical model described in Section 2.2 was implemented using 21 randomly generated cases with different values of solar radiation and energy demand. Additionally, planning horizon was varied, ranging from one to seven days, generating three cases for each of these values. For this experiment, there are four batteries (| | 4), limit of state of charge of of = 288 Wh (50%), battery capacity of = 576 Wh and 16 loads | | 16. Table 2 shows information from both the mathematical model (i.e., exact solution approach) and the heuristic. Unsurprisingly, by using an exact solution approach (the mathematical model in Table 2), it was not possible to solve to optimality in any of the 21 cases. As we discussed in Section 1, this is an non-deterministic polynomial-time hardness (NP-Hard) problem, and recent research has been devoted to designing heuristic and metaheuristic approaches [21]. Table 2 shows the required time to find the first integer solution and the obtained value of the OF with a maximum running time of 6 h. As can be seen in Table 2, it was not possible to find integer solutions for the cases in which the planning horizon was seven days. This was also true for two additional cases with a planning horizon of 5 and 6 days.
Since it is expected that the dispatching decision of one day affects the parameters of the next one, this result highlights the need to use alternative solution approaches to solve the dispatching phase for at least one week. Table 2 also shows that the heuristic was able to find good solutions in reduced computational times. Using the heuristic approach, not only were all the cases solved but also, the weighted demand satisfaction improved between 2.5% and 15.2%. Additionally, computational times were less than two seconds for all instances. As we discussed in Section 3, the optimization procedure will be performed for each iteration of the Monte Carlo simulation when solving the stochastic version of the problem. Therefore, this reduction in computational times is a desirable feature of the deterministic solution.

Including Variability
To quantify the impact of using our simulation-optimization approach, we assessed the performance of both the deterministic and the stochastic solutions when the solar radiation follows a beta probability function. Therefore, the 21 randomly generated cases were adapted, including two different levels of variability: medium and high. While the deterministic version of the problem was solved using the expected value of solar radiation, the stochastic version was tackled following the algorithm described in Section 3.1. Then, a new simulation process was run, and a confidence interval was built for two performance indicators with a confidence level of 95%. The number of iterations for the simulation was calculated by setting a desirable amplitude of the confidence interval for both indicators.
In this context, performance was assessed using (i) service level, defined as the percentage of iterations in which it was possible to meet the energy demand and (ii) missing of fulfilment, the average percentage of non-fulfilment.
For example, Figure 4 shows that for case number one, using the deterministic solution, it was possible to meet the energy demand only in 35% of the iterations. At the same time, in those cases in which demand was not fulfilled, an extra 5% of energy was needed. As can also be seen in Figure 4, stochastic solutions show better service levels the two variability scenarios. While with medium variability, this performance indicator is around 99% on average, when variability is high, it decreases for larger cases (i.e., cases with longer planning horizons) reaching levels of 70%. Likewise, stochastic solutions show lower levels of missing of fulfilment, for all the cases and variability levels.
One practical implication of an improvement on performance relates to customer satisfaction. Neglecting solar radiation variability could result in an overestimation of the system capacity during the dispatching phase and, therefore, mislead unit commitment decisions. Consequently, it might become impossible to use high-priority loads when needed. In this context, for non-interconnected zones, operational planning could play a key role in tackling energy poverty and improving user experience.

Operational Energy Planning at "Volcanes": A Case Study in Colombia
As mentioned before, the case study was carried out based on the school "Volcanes", located at the municipality of Caparrapí, Colombia. Table 3 shows the load demands, specifying the individual demand for the system. The code of each device identifies the consumption profile. The system consists of four solar panels and a set of batteries in which energy is stored. Solar panels can produce 260 W when receiving 1000 W/m 2 . At the maximum power point, they have a voltage of 31 V with a current of 8.39 A. The total area of each panel is 1645 m 2 and the batteries are lead-acid type with a nominal voltage of 12 V and capacity of 48 Ah. The priority of the loads is given by their consumption, taking into account the needs of the school and validated with the community. For example, during the week and between 7 am and 3 pm, computers (PC 1-5) should have higher priority than bulbs (BSC 1-4). During the weekend, the electrical outlet (TA) should have higher priority than external bulbs (EX). Figure 5 shows the consumption and nominal profile of the loads-this is the parameter , . Following the quantity of loads in Table 3, | | 16.  Historical data of a 39-week period (from August 2016 to June 2017) was used for input modeling [45]. Solar radiation was found to follow a beta4 distribution (Kolmogorov-Smirnov test [46] with 5% significance). As can be seen in Figure 6, in most weeks between 9:45 am and 2:00 pm, the maximum radiation values are obtained (1064-1549 Wh). Minimum radiation values are shown between 6:00 am and 7:00 am (0-90 Wh) and between 5:00 pm and 5:45 pm (0-87 Wh). It should be noted that there is no solar radiation between 12:00 am and 5:30 am and from 6:00 pm to 12:00 am. This is the energy supply , .

Consumption(Wh) Consumption(Wh)
Behavior of the radiation was analyzed month per month. An average of 145 Wh was obtained. The month in which more radiation was received is February 2017 with 185.49 Wh. This is 28% more than the average. The month with less radiation is March 2017 with 115.04 Wh. This is 30% less that the average. After the proposed simulation-optimization approach is applied, 1187 loads are switched on in a period of one week. It represents 55% of the total number of loads to be switched on. It was determined that 97.44% ± 4.96% of replicas can meet the demand. It means the solution is not able to turn on all the necessary loads in each period. However, loads with higher priority are turned on. This satisfies the basic needs of the school. Figure 7 shows the loads that can be turned on over a week. This solution can be replicated every week to take full advantage of the radiation level and the charge of the batteries installed. From Figure 7, it is observed that during the week the external bulbs (EXT1-4) are on between 5:00 and 5:45 am and between 7:00 and 8:30 pm. However, external bulb 1 (EXT1) will turn on additionally at 6:30 am or 6:45 am depending on the day. Regarding the classroom bulbs (BSC1-4), it is observed that BSC1 turns on from 7:00 am until 10:45 am every day of the week. The remaining classroom bulbs are on in a smaller proportion. This is because, during this period of time, there is sunlight, and it is not necessary to turn on all the classroom bulbs. It should be noted that on no-school days, only BSC1 is used. Of the five computers (PC1-5), PC1 is the one that turns on during the week. It turns on between 7:45 am and 2:45 pm. Other computers will be used during this time period, but should be turned off occasionally. The only computer that turns on is PC5 between 7:15 am and 2:45 pm. During the week, television (TV) is turned on day one from 9:15 am until 10:00 am. On no-school days, television is turned on from 9:00 am to 3:00 pm. Finally, the electrical outlet (TA) and the room bulb (BH) turn on between 5:00 am and 7:30 am and between 7:00 pm and 8:45 pm. It should be noted that each device has a priority of use that depends on the moment in time.

Concluding Remarks
For the deterministic version of the problem, an optimization model was formulated to support operational planning. Using exact solution methods, it was possible to obtain integer solutions only for 76.2% of a set of randomly generated cases. For the remaining 23.8%, no solution was found. Therefore, a heuristic procedure is designed improving objective function, on average, by 3.30% while decreasing computational time between 21998.28 and 21999.78 s. This is very important for the application since no big or expensive computers can be used in rural schools. In Colombia, these institutions have budget limitations for computational resources, partially explained by the lack of coordination between national and regional authorities [47]. A cheap solution like the one proposed here should also be easily scalable to institutions with a similar situation. Lastly, in order to tackle solar radiation variability, a simulation with optimization-based iterations was performed. Our results show that the stochastic solution improves the performance indicators of the system.
The approach presented in this article is promising for applications in which computational time is an issue and solutions need to be built with low use of resources. Using a heuristic procedure to solve the optimization problem implies a small sacrifice in the objective function compensated. However, this is compensated by the savings in computational effort. Moreover, the approach can handle random variables such as the uncertainty of the supplies. In this article, only solar radiation is considered. However, it can be easily extended to wind generation, which is another popular source of energy. Other applications for this method can be real-time control systems, in which optimization is important (model predictive control). Here, we tested the solution in a real rural school in Cundinamarca, Colombia, which means that it can be implemented in other rural school in the country or even in countries with similar cases. Funding: D.P. wants to thank to the Fondo de Ciencia, Tecnología e Innovación del Sistema General de Regalías (SGR), Gobernación de Cundinamarca (Colombia), convenio especial de cooperación No. SCTeI 016 de 2015 for its support to the Project entitled "Incremento en las estrategias sustentables en el uso del recurso de energía eléctrica para la población vulnerable en el depto. De Cundinamarca". Authors also thanks to Pontificia Universidad Javeriana and the project SIAP-9638.

Conflicts of Interest:
The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.