# The stickiness property
for antisymmetric nonlocal minimal graphs

Benjamin Baronowitz  ,  Serena Dipierro  and  Enrico Valdinoci  Department of
Mathematics and Statistics
University of Western Australia
35 Stirling Highway, WA 6009 Crawley, Australia<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract.

We show that arbitrarily small antisymmetric perturbations of the zero
function are sufficient to produce the stickiness phenomenon for planar
nonlocal minimal graphs (with the same quantitative bounds obtained for the
case of even symmetric perturbations, up to multiplicative constants).

In proving this result, one also establishes an odd symmetric version of the
maximum principle for nonlocal minimal graphs, according to which the odd
symmetric minimizer is positive in the direction of the positive bump and
negative in the direction of the negative bump.

###### Key words and phrases:

Nonlocal minimal surfaces, stickiness, qualitative and quantitative behavior.

###### 2020 Mathematics Subject Classification:

35R11, 49Q05

S. Dipierro and E. Valdinoci are members of AustMS. S. Dipierro is supported
by the Australian Research Council DECRA DE180100957 “PDEs, free boundaries
and applications”. E. Valdinoci is supported by the Australian Laureate
Fellowship FL190100081 “Minimal surfaces, free boundaries and partial
differential equations”.

##  1\. Introduction

Nonlocal minimal surfaces were introduced in [MR2675483] as minimizers of a
fractional perimeter functional with respect to some given external datum. As
established in [MR3516886], when the external datum is a graph in some
direction, so is the whole minimizer, hence it is also common to consider the
case of nonlocal minimal graphs, i.e. of nonlocal minimal surfaces which
possess a graphical structure.

Nonlocal minimal graphs exhibit a quite peculiar phenomenon, discovered in
[MR3596708], called “stickiness”. Roughly speaking, different from the
classical minimal surfaces, in the fractional setting remote interactions are
capable of producing boundary discontinuities, which in turn are complemented
by boundary divergence of the derivative of the nonlocal minimal graph.

We focus here on the case in which the ambient space is of dimension $2$: this
is indeed the simplest possible stage to detect interesting geometric patterns
and, differently from the classical case, it already provides a number of
difficulties since in the nonlocal framework segments are in general not the
boundary of nonlocal minimal objects. Also, in the plane the stickiness
phenomenon is known to be essentially “generic” with respect to the external
datum, see [MR4104542], in the sense that boundary discontinuities and
boundary singularities for the derivatives can be produced by arbitrarily
small perturbations of a given external datum. In particular, these behaviors
can be obtained by arbitrarily small perturbations of the flat case in which
the external datum is, say, the sublevel of the zero function.

See also [MR4178752] for an analysis of nonlocal minimal graphs in dimension
$3$, and [MR3926519, MR4184583, 2020arXiv201000798D] for other examples of
stickiness.

Up to now, all the examples of stickiness in the setting of graphs were
constructed by using “one side bumps” in the perturbation (for instance,
adding suitable positive bumps to the zero function, or to a given function
outside a vertical slab). In this paper, we construct examples of stickiness
in which the datum is antisymmetric. In a nutshell, we will consider
arbitrarily small perturbations of the zero function by bumps possessing odd
symmetry: in this case, in principle, one may fear that the effects of equal
and opposite bumps would cancel each other and prevent the stickiness
phenomenon to occur, but we will instead establish that the stickiness
phenomenon is persistent also in this class of antisymmetric perturbations of
the flat case.

Also, we provide quantitative bounds on the resulting boundary discontinuity
which turn out to be as good as the ones available for positive bumps (up to
multiplicative constants).

The mathematical notation used in this paper goes as follows. Given
$s\in(0,1)$, $a<b$ and $u_{0}\in C(\mathbb{R})$, we say that $u\in
L^{\infty}_{\rm loc}(\mathbb{R})$ is an $s$-minimal graph in $(a,b)$ if:

  * •

$u(x)=u_{0}(x)$ for every $x\in\mathbb{R}\setminus(a,b)$,

  * •

for every $L>\|u\|_{L^{\infty}((a,b))}$ and every measurable function
$v:\mathbb{R}\to\mathbb{R}$ such that $|v(x)|\leqslant L$ for all $x\in(a,b)$
we have that

$\mbox{\rm{Per}}_{s}(E_{u},\Omega_{L})\leqslant\mbox{\rm{Per}}_{s}(E_{v},\Omega_{L})$

where

$\displaystyle E_{u}:=\big{\\{}(x,y)\in\mathbb{R}^{2}{\mbox{ s.t.
}}y<u(x)\big{\\}},$ $\displaystyle
E_{v}:=\big{\\{}(x,y)\in\mathbb{R}^{2}{\mbox{ s.t. }}y<v(x)\big{\\}},$
$\displaystyle\Omega_{L}:=(a,b)\times(-L,L),$
$\displaystyle\mbox{\rm{Per}}_{s}(E,\Omega):={\mathcal{L}}_{s}(E\cap\Omega,E^{c}\cap\Omega)+{\mathcal{L}}_{s}(E\cap\Omega,E^{c}\cap\Omega^{c})+{\mathcal{L}}_{s}(E\cap\Omega^{c},E^{c}\cap\Omega),$
$\displaystyle\Omega^{c}:=\mathbb{R}^{2}\setminus\Omega$ and
$\displaystyle{\mathcal{L}}_{s}(A,B):=\iint_{A\times
B}\frac{dX\,dY}{|X-Y|^{2+s}}.$

See Figure 1 for a sketch of the interactions detected by the definition
above. See also [MR4279395] for further information about nonlocal minimal
graphs and in particular Theorem 1.11 for existence and uniqueness results for
nonlocal minimal graphs.

We also recall that, while the minimizers of the classical perimeter have
vanishing mean curvature, a minimizer $E$ in $\Omega$ of the fractional
perimeter has vanishing fractional mean curvature, in the sense that, for
every $P\in\Omega\cap(\partial E)$,

$\int_{\mathbb{R}^{2}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX=0.$

Here above and in the sequel, the singular integral is intended in the Cauchy
principal value sense and the above equation reads in the sense of viscosity,
see [MR2675483, Theorem 5.1] for full details.

Numerical examples of the stickiness properties have recently been provided by
[MR3982031, MR4294645]. The theory of nonlocal minimal surfaces also presents
a number of interesting offsprings, such as regularity theory [MR3090533,
MR3107529, MR3331523, MR3680376, MR3798717, MR3934589, MR3981295, MR4050198,
MR4116635], isoperimetric problems and the study of constant mean curvature
surfaces [MR2799577, MR3322379, MR3412379, MR3485130, MR3744919, MR3836150,
MR3881478], geometric evolution problems [MR2487027, MR3156889, MR3713894,
MR3778164, MR3951024], front propagation problems [MR2564467], phase
transition problems [MR3133422], etc.

Figure 1. The interactions contributing to the fractional perimeter.

Having formalized the setting in which we work, we can now state our result
about the stickiness phenomenon in the antisymmetric setting as an arbitrarily
small perturbation of the flat case:

######  Theorem 1.1.

Let $\epsilon_{0}>0$, $\bar{C}>0$, $d>d_{0}>0$ and $h>0$, with

(1.1)
$\frac{h^{2+s}}{(h^{2}+2)^{\frac{2+s}{2}}}>\frac{2\bar{C}}{(1+s){\left(d-d_{0}\right)^{1+s}}}.$

Then, there exist a constant $C>0$, depending only on $s$, $\epsilon_{0}$,
$\bar{C}$, $d$ and $h$, and an infinitesimal sequence of $\eta$’s for which
the following claim holds true.

Let $u_{0}\in C(\mathbb{R})$. Assume that

(1.2) $u_{0}(x)=0\quad{\mbox{ for all }}x\in(d,d+h),$ (1.3) $u_{0}(x)\leqslant
0\quad{\mbox{ for all }}x\in(d,+\infty),$ (1.4)
$\|u_{0}\|_{L^{\infty}(\mathbb{R})}\leqslant\bar{C}\eta,$ (1.5)
$u_{0}(x)=-u_{0}(-x)\quad{\mbox{ for all }}x\in(d,+\infty)$

and

(1.6)
$\int_{-\infty}^{-d-h}\frac{u_{0}(x)}{|d_{0}-d-x|^{2+s}}\,dx\geqslant\eta.$

Let $u$ be the $s$-minimal graph with $u=u_{0}$ in
$\mathbb{R}\setminus(-d,d)$.

Then,

(1.7) $u(x)=-u(-x)\quad{\mbox{ for all }}x\in(-\infty,+\infty),$ (1.8)
$u(x)\geqslant 0\quad{\mbox{ for all }}x\in(-\infty,0],$ (1.9) $u(x)\leqslant
0\quad{\mbox{ for all }}x\in[0,+\infty)$

and, for every $x\in(-d,0)$,

(1.10) $u(x)\geqslant\frac{\eta^{\frac{2+\epsilon_{0}}{1-s}}}{C}.$

Figure 2. The antisymmetric stickiness phenomenon, as given in Theorem 1.1.

See Figure 2 for a sketch of the geometric scenario described in Theorem 1.1.
Interestingly, the antisymmetric stickiness power $\frac{2+\epsilon_{0}}{1-s}$
in (1.10) is the same as the one detected in [MR3596708, Theorem 1.4] for the
even symmetric case (hence, the antisymmetric case appears to be also
quantitative in agreement with the even symmetric configurations, up to
multiplicative constants).

Notice also that condition (1.1) simply states that $d$ is sufficiently large
with respect to other structural parameters. For instance, one can take

$\displaystyle\bar{C}:=2^{2+s}(1+s),\qquad
d:=3^{\frac{2+s}{2(1+s)}}2^{\frac{3+s}{1+s}}+2\qquad d_{0}:=1\quad{\mbox{ and
}}\quad h:=1.$

Notice that with this choice of parameters, condition (1.1) is satisfied. We
also define

$d_{1}:=2\left(\left(\frac{10}{9}\right)^{\frac{1}{1+s}}-1\right)\quad{\mbox{
and }}\quad d_{2}:=2\left(10^{\frac{1}{1+s}}-1\right),$

and we take $u_{0}\in C(\mathbb{R})$ satisfying (1.2), (1.3), (1.4) and (1.5),
and such that

$u_{0}\geqslant\bar{C}\eta\chi_{(-d-h-d_{2},-d-h-d_{1})}\quad{\mbox{ in
}}(-\infty,-d).$

With this choice of $u_{0}$ we have that (1.6) is also satisfied. Indeed,
using the change of variable $z:=d_{0}-d-x$,

$\displaystyle\int_{-\infty}^{-d-h}\frac{u_{0}(x)}{|d_{0}-d-x|^{2+s}}\,dx\geqslant\int_{-d-h-
d_{2}}^{-d-h-d_{1}}\frac{u_{0}(x)}{|d_{0}-d-x|^{2+s}}\,dx$
$\displaystyle\qquad\qquad\geqslant\int_{-d-h-d_{2}}^{-d-h-
d_{1}}\frac{\bar{C}\eta}{|d_{0}-d-x|^{2+s}}\,dx=\bar{C}\eta\int_{d_{0}+h+d_{1}}^{d_{0}+h+d_{2}}\frac{dz}{z^{2+s}}$
$\displaystyle\qquad\qquad=\frac{\bar{C}\eta}{1+s}\left[\frac{1}{(d_{0}+h+d_{1})^{1+s}}-\frac{1}{(d_{0}+h+d_{2})^{1+s}}\right]$
$\displaystyle\qquad\qquad=2^{2+s}\eta\left(\frac{1}{2^{1+s}\frac{10}{9}}-\frac{1}{2^{1+s}10}\right)=2\eta\frac{8}{10}>\eta,$

as desired.

The gist of the proof of Theorem 1.1 is to employ the auxiliary function
constructed in [MR3596708, Corollary 7.2], which, by suitably scaling the
picture, provides a localized barrier, whose nonlocal mean curvature has
possibly a “wrong sign” somewhere, but this sign discrepancy is controlled by
a constant times delta. However, since in our framework the solution will
cross the horizontal axis, the previous barrier cannot be exploited as it is
and needs to be modified to stay below the “expected” negative regions of the
true solution.

To this end, we first need to localize a “safe region near the boundary” for
the true solution, that is an interval of well-determined length in which one
can be sure that the solution is positive. The detection of this safe region
will be also obtained by a barrier argument, with a different barrier based on
the analysis of smooth functions with small oscillations.

The objective then becomes to slide from below an appropriate barrier whose
positive portion occurs precisely in the above safe region. For this one has
to scale and modify the barrier in [MR3596708, Corollary 7.2] by accounting
for the external data’s bumps: roughly speaking, the positive bump on the left
will produce an advantageous term, while the negative bump on the right and
the modification needed to lower the barrier outside the safe region will
produce a disadvantageous term of the same order. Thus, playing around with
constants, one detects a natural structure assuring that the advantageous
terms is greater than the sum of the initial error on the nonlocal mean
curvature, the contribution of the negative bump and the error coming from the
barrier modification.

In view of these comments, to localize the safe region near the boundary and
perform the proof of Theorem 1.1, we establish a result of general interest,
which can be seen as a “maximum principle for antisymmetric $s$-minimal
graphs”. Roughly speaking, the classical maximum principle for $s$-minimal
graphs (see e.g. [MR2675483, MR3516886]) states that if the external data of
an $s$-minimal graph are positive (or negative) then so is the $s$-minimal
graph. Here we deal with an antisymmetric configuration, hence it is not
possible for the external data of an $s$-minimal graph to have a sign (except
in the trivial case of identically zero datum, which produces the identically
zero $s$-minimizer). Hence, the natural counterpart of the maximum principle
in the antisymmetric framework is to assume that the external datum “on one
side” has a sign, thus forcing the external datum on the other sign to have
the opposite sign. Under this assumption, we show that the corresponding
$s$-minimal graph maintains the antisymmetry and sign properties of the datum,
according to the following result:

######  Theorem 1.2.

Let $d$, $h>0$. Let $u:\mathbb{R}\to\mathbb{R}$ be an $s$-minimal graph in
$(-d,d)$, with $u\in C(-\infty,-d)\cap C^{1,\frac{1+s}{2}}(-d-h,-d)$. Assume
that

(1.11) $u(x)=-u(-x)\quad{\mbox{ for all }}x\in(d,+\infty)$

and that

(1.12) $u(x)\leqslant 0\quad{\mbox{ for all }}x\in(d,+\infty).$

Then,

(1.13) $u(x)=-u(-x)\quad{\mbox{ for all }}x\in(-\infty,+\infty),$ (1.14)
$u(x)\geqslant 0\quad{\mbox{ for all }}x\in(-\infty,0]$

and

(1.15) $u(x)\leqslant 0\quad{\mbox{ for all }}x\in[0,+\infty).$

Figure 3. The maximum principle for antisymmetric nonlocal minimal graphs, as
given in Theorem 1.2.

See Figure 3 for a sketch of the configuration detected in Theorem 1.2.

The proofs of Theorems 1.1 and 1.2 are contained in Sections 2 and 3,
respectively.

##  2\. Proof of Theorem 1.1

We give here the proof of Theorem 1.1. At this stage, we freely use Theorem
1.2, whose proof is postponed to Section 3.

We also use the notation $X=(x,y)$ to denote points in $\mathbb{R}^{2}$, with
$x$, $y\in\mathbb{R}$.

Let $\epsilon_{0}>0$. Given $\delta\in(0,1)$ sufficiently small, we apply
[MR3596708, Corollary 7.2] and we find a set $H_{\delta}\subset\mathbb{R}^{2}$
that contains the halfplane $\mathbb{R}\times(-\infty,0)$ and is contained in
the halfplane $\mathbb{R}\times(-\infty,\delta)$, such that

$\displaystyle\left\\{(x,y)\in H_{\delta}{\mbox{ s.t.
}}x\in(-\infty,-d)\cup\left(-d+d_{0},+\infty\right)\right\\}$
$\displaystyle\qquad=\left\\{(x,y)\in\mathbb{R}^{2}{\mbox{ s.t.
}}x\in(-\infty,-d)\cup\left(-d+d_{0},+\infty\right){\mbox{ and
}}y\in(-\infty,0)\right\\},$

with

(2.1) $\begin{split}&\left\\{(x,y)\in H_{\delta}{\mbox{ s.t.
}}x\in\left(-d,-d+d_{0}\right)\right\\}\\\
&\qquad\supseteq\left\\{(x,y)\in\mathbb{R}^{2}{\mbox{ s.t.
}}x\in\left(-d,-d+d_{0}\right)\times\left(-\infty,\delta^{\frac{2+\epsilon_{0}}{1-s}}\right)\right\\}\end{split}$

and satisfying, for each $P=(p,q)\in\partial H_{\delta}$ with
$p\in\left(-d,-d+d_{0}\right)$,

(2.2) $\int_{\mathbb{R}^{2}}\frac{\chi_{\mathbb{R}^{2}\setminus
H_{\delta}}(X)-\chi_{H_{\delta}}(X)}{|X-P|^{2+s}}\,dX\leqslant C\delta,$

for some constant $C>0$ depending only on $s$, $\epsilon_{0}$ and $d$.

The strategy will be to take

(2.3) $\eta:=C_{\star}\delta$

in the statement of Theorem 1.1, with $C_{\star}>0$ to be sufficiently large,
according to the following computation. We let

$\displaystyle{\mathcal{A}}:=\big{\\{}(x,y)\in\mathbb{R}^{2}{\mbox{ s.t.
}}x\in(-\infty,-d){\mbox{ and }}0<y<u_{0}(x)\big{\\}},$
$\displaystyle{\mathcal{B}}:=\big{\\{}(x,y)\in\mathbb{R}^{2}{\mbox{ s.t.
}}x\in(0,+\infty){\mbox{ and }}-\bar{C}\eta<y<0\big{\\}}$ and $\displaystyle
U_{\delta}:=(H_{\delta}\cup{\mathcal{A}})\setminus{\mathcal{B}},$

see Figure 4.

Figure 4. The set $U_{\delta}$ used in the proof of Theorem 1.1.

We observe that ${\mathcal{A}}\cap H_{\delta}=\varnothing$ and
${\mathcal{B}}\subset H_{\delta}$. Consequently,

$\displaystyle\chi_{U_{\delta}}=\chi_{H_{\delta}}+\chi_{{\mathcal{A}}}-\chi_{{\mathcal{B}}}$
and $\displaystyle\chi_{\mathbb{R}^{2}\setminus
U_{\delta}}=\chi_{\mathbb{R}^{2}\setminus(H_{\delta}\cup{\mathcal{A}})}+\chi_{{\mathcal{B}}}\leqslant\chi_{\mathbb{R}^{2}\setminus
H_{\delta}}+\chi_{{\mathcal{B}}}.$

From this and (2.2), for all $P=(p,q)\in\partial H_{\delta}$ with
$p\in\left(-d,-d+d_{0}\right)$,

(2.4) $\begin{split}&\int_{\mathbb{R}^{2}}\frac{\chi_{\mathbb{R}^{2}\setminus
U_{\delta}}(X)-\chi_{U_{\delta}}(X)}{|X-P|^{2+s}}\,dX\\\
&\qquad\leqslant\int_{\mathbb{R}^{2}}\frac{\chi_{\mathbb{R}^{2}\setminus
H_{\delta}}(X)+\chi_{{\mathcal{B}}}(X)-\big{(}\chi_{H_{\delta}}(X)+\chi_{{\mathcal{A}}}(X)-\chi_{{\mathcal{B}}(X)\big{)}}}{|X-P|^{2+s}}\,dX\\\
&\qquad\leqslant
C\delta+2\int_{{\mathcal{B}}}\frac{dX}{|X-P|^{2+s}}-\int_{{\mathcal{A}}}\frac{dX}{|X-P|^{2+s}}\end{split}$

We note that, for all $P=(p,q)\in\partial H_{\delta}$ with
$p\in\left(-d,-d+d_{0}\right)$,

$\displaystyle\int_{{\mathcal{B}}}\frac{dX}{|X-P|^{2+s}}\leqslant\int_{{\mathcal{B}}}\frac{dX}{|x-p|^{2+s}}=\bar{C}\eta\int_{0}^{+\infty}\frac{dx}{|x-p|^{2+s}}=\bar{C}\eta\int_{0}^{+\infty}\frac{dx}{(x-p)^{2+s}}$
$\displaystyle\qquad\leqslant\bar{C}\eta\int_{0}^{+\infty}\frac{dx}{\left(x+d-d_{0}\right)^{2+s}}=\frac{\bar{C}\eta}{(1+s){\left(d-d_{0}\right)^{1+s}}}=\frac{\bar{C}C_{\star}\delta}{(1+s){\left(d-d_{0}\right)^{1+s}}}.$

This and (2.4) yield that

(2.5) $\int_{\mathbb{R}^{2}}\frac{\chi_{\mathbb{R}^{2}\setminus
U_{\delta}}(X)-\chi_{U_{\delta}}(X)}{|X-P|^{2+s}}\,dX\leqslant
C\delta+\frac{2\bar{C}C_{\star}\delta}{(1+s){\left(d-d_{0}\right)^{1+s}}}-\int_{{\mathcal{A}}}\frac{dX}{|X-P|^{2+s}}.$

Also, by (1.4), if $X=(x,y)\in{{\mathcal{A}}}$ and $P=(p,q)\in\partial
H_{\delta}$ with $p\in\left(-d,-d+d_{0}\right)$,

$|y-q|\leqslant|y|+|q|\leqslant 1+\delta\leqslant 2.$

Moreover, owing to (1.2) and (1.5), for all $x\in(-d-h,-d)$ we have that
$u_{0}(x)=-u_{0}(-x)=0$.

Accordingly, if $X=(x,y)\in{{\mathcal{A}}}$ then $x\leqslant-d-h$. As a
result, if $X=(x,y)\in{{\mathcal{A}}}$ and $P=(p,q)\in\partial H_{\delta}$
with $p\in\left(-d,-d+d_{0}\right)$,

$\begin{split}&|x-p|\geqslant\min\big{\\{}|x-p|,\,|x+p|\big{\\}}\geqslant\min\big{\\{}p-x,\,-p-x\big{\\}}\\\
&\qquad\geqslant\min\left\\{-d-(-d-h),\,d-d_{0}-(-d-h)\right\\}=\min\left\\{h,2d-d_{0}+h\right\\}=h.\end{split}$

From these remarks we obtain that

$\displaystyle|X-P|=|x-p|\sqrt{1+\frac{|y-q|^{2}}{|x-p|^{2}}}\leqslant|x-p|\sqrt{1+\frac{2}{h^{2}}}=|x-p|\sqrt{\frac{h^{2}+2}{h^{2}}}.$

Hence,

(2.6)
$\begin{split}&\int_{{\mathcal{A}}}\frac{dX}{|X-P|^{2+s}}\geqslant\frac{h^{2+s}}{(h^{2}+2)^{\frac{2+s}{2}}}\int_{{\mathcal{A}}}\frac{dX}{|x-p|^{2+s}}\\\
&\qquad\qquad=\frac{h^{2+s}}{(h^{2}+2)^{\frac{2+s}{2}}}\int_{-\infty}^{-d}\frac{u_{0}(x)}{|x-p|^{2+s}}\,dx=\frac{h^{2+s}}{(h^{2}+2)^{\frac{2+s}{2}}}\int_{-\infty}^{-d-h}\frac{u_{0}(x)}{|x-p|^{2+s}}\,dx.\end{split}$

Since, if $p\in\left(-d,-d+d_{0}\right)$ and $x\in(-\infty,-d-h)$,

$|x-p|=p-x\leqslant-d+d_{0}-x=|d_{0}-d-x|,$

we deduce from (2.6) that

$\int_{{\mathcal{A}}}\frac{dX}{|X-P|^{2+s}}\geqslant\frac{h^{2+s}}{(h^{2}+2)^{\frac{2+s}{2}}}\int_{-\infty}^{-d-h}\frac{u_{0}(x)}{|d_{0}-d-x|^{2+s}}\,dx.$

This, (1.6) and (2.3) give that

$\int_{{\mathcal{A}}}\frac{dX}{|X-P|^{2+s}}\geqslant\frac{h^{2+s}\eta}{(h^{2}+2)^{\frac{2+s}{2}}}=\frac{C_{\star}h^{2+s}\delta}{(h^{2}+2)^{\frac{2+s}{2}}}.$

Plugging this information into (2.5) we infer that

(2.7) $\int_{\mathbb{R}^{2}}\frac{\chi_{\mathbb{R}^{2}\setminus
U_{\delta}}(X)-\chi_{U_{\delta}}(X)}{|X-P|^{2+s}}\,dX\leqslant
C\delta+\frac{2\bar{C}C_{\star}\delta}{(1+s){\left(d-d_{0}\right)^{1+s}}}-\frac{C_{\star}h^{2+s}\delta}{(h^{2}+2)^{\frac{2+s}{2}}}.$

Now we define

$\vartheta:=\frac{h^{2+s}}{(h^{2}+2)^{\frac{2+s}{2}}}-\frac{2\bar{C}}{(1+s){\left(d-d_{0}\right)^{1+s}}}$

and we observe that $\vartheta>0$, owing to (1.1). Then, we deduce from (2.7)
that

(2.8) $\int_{\mathbb{R}^{2}}\frac{\chi_{\mathbb{R}^{2}\setminus
U_{\delta}}(X)-\chi_{U_{\delta}}(X)}{|X-P|^{2+s}}\,dX\leqslant C\delta-
C_{\star}\vartheta\delta\leqslant-\frac{C_{\star}\vartheta\delta}{2}<0,$

as long as $C_{\star}$ is large enough.

We can now use $U_{\delta}$ as a barrier for the sliding method. Specifically,
by (1.4) and the maximum principle in [MR2675483], we know that
$|u(x)|\leqslant\bar{C}\eta$ for all $x\in\mathbb{R}$ and therefore a
downwards translation of $U_{\delta}$ with magnitude $1$, that we denote by
$U_{\delta}-(0,1)$, is completely contained in
$E_{u}:=\big{\\{}(x,y)\in\mathbb{R}^{2}{\mbox{ s.t. }}y<u(x)\big{\\}}$. We
then slide up till we reach a touching point: given $\tau\in(0,1]$, we notice
that by construction no touching can occur between $E_{u}$ and
$U_{\delta}-(0,\tau)$ at points with abscissa in $\mathbb{R}\setminus[-d,0]$.
But these touching points cannot occur in $\left(-d+d_{0},0\right]$ either,
thanks to (1.14). And they cannot occur in $\left[-d,-d+d_{0}\right]$, thanks
to (2.8) and the maximum principle in [MR2675483] (see also [MR4104542,
Theorem 1.4]).

As a consequence, we have that $U_{\delta}\subseteq E_{u}$. By virtue of (2.1)
and (2.3), this establishes (1.10).

Also, the claims in (1.7), (1.8) and (1.9) follow from (1.13), (1.14) and
(1.15). $\Box$

##  3\. Proof of Theorem 1.2

The claim in (1.13) about the odd symmetry of $u$ is a consequence of (1.11)
and [MR3596708, Lemma A.1].

Thus, to complete the proof of Theorem 1.2 it suffices to establish (1.14),
since (1.15) would then follow from (1.13) and (1.14).

To prove (1.14) we argue by contradiction and assume that (1.14) is violated.
We recall that $u$ is uniformly continuous in $(-d,d)$, owing to [MR3516886,
Theorem 1.1]. As a consequence, we can redefine $u$ at the extrema of $(-d,d)$
by setting

(3.1) $u(-d):=\lim_{x\searrow-d}u(x)\quad{\mbox{and}}\quad
u(d):=\lim_{x\nearrow d}u(x)$

and we have that

(3.2) $u\in C([-d,d]).$

Furthermore,

(3.3) if $x\in(-\infty,-d)$ then $u(x)=-u(-x)\geqslant 0$,

thanks to (1.11) and (1.13).

Figure 5. The geometry involved in the proof of Theorem 1.2 (Part 1:
detachment from the boundary).

This, (3.2) and our contradictory assumption give that there exists
$p\in[-d,0]$ such that

(3.4) $u(p)=\min_{[-d,0]}u<0.$

We claim that

(3.5) $p\in(-d,0).$

Indeed, since $u(0)=-u(0)$ by (1.13), we have that $u(0)=0$ and accordingly
$p\neq 0$. Consequently, to prove (3.5) it remains to show that $p\neq-d$. For
this, suppose that $p=-d$. This, together with (3.3), gives that

(3.6) $\lambda:=\lim_{x\nearrow-d}u(x)\geqslant 0>u(-d),$

see Figure 5.

Hence, by [MR4104542, Corollary 1.3(ii)] (see also [MR3532394] for related
results), we infer that there exists $\mu>0$ and a continuous function $v$ in
$\big{(}u(-d)-\mu,u(-d)+\mu\big{)}$ such that $v$ is the inverse function of
$u$ (in the above domain of definition, with the understanding that $v=-d$
along the jump discontinuity of $u$ detected in (3.6)).

Now we pick a sequence of points $q_{k}\in(-d,0)$ with $q_{k}\searrow-d$ as
$k\to+\infty$; in particular, by (3.1), we can assume that
$u(q_{k})\in\big{(}u(-d)-\mu,u(-d)+\mu\big{)}$ and we thus obtain that
$v(u(q_{k}))=q_{k}$. We have that

(3.7) $u(q_{k})<u(-d),$

otherwise, since $v(y)=-d$ for all
$y\in\big{[}u(-d),\min\\{\lambda,u(-d)+\mu\\}\big{]}$, we would have that
$-d=v(u(q_{k}))=q_{k}>-d$, which is a contradiction.

Having proved (3.7), we find that

$u(q_{k})<u(-d)=u(p)=\min_{[-d,0]}u\leqslant u(q_{k})$

and this is a contradiction that completes the proof of (3.5).

We now recall that $u\in C^{\infty}(-d,d)$, thanks to [MR3090533, MR3331523].
For this reason, we can write the Euler-Lagrange equation for nonlocal minimal
surfaces (see [MR2675483, Theorem 5.1]) in a pointwise sense in $(-d,d)$. In
particular, by (3.5),

(3.8) $\int_{\mathbb{R}^{2}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX=0,$

where $P:=(p,u(p))$ and

$E:=\big{\\{}(x,y)\in\mathbb{R}^{2}{\mbox{ s.t. }}y<u(x)\big{\\}}.$

We will now reach a contradiction by showing that (3.8) is violated, since
“the set $E$ is too large”. To this end, we consider the sets

(3.9) $\begin{split}&F_{-}:=\big{\\{}x\in\mathbb{R}{\mbox{ s.t.
}}u(x)<u(p)\big{\\}}\\\ {\mbox{and
}}\quad&F_{+}:=\big{\\{}x\in\mathbb{R}{\mbox{ s.t.
}}u(x)>-u(p)\big{\\}},\end{split}$

see Figure 6 (notice in particular that the sets $F_{-}$ and $F_{+}$ are given
by the bold intervals along the horizontal axis in Figure 6).

Figure 6. The geometry involved in the proof of Theorem 1.2 (Part 2: detection
of cancellations via isometric regions). In particular, the sets $F_{-}$ and
$F_{+}$ in (3.9) are given by the bold intervals along the horizontal axis.

We notice that, in light of (1.12), (3.3) and (3.4),

$F_{-}\subseteq(0,+\infty)\quad{\mbox{ and }}\quad F_{+}\subseteq(-\infty,0).$

Furthermore, by (1.13), we have that

(3.10) $x\in F_{-}$ if and only if $-x\in F_{+}$.

The strategy now is based on the detection of suitable cancellations via
isometric regions that correspond to the sets $F_{-}$ and $F_{+}$. The gist is
to get rid of the contributions arising from the complement of $E$ below the
line $\\{y=u(p)\\}$. For this, one has to use suitable transformations
inherited from the geometry of the problem, such as an odd reflection through
the origin, a vertical translation of magnitude $2u(p)$ and the reflection
along the line $\\{y=u(p)\\}$. The use of these transformations will aim, on
the one hand, at detecting isometric regions in $E$ which cancel the ones in
the complement of $E$ in the integral computations, and, on the other hand, at
maintaining a favorable control of the distance from the point $P$, since this
quantity appears at the denominator of the integrands involved.

To employ this strategy, we use the notation $X=(x,y)$ and the transformation

$\widetilde{X}=(\widetilde{x},\widetilde{y})=\widetilde{X}(X):=\big{(}-x,2u(p)-y\big{)}$

and we claim that

(3.11) if $x\geqslant 0$, then $|X-P|\geqslant|\widetilde{X}-P|$.

To check this, we observe that when $x\geqslant 0$ it follows that
$xp\leqslant 0$ and therefore

$\displaystyle|\widetilde{X}-P|^{2}-|X-P|^{2}$ $\displaystyle=$
$\displaystyle\Big{|}\big{(}-x,2u(p)-y\big{)}-\big{(}p,u(p)\big{)}\Big{|}^{2}-\Big{|}(x,y)-\big{(}p,u(p)\big{)}\Big{|}^{2}$
$\displaystyle=$
$\displaystyle\Big{|}\big{(}-x-p,u(p)-y\big{)}\Big{|}^{2}-\Big{|}\big{(}x-p,y-u(p)\big{)}\Big{|}^{2}$
$\displaystyle=$ $\displaystyle(-x-p)^{2}-(x-p)^{2}$ $\displaystyle=$
$\displaystyle 4xp$ $\displaystyle\leqslant$ $\displaystyle 0,$

which establishes (3.11).

Besides, we define

$\displaystyle G:=\Big{\\{}X=(x,y)\in\mathbb{R}^{2}\setminus E{\mbox{ s.t.
$x\in F_{-}$ and $y<2u(p)-u(x)$}}\Big{\\}}$

and we claim that

(3.12) if $X\in G$, then $\widetilde{X}\in E\cap\\{\widetilde{x}\in
F_{+}\\}\cap\\{-u(\widetilde{x})<\widetilde{y}<2u(p)+u(\widetilde{x})\\}$.

Indeed, if $X\in G$, by (1.13) we have that $y>u(x)=-u(-x)$ and accordingly

(3.13) $\widetilde{y}=2u(p)-y<2u(p)+u(-x)=2u(p)+u(\widetilde{x}).$

In particuar, $\widetilde{y}<u(\widetilde{x})$. From this, (3.10) and (3.13)
we obtain that

(3.14) $\widetilde{X}\in E\cap\\{\widetilde{x}\in
F_{+}\\}\cap\\{\widetilde{y}<2u(p)+u(\widetilde{x})\\}.$

Additionally,

$\widetilde{y}=2u(p)-y>2u(p)-\big{(}2u(p)-u(x)\big{)}=u(x)=u(-\widetilde{x})=-u(\widetilde{x}).$

Combining this with (3.14) we obtain (3.12), as desired.

Now, in light of (3.11) and (3.12), and observing that $d\widetilde{X}=dX$,

(3.15) $\begin{split}&\int_{\mathbb{R}^{2}\cap\\{x\in
F_{-}\\}\cap\\{y<2u(p)-u(x)\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)}{|X-P|^{2+s}}\,dX\\\
&\qquad=\int_{G}\frac{dX}{|X-P|^{2+s}}\leqslant\int_{G}\frac{dX}{|\widetilde{X}(X)-P|^{2+s}}\\\
&\qquad\leqslant\int_{E\cap\\{\widetilde{x}\in
F_{+}\\}\cap\\{-u(\widetilde{x})<\widetilde{y}<2u(p)+u(\widetilde{x})\\}}\frac{d\widetilde{X}}{|\widetilde{X}-P|^{2+s}}\\\
&\qquad=\int_{\mathbb{R}^{2}\cap\\{\widetilde{x}\in
F_{+}\\}\cap\\{-u(\widetilde{x})<\widetilde{y}<2u(p)+u(\widetilde{x})\\}}\frac{\chi_{E}(\widetilde{X})}{|\widetilde{X}-P|^{2+s}}\,d\widetilde{X}.\end{split}$

Now we use the even reflection across the line $\\{y=u(p)\\}$, that is we
define

(3.16) $X_{*}=(x_{*},y_{*})=X_{*}(X):=\big{(}x,2u(p)-y\big{)}.$

We notice that

(3.17)
$|X_{*}-P|=\Big{|}\big{(}x,2u(p)-y\big{)}-\big{(}p,u(p)\big{)}\Big{|}=\Big{|}\big{(}x-p,u(p)-y\big{)}\Big{|}=|X-P|.$

As a result,

$\displaystyle\int_{\mathbb{R}^{2}\cap\\{x\in
F_{-}\\}\cap\\{y>2u(p)-u(x)\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)}{|X-P|^{2+s}}\,dX\leqslant\int_{\mathbb{R}^{2}\cap\\{x\in
F_{-}\\}\cap\\{y>2u(p)-u(x)\\}}\frac{dX}{|X-P|^{2+s}}$
$\displaystyle\quad=\int_{\mathbb{R}^{2}\cap\\{x\in
F_{-}\\}\cap\\{y>2u(p)-u(x)\\}}\frac{dX}{|X_{*}(X)-P|^{2+s}}=\int_{\mathbb{R}^{2}\cap\\{x\in
F_{-}\\}\cap\\{y_{*}<u(x_{*})\\}}\frac{dX_{*}}{|X_{*}-P|^{2+s}}$
$\displaystyle\qquad=\int_{\mathbb{R}^{2}\cap\\{x\in
F_{-}\\}}\frac{\chi_{E}(X_{*})}{|X_{*}-P|^{2+s}}\,dX_{*}.$

From this and (3.15), changing the names of the integration variables, we
arrive at

$\displaystyle\int_{\mathbb{R}^{2}\cap\\{x\in
F_{-}\\}}\frac{\chi_{\mathbb{R}^{2}\setminus E}(X)}{|X-P|^{2+s}}\,dX$
$\displaystyle\qquad\leqslant\int_{\mathbb{R}^{2}\cap\\{x\in
F_{+}\\}\cap\\{-u(x)<y<2u(p)+u(x)\\}}\frac{\chi_{E}(X)}{|X-P|^{2+s}}\,dX+\int_{\mathbb{R}^{2}\cap\\{x\in
F_{-}\\}}\frac{\chi_{E}(X)}{|X-P|^{2+s}}\,dX,$

that is

$\int_{\mathbb{R}^{2}\cap\\{x\in F_{-}\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX\leqslant\int_{\mathbb{R}^{2}\cap\\{x\in
F_{+}\\}\cap\\{-u(x)<y<2u(p)+u(x)\\}}\frac{\chi_{E}(X)}{|X-P|^{2+s}}\,dX.$

Therefore,

(3.18) $\begin{split}&\int_{\mathbb{R}^{2}\cap\\{x\in F_{-}\cup
F_{+}\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX\\\
\leqslant&\int_{\mathbb{R}^{2}\cap\\{x\in
F_{+}\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX+\int_{\mathbb{R}^{2}\cap\\{x\in
F_{+}\\}\cap\\{-u(x)<y<2u(p)+u(x)\\}}\frac{\chi_{E}(X)}{|X-P|^{2+s}}\,dX\\\
=&\int_{\mathbb{R}^{2}\cap\\{x\in F_{+}\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)}{|X-P|^{2+s}}\,dX-\int_{\mathbb{R}^{2}\cap\\{x\in
F_{+}\\}\cap\\{y\in(-\infty,-u(x))\cup(2u(p)+u(x),+\infty)\\}}\frac{\chi_{E}(X)}{|X-P|^{2+s}}\,dX.\end{split}$

Now, recalling the notation in (3.16), we note that

(3.19) if $X\in(\mathbb{R}^{2}\setminus E)\cap\\{x\in F_{+}\\}$, then
$X_{*}\in E\cap\\{x_{*}\in F_{+}\\}\cap\\{y_{*}\in(-\infty,-u(x_{*}))\\}$.

Indeed, if $X$ is as above then, since $u(p)<0$,

(3.20) $y_{*}=2u(p)-y<2u(p)-u(x)<-u(x)=-u(x_{*}).$

Hence, since $u(x_{*})=u(x)>-u(p)>0$,

$y_{*}<-u(x_{*})<u(x_{*}).$

By gathering this and (3.20) we obtain the desired result in (3.19).

Thus, by (3.17) and (3.19) we deduce that

$\displaystyle\int_{\mathbb{R}^{2}\cap\\{x\in
F_{+}\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)}{|X-P|^{2+s}}\,dX=\int_{(\mathbb{R}^{2}\setminus E)\cap\\{x\in
F_{+}\\}}\frac{dX}{|X_{*}(X)-P|^{2+s}}$
$\displaystyle\qquad\leqslant\int_{E\cap\\{x_{*}\in
F_{+}\\}\cap\\{y_{*}\in(-\infty,-u(x_{*}))\\}}\frac{dX_{*}}{|X_{*}-P|^{2+s}}.$

We write this inequality in the form

$\displaystyle\int_{\mathbb{R}^{2}\cap\\{x\in
F_{+}\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)}{|X-P|^{2+s}}\,dX\leqslant\int_{\mathbb{R}^{2}\cap\\{x\in
F_{+}\\}\cap\\{y\in(-\infty,-u(x))\\}}\frac{\chi_{E}(X)}{|X-P|^{2+s}}\,dX.$

Comparing with (3.18) we thereby deduce that

(3.21) $\int_{\mathbb{R}^{2}\cap\\{x\in F_{-}\cup
F_{+}\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX\leqslant 0.$

We also remark that

$\mathbb{R}\setminus(F_{+}\cup F_{-})=\big{\\{}x\in\mathbb{R}{\mbox{ s.t.
}}|u(x)|\leqslant-u(p)\big{\\}}.$

From this observation, (3.8) and (3.21) we infer that

(3.22)
$0\leqslant\int_{\mathbb{R}^{2}\cap\\{|u(x)|\leqslant-u(p)\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX.$

It is now useful to remark that

if $|u(x)|\leqslant-u(p)$ and $y<u(p)$, then $X=(x,y)\in E$,

since in this setting $u(x)\geqslant u(p)>y$.

Consequently, the inequality in (3.22) gives that

$\displaystyle 0$ $\displaystyle\leqslant$
$\displaystyle\int_{\mathbb{R}^{2}\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y>u(p)\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX$
$\displaystyle\qquad+\int_{\mathbb{R}^{2}\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y<u(p)\\}}\frac{\chi_{\mathbb{R}^{2}\setminus
E}(X)-\chi_{E}(X)}{|X-P|^{2+s}}\,dX$ $\displaystyle=$
$\displaystyle\int_{\mathbb{R}^{2}\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y>u(p)\\}}\frac{dX}{|X-P|^{2+s}}-2\int_{E\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y>u(p)\\}}\frac{dX}{|X-P|^{2+s}}$
$\displaystyle\qquad-\int_{\mathbb{R}^{2}\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y<u(p)\\}}\frac{dX}{|X-P|^{2+s}}.$

Hence, recalling (3.17) and noticing that, if $y\geqslant u(p)$ then
$u(y_{*})\leqslant u(p)$, we find that

$\displaystyle 0$ $\displaystyle\leqslant$
$\displaystyle\int_{\mathbb{R}^{2}\cap\\{|u(x_{*})|\leqslant-u(p)\\}\cap\\{y_{*}<u(p)\\}}\frac{dX_{*}}{|X_{*}-P|^{2+s}}-2\int_{E\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y>u(p)\\}}\frac{dX}{|X-P|^{2+s}}$
$\displaystyle\qquad-\int_{\mathbb{R}^{2}\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y<u(p)\\}}\frac{dX}{|X-P|^{2+s}}$
$\displaystyle=$
$\displaystyle-2\int_{E\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y>u(p)\\}}\frac{dX}{|X-P|^{2+s}}.$

This entails that the set

(3.23) $E\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y>u(p)\\}$ is of null measure.

On the other hand, we set

$\eta:=\min\left\\{\frac{|u(p)|}{1+8\|u\|_{C^{1}(-d/2,d/2)}},\frac{d}{2}\right\\}$

and we claim that

(3.24)
$\left[-\eta,\eta\right]\times\left[\frac{u(p)}{2},\frac{u(p)}{4}\right]\subseteq
E\cap\\{|u(x)|\leqslant-u(p)\\}\cap\\{y>u(p)\\}.$

Indeed, if
$(x,y)\in\left[-\eta,\eta\right]\times\left[\frac{u(p)}{2},\frac{u(p)}{4}\right]$,
then

$\displaystyle
y-u(x)<\frac{u(p)}{4}-u(x)+u(0)\leqslant\frac{u(p)}{4}+\|u\|_{C^{1}(-d/2,d/2)}|x|$
$\displaystyle\qquad\leqslant\frac{u(p)}{4}+\frac{|u(p)|\,\|u\|_{C^{1}(-d/2,d/2)}}{1+8\|u\|_{C^{1}(-d/2,d/2)}}\leqslant\frac{u(p)}{4}+\frac{|u(p)|}{8}\leqslant\frac{u(p)}{8}<0,$

which gives that $(x,y)\in E$.

Moreover,

$\displaystyle|u(x)|=|u(x)-u(0)|\leqslant\|u\|_{C^{1}(-d/2,d/2)}|x|\leqslant\frac{|u(p)|\,\|u\|_{C^{1}(-d/2,d/2)}}{1+8\|u\|_{C^{1}(-d/2,d/2)}}\leqslant|u(p)|$

and $y\geqslant\frac{u(p)}{2}>u(p)$. These considerations prove (3.24).

From (3.23) and (3.24) we obtain the desired contradiction. The proof of
(1.14) is thereby complete. $\Box$

## References


# Feeder Microgrid Management on an Active Distribution System during a Severe
Outage

Valliappan Muthukaruppan, Ashwin Shirsat, Rongxing Hu, Victor Paduani, Bei Xu,
Yiyan Li, Mesut Baran, Ning Lu, David Lubkeman, and Wenyuan Tang Corresponding
author: Mesut Baran (baran@ncsu.edu)The authors are with Department of
Electrical and Computer Engineering, North Carolina State University, Raleigh,
NC 27695 USA (email: vmuthuk2, ashirsa, rhu5, vdaldeg, bxu8, yli257, nlu2,
dllubkem, wtang8, baran@ncsu.edu). This material is based upon work supported
by U.S. Department of Energy’s Office of Energy Efficiency and Renewable
Energy (EERE) under Solar Energy Technologies Office Award Number DE-
EE0008770.

###### Abstract

Forming a microgrid on a distribution system with large scale outage after a
severe weather event is emerging as a viable solution to improve resiliency at
the distribution level. This option becomes more attractive when the
distribution system has high levels of distributed PV. The management of such
feeder-level microgrid has however many challenges, such as limited resources
that can be deployed on the feeder quickly, and the limited real-time
monitoring and control on the distribution system. Effective use of the
distributed PV is also challenging as they are not monitored and controlled.
To handle these challenges, the paper proposes a 2-stage hierarchical energy
management scheme to securely operate these feeder level micorgrids. The first
stage of the scheme solves a sequential rolling optimization problem to
optimally schedule the main resources (such as a mobile diesel generator and
battery storage unit). The second stage adopts a dispatching scheme for the
main resources to adjust the stage-1 set-points closer to real-time. The
proposed scheme has unique features to assure that the scheme is robust under
highly varying operating conditions with limited system observability: (i) an
innovative PV forecast error adjustment and a dynamic reserve adjustment
scheme to handle the extreme uncertainty on PV power output, and (ii) an
intelligent fuel management scheme to assure that the resources are utilized
optimally over the multiple days of the restoration period. The proposed
algorithm is tested on sample system with real-time data. The results show
that the proposed scheme performs well in maximizing service to loads by
effective use of all the resources and by properly taking into account the
challenging operating conditions.

###### Index Terms:

distribution system restoration, feeder-level microgrid, energy management,
reserve management, forecast error correction.

## Nomenclature

$\overline{E}^{ES}_{i}$

kWh rating of ES $i\in\mathcal{N}^{ES}$

$\left\\{\underline{F},\overline{F}\right\\}_{i}$

Min/Max fuel limits of DG $i\in\mathcal{N}^{DG}$

$\overline{S}^{ES/DG}_{i}$

kVA rating of ES/DG $i\in\mathcal{N}^{ES/DG}$

$\left\\{\underline{SoC},\overline{SoC}\right\\}_{i}$

Min/Max SoC limits of ES $i\in\mathcal{N}^{ES}$

$\alpha_{i},\beta_{i}$

Fuel consumption coefficients of DG

$\gamma_{t}$

Reserve factor for GFM ES unit

$\alpha^{MSD}$

Minimum service duration of load groups

$\alpha^{up}$

Minimum up time for DG

$w_{i}$

Priority weight for load groups

$\theta_{i}$

Power factor angle for DG $i\in\mathcal{N}^{DG}$

$C^{up}$

Start cost for DG

$\lambda_{n}$

Forecast error normalization factor

$F_{f}$

Final desired fuel reserve for DG

$\\{P_{i,t,\phi},Q_{i,t,\phi}\\}^{D}$

Real and Reactive power load demand

$\\{P_{i,t,\phi},Q_{i,t,\phi}\\}^{ES}$

Real and Reactive power output of ES

$\\{P_{i,t},Q_{i,t}\\}^{DG}$

Real and Reactive power output of DG

$P_{i,t,\phi}^{PV}$

Real Power output of BTM PV

$x_{n,t}$

Status of switch connecting $n^{th}$ LG

$y_{i,t}$

Status of switch controlling $i^{th}$ DG

$C_{i,t}$

Start up cost of $i^{th}$ DG

$\\{P_{i,k,\phi},Q_{i,k,\phi}\\}^{D}$

Real and Reactive power load demand

$\\{P_{i,k,\phi},Q_{i,k,\phi}\\}^{ES}$

Real and Reactive power output of ES

$\\{P_{i,k},Q_{i,k}\\}^{DG}$

Real and Reactive power output of DG

$P_{i,k,\phi}^{PV}$

Real Power output of BTM PV

$x_{n,k}$

Status of switch connecting $n^{th}$ LG

$\hat{x}_{n,t}$

Stage-1 schedule of $n^{th}$ switch in $\Delta t$

$\hat{P}_{i,t}^{DG}$

Stage-1 DG schedule in $\Delta t$

$\hat{\epsilon}_{t}$

Forecast error correction factor in $\Delta t$

##  I Introduction

Recent increase in extreme weather events encourage utilities to take measures
to increase the resiliency especially at the distribution level in order to
provide service during such extreme events [1]. One promising technology is
the microgrids as it facilitates the operation of the healthy part of the
system during an extended outage caused by an extreme event [2, 3]. Formation
of a microgrid becomes more attractive especially on a distribution system
with high penetration of renewables, mainly distributed Photovoltaic (PV) [4].

Another promising technology for forming a microgrid on a distribution feeder
on demand is the mobile energy storage (MES) device as it provides the
flexibility to change location depending on changing grid or customer
conditions [5, 6, 7]. Hence, it becomes quite feasible to form a feeder level
microgrid during extended outage (day or more) [8, 9, 10]. If the feeder has
large penetration of distributed PV, then forming such a microgrid becomes
more attractive [5].

Forming a microgrid on a feeder has unique challenges [11, 12]. The main
feature is the mobile generation resources that can be brought to the
location. A typical approach is to bring a MES with a diesel generator (DG).
Both of these resources will have limited capacity and energy that needs to be
carefully managed [13]. Another issue is that the current distribution feeders
have only limited circuit switches for control, which makes it challenging to
ration the load when needed. Having a large amount of PV on the feeder can
help considerably with the microgrid operation. The challenge with utilizing
these resources is that most of these resources are not visible (i.e., not
monitored) and they are quite intermittent. These conditions make it
challenging to operate a feeder-level microgrid during an extended multi-day
outage. This paper focuses on this problem.

Existing literature address only some of the challenges identified. In [14,
15] authors focus on black start of a feeder-level microgrid by picking up
outage load with limited resources and with cold load pickup consideration in
[16, 17]. In [18, 19] load restoration over limited duration is considered and
it uses direct load control to selectively restore critical loads. In [20] a
similar problem is considered with direct load control but behind-the-meter
(BTM) PV is not included.

In [21, 22] load restoration problem over limited duration is considered and a
multi-microgrid formation and reconfiguration is adopted by leveraging
Distributed Energy Resource (DER) flexibility. The DER flexibility is
evaluated and implemented using an aggregator which coordinates with upstream
distribution system operator for restoration. For realizing the flexibility of
DERs and loads, author assumes energy storage (ES) devices collocated with the
PV, thermostatically controllable loads like HVAC can be fully leveraged by
the Distribution System Operator (DSO) at each house level. Such coordination,
flexibility on demand side and concept of DER aggregators is yet to be
implemented in many utilities. This framework also requires significant
communication between the aggregators and the DSO. In [11, 12] realistic
existing feeder infrastructure is considered in restoration but they do not
include the DER uncertainty in the restoration process.

This paper aims at development of a comprehensive microgrid management scheme
for the emerging case of operating a feeder level microgrid on an active
distribution system to provide service during a severe outage. The main
contribution of the proposed method is that all the main issues/challenges
associated with such a microgrid is considered and a robust management scheme
is developed to address these issues. The contributions include:

  * •

Realistic distribution system operating conditions are considered: limited
load and PV visibility and controllability, and limited controllable switches
on the feeder.

  * •

To make best use of the main resources for the microgrid - a DG, and a MES - a
novel DG fuel management scheme is adopted to ensure service to loads during
peak load conditions over multiple days.

  * •

To address the main challenge associated with estimating BTM PV variability, a
new forecast error estimation and adjustment strategy to correct the high
forecast error has been introduced which significantly increases the amount of
PV utilization during restoration.

  * •

A dynamic reserve adjustment strategy for the MES is also introduced in order
to securely operated the microgrid under severe cloud cover events and to
minimize the unintentional load shedding events.

The rest of the paper is organized as follows: sec-II introduces the feeder
microgrid management problem and the proposed scheme, sec-III outlines the new
robustness enhancement strategies proposed to handle the issues associated
with the BTM PV and multi-day restoration using limited resources. Section-IV
illustrates the performance of proposed scheme with a case study which
includes the IEEE 123 node system and operating conditions based on field
data.

##  II Feeder Microgrid Management Scheme

###  II-A Operating Conditions

To ensure realistic operating conditions for feeder level microgrid are
included in the management scheme, the following conditions are considered:

  * •

Utility owned mobile devices MES and DG are quickly deployed at a proper
location on the distribution feeder which has lost power. The MES and DG are
the main resources and are designed to be connected to feeder at substation or
other proper location with necessary infrastructure. The microgrid controller
with the energy management scheme can monitor and control these resources.

  * •

The feeder has high level of BTM PV which will provide supplemental power for
the microgrid. These resources are not visible to the microgrid controller, as
BTM PV is typically not monitored. Hence, there is a need to forecast at least
the net load (actual load minus the power from PV) for the management of the
microgrid.

  * •

The distribution feeder has only a limited remotely controllable
sectionalizing switches. The microgrid controller can use these switches to
adjust the load on the feeder, as these switches divide the loads in to load
groups which can be disconnected as needed. The feeder has also some critical
loads and they should have higher priority in the restoration.

###  II-B Energy Management Scheme

We propose a 2-stage hierarchical energy management framework as presented in
figure-1. Stage-1 is the scheduler to schedule the microgrid resources for the
next period (such as half or one hour) by taking in to account future load and
PV forecast. Stage-2 is the short term dispatching stage which determines the
proper dispatch levels for MES and the DG for the next dispatch period (of 1
to 5 minutes). After every dispatch cycle ($\Delta k$) and scheduling cycle
($\Delta t$) the real-time (RT) measurements are used to update the initial
condition of next dispatch or rolling horizon cycle respectively. These two
stages are outlined below.

Figure 1: Proposed Energy Management Framework

####  II-B1 Stage-1: Scheduling Problem

The stage-1 problem which solves a receding horizon optimization problem on
last day of restoration and a rolling horizon problem on other days of
restoration [23], to optimally allocate the available resources and control
the amount of load to be picked up in each interval through load groups over
the considered horizon which is typically 24 hours to ensure the adequate
availability of resources over the entire restoration period. The important
decision variables in this stage are the switch status $x_{n,t}$ and the
scheduled DG output $P^{DG}_{i,t}$. The Load and PV at each node given by
$P^{D/PV}_{i,t,\phi}$ is the stochastic variable obtained from forecast
information. The objective function is shown in (1) where the first term
maximizes the total expected load groups to be served denoted by $x_{n,t}$
with higher priority to load groups with critical loads (which have higher
weights $w_{n}$) and second term minimizes the start up cost of diesel
generator.

$\max_{x}\quad\sum_{t\in\mathcal{T}}\left[\sum_{i\in\mathcal{N}_{n}}x_{n,t}w_{n}\sum_{\phi\in\Phi}P_{i,t,\phi}^{D}\Delta
t-\sum_{i\in\mathcal{N}^{DG}}C_{i,t}\right]$ (1)

The constraints are defined for $t\in\mathcal{T}$, $\phi\in\Phi$, and
$n\in\mathcal{N}^{LG}$ unless explicitly stated. Equation (2a) highlights the
real power balance in the network over entire time horizon $\mathcal{T}$. The
circuit switches $x_{n,t}$ are the binary decision variables. The reactive
power balance in the network is shows in (2b). It is assumed that all BTM PVs
operate at unity power factor with no reactive power injection.

$\displaystyle\sum_{i\in\mathcal{N}^{ES}}P_{i,t,\phi}^{ES}+\sum_{i\in\mathcal{N}^{DG}}\frac{P_{i,t}^{DG}}{3}$
$\displaystyle=\sum_{i\in\mathcal{N}_{n}}x_{n,t}\left(P_{i,t,\phi}^{D}-P_{i,t,\phi}^{PV}\right)$
(2a)
$\displaystyle\sum_{i\in\mathcal{N}^{ES}}Q_{i,t,\phi}^{ES}+\sum_{i\in\mathcal{N}^{DG}}\frac{Q_{i,t}^{DG}}{3}$
$\displaystyle=\sum_{i\in\mathcal{N}_{n}}x_{n,t}Q_{i,t,\phi}^{D}$ (2b)

Equations (3a)-(3c) represent the real, reactive, and inverter limits of the
energy storage devices $\forall i\in\mathcal{N}^{ES}$. In equation-(3c),
$\gamma_{t}\leq 1$ indicates the reserve factor imposed on the inverter which
reduces the actual limits of the inverter to account for forecast error and
other modeling errors. Further, (3c) is a quadratic constraint which is
linearized using the m-sided polygon technique with m = 6 as proposed in [24].

$\displaystyle-
x_{n,t}\underline{P}_{i}^{ES}\leq\sum_{\phi\in\Phi}P_{i,t,\phi}^{ES}$
$\displaystyle\leq x_{n,t}\overline{P}_{i}^{ES}$ (3a) $\displaystyle
0\leq\sum_{\phi\in\Phi}Q_{i,t,\phi}^{ES}$ $\displaystyle\leq
x_{n,t}\overline{Q}_{i}^{ES}$ (3b)
$\displaystyle\left(\sum_{\phi\in\Phi}P_{i,t,\phi}^{ES}\right)^{2}+\left(\sum_{\phi\in\Phi}Q_{i,t,\phi}^{ES}\right)^{2}$
$\displaystyle\leq x_{n,t}\left(\gamma_{t}\overline{S}^{ES}_{i}\right)^{2}$
(3c)

Equations (4a)-(4b) represent the SoC constraints of the battery $\forall
i\in\mathcal{N}^{ES}$. The temporal change in SoC based on the expected output
of ES is shown in (4a) and the absolute SoC limits is defined in (4b) where
$\eta_{i}$ is the efficiency of $i^{th}$ ES unit.

$\displaystyle SoC_{i,t}$
$\displaystyle=SoC_{i,t-1}-\frac{\sum_{\phi\in\Phi}P_{i,t,\phi}^{ES}}{\overline{E}_{i}^{ES}\eta_{i}^{ES}}\Delta
t$ (4a) $\displaystyle\underline{SoC}_{i}$ $\displaystyle\leq
SoC_{i,t}\leq\overline{SoC}_{i}$ (4b)

Equations (5a)-(5d) represent the diesel generator constraints $\forall
i\in\mathcal{N}^{DG}$. Equations (5a) represents the real power limits of the
DG which is dispatched in fixed power factor mode with power factor angle
$\theta_{i}$ as indicated in (5b). The temporal change of fuel based on
expected dispatch of diesel generator is shown in (5c) and the fuel limits are
imposed using (5d).

$\displaystyle x_{i,t}\underline{S}_{i}^{DG}\cos\theta_{i}$ $\displaystyle\leq
P_{i,t}^{DG}\leq x_{i,t}\overline{S}_{i}^{DG}\cos\theta_{i}$ (5a)
$\displaystyle Q_{i,t}^{DG}$ $\displaystyle=P_{i,t}^{DG}\tan\theta_{i}$ (5b)
$\displaystyle F_{i,t}$
$\displaystyle=F_{i,t-1}-\left(y_{i,t}\alpha_{i}+\beta_{i}P_{i,t}^{DG}\right)\Delta
t$ (5c) $\displaystyle\underline{F}_{i}$ $\displaystyle\leq
F_{i,t}\leq\overline{F}_{i}$ (5d)

Equation (6a) represents the minimum service duration (MSD) constraint
$\forall n\in\mathcal{N}^{LG}$, $t\in\mathcal{T}\setminus\\{1\\}$. When a load
group is decided to be picked up by the microgrid, it has to be supplied for a
minimum duration of $\alpha^{MSD}$. A constraint to ensure load group
connectivity sequence to maintain radiality is also included.

$\displaystyle\sum_{t^{\prime}=t}^{t+\alpha^{MSD}-1}x_{n,t^{\prime}}$
$\displaystyle\geq\alpha^{MSD}(x_{n,t}-x_{n,t-1})$ (6a)

Equation (7a) represents the minimum up time ($\alpha^{up}$) constraint and
(7b), (7c) together compute the startup cost of diesel generators $\forall
i\in\mathcal{N}^{DG}$. Every time the diesel generator is switched on, a start
up cost of $C^{up}$ is incurred.

$\displaystyle\sum_{t^{\prime}=t}^{t+\alpha^{up}-1}y_{i,t^{\prime}}$
$\displaystyle\leq\alpha^{up}(y_{i,t}-y_{i,t-1})$ (7a) $\displaystyle C_{i,t}$
$\displaystyle\geq C^{up}\left(y_{i,t}-y_{i,t-1}\right)$ (7b) $\displaystyle
C_{i,t}$ $\displaystyle\geq 0$ (7c)

####  II-B2 Stage-2: Dispatching Problem

This stage utilizes a short-term forecast which is generally more accurate
than stage-1 forecast for actual dispatching of the resources and switch
control closer to real-time operation. The objective function is shown in (8)
where the first term maximizes the amount of load served with higher priority
to critical loads. Since the switch control is re-evaluated, the second term
minimizes the deviation of load group status from stage-1 results
$\hat{x}_{n,t}$. The third term minimizes the deviation of diesel generator
dispatch from the stage-1 scheduled value $\hat{P}^{DG}_{i,t}$. Here,
$k\in\mathcal{K}_{t}$ is the dispatch interval as indicated in figure-1.

$\max_{x}\quad\sum_{n\in\mathcal{N}^{LG}}\left[w_{n}x_{n,k}\sum_{\phi\in\Phi}P_{i,k,\phi}^{D}-w_{n}^{sw}\left[\hat{x}_{n,t}-x_{n,k}\right]^{2}\right]\\\
-\sum_{i\in\mathcal{N}^{DG}}\left[\hat{P}_{i,t}^{DG}-P_{i,k}^{DG}\right]^{2}$
(8)

The real power balance constraint is similar to stage-1 equation (2a) but
defined over $k\in\mathcal{K}_{t}$. The additional term involving
$\lambda_{n}$ and $\hat{\epsilon}_{t}$ is the forecast correction term which
is explained in section-III-B.

$\sum_{i\in\mathcal{N}^{ES}}P_{i,k,\phi}^{ES}+\sum_{i\in\mathcal{N}^{DG}}\frac{P_{i,k}^{DG}}{3}\\\
=\sum_{n\in\mathcal{N}^{LG}}\left[\sum_{i\in\mathcal{N}_{n}}x_{n,k}\left(P_{i,k,\phi}^{D}-P_{i,k,\phi}^{PV}\right)-\lambda_{n}\frac{\hat{\epsilon}_{t}}{3}\right]$
(9)

The reactive power constraint in (10) is similar to real power constraint.
Note that there is no forecast error correction term here, as the major issue
in forecast comes from BTM PV (as we will see in section-III-B) and there is
no reactive power contribution from PV..

$\sum_{i\in\mathcal{N}^{ES}}Q_{i,k,\phi}^{ES}+\sum_{i\in\mathcal{N}^{DG}}\frac{Q_{i,k}^{DG}}{3}=\sum_{n\in\mathcal{N}^{LG}}\sum_{i\in\mathcal{N}_{n}}x_{n,k}Q_{i,k,\phi}^{D}$
(10)

The constraints (3), (4), and (5) can be imported directly from stage-1
problem by replacing subscripts $t,\Delta t,\mathcal{T}$ by $k,\Delta
k,\mathcal{K}_{t}$ to represent the MES and DG constraints.

##  III Enhancements for Robustness

The performance of the main proposed scheme depends heavily on the accuracy of
the forecast on the net load, and as indicated earlier this forecast can have
significant error. The other challenge is the effective scheduling of the
limited fuel of the DG to ensure service to loads over multiple days. To
address these challenges, and thus to enhance the robustness of the proposed
EMS we propose the following strategies:

  * •

Multi-Day Fuel Rationing

  * •

Learning-based Forecast Correction

  * •

Dynamic Reserve Management

###  III-A Multi-Day Fuel Rationing

For the generation mix of this feeder microgrid we consider a MES, a DG, and
many distributed BTM PVs. MES serves as the grid forming (GFM) resource due to
its fast dynamic response, and DG serves as the main energy source. As stated
earlier in sec-II-B1, we solve a receding horizon problem on final day of
restoration and rolling horizon problem on all other days with a horizon of 24
hours primarily because smart meter data is available real-time at 15-min
interval. Hence, we can forecast the load in a rolling fashion updated every
30-min to 1hr but PV forecast requires weather information which is impossible
to obtain in a rolling fashion by the utility. Also it is very challenging to
accurately forecast PV irradiance beyond 24 hour period. For this reason we
can only solve a rolling horizon problem with fixed horizon of 24 hours
initially. Since the resources need to be rationed over multiple-days we
propose a multi-day fuel rationing framework that ensures availability of DG
over all days of a restoration process even under significant uncertainty in
load and PV. We develop a fuel management algorithm that can efficiently
ration the fuel over multiple days.

$F_{i,r}=\left.\begin{cases}F_{i,t-1}-\frac{(F_{i,t-1}-F_{f})(\mid\mathcal{T}\mid
d+t-1)}{\mid\mathcal{T}\mid\cdot\mid\mathcal{D}\mid},&\text{for
}d<\mid\mathcal{D}\mid\\\ F_{f}&\text{for }d=\mid\mathcal{D}\mid\\\
\end{cases}\right.$ (11) $F_{i,\mid\mathcal{T}\mid}\leq F_{i,r}$ (12)

Equation-(11) highlights a piece-wise linear function defined $\forall
i\in\mathcal{N}^{DG}$, $t\in\mathcal{T}$. $F_{i,r}$ is the minimum fuel
reserve target at the end of each control cycle in stage-1 scheduling problem,
which is enforced as a constraint in (12). $F_{i,t-1}$ is the initial fuel at
start of each control cycle in stage-1 problem, $F_{f}$ is the final reserve
desired at the end of multi-day restoration, which can be set to a non-zero
value to ration some emergency fuel at the end of restoration to supply only
critical loads or to resynchronize the distribution system with main feeder by
picking up all the loads in the system. This framework can be extended to
energy storage devices acting as grid following (GFL) resources by
implementing the same reserve function on SoC of the batteries.

###  III-B Learning-based Forecast Correction

The major issue with BTM PV is lack of real-time data. When individual houses
are net metered, it is easier to forecast the net load but dis-aggregating PV
from load and forecasting just the PV component becomes challenging. For the
EMS scheme, accurate forecast of the net load variation for each load group is
needed. The actual load patterns is usually quite predictable, it is the PV
variability that is challenging to forecast for each load group. We have
adopted an effective PV forecast method proposed in [25] for this purpose. The
method uses the irradiance forecast and the PV rating at each house to
estimate the BTM PV output. The low accuracy on individual house forecast can
become large when aggregated to feeder-level [26]. Handling such large errors
becomes very challenging since we have a MES with limited capacity.

Figure 2: High forecast error in short-term and day-head forecast of
aggregated BTM PV. (top) Real-time net load measurement, stage-1, and stage-2
forecast. (bottom) Total net load forecast error in stage-1 and stage-2
forecast. Zoomed plot shows stage-2 forecast more accurate than stage-1
forecast during no PV duration and error is close to zero.

Figure 2 shows the error in aggregated netload forecast which increases
significantly during the presence of PV (6:00 - 18:00 on day-1 and day-2).
There are two characteristic to the forecast error introduced by the PV. One
is the average forecast error due to error in PV rating assumed, unaccounted
rooftop PV systems and modeling errors. The other characteristic is the
instantaneous random errors introduced by cloud cover events. The errors can
be as high as 2 MW on a 3.5 MW circuit.

To address the average forecast error issue we introduce the learning-based
forecast correction scheme. The scheme uses the battery SOC to estimate the
forecast error in PV. This is based on the following relationship between the
two quantities:

$\Delta P^{D}_{t}\approx\kappa\Delta SoC_{j,t+1}$ (13)

######  Claim III.1.

In a islanded microgrid with battery as a slack generator and when the
proposed rolling optimization framework is used for energy management, then
the difference between the SOC computed at stage-1 and real-time SOC
measurement at each dispatch point given by $\Delta SoC_{j,t+1}$ is
approximately equal to the sum of average forecast error, average power loss
and average modeling error given by $\Delta P^{D}_{t}$ in the dispatch
interval $\Delta t$ with a constant
$\kappa=\frac{\overline{E}^{ES}\eta}{\Delta t}$.

This relationship is derived in the appendix. Hence, using this relationship,
we can use this estimated forecast error component from previous iteration to
adjust the stage-2 forecast error in the next control cycle as shown in (14).
Here, $K$ is the historical window length over which the average is taken.

$\hat{\epsilon}_{t}=\frac{\sum_{i=t-K}^{t-1}\Delta P^{D}_{i}}{K}$ (14)

This average value $\hat{\epsilon}_{t}$ is then normalized to each load group
based on their connected load (given by $\lambda_{n}$) and subtracted from the
stage-2 forecast as shown in (9).

###  III-C Dynamic Reserve Management

Apart from the large average forecast error during PV in fig-2, the other
issue is the large fast variations caused by cloud cover events which is
difficult to forecast even with accurate weather information. Since these are
random events, we cannot modify the forecast values to account for these
errors as we did in sec III-B but instead we can prepare the system to be
robust against such events.

A common method to account for unknown errors in the system while dispatching
is to use a power reserve on the generation resource, which is the MES in this
case. The reserve provides a head room to supply unknown forecast errors and
modeling errors that are not accounted for in the dispatch. This is
incorporated as $\gamma_{t}$ in (3c) where, $\gamma_{t}\in[0,1]$ is the
reserve factor for the MES which is responsible for power balance in the
microgrid. Even though reserves are efficient in operating the microgrid
securely, a large reserve can significantly affect the number of customers
served.

To handle the forecast errors introduced by cloud cover events optimally we
propose a dynamic reserve management strategy. The first step is to estimate
the random errors from available data. We use the previous historical forecast
errors estimates $\Delta P^{D}_{[t-q,t]}$ to estimate the random spikes using
a moving average (MA) model as shown in (15). The parameters of MA model can
be estimated from forecast error estimates obtained from previous few hours.

$\hat{\Delta P^{D}_{t+1}}=\mu+\sum_{i=1}^{q}\theta_{i}\Delta
P_{t-i}^{D}+\Delta P_{t}^{D}$ (15)

The algorithm to estimate the required reserve given in alg-1 is based on
comparing previous interval net load forecast $\hat{P}^{net}_{t}$ given in
(16) to the actual power measurement at PCC given by $P^{net}_{t}$ which is
the actual total net load in the circuit.

$\hat{P}^{net}_{t}=\sum_{i\in\mathcal{N}}\left(P^{D}_{i,t}-P^{PV}_{i,t}\right)$
(16)

High reserve requirements are required during two conditions, when actual load
is more than expected and when actual PV is lower than expected. In both these
conditions net load of microgrid is more than expected. When net load is lower
than expected a high reserve will hinder the optimality of the solution. This
is further validated in our case study in sec IV-A. So it is critical to
dynamically modify the reserve to securely and optimally restore the loads. By
comparing the forecasted net load and the actual measured total net load from
previous control cycle, we can decide the level of reserve requirement as
shown in algorithm 1. This is based on the assumption that net load behavior
in the next control cycle $\Delta t$ will be approximately same as previous
control cycle $\Delta t$ following a trend. When the trend is higher net load
than expected then the reserve is dynamically modified depending on expected
forecast error and left at minimum value $\alpha$ in other scenarios.

Even though this method will address the dynamic forecast errors to some
extent, we need a back up protection like under frequency load trip or other
load shedding algorithm to securely operate the microgrid since it is very
difficult to predict the instantaneous cloud cover events from average power
measurements.

Algorithm 1 Dynamic Reserve Management

0: Net load forecast $\hat{P}^{net}_{t}$, Net load measurements $P^{net}_{t}$,
and Forecast error estimate $\Delta\hat{P}^{D}_{t+1}$

0: Power Reserve Factor $\gamma_{t}$

1: for $t=k$ to $\mid\mathcal{H}\mid$ do

2: if ($P^{net}_{t}>\hat{P}^{net}_{t}$) then

3: $\gamma_{t}=(1-\Delta\hat{P}^{D}_{t+1})$

4: else

5: $\gamma_{t}=1-\alpha$, where $\alpha$ is minimum reserve desired.

6: end if

7: $\gamma_{t}=\mathbb{P}[\gamma_{t}]$ where $\mathbb{P}$ operator thresholds
any input to the set $[\underline{\gamma_{t}},\overline{\gamma_{t}}]$

8: return $\gamma_{t}$

9: end for

##  IV Case Study and Performance Evaluation

To illustrate the performance of the proposed energy management scheme a case
study is conducted. The sample system is the IEEE 123 node system shown in
fig-3. The power to the system is assumed to be unavailable due to outage and
to supply power to the loads during outage, utility brings in a mobile MES and
a DG to form the microgrid. Distribution feeder has 5 controllable switches
and they divide the feeder into 5 load groups as indicated on the figure.
Really high PV penetration of 80-100% is considered in each load zone as shown
in fig-3. The system resources are summarized in table-I and the simulation
parameters are given in table-II.

TABLE I: Location of Resources Resource | Location | Rating (kW) |  | Rating
---
(L/kWh)
| Critical
---
Load Nodes
48, 65, 76 | 210, 140, 245 kW | -
Battery | Substation | 2000 kW | 8000 kWh
| Diesel
---
Generator
Substation | 4000 kW | 10000 L
Rooftop PV | All nodes | 19 - 111 kW | -
Figure 3: Modified IEEE 123 Node Test Feeder TABLE II: Simulation Parameters Parameter | Value | Parameter | Value
---|---|---|---
$\Delta t$, $\Delta k$, $\Delta h$ | 30, 5, 1 min |  $\mathcal{T}$, $\mathcal{K}_{t}$ | 24 hr, 30 min
$w_{1}$,$w_{5}$ | 0.01 | $\alpha^{up}$ | 1 hr
$w_{2}$,$w_{3}$,$w_{4}$ | 0.4, 0.3, 0.2 | $\alpha^{MSD}$ | 2 hr
$\alpha_{i}$, $\beta_{i}$ $\forall i\in\mathcal{N}^{DG}$ | 84.87, 0.20 | $C^{up}$ | 6

The outage is assumed to be for a 48 hour period starting at midnight. The
forecast of load and PV for 4 day period is shown in fig. 4. The forecast are
obtained using the algorithm proposed in [25] trained on pecan street summer
data for load and corresponding weather data from nearby weather station in
Texas for irradiance forecast. As indicated in the figure, stage-1 load
forecast are obtained in a rolling manner every 30-min using smart meter data
while stage-1 irradiance forecast is obtained in a day-ahead stage since there
is no real-time data to update the forecasts. The day ahead forecast is
repeated to obtain a rolling forecast. Stage-2 load and irradiance are short
term forecast obtained closer to real-time. The total installed capacity of PV
at each node is utilized to get the final stage-1 and stage-2 PV forecast from
irradiance information. The following inference can be made from the load and
PV forecast results in fig. 4:

  * •

The forecast error in load is small due to availability of real-time smart
meter data.

  * •

Stage-2 forecast of load is much closer to real-time data compared to stage-1
forecast. While, short-term stage-2 forecast in PV provides no improvement
from stage-1 forecast due to unavailability of real-time information to update
the forecast.

  * •

The PV is consistently over estimated since rating of PV might be smaller than
actual installed capacity which leads to significant error when multiplied
with irradiance forecast. Also, one irradiance profile is used to obtain the
output of all PV units which leads to further error in forecast.

  * •

The cloud cover events are difficult to forecast from weather data and
eventually shows up as instantaneous random forecast errors in net load.

Figure 4: Total Load Forecast (top), and PV Forecast (bottom) against real-
time data for 4 days.

The simulations are carried out using a PC with Intel Core i7-11700F CPU @ 4.8
GHz processor and 64 GB RAM. The prposed schemes are implemented and solved in
Matlab using the Yalmip environment with GUROBI 9.5 solver. The microgrid is
simulated using OpenDSS and it is interfaced with Matlab using the COM
Interface. The simulations are run for 48 hours starting at Day-1 in fig. 4.

To assess the performance of the proposed scheme, the following performance
metrics are used: $P^{CL}_{\%}$ and $P^{NCL}_{\%}$ are the percentage of
critical and non-critical loads (NCLs) served. $P^{PV}_{\%}$ is the percentage
of available PV utilized during restoration. $T^{CL}$ and $T^{NCL}$ are the
average served duration of CL and NCL. $N^{CL}$ is the number of interruptions
in serving critical loads (CLs) and $N^{\mu G}$ is the number of times
microgrid is shutdown. $I^{\mu G}$ is the total duration microgrid was
shutdown due to unavailability of resourcses.

###  IV-A Performance Evaluation

To demonstrate the effectiveness of the forecast error adjustment and dynamic
reserve management strategies, we considered 4 different cases:

  * •

Base Case: Fixed power reserve of 400kW ($\gamma_{t}=0.8$), no forecast error
correction is included to correct the stage-2 forecast, and fixed target fuel
reserve through out two day restoration $F_{i,r}=500L$.

  * •

Case-1: Fixed power reserve of 400kW ($\gamma_{t}=0.8$) through out the two
day restoration period. But, stage-2 forecast is corrected in every stage-1
control cycle based on historical estimation using a window length of 5 hrs.
Fuel reserve management as per eqn-(11).

  * •

Case-2: 20% of current stage-1 net load forecast is kept as reserve in
stage-2, no forecast error correction included to correct the stage-2
forecast, and fuel reserve management as per eqn-(11).

  * •

Case-3: Reserve is adjusted dynamically based on algorithm-1. Stage-2 forecast
is corrected similar to case-2. Fuel reserve management as per eqn-(11).

The 48 hour restoration is simulated using the proposed management scheme and
using the load and PV profiles shown in fig-4. The results are summarized in
table-IV. To evaluate the performance we introduce additional metrics such as
$N^{\mu G}_{UnSch}$ which is the number of unique unscheduled shutdown of
microgrid due to under frequency load shedding which is assumed to last for 30
min and the total duration of such unscheduled shutdown is given by $T^{\mu
G}_{UnSch}$. $T^{\mu G}_{Total}$ is the sum of total scheduled $T^{\mu
G}_{Sch}$ and unscheduled $T^{\mu G}_{UnSch}$ shutdown of the microgrid.

TABLE III: Different cases to evaluate forecast error correction and dynamic reserve management module Cases |  | Error
---
Correction
Fuel Management |  | Reserve
---
$1-\gamma_{t}$
Base Case | No | Fixed $F_{i,r}=500L$ | 400 kW
Case-1 | Yes | As per eqn (11) | 400 kW
Case-2 | No | As per eqn (11) | $0.2\sum(\hat{\boldsymbol{P}}^{D}-\hat{\boldsymbol{P}}^{PV})$
Case-3 | Yes | As per eqn (11) | as per algorithm-1

###  IV-B Fuel Management

Figure-5 shows the actual usage of fuel over a two day period for day-1
simulation between base case with fixed reserve and case-4 with proposed fuel
management and corresponding DG output. The DG is predominantly used to either
charge the battery while serving low load or to completely serve the peak load
in the circuit. From the figure it is evident that the proposed fuel scheme
ensures higher fuel reserve for second day of restoration compared to fixed
fuel reserve scheme. This directly leads to about 7% increase in critical load
served which is highlighted in table-IV. The duration of critical loads served
is also considerably increased. This shows the importance of mulit-day fuel
rationing when only a day-ahead forecast is available.

Figure 5: (top) Comparison of DG dispatch between different fuel management
schemes. (bottom) Comparison of DG Fuel usage over two day restoration with
different fuel management schemes.

From the table we see that critical loads are not served 100%, this doesn’t
mean the loads are under outage during the restoration process. Since this EMS
is from a utility perspective, utility will communicate the hours critical
loads will not be served as scheduled outages from scheduling framework and
the critical loads will use the local generation during these outage hours.
Nevertheless, we want to minimize this duration since the local generation at
critical loads can be limited or unavailable.

###  IV-C Forecast Error Adjustment and Reserve Management

In table-IV, both on day-1 and day-3 we see that Case-2 with forecast error
adjustment performs extremely well across all metrics compared to Case-1
without any error adjustment. There is both significant increase in load
served and PV utilized during restoration which highlights the importance of
forecast error adjustment. Both the cases suffer from unscheduled outages due
to violation of battery limits caused by high forecast error which leads to
about 2 to 3 hours of unscheduled shutdown on day-1. Unscheduled outages are a
inconvenience to the customer, especially to critical loads, and hence it is
desirable to reduce such unscheduled outages during restoration .

Figure 6: (top) Average Net Load Forecast Error Estimation from Battery SoC on Day-1 48 hr restoration. (bottom) Comparison of Stage-2 net load forecast error before and after correction using SoC based average forecast error estimation. TABLE IV: Case study of Energy Management with different robustness enhancement modules over multiple days Metric | Day-1 | Day-3
---|---|---
Base Case | Case-1 | Case-2 | Case-3 | Base Case | Case-1 | Case-2 | Case-3
$P^{CL}_{\%}$ | 72.73% | 78.13% | 76.78% |  79.2% | 71.69% | 75.59% | 72.59% |  76.44%
$P^{NCL}_{\%}$ | 71.1% | 74.63% | 72.88% |  75.244% | 67.76% | 70.55% | 68.45% |  71.38%
$P^{PV}_{\%}$ | 78.33% | 86.06% | 79.63% |  86.92% | 81.88% | 88.27% | 83.68% |  90.35%
$T^{\mu G}_{sch}$ | 2h 40m | 1h 30m | 2h 10m | 2h 30m | 4h 5m | 4h 5m | 3h 30m | 3h 30m
$N^{\mu G}_{sch}$ | 4 | 4 | 4 | 3 | 6 | 5 | 5 | 2
$T^{CL}$ | 37h 15m | 38h 40m | 38h 15m | 39h 20m | 34h 45m | 36h 5m | 35h 25m | 37h 5m
$N^{CL}$ | 8 | 8 | 8 | 9 | 10 | 8 | 9 | 7
$T^{\mu G}_{UnSch}$ | 1h | 30m | 2h | 0 | 1h 30m | 1h | 1h 30m | 30m
$N^{\mu G}_{UnSch}$ | 2 | 1 | 4 | 0 | 3 | 2 | 3 | 1
$T^{\mu G}_{Total}$ | 3h 40m | 2h 10m | 4h 10m | 2h 30m | 5h 35m | 5h 5m | 5h | 4h

The reason for better performance of Case-2 is the forecast error adjustment
strategy. Figure-6 shows the forecast error estimated from SoC measurements
compared against the actual forecasts. This figure highlights that the
estimation strategy is good at estimating the average forecast error but
cannot account for the instantaneous errors. Also, after adjustment the
stage-2 forecast error is close to zero compared to before correction. The
high average error during PV is captured in this case and by adjusting the
forecast we are able to absorb about 8% higher PV on day-1, which eventually
leads to increase in load served and better performance. The percentage
increase on day-3 results is lower because the actual PV available on day-4
being a cloudy day is low.

Case-3 is a simple reserve management strategy where we keep a percentage of
net-load forecast as reserve with an intuition that high reserve is required
at high PV and high load duration which are mutually exclusive. But Case-3
seems to perform much poorer than Case-2 which just has fixed reserve through
out the restoration. Having a high reserve hinders optimality which is visible
in the percentage of load served and PV absorbed but such high reserve in
case-3 must reduce the number of unscheduled outages which it fails to do so
compared to case-2. For this purpose we use the proposed dynamic reserve
management strategy (Case-4) which retains the benefits (day-1) and sometimes
performs better (day-2) than Case-2 in terms of percentage load served and PV
absorbed but the main benefit of case-4 is the reduction in number and
duration of unscheduled outages. It is important to note that proposed
strategy cannot completely eliminate the unscheduled outages because in some
instances the forecast error can be so high that limited resources cannot
handle such errors through reserves and we rely on under-frequency load
shedding during these events.

Figure 7: Dynamic reserve adjustment (top) compared against actual forecast
error (bottom)

The reserves of case-3 and case-4 are shown in fig-7. Case-3 keeps a high
reserve even during low forecast error regions between 18:00-6:00 which
hinders the performance of the algorithm. Also, during cloud cover events with
high spikes the reserve is not high enough which leads to unscheduled outages.
In contrast case-4 keeps the reserve to minimum during low forecast error
regions (no PV) and a high reserve only during significant cloud cover events
that happens during net load zero region which is the main cause of battery
violations since diesel generator is not dispatched during this low net load
conditions. Hence, the algorithm intelligently keeps a high reserve around
9:00 and 16:00 where net load is zero. The other cloud cover events around
12:00 noon is not significant since the actual forecast of PV is higher than
load and keeping high reserve will only hinder optimality.

##  V Conculsion

This paper proposes a two-stage intelligent energy management strategy for a
feeder microgrid to be deployed on a distribution feeder to provide service
during a long outage. To assure that the proposed method has robust
performance, three special schemes have been introduced to addresses the
unique challenges associated with management of these microgrids: The multi-
day fuel rationing scheme assures that the load in the circuit is served over
multiple days of the restoration which increases the overall service to
critical loads and improves PV utilization. The novel state-of-charge based
forecast error estimation strategy which learns the forecast error online and
adjusts the near real-time forecast significantly increases the BTM PV
utilization during the restoration. The dynamic reserve management strategy
introduced to handle sudden cloud cover events that can lead to unscheduled
outages of the microgrid. A case study on a sample system illustrates the
effectiveness of the proposed scheme, as it shows that the proposed error
correction based reserve adjustment is quite effective even when there is high
variability in PV and considerably reduces the number of unscheduled outages
during the restoration process.

## References

  * [1] “Resilience of the U.S. Electricity System: A Multi-Hazard Perspective,” https://www.energy.gov/policy/downloads/resilience-us-electricity-system-multi-hazard-perspective.
  * [2] A. Hussain, V.-H. Bui, and H.-M. Kim, “Microgrids as a resilience resource and strategies used by microgrids for enhancing resilience,” _Applied Energy_ , vol. 240, pp. 56–72, Apr. 2019.
  * [3] A. Hirsch, Y. Parag, and J. Guerrero, “Microgrids: A review of technologies, key drivers, and outstanding issues,” _Renewable and Sustainable Energy Reviews_ , vol. 90, pp. 402–411, Jul. 2018.
  * [4] S. Booth, “Microgrid-Ready Solar PV - Planning for Resiliency,” _NREL/FS-7A40-70122_ , p. 2, 2017.
  * [5] “Mobile Energy Storage Study — Mass.gov,” https://www.mass.gov/service-details/mobile-energy-storage-study.
  * [6] H. Abdeltawab and Y. A. I. Mohamed, “Mobile energy storage sizing and allocation for multi-services in power distribution systems,” _IEEE Access_ , vol. 7, pp. 176 613–176 623, 2019.
  * [7] J. Kim and Y. Dvorkin, “Enhancing distribution system resilience with mobile energy storage and microgrids,” _IEEE Transactions on Smart Grid_ , 2019\.
  * [8] L. Che, M. Shahidehpour, and M. Shahidehpour, “Adaptive formation of microgrids with mobile emergency resources for critical service restoration in extreme conditions,” _IEEE Transactions on Power Systems_ , 2019.
  * [9] S. Lei, J. Wang, C. Chen, and Y. Hou, “Mobile emergency generator pre-positioning and real-time allocation for resilient response to natural disasters,” _IEEE Transactions on Smart Grid_ , 2018.
  * [10] S. Lei, C. Chen, H. Zhou, and Y. Hou, “Routing and scheduling of mobile power sources for distribution system resilience enhancement,” _IEEE Transactions on Smart Grid_ , 2019.
  * [11] N. C. Koutsoukis, P. S. Georgilakis, and N. D. Hatziargyriou, “Service restoration of active distribution systems with increasing penetration of renewable distributed generation,” _IET Generation, Transmission & Distribution_, vol. 13, no. 14, pp. 3177–3187, 2019.
  * [12] M.-G. Choi, J.-H. Choi, S.-Y. Yun, and S.-J. Ahn, “MILP-Based Service Restoration Method Utilizing Both Existing Infrastructure and DERs in Active Distribution Networks,” _IEEE Access_ , vol. 10, pp. 36 477–36 489, 2022.
  * [13] J. Dugan, S. Mohagheghi, and B. Kroposki, “Application of Mobile Energy Storage for Enhancing Power Grid Resilience: A Review,” _Energies_ , vol. 14, no. 20, p. 6476, Oct. 2021.
  * [14] C. Chen, J. Wang, F. Qiu, and D. Zhao, “Resilient Distribution System by Microgrids Formation After Natural Disasters,” _IEEE Transactions on Smart Grid_ , vol. 7, no. 2, pp. 958–966, Mar. 2016.
  * [15] S. Meng, R. Roofegari nejad, and W. Sun, “Robust Distribution System Load Restoration with Time-Dependent Cold Load Pickup,” _IEEE Transactions on Power Systems_ , pp. 1–1, 2020.
  * [16] C. Chen, J. Wang, and D. Ton, “Modernizing Distribution System Restoration to Achieve Grid Resiliency Against Extreme Weather Events: An Integrated Solution,” _Proceedings of the IEEE_ , vol. 105, no. 7, pp. 1267–1288, Jul. 2017.
  * [17] W. Yuan, J. Wang, F. Qiu, C. Chen, C. Kang, and B. Zeng, “Robust Optimization-Based Resilient Distribution Network Planning Against Natural Disasters,” _IEEE Transactions on Smart Grid_ , vol. 7, no. 6, pp. 2817–2826, Nov. 2016.
  * [18] H. Momen, A. Abessi, S. Jadid, M. Shafie-khah, and J. P. S. Catalão, “Load restoration and energy management of a microgrid with distributed energy resources and electric vehicles participation under a two-stage stochastic framework,” _International Journal of Electrical Power & Energy Systems_, vol. 133, p. 107320, Dec. 2021.
  * [19] S. Poudel and A. Dubey, “Critical Load Restoration Using Distributed Energy Resources for Resilient Power Distribution System,” _IEEE Transactions on Power Systems_ , vol. 34, no. 1, pp. 52–63, Jan. 2019.
  * [20] M. R. Kleinberg, K. Miu, and H.-D. Chiang, “Improving Service Restoration of Power Distribution Systems Through Load Curtailment of In-Service Customers,” _IEEE Transactions on Power Systems_ , vol. 26, no. 3, pp. 1110–1117, Aug. 2011.
  * [21] W. Liu and F. Ding, “Collaborative Distribution System Restoration Planning and Real-Time Dispatch Considering Behind-the-Meter DERS,” _IEEE Transactions on Power Systems_ , vol. 36, no. 4, pp. 3629–3644, Jul. 2021.
  * [22] ——, “Hierarchical Distribution System Adaptive Restoration With Diverse Distributed Energy Resources,” _IEEE Transactions on Sustainable Energy_ , vol. 12, no. 2, pp. 1347–1359, Apr. 2021.
  * [23] L. Glomb, F. Liers, and F. Rösel, “A rolling-horizon approach for multi-period optimization,” _European Journal of Operational Research_ , vol. 300, no. 1, pp. 189–206, Jul. 2022.
  * [24] H. Ahmadi and J. R. Martı´, “Linear Current Flow Equations With Application to Distribution Systems Reconfiguration,” _IEEE Transactions on Power Systems_ , vol. 30, no. 4, pp. 2073–2080, Jul. 2015.
  * [25] Y. Li, S. Zhang, R. Hu, and N. Lu, “A meta-learning based distribution system load forecasting model selection framework,” _Applied Energy_ , vol. 294, p. 116991, Jul. 2021.
  * [26] B. C. Erdener, C. Feng, K. Doubleday, A. Florita, and B.-M. Hodge, “A review of behind-the-meter solar forecasting,” _Renewable and Sustainable Energy Reviews_ , vol. 160, p. 112224, May 2022.

## Appendix

###  V-A Proof of Claim-III.1

###### Proof.

The estimated SoC from stage-1 problem $\forall j\in\mathcal{N}^{ES}_{GFM}$ is
given by (17) where, $SoC_{j,t}$ is the actual SoC of battery obtained from
previous feedback from system, $\hat{P}^{D}$ is the forecasted net load demand
in the interval $[t,t+\Delta t]$.

$\hat{SoC}_{j,t+\Delta t}=SoC_{j,t}+\\\
\frac{\hat{P}^{D}_{t}-\sum_{i\in\mathcal{N}^{DG}}P_{i,t}^{DG}-\sum_{i\in\mathcal{N}^{ES}_{GFL}}P_{i,t}^{ES}}{\overline{E}_{j}^{ES}\eta_{j}}\Delta
t$ (17)

The actual measured SoC of battery $\forall j\in\mathcal{N}^{ES}_{GFM}$ can be
formulated as shown in (18).

$SoC_{j,t+\Delta t}=SoC_{j,t}-\\\ \frac{\int_{t}^{t+\Delta
t}\left[P^{total}(t)-\sum_{i\in\mathcal{N}^{DG}}P_{i}^{DG}(t)-\sum_{i\in\mathcal{N}^{ES}_{GFL}}P_{i}^{ES}(t)\right]dt}{\overline{E}_{j}^{ES}\eta_{j}}$
(18)

Where, $P^{total}(t)$ is total net demand in the system with power loss.
Difference in SoC between measured and stage-1 estimation given by $\Delta
SoC_{j,t+\Delta t}$ can be simplified as shown in eqn-19.

$\Delta SoC_{j,t+\Delta t}=SoC_{t+\Delta t}-\hat{SoC}_{j,t+\Delta t}$ $\Delta
SoC_{j,t+\Delta t}=\left[\frac{\hat{P}^{D}\Delta t-\int_{t}^{t+\Delta
t}P^{total}(t)dt}{\overline{E}_{j}^{ES}\eta_{j}}\right]\\\
-\left[\frac{\sum_{i\in\mathcal{N}^{DG}\cup\mathcal{N}^{ES}_{GFL}}\left(\int_{t}^{t+1}P_{i}(t)dt-\hat{P}_{i,t}\Delta
t\right)}{\overline{E}_{j}^{ES}\eta_{j}}\right]$ (19)

Now, $\int_{t}^{t+\Delta t}P^{total}(t)dt$ can be approximated by the average
value $\overline{P}^{total}\Delta t$, where $\overline{P}^{total}$ is the
average demand value in the interval $[t,t+\Delta t]$. Further, the grid-
following resources dispatch in stage-2 over the horizon $\mathcal{K}$,
$P^{G}_{i,k}$ is kept close to the scheduling value $\hat{P}^{G}_{i,t}$
$\forall i\in\mathcal{N}^{G}_{GFL},k\in\mathcal{K},t\in\mathcal{T}$ through
stage-2 objective as shown in (8). This ensures the stage-2 dispatch is close
to scheduled value ($P_{i}(t)\approx\hat{P}_{i,t}$) which cancels out second
term in (19).

$\displaystyle\Delta SoC_{j,t+\Delta t}\overline{E}^{ES}_{j}\eta_{j}$
$\displaystyle\approx(\hat{P}^{D}-\overline{P}^{total})\Delta
t-0\quad\text{since,}\quad P_{i}(t)\approx\hat{P}_{i,t}$ (20)
$\displaystyle\text{let, }\Delta P^{D}_{t}$
$\displaystyle\approx(\hat{P}^{D}-\overline{P}^{total})$ $\displaystyle\Delta
P^{D}_{t}$ $\displaystyle\approx\frac{\Delta
SoC_{j,t+1}\overline{E}^{ES}_{j}\eta_{j}}{\Delta t}$ $\displaystyle\Delta
P^{D}_{t}$ $\displaystyle\approx\kappa\Delta SoC_{j,t+1}$

Therefore, the average forecast and modeling error in the interval
$[t,t+\Delta t]$ is approximately equal to the difference in SoC computed from
feedback measurements and stage-1 results with a constant factor
$\kappa=\frac{\overline{E}^{ES}_{j}\eta_{j}}{\Delta t}$. ∎