text:
\subsection{Dynamical Model}
Consider a linear model of CPS given by:
\begin{equation}\label{equ:system_model}
\begin{aligned}
    \mathbf{x}_{i+1} &= A\mathbf{x}_i \\
    \mathbf{y}_i &= C\mathbf{x}_i + \mathbf{e}_i,
    \end{aligned}
\end{equation}
where, $\mathbf{x}_i \in \mathbb{R}^n, \mathbf{y}_i \in \mathbb{R}^m$, with $m >n$, denote state vector and measurement vector at time $i$ respectively, $\mathbf{e}_i \in \mathscr{S}_k^m$ denotes the sparse attack vector, $A,C$ are system dynamic parameters. A control input may be included in the model above. However, since the control input is generally irrelevant to state estimation problems, we suppress in the model considered here. The following assumption is made regarding the CPS model above:
\begin{assumption}
The pair $(A,C)$ is full observable.
\end{assumption}
By iterating the system model \eqref{equ:system_model} $T$ time steps backwards, the $T$ horizon observation model is given by
\begin{equation}\label{equ:measurement_model}
    \mathbf{y}_T = H\mathbf{x}_{i-T+1}+\mathbf{e}_T,
\end{equation}
where $\mathbf{y}_T=[\mathbf{y}_i^{\top} \hspace{2mm} \mathbf{y}_{i-1}^{\top} \cdots \mathbf{y}_{i-T+1}^{\top}]^{\top} \in \mathbb{R}^{Tm}$ is a sequence of observation in the moving window $[i-T+1 \hspace{2mm} i]$, $\mathbf{x}_{i-T+1}\in \mathbb{R}^n$ is the state vector at time $i-T+1$, $\mathbf{e}_T = [\mathbf{e}_i^{\top} \hspace{2mm} \mathbf{e}_{i-1}^{\top} \cdots \mathbf{e}_{i-T+1}^{\top}]^{\top}$ is the sequence of attack vectors in the same moving window with $\mathbf{e}_i \in \mathscr{S}_{k}^m, \forall i \in [i-T+1 \hspace{2mm} i]$ and
$$ H = \begin{bmatrix} C A^{T-1} \\ \vdots \\ C A \\ C \end{bmatrix} =U \begin{bmatrix} \Sigma_1 \\0\end{bmatrix}V^{\top}
$$
where, $U = [U_1 \hspace{2mm} U_2]$, with $U_1 \in {\mathbb R}^{Tm \times n}, U_2 \in {\mathbb R}^{Tm \times Tm-n}$, is a matrix of left singular vectors of $H$, $V \in \mathbb{R}^{n \times n}$ is the corresponding matrix of right singular vectors, and $\Sigma_1 \in \mathbb{R}^{n \times n}$ is a diagonal matrix of the singular values of $H$, which are non-zero since $H$ is full rank.

\subsection{Threat Model}
We begin by defining a state decoder operation and a corresponding residual-based attack detector. Then we give a class of FDIA with guaranteed success against the defined decoder-detector pair. The following assumption is made concerning the attack vector:
\begin{assumption}
The attacker has knowledge of the system dynamics in \eqref{equ:system_model}.
\end{assumption}
\begin{definition}[Decoder]
Given a sequence of observation $\mathbf{y}_T \in \mathbb{R}^{Tm}$, a decoder for the measurement model in \eqref{equ:measurement_model} is a mapping of the form $\mathscr{D}:\mathbb{R}^{Tm} \mapsto \mathbb{R}^n$ given by
\begin{equation}\label{equ:decoder}
    \hat{\mathbf{x}} = \mathscr{D}(\mathbf{y}_T) =V\Sigma_1^{-1}\argmin_{\mathbf{z}}\|\mathbf{y}_T-U_1\mathbf{z}\|_{1}
\end{equation}
\end{definition}

\begin{definition}[Detector]
Given a positive threshold parameter $\epsilon$, a detector for the decoder $\mathscr{D}$ and the measurement model in \eqref{equ:measurement_model} is a binary classifier of the form
\begin{equation}\label{equ:detector}
    \mathscr{D}_{\epsilon}(\mathbf{y}_T) = \left\{
     \begin{array}{lr}
     1 &\quad\text{if}\hspace{0.2cm} \|\mathbf{y}_T - H\mathscr{D}(\mathbf{y}_T)\|_{1} > \epsilon \\
     0 &\quad\text{otherwise}
     \end{array}
     \right.
\end{equation}
\end{definition}

\begin{remark}
The positive class for $\mathscr{D}_{\epsilon}(\cdot)$ is regarded as the set of attacked (unsafe) measurements, while the negative class is regarded as the set of safe measurements. However, the negative class still contains unsafe measurements if the detector is compromised.
\end{remark}

\begin{definition}[Successful FDIA \cite{mo2010false}]\label{Def:Successful_FDIA}
Consider the CPS in \eqref{equ:system_model} and the corresponding measurement model \eqref{equ:measurement_model}, the attack sequence $\mathbf{e}_T\in\mathscr{S}^{Tm}_{km}$ is said to be $(\epsilon, \alpha)$-successful against the decoder-detector pair $\{\mathscr{D}, \mathscr{D}_{\epsilon}\}$ if
\begin{equation}\label{equ:def_success_FDIA}
     \|\mathbf{x}^{\star}-\mathscr{D}(\mathbf{y}_T)\|_2 \geq \alpha, \hspace{2mm} \mathscr{D}_{\epsilon}(\mathbf{y}_T)=0,
\end{equation}
where $\mathbf{y}_T = \mathbf{y}^{\star}_T+\mathbf{e}_T$ with $\mathbf{y}^*_T\in\mathbb{R}^{N}$ the true measurement vector, and $\mathbf{x}^{\star}$ is the true state vector.
\end{definition}
The following theorem gives a mechanism for constructing such $(\epsilon,\alpha)$-successful FDIA if the attack support is pre-determined.

\begin{theorem}\label{Thm:FDIA}
Given the support sequence $\mathscr{T} = \{\mathscr{T}_i \hspace{2mm} \mathscr{T}_{i-1} \cdots \mathscr{T}_{i-T+1}\}$ with $|\mathscr{T}_i|\leq k$. Let $\mathbf{z}_e$ be an optimal solution of the optimization program
\begin{equation}\label{equ:FDIA_optimization}
\begin{aligned}
    \textsf{Maxmize}&: \hspace{0.2cm} \|U_1 \mathbf{z}\|_{2}, \\
    \textsf{Subject to}&: \hspace{0.2cm} \| U_{1,\mathscr{T}^c}\mathbf{z}\|_{2} \leq \frac{\epsilon}{\sqrt{Tm-|\mathscr{T}|}}.
\end{aligned}
\end{equation}
If
$
\|U_{1,\mathscr{T}^c}\|_2 < \frac{1}{2\sqrt{Tm-|\mathscr{T}|}},
$
  then the FDIA
\begin{equation}\label{equ:FDIA}
    \mathbf{e}_{\mathscr{T}} = (U_{1,\mathscr{T}}) \mathbf{z}_e, \hspace{2mm} \mathbf{e}_{\mathscr{T}^c} = \mathbf{0}
\end{equation}
is $(\epsilon, \alpha)$-successful against the decoder-detector pair $\{\mathscr{D},\mathscr{D}_\epsilon\}$ for all
$$ \alpha \leq \frac{\epsilon}{2\sqrt{Tm}\overline{\sigma}}\left(\frac{1}{\overline{\sigma}_{\mathscr{T}^c}\sqrt{Tm-|\mathscr{T}|}} - 2\right),
$$
where $\overline{\sigma}_{\mathscr{T}^c}$ is the largest singular value of $U_{1,\mathscr{T}^c}$ and $\overline{\sigma}$ is the biggest nonzero singular value of $H$.
\end{theorem}
\begin{proof}
Since
$$\begin{aligned}
    \left\{\mathbf{z}\,\middle\vert\ \|\mathbf{z}\|_{2}\leq \frac{\epsilon}{\overline{\sigma}_{\mathscr{T}^c}\sqrt{Tm-|\mathscr{T}|}} \triangleq a\right\} \subset \\
\left\{\mathbf{z} \,\middle\vert\ \|U_{1,\mathscr{T}^c}\mathbf{z}\|_{2} \leq \frac{\epsilon}{\sqrt{Tm-|\mathscr{T}|}}\right\},
\end{aligned}
$$
then
$$\begin{aligned}\|U_1 \mathbf{z}_e\|_{2} \geq \max_{\|\mathbf{z}\|_{2} \leq a} \| \mathbf{z}\|_{2} \geq \frac{\epsilon}{\overline{\sigma}_{\mathscr{T}^c}\sqrt{Tm-|\mathscr{T}|}}.
\end{aligned}$$
Also
$$
    \mathbf{y}_T=\mathbf{y}^{\star}_T+P\begin{bmatrix}U_{1,\mathscr{T}} \\ 0\end{bmatrix}\mathbf{z}_e,
$$
with an appropriate permutation matrix $P$ satisfying $U_1 = P\begin{bmatrix}U_{1,\mathscr{T}} \\ U_{1,\mathscr{T}^c}\end{bmatrix}$. Hence, by decoder in \eqref{equ:decoder},
\begin{equation}\label{equ:solution_decoder}
    \begin{aligned}
    \hat{\mathbf{x}} &= V \Sigma_1^{-1} \argmin_{\mathbf{z}}\left\|U_1(\mathbf{z}^{\star}-\mathbf{z}) + P\begin{bmatrix}U_{1,\mathscr{T}}\\ 0\end{bmatrix}\mathbf{z}_e\right\|_{1} \\
&= V\Sigma_1^{-1}(\mathbf{z}^{\star}-\mathbf{z}_e^{\perp}),
\end{aligned}
\end{equation}
where $\mathbf{z}^{\star}=\Sigma_1 V^{\top}\mathbf{x}^{\star}$ corresponds to the true state and $\mathbf{z}_e^{\perp}$ is the projection given by
$$\mathbf{z}_e^{\perp} =\argmin_{\mathbf{z}}\left\|U_1\mathbf{z} + P\begin{bmatrix}U_{1,\mathscr{T}}\\ 0\end{bmatrix}\mathbf{z}_e\right\|_{1}.
$$
Then\footnote{Let $f(\mathbf{z}) = \left\|U_1\mathbf{z} - P\begin{bmatrix}0\\ U_{1,\mathscr{T}^c}\end{bmatrix}\mathbf{z}_e\right\|_{1}$ which is a convex function.The unique minimizer $\mathbf{z}_e^{\perp}$ satisfies $f(\mathbf{z}_e^{\perp}) \le f(\mathbf{z}_e^{\perp}-\mathbf{z}_e)$.},
$$\begin{aligned}\left\|U_1\mathbf{z}_e^{\perp} + P\begin{bmatrix}U_{1,\mathscr{T}}\\ 0\end{bmatrix}\mathbf{z}_e\right\|_{1} &\leq \left\|U_1\mathbf{z}_e^{\perp} - P\begin{bmatrix}0\\ U_{1,\mathscr{T}^c}\end{bmatrix}\mathbf{z}_e\right\|_1  \\
 \left\|U_{1,\mathscr{T}} \mathbf{z}_e\right\|_1-\|U_1\mathbf{z}_e^{\perp}\|_1  &\leq \|U_1\mathbf{z}_e^{\perp}\|_1+\left\|U_{1,\mathscr{T}^c} \mathbf{z}_e\right\|_1   \\
 2\|U_1\mathbf{z}_e^{\perp}\|_2 &\geq \left\|U_{1,\mathscr{T}} \mathbf{z}_e\right\|_1-\left\|U_{1,\mathscr{T}^c} \mathbf{z}_e\right\|_1 \\
\end{aligned}$$
Thus,
$$\begin{aligned}\|U_1\mathbf{z}_e^{\perp}\|_2 &\geq \frac{1}{2}\left(\left\|U_{1} \mathbf{z}_e\right\|_2-2\sqrt{Tm-|\mathscr{T}|}\left\|U_{1,\mathscr{T}^c} \mathbf{z}_e\right\|_2\right) \\
&\geq \frac{1}{2}\left(\frac{\epsilon}{\overline{\sigma}_{\mathscr{T}^c}\sqrt{Tm-|\mathscr{T}|}}-2\epsilon\right)
\end{aligned}$$
Next, according to the solution in \eqref{equ:solution_decoder}, it follows
$$\begin{aligned} \mathbf{x}^{\star} - \hat{\mathbf{x}} &= V \Sigma_1^{-1}\mathbf{z}_e^{\perp}\\
H(\mathbf{x}^{\star} - \hat{\mathbf{x}}) &= U_1 \mathbf{z}_e^{\perp}\\
\overline{\sigma} \|\mathbf{x}^{\star} - \hat{\mathbf{x}}\|_2 &\geq \frac{1}{\sqrt{Tm}}\|U_1 \mathbf{z}_e^{\perp}\|_1\\
\|\mathbf{x}^{\star} - \hat{\mathbf{x}}\|_2 &\geq \frac{1}{2\overline{\sigma}\sqrt{Tm}}\left(\frac{\epsilon}{\overline{\sigma}_{\mathscr{T}^c}\sqrt{Tm-|\mathscr{T}|}}-2\epsilon\right)
\end{aligned}$$
Since $
\|U_{1,\mathscr{T}^c}\|_2 < \frac{1}{2\sqrt{Tm-|\mathscr{T}|}},
$ the lower bound in the above inequality is positive.
Moreover,
$$\begin{aligned}
\|\mathbf{y}_T &-H\hat{\mathbf{x}}\|_{1} =\left\|U_1\mathbf{z}^{\star}+P\begin{bmatrix}U_{1,\mathscr{T}}\\ 0\end{bmatrix}\mathbf{z}_e-H\hat{\mathbf{x}}\right\|_{1} \\
&=\left\|U_1\mathbf{z}^{\star}+P\begin{bmatrix}U_{1,\mathscr{T}}\\ 0\end{bmatrix}\mathbf{z}_e-U_1(\mathbf{z}^{\star}-\mathbf{z}_e^{\perp})\right\|_{1}\\
& = \left\|U_1\mathbf{z}_e^{\perp}+P\begin{bmatrix}U_{1,\mathscr{T}}\\ 0\end{bmatrix}\mathbf{z}_e\right\|_{1} \\
&\leq\sqrt{Tm-|\mathscr{T}|}\left\|U_{1,\mathscr{T}^c} \mathbf{z}_e\right\|_2 \leq \epsilon
\end{aligned}$$
\end{proof}

\subsection{Prior Model}
The \emph{prior} information considered in this paper is an uncertain estimate of the support of the attack vector. There are many Machine Learning algorithms for estimating the location of attacks in a CPS \cite{ozay2015machine,mao2017least}. We will refer to such algorithm as \emph{localization algorithm} (or \emph{localization oracle}), and the resulting support estimate as \emph{support prior}.

Let $\mathscr{T}\subseteq\left\{1,2,\hdots, Tm\right\}$ be the actual support of the attacked nodes with the vector $\mathbf{q}\in\{0,\hspace{2mm}1\}^{Tm}$ the corresponding indicator
\begin{equation} \label{equ:support_indicator}
     \mathbf{q}_{i}=\left\{
     \begin{array}{lr}
     0 &\quad\text{if}\hspace{0.2cm} i \in \mathscr{T} \\
     1 &\quad\text{otherwise}
     \end{array}
     \right.
 \end{equation}

 Let $\hat{\mathscr{T}}$, be an estimate of $\mathscr{T}$, with the corresponding indicator $\hat{\mathbf{q}}\in\{0,1\}^{Tm}$ defined similarly to \eqref{equ:support_indicator}. Then an uncertainty model is defined as
  \begin{equation} \label{equ:uncertainty_model}
     \mathbf{q}_{i}=\epsilon_{i}\hat{\mathbf{q}}_{i}+(1-\epsilon_{i})(1-\hat{\mathbf{q}}_{i})
 \end{equation}
 where $\epsilon_i \sim \mathscr{B}(1,\mathbf{p}_i)$, with known $\mathbf{p}_i \in (0,1]$.

\begin{definition}[Positive Prediction Value, Precision, $\textsf{PPV}$]\label{Def:PPV_def}
Given the indicator vector estimate $\hat{\mathbf{q}}\in\{0,1\}^{Tm}$ of the unknown attack support indicator $\mathbf{q}\in\{0,1\}^{Tm}$, $\textsf{PPV}$ is the proportion of $\mathbf{q}$ that is correctly identified in $\hat{\mathbf{q}}$. It is given by
\begin{equation}\label{equ:PPV_def}
\textsf{PPV} =  \frac{\|\mathbf{q} \circ \hat{\mathbf{q}} \|_{0}}{\|\hat{\mathbf{q}}\|_{0}}
\end{equation}
\end{definition}

\subsection{Prior Pruning}
Based on the knowledge of uncertainty of \emph{prior}, a \emph{pruned support prior} $\hat{\mathscr{T}}_{\eta}^c$ is generated in two steps; first, an estimated safe nodes support vector is returned by the below \emph{offline} optimization program, given a reliability level $\eta \in (0,1)$,
\begin{equation}
    \begin{aligned}
    \textsf{Maximize} &\hspace{2mm} |\hat{\mathscr{I}}| \\
    \textsf{Subject to}& \hspace{2mm} \prod\limits_{i\in \hat{\mathscr{I}}}\mathbf{p}_i \geq \eta\\
    &\hspace{2mm}\mathscr{I}\in\left\{1,2,\hdots,Tm\right\}.
    \end{aligned}
\end{equation}

Then a \emph{pruned support prior} is obtained \emph{online} through a robust extraction:
\begin{align}\label{equ:pruning_alg}
    \hat{\mathscr{T}}_\eta^c = \mathscr{I}\cap\hat{\mathscr{T}}^c.
\end{align}

\begin{lemma}\label{Lem:PPV_1}
Given a \emph{support prior} $\hat{\mathscr{T}}$ generated by any underlying localization algorithm with associated uncertainty model in \eqref{equ:uncertainty_model}, if a \emph{pruned support prior} $\hat{\mathscr{T}}_{\eta}^c$ is generated by \eqref{equ:pruning_alg}, then it satisfies
\begin{equation}
    \textsf{Pr}\{\textsf{PPV}_{\eta} = 1\} \geq \eta,
\end{equation}
where $\textsf{PPV}_{\eta}$ is the precision of $\hat{\mathscr{T}}_{\eta}^c$ defined by \eqref{equ:PPV_def}.
\end{lemma}
\begin{proof}
From \eqref{equ:uncertainty_model}, it follows
$\mathbf{q}_i \hat{\mathbf{q}}_i = \epsilon_i \hat{\mathbf{q}}_i,
$
which implies
$$
\textsf{PPV} =  \frac{\|\mathbf{q} \circ \hat{\mathbf{q}} \|_{0}}{\|\hat{\mathbf{q}}\|_{0}} = \frac{\sum_{i=1}^N \mathbf{q}_i \hat{\mathbf{q}}_i}{\sum_{i=1}^N \hat{\mathbf{q}}_i}
=\frac{1}{|\hat{\mathscr{T}}^c|}\sum_{i=1}^N \epsilon_i \hat{\mathbf{q}}_i = \frac{1}{|\hat{\mathscr{T}}^c|}\sum_{i \in \hat{\mathscr{T}}^c}\epsilon_i
$$
Similarly,
$$ \textsf{PPV}_{\eta} = \frac{1}{|\hat{\mathscr{T}}^c_{\eta}|}\sum_{i \in \hat{\mathscr{T}}^c_{\eta}}\epsilon_i.
$$
Then, since $\hat{\mathscr{T}}_\eta^c\subseteq\mathscr{I}$ and $0<\mathbf{p}_i\le1$,
$$\textsf{Pr}\{\textsf{PPV}_{\eta} = 1\} = \textsf{Pr}\left\{\sum_{i \in \hat{\mathscr{T}}^c_{\eta}}\epsilon_i =|\hat{\mathscr{T}}^c_{\eta}|\right\} = \prod_{i\in \hat{\mathscr{T}}_{\eta}^c}\mathbf{p}_i \geq \eta
$$
\end{proof}

\subsection{Weighted $\ell_1$ Observer}
Consider the CPS in \eqref{equ:system_model}. Given a \emph{pruned support prior} $\hat{\mathscr{T}}^{c}_{\eta}$ in \eqref{equ:pruning_alg}, the weighted $\ell_1$ observer is given by the following moving horizon estimation problem:
\begin{equation}\label{equ:MHE_weighted_l1_minimization}
    \begin{aligned}
    \textsf{Minimize} &\hspace{2mm} \sum_{p=i-T+1}^i \|\mathbf{y}_p-C\mathbf{z}_p\|_{1,w}\\
    \textsf{Subject to}& \hspace{2mm} \mathbf{z}_{p+1} -A\mathbf{z}_p =0, \\
    &\hspace{10mm} p = i-T+1, \cdots,i-1 \\
    &\hspace{2mm} w = \left\{\begin{array}{ll}
1, & j \in \hat{\mathscr{T}}^{c}_{\eta} \\
\omega, & j \in \hat{\mathscr{T}}_{\eta}
\end{array}\right.
    \end{aligned}
\end{equation}
where, $0\leq \omega \leq 1$, $\|\mathbf{z}\|_{1, \mathrm{w}} \triangleq \sum_j \mathrm{w}_j |\mathbf{z}_j|$ is the weighted $1$-norm.

\begin{remark}
Consider the corresponding measurement model \eqref{equ:measurement_model}, then an equivalent form of \eqref{equ:MHE_weighted_l1_minimization} is given by
 \begin{equation}\label{equ:weighted_l1_minimization}
\begin{aligned}
\Minimize_{\mathbf{z}\in \mathbb{R}^n}\hspace{1mm} \|\mathbf{y}_T - H \mathbf{z}\|_{1, \mathrm{w}}, \hspace{1mm} \text{with} \hspace{1mm}w = \left\{\begin{array}{ll}
1, & j \in \hat{\mathscr{T}}^{c}_{\eta} \\
\omega, & j \in \hat{\mathscr{T}}_{\eta}
\end{array}\right.
\end{aligned}
\end{equation}
\end{remark}

\begin{theorem}\label{Thm:recovery_Withpruning}
Consider the CPS measurement model in \eqref{equ:measurement_model}, such that $U_2^{\top}$ satisfies the RIP condition
\begin{align}\label{equ: RIP}
   (1-\delta_{km})\|\mathbf{h}\|_{2}^2 \leq \|U_{2}^{\top}\mathbf{h}\|_{2}^2 \leq (1+\delta_{km})\|\mathbf{h}\|_{2}^2
 \end{align}
for all sparse vector $\mathbf{h} \in \mathscr{S}^{Tm}_{km}$.
Given a \emph{support prior} $\hat{\mathscr{T}}^c$ generated by an underlying localization algorithm with associated uncertainty model in \eqref{equ:uncertainty_model}.
Let $\hat{\mathscr{T}}_{\eta}^c$ be the \emph{pruned support prior} given by \eqref{equ:pruning_alg} with a parameter $\eta \in (0,1)$. If there exists an integer $a\geq \max\left\{\rho-1,1\right\}$, where $\rho = \frac{|\hat{\mathscr{T}}_{\eta}|}{km}$, such that
\begin{equation}
    \delta_{akm}+C\delta_{(a+1)km} \leq C-1
\end{equation}
then, with a probability of at least $\eta$, the estimation error due to the resilient observer in \eqref{equ:weighted_l1_minimization} can be upper bounded as
\begin{equation} \label{equ:error_bound3}
    \left\|\hat{\mathbf{x}}-\mathbf{x}^*\right\|_{2} \leq
   \frac{C_{1}}{\underline{\sigma}\sqrt{km}}\left(\omega\sigma_{km}(\mathbf{e})+(1-\omega)\|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}^c} \|_{1}\right),
\end{equation}
where $\underline{\sigma}$ is the smallest singular value of $H$, and
$$\begin{aligned} C_1 &= \frac{2a^{-\frac{1}{2}}(\sqrt{1-\delta_{(a+1)km}}+\sqrt{1+\delta_{a km}})}{\sqrt{1-\delta_{(a+1)km}}-\frac{1}{C}\sqrt{1+\delta_{akm}}}, \\  C &= \frac{a}{(\omega+(1-\omega)\sqrt{\rho -1})^2}.
\end{aligned}$$
\end{theorem}

\begin{proof}
To avoid repetition, we state upfront that all claims made in this rest of proof holds with a probability of at least $\eta$. Based on \emph{Lemma \ref{Lem:PPV_1}}, we have $\hat{\mathscr{T}}^{c}_{\eta} \subset \mathscr{T}^c$, $\mathscr{T}\subset \hat{\mathscr{T}}_{\eta}$ under this probability.

Let $\hat{\mathbf{x}}$ be a solution to \eqref{equ:weighted_l1_minimization} with $H\hat{\mathbf{x}} = H\mathbf{x}^*-\mathbf{h}$, then the corresponding minimizer $\hat{\mathbf{e}} = \mathbf{e}+\mathbf{h}$, thus $\|\hat{\mathbf{e}}\|_{1,\omega} \leq \|\mathbf{e}\|_{1,\omega}$. Following the definition of weight $w$ in \eqref{equ:weighted_l1_minimization},
$$\omega \|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}}+\mathbf{h}_{\hat{\mathscr{T}}_{\eta}} \|_{1} + \|\mathbf{e}_{\hat{\mathscr{T}}^{c}_{\eta}}+\mathbf{h}_{\hat{\mathscr{T}}^{c}_{\eta}}\|_{1} \leq \omega \|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}}\|_{1}
+ \|\mathbf{e}_{\hat{\mathscr{T}}^c_{\eta}}\|_{1}
$$
Since $1$-norm is decomposable for disjoint sets, such for $\mathscr{T}$ and $\mathscr{T}^c$, and $\hat{\mathscr{T}}^{c}_{\eta} \cap\mathscr{T} =\emptyset $, $\hat{\mathscr{T}}^{c}_{\eta} \cap \mathscr{T}^c =\hat{\mathscr{T}}^{c}_{\eta}$, $\mathscr{T}\cap \hat{\mathscr{T}}_{\eta}=\mathscr{T}$, it follows
$$\begin{aligned} \omega \|\mathbf{e}_{\mathscr{T}}+\mathbf{h}_{\mathscr{T}}\|_{1} +\omega\|\mathbf{e}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}+\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{1} + \|\mathbf{e}_{\hat{\mathscr{T}}^{c}_{\eta}}+\mathbf{h}_{\hat{\mathscr{T}}^{c}_{\eta}}\|_{1} \\
\leq \omega\|\mathbf{e}_{\mathscr{T}}\|_{1} +\omega \|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}\cap \mathscr{T}^c}\|_{1} +\|\mathbf{e}_{\hat{\mathscr{T}}^c_{\eta}}\|_{1}
\end{aligned}$$
By triangle inequality, it follows
$$ \omega\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{1}+\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta}^c}\|_{1} \leq \omega \|\mathbf{h}_{\mathscr{T}}\|_{1} +2\left(\|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}^c}\|_{1}+\omega\|\mathbf{e}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{1}\right)
$$
Adding and subtracting $\omega\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta}^c \cap \mathscr{T}^c}\|_{1} (= \omega\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta}^c}\|_{1})$ on LHS, $2\omega\|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}^c \cap \mathscr{T}^c}\|_{1} (= 2\omega\|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}^c}\|_{1})$ on RHS, it yields
$$\omega\|\mathbf{h}_{\mathscr{T}^c}\|_{1}+(1-\omega)\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta}^c}\|_{1} \leq  \omega\|\mathbf{h}_{\mathscr{T}}\|_{1} + \beta
$$
where, $\beta =2\left(\omega\|\mathbf{e}_{\mathscr{T}^c}\|_{1}+(1-\omega)\|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}^c}\|_{1}\right)$.
And since $\|\mathbf{h}_{\mathscr{T}^c}\|_{1} = \omega\|\mathbf{h}_{\mathscr{T}^c}\|_{1}+(1-\omega)\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta}^c}\|_{1}+ (1-\omega)\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{1}$, it yields
\begin{equation}\label{equ2}
\|\mathbf{h}_{\mathscr{T}^c}\|_{1} \leq \omega\|\mathbf{h}_{\mathscr{T}}\|_{1} + (1-\omega)\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{1}+ \beta
\end{equation}

Next, sort the coefficients of $\mathbf{h}_{\mathscr{T}^c}$ in descending order, and let $\mathscr{T}_j, j\in \{1,2,\cdots\}$ denote $j$th support in $\mathbf{h}_{\mathscr{T}^c}$ with size $akm \in \mathbb{Z}$, where $a>1$. Then $\|\mathbf{h}_{\mathscr{T}_{j}}\|_{2} \leq (akm)^{-1/2}\|\mathbf{h}_{\mathscr{T}_{j}}\|_{\infty}  \leq (akm)^{-1/2} \|\mathbf{h}_{\mathscr{T}_{j-1}}\|_{1}$, let $\mathscr{T}_0 = \mathscr{T} \cup \mathscr{T}_1$, we have
\begin{equation}\label{equ3}
\begin{aligned}
    \|\mathbf{h}_{\mathscr{T}_0^c}\|_{2} &\leq \sum_{j\geq 2}\|\mathbf{h}_{\mathscr{T}_j}\|_{2}\\ &\leq (akm)^{-1/2}\sum_{j\geq 1}\|\mathbf{h}_{\mathscr{T}_j}\|_{1} \\
    &\le  (akm)^{-1/2}\|\mathbf{h}_{\mathscr{T}^c}\|_{1}
    \end{aligned}
\end{equation}
Note that $|\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c|=(\rho-1)km $ since $\textsf{PPV}_{\eta}=1$, and $|\mathscr{T}| =km$, then $\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{1} \leq \sqrt{(\rho-1)km}\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{2}$ , $\|\mathbf{h}_{\mathscr{T}}\|_{1} \leq \sqrt{km}\|\mathbf{h}_{\mathscr{T}}\|_{2} \leq \sqrt{km}\|\mathbf{h}_{\mathscr{T}_0}\|_{2}$. Since $a \geq \rho-1$, we obtain $|\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c|=(\rho-1)km \leq akm =|\mathscr{T}_1|$,
then $\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{2} = \|\mathbf{h}_{\mathscr{T} \cup (\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c)}\|_{2} \leq \|\mathbf{h}_{\mathscr{T}_0}\|_{2}$.
Then combine with \eqref{equ2} and \eqref{equ3},
\begin{equation}\label{equ4}
    \|\mathbf{h}_{\mathscr{T}_0^c}\|_{2} \leq  \frac{\omega+(1-\omega)\sqrt{\rho-1}}{\sqrt{a}}\|\mathbf{h}_{\mathscr{T}_0}\|_{2} +\frac{\beta}{\sqrt{akm}}
\end{equation}

Since $ \|U_2^{\top}\mathbf{h}\|_{2}=\|(U_2^{\top}\mathbf{e}^{\star}-U_2^{\top}\mathbf{y}_T)-(U_2^{\top}\mathbf{e}-U_2^{\top}\mathbf{y}_T)\|_{2} =0$, it follows, based on triangle inequality and RIP condition,
$$\begin{aligned} \sqrt{1-\delta_{(a+1)km}}\|\mathbf{h}_{\mathscr{T}_0}\|_{2} \leq \|U_2^{\top}\mathbf{h}_{\mathscr{T}_0}\|_{2} \leq \|U_2^{\top}\mathbf{h}_{\mathscr{T}_0^c}\|_{2},
\end{aligned}$$
and
$$\begin{aligned} \|U_2^{\top}\mathbf{h}_{\mathscr{T}_0^c}\|_{2} \leq
\sum_{j\geq 2} \|U_2^{\top}\mathbf{h}_{\mathscr{T}_j}\|_{2}
\leq  \sqrt{1+\delta_{akm}}\sum_{j\geq 2} \|\mathbf{h}_{\mathscr{T}_j}\|_{2}  \\ \leq \frac{\sqrt{1+\delta_{akm}}}{\sqrt{akm}}\sum_{j\geq 1} \|\mathbf{h}_{\mathscr{T}_j}\|_{1} = \frac{\sqrt{1+\delta_{akm}}}{\sqrt{akm}} \|\mathbf{h}_{\mathscr{T}^c}\|_{1},
\end{aligned}$$
then combining with \eqref{equ2} yields
$$\begin{aligned} \sqrt{1-\delta_{(a+1)km}}\|\mathbf{h}_{\mathscr{T}_0}\|_{2} \leq  \omega\frac{\sqrt{1+\delta_{akm}}}{\sqrt{akm}} \|\mathbf{h}_{\mathscr{T}}\|_{1} + \\ (1-\omega)\frac{\sqrt{1+\delta_{akm}}}{\sqrt{akm}} \|\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{1} +\frac{\sqrt{1+\delta_{akm}}}{\sqrt{akm}}\beta.
\end{aligned}
$$
Notice, we have $\|\mathbf{h}_{\mathscr{T}}\|_{1}\leq \sqrt{km}\|\mathbf{h}_{\mathscr{T}_0}\|_{2} $ and $\|\mathbf{h}_{\hat{\mathscr{T}}_{\eta} \cap \mathscr{T}^c}\|_{1} \leq \sqrt{(\rho-1)km}\|\mathbf{h}_{\mathscr{T}_0}\|_{2}$, then manipulating the above inequality yields
$$ \|\mathbf{h}_{\mathscr{T}_0}\|_{2} \leq \frac{\frac{\sqrt{1+\delta_{akm}}}{\sqrt{akm}} \beta}{\sqrt{1-\delta_{(a+1)km}}-\frac{\omega+(1-\omega)\sqrt{\rho-1}}{\sqrt{a}}\sqrt{1+\delta_{akm}}}
$$

Since $\|\mathbf{h}\|_{2}\leq \|\mathbf{h}_{\mathscr{T}_0}\|_{2} +\|\mathbf{h}_{\mathscr{T}_0^c}\|_{2} $, combining the above inequality with \eqref{equ4} yields
$$ \|\mathbf{h}\|_{2}\leq \frac{2\frac{\sqrt{1-\delta_{(a+1)km}}+\sqrt{1+\delta_{akm}}}{\sqrt{akm}}\left(\omega\|\mathbf{e}_{\mathscr{T}^c}\|_{1}+(1-\omega)\|\mathbf{e}_{\hat{\mathscr{T}}_{\eta}^c}\|_{1}\right)}{\sqrt{1-\delta_{(a+1)km}}-\frac{\omega+(1-\omega)\sqrt{\rho-1}}{\sqrt{a}}\sqrt{1+\delta_{akm}}}
$$
where $\|\mathbf{e}_{\mathscr{T}^c}\|_{1} = \sigma_{km}(\mathbf{e})$ is best $k$ sparse approximation error of $\mathbf{e}$ defined in \eqref{equ:best_k_sparse_estimation_error}, and the estimation error $\|\hat{\mathbf{x}}-\mathbf{x}^*\|_{2} \leq \underline{\sigma}^{-1}\|\mathbf{h}\|_{2}$. Moreover, in order to ensure that the denominator is positive, the following condition is imposed:
$$ \delta_{akm} + C\delta_{(a+1)km} < C-1,
$$
where $C = \frac{a}{(\omega+(1-\omega)\sqrt{\rho-1})^2}$.
\end{proof}

\subsection{Numerical Simulation}
In this simulation, we compare the resiliency of three estimation scheme: $\ell_1$ observer without prior, weighted $\ell_1$ observer with prior generated an underlying attack localization algorithm with a true rate of $T_r=0.6$, and weighted $\ell_1$ observer with the pruned prior.
The system dimension is set as $m=20,n=10$, and then a full observable system is generated with random pair $(A,C)$ of independent Gaussian entries \cite{candes2005decoding}. For different attack percentage $P_A$, the FDIA is designed by \eqref{equ:FDIA} on random support $\mathscr{T}$.
By defining "success" as that the estimation error is less than $0.1\%$ of the real state, the success percentage is calculated from $1000$ trials. In Figure~\ref{Sucess_Attack_Tr60_m20_n10}, a performance comparison is presented for varying attack percentages.
As proved in literature \cite{fawzi2014secure}, \emph{$\ell_1$ observer without prior} cannot work when attack percentage is larger than $1/2$. The prior information can improve the resiliency of the estimator, but the improvement is limited because the precision of the prior information is not enough. By including pruning method, the resiliency is improved a lot.

\begin{figure}[h]
\begin{center}
\includegraphics[scale = 0.35]{Sucess_Attack_Tr60_m20_n10}
\end{center}
\caption{A comparison of estimation performance under different attack percentage $P_A \in [0,1]$ between \emph{$\ell_1$ Observer without Prior}, \emph{Weighted $\ell_1$ Observer with Prior}, and \emph{Weighted $\ell_1$ Observer with Prior Pruning}.}
\label{Sucess_Attack_Tr60_m20_n10}
\end{figure}

\subsection{Application Example}
In this subsection, an additional simulation is carried out for a more realistic IEEE Bus 14 system with a similar setup in \cite{anubi2020multi}. Here, $30\%$ of the measurement nodes are attacked. A comparison of the resulting phase angles estimation for an \emph{$\ell_1$ observer without prior} and the proposed \emph{weighted $\ell_1$ observer with prior pruning} is shown in Figure~\ref{angle_estimation}. Moreover, two comparison metrics for the estimation errors are given in Table~\ref{table_of_error_metric}. The first one is the \emph{root-mean-square} (RMS) value of error, the second one is the \emph{maximum absolute value} of the errors. Since the attack percentage is small, \emph{$\ell_1$ observer without prior} works well. However, there are still notable spikes that may induce closed loop instability. As shown in the Figure, the proposed \emph{weighted $\ell_1$ observer with prior pruning} completely eliminates the spikes!
\begin{table}[h]
\centering
 \caption{Error Metric Values\label{table_of_error_metric}}
 \begin{tabular}{||c | c | c | c | c | c | c||}
\hline
 &\multicolumn{3}{|c|}{RMS Metric}&\multicolumn{3}{|c||}{Max. Ans. Metric} \\
 \hline\hline
 & \textbf{LO} & \textbf{L1O} & \textbf{WL1P} & \textbf{LO} & \textbf{L1O} & \textbf{WL1P}\\
 \hline\hline
 $\delta_1$ & $2.5359$ & $0.0002$ & $0.00004$ & $5.7480$ & $0.0034$ & $0.0005$\\

$\delta_2$ & $2.3917$ & $0.0002$ & $0.0001$ & $5.7000$ & $0.0016$ & $0.0016$ \\

 $\delta_3$ & $2.6353$ & $0.0012$ & $0.0001$ & $7.5232$ & $0.0215$ & $0.0013$ \\

 $\delta_4$ & $2.3685$ & $0.0006$ & $0.0005$ & $5.7236$ & $0.0063$ & $0.0052$ \\

 $\delta_5$ & $2.6638$ & $0.0007$ & $0.0001$ & $7.8757$ & $0.0085$ & $0.0016$ \\
 \hline \hline
 \multicolumn{7}{||c||}{\textbf{LO}: Luenberger observer, \textbf{L1O}: $\ell_1$ observer without prior} \\
\multicolumn{7}{||c||}{\textbf{WL1P}: Weighted $\ell_1$ observer with prior pruning} \\
\hline
\end{tabular}
\end{table}

\begin{figure}[t]
\begin{center}
\includegraphics[scale = 0.3]{angle_estimation}
\end{center}
\caption{Angle estimation for IEEE-14 bus system under FDIA, "without Prior" means \emph{$\ell_1$ Observer without Prior}, "Prior with Pruning" means \emph{Weighted $\ell_1$ Observer with Prior Pruning}.}
\label{angle_estimation}
\end{figure}

\section*{NOTATION} \label{Sec: Notation}
The following notations and definitions are used throughout the whole paper: ${\mathbb R}, {\mathbb R}^n, {\mathbb R}^{n \times m}$ denote the space of real numbers, real vectors of length $n$ and real matrices of $n$ rows and $m$ columns respectively. ${\mathbb R}_{+}$ denotes positive real numbers. Normal-face lower-case letters $(e.g.\hspace{1mm} x \in {\mathbb R})$ are used to represent real scalars, bold-face lower-case letter $(e.g.\hspace{1mm} \mathbf{x} \in {\mathbb R}^n)$ represents vectors, while normal-face upper case $(e.g.\hspace{1mm} X \in {\mathbb R}^{n \times m})$ represents matrices. Let $\mathscr{T} \subseteq \{1,\dots, n\}$, then for a matrix $X \in {\mathbb R}^{m \times n}$, $X_{\mathscr{T}} \in {\mathbb R}^{|\mathscr{T}| \times n}$ is the sub-matrix obtained by extracting the rows of $X$ corresponding to the indices in $\mathscr{T}$. $\mathscr{T}^c$ denotes the complement of a set $\mathscr{T}$ and the universal set on which it is defined will be clear from the context. The support of a vector $\mathbf{x}$ is denoted by $\textsf{supp}(\mathbf{x}) \triangleq \{i|\mathbf{x}_{i} \neq 0\}$.  The set $\mathscr{S}_k^m \subset \mathbb{R}^n$ denotes the set of all vectors $\mathbf{v} \in \mathbb{R}^m$ such that $|\textsf{supp}(\mathbf{v})|\leq k$ (i.e the subset of $k$-sparse vectors). The \emph{best $k$-sparse approximation error} for $\mathbf{e} \in \mathscr{S}^m_{k}$ is given by
\begin{equation}  \label{equ:best_k_sparse_estimation_error}
    \sigma_k(\mathbf{e}) \triangleq \min_{\mathbf{z} \in \mathscr{S}^n_k} \|\mathbf{e}-\mathbf{z}\|_1 =\|\mathbf{e}_{\mathscr{T}^c}\|_1,
\end{equation}
where $\mathscr{T}$ is the support of $\mathbf{e}_i$ with first $k$ largest magnitude. The symbol $*$ denotes the convolution operator for vectors. The symbol $\circ$ denotes element-wise multiplication of two vectors and is defined as $\mathbf{z} = \mathbf{x} \circ \mathbf{y}$, where $\mathbf{z}_{i} = \mathbf{x}_{i} \cdot \mathbf{y}_{i}$. $\mathscr{B}(1,\mathbf{p}_i)$ denotes Bernoulli distributed variables with known independent probability $\mathbf{p}_i$.

\if false
and the probability mass function of sum of these variables is given in next Lemma.
\begin{lemma}(\cite{fernandez2010closed})  \label{Lem:PBD_CDF}
Given mutually independent Bernoulli random variables $\epsilon_i \sim \mathscr{B}(1,\mathbf{p}_i)$, $i=1,\cdots,N$, then
  \begin{equation}\label{equ:PMF}
      \textsf{Pr}\left\{\sum_{i=1}^{N} \epsilon_{i} = k-1 \right\}=\mathbf{r}(k), k=1,\cdots,N+1
  \end{equation}
  where,
  $$ \mathbf{r} = \beta \cdot \begin{bmatrix}-\mathbf{s}_1\\1 \end{bmatrix} * \begin{bmatrix}-\mathbf{s}_2\\1 \end{bmatrix} * \cdots * \begin{bmatrix}-\mathbf{s}_{N}\\1 \end{bmatrix}
  $$
  with
  $ \beta =\displaystyle\prod_{i=1}^{N} \mathbf{p}_{i}, \hspace{2mm} \mathbf{s}_i = -\frac{1-\mathbf{p}_{i}}{\mathbf{p}_{i}}
  $
\end{lemma}
\fi

\section{INTRODUCTION}\label{Sec: Introduction}
\input{1-introduction}

\section{MODEL DEVELOPMENT}\label{Sec: Model_Development}
\input{2-Model_Development}

\section{RESILIENT OBSERVER DESIGN} \label{Sec: Resilient_Observer_Design}
\input{3-Resilient_Observer_Design}

\section{NUMERICAL SIMULATION}\label{Sec: Simulation}
\input{4-Simulation}

\section{CONCLUSION}\label{Sec: Conclusion}
\input{5-Conclusion}

\medskip
\printbibliography

\end{document}

redpajama_set_name: RedPajamaArXiv


Get your perfect Call of Duty: WWII package.
The complete base game, including the campaign, multiplayer and Nazi Zombies game modes, The Resistance DLC and 1,000 COD WWII points.
The complete base game along with additional launch content and Season Pass, entitling players to additional post-launch content.


Prevail® Air Permeable Underpads for Low-Air-Loss Mattress Systems locks in the fluid, maximizing the air therapy benefit.
Prevail® Super Absorbent Breathable (Air Permeable) Underpads specifically designed for use with pressure relieving therapy beds.
Absorbent core locks in the fluid, maximizing the air therapy benefit.
Features an Integra Mat or bonded construction that reduces top sheet separation and clumping.
Prevail underpads are a disposable absorbent product designed to provide secondary leakage protection to protect beds, chairs and other surfaces.


At The Small Schools Multi Academy Trust Religious Education (RE) plays an important role in defining the Academies distinctive Christian Character. Religious education makes a unique contribution to the ethos of our schools and is reflected in our motto, "Believe and Achieve". As a Multi Academy Trust we recognise that spirituality lies at the heart of the curriculum. All members of the Multi Academy Trust experience Christianity through the life of the Trust, as well as taught through the curriculum. The skills, values and attitudes developed in RE are at the heart of a cohesive community where the children are encouraged to have a positive sense of identity and belonging. We aim to ensure that the RE curriculum is challenging, dynamic and relevant to pupils of all ages and encourage the children to ask questions to develop their critical thinking and understanding of the religions and beliefs which form our society, whilst promoting respect and diversity.
RE teaching is in line with the curriculum and the recommendations of the Statement of Entitlement for Church Schools, published by the National Society and adopted by the Lichfield Diocesan Board of Education (2012). The majority of our curriculum is based around developing the children's awareness of the traditions and beliefs of Christianity, particularly in relation to the Anglican Church. Christian morals and values are central to the ethos of our schools and we aim to give children a firm foundation of social values based on Christian principles. The encounter must be an open one which stems from and instils respect for different views and interpretations.
RE has a high profile within curriculum across the academy and learning activities provide fully for the needs of all learners. RE is taught discretely each week; our aim is to ensure that learners are inspired by the subject and develop a wide range of higher level skills such as enquiry, analysis, interpretation, evaluation and reflection to deepen their understanding of the impact of religion on the world.
Our curriculum is focussed around six broad themes; practices and ways of life, expressing meaning, diversity and belonging, responding to the natural world, values and commitments and beliefs and teaching. These broad themes allow children to explore multiple view points and religions and build comparison to identify similarities and differences, In addition they allow children to investigate their own beliefs and values and develop respect and sensitivity to those whose faiths and beliefs are different from their own.
A central element to our RE curriculum is ensuring children learn both about and from religion, through providing hands on experiences such as visiting places of worship, interacting with leaders of other faiths or handling religious artefacts we aim to ensure our RE is exciting and engaging for all children and creates a sense of awe and wonder. Children are assessed against criteria from both elements; learning from and learning about and assessment opportunities are catered to a range of learning styles to ensure progress for all children.
In addition, the subject contributes to other areas of education and human experience and plays an important part of the wider programme of spiritual, moral, social and cultural development.
context of a growing awareness and understanding of God.
Moral - helping each pupil develop their own informed values.
Social - helping pupils understand what shapes the values in our society.
Cultural – aiding pupils in exploring aspects of their own cultural heritage, and in developing positive attitudes towards diversity.
In our Worship experiences, we aim to promote the spiritual, moral, cultural and social development of the children. A whole school Act of Worship takes place daily, led by the Executive Head, Head of School or other members of staff. The Vicar leads the Act of Worship every fortnight and is regularly involved in the teaching of specific aspects of the Christian faith and special days, people and celebrations in the Church's year.
Special services to celebrate important times and festivals in the Christian Church's year take place at the local Church. Everyone in the school family, community and other friends and relatives are very welcome to join with us at these special times: Harvest, Christingle, Christmas, Mothering Sunday, Easter and also at the Leavers' Service at the end of the summer term. The children are actively involved in producing and leading the services showing their involvement with worship in the community.
Across our schools we aim to promote opportunities for personal reflection and spiritual development, this is encouraged across the curriculum specifically in our acts of worship and provides a safe environment for children to develop their own values and beliefs. Within the acts of worship, the children are actively involved both with discussions, stories, questions and prayers. The academy also encourages active engagement through developing the children as worship leaders allowing them to plan and deliver worship to the whole school.
At the heart of our multi academy trust are a set of core Christian Values; Love, Respect and Friendship, the development and understanding of these are rooted in all worships and curriculum opportunities. In addition our Worships are based around different Christian Values and are rotated half termly to allow for deeper understanding and application. These underpin our whole ethos, curriculum and school life and aim to ensure our pupils leave the academy as happy, secure, confident, caring and responsible citizens. The values we promote are important within the school community and are values which will be important throughout life. In all that we do we strive to share, encourage and demonstrate positive values with our children; in order that they understand them and strive to live by them.
Parents have the right to withdraw their children from religious education and collective worship. Children withdrawn will be supervised by a member of staff.


Fresh cow milk production reached 437 kl in 2014 in Mali, according to Faostat. This is 68 % more than in the previous year.
Historically, fresh cow milk production in Mali reached an all time high of 596 kl in 2010 and an all time low of 86.0 kl in 1961. When compared to Mali's main peers, fresh cow milk production in Algeria amounted to 2,167 kl, 132 kl in Burkina Faso, 126 kl in Mauritania and 487 kl in Niger in 2014.
Mali has been ranked 96th within the group of 172 countries we follow in terms of fresh cow milk production, 1 places behind the position seen 10 years ago.


text: Tarskavaig Tarscabhaig(Cod Bay) When crofting was introduced on Lord Macdonald's estates in 1811, 31 small crofts were created in Tarskavaig. They were small because the tenants were expected to supplement their income by fishing. A postal service from Broadford was introduced in 1867, but the first post office in the village was not opened till 1898. A Free Church chapel was built close to the shore, to the south of the township. The village of Tarskavaig has links to a well-known yacht-building firm in Sandbank Argyll, Alexander Robertson & Sons. Alexander Robertson, the company founder, was born in Inverkip in 1851, but his father, another Alexander, was born in Tarskavaig. Alexander senior was the son of the Gaelic bard, Donald Robertson, of No. 30 Tarskavaig. More information about Alexander Robertson & Sons (yachtbuilders) can be found on wikipedia. David Hutcheson, a descendant of the Robertsons of Sandbank, has written a 125 page history of Tarskavaig. It can be viewed and downloaded from The Islands Book Trust website https://islandsbooktrust.org/collections/digital-goods-vat-tax/products/tarskavaig-history

redpajama_set_name: RedPajamaCommonCrawl


"Home is where it all starts. It's where our children are raised and they learn how to go out into the world. It's home that gives people the grounding they need. Habitat creates this stability," says Linda Mealing, Habitat for Humanity of Huron Valley volunteer.
Habitat for Humanity of Huron Valley is committed to creating affordable homes for low-income families throughout Washtenaw County. With the attainment of a Habitat home, families find greater opportunities to improve their lives and children find more success in school and life.
Mealing, 68 of Ypsilanti, first heard about Habitat for Humanity's work when her two nieces received Habitat homes. She began volunteering in 2009, after a member of her church suggested it, knowing she wanted to do more in her retirement.
Mealing volunteers as a Family Support Partner, serving as a voice for the families on a one-on-one basis. Every family that is accepted into Habitat's program is assigned a partner who helps them stay on track, answers questions, see the bigger picture and makes sure they don't get overwhelmed. Mealing goes through the entire journey, start to finish, with her families.
Habitat for Humanity of Huron Valley's Executive Director, Sarah Stanton, believes Mealing's work is a key piece of the program.
Mealing also likes to get her hammer out and work on the homes as much as she can. She says she's learned a lot of useful skills and so do the families that work on their own homes.
In addition to providing new, affordable homes to families in need, Habitat for Humanity of Huron Valley also works to revitalize neighborhoods and partners with companies like DTE Energy to fulfill home appliance needs.
Mealing recalls Habitat restored homes in her neighborhood for the elderly who could no longer financially afford to keep their homes safe and up to date.
On a daily basis, Mealing gets to assist families in owning their first ever home. She gets to see the light in their eyes, their dreams coming true and their needs being met.


text:
\section*{Introduction}
Mechanical metamaterials can control softness via
a balance between the number of degrees of freedom of their components or nodes
and the number of constraints due to connections or links~\cite{Bertoldi2010,Sun2012,Kane2014,Chen2014,Lubensky2015, Meeussen2017, Rocks2017, Yan2017, Coulais2017}.
This balance, first noted by Maxwell~\cite{Maxwell1864} and later explored by Calladine~\cite{Calladine1978}, is termed \emph{isostaticity}.
In isostatic materials, softness can manifest itself via large-scale deformations, for example
as uniform Guest-Hutchinson modes~\cite{Guest2003, Kapko2009} or via periodic soft deformations, corresponding to so-called Weyl modes~\cite{Rocklin2016, Stenull2016}.
Uniform softness can be exploited to create extraordinary mechanical response~\cite{Florijn2014}, such as materials with a negative Poisson's ratio~\cite{Lakes1987,Bertoldi2010,Sun2012}.
Alternatively, localized softness has been programmed into isostatic materials in one and two dimensions
via a topological invariant called the polarization~\cite{Kane2014,Lubensky2015} that controls mechanical response and stress localization~\cite{Rocklin2016a} at an edge [including the edge of a disordered sample~\cite{Sussman2016}], an interface, or bound to a moving soliton~\cite{Chen2014}.
These mechanical~\cite{Nash2015, Susstrunk2015, Kariyado2015, Mousavi2015, Khanikaev2015,Yang2015,Huber2016} examples of topological metamaterials~\cite{Rechtsman2013, Lu2014, Wang2015} exhibit a general feature of topological matter~\cite{Hasan2010, Bernevig2013}: a correspondence between integer invariants in the bulk and response at a boundary.
Large-scale and localized deformations are deeply intertwined, as can be seen in demonstrations in which topological edge softness is created or destroyed by applying large uniform strains~\cite{Rocklin2017}.
A combination of topological polarization and localized defects can be used to program in softness or failure at a specified region in the material~\cite{Paulose2015,Paulose2015a}.

Although there are a number of examples of isostatic periodic structures in one, two, and three dimensions, the three-dimensional case is unique because all prior realizations of three-dimensional isostatic lattices include large-scale periodic deformations along continuous lines in momentum space~\cite{Paulose2015a, Stenull2016, Bilal2017}.
These Weyl lines define families of periodic soft modes in the material bulk
and contain a number of modes, which scales with the linear size of the structure.
As Ref.~\cite{Bilal2017} explores, Weyl lines can be useful to create a metamaterial surface with anisotropic elasticity, but in order to create a material whose top surface is much softer than the bottom, it proves necessary to collapse two Weyl lines on top of each other.
An alternative would be to find a metamaterial without Weyl lines.
However, these Weyl lines are generic and have a topological character which ensures that they cannot be annihilated locally---a single Weyl loop can only be destroyed by shrinking it to a point.
This presents a challenge in three-dimensional isostatic metamaterial design:
to achieve a gapped topological material analogous to those in two dimensions~\cite{Kane2014,Lubensky2015, Rocklin2017}, in which softness can be controlled and localized.

In this work, we design gapped topological materials by exploring the parameter space of the generalized stacked kagome lattice.
The mechanical response of this gapped structure is characterized by a nonzero topological polarization oriented along the $z$-axis.
This topological polarization $\bm P$ can be exploited
to localize soft modes in the material bulk by introducing topological defects within the lattice structure called \emph{dislocation loops}.
These dislocations are characterized by a topological invariant called the Burgers vector $\bm b$.
Along the dislocation, we show that the topological charge characterizing the softness or rigidity of the lattice (with unit cell volume $V_{cell}$) depends on the orientation $\bm{\hat\ell}$ of the dislocation line and is given by $\bm P\cdot(\bm b\times\bm{\hat\ell})/V_{cell}$ per unit length.

\begin{figure}[th!]
\centerline{\includegraphics{Fig1}}
\caption{\label{Fig1} Architecture of a gapped 3D stacked kagome lattice.
(A) We consider the doubled unit cell of a vertically stacked kagome lattice and rotate
the top constituent triangle around three perpendicular axes by angles $(\phi, \theta, \psi)$ [the bottom triangle is rotated by $(- \phi, - \theta, - \psi)$].
(B) For the choice $(\phi, \theta, \psi) = (\pi/3, 0, \pi/3)$ of these parameters, we obtain a unit cell that corresponds to a (gapped) 3D lattice with a topological polarization $\bm{P}$ and without Weyl lines.
(C) Architecture of the gapped metamaterial with one unit cell emphasized in the lower left and a 2D kagome sub-lattice surface highlighted in green on the top. }
\end{figure}

\section*{Gapping the stacked kagome lattice}
We examine the mechanics of metamaterial structures by using a lattice model for the displacements of nodes and strains of the links. Within this set-up, softness corresponds to displacements with small strains, and rigidity to strains that induce only small displacements.
We place a particle (or ball) at each vertex and connect the neighboring particles by linear springs. Such models capture the small-strain response of realistic structures that are either 3D-printed from soft polymers~\cite{Bilal2017, Paulose2015a} or assembled from construction sets~\cite{Chen2014} or laser-cut components~\cite{Paulose2015}.

The mechanics of this ball-spring model are captured via the linearized equation of motion
$\bm{\ddot X} = - D \bm{X}$, where $\bm{X}\equiv(\bm{x}_1,\dots,\bm{x}_N)$ is a $d \times N$-dimensional vector containing the displacements of all $N$ particles in $d$-dimensions relative to their equilibrium positions.
For a given lattice geometry, we calculate the dynamical matrix $D$, which relates the forces exerted by springs to displacements of the particles. For simplicity, we work in units in which all particle masses and spring constants are one.
In the linear regime, we can find the dynamical matrix by first relating the $N_B$-dimensional vector of spring extensions $\bm{S}\equiv(s_1,\dots,s_{N_B})$ to the displacements via $\bm{S}=R\bm{X}$.
The rigidity matrix $R$ contains $dN\times N_B$ entries determined by the equilibrium positions of the particles and the connectivity of the lattice.
In combination with Hooke's law, the matrix $R$ lets us calculate the dynamical matrix $D$ via the relation
$D = R^\mathsf{T}R $~\cite{Calladine1978,Sun2012}. Although the rigidity theory in terms of the matrix $R$ offers an equivalent description to the dynamical matrix, it contains additional physical insight.

The null spaces (kernels) of the rigidity matrix and its transpose give us information about lattice mechanics. From the relation $\bm{S}=R\bm{X}$, we note that $\ker R$ contains \emph{soft modes}, i.e., collective displacements that (to lowest order) do not stretch or compress any of the springs. On the other hand, the matrix $R^\mathsf{T}$ relates forces on particles $\bm{\ddot X}$ to spring strains via $\bm{\ddot X} = - R^\mathsf{T}\bm{S}$. From this relation, one notes that $\ker R^\mathsf{T}$ contains combinations of spring tensions that do not give rise to particle forces. These configurations are dubbed \emph{states of self-stress} because they define load-bearing states, which put the lattice under a static tension~\cite{Calladine1978}. Localized states of self-stress can concentrate uniformly-applied external loads, which can lead to selective buckling in pre-programmed regions~\cite{Paulose2015a}.

To calculate soft modes (states of self-stress), we need to numerically solve the configuration-dependent equation $R\bm{X} = \bm{0}$ ($R^\mathsf{T}\bm{S} = \bm{0}$).
However, a mathematical result called the rank-nullity theorem lets us compute the difference between the number of zero modes $N_0 \equiv \mathrm{null} R$ and the number of states of self-stress $N_B \equiv \mathrm{null} R^\mathsf{T}$~\cite{Calladine1978,Sun2012,Lubensky2015}. This difference, the rigidity charge $\nu \equiv N_0 - N_B$, is given by the dimensionality of the $R$ matrix: $\nu = dN - N_B$,
and can only change if balls or springs are either added or removed.
We focus on the special case of isostatic, or Maxwell, lattices defined by $\nu = 0$ when the system is considered under periodic boundary conditions. These lattices are marginally rigid and exhibit a symmetry between zero modes and states of self-stress.

To design metamaterials based on simple, repeated patterns, we focus on periodic structures. Periodicity allows us to explicitly calculate zero modes and states of self-stress in a large sample. We begin with a highly symmetric, ``undeformed'' lattice and explore its configuration space by changing the positions (but not the connectivity) of the nodes.
The specific geometry that we consider is illustrated in Fig.~1: this stacked kagome lattice has a coordination number $z\equiv 2N_B/N = 2d= 6$ and $\nu = 0$: the lattice is isostatic.
This lattice is based on the two-triangle unit cell shown in Fig.~\ref{Fig1}A, which corresponds to two unit cells of the kagome lattice stacked on top of each other. We deform the lattice via the orientations of the triangles, which are governed by the three angles $(\phi,\theta,\psi)$ of rotation around the $(x,y,z)$-axes, respectively.
Within this structure, we fix the connectivity and explore a range for each of the three angles from $-\pi/2$ to $\pi/2$. In Fig.~\ref{Fig1}B, we show the deformed unit cell for a particular structure in this region, which corresponds to the choice $(\phi,\theta,\psi) = (\pi/3,0,\pi/3)$. This unit cell builds the periodic geometry in Fig.~\ref{Fig1}C. We proceed to quantitatively demonstrate that this metamaterial is gapped.

We calculate the spectrum of normal modes by considering plane-wave solutions of the form $x_{\alpha,\bm{\ell}}=e^{i\bm \ell\cdot\bm{k}}x_{\alpha}$, where the index $\alpha$ refers to a particle within one unit cell and the lattice index $\bm\ell$ enumerates different unit cells within a lattice.
The three-component wavevector $\bm{k}$ is periodic in each component with $-\pi\leq k_x, k_y,k_z\leq\pi$. This collection of points forms the first Brillouin zone of the lattice.
The Bloch representation of the rigidity matrix in this plane-wave basis is $R_{\alpha\beta}(\bm{k})=R_{\alpha\beta}e^{i(\bm\ell_\alpha-\bm\ell_\beta)\cdot\bm{k}}$. A zero mode at wave-number $\bm{k}$ is a vector of unit-cell displacements $x_\alpha$ that solves $R_{\beta\alpha}(\bm{k})x_\alpha=0$. Thus, within this setting, zero modes correspond to the zeros of the complex function $\det R(\bm{k})$. By examining this function, we find the zero modes for different configurations of the lattice.

\begin{figure*}[th!]
\centerline{\includegraphics{Fig2}}
\caption{\label{Fig2} Phase diagram of stacked kagome lattices.
(A) The number of soft directions in $\phi$-$\psi$ parameter space for $\theta = 0$ (see Fig.~\ref{Fig1}A).
The region in dark blue corresponds to a gapped topological lattice, whereas the rest corresponds to lattices containing two (cyan), four (green) or six (orange) Weyl loops.
(B) Weyl loops along the line $(\theta,\psi)=(0,\pi/3)$, corresponding to magnified points in (A). As the gapped region is approached we observe Weyl loops shrinking towards the origin until they disappear.
(C) An example of the soft-mode structure in reciprocal space for a lattice with a single two-Weyl-loop structure (in blue), which corresponds to the cyan point $(\phi,\theta,\psi)=(5 \pi/12,0,\pi/4)$ in (A).
The integration of a phase around any contour (e.g., solid red circle) enclosing the Weyl loop
gives an integer winding number that we use to check the topological nature of the Weyl loop.
Note that this contour can be deformed into a pair of contours ($C_z$, dashed red lines) traversed in opposite directions: this illustrates how the winding number across the Brillouin zone changes as the Weyl loop is traversed.
Near the origin $(k_x,k_y,k_z) = (0,0,0)$, the loop manifests itself via soft directions (teal).
(D) The Brillouin zone of a lattice corresponding to the blue point $(\phi,\theta,\psi)=(\pi/3,0,\pi/3)$ in (A). Far away from the origin, the lattice is checked to be gapped. At the origin, the softest regions (blue region) align themselves with the soft directions of lattices with $(\theta,\psi)=(0,\pi/3)$ (red) or $(\phi,\theta)=(\pi/3,0)$ (green).}
\end{figure*}

\clearpage

We numerically evaluate zero modes for the stacked kagome lattice (see Supporting Information [SI] for details) and for most values of the angles $(\phi,\theta,\psi)$ find collections of zero modes along compact loops in $\bm{k}$-space [Fig.~\ref{Fig2}A]. These Weyl loops, shown in Fig.~\ref{Fig2}B--C,
are analogous to one-dimensional nodal lines~\cite{Volovik2007, Burkov2011, Lu2013} and experimentally observed zero-dimensional Weyl points~\cite{Xu2015,Lu2015} in three-dimensional electronic semimetals and photonic crystals.
In isostatic lattices, the number of Weyl loops is always even, because the materials' time-reversal symmetry maintains $\bm{k}\to-\bm{k}$ reflection symmetry in the Brillouin zone---each Weyl loop comes in a pair with its reflected partner~\cite{Stenull2016}.
Furthermore, these loops attach to the origin of the Brillouin zone along soft directions, i.e., the lines tangent to the Weyl loops at the origin. We count the number of loops by looking at soft directions in the neighborhood of the origin and, in Fig.~\ref{Fig2}A, plot this count in a slice of parameter space. In this phase diagram, we note regions with up to six different loops.
Strikingly, the middle of the diagram displays a region in which no Weyl loops exist [Fig.~\ref{Fig2}A,D].

We conclude that although Weyl loops are generic, the stacked kagome lattice also exhibits \emph{gapped} configurations which contain no Weyl loops.
In lattices with Weyl loops, the number of soft modes scales as the linear size of the system, whereas gapped lattices have only the three uniform Guest-Hutchinson modes~\cite{Guest2003, Kapko2009}.
We are not aware of any other realizations of a three-dimensional isostatic lattice which is gapped.
In the next section, we address the implications of the existence of a gap for the topological characterization of mechanical networks.

\section*{Topological rigidity in three dimensions}
For the gapped isostatic lattice, an integer topological invariant called the polarization can be computed from the bulk phonon spectrum~\cite{Kane2014}. This winding-number invariant is only well defined in the gapped case: in Weyl lattices, the closure of the gap prevents a consistent definition. That is, when Weyl loops are present, the polarization changes depending on the choice of contour in the Brillouin zone. The mechanical consequences of polarization are apparent by looking at the spectrum of the material's soft surface waves. The polarization controls which surfaces have more soft modes: this is the mechanical version of the bulk-boundary correspondence principle.

\begin{figure}[th!]
\centerline{\includegraphics{Fig3}}
\caption{\label{Fig3}  Surface modes for a gapped polarized lattice (C--D) compared to
a lattice with Weyl lines (A--B).
A surface mode in a plane perpendicular to the $z$-axis (see, e.g., Fig.~\ref{Fig1}C) decays exponentially into the material bulk with
amplitude $e^{- \kappa z/|\bm{a}_3|}$, where $\kappa$ is the inverse penetration depth.
Lattice with Weyl lines $(\phi,\theta,\psi)= (5 \pi/12, 0, \pi/4)$: (A) Inverse penetration depth for three soft surface modes over the entire surface Brillouin zone (sBZ).
For small wavevectors $(k_x, k_y)$, these modes correspond to the three acoustic phonons.
(B) Zoom on the middle mode in (A): the part of wavevector space for which
the mode lives on the bottom side is colored in red, whereas the part corresponding
to the top side is in blue. The projection of the Weyl line from the bulk
onto the sBZ corresponds to a line with infinite penetration depth ($\kappa = 0$) and separates the red and blue regions.
Gapped polarized lattice $(\phi,\theta,\psi)= (\pi/3, 0, \pi/3)$: (C) Inverse penetration depth for the same three modes as in (A), but for a different choice of parameters,
corresponding to the blue region in Fig.~\ref{Fig2}A.
(D) Zoom-in on the middle mode in (C), which [in contrast to the mode in (B)] is localized on the top of the lattice only:
it is blue throughout. Significantly, the bulk ($\kappa = 0$) modes that
correspond to the Weyl loop are absent for this lattice.}
\end{figure}

\clearpage

For the gapped polarized lattice, a polarization vector $\bm P$ can be computed via $\bm P=\sum_i m_i\bm a_i - \bm{d}_0$, where $\bm{a}_i$ are the three lattice vectors
and, in our case, the unit-cell dipole $\bm{d}_0 = -3 \bm a_3/2$. This dipole is computed via the expression $\bm{d}_0 \equiv - \sum_{b}\bm{r}_b$, where $\bm{r}_b$
are the positions of the bond centers relative to the center of mass of the unit cell~\cite{Kane2014}.
The three coefficients $m_i$ are winding numbers computed by integrating the phase change of the complex function $\det R(\bm k)$ across a straight-line contour crossing the Brillouin zone:
\begin{equation}
m_i \equiv - \frac{1}{2\pi i}\oint_{C_i}d k_i \frac{d}{d k_i}\ln\det R(\bm k), \label{eq:m}
\end{equation}
where $C_i$ is a closed contour in the $\bm{\hat k_i}$ direction~\cite{Kane2014}. This integral is well defined for contours along which $\det R(\bm k) \neq 0$. The integers $m_i$ depend smoothly on the choice of contour and are constant. In gapped lattices, the determinant is non-zero everywhere except for the origin $\bm k=0$. As a consequence, the winding numbers are independent of the chosen contour and the polarization is a topological invariant. On the other hand, Weyl loops partition the space of straight-line contours into lines going through the inside (outside) of the loop. Because $\det R(\bm k)=0$ along the loops, contours on either side of the loop can have different winding numbers $m_i$. Note that the combination of any two such contours taken in opposite directions can be smoothly deformed (without intersecting the Weyl line) into a small circle enclosing the Weyl line, as shown in Fig.~\ref{Fig2}C. The winding number $m_\mathrm{W}$ around this small circle is an invariant and equal to the difference between $m_i$ for $C_i$ on two sides of the loop. This topological protection guarantees that a single Weyl loop cannot be destroyed, which explains why Weyl loops only vanish by shrinking to the origin within the phase diagram in Fig.~\ref{Fig2}A. In summary, Weyl lines are protected by a winding number, and gapped lattices can have a well-defined topological polarization.

Bulk-boundary correspondence states that the topological invariants computed in the bulk can have significant effects on the mechanics of a sample with boundaries. We demonstrate this correspondence by computing the topological invariants and the spectrum of soft edge states in the stacked kagome lattice.
In the bulk of the Weyl lattice corresponding to $(\phi,\theta,\psi)=(5\pi/12,0,\pi/4)$ (c.f., Fig.~\ref{Fig2}A,C), the winding number around the loop is $m_\mathrm{W}=- 1$. The $m_3$ for contours inside and outside of the loop differ by one. To see the boundary counter-part of the correspondence, we compute the spectrum of soft edge states for this lattice with a stress-free surface parallel to the $xy$-plane (see Fig.~\ref{Fig1}C). In Fig.~\ref{Fig3}A--B,
we show the (signed) inverse penetration depth $\kappa$ for soft surface modes in the two-dimensional Brillouin zone for wavevector $(k_x,k_y)$. These plane-wave solutions of $\det R = 0$ correspond to soft modes of the form $e^{i\bm\ell\cdot(\bm{k}+i\bm{\kappa})}\bm{x}_\alpha$. The middle mode of Fig.~\ref{Fig3}A (see zoom in Fig.~\ref{Fig3}B) changes sides: for $\kappa > 0$ (red region), the soft mode is attached to the bottom surface and decays upwards into the bulk, whereas for $\kappa < 0$ (blue region), the mode is localized on the top surface.
For the line corresponding to the projection of the Weyl loop onto the two-dimensional Brillouin zone, the penetration depth is infinite and $\kappa = 0$, because the Weyl lines are soft modes in the bulk.
The difference in $m_3$ between inside and outside of the Weyl loops is  the bulk invariant $m_\mathrm{W} = - 1$, which corresponds to the difference between the number of soft surface modes across the projected Weyl loop. This connection follows from Cauchy's argument principle for Eq.~(\ref{eq:m}), which states that the third winding numbers $m_3$ count, up to a constant, the number of modes with zero energy ($\det R = 0$) at a boundary. These observations confirm bulk-boundary correspondence.

A similar correspondence exists in the gapped lattice. There, the invariant polarization pointing along the $z$-axis is given by $\bm P = \bm a_3/2$, which we computed for parameters $(\phi,\theta,\psi)=(\pi/3,0,\pi/3)$.
To understand the effect of $\bm P$ on the boundary, note in analogy with electromagnetism that the rigidity charge $\nu^S$ in region $S$ is related to the flux of polarization $\bm P$ through the region's boundary $\partial S$. In a nearly-uniform, gapped lattice~\cite{Kane2014}
\begin{equation}
\nu^S{=}\oint_{\partial S}\frac{dA}{V_{cell}}\bm{\hat n}\cdot \bm P,
\label{eq:nu}
\end{equation}
where $\bm{\hat n}$ is the boundary's inward normal and $V_{cell}=\det(\bm a_1,\bm a_2,\bm a_3)$ is the volume of a unit cell.
The difference between the rigidity charges at the top versus the bottom is determined by $\bm P$: for the stacked kagome, we expect one more band of soft modes along the top surface relative to the bottom.
The total number of soft surface modes is three as a result of three bonds being cut per unit cell~\cite{Kane2014}.
We plot the inverse penetration depth $\kappa$ for this extra band in Fig.~\ref{Fig3}C--D within the two-dimensional Brillouin zone.
Unlike the Weyl lattice, the whole middle band in Fig.~\ref{Fig3}C is blue: the top surface has more soft modes and is, therefore, softer.

Equation~(\ref{eq:nu}) shows that topological polarization acts analogously to an electric polarization. Inside a homogeneous polarized material, the charge is zero. However, if homogeneity is broken by the presence of boundary or \emph{defects}, charge can accumulate at these spots.

\begin{figure*}[h!]
\centerline{\includegraphics{Fig4}}
\caption{\label{Fig4} Dislocation loops can act as dipoles of topological charge.
(A)~Architectures of a dislocation loop in a periodic lattice:
in a screw dislocation, the burgers vector $\bm{b}$ is parallel to the dislocation line
whereas in an edge dislocation, the burgers vector and the dislocation line are perpendicular.
 A dislocation loop can combine edge and screw dislocations along its contour. In a gapped polarized lattice (with polarization $\bm{P}$), a dislocation line separates in space an edge dislocation segment that carries zero modes (ZM, blue) from an edge dislocation segment that carries states of self-stress (SSS, red), for example via screw dislocations (green) which carry no such charges. The net topological charge, defined as the difference between ZM and SSS, is zero when summed over a dislocation loop contour. Nevertheless, the dislocation loop carries a topological charge dipole, which is in this example is parallel to the Burgers vector $\bm{b}$.
(B)~Geometry of the dislocation loop. Each prism represents a unit cell with triangles oriented according to $(\phi,\theta,\psi)=(\pi/3,0,\pi/3)$ as shown in the inset. (C)~Numerical results for the four lowest-frequency phonons (excluding the trivial translations) in a large ($36\times36\times18$ unit cells) gapped lattice that has a dislocation loop: warmer color signifies a larger displacement within that unit cell. Plotted are only those unit cells that have displacements above a cutoff of 20\% of the maximum.
Note the localization of the softest modes to the near side of the loop, in accordance with part A and the polarization $\bm{P}$.
(Caption continued on next page.)
}
\end{figure*}
\clearpage
\addtocounter{figure}{-1}
\begin{figure*}[th!]
\caption{(continued)
(D)~The common logarithm of lowest mode frequencies $\log \omega^2$ for the $N$-th lowest mode, plotted versus $N$, comparing two large samples (same size as C).
The mode frequencies are significantly lower for the dislocated lattice (black), than in the non-dislocated case (gray).
Without a dislocation, the lowest eigenmodes are the extended plane-wave acoustic phonons, whereas
with the dislocation, these modes include both the acoustic phonons and the modes localized along the dislocation loop (in blue: $4$ lowest modes, whose eigenvectors are plotted in C).
}
\end{figure*}

\section*{Local rigidity and softness at dislocations}
For many applications such as cushioning~\cite{Bilal2017}, programmed assembly~\cite{Hawkes2010}, or controlled failure~\cite{Paulose2015a}, materials need to be designed with build-in softness or rigidity. The previous section illustrates how a topological invariant can be programmed into a material to design softness at boundaries. We now proceed to show how a flux of this topological polarization can be harnessed to create soft regions in the bulk of the material. Equation~(\ref{eq:nu}) suggests that a nonzero polarization flux can be achieved by considering an inhomogeneous material.
We choose dislocations to provide this inhomogeneity.

The natural defects in three-dimensional crystals are line dislocations: displacements of unit cells along straight lines which terminate at material boundaries. In many crystalline solids, such defects control mechanical deformations and plasticity. As these lines are topologically protected, the only way to confine a line dislocation to a finite region is to form a closed \emph{dislocation loop} (Fig.~\ref{Fig4}A--B).  The topological invariant characterizing a dislocation is the Burgers vector $\bm b$, which expresses the effect of a dislocation on the surrounding lattice. If a particle at point $\bm x$ is displaced to $ \bm x+\bm u(\bm x)$ according to displacement $\bm u(\bm x)$, then the Burgers vector is $\bm b=\oint_C d\bm u$, where $C$ is any closed contour surrounding the dislocation line.

Line dislocations come in two primary types called \emph{edge} and \emph{screw} dislocations, which are distinguished by the orientation of the Burgers vector $\bm b$ relative to the dislocation line direction $\hat{\bm\ell}$ (see Fig.~\ref{Fig4}A). For edge dislocations, these vectors are orthogonal: the displacement $\bm u$ pushes unit cells apart in order to insert a half-plane of unit cells that extends from the dislocation line in the direction $\hat{\bm\ell}\times\bm b$. In this way, edge dislocations are three-dimensional generalizations of two-dimensional point dislocations.
By contrast, for screw dislocations, the $\bm b$ and $\hat{\bm\ell}$ vectors are parallel: the displacement $\bm u$ pushes neighboring cells apart along the loop. In this way, screw dislocations give rise to an inherently three-dimensional spiral structure. Along a dislocation loop, the Burgers vector $\bm b$ is constant, but the line direction $\hat{\bm\ell}$ changes. As Fig.~\ref{Fig4}A shows, a dislocation loop can contain both edges and screws.

The interplay between line properties $(\bm{\hat{\ell}},\bm{b})$ and the topological polarization $\bm P$ of the lattice leads to a rigidity charge that can be expressed in a simple formula (derived in the SI) for the charge line density $\rho^L$ localized at the dislocation line:
\begin{equation}
\rho^L = \frac{1}{V_{cell}} \bm P\cdot(\bm b\times\bm{\hat\ell}),
\label{eq:rhol}
\end{equation}
where  $V_{cell}$ is the unit cell volume of the uniform lattice.
Equivalently, we could write $\rho^L = V_{cell}^{-1}\det (\bm P,\bm b,\bm{\hat\ell})$. Consequently, for any configuration in which these three vectors are not linearly independent, the charge density is zero. In particular, screw dislocations do not carry charge---their Burgers vector points along the line direction. For edge dislocations, the charge is negative (positive) if $\bm P$ points along (against) the half-plane direction $\hat{\bm\ell}\times\bm b$: intuitively, these rows carry extra polarization out or in.
Note that for any dislocation loop, the net charge $\oint_L \rho^L$ is zero.
However, dislocations loops do separate rigidity charges in space, resulting in a topological charge dipole, as shown in Fig.~\ref{Fig4}A.
This dipole is oriented within the $\bm{b}\,$-$\bm{P}$ plane and can be quantified by a dipole moment, which captures the amount of charge separated and the distance of separation.
The dipole's $\bm{b}$ ($\bm{P}$) component is the area of the dislocation loop after it is projected onto the plane formed by the two vectors $\bm{b}$ ($\bm{P}$) and $\bm{b}\times\bm{P}$ (see SI for derivation). For the loop composed of edge and screw dislocations
shown in Fig.~\ref{Fig4}A, all of the charges are localized along the edge dislocations and the dipole moment lies along $\bm b$.

A positive charge density corresponds to soft modes localized along part of the dislocation loop. We investigate this localized softness within a polarized lattice using the configuration shown schematically in Fig.~\ref{Fig4}A and plotted for a small sample in Fig.~\ref{Fig4}B. In Fig.~\ref{Fig4}C, we demonstrate that for the softest (i.e., lowest-frequency non-translational) modes of this lattice, the unit cells with the largest displacements are localized along the near side of the loop, in agreement with Eq.~(\ref{eq:rhol}).
This can be contrasted with the lowest-energy modes of a sample without a dislocation: in Fig.~\ref{Fig4}D, we show that the dislocated lattice has soft modes at lower frequencies.
Whereas the lowest modes in the non-dislocated sample are the largest-wavelength acoustic phonons that fit within the periodic box, the lowest modes of the dislocated lattice are a combination of these acoustic phonons
and the many localized modes such as the ones shown in Fig.~\ref{Fig4}C.

\section*{Conclusions and outlook}
We demonstrated that dislocation loops in polarized three-dimensional lattices
can localize soft modes at one part of the loop. Similarly, the opposite part of the loop localizes states of self-stress. Although the total topological charge along the loop is zero, the loop acts as a topological charge dipole. The loop can be used to localize topologically protected soft modes in one region of the material, inside the material bulk. Furthermore, because the lattice we have designed is gapped, these localized soft modes do not compete with the many extended, periodic soft modes characteristic of Weyl loops.

Let us conclude by discussing applications for both three-dimensional gapped lattices and for localized softness. Three-dimensional topological materials may be useful in the design of devices for cushioning, including safety devices such as helmets, in which a hard side necessary to resist stress exists in combination with a soft side necessary to cushion the body~\cite{Bilal2017}.
Localized softness can be used to pre-program large displacements and to isolate the rest of the material from strain. In active mechanical metamaterials~\cite{Khanikaev2015,Yang2015, Souslov2017}, for example in self-folding origami~\cite{Hawkes2010, Silverberg2014, Chen2016}, these mechanisms could be actuated using motors to strain the material in a well-controlled way~\cite{Paulose2015}. Furthermore, localized states of self-stress can be used to control material failure via either buckling or fracture, which then isolates the rest of the material from failure even under a large external load.

\section*{Acknowledgments}
We thank Paul Baireuther and Bryan G.~Chen for fruitful discussions and T.~C.~Lubensky for a critical reading of the manuscript.
We gratefully acknowledge funding from FOM, NWO, and Delta Institute for theoretical physics.

\clearpage

\makeatletter
\def\@seccntformat#1{%
  \expandafter\ifx\csname c@#1\endcsname\c@section\else
  \csname the#1\endcsname\quad
  \fi}
\makeatother

\setcounter{equation}{0}
\setcounter{figure}{0}

\section{Supporting Information for ``Localizing softness and stress along loops in three-dimensional topological metamaterials''}

\vspace{2em}

\renewcommand{\theequation}{S\arabic{equation}}
\renewcommand{\figurename}{{\bf Fig. }}
\renewcommand{\thefigure}{{\bf S\arabic{figure}}}

\section{Generalized stacked kagome lattice}

The generalized stacked kagome lattice is built by repeating a given unit cell along the following lattice vectors
\begin{equation}\label{eq:PrimVectors}
\begin{array}{ccc}
\bm a_1=a(1,0,0)
&\bm a_2=a(1/2,\sqrt3/2,0)
&\bm a_3=a(0,0,1)
\end{array},
\end{equation}
where $a$ defines the lattice spacing. In the undeformed kagome lattice the equilibrium positions of the particles in the unit cell is given by the triangular prism
\begin{equation}
\begin{array}{ccc}
\bm r(1)=-\tfrac{1}{4}\bm a_3
&\bm r(2)=\tfrac{1}{2}\bm a_1-\tfrac{1}{4}\bm a_3
&\bm r(3)=\tfrac{1}{2}\bm a_2-\tfrac{1}{4}\bm a_3
\vspace{4pt}\\
\bm r(4)=+\tfrac{1}{4}\bm a_3
&\bm r(5)=\tfrac{1}{2}\bm a_1+\tfrac{1}{4}\bm a_3
&\bm r(6)=\tfrac{1}{2}\bm a_2+\tfrac{1}{4}\bm a_3,
\end{array}
\end{equation}
where particles $s=1,2,3$ ($s=4,5,6$) constitute the bottom (top) triangle of the prism.

Within the unit cell we connect the particle pairs (2,3), (3,1), (1,2), (5,6), (6,4), (4,5), (1,4), (2,5), (3,6) by springs. Horizontally between unit cells we connect pairs (3,2), (1,3), (2,1), (6,5), (4,6), (5,4) according to the nearest neighbor principle. Finally we connect (4,1), (5,2) and (6,3) between the unit cells in the positive vertical direction.

We let a deformation act on the unit cell by rotating the triangles around three orthogonal axes attached to the triangles. The first two axes are coplanar to the triangle, with the first one bisecting one of its angles. The third axis is perpendicular to the triangle. The three axes meet in the center of the triangle. We then define our deformation by rotating the top triangle around the three axis by angles $\phi$, $\theta$, and $\psi$ respectively, while rotating the bottom triangle by $-\phi$, $-\theta$, and $-\psi$. Note that by taking our rotation axes relative to the triangles, our rotation group is commutative and we do not need to worry about the order of rotations. Explicitly, the equilibrium positions of the particles in the generalized stacked kagome can be given as
\begin{widetext}
\begin{align}
\bm r(1)&=\frac{a}{2\sqrt{3}}
\begin{pmatrix}
s_\psi c_\theta-\sqrt{3}c_\psi c_\phi-c_\psi s_\phi s_\theta
\\-c_\psi c_\theta-\sqrt{3} s_\psi c_\phi- s_\psi s_\phi s_\theta
\\- c_\phi s_\theta+ \sqrt{3}s_\phi+ \sqrt{3}
\end{pmatrix}+\frac{ar}{2}\begin{pmatrix}c_\alpha\\ s_\alpha\\0\end{pmatrix}\nonumber
\\\bm r(2)&=\frac{a}{2\sqrt{3}}
\begin{pmatrix}
s_\psi c_\theta+\sqrt{3}c_\psi c_\phi-c_\psi s_\phi s_\theta
\\-c_\psi c_\theta+ \sqrt{3}s_\psi c_\phi- s_\phi s_\phi s_\theta
\\- c_\phi s_\theta- \sqrt{3}s_\phi+ \sqrt{3}
\end{pmatrix}+\frac{ar}{2}\begin{pmatrix}c_\alpha\\ s_\alpha\\0\end{pmatrix}\nonumber
\\\bm r(3)&=\hspace{5pt}\frac{a}{\sqrt{3}}
\begin{pmatrix}
- s_\psi c_\theta+c_\psi s_\phi s_\theta
\\c_\psi c_\theta+ s_\psi s_\phi s_\theta
\\ c_\phi s_\theta + \sqrt{3}/2
\end{pmatrix}+\frac{ar}{2}\begin{pmatrix}c_\alpha\\ s_\alpha\\0\end{pmatrix}\nonumber
\\\bm r(4)&=\frac{a}{2\sqrt{3}}
\begin{pmatrix}
- s_\psi c_\theta-\sqrt{3}c_\psi c_\phi-c_\psi s_\phi s_\theta
\\-c_\psi c_\theta+\sqrt{3} s_\psi c_\phi+ s_\psi s_\phi s_\theta
\\ c_\phi s_\theta- \sqrt{3}s_\phi- \sqrt{3}
\end{pmatrix}-\frac{ar}{2}\begin{pmatrix}c_\alpha\\ s_\alpha\\0\end{pmatrix}\nonumber
\\\bm r(5)&=\frac{a}{2\sqrt{3}}
\begin{pmatrix}
- s_\psi c_\theta+\sqrt{3}c_\psi c_\phi-c_\psi s_\phi s_\theta
\\-c_\psi c_\theta- \sqrt{3}s_\psi c_\phi+ s_\phi s_\phi s_\theta
\\ c_\phi s_\theta+ \sqrt{3}s_\phi- \sqrt{3}
\end{pmatrix}-\frac{ar}{2}\begin{pmatrix}c_\alpha\\ s_\alpha\\0\end{pmatrix}\nonumber
\\\bm r(6)&=\hspace{5pt}\frac{a}{\sqrt{3}}
\begin{pmatrix}
 s_\psi c_\theta+c_\psi s_\phi s_\theta
\\c_\psi c_\theta- s_\psi s_\phi s_\theta
\\- c_\phi s_\theta - \sqrt{3}/2
\end{pmatrix}-\frac{ar}{2}\begin{pmatrix}c_\alpha\\ s_\alpha\\0\end{pmatrix},
\end{align}
\end{widetext}
where $s_\cdot=\sin(\cdot)$ and $c_\cdot=\cos(\cdot)$.

\clearpage
\begin{widetext}
Following the above order of bond enumeration, the compatibility matrix takes the form:
\begin{equation}
\label{eq:RigMat}
R^\dagger(\bm k)=\begin{pmatrix}
\begin{array}{c}
\bm{\hat b}^\mathsf{T}_1\\\bm{\hat b}^\mathsf{T}_2\\\bm{\hat b}^\mathsf{T}_3\\
\bm{\hat b}^\mathsf{T}_4\\\bm{\hat b}^\mathsf{T}_5\\\bm{\hat b}^\mathsf{T}_6\\
\bm{\hat b}^\mathsf{T}_7\\\bm{\hat b}^\mathsf{T}_8\\\bm{\hat b}^\mathsf{T}_9\\
\bm{\hat b}^\mathsf{T}_{10}\\\bm{\hat b}^\mathsf{T}_{11}\\\bm{\hat b}^\mathsf{T}_{12}\\
\bm{\hat b}^\mathsf{T}_{13}\\\bm{\hat b}^\mathsf{T}_{14}\\\bm{\hat b}^\mathsf{T}_{15}\\
\bm{\hat b}^\mathsf{T}_{16}\\\bm{\hat b}^\mathsf{T}_{17}\\\bm{\hat b}^\mathsf{T}_{18}
\end{array}&
\begin{array}{ccc|ccc}
0&-1&1&\phantom{\hat b^\mathsf{T}_1}\\
1&0&-1&\phantom{\hat b^\mathsf{T}_1}&O\\
-1&1&0&\phantom{\hat b^\mathsf{T}_1}\\\hline
\phantom{\hat b^\mathsf{T}_1}&&&0&-1&1\\
\phantom{\hat b^\mathsf{T}_1}&O&&1&0&-1\\
\phantom{\hat b^\mathsf{T}_1}&&&-1&1&0\\\hline
\phantom{\hat b^\mathsf{T}_1}&&&&&\\
\phantom{\hat b^\mathsf{T}_1}&-I&&&I&\\
\phantom{\hat b^\mathsf{T}_1}&&&&&\\\hline
0&e^{-i k_1}&-1&\phantom{\hat b^\mathsf{T}_1}\\
-1&0&e^{-i k_2}&\phantom{\hat b^\mathsf{T}_1}&O\\
e^{-i (k_1-k_2)}&-1&0&\phantom{\hat b^\mathsf{T}_1}\\\hline
\phantom{\hat b^\mathsf{T}_1}&&&0&e^{-i k_1}&-1\\
\phantom{\hat b^\mathsf{T}_1}&O&&-1&0&e^{-i k_2}\\
\phantom{\hat b^\mathsf{T}_1}&&&e^{-i(k_1-k_2)}&-1&0\\\hline
\phantom{\hat b^\mathsf{T}_1}&&&&&\\
\phantom{\hat b^\mathsf{T}_1}&e^{-i k_3}I&&&-I&\\
\phantom{\hat b^\mathsf{T}_1}&&&&&\\
\end{array}
\end{pmatrix}
\end{equation}
\end{widetext}
where the initial spring direction $\bm b_i$ is defined as the normalization of $\bm r(s')-\bm r(s)$, where $(s,s')$ is the $i$-th spring. In this notation, each number in the matrix grid represents a multiple of the triple on the far left in the same row. This turns the matrix into an $18\times18$ square matrix. As the initial bond directions depend smoothly on the initial particle positions, each element of $Q$ depends smoothly on the deformation parameters $\phi$, $\theta$, and $\psi$.

\section{Soft directions}

Weyl lines are formed by collections of points satisfying $\det R(\bm k)=0$. This equation is often quite hard to solve. However, with observational evidence suggesting that the Weyl loops generally connect to the origin, we can formally expand the determinant around $k\equiv||\bm k||=0$ as
\begin{equation}
\det R(\bm k)
=\sum_{\bm{n}\in \mathbb{Z}^3}c_{\bm n}e^{-i\bm n\cdot\bm k}
=i k^3\mathcal P_3(\hat{\bm k})+k^4\mathcal P_4(\hat{\bm k})+\mathcal O(k^5),
\end{equation}
for finitely many nonzero coefficients $c_{\bm n} \in \mathbb{R}$, where $\mathcal P_m$ is an $m$-th order homogeneous polynomial with real coefficients in three variables over the unit sphere $||\hat{\bm k}||=1$. The vanishing constant, linear, and quadratic polynomials reflect the triple multiplicity of the zero at $\bm k=0$.

A Weyl line can only attach to the origin at a direction $\hat{\bm k}$ where both $\mathcal P_3(\hat{\bm k})$ and $\mathcal P_4(\hat{\bm k})$ vanish.
Suppose, for instance, that $\mathcal P_4(\hat{\bm k}) \ne 0$ and note that the error term is dominated by $M k^5$ for all $k<\delta$, where $M$ and $\delta$ are fixed positive numbers.
As a consequence, the relation
\begin{equation}
|\det R(\bm k)|> k^4\left(|\mathcal P_4(\hat{\bm k})|-M k\right)
\end{equation}
holds for all $k<\delta$. In particular, this means $|\det R(\bm k)|>0$ for all $k<\min\left(\delta,\frac{1}{M}|\mathcal P_4(\hat{\bm k})|\right)$. We call these common zeroes of $\mathcal P_3(\hat{\bm k})$ and $\mathcal P_4(\hat{\bm k})$ over the unit sphere the soft directions of $\det R$.

As any Weyl loop connects to the origin along two directions, we find that a configuration with $2n$ soft directions can accommodate no more than $n$ pairs of Weyl loops.
In the main text we show the $\theta=0$ slice of soft-mode phase space, containing a region where the lattice does not exhibit any soft directions. This absence of soft directions is rare. Sampling of the $0<\phi,\theta,\psi<\pi/2$ phase space with step size $\pi/60$ we obtain the following data for the prevalence of different numbers of soft directions:
\begin{table}[ht!]\centering
\begin{tabular}{lr|r|r|r|r|r|r}
Number of soft directions:&0&2&4&6&8&10&12\\\hline
Percentage of phase space:\phantom{$\Big|$}&$<$0.1&20.0&53.5&17.3&6.9&2.2&$<$0.1
\end{tabular}
\end{table}

\section{Rigidity charge bound to line dislocations}

The dislocation acts on the polarization and the unit cell volume through the lattice vectors $\bm a_i$. Letting the displacement function $\bm u$ act on both the initial and terminal point of the vector $\bm a_i$, we obtain the lattice vector field in the dislocated lattice:
\begin{equation}
\tilde{\bm a}_i(\bm x)
=[\bm x+\bm a_i+\bm u(\bm x+\bm a_i)]-[\bm x+\bm u(\bm x)]
=\bm a_i+(\bm a_i\cdot\bm\nabla)\bm u(\bm x),
\end{equation}
where we assume that higher order differentials of $\bm u$ are negligible.
The polarization changes as
\begin{equation}
\tilde{\bm P}(\bm x)
=\sum_i m_i\tilde{\bm a}_i
=\bm P+(\bm P\cdot\bm\nabla)\bm u(\bm x).
\end{equation}
The unit cell volume changes as
\begin{equation}
\tilde{V}_{cell}(\bm x)
=\det( \tilde{\bm a}_1,\tilde{\bm a}_2,\tilde{\bm a}_3)
=V_{cell}(1+\bm\nabla\bm u(\bm x)).
\end{equation}
We start from Eq.~(1) in the main text. In a canonical basis with coordinates $(x^1,x^2,x^3)$, a surface element is given by $dx^\ell\wedge dx^k=\epsilon^{i\ell k}dS_i$, where $d\bm S$ is positively oriented. In this notation, the change in rigidity charge in the region $S$ from the $\nu^S=0$ periodic case is given by
\begin{align}
\tilde{\nu}^S-\nu^S
&=\int_{\partial S}(\tilde{V}_{cell}^{-1}\,\tilde{\bm P}-V_{cell}^{-1}\,\bm P)\cdot(-d\bm S)\nonumber\\&
=V_{cell}^{-1}\int_{\partial S} P_i\partial_ju_j-P_j\partial_ju_i\,dS_i\nonumber\\&
=V_{cell}^{-1}\int_{\partial S}P_i \epsilon_{ijk}\partial_\ell u_j\,\epsilon^{m\ell k}dS_m\nonumber\\&
=V_{cell}^{-1}\int_{\partial S}\big(\bm P\times d\bm u\big)_k\wedge dx^k. \label{eq:dislcharge}
\end{align}

Consider a straight dislocation line segment $L$ in direction $\hat{\bm\ell}$ and let $S$ be the cylinder with axis $L$ and radius $\mathcal R$. Without loss of generality, we choose canonical coordinates such that $L=\{t\hat{\bm\ell}: 0<t<T\}$ with $\hat{\bm\ell}$ pointing along the third basis vector. Furthermore, we take $T\hat{\bm\ell}$ to be an integer combination of lattice vectors $\bm a_i$. If the displacement function respects the periodicity of the lattice, the top ($t=T$) and bottom ($t=0$) surfaces of $S$ will cancel by orientation. It remains to calculate Eq.~(\ref{eq:dislcharge}) along the vertical boundary of $S$. First consider the $k=1$ term:
\begin{equation}
\int_{\partial S}\big(\bm P\times d\bm u\big)_1\wedge dx^1=\int_0^{2\pi}\left(\bm P\times\textstyle\int_{L_\theta}d\bm u\right)_1\,d\mathcal R\cos\theta=0,
\end{equation}
where $\partial S=\bigcup_\theta L_\theta$ is a decomposition of the cylindrical boundary in vertical lines at different angle $\theta$. The last equality holds since $\int_{L_\theta}d\bm u$ is independent of $\theta$ as can be shown by contour integration
\begin{equation}
\int_{L_\theta}d\bm u-\int_{L_{\theta'}}d\bm u=\int_\theta^{\theta'}d\bm u\big|^{t=T}-d\bm u\big|^{t=0}=0.
\end{equation}
The $k=2$ term vanishes similarly. Hence we find that the rigidity charge is given entirely by the term corresponding to the line direction. Let $\partial S=\bigcup_t C_t$ be the decomposition of the cylindrical boundary in horizontal circles at different height $t$. Recalling the definition of the Burgers vector, we then obtain
\begin{equation}
\tilde{\nu}^S
=V_{cell}^{-1}\int_0^T\left(\bm P\times\textstyle\int_{C_t}d\bm u\right)_3\, dt=\displaystyle\frac{1}{V_{cell}}(\bm P\times\bm b)\cdot T\hat{\bm\ell}.
\end{equation}
The rigidity charge density along the line can then be given by $\rho_L\equiv\nu_L/T=V_{cell}^{-1}\det(\bm P,\bm b,\hat{\bm\ell})$, where the three vectors $\bm P$, $\bm b$ and $\hat{\bm\ell}$ form the columns of a $3\times3$ matrix.

\section{Dislocation loop dipole moment}

Let $\bm r(s)$ be an arc length parametrization of the loop $L$. The total charge along the loop is given by
\begin{equation}
\nu^L=\oint_L\rho^L(\bm r)ds\propto\oint_L \bm{r}'(s)ds=0.
\end{equation}
As the loop separates charges of different sign, it is sensible to consider the associated dipole moment
\begin{equation}
\bm d^L\equiv\oint_L \bm r\rho^L(\bm r)ds=\oint_L \bm r(s)\frac{(\bm P\times\bm b)\cdot \bm r'(s)}{V_{cell}}ds.
\end{equation}
Without loss of generality we can choose our basis $(\hat{\bm e}_1,\hat{\bm e}_2,\hat{\bm e}_3)$ such that $\bm P\times\bm b$ points in the direction of $\hat{\bm e}_3$. We then obtain:
\begin{equation}
\bm{d}^L=\frac{||\bm P\times\bm b||}{V_{cell}}\oint_L \bm{r} \, dr_3
=\frac{||\bm P\times\bm b||}{V_{cell}}(\mathcal A_{13},\mathcal A_{23}, 0),
\end{equation}
where $\mathcal A_{ij}$ is the signed area of the projection of the dislocation loop on the plane generated by the basis $(\hat{\bm e}_i,\hat{\bm e}_j)$.

\clearpage