**3. Problem Formulation**

The objective of the method presented is to provide the structure of the RES and ES allocation in the distribution system along with its operation, ensuring minimal costs of the distribution system development without violating the technical standard of the grid operation. Three types of RES (wind turbine (WT), photovoltaic (PV) and biogas (BG)) can be connected in each node with rated power depending on the value of typical units. RES units are represented by representative generation profiles, reflecting typical operating conditions. Due to the grid, overloads can be identified and eliminated by proper RES allocation and curtailment and energy storing. The optimization is performed, preserving technical constraints. Technical constraints refer to the network transmission capacity, nodal voltage standards and power exchange with the transmission system. Power flow in each line is calculated for each time step. The generation for RES is calculated as a product of the installed capacity and

generation profile and load consumption is calculated similarly. Power exchange with the transmission system depends on the power of each transformer which connects these two systems.

The optimization model, which is described in the next section, consists of general mathematical formulations for mixed-integer linear problems which can be used by other optimization software.

#### *3.1. Optimization Model*

The minimization of the total costs of the distribution system development in relation to one year (installation and operating RES and energy storages) is the objective function. Total costs consist of three elements:

Fixed costs (CAPEX and OPEXfix) of new RES,

Fixed costs (CAPEX and OPEXfix) of new energy storages,

Operating costs (OPEXvar).

CAPEX there are investment costs connected with equipment and construction works costs. OPEXfix there are costs intended for taxes, maintenances and salaries. OPEXvar are the costs of fuels.

$$\min \left\{ \sum\_{n \in \mathcal{N}} \left( \sum\_{d \in D} \left( \sum\_{r \in \mathbb{R}} (p\_{n,d,r} \* P\_{d,r} \* \mathcal{F} \mathcal{C}\_{d,r}^{\mathrm{RES}}) \right) + p\_n^{\mathrm{ES}} \* \mathcal{C}^{\mathrm{ES}} \* \mathcal{F} \mathcal{C}^{\mathrm{ES}} + \sum\_{d \in D} \left( \mathcal{E}\_{n,d} \* V \mathcal{C}\_{d,r}^{\mathrm{RES}} \right) \right) \right\} \tag{1}$$

subject to

$$
\begin{array}{ccccc}
\vee & \vee & \vee \\
 n \in N & d \in D & r \in R
\end{array}
\quad p\_{n,d,r} \in N \tag{2}
$$

The DG generation in each node for each type of RES (*En,d*) is calculated as a sum of energy generation of each DG type in one node for the entire year.

$$\forall n \in N \qquad d \in D \qquad \qquad E\_{n,d} = \sum\_{t \in T} E\_{n,d,t}^{GEN\\_RES} \tag{3}$$

The DG generation in each time period (*En,d,tGEN-RES*) is a function of the optimization variable *pn,d,r* and is calculated a product of the RES rated power (*Pd,r*), current power utilization level (according to the generation profile *Pd,ravail*), optimization variable *pn,d,r* and minus energy curtailment (*en,d,tcurt*). The *pn,d,r* is a three-dimension variable, where dimension *n* refers to the number of the node, dimension *d* refers to the RES technology type (i.e., PV) and dimension *r* refers to the rated power (2). Energy curtailment (*en,d,tcurt*) is a variable which curtails the generation from the reference generation profiles for each node (*n*), RES type (*d*) and each time period (*t*).

$$\forall n \in N \quad d \in D \quad t \in T \quad E\_{n,d,t}^{GEN\\_RES} = \sum\_{r \in R} \left( \left( P\_{d,r} \cdot p\_{n,d,r} \right) \cdot P\_{d,t}^{avail} \right) - \epsilon\_{n,d,t}^{curt} \tag{4}$$

Another important parameter is the generation in the entire node (*En,tGEN\_node)* which depends on the sum of DG generation for each RES type in single node and energy exchange with energy storages that are installed in the node.

$$\begin{array}{cccc} \mathsf{V} & \mathsf{V} & \mathsf{V} \\ m \in \mathsf{N} & t \in T & \end{array} \\ E\_{n,t}^{\text{GEN\\_mode}} = \left( \sum\_{d \in D} \left( E\_{n,d,t}^{\text{GEN\\_RES}} \right) - e\_{n,t}^{to \text{ES}} + e\_{n,t}^{fron \text{ES}} \right) \tag{5}$$

The typical RES sizes of each renewable technology are predefined and are included in the two-dimension matrix *Pd,r* (4). As a result of multiplication of variable matrix *pn,d,r* and parameters matrix *Pd,r* the installed power of all RES in every node is calculated.

$$\begin{array}{ccccc} & & & \text{dimension } D\\ \vee & \vee & & & \vee\\ d \in D & r \in R & & & \\ & d \in D & r \in R & & \end{array} \quad P\_{d,r} = \left[ \begin{array}{ccccc} & & & & \\ & P\_{1,d\_1} & P\_{1,d\_2} & P\_{1,d\_3} \\ & P\_{2,d\_1} & P\_{2,d\_2} & P\_{2,d\_3} \\ & \vdots & & \vdots & & \\ \vdots & & \vdots & & \vdots \\ & P\_{r,d\_1} & P\_{r,d\_2} & P\_{r,d\_3} \end{array} \right] \quad \text{dimension } R \end{array} \tag{6}$$

The optimization model in general is formulated as for an integer programming, however, it includes some non-linear elements such as power losses or voltage drops. In order to improve computation time and reformulate the original problem to the Mixed Integer-Linear Programming (MILP), non-linear components are linearized. The quadratic power losses are linearized by the spline of five linear functions—Figure 1. Each linear function was created for different values of power flow in relation to power line capacity. The first function was created for power flow in rage 0–100% of line capacity. The rest of the functions were designed for power flow in range: *f2*: 20%–100%, *f3*: 40%–100%, *f4*: 60%–100% and *f5*: 80%–100% of line capacity. The final power losses were created as a sum of all mentioned linear functions.

**Figure 1.** Power losses linearization—spline function consisting of five linear functions.

Each linearising function refers to power line capacity, this is applicable for every power line of the grid model (7)–(11).

$$\begin{array}{cccc} \vee & \vee & \vee\\ t \in T & n, w \in N & f \mathbf{1}\_{n,w,t} = 0.2 \ast \frac{\text{Line.const}\_{\mathbf{u},\mathbf{z}} \ast R\_{\mathbf{u},\mathbf{w}}}{\text{U}\_{\mathbf{u},t}^2} \ast PF\_{\mathbf{u},\mathbf{w},t} \end{array} \tag{7}$$

∀ *t* ∈ *T* ∀ *<sup>n</sup>*, *<sup>w</sup>* <sup>∈</sup> *<sup>N</sup> <sup>f</sup>* <sup>2</sup>*n*,*w*,*<sup>t</sup>* <sup>=</sup> 0.4 <sup>∗</sup> *Line*.*constn*,*<sup>w</sup>* ∗ *Rn*,*<sup>w</sup> Un*,*t* <sup>2</sup> <sup>∗</sup> *PFn*,*w*,*<sup>t</sup>* <sup>−</sup> 0.08 <sup>∗</sup> *Line*.*const*<sup>2</sup> *<sup>n</sup>*,*<sup>w</sup>* ∗ *Rn*,*<sup>w</sup> Un*,*t* <sup>2</sup> (8)

$$\begin{array}{cccc} \vee & \vee & \vee\\ t \in T & n, w \in N & f \mathbf{3}\_{n,w,t} = 0.4 \ast \frac{\text{Line.const}\_{n,w} \ast R\_{\text{tr},w}}{\text{U}\_{\text{tr},t}^2} \ast P \mathbf{F}\_{w,k,t} - 0.16 \ast \frac{\text{Line.const}\_{n,w}^2 \ast R\_{\text{tr},w}}{\text{U}\_{\text{tr},t}^2} &\\ & & & (9) & \end{array} \tag{9}$$

$$\begin{array}{cccc} \vee & \vee & \vee\\ t \in T & n, w \in \mathbb{N} & f \mathbf{4}\_{n,w,t} = 0.4 \ast \frac{\text{Line.const}\_{n,w} \ast \text{R}\_{n,w}}{\text{U}\_{n,t}^{2}} \ast \text{PF}\_{w,k,t} - 0.24 \ast \frac{\text{Line.const}\_{n,w}^{2} \ast \text{R}\_{n,w}}{\text{U}\_{n,t}^{2}} &\\ & & & (10) \end{array} \tag{10}$$

∀ *t* ∈ *T* ∀ *<sup>n</sup>*, *<sup>w</sup>* <sup>∈</sup> *<sup>N</sup> <sup>f</sup>* <sup>5</sup>*n*,*w*,*<sup>t</sup>* <sup>=</sup> 0.4 <sup>∗</sup> *Line*.*constn*,*<sup>w</sup>* ∗ *Rn*,*<sup>w</sup> Un*,*t* <sup>2</sup> <sup>∗</sup> *PFw*,*k*,*<sup>t</sup>* <sup>−</sup> 0.32 <sup>∗</sup> *Line*.*const*<sup>2</sup> *<sup>n</sup>*,*<sup>w</sup>* ∗ *Rn*,*<sup>w</sup> Un*,*t* <sup>2</sup> (11)

Consequently, the energy lost on the power lines in each time step is as in (12),

$$\begin{array}{cccc} \vee & \vee & \vee \\ t \in T & n, w \in N \end{array} \begin{array}{c} P^{\text{losses}}\_{n,w,t} = f \mathbf{1}\_{n,w,t} + f \mathbf{2}\_{n,w,t} + f \mathbf{3}\_{n,w,t} + f \mathbf{4}\_{n,w,t} + f \mathbf{5}\_{n,w,t} \end{array} \tag{12}$$

while the total power losses in the entire time horizon are as (13).

$$\bigcup\_{t \in T}^{\mathsf{V}} P\_t^{\text{losses}} = \sum\_{n, \mathbf{w} \in \mathsf{N}} \left( P\_{n, \mathbf{w}, t}^{\text{losses}} \right) \tag{13}$$

The same assumption of constant nodal voltages and omission of reactive power and reactance allows for the linearization of voltage drops (14).

$$\begin{array}{cccc}\vee & \vee & \vee\\\text{it}\in T & n,w\in N\end{array} \quad \Delta l I\_{n,w,t} = \frac{PF\_{n,w,t} \* R\_{n,w}}{\mathcal{U}\_{n,t}}\tag{14}$$

Power flow in a power line, as well as nodal voltage, depends on power balance in each node depending on the grid structure and temporary power generation and demand. While the generation depends on weather conditions and generation structure, the demand for power in each time period depends on the load size and its type (referring to user behaviors).

$$\forall n \in N \qquad t \in T \quad \sum\_{n \in N} (PF\_{n, \nu, t}) = E\_{n, t}^{GEN\\_node} - \sum\_{l \in L} (Load\_{n, l} \* P\_{l, t}^{DEM}) \tag{15}$$
