R Rocchetta

Y F Li

E Zio

Risk Assessment and Risk-Cost Optimization of Distributed Power Generation Systems Considering Extreme Weather Conditions

Keywords: Distributed generation, AC power flow, extreme weather conditions, probabilistic risk assessment, Monte Carlo simulation, weather modelling

Risk assessment and risk-cost optimization of distributed power generation systems considering extreme weather conditions

R. Rocchetta, Yanfu Li, Enrico

Introduction

Existing power grids have been developed to meet the requirements of conventional single direction power delivery from centralized high-capacity generation units (e.g. thermal plants, nuclear power plants, etc) to various end-user loads (e.g. industry, commerce, residence, etc). The energy challenges faced by Europe and the rest of the world are changing the landscape of power systems.

Renewable energy resources, often geographically separated from the traditional power sources, are increasingly integrated into the distribution network in the form of distributed generators (DGs), such as photovoltaic panels and wind turbines. Owing to the random nature of these resources, DGs behave quite differently from conventional generators and they inject considerable amounts of uncertainty into power system operation; this uncertainty puts pressure on decision makers to properly assess the risk of the modern distribution networks integrated with DGs.

Unlike power system reliability assessments that focus on the evaluation of quantities such as system average interruption duration index (SAIDI), system average interruption frequency index (SAIFI) and expected energy not supplied (EENS) [START_REF] Brown | Electric power distribution reliability[END_REF] to reflect the ability to supply adequate electric service over the long term [START_REF] Kundur | Definition and classification of power system stability IEEE/CIGRE joint task force on stability terms and definitions[END_REF], probabilistic risk assessment (PRA) aims to estimate the probability (or frequency) of disturbances to system operation and their consequences [START_REF] Mccalley | Probabilistic security assessment for power system operations[END_REF]: these two elements are the constituents of the risk. Extreme weather conditions (e.g. high wind, thunderstorm, heavy snow, etc) can significantly affect system risk by increasing the frequency of failures of the power components and/or inducing severe damage [START_REF] Mccalley | Computation of contingency probabilities for electric transmission decision problems[END_REF].

In the past decades, many research works have been devoted to the risk assessment of power systems [START_REF] Mccalley | Probabilistic security assessment for power system operations[END_REF][START_REF] Mccalley | Computation of contingency probabilities for electric transmission decision problems[END_REF][START_REF] Gabbar | Risk-based performance analysis of microgrid topology with distributed energy generation[END_REF][START_REF] Javadian | Analysis of protection system's risk in distribution networks with DG[END_REF][START_REF] Arroyo | A risk-based approach for transmission network expansion planning under deliberate outages[END_REF][START_REF] Xiao | Power system risk assessment and control in a multiobjective framework[END_REF][START_REF] Alizadeh Mousavi | Blackouts risk evaluation by Monte Carlo Simulation regarding cascading outages and system frequency deviation[END_REF][START_REF] Guikema | Prestorm estimation of hurricane damage to electric power distribution systems[END_REF][START_REF] Javadian | Adaptive centralized protection scheme for distribution systems with DG using risk analysis for protective devices placement[END_REF][START_REF] Nateghi | Comparison and validation of statistical methods for predicting power outage durations in the event of hurricanes[END_REF][START_REF] Volkanovski | Application of the fault tree analysis for assessment of power system reliability[END_REF]. A number of studies have focused on transmission systems [7-9, 13, 31, 32, 39]; distribution network risk analysis [START_REF] Javadian | Analysis of protection system's risk in distribution networks with DG[END_REF][START_REF] Javadian | Adaptive centralized protection scheme for distribution systems with DG using risk analysis for protective devices placement[END_REF] has also been performed to analyze the response of protection devices/systems. Volkanovski, Cepin and Mavko [START_REF] Volkanovski | Application of the fault tree analysis for assessment of power system reliability[END_REF] have studied power grid reliability by fault tree analysis, considering voltage drop and power flow. Differently, Guikema et al. [START_REF] Guikema | Prestorm estimation of hurricane damage to electric power distribution systems[END_REF] and Nateghi et al. [START_REF] Nateghi | Comparison and validation of statistical methods for predicting power outage durations in the event of hurricanes[END_REF] have focused on the estimation of hurricane damage on distribution networks, using statistical tools to account for historical data. More recently, Gabbar et al. [START_REF] Gabbar | Risk-based performance analysis of microgrid topology with distributed energy generation[END_REF] have proposed an integrated framework for risk-based performance analysis of microgrids with DGs installed.

To the authors' knowledge, none of the existing works have considered the impact of extreme weather conditions within the framework of the risk assessment of distribution networks or the optimization of the DGs nominal power considering these conditions. Recently, Alvehag and Söder [START_REF] Alvehag | A reliability model for distribution systems incorporating seasonal variations in severe weather[END_REF] have conducted a reliability assessment for a distribution system considering the influence of extreme weather events (e.g. high wind and lightning). More specifically, in their model the extreme weather events affect the system by causing the overhead lines to fail. Kirschen and Jayaweera also remarked that line performance can be significantly affected by weather conditions [START_REF] Kirschen | Comparison of risk-based and deterministic security assessments[END_REF].

In this paper, we originally develop a simulation-based probabilistic risk assessment framework of DG systems that also considers severe weather. Based on the indications found in literature, we consider high wind and lightning as two major threats that can significantly increase the failure rates of distribution lines.

Furthermore, it we consider the optimal integration of DG within the power grid, which can provide several benefits (e.g. reduced power losses and improved voltage profile) [START_REF] Alarcon-Rodriguez | Multi-objective planning of distributed energy resources: A review of the state-of-the-art[END_REF]. Optimal integration of DGs needs to consider multiple conflicting objectives on which the decision makers must find satisfactory trade-off solutions. Mena et al. [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF] optimize the allocation of DGs in a reliability-cost bi-objective framework of simulation and optimization. Niknam et al. [START_REF] Niknam | A modified honey bee mating optimization algorithm for multiobjective placement of renewable energy resources[END_REF] optimize the size and allocation of DGs considering objective like minimizing costs, emissions and losses. In this paper, we propose an innovative risk-cost optimization for DG sizing, with the bi-objective of minimizing risk considering normal and extreme weather events, and the system investment and operative costs.

The rest of the paper is organized as follows. Section 2 presents the risk definition, the severity functions and the distribution line failure probability models, taking into consideration the two environmental threats of high wind and lightning. Section 3 describes the weather modeling and the power component modeling. Section 4 presents the Monte Carlo (MC) simulation procedure for risk estimation and the risk-cost bi-objective framework for optimal DG sizing. Section 5 describes the case study o a relatively complete DG system exposed to extreme weather conditions. Section 6 presents the DG system risk assessment and optimization results, and their analysis. Conclusions are presented in Section 7.

Risk Concepts

Definition of risk

We adopt a quantitative definition of risk as the product of the probability of occurrence of the undesired event (i.e. contingency) and the related consequence (i.e. severity) [START_REF] Zio | An introduction to the basics of reliability and risk analysis[END_REF][START_REF] Čepin | Assessment of Power System Reliability: Methods and Applications[END_REF]. To take into account more than one undesired event [START_REF] Čepin | Advantages and difficulties with the application of methods of probabilistic safety assessment to the power systems reliability[END_REF], the definition is extended by summing all contributions as:

(1), where p( ) is the probability of occurrence of the undesired event and is the severity of the related consequences. In probabilistic risk assessment (PRA), contingencies frequencies are used as probabilities and severities functions as consequences [START_REF] Mccalley | Probabilistic security assessment for power system operations[END_REF]. In the context of power systems, contingency is defined as the unexpected loss of one or more elements (e.g. distribution line, transformer or generator) comprising the power system [START_REF] Mccalley | Computation of contingency probabilities for electric transmission decision problems[END_REF]. Over-load, related with the feeders thermal limits, and bus voltage magnitude, related with frequency and system balance, are both indicators of power system stress and are used to represent the consequences for the risk calculation [START_REF] Xiao | Power system risk assessment and control in a multiobjective framework[END_REF]. Thus, the risk index associated with one contingency can be expressed as follows for the whole power network:

(2), where χ is the set of all operational and environmental conditions (e.g. wind speed, ground strike density, solar irradiation, temperature), C i is the i-th contingency, is the overload severity for line k in the conditions of C i and χ, is the low voltage severity for the node (or bus) b, is the risk associated with overload, is the risk associated with low voltage, L is the total number of lines in the system and B is the total number of nodes in the system. The composite risk due to all contingencies is, then, obtained as:

(3),

where N is the total number of contingencies. The severity functions and the probability models adopted are illustrated in the subsequent Sections 2.2 and 2.3, respectively.

Severity functions

The low voltage severity function measures the extent of a violation in terms of voltage magnitude drop at one node. There are three types of severity functions: continuous, percentage and discrete [START_REF] Ni | Online risk-based security assessment[END_REF]. The one selected for our study is the continuous function, because it measures the extent of the violation by reflecting the realistic sense that a performance close to, but within a performance limit, is, in fact, risky [START_REF] Ni | Online risk-based security assessment[END_REF].

The continous low voltage severity function adopted is as follows [START_REF] Mccalley | Probabilistic security assessment for power system operations[END_REF] : The severity function for overload is specifically defined for each circuit (distribution lines and transformers) and it measures the extent of violation in terms of excessive power flow as the percentage of rating (PR). The mathematical expression for this severity in the line k is presented as follows :

(4), (5)
(6).

In line with [START_REF] Ni | Online risk-based security assessment[END_REF], the deterministic limit for violation is , the near violation region is , the value under 0.9 is regarded as safe, , d = 10 and .

Probability of contingency considering extreme weather conditions

The probability of each contingency takes into account the failure of a distribution line by a Poisson distribution function [START_REF] Ni | Online risk-based security assessment[END_REF]: [START_REF] Arroyo | A risk-based approach for transmission network expansion planning under deliberate outages[END_REF], where is the probability of line contingency k in the next 1 hour, is the failure rate of the distribution line [failure/(h*km)], assumed constant, and is the length of the k-th line. Different from the existing PRA studies and in line with [START_REF] Alvehag | A reliability model for distribution systems incorporating seasonal variations in severe weather[END_REF], in our study account also for extreme wind speed and lightning events: For , the following expression is adopted [START_REF] Alvehag | A reliability model for distribution systems incorporating seasonal variations in severe weather[END_REF] : [START_REF] Alizadeh Mousavi | Blackouts risk evaluation by Monte Carlo Simulation regarding cascading outages and system frequency deviation[END_REF], where ɤ1, ɤ2 [s/m] and ɤ3 are fitting parameters and [m/s] is the critical wind speed value. If the wind speed is higher than , the failure rate of the line grows in an exponential way, whereas if the wind speed is below this value the contribution due to the wind is null.

Lightning activities generate a 'lightning failure rate' due to the increment of ground strike density per square km and hour. The relation is statistically linear and, like for the wind speed, if the value is lower than a defined ground strike density threshold, in our case set to 0, the lightning activity does not contribute to the failure rate, i.e. [START_REF] Alvehag | A reliability model for distribution systems incorporating seasonal variations in severe weather[END_REF]: [START_REF] Guikema | Prestorm estimation of hurricane damage to electric power distribution systems[END_REF], where is the lightning flash density, i.e. the number of ground strikes per square km per hour and is the fitting parameter.

For , we take [START_REF] Alvehag | A reliability model for distribution systems incorporating seasonal variations in severe weather[END_REF]: [START_REF] Javadian | Adaptive centralized protection scheme for distribution systems with DG using risk analysis for protective devices placement[END_REF],

where is a constant value.

Weather and Power Components Modelling

Weather modelling

To consider distributed generation by renewable sources and the effects of extreme weather conditions on the failure behaviour of the system components, the factors to be modelled in our work include wind, solar irradiance and lightning. As often done in the literature, the Weibull distribution is used to model the random wind speed at time t; the probability density function (PDF) of wind speed is [START_REF] Li | Uncertainty analysis of the adequacy assessment model of a distributed generation system[END_REF][START_REF] Li | A multi-state model for the reliability assessment of a distributed generation system via universal generating function[END_REF]: [START_REF] Nateghi | Comparison and validation of statistical methods for predicting power outage durations in the event of hurricanes[END_REF],

where and are the shape parameter and the scale parameter, respectively, which can be obtained by fitting eq. ( 12) to historical wind speed data at time t.

The beta PDF used to model the solar irradiance is [START_REF] Li | Uncertainty analysis of the adequacy assessment model of a distributed generation system[END_REF][START_REF] Li | A multi-state model for the reliability assessment of a distributed generation system via universal generating function[END_REF] :

(13), (14),
where is the solar irradiance [kW/m 2 ] at time t, and and are fitting parameters.

The model for random lightning occurrence is obtained from historical data of mean occurrence hours per year (T light ) and mean ground flash densities (Ng). In the Monte Carlo procedure here adopted, the lightning occurrence is sampled from the occurrence probability obtained by dividing T light by the total hours in the year. When the lightning is sampled to occur, Ng is randomized uniformly as in [START_REF] Alvehag | A reliability model for distribution systems incorporating seasonal variations in severe weather[END_REF]. Thunderstorm events are obtained when both the wind speed is higher than and a lightning event is sampled.

Power components modelling

This Section presents the models of the generation and consumption components in the DG system, as described in previous work by some of the co-authors [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF].

Photovoltaic panel (PV)

The power output from a single solar cell is obtained using the following equations and the randomized irradiance [START_REF] Li | Uncertainty analysis of the adequacy assessment model of a distributed generation system[END_REF][START_REF] Li | A multi-state model for the reliability assessment of a distributed generation system via universal generating function[END_REF]:

(15), (16), (17),

Wind turbine (WT)

The WT uses the kinetic wind energy to generate electricity. The power output of one WT can be modelled as a function of wind speed as [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF]: [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF], where is the cut-in wind speed [m/s], is the rated wind speed [m/s], is the cut-out wind speed [m/s] and is the rated power [kW].

Electric vehicle (EV)

In this work, EVs are considered as battery EVs with three possible operating states: charging, discharging and disconnected [START_REF] Clement-Nyns | The impact of vehicle-to-grid on the distribution grid[END_REF]. The power output of one EV block, made by a number of individual EVs considered as a whole, [START_REF] Li | Uncertainty analysis of the adequacy assessment model of a distributed generation system[END_REF][START_REF] Li | A multi-state model for the reliability assessment of a distributed generation system via universal generating function[END_REF] is calculated using:

(21), (22),
where t is the hour of the day [h], the residence time interval for operating state op [h],

is the rated power [kW], is the discharging probability at time t, is the charging probability at t and is the probability of disconnection at time t.

Storage devices (ST)

The behaviour of the STs is similar to that of the EVs with only charge and discharge states. The level of charge is assumed to follow a uniform distribution and for the sake of simplicity we suppose that only discharge is allowed as the operative state [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF] :

(23), (24), (25)
, where is the level of charge in the battery [kJ], SE is the specific energy of the active chemical [kJ/kg], MT is the total mass of the active chemical in the battery [kg], is the rated power [kW] and is the discharging time interval [h].

Load demand (L)

Demand of power, overall and nodal, is described by the daily power demand curve (Fig. 2). At each hour, there is a specific mean power demand and variance (obtained from historical data). The uncertainty for the load profile is modelled by a normal distribution [START_REF] Atwa | Optimal renewable resources mix for distribution system energy loss minimization[END_REF].

Fig.2

The nodal demand of power is deducted from the overall demand in the network and is modelled as: [START_REF] Zimmerman | MATPOWER: Steady-state operations, planning, and analysis tools for power systems research and education[END_REF], where is the mean load at node b at time t, is the standard deviation at node b at time t, is the power demand (Load) at node b [kW], is the normal PDF of power demand in node b at time t, (-) is the cumulative distribution function (CDF).

Main supply (MS)

The MS spots in the distribution network are the power stations connected to the transmission system. The distribution transformers are located on these spots and provide the voltage level of the customers. In this study, the MS is modelled like a slack bus, as the generator used for balancing the production and the demand by injecting power into the DG system.

Simulation Procedures for the Risk Assessment and Multiobjective

optimization of DG systems exposed to extreme weather conditions

Non-sequential Monte Carlo

Non-sequential MC [START_REF] Li | Uncertainty analysis of the adequacy assessment model of a distributed generation system[END_REF][START_REF] Johansson | Reliability and vulnerability analyses of critical infrastructures: Comparing two approaches in the context of power systems[END_REF] is used to evaluate the DG system model and compute the composite risk (3), accounting for the variability of the environmental and operational states, . At each simulation run, first an hour of the day t is randomly sampled from a uniform distribution U [START_REF] Brown | Electric power distribution reliability[END_REF][START_REF] Subrahmanyam | A Simple Approach of Three phase Distribution System Modeling for Power Flow Calculations[END_REF]; then, the environment vector χ, including wind speed, lightning and irradiance, is generated in correspondence of that t. These actions consist of the first step in the flow chart shown in Fig. 3, which refers to a single MC iteration. At the second step, the DGs power outputs, the distribution lines failure rates and contingencies probabilities (eqs. ( 7)-( 8)) are obtained for the sampled χ. At the third step, one contingency is selected, the AC power flow is performed under the scenario induced by this contingency and χ, and the system performance is evaluated. At the fourth step, the probability of and the severity values relevant to are computed, and eq. ( 2) is used to compute . At the fifth step, is stored for . At the last (sixth) step, is obtained according to eq. ( 3). The procedure is repeated for a large number of times and the results are the contingencies probabilities, severities and risk values for different sampled scenarios.

Fig. 3 4.

AC power flow

Compared to DC power flow, which can be fast but less accurate, AC power flow provides a closer approximation to the real state of the system for the evaluation of the severity functions presented in Section 2.2. The single-phase AC power model is widely accepted as a high-quality approximation of the steady-state behaviour of real-world power flow [START_REF] Bent | Approximating Line Losses and Apparent Power in AC Power Flow Linearizations[END_REF]; more accurate but complex models, like AC three phase complete power flow [START_REF] Subrahmanyam | A Simple Approach of Three phase Distribution System Modeling for Power Flow Calculations[END_REF], can also be considered. Taking into account the balance between computation efficiency and approximation accuracy, this study adopts the singlephase AC power model.

Multiobjective optimization (MOO) of optimal power sizing

MOO problems are commonly defined for optimal DG allocation and sizing [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF]. In particular, riskbased multiobjective optimization (RBMOO) offers the capability of controlling system risk via trade-offs with the other values (e.g. economical) of the solutions, as explained in [START_REF] Xiao | Power system risk assessment and control in a multiobjective framework[END_REF]. We embed the proposed simulation-based PRA including extreme weather events within the risk-based multiobjective optimization of the sizing of DG installation in a distribution grid. The decision variables are the nominal power of the sources to be installed on each node of the grid, grouped in the decision vector:

DG system description

An 11-nodes, DG and radial system [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF], modified from the IEEE 13 nodes radial distribution network [START_REF] Kersting | Radial distribution test feeders[END_REF], is considered (Fig. 5 and Tab. 1) as example case study for the PRA modelling framework proposed. The spatial structure of the IEEE 13 network has not been changed, but we neglect the regulator, capacitor and switch, remove the feeders of zero length and we consider all the lines has overhead lines. The modifications are made so that it becomes of interest to consider the integration of renewable DG units and perform a PRA including severe weather conditions. In spite of its small scale, the systems are relatively complete. We assume that weather condition is homogeneous across all its components. MS is installed on the 1 st node and the configuration of DG (electric vehicles EVs, photovoltaic cells PVs, storages STs, and wind turbines WTs) follows an optimal solution found in [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF] by maximizing the reliability of power supply and minimizing the global cost of the system. The details are presented in Table . 1. For example, the entry at 2 nd row and 3 rd column indicates that there are 62 PVs being installed on node 1.

Fig.5

Table 1 Table 2 presents the parameters of the DG models, of the weather models and of the failure rates.

Table 2

Fig. 6 shows the hourly probability data for the EV operation states.

Fig .6

The parameters of the distribution lines relevant to the AC power flow computation are presented in Table 3.

Tab. 3

An RBMOO has, then, been defined to aim at finding the optimal size and location of the renewable technologies considered, comparisons and result analysis are presented in Section 6.

Simulation setting

The main supply (MS) situated in node 1 is the slack bus (i.e., swing or reference bus [START_REF] Zimmerman | MATPOWER: Steady-state operations, planning, and analysis tools for power systems research and education[END_REF]) used for the AC power flow. It is typical to assign given reference voltage magnitude and voltage angle to this bus: in our case,

[kV] which equals to the voltage magnitude of the original case and , which is commonly set for the slack bus. The base apparent power and base voltage used for the per-unit system are 100 [MVA] and 4.16 [kV], respectively. Based on the base voltage, p.u. is used, following [START_REF] Zou | Voltage support by distributed generation units and shunt capacitors in distribution systems[END_REF]. Note that the system does not have a protection device that allows an islanding mode, which means that a distribution line failure can result to part of the system being in black out. For example, the failure of line 24 can produce the black outs in nodes 4 and 5. The portion of the system in black out state has zero power flow and node voltage.

The Monte Carlo simulation framework proposed for the PRA has been applied to the network, one time with DG installed as explained in [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF] (Fig. 5) and another time without DG (i.e. the original radial distribution network). Comparisons of the results obtained has been made considering or not the effects of severe weather conditions in the failure model. The PRA for the case without taking into account the extreme weather conditions influence (equations ( 9) and ( 10) in Section 2 neglected) has been performed setting the failure rate at a constant value equal to .

The contingency list considered includes the 'N-1' overhead lines contingencies, similar to other assessments of literature [START_REF] Marsadek | Assessment and classification of line overload risk in power systems considering different types of severity functions[END_REF] and also considering the criticality of overhead feeders in stormy or windy environments [START_REF] Bent | Approximating Line Losses and Apparent Power in AC Power Flow Linearizations[END_REF]. The proposed model is flexible and allows the introduction of various contingencies for various components, generators included, but an exhaustive list would be out of the scope of this paper.

The NSGA-II was, then, applied to solve the DG sizing and allocation MOO problem introduced Section 5.1. Optimization parameters were set as follows: starting population size N=100, number of generations gen=50, MC runs equal to 200, crossover probability pc = 0.9 and mutation probability pm = 1/44, where 44 is the number of decision variables (power sizes in each node). The scenario time has been set equal to 1 hour and equals to 10 years. Upper power size limits and cost data are presented in Table 4. The lower limits for power size are fixed to zero.

Table 4

Results and Analysis

In this Section, we first compare the PRA results obtained considering extreme environments with those obtained under normal weather conditions (subsection 6.1). This is done for both case of DG or not DG systems, to explore the effects of extreme weather conditions on the system risk for the two system configurations. Secondly, we show the contingency ranking for the DG system considering extreme weather conditions (subsection 6.2), in order to identify the worst contingencies. After this, we compare with and without DG systems under both extreme and normal weather states (subsection 6.3), and this will show the lower risk and severity associated to the DG system. We indicate with the notation the unexpected failure in the connection between nodes b and nodes b'. Finally, we present the Pareto front and three optimal solutions for the generator sizing obtained by solving the RBMOO problem with NSGA-II (section 6.4).

PRA: Effects of extreme weather conditions

The expected risk values for the no-DG and the DG systems are shown in Figs. 7 and8, respectively. The risk is computed as the sum of all the recorded s at each iteration, divided by the number of scenarios sampled at the iteration. An increment in the expected risk could be attributed to some heavy environmental conditions. It is seen that the average risk has a sensibly higher value if lightning and high wind are included in the model, for both systems without or with DG: at convergence, the ratios are and respectively.

Fig. 7 Fig. 8 Dotted lines in the plots refer to the expected risk if only the lightning is considered; influence of the high wind can be derived by looking at the difference between this line and the dashed line (high wind and lightning).

Fig. 9(a) show the classical probability-low voltage severity plot for the contingencies of the DG system. Each dot in this plot is a risk value computed in a single MC run (step 5). It is observed that, as expected, contingencies C 12 and C 26 result in higher risks, due to their higher severities and larger probabilities. This is intuitive by taking into account the topology of the system: the lines 1-2 and 2-6 are two central lines and the failure of either will cause a black out of a considerable amount of nodes, resulting in a large value of severity; additionally, these two lines are the longest ones in the system, so that their failure rates [failure/h] have values higher than the others.

Fig. 9(b) displays the probability-overflow severity plot. We can observe that the severity values are mainly lower than those in Fig. 9(a), while the probability values are the same in the two Figures.

The differences in severity values can be understood by analyzing the failure of one branch in the system, e.g. line 2-6: if it fails, there is a high low-voltage consequence due to the large number of nodes with zero voltage (in black out), but in that part of the system the active and reactive power flows are zero, and hence the overflow severity is lower. 10 that a large number of violation scenarios occur. These are obviously the ones with high load demands and, hence, high power flows in the lines. Note that like some other contingencies, has scenarios with a zero overload severity due to the low power demand by the customers (e.g. in night hours). These scenarios are located inside the secure zone. Fig. 11 shows that each violation scenario leads to a DV: this is due to the fact that our model does not include protection schemes and no islanding mode is available. Hence, if we consider a contingency in line 2-4, it will cause a DV violation in nodes 4 and 5.

PRA: Contingency ranking

Contingency ranking is usually performed to identify the contingencies that cause the worst problems in the system [START_REF] Singh | Fast voltage contingency screening and ranking using cascade neural network[END_REF]. The proposed MC simulation allows retrieving a contingency ranking.

We use the percentage risk contribution as the ranking criterion and show its values under various load demand levels for each contingency. The x-axis, in Figs. 12(a-c) and 13(a-c), represents the power demand in [kW] and the y-axis shows . Figs. 13(a-c) displays the quadratic fit for the percentage risk contributions. It is observed that the rank changes with different load levels. For example, under stressed condition (high load demand), the risk associated to C 810 is larger than that associated to C 611 and C 68 . If lower power demand is considered, higher risks are due to C 611 and C 68 than to C 810 . For the majority of sampled scenarios, is the third higher contributor to the risk, with value around 3% as shown in Fig. 13(a). The first and second worst contingencies are ≈52%, as shown in Fig. 13(c) and ≈ 31.5%, as shown in Fig. 13(b). As to C 68 , who has the third highest low voltage severity (in Fig. 9), lower relevance is partially due to the low probability of occurrence.

PRA: Comparisons between the without and with DG systems

As complement to Section 6.1, we compare the risks of the two systems under different weather conditions. Under normal weather conditions, risk is very low for both systems (Fig. 14): the results

show that the risk of the network without DG ( ) is slightly higher (~1%) than that of the system with DG ( ). We also analyze the difference between the two systems under extreme weather conditions. Fig. 15 shows that the risk value of the system without DG ( ) is approximately 1.3% higher than that of the with-DG system ( ),

if extreme weather conditions are taken into account. The risk reduction in the with-DG system is mainly attributed to the improved mean voltage profile, as displayed in Fig. 16. In general the risks for both systems are much higher under the extreme weather conditions than under the normal conditions. The comparisons between the without-and with-DG systems reveal that the DG installation has reduced of about 1.3% the risk under extreme weather conditions and of about 1% the risk under normal conditions. extreme weather conditions. Fig. 18 shows three solutions within the optimal front, with three levels of risk and cost values: is the solution of highest expected risk and lowest expected cost, is the best compromise solution obtained using the min-max approach [START_REF] Belegundu | Optimization concepts and applications in engineering[END_REF] to compromise between expected risk and expected cost, and is the solution of lowest expected risk and highest expected cost.

Figs. [START_REF] Li | A multi-state model for the reliability assessment of a distributed generation system via universal generating function[END_REF][START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF][START_REF] Clement-Nyns | The impact of vehicle-to-grid on the distribution grid[END_REF][START_REF] Atwa | Optimal renewable resources mix for distribution system energy loss minimization[END_REF] show the results for the sizing of different DG types on each node. The y axis contains the node number and the x axis represents the related power size . We also show in Fig. 23 and Table 5 that, coherently with the model, the lowest cost solution corresponds to one with lowest power installed and the highest solution corresponds to the one with higher amount of power installed, whereas higher risk correspond to lower amount of power installed. The x axis in Fig. 23 represents the global nominal power installed in the network and the y axis reports the DG technology type k. The amount of PV power installed is lower than other technologies: this is explainable if one considers the investment cost , which is one order of magnitude higher than for the other technologies.

Fig.5

Table 5 shows the global power installed in each node, if we compare the less risky solution with the more risky one, the tendency is a global increment of the power installed in each node. We can also observe that a less risky solution allocates more power in a central node for the system (node 6)

improving the voltage profile and reducing the flow of power in the central lines. The RBMOO solutions appear to be more expensive but less risky than the original ones assessed in the preceding sections. Comparison values are reported in Tables 6 and7. Table 6 contains the expected risk and expected cost values. In Table 7, we show the risk reduction for the four systems with-DG compared to the without-DG system. We can observe that the optimal solutions, accounting for the extreme events, further reduce the risk than the one for which extreme events are not considered. Compared to the without-DG case, the cheapest solution reduces risk of almost 1.8 % if extreme weather conditions are considered (around 1.5 %, without accounting for the extreme weather events) whereas the more expensive solution reduces it more than ~2.56 % (~2.05 % without considering the extreme weather events). The comparisons between the RBMOO optimal solutions and the system in Table 1 reveal that they are actually non-dominated to each other with respect to the expected risk and the expected cost.

To investigate the effects of DG installations significance tests have been performed to compare the risk values of the three Pareto optimal DG systems and the radial system. The Wilcoxon rank sum tests are used in place of the standard t-tests because the distributions of risk values are nonnormal [START_REF] Wilcoxon | Individual comparisons by ranking methods[END_REF].The results in Table 8 show rejections of null hypothesis of equal medians at 0.1% significance level, so that we can accept the alternative hypothesis that the DG installations in the three Pareto solutions do help to reduce the system risk. Note that we do not investigate the effects of severe weather conditions onto the system risk because one goal of this work is to include the severe weather into the system simulation model for closer adherence to reality.

Tab. 6

Tab. 7

Tab.8

Discussions

The results presented in the previous sections show percentage reductions for the DG systems risks, even if severe weather conditions are accounted for. These reductions bear intuitive economichuman relevance: the (lines) contingencies can lead to blackouts of one or several nodes, so that the customer power demands are not satisfied, thereby creating problems in both economic and human terms. We do not quantitatively express this , because it would require considering also greenhouse gas reduction, line power loss reduction, economic return due to the usage of renewable energy, and several others aspects, which would lead us out of the scope of the work.

Note that reduction of 1% or 1.3% is numerically small, but not irrelevant for the networks considered. For a significant number of scenarios, the system stress due to severe weather conditions is much lower in the DG system than in the radial system. This implies that the incorporation of DGs, with stochastic power production dependent on weather conditions, does not deteriorate the system even if severe weather conditions are accounted for. On the contrary, it slightly reduces risk.

Additionally, we should take into account also the specific setting of the system, with weakly connected nodes that cannot lead to numerically big reductions of the risk without also accounting for security schemes such as the islanding mode and different grid topologies. Such considerations are part of our intended future developments.

CONCLUSIONS

In the work presented in this paper, a simulation-based framework for the probabilistic risk assessment of DG systems, considering also extreme weather conditions. MC non-sequential algorithm has been implemented, for accounting of extreme wind speed and lightning events, whose effects have been originally included onto the probabilistic failure models of the system components. Application to a case study shows that under extreme weather conditions, there is an increment in the expected system risk, as expected. Severity-probability visualizations are presented and contingency ranking is performed, under various operational conditions. Comparison with a without-DG system confirms the benefits of DG installation in terms of bus voltage, line flows and post-contingency severities, particularly under extreme weather conditions.

The probabilistic risk assessment framework has been embedded within a multi-objective optimal power generation sizing framework to obtain the best nominal power size of different technologies to be allocated in each node of the DG system. The optimal solutions turn out to be best in terms of expected cost and risk, as computed according to the probabilistic risk assessment model proposed to incorporate severe weather.

The computational effort added is, compensated by the insights that can be gained from the simulation and optimization results. In this sense, a decision maker can be provided with riskinformed knowledge on optimal DGs allocation and the benefits not only of emission reduction etc., but also of system performance under severe weather conditions.

In synthesis, the original contributions of the work are:

PRA simulation framework:

-Inclusion of extreme weather conditions and their effects on the system components failure;

-Severity-probability visualization -Contingency ranking within for the assessment;

RBMOO framework:

-Inclusion of the extreme weather conditions by embedding the PRA within the MOO framework; Fig. 2 illustrative example of overall mean power demand [START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF] Fig. Power demand

[kW] t[h]

Step 1: operational and environmental state sampling (χ)

χ =[t, , ,Ng, , , ]
Step 2: using the formulas to obtain i=1

Step 6: Computing End of contingency inside the list (i=N+1)?

yes no

Step 3: Solving AC power flow in contingency ( ) state ( ,..., ,..... );( ,... ,... ); i=i+1

Step 4: Computing

Step 5: Recording

Fig.4 Flow chart of RBMOO-NSGA-II

Step 2: Randomize initial population p0 p 0 = ; the i-th element of the population is the power size of the generators in each node: with .

Step 1: Set Generation number gen, starting population size N and Monte Carlo loops.

Set Generation number gen, starting population size N and Monte Carlo runs.

Step 3: Perform MC steps from 1 to 6, as presented in Fig 3, to obtain the expected objectives for each chromosome of the initial population p 0.

Step 4: The initialized population p0 is sorted based on non-domination with respect to the objective functions; crowding distance is assigned to each element. Ranked non-dominated fronts are obtained

, where is the best front and the less good front.

Get g=0 and start evolutionary loop.

Step 5: The parents individuals are selected by using a binary tournament selection based on the crowed-comparison-operator: an intermediate population p g ' is obtained.

Step 6: Polynomial mutation and binary crossover are applied to reproduce the population p g ' and create an offspring population o g . g=g+1

Step 7: Objective functions are evaluated o g.

Step 8: o g is combined with the current generation population p g to form a union population u g = p g o g. Chromosomes in u g are sorted based on non-domination and crowding distance assigned.

Step 9: Select the N best solutions from u g to create the next parent population p g+1

Return p g+1 and the best front as the optimal set of solutions. Table 1. Allocation of the different power sources on the generation nodes of the distribution network of Fig. 5.

  , where is the deterministic limit (DL) of the voltage, is the reference voltage and is the voltage magnitude in per-unit (p.u.) in the node or bus b. In this study, we set p.u. and p.u., following [27]. Fig 1 illustrates eq. (4), where the deterministic violation region (DV) contains the values satisfying and the near violation region (NV) contains the values satisfying .

Fig. 1

 1 Fig.1

  , where ambient temperature [ºC ], is the nominal cell operating temperature [ºC ], cell temperature [ºC ], short circuit current [A], current temperature coefficient [mA/º C], open circuit voltage [V], voltage temperature coefficient [mV/º C], voltage at maximum power [V], current at maximum power [A]. FF fill factor, number of photovoltaic cells, PV power output [W].

  power size of generator type k (e.g. wind turbine, EV, etc) on node b. Given , the simulation model is run under various scenarios to evaluate the objective functions, which are the expected risk , obtained as the mean of the risk values over all MC runs, and the expected global cost obtained by averaging the values of eq. (28) over all MC runs. The global cost sums up the operation and maintenance costs and the investment cost related to the specific DG size, similarly to what proposed by Mena and co-authors [20]: of the scenario [h], is the horizon of analysis [h], are the operation and maintenance costs of the total power supply and generation [$], is the investment cost of the k-th generator technology [$/kWh], are the operation and maintenance costs of the power source (with = , , , ) [$/kWh], is the global investment cost of all the technology installed on node b [$] and is the global cost [$]; 4.4 Solving procedure for the MOO For the MOO solution, we adopt the Non-Dominated Sorted Genetic Algorithm version two (NSGA-II), because it has been shown to be very effective through the incorporation of elitism in the search and with no need for the sharing parameter to be chosen a priori [35]. Briefly the NSGA-II applied to a MOO problem proceeds as follows: at step 1, we set the generation number gen, starting population size N and MC runs number; at step 2, the starting population p 0 is randomly created, each one of the N chromosomes inside p 0 is associated a decision vector and each element is constrained between preselected lower and upper power limits (in fact, the constraints are incorporated into the global cost objective function. The unfeasible solutions are subject to the death penalty, i.e. their fitness values are set to positive infinite.); the 3rd step consists of the MC simulation for evaluating the objectives for each chromosome in p 0 ; at step 4, the population is sorted and ranked based on non-domination and crowding distance; the 5 th step starts the evolutionary loop for the search of optimal (dominant) solutions: the generation number g in the first loop is set to zero, the parent chromosomes are selected and an intermediate population p g ' is obtained; at step 6, the polynomial mutation and binary crossover (having probabilities pm and pc, respectively) are applied to the intermediate population p g ' and the result is the offspring population o g ; the 7 th step consists of the MC simulations used to evaluate the objectives for each one of the vectors inside o g ; at step 8, o g is combined with p g to form the union population u g = p g o g , whose chromosomes are, then, sorted based on non-domination and crowding distance; at step 9, the first N elements of the sorted u g are selected to create the new population p g+1 ; finally, the generation number is updated g=g+1 and the evolutionary process is repeated until g is equal to gen runs. Fig 4 sums up the procedure to solve the MOO by NSGA-II.

  Fig 4

Fig. 9 Figs. 10

 910 Fig.9(a)

Fig. 10 Fig. 11

 1011 Fig.10 Fig.11 It can be observed in Fig. 10 that a large number of violation scenarios occur. These are obviously

Fig. 12

 12 Fig. 12(a) Fig. 12(b) Fig. 12(c)

Fig. 13

 13 Fig.13(a)

Fig. 14

 14 Fig.14

Fig. 15

 15 Fig.15

Fig. 16 6. 4

 164 Fig.16

Fig. 19 - 22 Fig. 23

 192223 Fig.19-22

Fig. 5

 5 Fig.5 Modification from IEEE 13 nodes radial distribution network to the 11 nodes DG system.

Fig. 6

 6 Fig.6Hourly per day probability data of EV operating states[START_REF] Mena | Optimal sizing and allocation of distributed generation for reliable energy distribution accounting for uncertainty[END_REF]

Fig. 7

 7 Fig.7 Expected risk index as a function of MC iterations, in the DG system

Fig. 12 (

 12 Fig.10 and visualization. Fig.11 and visualization. Dots= severity values; Dots= severity values; Solid line=set of results in DV,NV or secure zones Solid line=sets of results in DV zone; Dashed line= set of results under normal weather

Figs. 19 -

 19 Figs.19-22 Optimal power size in each node and for each technology

  3 Flow chart of a single MC run

	2900
	2700
	2500
	2300
	2100
	1900
	1700
	1500
	0	2	4	6	8	10	12	14	16	18	20	22	24

Table 2 .

 2 Parameters of DG models, weather models and failure rates

		DG model parameters
	WT	PV	ST	EV
		Voc=55.5 [V]
	[m/s]	Ta=30 [°C]
	[m/s]
		Weather model parameters
	Wind speed	Solar irradiance	Lightning	High wind
	2.3			Vcritic = 9.86[m/s]

Table 3 .

 3 Distribution line parameters

	1,2
	1									1	Discharging
						0,9				0,96	Charging
	0,8				0,7	0,7	0,73	0,69	0,74	0,77	0,81	Disconnected
	0,6					0,59				0,59
		0,45 0,47	0,45 0,43 0,44 0,44						0,51	0,42	0,43 0,44 0,45
	0,2 0,4	0,36	0,19 0,2 0,19 0,18 0,19 0,19 0,17 0,4 0,21 0,07 0,33 0,36 0,39 0,37 0,37 0,43	0,12 0,09 0,08 0,09 0,08 0,29 0,18 0,16 0,21 0,19 0,22	0,03 0,13	0,12 0,29 0,34 0,15 0,17 0,19 0,19 0,19 0,38 0,37 0,36 0,41
	0				0,09	0,03			0,07 0,06	0,01	0	0
		1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
				node b	node b'	R[ /km] X[ /km] Ampacity[A] length[km]
				1	2	0.116	0.371	365	0.610
				2	3	0.368	0.472	170	0.152
				2	4	0.696	0.555	115	0.152
				2	6	0.116	0.371	365	0.610
				4	5	0.696	0.555	115	0.091
				6	7	0.303	0.252	165	0.152
				6	8	0.696	0.555	115	0.091
				6	11	0.116	0.371	365	0.305
				8	9	0.696	0.555	115	0.091
				8	10	0.377	0.318	115	0.244

Table . 4

 . Nodal power size limit for the technology, equal in each node.

	k	[kW]	[$/kWh]	[$/kWh]
	PV	10	0.41	0.0000376
	W	100	0.026	0.0390000
	EV	33.65	0.028	0.0220000
	ST	100	0.0771	0.0000462

Table 5

 5 Nominal power in each node for the three solutions.

		[kW]	[kW]	[kW]
	Node1	49.24	84.93	177.61
	Node2	69.68	93.24	161.72
	Node3	143.44	55.63	122.49
	Node4	99.56	86.24	158.81
	Node5	109.35	66.70	129.99
	Node6	28.64	175.99	197.77
	Node7	48.85	114.73	134.07

Table 6

 6 Risk index results for the different systems analyzed

	Network			E(Cg) [$]
	configuration	Considering extreme events	No extreme events
	Radial	0.02284	0.002237	0
	Allocation Table 1	0.02256	0.002215	~31
		0.02244	0.002204	~50
		0.02238	0.0022	~64
		0.02227	0.002193	~94

Table 7

 7 Risk reduction with respect to the radial case (in percentage)

	Network	Percentage of risk reduction	Percentage of risk reduction
	configuration	with respect to the radial case	with respect to the radial
		(considering extreme events )	case (neglecting extreme
			events)
	Allocation Table 1	-1.24	-1.01
		-1.78	-1.49
		-2.06	-1.73
		-2.56	-2.00

Table 8

 8 Results of the Wilcoxon rank sum tests in comparing risk distributions

	Paired comparisons	p-value
	vs. radial	<0.001
	vs. radial	<0.001
	vs. radial	<0.001