3.2.1. Overloading Cost of the Transformer

Overloading of a transformer occurs due to higher loads at connection points such as, for example, with charging a large number of EVs or the operation of domestic HPs. Thermal overloading generally affects the insulation of the transformer windings; however, the involved dynamics in the loading enable the transformer to be overloaded for some time. Consequently, the amount of flexibility that needs to be procured should be aligned with the thermal status and the respective cost of overloading of the transformer.

Based on the provided load-profiles and historical values of the transformer load, the TA generates a set of probable loading scenarios to tackle the inherent uncertainties. For each of these scenarios, it estimates the imminent loss-of-life and respective overloading cost of the transformers. To do so, the following steps are followed:

Step 1: In this step, the TA generates the set of scenarios for probable loading. This is to counter the uncertainties associated with real-time loading, deviations from day-ahead schedule and forecast errors of the local generation technologies. In reality, correlations exist between the loads of consecutive hours. Therefore, instead of calculating single values for each time step, pseudo-random profiles of loading are generated based on Copula theory [33].

A copula is a multivariate probability distribution where the marginal-distributions are uniform. According to the Sklar's theorem, if those marginal distributions are continuous, the copula is unique. Let F be a 24-dimensional distribution function with continuous margins, (*F*1, *F*<sup>2</sup> ...... , *F*24). A copula

can be fitted, using the cumulative distribution functions (CDFs) of the historical load and the expected load profile provided by the aggregator. The marginal distribution here is the cumulative distribution functions generated by a Kernel Density Estimator (KDE). Then, the CDFs are used as inputs to generate the copula reflecting the correlation between the loads at different time steps. Next, the synthetic profiles generated from the fitted copula are transformed back to the original scale by applying the inverse cumulative distribution function. The whole process of scenario generation can be summarized by Figure 3.

**Figure 3.** Generating day-ahead loading scenarios using copula.

Step 2: The expected load is converted to the resulting hottest-spot temperature, θ*H*, of the transformer. According to IEEE Std C57.12.00-1993, the hottest-spot temperature is defined as the highest temperature of the winding at the operating condition, and is the main element for calculating the expected life of a transformer [34]. To determine this, a ratio called the load multiplex, *K*, is calculated and then used to calculate the top oil temperature rise, Δθ*TO*, and the hottest-spot temperature rise, Δθ*H*, as:

$$K = \frac{\text{Expected load}}{\text{Rated load}}\tag{10}$$

$$
\Delta\theta\_{TO} = \Delta\theta\_{TO,R} \left[ \frac{K^2 R + 1}{R + 1} \right] \tag{11}
$$

$$
\Delta\theta\_H = \Delta\theta\_{H,\mathcal{R}} \mathcal{K}^2 \tag{12}
$$

where Δθ*TO*,*<sup>R</sup>* is the top oil temperature rise at rated load, Δθ*H*,*<sup>R</sup>* is hottest-spot temperature at the rated load and *R* is the ratio of load loss at rated load to no-load loss at rated load. Then, θ*<sup>H</sup>* is found by summing the ambient temperature θ*<sup>A</sup>* with the above mentioned temperature rises.

$$
\Delta\theta\_H = \theta\_A + \Delta\theta\_{TO} + \Delta\theta\_H \tag{13}
$$

Step 3: In this step, the ageing acceleration factor (*FAA*) and equivalent ageing factor (*Feqv*) are calculated for the particular combination of loads and temperature for a duration of thirty min. This is due to the fact that degradation of the insulation is realized when the transformer is generally overloaded for half-an hour.

$$F\_{AA} = e^{\left[\frac{15000}{383} - \frac{15000}{\vartheta\_H + 273}\right]} \tag{14}$$

$$F\_{eqv} = \frac{\sum\_{n=1}^{N} F\_{AA} \Delta t\_n}{\sum\_{n=1}^{N} \Delta t\_n} \tag{15}$$

where *N* is the number of time intervals, the duration of each is Δ*tn* hours and normal insulation life of the transformer is *Tinl* hours. According to IEEE standards, the normal insulation life of a well dried, oxygen free distribution transformer is 180,000 h or 20.55 years [2,16].

Step 4: Based on the calculated equivalent factor, per unit loss of life of the transformer *Tlol* is determined.

$$T\_{lol} = \frac{F\_{eqv}t}{T\_{inl}}\tag{16}$$

The aging cost *Cag* can be determined by the loss of life with the total owning cost (TOC) *Co* of the transformer. The TOC method is considered to be one of the most cost and resource efficient methods for economic analysis of a transformer [35]. In addition to the initial cost of the transformer, TOC considers the operation and maintenance cost of the transformer and is calculated over the life span of the asset. The TOC can be determined from purchase cost *CP*, cost of no-load loss *CNL* and cost of load loss *CLL* of the transformer [35].

$$\mathbf{C}\_{\mathcal{o}} = \mathbf{C}\_{P} + \mathbf{C}\_{\mathrm{NL}} + \mathbf{C}\_{\mathrm{LL}} \tag{17}$$

$$\mathbb{C}\_{q\text{g}} = T\_{lol}\mathbb{C}\_{o} \tag{18}$$

If *Cag* is greater than the aging cost at nominal rating of the transformer *Cag*,*R*, then the overloading cost *COL* is determined from the arithmetic difference. Otherwise, the overloading cost is assumed to be zero.

$$\mathcal{C}^{OL} = \begin{cases} \mathcal{C}\_{\text{ag}} - \mathcal{C}\_{\text{ag},R} & \text{when } \mathcal{C}\_{\text{ag}} > \mathcal{C}\_{\text{ag},R} \\ 0 & \text{otherwise} \end{cases} \tag{19}$$

## 3.2.2. Calculation of Dynamic Tariff

The TA calculates the dynamic tariff that satisfies the thermal constraints for each of the loading scenarios *s* ∈ *S*. This can be mathematically presented as the following optimization problem:

$$\min \sum\_{s=1, s \in S}^{s=s\_N} \sum\_{t=1, t \in T}^{t=t\_N} (p\_t^{DSO} - p^{DSO}) P\_t^s \text{ } \forall s \in S \text{ if } (C\_{t,s}^{OL} > 0) \tag{20}$$

subject to,

$$
\overline{p\_t^{DSO}} = \overline{p^{DSO}}\tag{21}
$$

$$p\_t^{DSO} = \frac{1}{5} \sum\_{t\_f=t-2}^{t+2} p\_t^{DSO} + p\_{adj} \text{ if } \mathbb{C}\_{t\_f, s}^{OL} > 0 \tag{22}$$

$$p\_{adj,t} \le \max\limits(\mathbb{C}^{OL}\_{t,s}).\tag{23}$$

*pDSO <sup>t</sup>* is the decision variable and denotes the dynamic network tariff at each time step. *Ps <sup>t</sup>* is the load at time *t* ∈ *T* in scenario, *s*. Constraints in Equation (21) dictate that the average network tariff should be the same for both normal operations and tariff-based DR cases. Constraints in Equations (22) and (23) limit the range of tariff adjustment in time and monetary values, respectively.

The resulting dynamic network tariff is sent to the Aggregator for updating the price for the end-users. This process is continued as an iterative process until the issues are resolved. In this work, we have limited the process to up to ten iterations. The whole process can be represented by the flowchart as shown in Figure 4.

**Figure 4.** Flowchart depicting the processes.

#### **4. Simulation Setup**
