On the Sensitivity of Local Flexibility Markets to Forecast Error: A Bi-Level Optimization Approach

: The large-scale integration of intermittent distributed energy resources has led to increased uncertainty in the planning and operation of distribution networks. The optimal ﬂexibility dispatch is a recently introduced, power ﬂow-based method that a distribution system operator can use to effectively determine the amount of ﬂexibility it needs to procure from the controllable resources available on the demand side. However, the drawback of this method is that the optimal ﬂexibility dispatch is inexact due to the relaxation error inherent in the second-order cone formulation. In this paper we propose a novel bi-level optimization problem, where the upper level problem seeks to minimize the relaxation error and the lower level solves the earlier introduced convex second-order cone optimal ﬂexibility dispatch (SOC-OFD) problem. To make the problem tractable, we introduce an innovative reformulation to recast the bi-level problem as a non-linear, single level optimization problem which results in no loss of accuracy. We subsequently investigate the sensitivity of the optimal ﬂexibility schedules and the locational ﬂexibility prices with respect to uncertainty in load forecast and ﬂexibility ranges of the demand response providers which are input parameters to the problem. The sensitivity analysis is performed based on the perturbed Karush–Kuhn–Tucker (KKT) conditions. We investigate the feasibility and scalability of the proposed method in three case studies of standardized 9-bus, 30-bus, and 300-bus test systems. Simulation results in terms of local ﬂexibility prices are interpreted in economic terms and show the effectiveness of the proposed approach.


Background and Motivation
The large-scale integration of intermittent renewable energy sources (RES) and the change in the role of end users, from energy consumers to prosumers with new technologies (e.g., distributed energy resources (DER) and storage facilities), have increased the uncertainty in the operation of the power systems at all levels. One way to overcome the uncertainty is to increase flexibility in the distribution system by harnessing the flexibility that is available to the end users by implementing demand response (DR) programs.
In general, the goal of a DR program is to motivate end users to adapt changes in their local production/consumption profiles in response to a price (i.e., through a market-based mechanism) or a command (i.e., originated by a non-market-based mechanism) signal [1,2]. Regardless of the mechanism considered in place, the distribution system operator (DSO) needs to become active in network management by taking a coordinating role through implementing a proper policy that provides an adequate price or command signal. This for one, would enable the DSO to harness the flexibility from DR providers and overcome the challenges related to uncertainty.
To do so, the DSO needs to determine the amount of flexibility that is needed to maintain the secure and reliable operation of the network. The DSO also needs to determine a reliable projection (i.e., prediction) of the amount of flexibility that will potentially be available at its disposal, at different times and locations during the real-time operation. Once this information is available, the DSO can implement the mechanism to determine (i.e., estimate or calculate) a proper price or command signal such that it encourages DR providers to engage in network management by offering their flexibility in an efficient manner.

Literature Review
This subsection briefly reviews relevant academic literature.

Harnessing Demand Side Flexibility in the Distribution Network
Harnessing the flexibility potential of the consumers has been investigated in the literature. In [3][4][5][6] authors investigate mechanisms that can provide a price or command signal to DR providers to harness their flexibility and maintain power balance in power system, neglecting the network constraints. Subsequently, references [7][8][9] introduce new algorithms that enable an aggregator to estimate the flexibility that is potentially available to DR providers under certain assumptions. These works provide valuable insight into the existing network management potential through harnessing end users' flexibility. What is missing in these studies is to account for the impact of including the topology and operational constraints of the distribution grid (e.g., voltage limits). The grid constraints play a crucial role in determining the accurate amount of flexibility that is required to maintain a secure and reliable operation of the network.
A common method to account for network constraints in the analysis of the operation of an electricity network (including flexibility analysis) is to use optimal power flow (OPF)-based models [10][11][12][13][14][15][16]. In this regard, Capitanescu et al. [13] proposed an OPF model that seeks to minimize the curtailment of the DERs while maintaining voltage constraints in a distribution grid. While the existing works are mostly concerned with the technical aspects of the grid, to ensure a successful network management, the economic aspects also need to be considered properly [17].

Pricing Flexibility in Distribution Network
References [18][19][20] introduce OPF-based models to analyze the relation of power flows, network congestion and voltage constraints on the formation of distribution locational marginal prices (DLMP) in local energy markets. However, a limitation of AC-OPF is that the non-linearities of AC power flows introduce computational complexities. Li et al. [18] use DC-OPF which is a linear approximation of the full AC-OPF to investigate the day-ahead DLMPs. Likewise, Yuan et al. [19] propose a new model for pricing the electric power at the distribution level, considering voltage constraints, using linearized AC power flows. Papavasiliou [20] uses second-order cone (SOC) relaxation of AC power flows and proposes a novel framework to study the formation of DLMPs in distribution network, considering congestion, voltage constrains, real and reactive power limits and losses, as well as the non-linearities of the AC power flow in the distribution grid.
In recent years in parallel with higher penetration of DERs, with their variable and unpredictable production, there has been a growing interest in flexibility (defined as in [21]) as the service to be traded, instead of absolute energy, in local energy communities [22][23][24][25][26]. To this end, authors in [27] introduce the optimal flexibility dispatch (OFD) as an optimization problem. The OFD is an auction-based platform for trading flexibility in distribution grids. The OFD seeks to minimize the costs of curtailing RES and, as a side product, determines the quantity and price of active and reactive power flexibility that the DSO needs to procure from different flexibility providers at different locations in the network.

Second-Order Cone Relaxation Error
The OFD framework proposed in [27] considers network constraints, end users' flexibility limits as well as the nonlinearity of AC power flow equations. The problem is that considering such non-linearities in the OFD makes the problem non-convex and non-tractable for large-scale systems. To make the AC-OFD problem tractable, authors in [27] use the second-order cone (SOC) relaxation of the AC power flow formulation. However, the exactness of the solution of the SOC-relaxed problem is an issue of concern, as it can only be guaranteed under certain simplifying assumptions [28,29]. It should be noted that inexact solutions do not coincide with Kirchhoff's laws and therefore, lack a physical interpretation.
SOC-relaxation error minimization can be done in various ways. One way would be to modify the problem formulation by applying the penalty function method, which is to include the relaxation error as an auxiliary term in the objective function of the SOC-OFD problem and penalize it by a positive weighting parameter [30,31]. The problem is that such modification yields solutions (e.g., the flexibility prices, power flows) which are different from the original SOC-OFD problem. In addition, one would see that the results of the modified SOC-OFD problem are sensitive to the weighting parameter. This itself would raise questions regarding appropriate values for the weighting parameter [30][31][32][33].

Research Gap
In summary, most works on harnessing flexibility from the distribution network first consider energy instead of flexibility as the commodity to be traded. Secondly, the existing works either follow a purely economic approach neglecting technical constraints, including power flow constraints, or vice versa. In addition, thirdly, even if OPF-based models are used (through which technical constraints are accounted for), they are computationally costly, due to non-linearities of the AC power flow, for large-scale problems, or yield solutions for which the feasibility of the results cannot be guaranteed due to errors inherent in the approximation and/or relaxation techniques applied. Fourthly, none of the studies above account for the impact of inaccurate input data, e.g., forecast errors, on the optimal energy or flexibility programs.
Please note that the accuracy of the input data is an important factor because, in real world applications, the amount of flexibility that the DSO can have at its disposal in real time is uncertain and is not fully known beforehand. As a result, estimates of the available flexibility values shall be used in practice. This can affect the operational decisions. Therefore, it is important to study the impact of uncertainty in the potential flexibility that is available to the DSO on the flexibility schedules and flexibility prices.

Contributions
This paper is follow-up of the work by Torbaghan et. al. in [27], where a SOC relaxation of the AC power flow is used to solve the OFD problem. In this work we propose a bi-level optimization problem formulation to minimize the relaxation error in the SOC-OFD problem in order to address the above-mentioned shortcomings. The major contributions of this paper are as follows: • proposes a new formulation for the SOC-OFD problem that seeks to find an optimal flexibility schedule (i.e., amount and price of flexibility to be dispatched at each bus) with zero relaxation error. This is done by formulating a bi-level optimization model that seeks to solve the SOC-OFD problem while minimizing the SOC-relaxation error, • introduces an effective reduction/relaxation technique that makes the bi-level problem tractable both analytically and numerically, and • provides an analytical framework for performing local sensitivity analysis on the reduced OFD problem, which determines the robustness of the resulted optimal flexibility dispatch against imperfect forecasts of input data.
The bi-level OFD framework introduced here has two key characteristics. First, it results in a solution with zero relaxation error. Secondly, it returns 'locational flexibility prices' (LFPs) as explicit optimization variables that are otherwise only available as the dual of the power balance constraint at every bus, and thus only available through specific solvers.

Outline
The remainder of this paper is organized as follows. Section 2 first provides a brief overview of the SOC-OFD problem. It then dives into the mathematical formulation of the bi-level optimization model that seeks to minimize the relaxation error while solving the SOC-OFD problem. Section 3 presents the numerical results for the case studies. Sections 4 and 5 conclude the paper and discuss directions for future research.

Problem Formulation
The OFD is defined as an optimization problem that seeks to minimize the curtailment of DERs while keeping the distribution grid congestion free. Inspired by [27] we use the SOC formulation of the AC power flow as a practical approach to model the network constraints while accounting for the non-linearities of AC power flow. To keep the SOC-relaxation error as small as possible, we propose a bi-level structure in which the upper level (UL) problem seeks to minimize the SOC-relaxation error, while the lower-level (LL) problem solves the SOC-OFD problem. As the lower-level problem is convex, using Karush-Kuhn-Tucker (KKT) conditions, we recast the bi-level structure into a mathematical problem with complementarity constraints (MPCC) [34].
An MPCC problem is in general very difficult to solve. Thus, to make the problem tractable, we apply an innovative relaxation method to relax the complementarity constraints and reformulate the problem as a non-linear problem (NLP), without losing accuracy. In this sense, the proposed methodology establishes a good trade-off between accuracy and tractability. Solving the proposed bi-level optimization problem enables the DSO to determine the quantity of the active and reactive power flexibility, as well as the distribution locational marginal prices of flexibility, with minimum error. Please note that the convexity of SOC-OFD allows one to apply local sensitivity analysis to determine the sensitivity of the optimal flexibility dispatch program with respect to perturbation in input data (e.g., available flexibility) which mimics imperfect forecasts [35].

Assumptions
Consider a distribution grid with N buses. We use indices i or j and i, j ∈ Ω N to refer to buses where Ω N is the set of buses. Assume that every bus can comprise of a number of households with flexible and non-flexible units. A flexible unit can be a generator (e.g., DERs) or a load (e.g., DR). For the sake of simplicity, we assume that all generators are flexible. Loads can be flexible or non-flexible. Suppose that each bus is operated by an aggregator. Each aggregator receives the "base-power program". It entails supply and demand forecast of the flexible and non-flexible units. The aggregator also receives the "available flexibility" forecast of the flexible units. The aggregator processes this information and provides them to the DSO. For every bus, the aggregator provides the base active and reactive program (i.e., p 0 u , q 0 u ) where p and q denote respectively active and reactive power of unit u ∈ Ω k i in bus i. Index k determines whether unit u in bus i is a generator (i.e., {k = G|u ∈ Ω G i }) or load (i.e., {k = D|u ∈ Ω D i }). Ω G i and Ω D i denote respectively the set of generation or consumption units in bus i. The aggregator also provides the upper and lower bounds of active and reactive flexibility of the flexible units (∆p u , ∆p u , ∆q u , ∆q u ) where ∆p u and ∆q u denote respectively the active and reactive flexibility of unit u ∈ Ω k i . Please note that x, x denote respectively the lower and upper bound values of an arbitrary variable x. The operational limits include power flow limits of the lines and bus voltage limits. The power flow limits are denoted as p i,j , p i,j , q i,j , q i,j , where p i,j and q i,j denote the active and reactive power flow of the line (i, j) ∈ Ω E connecting bus i to bus j. Ω E is the set of power flows when power flows from the sending bus (i.e., bus i) to the receiving end (i.e., bus j). Please note that Ω E r will be used to refer to the set of power flows in the reverse direction (i.e., for line i, j, and the convention of bus i being the send bus and j being the receiving bus, we have . In this formulation |v i | and |v i | are the upper and lower bound of the bus voltage, and θ i,j is the phase angle difference between buses i and j.

Second-Order Cone Relaxed OFD (SOC-OFD)
The SOC-OFD problem, as proposed in [27], is formulated as follows: ∀i ∈ Ω N where τ DER is the fixed tariff of DER curtailment in A C/MWh. g sh i and b sh i are the shunt conductance and admittance of bus i. Likewise, g s i,j and b s i,j are series conductance and admittance of the line (i, j).
} is the set of independent optimization variables and Γ LL is the set of Lagrange dual variables. The Greek letters between parenthesis denote Lagrangian dual variables for the corresponding constraint.
The objective function (1) states that the DSO seeks to minimize the cost of curtailing the production of DERs. This is a weak assumption. One could consider other objectives as long as the objective function is convex. Constraints (2)-(4) enforce lower and upper bounds of the square voltage and voltage products. Constraint (5) reflects the SOC-relaxation of the AC power flow. Constraints (6) and (7) (17) and (18) are the Kirchhoff's current law (KCL) which enforce real and reactive power balance.
The SOC-OFD problem can be inexact due to relaxation error. Reference [27] showed that the relaxation error can be larger when flexibility is scarce. The non-exact solution lacks physical interpretation and at best, can be treated as an approximation of the exact solution.
In the following subsection we propose a novel approach to reduce the relaxation error.

Bi-Level Model for the OFD Problem
One can define the SOC-relaxation error associated with (5) as follows: If SOC i,j = 0, that is if the inequality (5) becomes an equality, then the relaxation error is zero and the solution of the SOC-OFD coincides with of the full-AC-based OFD. However, if the error is non-zero, then the results of the OFD problem are inexact and have no physical interpretation. To minimize the SOC-relaxation error, we introduce a novel alternative approach. The idea is to consider the SOC-relaxation error minimization as a separated optimization problem that should be solved simultaneously with the SOC-OFD problem. This can be done by formulating a bi-level optimization model in which the UL problem seeks to minimize the SOC-relaxation error as defined in (19). The LL problem is the original SOC-OFD problem (i.e., problem (1)-(18)). The bi-level OFD problem takes on the following form: subject to where (1)-(18) is the original SOC-OFD problem. The objective (20) is to minimize the aggregated SOC-relaxation error. The SOC-relaxation error minimizing problem (i.e., UL problem) is constrained by the SOC-OFD problem (i.e., LL problem). The decision variables of the UL problem includes the decision variables of the LL problem (i.e., Ξ LL ) in addition to the dual variables associated with constraints of the LL problem (i.e., Γ LL ).
Since the constraining LL problem is convex, it can be replaced by its corresponding KKT conditions. This results in a single level optimization problem [34]. which is a mathematical program with complementarity constraints (MPCC) [34]. The resulting MPCC is a non-linear and non-convex problem, and therefore is NP-hard.
There are two sources of non-linearities in the resulting MPCC. The first source is the quadratic terms, i.e., (w i i,j ) 2 and (w r i,j ) 2 , and the bi-linear products, i.e., w i w j , in the objective function (20) with SOC-relaxation error given by (19) as well as in the SOC-relaxed power flow constraint in the LL problem (5). The second source is the complementarity constraints in the KKT optimality conditions of the LL problem. The complementarity constraints are of the form 0 ≤ a ⊥ b ≥ 0. Each complementarity constraint can be written equivalently as a set of three constraints, which are 0 ≤ a, 0 ≤ b, and a · b = 0, where the last constraint is non-linear and extremely non-convex. The MPCC problem must be simplified to become tractable.
A possible simplification is to relax equality a · b = 0 with an inequality constraint as a · b ≤ δ, where δ is a non-negative variable to be minimized, as discussed in [36]. Such relaxation transforms the MPCC into an NLP problem, which can be solved using available solvers especially for small problems. For large-scale problems, however, the proposed relaxation can be problematic [37].
To reduce the size and complexity of the problem, and also to improve the accuracy of the results, in the following we propose a new relaxation technique which is to first reconstruct the bi-level formulation of the problem, then transform the resulted MPCC to a NLP problem using the relaxation discussed above.
Now consider shifting constraints (23) and (25) to the UL and reformulating the problem as follows: subject to subject to As suggested in [38], one can show that problem (21)-(26) is equivalent to problem (27)-(32), when the dual variables l 1 , l 2 = 0.
In a similar way, we reformulate the bi-level OFD problem by shifting all constraints from the LL to the UL except for the active power flexibility limit (6) and active power balance (17). The proposed reformulation reduces the LL problem to a new optimization problem that seeks to minimize curtailment of DERs subject to determine the active power flexibility only. As the constraints of the UL problem are transparent to the LL problem, the UL induces a proper choice of voltages, active and reactive power flows, and reactive power flexibility in a way that a minimum of the sum of relaxation errors of all SOC constraints is reached. It is apparent that if the shifted constraints are not binding (i.e., their associated Lagrange multipliers are zero) the reduced bi-level OFD problem is equivalent to the original bi-level problem. As a result, an optimal solution of the reduced problem can be matched to a feasible solution of the original bi-level problem if the constraints that are shifted from the LL to UL are not active at the optimal solution. This point is further investigated in the next section.
Once the LL problem is reduced by shifting the constraints to the UL problem, the reduced LL problem can be replaced by its corresponding KKT conditions. The resulted complementarity conditions in the reduced MPCC can be further relaxed using the technique discussed in Section 2.3, as follows: subject to (1)- (16), (18), where is the set of decision variables. Constraints (1)- (16), and (18) include the objective function of the SOC-OFD problem, the square voltage and voltage product limits, the SOC-relaxation of AC power flows, reactive power flexibility limit, and reactive power balance, which are shifted to the UL problem. Constraints (34) and (35) enforce a threshold on complementarity relaxation, where parameter δ max is the maximum allowed constraint relaxation. Constraints (36)-(44) are the KKT conditions associated with the reduced LL problem (i.e., reduced SOC-OFD). Equalities (36) and (37) represent the derivatives of the Lagrangian of the LL problem with respect to its decision variable ∆p k u (i.e., ∂L LL /∂∆p k u = 0), for DER and load units, respectively.
Problems (1)-(16), (18), (33)-(44) (referred to in this paper as bi-level-OFD) is a NLP which can be solved in a feasible time by available solvers. The proposed reformulation effectively reduces the computation complexity of the problem with the lowest loss of accuracy possible.

Derivation of the Sensitivity Matrix in the OFD Problem
In this paper, we also aim to investigate the sensitivity of the optimal solution of the OFD problem to perturbations in input data. One can study the sensitivity by deriving a linear sensitivity matrix based on the perturbed KKT equation proposed in [35]. We induce a perturbation in the KKT conditions and the objective function of problem (1)-(18) (i.e., Ψ LL ), with respect to the decision variables Ξ LL , the Lagrange multipliers, and the input data parameters Then, as proposed in [35] we obtain the sensitivity matrix M = dX/da, which determines the sensitivities of the optimal solution along with the Langrange multipliers to changes in data. Suppose that inactive constraints are removed from the analysis, and that the active inequality constraints remain active after the perturbation. As the LL problem represents the (reduced) original SOC-OFD problem, one can see that an optimal solution of the reduced-bi-level problem includes the optimal solution of the SOC-OFD. Therefore, at the optimal solution, the sensitivities of the SOC-OFD problem can be evaluated at the optimal solution of the reduced bi-level-OFD problem.

Numerical Results
In this section, we present the numerical results obtained using the SOC-OFD, i.e., problem (1)-(18), and the bi-level-OFD, i.e., problem (33)- (39), for three case studies. We compare the results in terms of (i) the SOC-relaxation error [ SOC , as defined in Equation (19)], (ii) the DER curtailment [Ψ LL as defined in Equation (1)], and (iii) locational flexibility prices (λ p i ) which are derived as dual variables of active power balance constraint. Furthermore, we numerically analyze the sensitivity of the OFD problem to perturbations in input parameters (i.e., ∆p u , ∆p u , p 0 u ).

Case Study
In power system research it is common practice to use synthetic electric grid models as case studies. The reason is, they provide a standard benchmark that enable researchers to validate newly developed analytic methods against the existing ones, without using any confidential information of a critical energy infrastructure [39].
In this study, to investigate the feasibility and scalability of the proposed framework, three power grids from MATPOWER, with 9-bus, 30-bus, and 300-bus, are considered. The grid data and topology of the three case studies are adopted from a subset of automated tests prepared for MATPOWER [40][41][42].
MATPOWER is an open-source MATLAB-based power system simulation package. It provides a set of (optimal) power flow-based tools that are specifically designed for research purposes [43]. The network parameters are customized to make the technical characteristics of the test networks consistent with those of distribution networks. Table 1 provides an overview of the grid parameters, including the number of buses with DER resource, number of branches, average per unit (p.u.) values of the series resistance (R), inductance (X), and shunt susceptance (B) for each case study.  Table 2 provides an overview of the three case studies used in this work. The aggregated supply and demand represent the total generation and total demand of the DER and flexible load buses in p.u. values, respectively. The aggregated available flexibility of the DER and the load buses represents the total available active power that the DERs and the flexible load units can increase/decrease. The selected cases allow for investigating the scalability of the proposed framework.

1.
Each bus represents a cluster of aggregated end users which can be either a net flexible generator or net flexible load.

2.
The OFD problem is solved per single time-step. Inter-temporal linkages are neglected.

3.
A fixed cost for curtailing DER production is considered, and τ DER cur = 0.05 A C /kWh.

4.
Cost of load curtailment is assumed to be zero.

5.
Load is modelled as negative generation, i.e., negative flexibility for a net flexible load bus indicates an increase in the total consumption of that bus.
Please note that although the results provided in the context of this paper correspond to a specific objective for the DSO (i.e., minimizing cost of curtailment), the proposed model and analytical discussion is entirely general and can be applied when considering other objectives for the DSO as long as the objective function remains linear.

Simulation Results
The SOC-OFD problem and the bi-level-OFD problems are implemented in MATLAB 2017b, and solved by optimization toolbox YALMIP [44]. The computer used ran Windows 7 Enterprise 64-bit, with an Intel Xeon CPU E5-1650 processors running at 3.20 GHz and 24 GB memory. The bi-level-OFD problem is solved for δ max = 1 × 10 −6 . )   Table 3 presents the minimum, maximum, and average SOC-relaxation error (i.e., SOC = ∑ i,j∈Ω E SOC ij ) for each case study, calculated by solving the SOC-OFD and bi-level-OFD models. One can observe that for Case 1 and Case 2, SOC is minimized to zero when the bi-level-OFD model is used. For Case 3, this figure is reduced by an average of 89% when compared with the results of SOC-OFD model. This shows that the bi-level-OFD model effectively eliminates the SOC-relaxation error for the 9-Bus and 30-Bus test systems, and reduces it considerably for the 300-Bus system. It also suggests that solving the bi-level-OFD model indeed results in solutions that have physical interpretation, as outlined in the previous section. Table 3. The maximum, minimum, and average SOC-relaxation error SOC of all line flows, determined by SOC-OFD and bi-level-OFD models for the three cases. Now consider Figure 1 which depicts the distribution of bus voltages calculated by both models for each case. The horizontal red lines inside the boxes indicate the median of the voltages in the network. One can observe that under each case the bus voltages of the majority of buses are lower when calculated by the bi-level-OFD model as compared with the SOC-OFD model.  Figure 2 presents the aggregated production of DERs calculated by the SOC-OFD model (blue bar) and the bi-level-OFD model (red bar). One can see that the bi-level-OFD model yields the same aggregated DER production as that of the SOC-OFD model. This implies that the bi-level model finds a solution in which the flexible units are dispatched in a such way that the aggregated DER production (and thus DER curtailment) remains the same as for the SOC-model, but with no SOC-relaxation error.   One can see that the LFP of every DER bus (i.e., the red dots in Figure 3) is a positive value which is lower than or equal to the DER curtailment tariff set by the regulators (i.e., τ DER = 5 A Ccent/kWh). Likewise, one can see in Figure 4 that for the three cases, the flexibility dispatched from the DER buses is positive. This implies that the DSO pays the DERs to increase their production (i.e., decrease curtailment). Now consider the flexible load buses (i.e., blue dots in Figure 3). For Case 1 and Case 2 the LFPs at load buses are very close to zero. In Case 3, however, there are certain load buses with a negative LFP. From (37), (41)-(42) one can observe that for load buses, λ p i gets a non-zero value only when the boundaries of active flexibility are reached.

DER Curtailment (Ψ)
From Figure 4, one can also observe that in Case 1, the bi-level model yields negative flexibility values for load buses. This implies that the model steers the load buses to increase their load to match DER production and reduce DER curtailment. This contrasts with Case 2 and Case 3, in which some load buses get a positive flexibility value. A positive value means that the model sets the price such that it enables the DSO to steer the load buses to reduce their load. This observation can also be explained through curtailment. In Case 2 and Case 3, DERs are producing at their maximum capacity (the zero curtailment levels in Figure 2), therefore, as the total production has reached its maximum limit, the model tends to decrease the load in certain buses such that the total demand equals the total (maximum) supply.

Sensitivity of LFP to Data Perturbation
For Case 1 and Case 2, Figures 5 and 6 respectively present the aggregated sensitivity of LFPs of all buses (j), to perturbation in the lower-or upper bound of active flexibility ∆p i , ∆p i , as well as the base-power program p 0 i of a certain bus i, that is:

Discussion
The results of the bus voltages calculated by the SOC-OFD and bi-level-OFD show that in general, the bus voltages are lower when calculated by the bi-level-OFD model. This observation can be explained via the relaxation error. As outlined in [27], the relaxation error can be interpreted as an artificially increased active and reactive loss over the lines. Now consider the results of the SOC-OFD model (in which no boundary is superimposed on increasing the relaxation error). One can observe that in the SOC-OFD model, the solver finds solutions with a higher relaxation error, and therefore higher (artificial) losses. This leads to higher voltage differences as well as higher network voltage levels, especially at the DER buses. The opposite is true for the bi-level-OFD model, in which the relaxation error is minimized, and as a result, the bus voltages obtain lower values.
The results of the flexibility prices show that if flexibility is dispatched such that consumers' flexibility boundaries are binding, the LFP of the respective bus equals the penalty tariff. This implies that the flexible loads have incentive to provide flexibility services such that they reach their full flexibility limits. Please note that other settings could also be possible, for instance, those in which flexible loads are also remunerated when deviations remain within stated boundaries, which is beyond the scope of this work.
The local sensitivity analysis and the results provided are dependent on the specific network configuration. However, the methodology is generic. In this light, we argue that the DSO can use the local sensitivity analysis to identify the relative importance of (uncertain) input data (e.g., forecast on flexibility limits and load) on the solution of OFD problem. This includes the optimal flexibility dispatch and the formation of LFPs.
In the context of this work we assumed flexibility providers and DERs are competitive but are not able to exercise market power. However, the absence of an adequate regulatory framework, especially in the short run, in the context of local energy and flexibility trading platforms, may encourage opportunistic participants to engage in strategic behavior. Therefore, regardless of the market setting considered, such behaviors could significantly influence the operation of the local energy market and from there, the decision making process of the DSO. Different mathematical models would be necessary for capturing strategic behavior of flexibility providers.

Conclusions
The OFD is an optimal power flow-based method that enables the DSO to determine the amount of flexibility (and associated prices) that it needs to dispatch from flexible resources to meet its operational constraints. Due to non-linearities in the power flow, the SOC relaxation of the AC power flow is typically used in the OFD problem, i.e., the SOC-OFD model. However, a solution to the SOC-OFD model can be inaccurate and thus will lack physical interpretation. When formalized, the SOC-OFD model takes certain input parameters that neither the DSO nor the flexibility providers are entirely aware of, and therefore, are highly uncertain (e.g., forecasts on available flexibility in the day-ahead).
In this work, we introduced a bi-level optimization problem in which the UL seeks to minimize the SOC-relaxation error and the LL seeks to solve the SOC-OFD. The objective of the LL is set to minimize the curtailment of DERs subject to operational limits. Next, we derived a linear sensitivity matrix based on the perturbed KKT conditions of the SOC-OFD problem. We use the sensitivity matrix to investigate the sensitivity of the optimal flexibility dispatch solution (i.e., LFP and ∆p i ) to forecast errors in input data (i.e., ∆p i , ∆p i , p 0 i ) . To solve the proposed bi-level problem, we first reduced the complexity by shifting certain constraints from the LL to the UL. We showed that such reformulation causes no loss of accuracy as long as the shifted constraints are not binding. Secondly, we recast the reduced bi-level problem as a MPCC by replacing the LL with its corresponding KKT conditions. Finally, we relaxed the complementarity constraint and reformulated the problem as a NLP which is computationally tractable, even for large-scale problems. The proposed framework is general and can be applied regardless of network topology, as long as the objective of the lower-level problem is linear.
We compared the results of our bi-level-OFD model for three test networks with that of the SOC-OFD model. The analysis shows that for the same level of DER production, the bi-level-OFD model can effectively reduce the SOC-relaxation error.
Our results show that the DSO pays the DERs when they provide flexibility by increasing their production. However, flexible loads have no incentive to provide flexibility services unless they reach their operational flexibility boundaries. This is due to the way the OFD problem is formulated in this work. Therefore, one avenue for research in the future is to investigate other settings that provide an incentive by reimbursing flexible loads even before hitting the boundaries of their flexibility range.
A practical implication of our proposed sensitivity analysis framework is that the DSO can identify the forecast error of which bus(es) has the most prominent impact on the optimal flexibility dispatch and the LFPs of the other buses in the network. This would allow the DSO to take adequate preventive actions to restrain unintentional deviations from the optimal values, and therefore, ensure the economic efficiency of the resulting flexibility dispatch in practice.
A direct extension of the current work would be to investigate sensitivity of an optimal flexibility schedule to uncertainty associated with exercising such strategic behavior. Additionally, along with the recent development in peer-to-peer distributed energy trading [45], there is a growing need for determining the hosting capacity of the distribution network for enabling such trades in advance. In this regard, one application of the proposed methodology in practice is for the DSO to estimate the hosting capacity of its network, considering the uncertainties associated with imperfect forecast and/or, the (un)intentional non-compliance of peers from the scheduled trades in the real-time operation.
In conclusion, the methods and resulting insights of this paper can support the DSO in its many activities related to planning and operation of the distribution grid with increasing amounts of (intermittent) distributed resources.

Conflicts of Interest:
The authors declare no conflict of interest.