text:
\begin{document}
\allowdisplaybreaks
\def\R{\mathbb R}  \def\ff{\frac} \def\ss{\sqrt} \def\B{\mathbf
B} \def\W{\mathbb W}
\def\N{\mathbb N} \def\kk{\kappa} \def\m{{\bf m}}
\def\ee{\varepsilon}\def\ddd{D^*}
\def\dd{\delta} \def\DD{\Delta} \def\vv{\varepsilon} \def\rr{\rho}
\def\<{\langle} \def\>{\rangle} \def\GG{\Gamma} \def\gg{\gamma}
  \def\nn{\nabla} \def\pp{\partial} \def\E{\mathbb E}
\def\d{\text{\rm{d}}} \def\bb{\beta} \def\aa{\alpha} \def\D{\scr D}
  \def\si{\sigma} \def\ess{\text{\rm{ess}}}
\def\beg{\begin} \def\beq{\begin{equation}}  \def\F{\scr F}
\def\Ric{\text{\rm{Ric}}} \def\Hess{\text{\rm{Hess}}}
\def\e{\text{\rm{e}}} \def\ua{\underline a} \def\OO{\Omega}  \def\oo{\omega}
 \def\tt{\tilde} \def\Ric{\text{\rm{Ric}}}
\def\cut{\text{\rm{cut}}} \def\P{\mathbb P} \def\ifn{I_n(f^{\bigotimes n})}
\def\C{\scr C}      \def\aaa{\mathbf{r}}     \def\r{r}
\def\gap{\text{\rm{gap}}} \def\prr{\pi_{{\bf m},\varrho}}  \def\r{\mathbf r}
\def\Z{\mathbb Z} \def\vrr{\varrho} \def\ll{\lambda}
\def\L{\scr L}\def\Tt{\tt} \def\TT{\tt}\def\II{\mathbb I}
\def\i{{\rm in}}\def\Sect{{\rm Sect}}  \def\H{\mathbb H}
\def\M{\scr M}\def\Q{\mathbb Q} \def\texto{\text{o}} \def\LL{\Lambda}
\def\Rank{{\rm Rank}} \def\B{\scr B} \def\i{{\rm i}} \def\HR{\hat{\R}^d}
\def\to{\rightarrow}\def\l{\ell}\def\iint{\int}
\def\EE{\scr E}\def\Cut{{\rm Cut}}
\def\A{\scr A} \def\Lip{{\rm Lip}}
\def\BB{\scr B}\def\Ent{{\rm Ent}}\def\L{\scr L}
\def\R{\mathbb R}  \def\ff{\frac} \def\ss{\sqrt} \def\B{\mathbf
B}
\def\N{\mathbb N} \def\kk{\kappa} \def\m{{\bf m}}
\def\dd{\delta} \def\DD{\Delta} \def\vv{\varepsilon} \def\rr{\rho}
\def\<{\langle} \def\>{\rangle} \def\GG{\Gamma} \def\gg{\gamma}
  \def\nn{\nabla} \def\pp{\partial} \def\E{\mathbb E}
\def\d{\text{\rm{d}}} \def\bb{\beta} \def\aa{\alpha} \def\D{\scr D}
  \def\si{\sigma} \def\ess{\text{\rm{ess}}}
\def\beg{\begin} \def\beq{\begin{equation}}  \def\F{\scr F}
\def\Ric{\text{\rm{Ric}}} \def\Hess{\text{\rm{Hess}}}
\def\e{\text{\rm{e}}} \def\ua{\underline a} \def\OO{\Omega}  \def\oo{\omega}
 \def\tt{\tilde} \def\Ric{\text{\rm{Ric}}}
\def\cut{\text{\rm{cut}}} \def\P{\mathbb P} \def\ifn{I_n(f^{\bigotimes n})}
\def\C{\scr C}      \def\aaa{\mathbf{r}}     \def\r{r}
\def\gap{\text{\rm{gap}}} \def\prr{\pi_{{\bf m},\varrho}}  \def\r{\mathbf r}
\def\Z{\mathbb Z} \def\vrr{\varrho} \def\ll{\lambda}
\def\L{\scr L}\def\Tt{\tt} \def\TT{\tt}\def\II{\mathbb I}
\def\i{{\rm in}}\def\Sect{{\rm Sect}}  \def\H{\mathbb H}
\def\M{\scr M}\def\Q{\mathbb Q} \def\texto{\text{o}} \def\LL{\Lambda}
\def\Rank{{\rm Rank}} \def\B{\scr B} \def\i{{\rm i}} \def\HR{\hat{\R}^d}
\def\to{\rightarrow}\def\l{\ell}
\def\8{\infty}\def\I{1}\def\U{\scr U}
\maketitle

\begin{abstract} By using the heat kernel parameter expansion   with respect to the frozen SDEs, the intrinsic derivative is estimated for the law of
 Mckean-Vlasov SDEs with respect to the initial distribution.  As an application,
 the total variation distance between the laws of two solutions  is bounded by the Wasserstein distance for initial distributions. These extend some recent  results
 proved  for  distribution-free noise     by using  the coupling  method  and Malliavin calculus.
 \end{abstract} \noindent
 AMS subject Classification:\  60H1075, 60G44.   \\
\noindent
 Keywords: Mckean-Vlasov SDEs, intrinsic derivative, $L$-derivative,  heat kernel parameter expansion.
 \vskip 2cm

\section{Introduction}
Let $\scr P_2$ be the set of all probability measures on $\R^d$ with finite second moment, which is called the Wasserstein space under the metric
$$\W_2(\mu,\nu):= \inf_{\pi\in \C(\mu,\nu)} \bigg(\int_{\R^d\times\R^d} |x-y|^2 \pi(\d x,\d y)\bigg)^{\ff 1 {2}},\ \ \mu,\nu\in \scr P_{2},$$ where $\C(\mu,\nu)$ is the set of all couplings of $\mu$ and $\nu$.   Consider the following distribution dependent SDE on $\R^d$:
\beq\label{E1} \d X_t^\mu= b_t(X_t^\mu , \L_{X_t^\mu})\d t+\si_t(X_t^\mu,\L_{X_t^\mu})\d W_{t},\ \ \L_{X_0^\mu}=\mu\in\scr P_2,\end{equation}
where $W_t$ is an $m$-dimensional Brownian motion on a complete filtration probability space $(\OO,\{\F_t\}_{t\ge 0},\P)$, $\L_{X_t}$ is the law of $X_t$, and
\begin{align*} &b: \R_+\times\R^d\times \scr P_2\to \R^d,\ \ \si: \R_+\times \R^d\times \scr P_2\to \R^d\otimes\R^m
 \end{align*}
 are measurable.  This type equations, known as  Mckean-Vlasov or mean field SDEs, have been intensively investigated and applied, see for instance the monograph \cite{CDLL19}
 and references therein.

 To characterize  the regularity of   the law $\L_{X_t^\mu}$ with respect to the initial distribution $\mu$,  we  investigate  the derivative estimate of the functions
 $$\scr P_2\ni\mu\mapsto P_tf(\mu):= \E f(X_t^\mu),\ \ f\in \B_b(\R^d), t>0.$$
 When  the noise coefficient $\si_t(x,\mu)$ does not depend on $\mu$, the Harnack inequality and derivative formula have been established in \cite{FYW1, RW}
 for   $P_tf $    by using the coupling by change of measures and Malliavin calculus respectively. See also \cite{BRW20, HRW, HW, Song} for extensions to distribution-path dependent SDEs/SPDEs,  singular distribution dependent SDEs, and  distribution dependent SDEs with jumps, where in \cite{Song} allows the noise to be also  distribution dependent and establishes the gradient estimate on $P_t f(x):=(P_tf)(\dd_x)$ when the initial distribution is a Dirac measure.  In this paper, we   estimate the derivative
 of $P_tf(\mu)$  in $\mu$ by using the heat kernel parameter expansion  with respect to the frozen SDE
 \beq\label{E1'} \d X_t^{z,\mu} = b_t(z,  \mu_t)\d t+\si_t(z,\mu_t)\d W_{t}\end{equation}
 for fixed $(z,\mu)\in\R^d\times \scr P_2(\R^d),$ where  $\mu_t:= \L_{X_t^\mu}.$  Since this SDE has constant coefficients, the solution has a Gaussian heat kernel which can be easily analyzed.

Before introducing the main result, we first recall the intrinsic derivative and $L$-derivative for functions on $\scr P_2$ which go back to \cite{AKR} where the intrinsic derivative on the configurations space is introduced, see \cite{RW19} for the link of different derivatives for measures.

\begin{defn}\label{defLio} Let $f:\scr P_2\to \mathbb{R}$ and $g:\mathbb{R}^d\times \scr P_2\to \mathbb{R}$.
\begin{enumerate}
\item[(1)] $f$ is called intrinsically differentiable, if for any $\mu\in \scr P_2$,
$$L^2(\R^d\to\R^d;\mu)\ni\phi\mapsto D_\phi^Lf(\mu):=\lim_{\vv\downarrow 0} \ff{f(\mu\circ({\rm Id}+\vv\phi)^{-1})- f(\mu)}\vv\in\R$$
is a well defined bounded linear functional. In this case,  the    unique map
$$\scr P_2\ni\mu\mapsto D^Lf(\mu)\in L^2(\R^d\to\R^d;\mu)$$
such that $D_\phi^Lf(\mu)= \<\phi, D^Lf(\mu)\>_{L^2(\mu)}$ holds for any $\mu\in \scr P_2$  and $\phi\in L^2(\R^d\to\R^d;\mu)$ is called the intrinsic derivative of $f$, and we denote
$$\|D^Lf(\mu)\|:= \|D^Lf(\mu)(\cdot)\|_{L^2(\mu)},\ \ \mu\in \scr P_2.$$
If moreover
$$\lim_{\mu(|\phi|^2)\to 0}\frac{f(\mu\circ({\rm Id}+ \phi)^{-1})-f(\mu)-D_\phi^Lf(\mu)}{\sqrt{\mu(|\phi|^2)}}=0,\ \ \mu\in \scr P_2,$$ we call $f$ $L$-differentiable, and in this case $D^Lf$ is also called the $L$-derivative of $f$.
\item[(2)] We denote $f\in C^1(\scr P_2)$, if $f$ is $L$-differentiable and  its $L$-derivative has a version $D^L f(\mu)(x)$ jointly continuous in $(x,\mu)\in \R^d\times \scr P_2.$

 \item[(3)] We denote $g\in C^{1,1}(\R^d\times\scr P_2)$, if $g(x,\cdot)\in C^1(\scr P_2)$ for $x\in \R^d$, $g(\cdot,\mu)\in C^1(\R^d)$ for $\mu\in \scr P_2$,
 $g(x,\mu), \nn g(\cdot,\mu)(x)$ are jointly continuous in $(x,\mu)\in \R^d\times\scr P_2$,  and $D^Lg(x,\cdot)(\mu)(y)$  has a version jointly continuous in $(x,y,\mu)\in \R^d\times\R^d\times\scr P_2.$
 \item[(4)] A vector- or matrix-valued function is said in a class defined above, if so are its component functions.
 \end{enumerate}
\end{defn}

To estimate the intrinsic derivative of $P_tf(\mu)$, we need the following condition. Let $|\cdot|$ and $\|\cdot\|$ denote the norm in $\R^d$ and the operator norm for linear operators repsectively.

 \beg{enumerate} \item[{\bf (H)}] For any $t\ge 0$, $b_t,\si_t\in C^{1,1}(\R^d\times\scr P_2),$  and  there exists an increasing function $K: [0,\infty)\to [0,\infty)$ such that
 for any  $ t\ge 0, x,y\in\R^d$ and $ \mu\in \scr P_2(\R^d),$
$$  K_t^{-1}  {\rm Id} \le  (\si_t\si_t^*) (x,\mu)\le  K_t {\rm Id},$$
\beg{align*}  &\ |b_t(x,\mu)|+ \|\nn b_t(\cdot,\mu)(x)\|+ \|D^L \{b_t(x,\cdot)\}(\mu)\| \\
&\qquad  +  \|\nn \{\si_t(\cdot,\mu)\}(x)\|^2+  \|D^L\{\si_t(x,\cdot)\}(\mu)\|^2
  \le K_t,\end{align*}
\beg{align*}   &\|D^L \{b_t(x,\cdot)\}(\mu)- D^L \{b_t(y,\cdot)\}(\mu)\|+ \|D^L \{\si_t(x,\cdot)\}(\mu)- D^L\{ \si_t(y,\cdot)\}(\mu)\|\\
&\le K_t |x-y|^2.
  \end{align*}
\end{enumerate}

 It is well known that    SDE \eqref{E1} is well-posed under the assumption {\bf (H)}, so that $P_tf$ is  well defined on $\scr P_2$ for any $t\ge 0$ and $f\in \B_b(\R^d).$
 In general, for any $s\ge 0$ and $X_{s,s}^\mu\in L^2(\OO\to\R^d, \F_s,\P)$ with $\L_{X_{s,s}^\mu} =\mu$, let $X_{s,t}^\mu$ be the unique solution of \eqref{E1} for $t\ge s$:
 \beq\label{E11} \d X_{s,t}^\mu= b_t(X_{s,t}^\mu , \L_{X_{s,t}^\mu})\d t+\si_t(X_{s,t}^\mu,\L_{X_{s,t}^\mu})\d W_{t},\ \ t\ge s, \L_{X_{s,s}^\mu}=\mu\in\scr P_2.\end{equation}
 We denote  $P_{s,t}^*\mu=\L_{X_{s,t}^\mu}$ and investigate the regularity of
 $$P_{s,t}f(\mu):= \E f(X_{s,t}^\mu)= \int_{\R^d} f\d (P_{s,t}^*\mu),\ \ f\in \B_b(\R^d).$$
 By the uniqueness, we have the flow property
 $$P_{s,t}^*= P_{r,t}^*P_{s,r}^*,\ \ 0\le s\le r\le t.$$ However,  due to the distribution dependence,  $P_{s,t}$ is no-longer a semigroup, i.e. in general
 $P_{s,t}\ne P_{r,t}P_{s,r}$  and
 $$P_tf(\mu)\ne \int_{\R^d} P_tf(x)\mu(\d x),$$ so that the regularity of $P_tf(\mu)$ in $\mu\in\scr P_2$ can not be  deduced from that of $P_tf(x):= P_tf(\dd_x)$ for $x\in \R^d$, see for instance \cite{FYW1} for details.

We now state  the main result of the paper as follows.

 \beg{thm}\label{T1} Assume {\bf (H)}.  Then for any $t>s$ and $f\in \B_b(\R^d)$, $P_{s,t}f$ is
  $L$-differentiable, and there exists an increasing function $C: [0,\infty)\to (0,\infty)$ such that
\beq\label{EST1} \|D^L P_{s,t} f(\mu)\|\le \ff{C_t \|f\|_\infty}{\ss {t-s}},\ \ t>s\ge 0, f\in \B_b(\R^d).\end{equation}
 Consequently, for any $t>s\ge 0,   \mu,\nu\in \scr P_2,$
\beq\label{EST2} \|P_{s,t}^*\mu- P_{s,t}^*\nu\|_{var}:= \sup_{\|f\|_\infty\le 1} |P_{s,t}f(\mu)-P_{s,t}f(\nu)|\le \ff{C_t\|f\|_\infty}{\ss {t-s}} \W_2(\mu,\nu).\end{equation}
\end{thm}

\paragraph{Remark 1.1.} We may also apply Malliavin calculus to establish a derivative formula for $D^L P_{s,t} f(\mu)$  as in \cite{Song}, where the usual derivative in initial points
(rather than in initial distributions) are studied. However, in this way we need stronger conditions on the coefficients,  i.e.    $b_t(x,\mu)$ and $\si_t(x,\mu)$ also have bounded second order derivatives in $x$. Let us explain this in more details.

Firstly, under  {\bf (H)},  the Malliavin matrix
  $$M_{s,t}:= \big\{\<D (X_{s,t}^\mu)_i, D (X_{s,t}^\mu)_j\>_\H\big\}_{1\le i\le j}$$
  is invertible with $\E\|M_{s,t}^{-1}\|^2<\infty$ for $t>s\ge 0$, where $D$ is the Malliavin gradient,   $\H$ is the Cameron-Martin space in Malliavin calculus, and $(X_{s,t}^\mu)_i$ is   the $i$-th component of $X_{s,t}^\mu$.

  Next, for any $\phi\in L^2(\R^d\to\R^d;\mu)$, let $v_{s,t}^\phi= D_\phi^L X_{s,t}^\mu$, which exists  in $L^2(\P)$ and satisfies
  $$\E|v_{s,t}^\phi|^2\le c(t) \mu(|\phi|^2)$$ for some constant $c(t)>0$, see \cite[Proposition 3.2]{RW}.

  Then for any $f\in C_b^1(\R^d)$, by the chain rule and the integration by parts formula for
  the Malliavin gradient $D$, we have
  \beg{align*} & D_\phi^LP_{s,t}f(\mu)= \E\<\nn f(X_{s,t}^\mu), v_{s,t}^\phi\> = \sum_{i=1}^d \E[\pp_i f(X_{s,t}^\mu) (v_{s,t}^\phi)_i] \\
  &= \sum_{i,j,k=1}^d\E\big[ \pp_i f(X_{s,t}^\mu) (M_{s,t})_{ij} (M_{s,t}^{-1})_{jk} (v_{s,t}^\phi)_k\big]\\
  &=  \sum_{i,j,k=1}^d\E\big[\<D f(X_{s,t}^\mu), D(X_{s,t}^\mu)_j\>_\H    (M_{s,t}^{-1})_{jk} (v_{s,t}^\phi)_k\big]\\
  &= \sum_{i,j,k=1}^d \E\big[f(X_{s,t}^\mu) D^*\big\{ (M_{s,t}^{-1})_{jk} (v_{s,t}^\phi)_k D(X_{s,t}^\mu)_j\big\}\big],\end{align*}
  where $D^*$ is the Malliavin divergence. To make the above calculations meaningful, we need to verify that $(M_{s,t}^{-1})_{jk} (v_{s,t}^\phi)_k D(X_{s,t}^\mu)_j$ belongs to the domain of   $D^*$, for which the second order derivatives of coefficients will be involved. For instance, as shown in \cite[Proposition 3.2]{RW} that $v_{s,t}^\phi$ solves an SDE
  involving in the first order derivatives of $b$ and $\si$,  making Malliavin derivative to this SDE we see that $D v_{s,t}^\phi$  solves an SDE containing  the second order derivatives of coefficients.

  \

The remainder of the paper is organized as follows. In Section 2,  we formulate $P_{s,t}f(\mu)$ using classical SDEs with parameter $\mu$ and the parameter expansion of heat kernels with respect to the  frozen SDE \eqref{E1'},  and estimate the $L$-derivative for functions of $P_{s,t}^*\mu$.   With these preparations, we prove Theorem \ref{T1} in Section 3.

 \section{Preparations}

 We first represent $P_{s,t}f(\mu)$ by using a Markov semigroup $P_{s,t}^\mu$ with parameter $\mu$, then introduce the heat kernel expansion of $P_{s,t}^\mu$ with respect to the frozen SDEs. Since the frozen SDE has explicit Gaussian heat kernel, this enables us to calculate the intrinsic derivative of $P_tf(\mu)$ with respect to $\mu$.

\subsection{A representation of $P_{s,t}$}

   For any $s\ge 0, x\in\R^d$ and $\mu\in \scr P_2$, consider the decoupled SDE
   \beq\label{ES0} \d X_{s,t}^{x,\mu}= b_t(X_{s,t}^{x,\mu}, P_{s,t}^*\mu)\d t + \sigma_t(X_{s,t}^{x,\mu}, P_{s,t}^*\mu)\d W_t, \ \ X_{s,s}^{x,\mu}=x, t\ge s.\end{equation}
   In this SDE,  the measure variable $P_{s,t}^*\mu$ is fixed, so that it reduces to the classical time inhomogeneous SDE.  Let $P_{s,t}^\mu$ be the associated Markov semigroup, i.e.
   $$P_{s,t}^\mu f(x)= \E f(X_{s,t}^{x,\mu}),\ \ t\ge s, f\in \B_b(\R^d), x\in \R^d.$$
Since  $X_{s,t}^\mu$ solves \eqref{ES0} with the random initial value    $X_{s,s}^\mu$ replacing $x$,  and since $\L_{X_{s,s}^\mu}=\mu,$  by the standard Markov property of solutions to \eqref{ES0}, we have
\beq\label{ES1} P_{s,t}f(\mu):=\E f(X_{s,t}^\mu)= \int_{\R^d} P_{s,t}^\mu f(x) \mu(\d x),\ \ t\ge s, f\in \B_b(\R^d), \mu\in \scr P_2.\end{equation}
Since for any $g\in C_b^1(\R^d)$ the function
$\mu\mapsto \mu(g):= \int_{\R^d}g\d\mu$ is $L$-differentiable with
$D^L \mu(g)= \nn g$,  we first study the derivative of $P_{s,t}^\mu f(x)$ in $x$.

\beg{lem}\label{L0} Assume {\bf (H)}. Then for any $f\in \B_b(\R^d)$ and $t>s\ge 0$,
we have $P_{s,t}^\mu f\in C^1(\R^d)$ such that $(\nn P_{s,t}^\mu f)(x)$  is continuous in $(x,\mu)\in \R^d\times\scr P_2$, and
\beq\label{ES2} \|\nn P_{s,t}^\mu f\|_\infty \le \ff{CK_t \|f\|_\infty }{\ss{t-s}} \e^{C K_t},\ \ t>s, f\in \B_b(\R^d), \mu\in \scr P_2\end{equation}
holds for some constant $C>0.$\end{lem}

\beg{proof} Since {\bf (H)} implies that $P_{s,t}^*\mu$ is Lipschitz continuous in $\mu\in \scr P_2$, see for instance \cite{FYW1}, the desired assertions follow from    {\bf (H)} and the   Bismut formula
\beq\label{BS} \nn_v P_{s,t}^\mu f(x)  =  \E\bigg [\ff{f(X_{s,t}^{x,\mu})}{t-s}\int_s^t \big\<\{\si_r (\si_r\si_r^*)^{-1}\}(X_{s,r}^{x,\mu}, P_{s,r}^*\mu)  v_{s,r}^{x,\mu}, \d W_r\big\>\bigg],\ \ v\in\R^d\end{equation} for $ f\in \B_b(\R^d)$,
where $v_{s,t}^{x,\mu}:=\ff{\d}{\d \vv} X_{s,t}^{x+\vv v, \mu}|_{\vv=0} $ solves the linear SDE
\beq\label{LL}\d v_{s,t}^{x,\mu}= \{\nn_{v_{s,t}^{x,\mu}}b_t(\cdot, P_{s,t}^*\mu) \}(X_{s,t}^{x,\mu} )\d t +\{\nn_{v_{s,t}^{x,\mu}} \si_t(\cdot,P_{s,t}^*\mu)\}(X_{s,t}^{x,\mu}) \d W_t,
\ \ t\ge s, v_{s,s}^{x,\mu}=v.\end{equation}  By {\bf (H)},   $v_{s,t}^{x,\mu}$ is continuous in $(x,\mu)\in \R^d\times\scr P_2$ and
$$\E |v_{s,t}^{x,\mu}|^2\le |v|^2 \e^{C K_t},\ \ t\ge s, v\in \R^d$$ holds for some constant $C>0$, so that \eqref{BS} implies that $(\nn P_{s,t}^\mu f)(x)$  is continuous in $(x,\mu)\in \R^d\times\scr P_2$  and satisfies \eqref{ES2}.

To prove \eqref{BS}, for fixed $t>s$, take
$$h_u=\int_s^u \{\si_r^* (\si_r\si_r^*)^{-1}\}(X_{s,r}^{x,\mu}, P_{s,r}^*\mu)  v_{s,r}^{x,\mu} \d r,\ \ u\in [s,t].$$   Then the Malliavin derivative $w_r^{x,\mu}:= D_h X_{s,r}^{x,\mu}$ along $h$ solves the SDE
\beg{align*} &\d w_{s,r}^{x,\mu} = \Big[\{\nn_{w_{s,r}^{x,\mu}}b_r(\cdot, P_{s,r}^*\mu) \}(X_{s,r}^{x,\mu} )+\si_r(X_{s,r}^{x,\mu}, P_{s,r}^*\mu) h_r'\Big] \d r +\{\nn_{w_{s,r}^{x,\mu}} \si_r(\cdot,P_{s,r}^*\mu)\}(X_{s,r}^{x,\mu}) \d W_r\\
&=\Big[\{\nn_{w_{s,r}^{x,\mu}}b_r(\cdot, P_{s,r}^*\mu) \}(X_{s,r}^{x,\mu} )+ v_{s,r}^{x,\mu} \Big] \d r +\{\nn_{w_{s,r}^{x,\mu}} \si_r(\cdot,P_{s,r}^*\mu)\}(X_{s,r}^{x,\mu}) \d W_r ,\ \ r\in [s,t],  w_{s,s}^{x,\mu}=0,\end{align*} see for instance \cite[Proposition 3.5]{RW}.
It is easy to see from \eqref{LL} that $\bar v_r:= (r-s)v_{s,r}^{x,\mu}$ solves the same equation. By the uniqueness we obtain
$(t-s)v_{s,t}^{x,\mu}= D_h X_{s,t}^{x,\mu}$, so that the chain rule and the integration by parts formula yield
\beg{align*} &\nn_v P_{s,t}^\mu f(x)= \E\<\nn f(X_{s,t}^{x,\mu}), v_{s,t}^{x,\mu}\>=\ff 1 {t-s} \E\<\nn f(X_{s,t}^{x,\mu}), D_h X_{s,t}^{x,\mu}\>\\
&= \ff 1 {t-s} \E D_h\{ f(X_{s,t}^{x,\mu})\} =   \E\bigg[ \ff{f(X_{s,t}^{x,\mu})}{t-s} \int_s^t \big\<\{\si_r^\ast (\si_r\si_r^*)^{-1}\}(X_{s,r}^{x,\mu}, P_{s,r}^*\mu)  v_{s,r}^{x,\mu}, \d W_r\big\>\bigg].\end{align*}
\end{proof}

Combining \eqref{ES1} with Lemma \ref{L0}, we have the following result.

\beg{lem}\label{L1} Assume {\bf (H)}. Let $t>s$ and $f\in \B_b(\R^d)$. If for any $x\in\R^d$, the function $\mu\mapsto P_{s,t}^\mu f(x)$ is $L$-differentiable with
\beq\label{ES3}  \sup_{x\in \R^d} \big\|D^L \{P_{s,t}^\cdot f(x)\}(\mu)\big\|<\infty,\end{equation}
then $P_{s,t}f(\mu)$ is $L$-differentiable in $\mu$ with
\beq\label{ES4} D^L P_{s,t}f(\mu)= \nn P_{s,t}^\mu f + \int_{\R^d} D^L \{P_{s,t}^\cdot f(x)\}(\mu) \mu(\d x).\end{equation}
Consequently,  there exists a constant $C>0$ such that for any $f\in \B_b(\R^d)$ and $\mu\in\scr P_2,$
\beq\label{ES5} \big\|D^L P_{s,t}f(\mu)\big\|\le \ff{CK_t \|f\|_\infty }{\ss{t-s}} \e^{C K_t}+ \sup_{x\in \R^d} \big\|D^L \{P_{s,t}^\cdot f(x)\}(\mu)\big\|, \ \ t>s\ge 0.\end{equation}
\end{lem}

\beg{proof} Obviously, \eqref{ES5} is implied by \eqref{ES2} and \eqref{ES4}. So, we only need to prove that $P_{s,t}f(\mu)$ is $L$-differentiable and satisfies \eqref{ES4}.

(1) We first prove that  $P_{s,t}f(\mu)$ is  intrinsically differentiable and satisfies \eqref{ES4}.   For any $g\in C_b^1(\R^d)$, the function
$\mu\mapsto \mu(g):= \int_{\R^d}g\d\mu$ is $L$-differentiable with
$D^L \mu(g)= \nn g$. So, for any $\mu\in \scr P_2$, the function
\begin{align}\label{Cy}\scr P_2\ni\nu\mapsto  P_{s,t}^\mu f(\nu):= \int_{\R^d} P_{s,t}^\mu f\d\nu
 \end{align}
 is $L$-differentiable with
$D^L (P_{s,t}^\mu f)(\nu)= \nn P_{s,t}^\mu f,\ \nu\in\scr P_2.$ Combining this with \eqref{ES1}, \eqref{ES2} and \eqref{ES3},  and using the dominated convergence theorem, we conclude that the map
$$L^2(\R^d\to\R^d;\mu)\ni \phi\mapsto  D^L_\phi P_{s,t}f(\mu)= \<\nn P_{s,t}^\mu f,\phi\>_{L^2(\mu)} + \int_{\R^d} D^L_\phi \{P_{s,t}^\cdot f(x)\}(\mu) \mu(\d x)$$
is a bounded linear functional, so that by definition,
$P_{s,t}f(\mu)$ is  intrinsically differentiable in $\mu\in \scr P_2$, and the formula \eqref{ES4} holds true.

(2) By \eqref{ES4}, for any $\phi\in L^2(\R^d\to\R^d;\mu)$, we have
\beg{align*}& P_{s,t}f (\mu\circ({\rm Id}+\phi)^{-1} )- P_{s,t}f (\mu )- D_\phi^L P_{s,t}f(\mu)\\
&= \int_{\R^d} \big\{P_{s,t}^{\mu\circ ({\rm Id}+\phi)^{-1}} f(x+\phi(x)) - P_{s,t}^{\mu\circ ({\rm Id}+\phi)^{-1}} f(x) - \<\nn P_{s,t}^{\mu\circ ({\rm Id}+\phi)^{-1}} f(x), \phi(x)\>\big\}\mu(\d x)\\
& \quad + \int_{\R^d} \big\{P_{s,t}^{\mu\circ ({\rm Id}+\phi)^{-1}} f(x)- P_{s,t}^\mu f(x) - D_\phi^L[P_{s,t}^\cdot f(x)](\mu) \big\}\mu(\d x)\\
&\quad +\int_{\R^d} \<\nn P_{s,t}^{\mu\circ ({\rm Id}+\phi)^{-1}} f(x)- \nn P_{s,t}^\mu f(x), \phi(x)\>\mu(\d x).\end{align*}
Combining this with Lemma \ref{L0}, \eqref{ES3}, and the $L$-differentiability of $P_{s,t}^\mu f(x)$  in $\mu$, we may apply the dominated convergence theorem to derive
$$\lim_{\|\phi\|_{L^2(\mu)}\downarrow 0} \ff{|P_{s,t}f (\mu\circ({\rm Id}+\phi)^{-1} )- P_{s,t}f (\mu )-D_\phi^L P_{s,t}f(\mu)|}{\|\phi\|_{L^2(\mu)}} =0,$$
that is, $P_{s,t}f(\mu)$ is $L$-differentiable.
\end{proof}

 According to Lemma \ref{L1},  to estimate $\|D^L P_{s,t} f(\mu)\|,$ it remains to investigate the  $L$-derivative of $P_{s,t}^\mu f(x)$ in $\mu$. To this end, we let $p_{s,t}^\mu(x,y)$ be the  heat kernel of $P_{s,t}^\mu$ for $t>s$, which exists and is differentiable in $x$ and $y$  under conditions {\bf (H)}. We have
\beq\label{MA0} P_{s,t}^\mu f(x)= \int_{\R^d} p_{s,t}^\mu (x,y) f(y) \d y,\ \ f\in \B_b(\R^d), t>s, x\in \R^d.\end{equation}
 So, to investigate the   $L$-derivative of $P_{s,t}^\mu f(x),$ we need to study  that of $p_{s,t}^\mu(x,y)$, for which we will use the heat kernel parameter expansion.

 \subsection{Parameter expansion for $p_{s,t}^\mu$}

Since  heat kernel  $p_{s,t}^\mu$ is less explicit, we make use of its parameter expansion with respect to the heat kernel of the Gaussian process
 $$X_{s,r,t}^{x,\mu,z}= x + \int_r^t b_u(z, P_{s,u}^*\mu)\d u + \int_r^t \si_u(z, P_{s,u}^*\mu)\d W_u,\ \ t\ge r\ge s\ge 0, x\in\R^d$$
 for fixed $z\in \R^d$ and $\mu\in \scr P_2$.  For any $t\ge r\ge s\ge 0$, let
\beq\label{MA}\beg{split} & m_{s,r,t}^{\mu,z} := \int_r^t b_u(z, P_{s,u}^*\mu)\d u,\ \ m^{\mu,z}_{s,t}:=m^{\mu,z}_{s,s,t},\\
& \ a_{s,r,t}^{\mu,z} :=  \int_r^t (\si_u\si_u^*)(z, P_{s,u}^*\mu)\d u, \ \ a^{\mu,z}_{s,t}:=a^{\mu,z}_{s,s,t}.\end{split}\end{equation}
By {\bf (H)}, we have
\begin{align}\label{mv}
|m_{s,r,t}^{\mu,z}|+|a_{s,r,t}^{\mu,z}|\leq (t-r)K_t,\ \ t\ge r\ge s\ge 0.
\end{align}
Obviously,  the law of $X_{s,r,t}^{x,\mu,z}$ is the $d$-dimensional normal distribution entered at $x+m_{s,r,t}^{\mu,z}$ with covariance matrix $a_{s,r,t}^{\mu,z},$ i.e.  the distribution density function   is
\beq\label{ES6} p_{s,r,t}^{\mu,z}(x,y)= \ff{\exp[-\ff 1 2 \<(a_{s,r,t}^{\mu,z})^{-1}(y-x-m_{s,r,t}^{\mu,z}),  y-x-m_{s,r,t}^{\mu,z}\>]}{(2\pi)^{\ff d 2} ({\rm det} \{a_{s,r,t}^{\mu,z}\})^{\ff 1 2}},\ \ y\in \R^d, t>r\ge s.\end{equation}
When $r=s$, we simply denote $p_{s,t}^{\mu,z}=p_{s,s,t}^{\mu,z}$, so that
\beq\label{ES6'} p_{s,t}^{\mu,z}(x,y)= \ff{\exp[-\ff 1 2 \<(a_{s,t}^{\mu,z})^{-1}(y-x-m_{s,t}^{\mu,z}),  y-x-m_{s,t}^{\mu,z}\>]}{(2\pi)^{\ff d 2} ({\rm det} \{a_{s,t}^{\mu,z}\})^{\ff 1 2}},\ \ y\in \R^d, t>  s.\end{equation}
For any $0\le s\le r<t$ and $y,z\in\R^d$, let
\beq\label{ES7}\beg{split}  H_{s,r,t}^{\mu} (y,z):=  & \left\<b_r (z, P_{s,r}^*\mu) -b_{r}(y, P_{s,r}^*\mu),\nabla p_{s,r,t}^{\mu,z}(\cdot,z)(y)\right\>\\
&+\frac{1}{2}\mathrm{tr}\left[\left\{(\sigma_{r}\sigma^\ast_{r})(z, P_{s,r}^*\mu)-(\sigma_{r}\sigma^\ast_{r})(y, P_{s,r}^*\mu)\right\}\nabla^2 p_{s,r,t}^{\mu,z} (\cdot,z)(y)\right].
\end{split} \end{equation}
By the parameter expansion, see for instance \cite[Lemma 3.1]{KM00},  we have
\beq\label{ES8} p_{s,t}^\mu(x,z) = p_{s,t}^{\mu,z}(x,z) + \sum_{m=1}^\infty \int_s^t \d r \int_{\R^d} H_{s,r,t}^{\mu,m}  (y,z) p_{s,r}^{\mu,z} (x,y) \d y,\end{equation}
where   $H_{s,r,t}^{\mu,m}$ for $m\in \mathbb N$ are defined by
\beq\label{HU} \beg{split}& H_{s,r,t}^{\mu, 1}:= H_{s,r,t}^{\mu},\\
&H_{s,r,t}^{\mu, m} (y,z) := \int_r^t\d u \int_{\R^d}  H_{s,u,t}^{\mu, m-1}(z',z) H_{s,r,u}^{ \mu}(y,z') \d z',\ m\ge 2.\end{split}\end{equation}
Combining  \eqref{ES8}  with \eqref{MA}, \eqref{ES6} and \eqref{ES6'},  to estimate $D^L P_{s,t}^\mu f$,  it suffices to study the $L$-derivative of  $b_r(y, P_{u_1,u_2}^* \mu)$ and
$ (\si_r\si_r^*)(y, P_{u_1,u_2}^* \mu)$ in $\mu$ for $r\ge 0$ and $u_2\ge u_1\ge 0$.
So, we present the following lemma.

\beg{lem}\label{L2}  Assume {\bf (H)} and let $t>s\ge 0$. Then for any $F\in C^1(\scr P_2)$ with bounded $\|D^L F\|$,  $F(P_{s,t}^*\mu)$ is $L$-differentiable in $\mu$ such that
\beq\label{LL1} \|D^L F(P_{s,t}^*\cdot)(\mu)\|\le \|D^L F\|_\infty \e^{4 K_t (t-s)}.\end{equation} Consequently,
for any $r\ge 0, t\ge s\ge 0$ and $y\in \R^d$, $b_r(y, P_{s,t}^* \mu)$ and $(\si_r\si_r^*)(y,P_{s,t}^*\mu)$ are   $L$-differentiable in $\mu$, and
$$\max\Big\{\|D^L b_r(y, P_{s,t}^*\cdot)(\mu)\|,  \|D^L (\si_r\si_r^*)(y, P_{s,t}^*\cdot)(\mu)\|\Big\}\le K_r\e^{4K_t(t-s)}, \ \  \mu\in \scr P_2.$$
\end{lem}

\beg{proof} It  suffices  to prove the first assertion.
We first prove the intrinsic differentiability.  Let $\mu\in \scr P_2$ and $\phi\in L^2(\R^d\to\R^d;\mu).$  Since  $\L_{X_{s,s}^\mu}=\mu$ implies
$$\L_{X_{s,s}^\mu+\vv \phi(X_{s,s}^\mu)}=\mu\circ({\rm Id}+\vv\phi)^{-1},\ \  \vv\ge 0,$$ we have $\L_{X_{s,t}^\vv}= P_{s,t}^*(\mu\circ({\rm Id}+\vv\phi)^{-1})$ for $X_{s,t}^\vv$ solving \eqref{E11}   with initial value  $X_{s,s}^\vv= X_{s,s}^\mu+\vv \phi(X_{s,s}^\mu).$
By \cite[Proposition 3.1]{RW} for $\eta=\phi(X_0^\mu)$ and \cite[(4.21)]{RW} for time $s$ replacing $0$,  for any $\dd\ge 0$,
$$v_{s,t}^{\phi,\dd}:= D_\phi^L X_{s,t}^\dd= \lim_{\vv\downarrow 0} \ff{X_{s,t}^{\dd+\vv}-X_{s,t}^\dd}\vv,\ \ t\ge s$$
exists in  $L^2(\OO\to C([s,T];\R^d);\P)$ for any $T>0$, and solves the linear SDEs:
\begin{equation}\begin{split}\label{vp}
\d v_{s,t}^{\phi,\dd}&= \Big[\nabla_{v_{s,t}^{\phi,\dd}}b_t(X^\dd_{s,t}, \L_{X^\dd_{s,t}})+\E \big\{\<D^Lb_t(z,\cdot)( \L_{X^\dd_{s,t}})(X^\dd_{s,t}), v_{s,t}^{\phi,\delta}\>\big\}\big|_{z=X^\dd_{s,t}}\Big]\d t\\
&\quad +\Big[\nabla_{v_{s,t}^{\phi,\dd}}\sigma_t(X^\dd_{s,t}, \L_{X^\dd_{s,t}})+\E \big\{D^L\sigma_t(z,\cdot)( \L_{X^\dd_{s,t}})(X^\dd_{s,t}) v_{s,t}^{\phi,\dd}\big\}\big|_{z=X^\dd_{s,t}}\Big]\d W_{t},\\
&\quad    v_{s,s}^{\phi,\dd}=\phi(X_0),\ \ t\ge s.
\end{split}\end{equation}
Fromm  {\bf (H)} we see that $v_{s,t}^{\phi,\vv}$ is continuous in $\vv$ and
\beq\label{ES9} \E |v_{s,t}^{\phi,\dd}|^2\le \mu(|\phi|^2) \e^{8 (t-s) K_t},\ \ t\ge s, \phi\in L^2(\R^d\times\R^d;\mu).\end{equation}

By the chain rule, see for instance \cite[Proposition 3.1]{RW}, we have
\beq\label{LL2} D^L_\phi F( P_{s,t}^*\cdot)(\mu) =\ff{\d}{\d\vv} F(\L_{X_{s,t}^\vv})\Big|_{\vv=0}= \E\<(D^L F)(P_{s,t}^*\mu)(X_{s,t}^\mu), v_{s,t}^{\phi,0}\>.\end{equation}
Combining this with {\bf(H)} and \eqref{ES9},  we obtain
\beg{align*} &|D^L_\phi  F( P_{s,t}^*\cdot)(\mu)|\le \|(D^L F)(P_{s,t}^*\mu)\|\ss{ \E |v_{s,t}^{\phi,0}|^2}\\
&\le \|\phi\|_{L^2(\mu)}  \|D^LF\|_\infty  \e^{4(t-s) K_t},\ \ \phi\in L^2(\R^d\to\R^d;\mu).\end{align*}
Therefore, $ F(P_{s,t}^* \mu)$ is intrinsically differentiable in $\mu$ such that \eqref{LL1} holds.

It remains to verify the $L$-differentiability.  By the chain rule and \eqref{LL2}, we obtain
\beg{align*} &F(P_{s,t}^* \mu\circ({\rm Id} +\phi)^{-1}) - F(P_{s,t}^*\mu)- D_\phi^L F(P_{s,t}^*\cdot)(\mu)
= \int_0^1\ff{\d}{\d\vv} F(\L_{X_{s,t}^\vv})\d\vv - D_\phi^L F(P_{s,t}^*\cdot)(\mu) \\
&= \int_0^1 \big\{\E\<(D^L F)(P_{s,t}^*\mu\circ({\rm Id}+\vv\phi)^{-1}) (X_{s,t}^\vv), v_{s,t}^{\phi,\vv}\>-  \E\<(D^L F)(P_{s,t}^*\mu)(X_{s,t}^\mu), v_{s,t}^{\phi,0}\>\big\}\d\vv.\end{align*}
Combining this with $F\in C^1(\scr P_2)$ with bounded $\|D^LF\|$,   the continuity of $v_{s,t}^{\phi,\vv}$ in $\vv$, \eqref{ES9}, and that $X_{s,t}^\vv\to X_{s,t}^\mu$ when $\|\phi\|_{L^2(\mu)}\to 0$, by the dominated theorem we prove
$$\lim_{\|\phi\|_{L^2(\mu)} \downarrow 0} \ff{|F(P_{s,t}^* \mu\circ({\rm Id} +\phi)^{-1}) - F(P_{s,t}^*\mu)- D_\phi^L F(P_{s,t}^*\cdot)(\mu)|}{\|\phi\|_{L^2(\mu)}}=0,$$
thus, $F(P_{s,t}^*\mu)$ is $L$-differentiable in $\mu$.

\end{proof}

\section{Proof of Theorem \ref{T1} }
According to Lemma \ref{L1}, \eqref{MA0} and \eqref{ES8}, to estimate $\|D^L P_{s,t}f(\mu)\|$, it suffices to handle
the   derivative of $ p_{s,t}^\mu$ and $H_{s,r,t}^{\mu,m}$ in $\mu$.  To this end, for fixed $T>0$, we introduce the Gaussian  heat kernel
 \beq\label{ES10} h_{T}(s,y)=\ff{\exp[-\ff{|y|^2}{8sK_T} ]}{(8\pi s K_T)^{\ff d 2}},\ \ y \in \R^d, s>0,\end{equation} which satisfies the Chapman-Kolmogorov equation
\beq\label{CK} \int_{\R^d} h_T(s_1, y-z) h_T(s_2, z)\d z= h_T(s_1+s_2, y),\ \ s_1,s_2>0, y\in \R^d.\end{equation}
By {\bf(H)}, there exists a constant $K_1(T)$, which increases in $T$, such that
 \begin{align*}  p_{s,r,t}^{\mu,z}(y,z)
&\le K_1(T) h_{T}(t-r, y-z) \e^{-\ff{|y-z|^2}{8(t-r)K_T}},\ \ y,z\in\R^d, 0\le s\le r<t\le T,\mu\in \scr P_2.\end{align*}
Consequently, there exists a constant $K_2(T)$, which increases in $T$, such that
\beq\label{G2}\beg{split} &p_{s,r,t}^{\mu,z}(y,z) \Big(1+ \ff{|y-z|^2}{t-r} +\ff{|y-z|}{(t-r)^{\ff 1 2}}\Big)\\
& \le K_2(T) h_{T}(t-r, y-z),\ \ y,z\in\R^d, 0\le s\le r<t\le T,\mu\in \scr P_2.\end{split}\end{equation}

\beg{lem}\label{L3} Assume {\bf (H)}.  There exists a constant  $\bar{K}_T>0$ which increases in $T>0$, such that  for any $ 0\le s\le r<t\le T, y,z\in\R^d$ and $m\ge 1$,
$p_{s,r,t}^{\mu,z}(y,z)$ and $H_{s,r,t}^{\mu,m} $ are $L$-differentiable in $\mu\in \scr P_2$ satisfying
\beg{align} & \label{Ta1} \|D^L\{ p_{s,r,t}^{\cdot,z}(y,z) \}(\mu) \|\le \bar{K}_T  h_T (t-r, y-z),\\
&\label{Ta2}  | H_{s,r,t}^{\mu,m}(y,z) |  \le \ff{\bar{K}_T^m (t-r)^{\ff m 2 -1} }{\GG(\ff m 2)}  h_T  (t-r, y-z),\ \ m\ge 1,\\
& \label{Ta3} \|D^L \{H_{s,r,t}^{\cdot,m}(y,z)\} (\mu) \|\le   \ff{m \bar{K}_T^m (t-r)^{\ff m 2 -1} }{\GG(\ff m 2)}  h_T (t-r, y-z),\ \ m\ge 1.\end{align}

\end{lem}
\beg{proof}   By {\bf (H)}, we have $|m_{s,r,t}^{\mu,z}|\le (t-r)K_T$, so that   \eqref{G2}  yields
\beq\label{ES14}\begin{split} & p_{s,r,t}^{\mu,z} (y,z) \Big(1+\ff{|y-z-m_{s,r,t}^{\mu,z}|^2} {t-r}+ \ff{|y-z-m_{s,r,t}^{\mu,z}|}{(t-r)^{\ff 1 2}} \Big) \le C_1(T)   h_T(t-r, y-z)  \end{split}\end{equation}
for some constant $C_1(T)>0$ increasing in $T$, and all $0\le s\le r<t\le T, \mu \in \scr P_2$ and $y,z\in \R^d$.
 Combining this with {\bf (H)},  \eqref{ES6}, \eqref{ES14} and applying Lemma \ref{L2}, we   prove the $L$-differentiability of $p_{s,r,t}^{\mu,z}(y,z)$ in $\mu\in\scr P_2$ and  the estimate \eqref{Ta1}.

Next,  by {\bf (H)}, \eqref{ES6}, \eqref{ES7} and \eqref{ES14}, we find   constants $C_2(T), C_3(T)>0$ increasing in $T>0$ such that for any  $0\le s\le r<t\le T, \mu \in \scr P_2$ and $y,z\in \R^d$,
\beq\label{A1} \beg{split} |H_{s,r.t}^\mu(y,z)|&\le C_2(T)p_{s,r,t}^{\mu,z} (y,z) |y-z|\Big(\ff {1}{t-r}+\ff{|y-z-m_{s,r,t}^{\mu,z}|^2}{ (t-r)^2}+ \ff{|y-z-m_{s,r,t}^{\mu,z}|}{ t-r} \Big) \\
&\le  C_3(T) (t-r)^{-\ff 1 2}  h_T(t-r, y-z).   \end{split}\end{equation}
Assume that for some $k\ge 1$ we have
$$ | H_{s,r,t}^{\mu,k}(y,z) |   \le C_3(T)^k    (t-r)^{\ff k 2 -1} \Big(\prod_{i=1}^{k-1} \bb\Big(\ff i 2, \ff 1 2\Big) \Big) h_T  (t-r, y-z).$$
Combining this with \eqref{HU}, \eqref{CK},  and \eqref{A1},  we derive
\beg{align*} & | H_{s,r,t}^{\mu,k+1}(y,z) | \le  \int_r^t\d u \int_{\R^d}  |H_{s,u,t}^{\mu, k}(z',z) H_{s,r,u}^{\mu}(y,z') | \d z'\\
&\le  C_3(T)^{k+1}  h_T(t-r, y-z)   \Big(\prod_{i=1}^{k-1} \bb\Big(\ff i 2, \ff 1 2\Big) \Big)  \int_r^t (t-u)^{\ff k 2-1} (u-r)^{-\ff 1 2}\d u\\
&=  C_3(T)^{k+1} (t-r)^{\ff{k+1}{2}-1}  h_T(t-r, y-z)   \Big(\prod_{i=1}^{k} \bb\Big(\ff i 2, \ff 1 2\Big) \Big).\end{align*}
In conclusion, for any $m\ge 1$, we have
$$| H_{s,r,t}^{\mu,m}(y,z) |  \le C_3(T)^m  (t-r)^{\ff m 2 -1}   \Big(\prod_{i=1}^{m-1} \bb\Big(\ff i 2, \ff 1 2\Big) \Big) h_T  (t-r, y-z),$$
which implies \eqref{Ta2}  for $\bar{K}_T= C_3(T)\GG(\ff 1 2)$,  since
\beq\label{AU} \prod_{i=1}^{m-1} \bb\Big(\ff i 2, \ff 1 2\Big)= \prod_{i=1}^{m-1} \ff{\GG(\ff i 2)\GG(\ff 1 2)}{\GG(\ff {i+1}2)} = \ff{\GG(\ff 1 2)^m}{\GG(\ff m 2)}.\end{equation}

Finally, by {\bf (H)}, \eqref{ES6'}, \eqref{ES7}, Lemma \ref{L2} and \eqref{A1}, we see that $H_{s,r,t}^{\mu,m}$ is $L$-differentiable in $\mu$, and there exist  constants  $C_4(T), C_5(T)\ge C_3(T)$ increasing in $T>0$ such that
\beq\label{A2} \beg{split}   &\| D^L\{H_{s,r,t}^\cdot (y,z)\}(\mu)\|\\
&\le C_4(T)p_{s,r,t}^{\mu,z} (y,z) |y-z|\Big(\ff {1}{t-r}+\ff{|y-z-m_{s,r,t}^{\mu,z}|^2}{ (t-r)^2}+ \ff{|y-z-m_{s,r,t}^{\mu,z}|}{ t-r} \Big) \\
&\le  C_5(T) (t-r)^{-\ff 1 2}  h_T(t-r, y-z).   \end{split}\end{equation}
Assume that for some $k\ge 1$ we have
$$\| D^L \{H_{s,r,t}^{\cdot,k}(y,z)\}(\mu)\|   \le k C_5(T)^k    (t-r)^{\ff k 2 -1} \Big(\prod_{i=1}^{k-1} \bb\Big(\ff i 2, \ff 1 2\Big) \Big) h_T  (t-r, y-z).$$
Combining this with \eqref{HU}, \eqref{CK},  and \eqref{A2},   we derive
\beg{align*} & \| D^L \{H_{s,r,t}^{\cdot,k+1}(y,z) \}(\mu)\|  \\
& \le   \int_r^t\d u \int_{\R^d} \Big\{\|D^L\{H_{s,u,t}^{\cdot, k}(z',z)\}(\mu)\|\cdot | H_{s,r,u}^{\mu}(y,z') |  \\
&\qquad\qquad  \qquad \quad+ |H_{s,u,t}^{\mu, k}(z',z)  | \cdot\| D^L\{H_{s,r,u}^{\cdot}(y,z')\}(\mu)\|\Big\}\d z'\\
&\le  (k+1) C_5(T)^{k+1}  h_T(t-r, y-z)   \Big(\prod_{i=1}^{k-1} \bb\Big(\ff i 2, \ff 1 2\Big) \Big)  \int_r^t (t-u)^{\ff k 2-1} (u-r)^{-\ff 1 2}\d u\\
&=  (k+1)C_5(T)^{k+1} (t-r)^{\ff{k+1}2-1}  h_T(t-r, y-z)   \Big(\prod_{i=1}^{k} \bb\Big(\ff i 2, \ff 1 2\Big) \Big).\end{align*}
This together with \eqref{AU} implies   \eqref{Ta3} for $\bar{K}_T=C_5(T)\GG(\ff 1 2). $
\end{proof}

We are now ready to prove the main result.

\beg{proof} [Proof of Theorem \ref{T1}] By Lemma \ref{L3} with \eqref{ES8} and \eqref{CK},   $p_{s,t}^\mu(x,z)$ is
$L$-differentiable in $\mu$ for $t>s$, and
there exists a constant $\dd_T>0$ increasing in $T>0$ such that
  \beg{equation}\begin{split} \label{Dp} & \|D^L \{p_{s,t}^\cdot(x,z)\}(\mu)\|\le \bar{K}_T h_T(t-s, x-z)\\&\qquad + \sum_{m=1}^\infty \ff{(m+1) \bar{K}_T^{m+1}}{\GG(\ff m 2)} \int_s^t (t-s)^{\ff m 2 -1} \d r \int_{\R^d} h_T(t-r, y-z) h_T(r-s, x-y)\d y\\  &\le \dd_T h_T(t-s, x-z). \ \ \  0\le s<t\le T, x,z\in \R^d, \mu\in \scr P_2.
  \end{split}\end{equation}
 This and \eqref{ES1} imply that $P_{s,t}^\mu f(x)$ is $L$-differentiable in $\mu$  such that
 $$\|D^L\{P_{s,t}^\cdot f(x)\}(\mu)\|\le \|f\|_\infty\int_{\R^d} \|D^L \{p_{s,t}^\cdot (x,z)\}(\mu)\|\d z\le \dd_T\|f\|_\infty$$
 holds for all $0\le s<t\le T, f\in \B_b(\R^d)$ and $\mu\in \scr P_2$. Combining this with Lemma \ref{L1}, we prove that $P_{s,t}f(\mu)$ is $L$-differentiable in $\mu$ and \eqref{EST1} holds for some increasing $C: [0,\infty)\to (0,\infty).$
According to the proof of \cite[Corollary 2.2(2)]{RW}, we can show that \eqref{EST1} implies \eqref{EST2}. We include below  a simple  proof   for completeness.

Since $C_b^1(\R^d)$ is dense in $L^1(P_{s,t}^*\mu+P_{s,t}^*\nu)$, \eqref{EST2} is equivalent to
 \beq\label{EST2'} |P_{s,t} f(\mu)- P_{s,t} f(\nu)|\le \ff{C_t\|f\|_\infty}{\sqrt{t-s}}\W_2(\mu,\nu),\ \ t>s, f\in C_b^1(\R^d), \mu,\nu\in \scr P_2.\end{equation}
Let $f\in C_b^1(\R^d)$ be fixed. We first prove this inequality  for     $\mu,\nu$ with compact supports.

 Let $\xi,\eta$ be two bounded random variables such that $\L_\xi=\mu, \L_\eta=\nu$ and
 $$\E |\xi-\eta|^2= \W_2(\mu,\nu)^2.$$ By  Proposition 3.1 in \cite{RW} and \eqref{EST1}, we obtain
\beg{align*}  &|P_{s,t}f(\mu)-P_{s,t}f(\nu)|= \bigg|\int_0^1 \ff{\d }{\d r} P_{s,t} f(\L_{r \xi+(1-r)\eta})  \d r\bigg| \\
&\le \int_0^1 \big|\E\<D^L P_{s,t}f(\L_{r \xi+(1-r)\eta})  (r \xi+(1-r)\eta), \xi-\eta\> \big|\d r \le \ff{C_t \|f\|_\infty}{\ss{t-s}} \W_2(\mu,\nu).\end{align*}
So, \eqref{EST2'} holds.

Next, for any $\mu,\nu\in \scr P_2$,  we choose   $\{\mu_n,\nu_n\}_{n\ge 1}\subset \scr P_2$ with compact supports such that
$$\lim_{n\to\infty}\big\{\W_2(\mu,\mu_n)+\W_2(\nu,\nu_n)\big\}=0.$$
Then  by the last step,
\beq\label{LST2} |P_{s,t}f(\mu_n)-P_{s,t}f(\nu_n)|\le \ff{C_t \|f\|_\infty}{\ss{t-s}} \W_2(\mu_n,\nu_n),\ \ n\ge 1.\end{equation} If $P_{s,t}f(\gg)$ is continuous in $\gg\in \scr P_2$, then by letting $n\to\infty$ we obtain the desired estimate \eqref{EST2'}.  To prove the continuity,
for any $\gg_1,\gg_2\in \scr P_2$,  let $\xi_1,\xi_2$ be $\F_0$-measurable random variables  such that $\L_{\xi_i}=\gg_i, i=1,2,$ and
$$\W_2(\gg_1,\gg_2)^2 =\E|\xi_1-\xi_2|^2.$$  For any   $\vv\in [0,1]$, let $X_{s,t}^\vv$ solve \eqref{E11} with initial value $X_{s,s}^\vv:= \vv\xi_1+(1-\vv)\xi_2$.
By \cite[Proposition 3.2 and (4.2)]{RW},
$$\nn_{\xi_1-\xi_2} X_{s,t}^\vv:=\ff{\d}{\d\vv} X_{s,t}^\vv$$ exists in $L^2(\P)$ with
$$\E|\nn_{\xi_1-\xi_2} X_{s,t}^\vv|^2\le c(t)\E|\xi_1-\xi_2|^2= c(t) \W_2(\gg_1,\gg_2)^2$$ for some constant $c(t)>0$.  Then
\beg{align*} &|P_{s,t}f(\gg_1)-P_{s,t}f(\gg_2)| =|\E f(X_{s,t}^1)-\E f(X_{s,t}^0)|  = \bigg|\int_0^1 \ff{\d }{\d \vv} \E f(X_{s,t}^\vv)    \d \vv \bigg| \\
&\le \int_0^1 \big|\E\<\nn f ( X_{s,t}^\vv), \nn_{\xi_1-\xi_2} X_{s,t}^\vv\> \big|\d s \le \ss{ c(t) }\|\nn f\|_\infty  \W_2(\gg_1,\gg_2).\end{align*}
Therefore, $P_{s,t}f(\gg)$ is continuous in $\gg\in\scr P_2$ and the proof is then finished.  \end{proof}


text:
\begin{document}
\tikzstyle{block} = [draw, fill=white!20, rectangle,
    minimum height=3em, minimum width=6em]
\tikzstyle{sum} = [draw, fill=blue!20, circle, node distance=1cm]
\tikzstyle{input} = [coordinate]
\tikzstyle{output} = [coordinate]
\tikzstyle{pinstyle} = [pin edge={to-,thin,black}]

\tikzstyle{int}=[draw, fill=blue!20, minimum size=2em]
\tikzstyle{init} = [pin edge={to-,thin,black}]

\maketitle

\begin{abstract}
\kas{\icl{We consider the problem of distributed secondary frequency regulation in power networks such that stability and an optimal power allocation are attained.
This is a problem that has been widely studied in the literature, and two main control schemes have been proposed, usually referred to as 'primal-dual' and
'distributed averaging proportional-integral (DAPI)' respectively.}
 However, each has its limitations, with the former requiring knowledge of uncontrollable demand, which can be difficult to obtain \icl{in real time}, and \icl{with} the existing literature on the latter \icl{being based on static models for} generation and demand.
We propose a novel control scheme that overcomes these issues by making use of generation measurements in \icl{the control policy.} In particular, our analysis
allows distributed stability and optimality guarantees to be deduced with practical measurement requirements and permits a broad range of
linear generation dynamics, \icl{that can be of} higher order, to be incorporated in the power network.

\icl{We show how the controller parameters can be selected in} a computationally efficient way by solving appropriate linear matrix inequalities (LMIs).
 Furthermore, we demonstrate how the proposed analysis applies \ak{to several} examples of turbine governor models.}

 The practicality of our analysis is demonstrated with simulations on the Northeast Power Coordinating Council (NPCC) 140-bus system that verify that our proposed controller achieves convergence to \icl{the} nominal frequency and an economically optimal power allocation.

\end{abstract}

\section{Introduction}\label{sec:Introduction}

\textbf{Motivation:}
\il{There is currently a  growing attention} on renewable sources of generation as a result of environmental concerns,
with their penetration in power \il{networks expected to} grow over the next years \cite{lund2006large, ipakchi2009grid}.
This will greatly increase the amount of active elements in the power network making its electromechanical behaviour less predictable and centralized control approaches increasingly difficult to implement. This highlights the importance of investigating distributed schemes that will guarantee power network stability when such devices are included. Over the past few years, these concerns  have motivated research on distributed schemes with applications on both primary \cite{kasis2016primary, molina2011decentralized} and secondary frequency control \cite{trip2016internal, P30_Low}.

The introduction of highly distributed schemes for frequency control raises  an issue of economic optimality in the power allocation. Attempts to resolve this issue in \ak{the literature resulted in} devising appropriately constructed optimization problems that ensure economic optimality and \ak{asking} for the system equilibria to be solutions of these problems. \ak{It is evident in \il{the} literature that a synchronising variable is useful for optimality to be achieved.
 While in primary control studies frequency is used as the synchronizing signal \cite{kasis2016primary, zhaotopcu, devane2016primary}, in secondary control studies some other signal, resulting from a suitably designed controller, has to be employed (e.g. \cite{trip2016internal, kasis2017secondary, dorfler}).
Therefore, the study of how distributed controllers  that achieve optimality for secondary frequency regulation should be designed is an interesting \icl{problem of practical relevance}.}

\textbf{Literature survey:}
There are many recent studies associated with  \icl{stability and optimality in distributed} secondary frequency control. Many of those, involve control schemes with dynamics that follow from a primal/dual algorithm associated with some optimal power allocation optimisation problem \cite{P30_Low},  \cite{kasis2017secondary}, \cite{low2014distributed}, \cite{zhao2015distributed}, \cite{li2016connecting}.
This approach allows to take into account economic considerations along \il{with} the objectives of secondary frequency control. Furthermore, it allows for stability guarantees when high order and non linear generation dynamics are considered. However, such schemes require knowledge of demand, which can in some cases limit their practicality.

A different approach for optimal distributed secondary frequency regulation involves the use of distributed averaging proportional integral (DAPI) controllers \cite{trip2016internal, dorfler, dorfler2014plug, simpson2013synchronization, andreasson2014distributed}.
DAPI controllers are simpler than primal/dual inspired ones, requiring only knowledge of local frequency and exchange a synchronization signal without requiring any generation or load measurements.
On the other hand, existing results in \icl{the literature incorporating DAPI controllers} do not accommodate

\icl{higher-order generation dynamics and restrict the anlysis to static generation.}

\kas{For a thorough survey of distributed approaches for stability and optimality in power systems, see \cite{molzahn2017survey}.}

\textbf{Main contributions:}
In this paper, we propose a \il{control scheme, that will be referred to as} distributed averaging dynamic output control (DADOC). \kas{\icl{A distinctive feature of this scheme is that it overcomes} the limitations of DAPI and Primal-Dual controllers allowing for stability and optimality guarantees when high order generation dynamics are considered, \icl{without} imposing load measurement requirements.}

\icl{In particular, the} proposed scheme dynamics achieve the synchronization of an exchange variable,  which is necessary for optimality, and contain frequency dependent terms that ensure that \il{the} frequency attains its nominal value at equilibrium.
Furthermore, DADOC schemes include feedback from generation output, which is a key part in \ak{their design}, enhancing their stability properties, and allowing the inclusion of a broad range of linear generation dynamics.

\icl{DADOC schemes have} the advantage over DAPI schemes that they allow \ank{the inclusion of} \icl{higher order} generation dynamics, by imposing only an additional condition for knowledge of generation output. Although primal dual inspired schemes are also able to provide stability guarantees \ank{when high order generation schemes are considered}, those impose a requirement for knowledge of demand, which can
\kas{be} restrictive. DADOC controllers share benefits from both schemes, with \ak{mild measurement} requirements (generation and frequency) and

\icl{allow the} inclusion of highly relevant generation dynamics that are not included in any of the aforementioned control schemes.

\kas{Our analysis provides \icl{also} conditions for the design of the controller gains such that stability and optimality are guaranteed.}
An important feature \icl{of} the proposed conditions, is that those can be efficiently verified by means of a linear matrix inequality (LMI). Various examples of relevant generation dynamics are provided to demonstrate the importance of our \kas{contribution.}

\textbf{Paper structure:} The remainder of the paper is structured as follows: Section \ref{sec: Preliminaries} provides some basic notation and in Section \ref{Problem_Formulation} we present the power network model and the generation dynamics. In Section \ref{sec:Power_Command_Dynamics} we present our proposed controller and provide conditions such that an optimal power allocation is achieved at steady state.
In Section \ref{sec:Convergence} we present the main stability result of this paper, \il{and} the relevant conditions imposed to guarantee it. In Section \ref{sec:Discussion} we discuss various applications of the proposed results and clarify the contribution of our analysis by comparing it with existing literature. Our results are numerically validated in Section \ref{Simulation_NPCC}, where convergence to nominal frequency and an economically optimal allocation are demonstrated at the presence of high order turbine governor dynamics. Finally, conclusions are drawn in Section \ref{Conclusions}.

\section{Notation}\label{sec: Preliminaries}

 Real numbers are denoted by $\mathbb{R}$, and the set of $n$-dimensional vectors with real entries is denoted by $\mathbb{R}^n$.
For a function $f(q)$, $f:\mathbb{R}\rightarrow \mathbb{R}$, we denote its first derivative by $f'(q) = \tfrac{d}{dq} f(q)${\color{black}, its inverse by $f^{-1}(.)$.
A function $f:\mathbb{R}^n \rightarrow \mathbb{R}$ is said to be positive definite if $f(0) = 0$ and $f(x) > 0$ for every $x\neq0$.
 We write $\vect{0}_n$ and $\vect{1}_n$ to denote $n \times 1$ vectors with all elements equal to $0$ and $1$ respectively. We denote by $\mathcal{I}_k$ for some $k>0$ the identity matrix of rank $k$.
For a discrete set $S$, the term $|S|$ denotes its cardinality.
A matrix $A$ is said to be Hurwitz if all its eigenvalues lie on the open left half plane \cite{zhou1998essentials}. \ak{A matrix $A$ is said to be positive definite (semi-definite) when $x^TAx > 0$ for all $x \neq 0$ (when $x^TAx \geq 0$ for all $x$).}
Finally, $\ak{\Ima(A)}$ denotes the range of a matrix $A \in \mathbb{R}^{m \times n}$.

\section{Power network dynamics} \label{Problem_Formulation}

In this section, we present a mathematical description of the power network and a general class of linear generation dynamics that is considered within the rest of the manuscript.
Furthermore, we discuss how stability considerations within \icl{the} primary frequency control timeframe can affect the choice of generation droop coefficients, following existing results \icl{in the} literature. The latter will be \ank{associated with the secondary frequency regulation analysis}
presented later in this manuscript.

\subsection{Network model}\label{Network_Model}

 We use a connected graph $(N,E)$ to describe the power network, where $N = \{1,2,\dots,|N|\}$ is the set of buses and $E \subseteq N \times N$ the set of transmission lines connecting the buses.
There are two types of buses in the network, buses with inertia and buses without inertia.
Since generators have inertia, it is reasonable to assume that only buses with inertia have non-trivial generation dynamics. We define $G = \{1,2,\dots,|G|\}$ and $L= \{|G|+1,\dots,|N|\}$ as the sets of buses with and without inertia respectively such that $|G| + |L| = |N|$.
Moreover, the term $(i,j)$ denotes the link connecting buses $i$ and $j$. The graph $(N,E)$ is assumed to be directed with an arbitrary direction, so that if $(i,j) \in E$ then $(j,i) \notin E$. Additionally, for each $j \in N$, we use $i:i\rightarrow j$ and $k:j\rightarrow k$ to denote the sets of buses that precede and succeed bus~$j$ respectively. It should be noted that the form of the dynamics in~\eqref{sys1}--\eqref{dynsys} below is not affected by changes in graph ordering, and our results are independent of the choice of direction.  We make the following assumptions for the network: \newline
1) Bus voltage magnitudes are $|V_j| = 1$ p.u. for all $j \in N$.

2) Lines $(i,j) \in E$ are lossless and characterized by their susceptances $B_{ij} = B_{ji} > 0$. \newline
3) Reactive power flows do not affect bus voltage phase angles and frequencies. \newline

Following the above, we can make use of swing equations to describe the rate of change of frequency at generation buses. Moreover, power must be conserved at each of the load buses.
 This motivates the following system dynamics (e.g.~\cite{Bergen_Vittal}),
\begin{subequations} \label{sys1}
\begin{equation}
\dot{\eta}_{ij} = \omega_i - \omega_j, \; (i,j) \in E, \label{sys1a}
\end{equation}
\begin{equation}
M_j \dot{\omega}_j = - p_j ^L + p_j^M - \Lambda_j\omega_j - \sum_{k:j\rightarrow k} p_{jk} + \sum_{i:i\rightarrow j} p_{ij}, \; j\in G, \label{sys1b}
\end{equation}
\begin{equation}
 0 = - p_j ^L -  \Lambda_j\omega_j - \sum_{k:j\rightarrow k} p_{jk} + \sum_{i:i\rightarrow j} p_{ij}, \; j\in L, \label{sys1c}
\end{equation}
\begin{equation}
 p_{ij}=B_{ij} \sin\eta _{ij}, \; (i,j) \in E. \label{sys1d}
\end{equation}
\end{subequations}
In system~\eqref{sys1}, the time-dependent variable $\omega_j$ represents the deviation of the frequency at bus $j$ from its nominal value, namely 50Hz (or 60Hz). The time dependent variables $p^M_j$ represent
the mechanical power injection to the generation bus $j$.
 The constant $\Lambda_j$ represents the frequency damping coefficient at bus $j$.
 The time-dependent variables $\eta_{ij}$ and $p_{ij}$ represent, respectively, the power angle difference\footnote{The quantities $\eta_{ij}$ represent the phase differences between buses $i$ and $j$, given by $\theta_i - \theta_j$, i.e. $\eta_{ij} = \theta_i - \theta_j$. The angles themselves must also satisfy $\dot{\theta}_j = \omega_j$ at all $j \in N$. This equation is omitted in \eqref{sys1} since the power transfers are functions of the phase differences only.}
 and the power transferred from bus $i$ to bus $j$.
  The constant $M_j > 0$ denotes the generator inertia.
  The response of the system~\eqref{sys1} will be studied, when a step change $p_j^L$,
  $j \in N$ {occurs in} the uncontrollable demand.

\subsection{Generation Dynamics}

To investigate broad classes of dynamics and control policies for generation systems,
we consider dynamical systems of the form
\vspace{-.1cm}
\begin{equation} \label{dynsys}
\begin{aligned}
&\dot{x}^M_j = A_jx^M_j + B_ju_j,\\
&p^M_j = C_jx^M_j + D_ju_j,
\end{aligned}, j \in G
\end{equation}
with input $u_j(t) \in \mathbb{R}$, state $x^M_j(t) \in \mathbb{R}^{n_j}$,
output $p^M_j(t) \in \mathbb{R}$ and corresponding matrices $A_j \in \mathbb{R}^{n_j \times n_j}, B_j \in \mathbb{R}^{n_j}_j, C_j \in \mathbb{R}^{1 \times n_j}$ and $D_j \in \mathbb{R}$.
We assume in~\eqref{dynsys} that $A_j$ is Hurwitz which implies
that given any constant input \ak{$u_j(t) = \bar{u}_j$}, there exists an
 asymptotically stable equilibrium point $\bar{x}^M_j \in \mathbb{R}^{n_j}$, such that $A_j\bar{x}^M_j + B_j\bar{u}_j = 0$.
 Correspondingly, there exists a constant
$K_j \in \mathbb{R}$, satisfying $K_j = -C_j A_j^{-1}B_j + D_j$, such that for any constant input $\bar{u}_j$ and corresponding state $\bar{x}_j$, the output $\bar{p}^M_j$ is given by
\begin{equation}
\bar{p}^M_j = C_j \bar{x}^M_j + D_j \bar{u}_j = K_j \bar{u}_j.
\label{iochar}
\end{equation}

Our aim in this paper is to provide conditions on the dynamics described in \eqref{dynsys} and the choice of input $u_j$ that allow for stability and optimality guarantees and ensure the satisfaction of  secondary frequency control objectives.

\subsection{Stability in primary frequency control}\label{Primary_control_convergence}

The main objective of primary frequency regulation is to balance generation and demand.
Within this subsection, we will consider conditions imposed in \icl{the} literature to provide decentralized stability guarantees for primary frequency regulation. These conditions will later be extended when dynamics associated with secondary frequency control are taken into account.

 We let the generation input $u_j$ be described by
\begin{equation}\label{primary_input}
u_j = - k_{d,j} \omega_j, j \in G
\end{equation}
where $k_{d,j}$ are positive constants related to the droop gains.
To obtain decentralized conditions on the design parameters $k_{d,j}$ and generation dynamics \eqref{dynsys} such that stability is guaranteed, we will make use of the results presented in \cite{kasis2016primary}.
From the analysis in \cite{kasis2016primary}, the following \ak{input strict passivity assumption} should hold for the system \eqref{dynsys},\eqref{primary_input} when some frequency damping term $\Lambda_j \omega_j$ is present\footnote{Note that the condition in \cite{kasis2016primary} also applies to non linear systems and is defined locally. However, here it is simplified to fit the purpose of this study.}.

\ank{
\begin{assumption}\label{assum_passivity}
Consider the systems with input $-\omega_j$ and output $p^M_j - \Lambda_j \omega_j$ described by \eqref{sys1},\eqref{dynsys},\eqref{primary_input}.
Then, for each of these systems, there exists a symmetric
positive definite matrix $P_j = P^T_j \in \mathbb{R}^{n_j \times n_j}$ such that the functions \ank{$V_j(x^M_j) = \frac{1}{2}(x^M_j)^T P_j x^M_j$} satisfy

\begin{equation}\label{passivity_condition}
\dot{V}_j \leq -\omega_j p_j^M  - \kas{\mu_j} (-\omega_j)^2,
\end{equation}
for some \kas{$\mu_j > 0$}, for all $j \in G$.
\end{assumption}
}

Assumption \ref{assum_passivity} follows directly from the conditions imposed in \cite{kasis2016primary}. It has the important property that it can be efficiently verified with appropriate LMI conditions that follow from the KYP lemma \cite{khalil1996nonlinear}. In particular, it can be verified that Assumption \ref{assum_passivity} holds for given gain $k_{d,j}$ when there exists a positive definite symmetric matrix $P_j$, such that
\begin{equation}\label{primary_LMI}
\hat{Q}_j = \begin{bmatrix}
\frac{P_jA_j + A_j^TP_j}{2} & \frac{k_{d,j}P_jB_j - C_j^T}{2} \\
\frac{k_{d,j}B^T_jP_j - C_j}{2}  & -\hat{\Lambda}_j - D_jk_{d,j}
\end{bmatrix} \leq 0
\end{equation}
holds for some $\hat{\Lambda}_j < \Lambda_j$.

We will see in the following sections how the condition in \eqref{primary_LMI} extends when dynamics related with secondary frequency regulation are also included.  This will demonstrate that some of the imposed decentralized conditions on stability for secondary frequency regulation can be seen as existing stability requirements  for primary frequency control.

\section{Distributed averaging dynamic output controller}\label{sec:Power_Command_Dynamics}

In this section we propose a novel \ank{secondary} control scheme, called distributed averaging dynamic output controller (DADOC) which, as discussed in Section \ref{sec: Comparison}, offers advantages over existing distributed schemes for secondary frequency regulation.

We consider a communication network described by a connected graph $(G, \tilde{E})$. We \ak{propose} the following DADOC scheme
\begin{equation}\label{p_command}
\gamma_{j} \dot{p}^c_j = p^M_j - K_j u_j - k_{f,j} \omega_j + \hspace{-3mm}\sum_{i: (i,j) \in \tilde{E}}\hspace{-3mm} \alpha_{ij} (p^c_i - p^c_j),  j \in G
\end{equation}
where $\alpha_{ij} = \alpha_{ji}, \gamma_j$ and $k_{f,j}$ are positive constant gains of the controller, and $p^c_j$ is a power command signal.
\icl{The generation input in \eqref{primary_input} is also extended to}
\begin{equation}\label{input_description}
u_j = k_{c,j}p^c_j -k_{d,j}\omega_j, j \in G
\end{equation}
where $k_{c,j}, k_{d,j}$ are positive \kas{design} constants.

 The controller dynamics in \eqref{p_command} ensure that power command variables synchronize at steady state, a necessary feature \ak{for} optimality interpretation.
Furthermore, the frequency dependent term in \eqref{p_command} ensures that when generation and power command variables reach steady state and power command variables synchronize, then \icl{the} steady state frequency must attain its nominal value. The term $p^M_j - K_j u_j$ in \eqref{p_command} does not affect the steady state of the power command variables. However, as shall be seen in the subsequent analysis, it has a pivotal role in providing stability guarantees when high order generation dynamics of the form \eqref{dynsys} are considered.

We choose the generation input $u_j$ to be a weighted sum of frequency and power command, allowing the weight coefficients \kas{to} be design parameters. In what follows,  we provide conditions on how these parameters should be chosen such that convergence to nominal frequency is achieved while also taking economic considerations into account.

\subsection{Equilibrium analysis}

We now describe what is meant by an equilibrium
of the interconnected system~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}.

\begin{definition} \label{eqbrdef}
The point $\beta^* = (\eta^*, \omega^*,x^{M,*}, p^{c,*})$ defines an equilibrium of the system~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description} if all time derivatives of ~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description} are equal to zero at this point.
\end{definition}

In our analysis, we shall consider conditions on  controller design variables such that the generation equilibrium values solve an economic optimisation problem, as will be described in Section \ref{sec:optim}.

Throughout the paper, it is assumed that there exists some equilibrium of~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}, denoted by $\beta^* = (\eta^*, \omega^*, x^{M,*}, p^{c,*})$, as defined in Definition~\ref{eqbrdef}.
Note that the study on the existence of equilibria is beyond the scope of this paper and the interested reader is referred to \cite{dorfler2013synchronization}.
Furthermore, we use $(p^*, p^{M,*}, u^*)$ to represent the equilibrium values of respective quantities
in~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}.

The following lemma, proven in the appendix, \ak{characterizes the equilibria of the system~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}. It demonstrates that \il{the} frequency attains its nominal value at steady state and that \il{the} power command variables synchronize, a useful property when an optimality interpretation of the equilibria is desired.}

 \begin{lemma}\label{eqbr_lemma}
Any equilibrium point $\beta^*$ given by Definition \ref{eqbrdef} satisfies $\omega^{*} = \vect{0}_{|N|}$ and $p^{c,*} \in \ak{\Ima}(1_{|G|})$.
 \end{lemma}

Furthermore, the power angle differences at the considered equilibrium are assumed to satisfy the following security constraint. It will be seen that Assumption \ref{assum_angles} is required for the convergence results presented in Section \ref{sec:Convergence}.

\begin{assumption} \label{assum_angles}
$| \eta^*_{ij} | < \tfrac{\pi}{2}$ for all $(i,j) \in E$.
\end{assumption}

The stability and optimality properties of such equilibria will be studied in the following sections.

\subsection{Optimality analysis} \label{sec:optim}

We aim to study how generation should be adjusted in order to meet the step change in frequency independent demand and simultaneously minimize the cost that comes from the deviation in the power generated. We now introduce an optimization problem, which we call the optimal generation regulation problem (OGR), that can be used to achieve this goal.

A quadratic cost is supposed to be incurred when \ank{the generation output at bus $j$ is $p_j^M$.}

Note that quadratic cost functions are frequently used in the literature, \cite{zhaomallada, lavaei2012zero}, motivated by the fact that a convex function can be locally approximated by a quadratic one.
 The problem is to find the vector $p^M$ that minimizes this total cost and simultaneously achieves power balance. More precisely, the  following optimization problem is considered

\begin{equation} \label{Problem_To_Min}
\begin{aligned}
&\hspace{-2em}\underline{\text{OGR:}} \\
&\min_{p^M} \sum\limits_{j\in G} \frac{1}{2} q_j (p_{j}^M)^2,  \hspace{-1.5em}\\
&\text{subject to } \sum\limits_{j\in  G} p_j^M = \sum\limits_{j\in  N}  p_j^L,
\end{aligned}
\end{equation}
\ank{where $q_j$ are positive cost coefficients associated with the generation cost at bus $j$.
Note that more general quadratic cost functions can be considered, following similar approaches as in relevant literature (e.g. \cite{trip2017distributed}). However, we opt for the cost functions in \eqref{Problem_To_Min} for simplicity.}
The equality constraint in~\eqref{Problem_To_Min} requires all the frequency-independent \il{loads} to be matched by the total generation.
This ensures that when system \eqref{sys1} is at equilibrium, then frequency will be at its nominal value, as follows from summing \eqref{sys1b}--\eqref{sys1c} at equilibrium over all buses.

Within the paper, we aim to specify properties on the control dynamics of $p^M$, described in \eqref{dynsys}, that ensure that those quantities converge to values at which optimality can be guaranteed for \eqref{Problem_To_Min}.
Below, we demonstrate how the controller gains in \eqref{p_command}--\eqref{input_description} can be chosen to ensure that the equilibrium points of the system are solutions\footnote{Note that an equilibrium point is a solution to the OGR problem when at that point \ank{the value of $p^M$ is optimal for \eqref{Problem_To_Min}.}} to the OGR problem \eqref{Problem_To_Min}. We will then demonstrate in the subsequent section, how convergence to optimality can be achieved.

\begin{proposition} \label{optthm}
Consider equilibria of \eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description},  characterized by Lemma \ref{eqbr_lemma}. Then, if the control dynamics in \eqref{input_description} are chosen such that
\begin{equation}\label{opt_cond}
k_{c,j} = \frac{1}{q_j K_{j}}, j \in G,
\end{equation}
holds, then the equilibrium values $p^{M,*}$ are optimal for the OGR problem \eqref{Problem_To_Min}.
\end{proposition}

\begin{remark}
\ak{Proposition} \ref{optthm} provides conditions on the choice of the design variable $k_{c,j}$ such that the equilibria of system \eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}, described by Lemma \ref{eqbr_lemma} are solutions to the OGR problem \eqref{Problem_To_Min}. Note that the design variables $k_{f,j}$ and $k_{d,j}$ do not appear in the optimisation problem, since those are gains on frequency deviation, which becomes zero at equilibrium. However, it will be seen in the following section that their values have a significant impact on the stability properties of the system.
\end{remark}

\section{Stability analysis}\label{sec:Convergence}

This section contains our main convergence results. We provide appropriate conditions on the choice of gains in \eqref{p_command}--\eqref{input_description}, applicable to highly relevant generation schemes,
and show that when those are satisfied, then convergence is guaranteed.

\subsection{Controller design conditions}

In this section, we impose a condition involving design constants $k_{c,j}, k_{d,j}$ and $k_{f,j}$, which is used in the convergence theorem presented in Section \ref{sec:Main_result} below. We then explain how this condition can be numerically tested in a computationally efficient way. The considered condition is presented below.

\begin{design condition}\label{assum_gains}
For each generation bus $j$, with dynamics described by \eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}, the controller parameters $k_{c,j}, k_{d,j}$ and $k_{f,j}$ are such \ak{that
\begin{equation}\label{strictness_condition}
\begin{bmatrix}

-K_jk_{c,j} + D_j k_{c,j} & r^T_j \\
r_j & \hat{Q}_j \\
\end{bmatrix}  \leq 0,
\end{equation}
where $r_j = \begin{bmatrix}
\frac{k_{c,j}B_j^TP_j + C_j}{2}  & \frac{k_{f,j}-k_{d,j}K_j + D_j k_{d,j} - D_jk_{c,j}}{2}
\end{bmatrix}^T$,}  holds for some $P_j = P_j^T > 0$ and some $\hat{\Lambda}_j < \Lambda_j$.

\end{design condition}

Design condition \ref{assum_gains} is the main stability condition imposed on this paper, and is feasible for a broad class of linear systems, as discussed in Section \ref{sec:Discussion}.

Design condition  \ref{assum_gains} can be interpreted as an extension of Assumption \ref{assum_passivity} to secondary frequency control. This is since Assumption \ref{assum_passivity} is a necessary condition for Design condition \ref{assum_gains} to hold, as follows from noting that the matrix $\hat{Q}_j$ in \eqref{primary_LMI} is a \ak{principal} submatrix of the one considered in \eqref{strictness_condition}.  Therefore, part of the stability conditions imposed for secondary frequency regulation can be seen as conditions for stability in primary frequency control.

\ank{
\begin{remark}
The inequality condition  in \eqref{strictness_condition} is an LMI with respect to the matrix $P_j$ and  design variables $k_{f,j}$ and can therefore be verified in a computationally efficient way. Furthermore, the condition in \eqref{strictness_condition} can be used to formulate various optimization problems that may make use of the flexibility in choosing the matrix $P_j$ and the design variable $k_{f,j}$. An example of such problem would be to obtain the minimum frequency damping $\Lambda_j$ that is required for \eqref{strictness_condition} to hold for particular generation dynamics. Hence, \eqref{strictness_condition} can be also useful in system design.
\end{remark}
}

\subsection{Main result}\label{sec:Main_result}

We are now in a position to state our main result, demonstrating convergence to an optimal point of \eqref{Problem_To_Min} where frequency attains its nominal value.

\begin{theorem}\label{combined_thm}
\ank{Consider an equilibrium of~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description} such that Assumption~\ref{assum_angles} holds and let Design condition \ref{assum_gains} and \eqref{opt_cond} be satisfied.
Then, there exists an open neighborhood of initial conditions about this equilibrium such that the solutions of~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}  asymptotically converge to a set of equilibria that solve the OGR problem \eqref{Problem_To_Min} with $\omega^* = \vect{0}_{|N|}$.}
\end{theorem}

\ak{Theorem \ref{combined_thm} demonstrates local convergence to an optimal solution of the OGR problem \eqref{Problem_To_Min} that, following Lemma \ref{eqbr_lemma}, satisfies $\omega^* = \vect{0}_{|N|}$.
The main conditions for stability are Assumption~\ref{assum_angles}, which is abundant in power literature, and Design condition \ref{assum_gains}.  In the following section, we demonstrate the relevance of Design condition \ref{assum_gains}, explaining how it applies to various generator models.
}

\section{Discussion}\label{sec:Discussion}

In this section we demonstrate the \ak{applicability} of our main results with examples of highly relevant generation schemes that fit within the considered analysis. Moreover, we explain the contribution of our results relative to the existing literature.

\subsection{Applications of main results}\label{sec:Applications}

To demonstrate the \ak{relevance} of our analysis, we provide examples of first and second order turbine governor dynamics and explain how \ak{our proposed conditions apply to them.}

Consider first order generation dynamics described by
\begin{equation}\label{first_order}
\tau_j \dot{p}^M_j = -p^M_j + K_j(k_{c,j}p^c_j - k_{d,j}\omega_j),
\end{equation}
for some constant $\tau_j > 0$, coupled with the controller \eqref{p_command}--\eqref{input_description} and some frequency damping $\Lambda_j$. For this system, it can be shown that for any positive values for $\tau_j, q_j, K_j, k_{d,j}$ and $\Lambda_j$
 there always exist positive constants $k_{f,j}, k_{c,j}$ and some $P_j = P^T_j > 0$  such that both Design condition \ref{assum_gains} and optimality condition \eqref{opt_cond} are satisfied.

An important aspect of the proposed framework is that it allows the inclusion of high order generation dynamics. A significant example of such, is the following second order model describing turbine governor dynamics (e.g. \cite{Bergen_Vittal}),
\begin{subequations}\label{second_order}
\begin{equation}
 \dot{\alpha}_j = -\frac{1}{\tau_{a,j}}(\alpha_j - K_j(k_{c,j}p^c_j - k_{d,j}\omega_j)),
\end{equation}
\begin{equation}
 \dot{p}^M_j = -\frac{1}{\tau_{p,j}}(\ank{p^M_j} - \alpha_j),
\end{equation}
\end{subequations}

where $\alpha_j$ is the internal state of the model and $\tau_{a,j}, \tau_{b,j} > 0$ are time constants associated with the generation dynamics. We considered the case where \eqref{second_order} is coupled with the power command dynamics described by \eqref{p_command} and some frequency damping $\Lambda_j$.
\ak{The following lemma, proven in the appendix, provides a sufficient condition for the frequency damping $\Lambda_j$ such that Design condition \ref{assum_gains} holds for the considered system.

\begin{lemma}\label{suff_cond_lemma}
Consider a bus $j$ with generation dynamics described by \eqref{second_order} coupled with the power command dynamics \eqref{p_command} and \ank{some frequency damping $\Lambda_j$. Then, Design condition \ref{assum_gains} holds for any positive values of $\tau_{a,j}, \tau_{p,j}$ if  \[\Lambda_j > \frac{K_j}{3k_{c,j}}(k_{c,j}^2 - k_{c,j}k_{d,j} + k_{d,j}^2))\] holds. }

\end{lemma}

Lemma \ref{suff_cond_lemma} provides a sufficient condition for the value of frequency damping $\Lambda_j$ such that Design condition \ref{assum_gains} holds at a bus with generation and power command dynamics described by \eqref{second_order} and \eqref{p_command} respectively. Note that this condition can \ank{be made less conservative} when particular values for $\tau_{a,j}, \tau_{p,j}$ are considered.}

Note that numerical analysis demonstrates that our proposed results also apply to higher than second order generation schemes. However, these results are omitted due to space constraints.

\subsection{Comparison with existing literature}\label{sec: Comparison}

The problem of addressing issues of stability and optimality for secondary frequency regulation in a distributed way has been widely studied over the last years. Most studies focused on two particular schemes which are briefly discussed below.

A trend in secondary frequency control literature is to consider a distributed averaging proportional integral (DAPI) controller \cite{trip2016internal, dorfler, dorfler2014plug, simpson2013synchronization,  andreasson2014distributed}. The benefits of this scheme lie in its simplicity, since further than exchanging a synchronizing variable, it only requires knowledge of the local frequency which is easily obtainable. However, existing stability results along this setting are limited to static generation models and quadratic cost functions.

Another approach in literature is to consider power command dynamics that follow from a primal/dual algorithm associated with some optimization problem \cite{P30_Low}, \cite{kasis2017secondary}, \cite{low2014distributed}, \cite{zhao2015distributed}, \cite{li2016connecting},  which ensures secondary frequency control objectives are met. For convenience, we shall refer to those
schemes as 'Primal-Dual'.
Primal-Dual schemes have the significant advantage that they allow for stability guarantees when high order and non linear generation dynamics are included. Furthermore, they allow for economic optimality when general convex cost functions are considered.
 However, these controllers are more complicated than DAPI schemes, requiring knowledge of demand, which can be difficult to obtain in many cases, and generation.

DADOC schemes share benefits of both mentioned control schemes, allowing for stability guarantees when high order generation dynamics are considered and requiring measurements of frequency and generation output, which are not restrictive to obtain.
In particular, DADOC schemes have the advantage over DAPI schemes that they allow the inclusion of high order generation dynamics.
This comes in the expense of an additional requirement for generation output measurements.
Both Primal-Dual and DADOC schemes allow for stability and optimality guarantees when high order generation dynamics and quadratic cost functions are considered.
However, Primal-Dual schemes additionally allow to incorporate non linear generation dynamics and general convex cost functions.
The main disadvantage of Primal-Dual schemes is their requirement for measurements of local demand, which can \kas{be} difficult to obtain.
On the other hand, DADOC schemes allow for stability guarantees for high order generation dynamics without requiring knowledge of demand.
It is important to also note that DADOC schemes allow the inclusion of \il{classes of}
generation dynamics that cannot be incorporated when a Primal-Dual controller is implemented. One such example is described by \eqref{second_order} where, unlike when Primal-Dual schemes are considered, no static dependence on power command is required.
Table \ref{table_comparison} summarizes the comparison between DADOC, DAPI and Primal-Dual schemes.

\kas{A further feature of DADOC schemes is that the required equilibrium condition for frequency, i.e. that $\omega_j^* = 0, j \in N$, follows from the dynamics at each controller. This is important when a controller is withdrawn from the network due to a failure, since the equilibrium condition $\omega^* = 0$ will still hold. By contrast, the equilibrium condition for frequency when Primal-Dual schemes are considered follows from the aggregate controller dynamics, and is not in general valid when a single controller fails. Hence, DADOC schemes have better robustness to controller failure over Primal-Dual schemes.}

\begin{table}
\resizebox{0.49\textwidth}{!}{
 \begin{tabular}
 {|c | c | c | c |}
 \hline
  & DAPI & DADOC & Primal-Dual \\ [0.5ex]
 \hline
 Allowable cost  & Quadratic & Quadratic & Convex \\
 function models & & & \\[0.5ex]
 \hline
 Allowable generation & Static & High Order & High Order \\
 dynamics & & & and Non linear \\[0.5ex]
 \hline
 Required & Frequency & Frequency  & Generation \\
 measurements & & Generation & Demand \\[0.5ex]
 \hline
\end{tabular}
}
 \caption{Comparison between our proposed (DADOC) controller and the two dominant schemes in literature on allowable cost function models and generation dynamics and required measurements for their implementation.}
 \label{table_comparison}
 \vspace{-4mm}
\end{table}

\ak{A further attempt to address distributed secondary frequency regulation issues has been made in \cite{trip2017distributed}, where the authors propose a controller which allows for stability guarantees when first or second order generation dynamics are considered, by making use of measurements of generators internal states. The main differences of DADOC schemes compared to those in \cite{trip2017distributed}, are the different measurement requirements for their implementation (generation output instead of internal states) and that the former  \ank{provides design conditions for power networks with} higher than second order generation dynamics.}

\section{Simulation on the NPCC 140-bus system} \label{Simulation_NPCC}

In this section we verify our analytic results with numerical simulations
on the Northeast Power Coordinating
Council (NPCC) 140-bus interconnection system, performed using the Power System Toolbox~\cite{cheung2009power}. This model is more detailed and realistic than our analytical one, including line resistances, a DC12 exciter model, a subtransient reactance generator model, and {higher} order turbine governor models\footnote{The details of the simulation models can be found in the Power System {Toolbox \cite{cheung2009power} manual and} data file datanp48.}.

The test system consists of 93 load buses serving different types of loads including constant active loads and 47 generation buses. The overall system has a total real power of 28.55GW. For our simulation, we added three loads on units 2, 9, and 17, each having a step increase of magnitude $1$ p.u. (base 100MVA) at $t=1$ second.

To demonstrate the applicability of the proposed controller, the dynamics in \eqref{p_command}--\eqref{input_description} where implemented on 16 generators with third, fourth and fifth order turbine governor dynamics.
Furthermore,  in order to verify our optimality analysis, we considered quadratic cost functions with cost coefficients equal to $K^{-1}_j$ that penalised deviations on generation outputs. The choice of cost coefficients relates high cost coefficients with small droop gains, and is consistent with current results in literature on optimal frequency regulation \cite{kasis2016primary}. However, it has been numerically verified that the  stability and optimality properties of the system, demonstrated below, are retained for a broad range of cost coefficient values.
Controller parameters were selected such that Design condition  \ref{assum_gains} and optimality codition \eqref{opt_cond} were satisfied.

\begin{figure}[t]
\centering
\includegraphics[trim = 1mm 0mm 10mm 0mm, height = 2.2in,width=3.115in,clip=true]{Frequency.eps}

\vspace{-1mm}
\caption{Frequency at bus 101 of the NPCC network.}
\label{Frequency}
\vspace{-3mm}
\end{figure}

 The frequency response of a randomly selected bus is depicted on Figure \ref{Frequency}. There, it is demonstrated that \il{the} frequency returns to its nominal value, hence numerically validating the analytic convergence results of Theorem \ref{combined_thm}.

\begin{figure}[t]
\centering
\includegraphics[trim = 1mm 0mm 10mm 0mm, height = 2.2in,width=3.115in,clip=true]{Marginal.eps}

\vspace{-1mm}
\caption{Marginal cost of generation buses contributing to secondary frequency control.}
\label{Marginal}
\vspace{-3mm}
\end{figure}

Furthermore, from Fig. \ref{Marginal}, it is observed
that the marginal costs at all 16 generators that contribute to secondary frequency control converge to the same value.  This illustrates the optimality in the power allocation
among generators, since equality in the marginal cost is necessary to solve \eqref{Problem_To_Min}.

\vspace{-2mm}
\section{Conclusion}\label{Conclusions}

We have considered the problem of designing distributed control schemes such that stability and an optimal power allocation can be guaranteed for secondary frequency control.
We proposed a distributed averaging dynamic output controller (DADOC) which ensures stability and optimality  when high order linear generation dynamics and quadratic cost functions are considered and also that frequency attains its nominal value at steady state.
DADOC controllers are advantageous compared to existing schemes since they allow for stability \ank{and optimality} guarantees for a broad class of generation dynamics and are also easy to implement, in terms of measurement requirements.
Moreover, we demonstrated the applicability of our results with examples \il{of
 high} order generation dynamics and explained how the design parameters can be selected in a computationally efficient way by appropriate LMI conditions. Our results have been numerically verified with simulations on the NPCC 140-bus system. Interesting potential extensions in the analysis include considering controllable loads \il{in} the proposed framework, incorporating non linear generation dynamics and general convex cost functions, and extending the analysis to include voltage dynamics.

\appendix

This appendix contains the proofs of all the results presented in this paper.

Throughout the proofs we will make use of the following equilibrium equations for the dynamics in \eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description},
\begin{subequations} \label{eqbr}
\begin{equation}
0 = \omega^*_i - \omega^*_j, \; (i,j) \in E, \label{eqbr1}
\end{equation}
\begin{equation}
0 = - p_j^L + p_j^{M,*} - \Lambda_j \omega^*_j - \sum_{k:j\rightarrow k} p^*_{jk} + \sum_{i:i\rightarrow j} p^*_{ij}, \; j\in G, \label{eqbr2}
\end{equation}
\begin{equation}
0 = - p_j^L - \Lambda_j \omega^*_j - \sum_{k:j\rightarrow k} p^*_{jk} + \sum_{i:i\rightarrow j} p^*_{ij}, \; j\in L, \label{eqbr3}
\end{equation}
\begin{equation}\label{eqbr4}
p^{*}_{ij} = B_{ij} \sin \eta_{ij}^*, \; (i,j) \in E,
\end{equation}
\begin{equation}
p^{M,*}_j = K_{j} u^*_j, \; j \in G, \label{eqbr5}
\end{equation}
\begin{equation}\label{eqbr6}
0 = p^{M,*}_j - K_{j} u^*_j - k_{f,j} \omega^*_j  + \hspace{-3mm}\sum_{i: (i,j) \in \tilde{E}}\hspace{-3mm} \alpha_{ij} (p^{c,*}_i - p^{c,*}_j), j \in G
\end{equation}
\begin{equation}\label{eqbr7}
u^*_j =
k_{c,j} p^{c,*}_j -k_{d,j} \omega^*_j, j \in G
\end{equation}
\end{subequations}

\emph{Proof of Lemma~\ref{eqbr_lemma}:}
From \eqref{eqbr1} it follows that \kas{$\omega^*_i = \omega^*_j$ for all $(i,j) \in E$, which results to $\omega^* \in \Ima(\vect{1}_{|N|})$}.
Then, summing \eqref{eqbr6} over all $j \in G$ results to
$
 \sum_{j \in G} p^{M,*}_j -  (K_ju^*_j + k_{f,j}\omega^*_j) = 0
$,
which by \eqref{eqbr5} and $k_{f,j} > 0$ implies that $\omega^* = \vect{0}_{|N|}$.
Since $\omega^*_j = 0, j \in N$, it follows by \eqref{eqbr5}, \eqref{eqbr6} and the fact that the communication graph is connected that $p^{c,*} \in \ak{\Ima}(\vect{1}_{|G|})$.
\hfill $\blacksquare$

\emph{Proof of \ak{Proposition}~\ref{optthm}:}
The OGR optimization problem~\eqref{Problem_To_Min} is convex and has a continuously differentiable cost function.
 Thus, a point $\bar{p}^M$ is a global minimum for~\eqref{Problem_To_Min} if and only if it satisfies the KKT conditions~\cite{boyd2004convex}
\begin{subequations} \label{kkt}
\begin{gather}
q_j \bar{p}^M_j = \nu, \; j \in G, \label{kkt1} \\
\sum\limits_{j\in  G} \bar{p}_j^M = \sum\limits_{j\in  N} p_j^L, \label{kkt2}
\end{gather}
\end{subequations}
for some constant $\nu \in \mathbb{R}$. It will be shown below that these conditions are satisfied by the equilibrium values $\bar{p}^M = p^{M,*}$ defined by equations~\eqref{eqbr5} and~\eqref{eqbr7} when \eqref{opt_cond} holds.

From Lemma \ref{eqbr_lemma}, it follows that $\omega^* = \vect{0}_{|N|}$ and $p^{c,*} \in \ak{\Ima}(\vect{1}_{|G|})$.
Then, let $\nu = p^{c,*}_j$ and note that is common at every bus since power command variables synchronize at steady state.
 Therefore, it follows that $(q_j)^{-1} \nu  = (q_j)^{-1} p^{c,*}_j = (q_j K_j k_{c,j})^{-1}p^{M,*}_j = p^{M,*}_j$, by $\omega^* = \vect{0}_{|N|}$  and equations~\eqref{eqbr5}, \eqref{eqbr7} and~\eqref{opt_cond}. Thus, the optimality condition~\eqref{kkt1} holds.

Summing equations~\eqref{eqbr2} and~\eqref{eqbr3} over all $j \in G$ and $j \in L$ respectively and using that $\omega^* = \vect{0}_{|N|}$ shows that~\eqref{kkt2} also holds.

Hence, the values $\bar{p}^M = p^{M,*}$ satisfy the KKT conditions~\eqref{kkt}. Therefore, the equilibrium values $p^{M,*}$ define a global minimum for~\eqref{Problem_To_Min}.
\hfill $\blacksquare$

\emph{Proof of Theorem~\ref{combined_thm}:}
We will use the dynamics in \eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description} and the matrices $P_j$ in Design condition \ref{assum_gains} to define a Lyapunov function for the system~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}.

Firstly, let $V_F (\omega^G) = \frac{1}{2}\sum_{j \in G} M_j (\omega_j - \omega^*_j)^2$. The time-derivative of $V_F$ along the trajectories of~\eqref{sys1}--\eqref{dynsys} is given by
\begin{align*}
\hspace{-0.5mm}\dot{V}_F &=\hspace{-1.5mm} \sum_{j \in N}\hspace{-0.5mm}(\omega_j - \omega^*_j) \Bigg(\hspace{-2.0mm}-p^L_j + p^M_j - \Lambda_j\omega_j -\hspace{-2.5mm} \sum_{k:j\rightarrow k} p_{jk} + \hspace{-2.0mm}\sum_{i:i\rightarrow j} p_{ij}\hspace{-0.5mm}\Bigg),
\end{align*}
by substituting~\eqref{sys1b} for $\dot{\omega}_j$ for $j \in G$ and adding extra terms for $j \in L$, which are equal to zero by~\eqref{sys1c}. Subtracting the product of $(\omega_j - \omega^*_j)$ with each term in~\eqref{eqbr2} and~\eqref{eqbr3}, this becomes
\begin{align}
\dot{V}_F =& \sum_{j \in G} \hspace{-0.5mm}(\omega_j - \omega^*_j) (p^M_j - p^{M,*}_j)
\hspace{-0.5mm} - \hspace{-0.2em} \sum_{j \in N}  \Lambda_j(\omega_j - \omega^*_j)^2\hspace{-0.5mm}  \nonumber \\
& +
\sum_{(i,j) \in E} (p_{ij} - p^*_{ij}) (\omega_j - \omega_i), \label{VFdiff}
\end{align}
using the equilibrium condition~\eqref{eqbr1} for the final term.

Furthermore, let $V_C(p^c) = \frac{1}{2}\sum_{j \in N} \gamma_j(p^c_j - p^{c,*}_j)^2$. Using \eqref{p_command} the time derivative of $V_C$ can be written as
\begin{multline}
\dot{V}_C = \sum_{j \in N}  (p^c_j - p^{c,*}_j) \Big((p^M_j - p^{M,*}_j) -K_j(u_j - u^*_j) \\
 - k_{f,j}(\omega_j - \omega^*_j)
+\hspace{-3mm}\sum_{i: (i,j) \in \tilde{E}}\hspace{-3mm} \alpha_{ij} [ (p^c_i-p^{c,*}_i) - (p^c_j - p^{c,*}_j)] \Big).
\end{multline}

Additionally, define $V_P(\eta) = \sum_{(i,j) \in E} B_{ij} \int_{\eta^*_{ij}}^{\eta_{ij}} ( \sin \theta - {\color{black}\sin \eta^*_{ij}} ) \, d\theta$. Using~\eqref{sys1a} and~\eqref{sys1d}, the time-derivative is given~by
\begin{align}
\dot{V}_P &= \sum_{(i,j) \in E} B_{ij} (\sin \eta_{ij} - \sin \eta^*_{ij}) (\omega_i - \omega_j) \nonumber \\
&= \sum_{(i,j) \in E} (p_{ij} - p^*_{ij}) (\omega_i - \omega_j). \label{VPdiff}
\end{align}

Furthermore, from Design condition \ref{assum_gains},
it follows that there exist gains $k_{f,j},k_{c,j},k_{d,j}$ and a positive definite matrix $P_j = P^T_j$ such that \eqref{strictness_condition} holds. Then, \ank{let $V^M_j(x^M_j) = \frac{1}{2}(x^M_j - x^{M,*}_j)^TP_j(x^M_j - x^{M,*}_j)$ and note that it is positive definite. Following \eqref{dynsys}, the time derivative of $V^M_j$ is given by

\begin{multline}\label{VM_derivative}
\dot{V}^M_j = \frac{1}{2}(x^M_j - x^{M,*}_j)^T(P_jA_j + A_j^TP_j) (x^M_j - x^{M,*}_j) \\
+ (x^M_j - x^{M,*}_j)^T B_j (u_j - u^*_j),
\end{multline}
where $u_j$ and $u^*_j$ are given by \eqref{input_description} and \eqref{eqbr7} respectively.}

Based on the above, we \ak{consider the Lyapunov candidate
\begin{equation}\label{Lyapunov_function}
V(\eta, \omega^G, x^M, p^c) = V_F + V_P
 + \sum_{j \in G} V^M_j + V_C.
\end{equation}
Using \eqref{VFdiff}  - \eqref{VM_derivative}, the time derivative of V is given by
\begin{align}\label{dot_V}
\dot{V} =& \sum_{j \in G} \bigl[(\omega_j - \omega^*_j) (p^M_j - p^{M,*}_j)
 + \dot{V}^M_j  \nonumber\\
&\hspace{-7.5mm} + (p^c_j - p^{c,*}_j) \Big((p^M_j - p^{M,*}_j) -K_j(u_j - u^*_j)  - k_{f,j}(\omega_j - \omega^*_j)\Big)]
\nonumber \\
&\hspace{-7.5mm}
 - \sum_{j \in N} \Lambda_j(\omega_j - \omega^*_j)^2  - \sum_{(i,j) \in \tilde{E}} \alpha_{ij}((p^c_i - p^{c,*}_i) - (p^c_j - p^{c,*}_j))^2.
\end{align}
Using~\eqref{VM_derivative}, it therefore holds that
\begin{equation}\label{dot_V_2}
\dot{V} =  - \hspace{-2.5mm} \sum_{(i,j)\in \tilde{E}} \hspace{-2.0mm}  \alpha_{ij}((p^c_i - p^{c,*}_i) - (p^c_j - p^{c,*}_j))^2  + \sum_{j\in G}z_j^T Q_j z_j,
\end{equation}
where \ank{$z_j = \begin{bmatrix}
(p^c_j - p^{c,*}_j) & (x^M_j - x^{M,*}_j)^T & -(\omega_j -\omega^*_j)
\end{bmatrix}^T$} and $Q_j$ is given by

\begin{equation}\label{Suff_cond_proof}
Q_j = \begin{bmatrix}
-K_jk_{c,j} + D_j k_{c,j} & r^T_j \\
r_j & \begin{smallmatrix}
\frac{P_jA_j + A_j^TP_j}{2} & \frac{k_{d,j}P_jB_j - C_j^T}{2} \\
\frac{k_{d,j}B^T_jP_j - C_j}{2}  & -\Lambda_j - D_jk_{d,j}
\end{smallmatrix}
\end{bmatrix}.
\end{equation}
Note that $Q_j$ in \eqref{Suff_cond_proof} is identical to the matrix in \eqref{strictness_condition} when $\Lambda_j$ is replaced by some $\hat{\Lambda}_j < \Lambda_j$ and using Design condition \ref{assum_gains} on}
 \eqref{dot_V_2}, it therefore holds that
\begin{align}\label{vdotineq}
& \dot{V} \le - \hspace{-1.5mm} \sum_{j \in N} \bar{\Lambda}_j(\omega_j - \omega^*_j)^2  - \hspace{-2.5mm} \sum_{(i,j)\in \tilde{E}} \hspace{-2.0mm}  \alpha_{ij}((p^c_i - p^{c,*}_i) - (p^c_j - p^{c,*}_j))^2
 \nonumber \\
& \le 0
\end{align}
for $\bar{\Lambda}_j = \Lambda_j - \hat{\Lambda}_j > 0, j \in N$.
Clearly $V_F$ and has a strict global minimum at $\omega^{G,*}$ and  $V^M_j$ have strict global minima at $x^{M,*}_j$ for all $j \in G$ respectively.
 Moreover, $V_C$  has a strict global minimum at $p^{c,*}$.
 Furthermore, Assumption~\ref{assum_angles} guarantees the existence of some neighborhood of each $\eta^*_{ij}$ in which $V_P$ is increasing. Since the integrand is zero at the lower limit of the integration, $\eta^*_{ij}$, this immediately implies that $V_P$ has a strict local minimum at $\eta^*$. Thus, $V$ has a strict local minimum at the point $\ak{\Gamma^*} := (\eta^*, \omega^{G,*}, x^{M,*}, p^{c,*})$.
\kas{Therefore there exists a \icl{connected} set $\Xi := \{(\eta,  \omega^G, x^M, p^c) \colon V \le \epsilon\}$ containing $\Gamma^*$ where, for sufficiently small $\epsilon > 0$, $V$ is a nonincreasing function of all the system states, as follows from ~\eqref{vdotineq}, and has a strict local minimum at $\ak{\Gamma^*}$.  Therefore, $\Xi$ contains $\Gamma^*$ and is compact and positively invariant for~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description}.}

Lasalle's Invariance Principle can now be applied with the function $V$ on the compact positively invariant set~$\Xi$. This guarantees that all solutions of~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description} with initial conditions $(\eta(0),  \omega^G(0), x^M(0), p^c(0)) \in \Xi$ converge to the largest invariant set within $\Xi \, \cap \, \{(\eta,  \omega^G, x^M, p^c) \colon \dot{V} = 0\}$. We now consider this invariant set. If $\dot{V} = 0$ holds within $\Xi$, then~\eqref{vdotineq} holds with equality, hence  we must have $\omega = \omega^*$ for all $j \in N$ and $(p^c - p^{c,*}) \in \ak{\Ima}(\vect{1}_{|G|})$.

 Then, summing \eqref{eqbr2}--\eqref{eqbr3} over all $j \in G$ and $j \in L$ respectively, it follows that $\sum_{j \in G} p^{M}_j =  \sum_{j \in N} (p^L_j + \Lambda_j \omega^*_j) = c_1$, where $c_1$ is constant.
Furthermore, summing \eqref{p_command} over all $j \in G$, it follows that $\sum_{j \in G} \dot{p}^c_j = \sum_{j \in G}  (p^M_j - K_j u_j - k_{f,j} \omega^*_j) = \sum_{j \in G}  (p^M_j - K_j k_{c,j}p^c_j + (K_jk_{d,j} - k_{f,j}) \omega^*_j) = c_2 -\sum_{j \in G}(K_j k_{c,j}p^c_j)$, for some constant $c_2$.
Letting $\hat{p}^c = \sum_{j \in G} p^c_j$, from $(p^c - p^{c,*}) \in \ak{\Ima}(\vect{1}_{|G|})$ it holds that $\dot{\hat{p}}^c = c_2 - c_3 \hat{p}^c$, where $c_3 = \sum_{j \in G}(K_j k_{c,j})$. \ak{Therefore, it is trivial to show that $\hat{p}^c$ converges to some constant value, and since $(p^c - p^{c,*}) \in \Ima(\vect{1}_{|G|})$, it follows that also $p^c$ converges to some constant value $\bar{p}^c$.
The convergence of $(\omega, p^c)$ to $(\omega^*, \bar{p}^{c})$ results to $u_j = \bar{u}_j$ and by the conditions imposed on \eqref{dynsys} to $x^{M}_j = \bar{x}^{M}_j$ for all $j \in G$, where $(\bar{u}_j, \bar{x}^M_j)$ are constants.
Finally,  the constancy of $(\omega, p^M)$ guarantees from \eqref{eqbr2}--\eqref{eqbr4} that $\eta$ and $p$ are also constant.  Furthermore, by summing \eqref{eqbr2}--\eqref{eqbr3} and using \eqref{eqbr7} and the synchronisation of $p^c_j$ variables, it follows that $\bar{p}^c$ is unique and therefore equal to $p^{c,*}$, which also implies that $(\bar{u}_j, \bar{x}^M_j)$ are equal to $(u^*_j, x^{M*}_j)$.}

Therefore, we conclude by Lasalle's Invariance Principle that all solutions of~\eqref{sys1}, \eqref{dynsys}, \eqref{p_command}, \eqref{input_description} with initial conditions $(\eta(0),  \omega^G(0), x^M(0),  p^c (0)) \in \Xi$ converge to the set of equilibrium points defined in Definition~\ref{eqbrdef}, characterized by Lemma \ref{eqbr_lemma}. Finally, choosing for $S$ any open neighborhood of $\ak{\Gamma^*}$ within $\Xi$ completes the convergence proof.
  Finally,  as follows from \ak{Proposition} \ref{optthm}, when \eqref{opt_cond} holds, then the described equilibria are solutions to \eqref{Problem_To_Min}.
\hfill $\blacksquare$

\ak{Proof of Lemma \ref{suff_cond_lemma}:
The proof follows by analytically evaluating the eigenvalues of the matrix in \eqref{strictness_condition} with the matrices $A_j, B_j, C_j$ and $D_j$ following from \eqref{second_order} and\footnote{Note that these choices for $P_j$ and $k_{f,j}$ where numerically seen to minimize the required amount of frequency damping $\Lambda_j$ such that Design condition \ref{assum_gains} was satisfied when generation dynamics described by \eqref{second_order} were considered. Also, note that Design condition \ref{assum_gains} is also feasible for other choices of $P_j$ and $k_{f,j}$.} selecting $P_j = \frac{1}{K_j k_{c,j}} \begin{bmatrix}
\tau_{a,j} & 0 \\
0  & \tau_{p,j}
\end{bmatrix}$ and  $k_{f,j} = \frac{K_j}{2}(k_{c,j} + k_{d,j})$. In particular three of the eigenvalues of the resulting matrix are \kas{non-positive} and the fourth is negative for all positive values of $\tau_{a,j}, \tau_{p,j}$ if $\Lambda_j > \frac{K_j}{3k_{c,j}}(k_{c,j}^2 - k_{c,j}k_{d,j} + k_{d,j}^2))$ holds.
\hfill $\blacksquare$
}

\balance

\bibliographystyle{IEEEtran}
\bibliography{andreas_bib}

\end{document}

config: arxiv


TITLE: Can two-point sets be Borel?
QUESTION [18 upvotes]: Recall that a two-point set is a subset of the plane which meets every line in exactly two points. Such a set was first constructed by Mazurkiewicz in 1914.
I wonder if the following question of Erdos is solved:

Question. Is there a two-point set which is Borel?

See On sets which meet each line in exactly two points, where Mauldin says that he “believes” he first heard of the problem from Erdos, who in turn
said that it had been around since he (Erdos) was a “baby.”
Remarks. A two-point set can not be $F_\sigma$ (by a result of Larman). Also it is known that if a two-point set is analytic, then it is Borel.

REPLY [15 votes]: A two-point set cannot be $F_\sigma$, as Mohammad mentions in his question. Also,

A two-point set cannot contain a dense $G_\delta$ subset of an arc.

This was proved by Gareth Davies in his thesis (Oxford, 2011), but I do not think he ever published this result. To my knowledge, no better results are known in the Borel-sets-shouldn't-work direction.
In the other direction, the best known results belong to Arnie Miller. He showed that

If $V=L$ then there is a co-analytic two-point set.

See this paper of his from 1989 where he writes about this question and other such things, like nicely-definable Hamel bases and MAD families. Also see this related paper of Zoltan Vidnyansky's, which appeared more recently and extended/streamlined some of Miller's work.
In an unpublished manuscript that he made available here on his website, Arnie Miller also showed that

It is consistent for the Axiom of Choice to fail badly, but still to have two-point sets (e.g., two-point sets can exist even when there is no well-ordering of the reals).

Ben Chad did quite a bit of work about 5-10 years ago in trying to eliminate the Axiom of Choice from constructions of two-point sets as much as possible. Some of this stuff made it into this paper, joint with Robin Knight and Rolf Suabedissen.


TITLE: distribution associated with a discontinuous function
QUESTION [1 upvotes]: Let $f\colon\mathbb{R}\to\mathbb{R}$ be such that, for every $n\in\mathbb{Z}$, $f$ is differentiable on $\left(n,n+1\right)$ and $n$ is a discontinuity of first kind of $f$. We define
$$T_f(\phi)=\int_{\mathbb{R}}f(x)\phi(x)dx,\quad\text{where }\phi\text{ is a test function}.$$
I understand that, when $f$ is continuously differentiable, $T_f$ is well defined and the derivative of $T_f$ is just $T_{f'}$. However, these two do not seem clear in this case.
My attempt: Since $n$ and $n+1$ are discontinuities of first kind of $f$, we can define a new function, say $g_{n}$ such that $g_{n}$ is continuous on $[n,n+1]$ and equal to $f$ a.e. Thus $f$ is (Lebesgue) integrable on $[n,n+1]$. We then conclude that $f$ is locally integrable on $\mathbb{R}$ and therefore $T_{f}$ is well defined.
My question: 1) Is the above reasoning (to show that $T_f$ is well defined) correct? And:
2) What would the derivative of $T_f$ look like in this case? I think it is not $T_{f'}$ as $f'$ is not defined at $n$.
Any help/hint is highly appreciated.

REPLY [1 votes]: If $f$ is absolutely continuous on each segment $[n, n + 1]$ when defined to take the values $f(n + 0)$ and $f((n + 1) - 0)$ at the endpoints, integration by parts shows that the distributional derivative is
$$f'(x) + \sum_n (f(n + 0) - f(n - 0)) \delta(x - n),$$
where $f'(x)$ is the ordinary derivative. If the conditions for integration by parts do not hold, the distributional derivative may not be representable in this form. Also, if the test functions do not have finite support, it is necessary to add a condition on how fast $f$ can grow.


TITLE: What would be good example of antisymmetry for this relation?
QUESTION [2 upvotes]: We have set
$X=\{a,b,c,d\}$
and relation
$R=\{(a,a),(a,b),(a,c),(a,d),(b,b),(c,c),(d,b),(d,c),(d,d)\}$.
It is obvious that it is not symmetric and I suppose that it is antisymmetric but I can't come up with some good example to show it.
I know the rule of antisymmetry $xRy \wedge yRx \rightarrow x = y$, but I am not sure how to apply it on this relation.

REPLY [0 votes]: You need to exhaustively check through most of $R$. The following algorithm will allow you to determine if it is antisymmetric or not.

isAntisymmetric $\leftarrow$ True
for each $(x,y) \in R$:
$\quad$ if $x<y$:
$\quad\quad$ if $(y,x) \in R$: //use that $R$ (in this case) is lexographically ordered to quickly determine this
$\quad\quad\quad$ isAntisymmetric $\leftarrow$ False
return isAntisymmetric


TITLE: Make $kt^2+(3k+1)t+4k+1$ constant?
QUESTION [3 upvotes]: Find $k$ such that $kt^2+(3k+1)t+4k+1=0$ is an identity (i.e. true for all $t$).
E.g. $k=t+1$ doesn't work since you end up with a third degree polynomial in $t$ which determines $t$, making $t$ constant, so it isn't true for all $t$. Similarly $k=f(t)$ where $f$ is a polynomial doesn't work for the same reasons.
So we probably need $k=f(t,u)$ where $u$ is a variable.
Is it possible? Thanks.
EDIT: Is it possible to do this while making $k\in \mathbb{Z}[x]$? (i.e. a polynomial with integer coefficients)

REPLY [1 votes]: Hints: First, solve for $k$ to get
$$
k=-(t+1)\frac{1}{t^3+3t+4}\cdot
$$
To write the above as polynomial in $t$ (with infinite order), find $a$ and $b$ such that $t^3+3t+4=(t-a)(t-b)$ (in this case $a$ and $b$ are complex). Then
$$
\frac{1}{t^2+3t+4}=\frac{1}{(t-a)(t-b)}=\frac{1}{b-a}\left(\frac{a^{-1}}{1-t a^{-1}}-\frac{b^{-1}}{1-tb^{-1}}\right)
$$
from which the $n$-th order term of $(t^2+3t+4)^{-1}$ is
$$
c_n\equiv\frac{1}{b-a}(a^{-(n+1)}-b^{-(n+1)})
$$
which is a real (actually rational) number for all $n\geq -1$.  Then for $n\geq 0$, the $n$-th order term for $k$ is
$$
d_n\equiv-(c_n+c_{n-1}),\quad n\geq 0.
$$
Of course, $d_n$'s obtained this way agree with the Taylor expansion on $k(t)$ around $t=0$.