text:
\begin{document}

\author{Franz Achleitner\thanks{Technische Universit\"at Wien, Institute for Analysis and Scientific Computing, Wiedner Hauptstraße 8-10 1040 Wien, Austria},~~Christian Kuehn\thanks{Technical University of Munich (TUM), Faculty of Mathematics, 85748 Garching bei M\"unchen, Germany},~~Jens M.~Melenk\footnotemark[1]~~and~~Alexander Rieder\thanks{University of Vienna, Faculty of Mathematics, Oskar-Morgenstern-Platz 1, 1090 Wien, Austria}
}

\title{Metastable Speeds in the Fractional Allen-Cahn Equation}

\maketitle

\begin{abstract}
We study numerically the one-dimensional Allen-Cahn equation with the spectral fractional Laplacian $(-\Delta)^{\alpha/2}$ on intervals with homogeneous Neumann boundary conditions.
In particular, we are interested in the speed of sharp interfaces approaching and annihilating each other.
This process is known to be exponentially slow in the case of the classical Laplacian.
Here we investigate how the width and speed of the interfaces change if we vary the exponent~$\alpha$ of the fractional Laplacian.
For the associated model on the real-line we derive asymptotic formulas for the interface speed and time--to--collision in terms of~$\alpha$ and a scaling parameter~$\varepsilon$.
We use a numerical approach via a finite-element method based upon extending the fractional Laplacian to a cylinder in the upper-half plane,
and compute the interface speed, time--to--collapse and interface width for~$\alpha\in(0.2,2]$.
A comparison shows that the asymptotic formulas for the interface speed and time--to--collision give a good approximation for large intervals.
\end{abstract}

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

In this work we study the fractional Allen-Cahn equation
\be
\label{eq:AC}
\partial_t u = -\eps^\alpha(-\Delta)^{\alpha/2}u+u(1-u^2)=:
-\eps^{\alpha}(-\Delta)^{\alpha/2}u-f(u),
\ee
for $u=u(x,t)$, where $(x,t)\in\Omega\times[0,T)$, $\Omega:=[-L,L]$ is the spatial
domain for $L>0$, $-(-\Delta)^{\alpha/2}$ is the fractional Laplacian for
$\alpha\in(0,2]$, and $0<\eps\ll 1$ is a small parameter.
We assume homogeneous Neumann boundary conditions and an initial condition
\benn
\partial_xu(-L,t)=0=\partial_x u(L,t),\qquad u(x,0)=u_0(x) \,,
\eenn
for some given initial datum $u_0 =u_0(x)$.
We consider the \textit{spectral} fractional Laplacian, which is defined as
\benn
(-\Delta)^{\alpha/2} u:=\sum_{n=0}^{\infty}{\lambda_{n}^{\alpha/2} (u,\varphi_n)_{L^2(\Omega)} \, \varphi_n}
\eenn
where $(\lambda_n,\varphi_n)$ are the eigenvalues and eigenfunctions of the Laplacian with homogeneous Neumann boundary conditions
scaled such that $\|\varphi_n\|_{L^2(\Omega)} = 1$.
Other possible definitions for fractional Laplacians on bounded domains with Neumann-type boundary conditions are discussed and compared in~\cite[\S 7]{DiRoVa:2017} and~\cite[\S 6]{Grubb:2016}.

\smallskip
The potential $F$ associated to $f(u) =u(u^2 -1)$ via $F'(u)=f(u)$ is given by
\benn
F(u):=\frac14 u^4-\frac12 u^2 \,.
\eenn
There are three homogeneous steady states for~\eqref{eq:AC} given by $u_*=-1,0,1$.
Using the linearized problem
\benn
\partial_t U=\underbrace{\left[-\eps^\alpha(-\Delta)^{\alpha/2}
-\txtD_u f(u_*)\right]}_{:=L(u_*)}U,\qquad U=U(x,t),
\eenn
one easily checks that $u_*=\pm1$ are locally asymptotically stable since the spectrum of $L(u_*)$ is contained in $\{z\in\C:\Re(z)<0\}$.

The state $u_*=0$ is unstable.
Indeed, one can also view $u_*=\pm 1$ as global minima and $u_*=0$ as a local maximum of the potential $F$.

For the case of $\alpha=2$, so that $-(-\Delta)^{\alpha/2}=\Delta$, even more global dynamics of~\eqref{eq:AC} is well understood; see Section~\ref{subsec:classical} for a more detailed technical review.
Here we just emphasize that the Allen-Cahn equation~\eqref{eq:AC} can exhibit metastability as shown in Figure~\ref{fig:alpha2}.
\begin{figure}[h]
  \begin{subfigure}[b]{0.32\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation/out_0.pdf}
    \caption{solution at time $t=0$.}
  \end{subfigure}
  \begin{subfigure}[b]{0.32\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation/out_1771.pdf}
    \caption{solution at time $t=1771$.}
  \end{subfigure}
  \begin{subfigure}[b]{0.32\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation/out_177216.pdf}
    \caption{solution at time $t=177216$.}
  \end{subfigure}

  \caption{Evolution of solution for the classical Allen-Cahn equation~\eqref{eq:AC}
    for~\eqref{eq:TL:simple} with $L=5$, $\eps=0.3$ and $\alpha=2$.}
  \label{fig:alpha2}
\end{figure}

We observe that suitable initial conditions $u_0$ get attracted very quickly to solutions composed of multiple, say $n$, very sharp interfaces.
Each interface is known to have width $\cO(\eps)$ for $\alpha=2$; see Section~\ref{subsec:classical}.
For a long time, this $n$-interface solution appears to be stationary, but it evolves on an exponentially long time scale $\cO(\txte^{K/\eps})$ for some constant $K>0$; see again Section~\ref{subsec:classical} for details.
On this exponentially long time scale the interfaces move towards each other and then annihilate; see Figure~\ref{fig:alpha2}.
This effect is the essence of \emph{metastable behavior}, i.e., apparent stationarity on extremely long, yet still transient, time scales.

\medskip
In this work, we are interested in the influence of the parameter $\alpha$ on metastable interface motion.
It is natural to ask how the order of the fractional Laplacian influences the interface speed as well as  the interface width.
For $\alpha=2$, these effects have been quantified, see, e.g.,~\cite{CarrPego}.
However, for $\alpha=(0,2)$ a precise quantification has not been carried out yet.

However, this problem has been studied for~\eqref{eq:AC} with the fractional Laplacian on the real line.
For $\alpha\in(1,2)$ the dependence of the interface speed has been quantified by heuristic arguments and verified numerically via a pseudo-spectral method~\cite{NecNepGol:2008}.
Moreover, equation~\eqref{eq:AC} with periodic potential $F$ and external forcing on the real line has been studied as a model in crystal dislocation dynamics, see~\cite{GoMo:2012, DiFiVa:2014, DiPaVa:2015, PaVa:2015, PaVa:2016, PaVa:2017}.
They characterize the long-time behavior of solutions which includes metastable scenarios and give rigorous proofs.

It is not obvious if these results carry over to the fractional Allen-Cahn equation~\eqref{eq:AC} with homogeneous Neumann boundary conditions.
For example, the fractional Laplacian on the real line is uniquely defined, e.g. see~\cite{Kwa:2017}, whereas on bounded intervals several distinct definitions exist, see~\cite{SerVal:2014, Grubb:2016, DiRoVa:2017, Lischke+etal:2020}.
However our numerical studies of the interface motion for~\eqref{eq:AC} with homogeneous Neumann boundary conditions indicate that the interface speed is governed by the same asymptotic dependence on~$\eps$ and~$\alpha$.
Additionally, we characterize the dependence of the interface width on these parameters.
In particular, our main observations are:
\begin{itemize}
 \item[(R1)] \textit{Interface speed:}
  For $\alpha\in(0.5,2)$, the magnitude of the interface speed in case of two interfaces is approximately
  \begin{equation} \label{interface.speed}
   |s (\eps,\alpha)|
   \approx \eps^{1+\alpha}\ \frac{4}{\alpha} \frac{2^{\alpha} \Gamma((1+\alpha)/2)}{\sqrt{\pi} |\Gamma(-\alpha/2)|}\ \gamma\ \frac{1}{|x_1 -x_2|^{\alpha}} \,.
  \end{equation}
  where $x_1$ and $x_2$ are the centers of the interfaces and $\gamma := (\int_\R (v'(x))^2 \txtd x)^{-1}$ is the inverse of the semi-norm of a basic layer solution~$v$ of~\eqref{eq:BLS:1}--\eqref{eq:BLS:2}.
 \item[(R2)] \textit{Interface width:}
 For $\alpha\in(0.2,2)$, the interface width is approximately
 \begin{equation} \label{interface.width}
  w(\eps,\alpha) = b(\alpha) \eps^{\kappa_1 \alpha^{-1} +  \kappa_2}\,,
 \end{equation}
where we numerically estimated $\kappa_1\approx -0.168298$ and $\kappa_2\approx 1.11709$. We note that the numerically estimates suggest the \emph{conjectures} $\kappa_1=1/6$ and $\kappa_2=10/9$.
\end{itemize}

The paper is structured as follows:
In Section~\ref{sec:analysis}, we briefly review analytical results for the Allen-Cahn equation involving the Laplacian ($\alpha=2$) with homogeneous Neumann boundary conditions and the fractional Laplacian ($\alpha\in(0,2)$) on the real line, respectively.
In Section~\ref{sec:numerics}, we describe the numerical setup used to simulate~\eqref{eq:AC} for $\alpha\in(0,2]$.
Moreover, we compare our numerical method with a spectral method in~$\S$\ref{subsec:comparision_methods}.
The main results are presented and discussed in Section~\ref{sec:results}.
A brief summary and an outlook to future open problems can be found in Section~\ref{sec:outlook}.

\section{Asymptotic analysis}
\label{sec:analysis}

\subsection{The Classical Case}
\label{subsec:classical}

In this section, we briefly review the known results for metastability of the
classical Allen-Cahn equation
\be
\label{eq:AC1}
\partial_t u = \eps^2\partial_{xx} u+u(1-u^2),
\ee
with Neumann boundary conditions posed on the interval $[0,1]$ as studied
by several authors; here we mainly follow~\cite{CarrPego}.
Note that considering the interval $[0,1]$ is equivalent to picking $[-L,L]$ up to shifting and
scaling the $x$-coordinate.
A key auxiliary tool to construct $n$-interface solutions is to first construct a single layer.
This requires the solution of the auxiliary two-point boundary value problem
\be
\label{eq:BVP}
\eps^2 \frac{\txtd^2 \phi}{\txtd x^2}=\phi(1-\phi^2),\qquad
\phi(-\ell/2)=0=\phi(\ell/2).
\ee
For $\eps>0$ sufficiently small, it can be proven that~\eqref{eq:BVP} has
a unique solution $\phi(x,\ell)$ for each given $\ell>0$, which is positive for
$|x|<\ell/2$.
The shape of $\phi(x,\ell)$ corresponds to a spike centered at $x=0$, i.e., an interior layer solution in the terminology of multiscale analysis~\cite{KuehnBook}.
The width of this spike is $\cO(\eps)$ and it can also be interpreted as a homoclinic orbit in $(\phi,\phi')$-coordinates to $(0,0)$.
Next, one wants to construct an interface (layer) with location $h=h(t)$ for~\eqref{eq:AC1}.
Consider a monotone function
\benn
\xi\in C^\I(\R,[0,1]),\qquad
\xi(x)=\left\{
\begin{array}{ll}
0\quad\text{for $x\leq -1$,}\\
1\quad\text{for $x\geq 1$,}\\
\end{array}
\right.
\eenn
and define the \emph{approximate metastable interface},
$u^h=u^h(x)$ for some $0<h<1$ by
\be
\label{eq:hlayer}
u^h(x):=-\left[1+\xi\left(\frac{x-h}{\eps}\right)\right]
\phi(x,2h)+\xi\left(\frac{x-h}{\eps}\right)\phi(x-1,2-2h)
\ee
so that $u^h(0)<0$, $u^h(h)=0$, $u^h(1)>0$, and the width of the interface is $\cO(\eps)$.
Obviously one can construct an approximate $n$-interface solution by a similar procedure.
Since any $n$-interface solution has $n$ natural coordinates given by the positions $h$ of each interface,
one can define an $n$-dimensional manifold $\cM$ in $H^1([0,1])$ parametrized by the positions.
It is shown in~\cite{CarrPego} that this manifold is locally attracting for large classes of initial data for~\eqref{eq:AC1}.
Once the solution is close to the approximately invariant manifold $\cM$, there exists an ODE describing the motion of the positions $h$ on $\cM$.
It turns out that all interface positions move at equal asymptotic speed in a generic setting away from any annihilation events.
Their speed in this regime can be determined by just looking at one interface.
It is highly non-trivial to prove that the ODE for its position is given by
\be
\label{eq:hmotion}
h'=K_0\eps \left[q\left(\frac{\eps}{2h}\right)
-q\left(\frac{\eps}{2-2h}\right)\right]+R_1,\qquad q(r):=F(\phi(0,\ell)),~
r:=\eps/\ell,
\ee
where $K_0>0$ is a computable order $\cO(1)$ constant, and $R_1$ turns out to be a
negligible remainder term~\cite{CarrPego}.
Metastability arises because one can also prove that
\benn
q\left(r\right)=K_1\txte^{-K_2 \ell/\eps}+R_2,
\eenn
for constants $K_{1,2}>0$ and $R_2$ turns out to be another higher-order terms irrelevant for computing the leading-order of the vector field in~\eqref{eq:hmotion}, i.e., the local interface speed is exponentially small for $\alpha=2$.

\subsection{Analysis of fractional Allen-Cahn equations}
For $\alpha\in(0,2)$, the results for the classical case lead to the following questions:
\begin{itemize}
 \item[(Q1)] Does the metastable interface speed $s=s(\eps,\alpha)$ change?
 If so, what is the graph of $s(\eps,\alpha)$ as a function of $\alpha$ and/or $\eps$?
 \item[(Q2)] Does the interface width $w=w(\eps,\alpha)$ change?
 If so, what is the graph of $w(\eps,\alpha)$ as a function of $\alpha$ and/or $\eps$?
\end{itemize}
The questions (Q1)-(Q2) are our main focus as they provide immediate information
on the relevant dynamical behavior.

We discuss the ambitious task to extend the results by Carr and Pego in the Section~\ref{sec:outlook}.
Next, we will present the analysis of a fractional Allen-Cahn equation on $\R$,
which provides approximations for the interface speed and time--to--collapse.
In Section~\ref{sec:results} we show by numerical simulations that these formulas also give a good approximation for our model on a bounded interval.

\subsubsection{Analysis of fractional Allen-Cahn equations on $\R$}
Equation~\eqref{eq:AC} with periodic potential $F$ and external forcing on the real line $x\in\R$ models the dynamics of crystal dislocations, see~\cite{GoMo:2012, DiFiVa:2014, DiPaVa:2015, PaVa:2015, PaVa:2016, PaVa:2017}.
First, the evolution of dislocations given as a superposition of transitions with the same orientation has been studied. It has been shown that these transitions repel each other~\cite{GoMo:2012, DiFiVa:2014, DiPaVa:2015}.
Patrizi and Valdinoci~\cite{PaVa:2015, PaVa:2016, PaVa:2017} considered also dislocations given as a superposition of transitions with arbitrary orientations and studied again the long-time behavior.
For example, two transitions of opposite orientation attract each other (if no external force is present), which shows the metastable behavior we are interested in.
In general, they study the long-time behavior of well--prepared initial data modeling an arbitrary (but finite) number of transitions with any order of orientations.
We report here their results for the equation
\be
\label{eq:AC:rescaled}
\partial_t u_\eps = \tfrac1{\eps} \Big(-(-\Delta)^{\alpha/2}u_\eps -\tfrac1{\eps^\alpha} F'(u_\eps)\Big),
\ee
on $x\in\R$ with (non-periodic) potential $F(u) =\tfrac{u^4}4 -\tfrac{u^2}2$ and superposition of transitions with alternating orientations (and note the minor modifications needed along the way):

\medskip
Equations~\eqref{eq:AC} and~\eqref{eq:AC:rescaled} are related via a rescaling of time $t\mapsto t/\eps^{1+\alpha}$.

\medskip
First, a well--prepared initial datum $u_0$ is constructed from basic layer solutions of~\eqref{eq:AC} (instead of initial layer solutions):
A basic layer solution~$v:\R\to[-1,1]$ is a stationary solution of~\eqref{eq:AC}, i.e. it solves
\begin{equation} \label{eq:BLS:1}
 -(-\Delta)^{\alpha/2}v -f(v) = 0 \,, \qquad x\in\R\,,
\end{equation}
and satisfies
\begin{equation} \label{eq:BLS:2}
 v'(x)>0 \text{ for all } x\in\R\,, \quad v(-\infty)=-1\,, \quad v(0)=0\,, \quad v(+\infty)=1\,.
\end{equation}
Here, with a slight abuse of notation $-(-\Delta)^{\alpha/2}$ denotes the fractional Laplacian acting on functions on the real line.
Note that the fractional Laplacian given on the real line is equivalent to a singular integral representation, see~\cite{CafSil:2007,Kwa:2017}.
The existence and uniqueness (due to $v(0)=0$) of a basic layer solution has been proved in~\cite[Thm. 2.4]{CabSir:2015}.

\medskip
Then, a well--prepared initial datum $u_0$ is constructed from shifted basic layer solutions~$v$ of alternating orientation.
For example, we study the evolution of a transition layer given by
\begin{equation} \label{eq:TL:simple}
 u_\eps^0(x)
 := v\bigg(\frac{x-x_1^0}{\eps}\bigg) + v\bigg(-\frac{x-x_2^0}{\eps}\bigg) +1
\end{equation}
for some $x_2^0 < x_1^0$.
More generally, we study the evolution of an initial datum $u_\eps^0$ given by
\begin{equation} \label{eq:TL}
 u_\eps^0(x) := \sum_{i=1}^{2K} v_{\eps,i}^0(x) +1
\end{equation}
where $K\in\N$, $v_{\eps,i}^0(x) = v\Big(\zeta_i \frac{x-x_i^0}{\eps}\Big)$ with $x_{2K}^0 < x_{2K-1}^0 < \ldots < x_2^0 < x_1^0$ and
\[ \zeta_i = \begin{cases} +1 &\text{if $i$ is odd,} \\ -1 &\text{if $i$ is even,} \end{cases}
     \qquad i\in\{1,2,\ldots,2K\} \,.
\]

It can be shown that, for sufficiently small $\eps>0$, the solution $u_\eps(x,t)$ of~\eqref{eq:AC:rescaled} with initial datum $u_\eps^0$ is approximately of the form
\begin{equation}
  u_\eps(x,t)
   \approx \sum_{i=1}^{2K} v\Big(\zeta_i \frac{x-x_i(t)}{\eps}\Big) +1 \,,
\end{equation}
where the functions $x_i(t)$ satisfy approximately (i.e., up to higher order terms in $\eps$) a system of ODEs
\begin{equation} \label{ODE:centers}
 \begin{cases}
   \frac{\txtd x_i}{\txtd t} = \frac{\eps^{1+\alpha}}{\alpha} 4 \frac{2^{\alpha} \Gamma((1+\alpha)/2)}{\sqrt{\pi} |\Gamma(-\alpha/2)|} \gamma \sum_{i\ne j} \zeta_i \zeta_j \frac{x_i -x_j}{|x_i -x_j|^{1+\alpha}} \,, \\
   x_i(0) = x_i^0 \,,
 \end{cases}
\end{equation}
for $i\in\{1,2,\ldots,2K\}$ and $\gamma := (\int_\R (v'(x))^2 \txtd x)^{-1}$.
The centers $x_i$ are initially ordered as $x_{2K}^0 < x_{2K-1}^0 < \ldots < x_2^0 < x_1^0$, then for sufficiently small times $t>0$ this order will persist for the solution of~\eqref{ODE:centers}, i.e.
\begin{equation} \label{ineq:centers}
  x_{2K}(t) < x_{2K-1}(t) < \ldots < x_2(t) < x_1(t) \,.
\end{equation}
System~\eqref{ODE:centers} is well-defined until the first collision at time $T_C$, which is defined as the time such that~\eqref{ineq:centers} holds for all $t\in[0,T_C)$ and
\begin{equation} \label{cond:TC}
 \exists i_C \in\{1,2,\ldots,2K-1\}\,: \qquad
 x_{i_{C} +1}(T_C) =x_{i_C}(T_C)\,.
\end{equation}

More precisely, it is proved that the solution $u_\eps$ approaches a superposition of sharp transitions with moving centers $x_i(t)$ in the following sense:
\begin{thm}[cf. {\cite[Thm. 1.1]{PaVa:2015}}] \label{thm:ueps}
Suppose $\alpha\in[1,2]$.
Let
\begin{equation}\label{eq:AL}
 v(t,x) = \sum_{i=1}^{2K} \sgn(\zeta_i(x-x_i(t))) +1\,,
\end{equation}
where $\sgn$ is the sign function and $(x_i(t))_{i=1,\ldots,2K}$ is the solution to~\eqref{ODE:centers}.
Then, for every $\eps>0$ there exists a unique solution $u_\eps$ of~\eqref{eq:AC:rescaled}.
Furthermore, as $\eps\to0^+$, the solution $u_\eps$ exhibits the following asymptotic behavior:
\begin{align}
 \label{est:above}
 \limsup_{\stackrel{(t',x')\to(t,x)}{\eps\to 0^+}} u_\eps(t',x')
 \leq \limsup_{(t',x')\to(t,x)} v(t',x')
\intertext{and}
 \label{est:below}
 \liminf_{\stackrel{(t',x')\to(t,x)}{\eps\to 0^+}} u_\eps(t',x')
 \geq \liminf_{(t',x')\to(t,x)} v(t',x')
\end{align}
for any $(t,x)\in[0,T_C)\times\R$.
\end{thm}
\begin{proof}[Sketch of proof]
 The initial conditions $u_\eps^0$ under consideration take values in $[-1,1]$.
 Due to a maximum principle for~\eqref{eq:AC:rescaled}, the corresponding solution $u_\eps$  of~\eqref{eq:AC:rescaled} with initial datum $u_\eps^0$ will again satisfy $-1\leq u_\eps (t,x)\leq 1$ for all $(t,x)$.
 We consider~\eqref{eq:AC:rescaled} with the periodic potential $F_{per}$ such that $F_{per}(u) =F(u)$ for all $u\in[-1,1]$.
 Then $F_{per}\in C^{2,\beta}(\R)$ is H\"older continuous.
 Now we can use the result of~\cite[Thm. 1.1]{PaVa:2015} for the modified problem.
 However, since $-1\leq u_\eps (t,x)\leq 1$ for all $(t,x)$, the statement carries over to the original problem unchanged.
\end{proof}
\begin{remark}
\begin{enumerate}
\item The restriction to $\alpha\in[1,2]$ seems to be due to technical conditions in the proofs. In case of $\alpha\in(0,1)$, the assumption $F\in C^{3,\beta}(\R)$ for some $\beta>0$ is needed, see~\cite{DiFiVa:2014}.
\item In the studies on crystal dislocation dynamics the stable stationary states are taken to be 0 and 1 compared to -1 and 1 as in our setting.
That means, we have to include e.g. an affine transformation $y \mapsto -2y+1$ which gives an additional multiplicative constant~$4$.
\item
Moreover, we have to include a normalization factor $\frac{2^{\alpha} \Gamma((1+\alpha)/2)}{\sqrt{\pi} |\Gamma(-\alpha/2)|}$, due to the singular integral representations used in~\cite[Lemma 5.1]{stinga_torrea} compared with~\cite[(1.2)]{PaVa:2015}.
\end{enumerate}
\end{remark}

In case of two initial transitions with opposite orientation the following estimate for the time--to--collision~$T_C$ holds
\begin{thm}[{cf. \cite[Thm. 1.2]{PaVa:2015}}] \label{thm:TC}
 Let $K=1$. Let $d_0 := x_1^0 -x_2^0 >0$.
 Then,
 \begin{equation} \label{eq:TC}
 T_C = \tfrac{C_\alpha d_0^{1 +\alpha}}{2(1 +\alpha)\gamma}
 \qquad\text{with }
  C_\alpha
  :=1/\bigg( \frac{4}{\alpha} \frac{2^{\alpha} \Gamma((1+\alpha)/2)}{\sqrt{\pi} |\Gamma(-\alpha/2)|}\bigg)
  =\frac{\Gamma(1-\alpha)}{2\sec(\alpha \frac{\pi}{2})}\,.
 \end{equation}
\end{thm}

Then,~\eqref{est:above} and~\eqref{est:below} imply that for any $x\ne x_c$ where $x_c:=x_1(T_C)=x_2(T_C)$, we have
\[ \lim_{t\to T_C^-} \lim_{\eps\to 0^+} u_\eps(t,x) =1 \,. \]
That means two dislocations annihilate each other after their collision.
However, the limit of $u_\eps(t,x)$ keeps a memory of them, in the sense that $u_\eps$ at the point $x_c$ does not approach the state~$1$ in the limit.
\begin{thm}[{\cite[Thm. 1.4]{PaVa:2015}}] \label{thm:xc}
 Let $K=1$. Let $u_\eps$ be the solution to~\eqref{eq:AC:rescaled}, then
 \begin{equation}\label{est:limsup:below}
 \liminf_{\stackrel{t\to T_C^{-}}{\eps\to 0^+}} u_\eps(t,x_c) \leq -1 \,.
 \end{equation}
\end{thm}
Nonetheless, after some time $\Teps$ -- which is only slightly larger than the collision time $T_C$ -- the solution $\ueps$ will become small like $\eps$. More precisely,
\begin{thm}[cf. {\cite[Thm. 1.1 \& Thm. 1.2]{PaVa:2016}}] \label{thm:Teps}
Let $K=1$.
Under the assumptions of Theorem~\ref{thm:ueps}, let $\ueps$ be a solution of~\eqref{eq:AC:rescaled}.
Then, there exists $\eps_0>0$ such that for all $\eps<\eps_0$ there exist $\Teps,\rho_\eps>0$ such that
\[ \Teps = T_C +o(1)\,, \quad \rho_\eps =o(1) \quad \text{as $\eps\to 0$} \]
and
\begin{equation}
 1-\ueps(\Teps,x) \leq \rho_\eps \quad \text{for all $x\in\R$.}
\end{equation}
Moreover, there exist $\tilde{\eps}_0>0$ and $c>0$ such that for any $\eps<\tilde{\eps}_0$ we have
\begin{equation}
 |1-\ueps(t,x)| \leq \rho_\eps \exp\big( c\tfrac{\Teps -t}{\eps^{1+\alpha}}\big)\,, \quad \text{for all $x\in\R$ and $t\geq\Teps$.}
\end{equation}
\end{thm}

\section{Numerical Methods}
\label{sec:numerics}

In recent years the numerical simulations of partial differential equations with fractional space derivatives has attracted a lot of attention.
Some references on (various forms of) fractional Allen-Cahn equations in one and multi-dimensional spatial settings are the following:~\cite{BurHalKay:2012,Zeng+etal:2014,StoYue:2015,AchleitnerKuehn2,SimYanMor:2015,SonXuKar:2016,ZhaWenFen:2016,HouTanYan:2017,AlzKha:2019,LeeLee:2019,LiuCheWan:2020,HePanHu:2020}.

In order to discretize \eqref{eq:AC}, we consider an alternative formulation of the fractional Laplacian, based on the Caffarelli-Silvestre extension~\cite{CafSil:2007}, as presented in \cite{stinga_torrea} for the Neumann problem.
Consider the solution $\UU: \R_+ \times \Omega  \times \R_+ \to \R$ to:
\begin{subequations}
  \label{eq:lifted_eq}
\begin{align}
  -\operatorname{div}\left( y^{1-\alpha} \UU(t) \right) &=0\,,  && \text{ in $\Omega \times \R_+$ },\ \forall t>0\,, \\
  d_\alpha \trace \UU_t + \eps^2 \partial_\nu^\alpha \UU   &=-d_\alpha f(\trace \UU )\,, && \text{ in $\Omega \times \{0\} \times \R_+$ }, \\
  \partial_\nu \UU(t) &=0\,,  &&\text{ on $\partial \Omega \times \R_+$},\ \forall t>0\,,\\
  \trace \, \UU(0)&=u_0\,, && \text{ in $\Omega$.}
\end{align}
\end{subequations}
Here $d_\alpha:=2^{1-\alpha}\Gamma(1-\alpha/2)/\Gamma(\alpha/2)$, $\partial_\nu$ denotes the
normal derivative with respect to the lateral boundary,
$\trace $ denotes the trace operator with respect to the artificial variable $y$ at $0$  and $\partial_\nu^\alpha$ is defined by
\begin{align*}
  \partial_\nu^{\alpha}  \UU:=  -\lim_{y\to 0^+}{y^{1-\alpha} \partial_y \, \UU(\cdot,y)}.
\end{align*}
It can then be proven that if $\UU$ solves (\ref{eq:lifted_eq}), then $\trace \UU$ solves~\eqref{eq:AC}.

The discretization is based on an semi-implicit BDF2 discretization of the time variable, and a $hp$-finite element method to discretize the $x$ and $y$ variables, similar to what was presented in~\cite{tensor_fem} for the linear stationary case and~\cite{pde_frac_parabolic} for the linear time-dependent fractional diffusion problem, but which is based on a low order
discretization in space. An analogous $hp$-type method for the linear fractional heat equation,
but using different time-stepping and Dirichlet boundary condition was analyzed in detail in
\cite{fractional_heat}.

The semi-implicit BDF2 (SBDF2) method for the ODE  $u_t = \mathscr{L} u + f(t,u)$ with stiff operator $\mathscr{L}$ is given by
\begin{align*}
u^{n+1}&=\frac{4}{3} u^n - \frac{1}{3} u^{n-1} + \frac{2}{3} \Delta t\left[2 f(t_n,u^n) - f(t_{n-1},u^{n-1})\right] + \frac{2}{3} \Delta t \mathscr{L} u^{n+1}.
\end{align*}
For the discretization in the auxiliary variable $y$, we fix $\sigma \in (0,1)$ and
consider a graded mesh $\cT_y:=\{y_0,\dots y_N\}$ with $y_0:=0$, $y_{j}:=\cY \sigma^{N-j}$. Here $\cY>0$ is
a cutoff parameter. In $\Omega:=(-L,L)$, we consider a uniform mesh with size $h>0$, denoted by $\cT_x$.
Denoting by $\cS^{p,1}(\cT)$ the space of continuous piecewise polynomials of degree $p$, on  a triangulation $\cT$ we solve in each time step
for $\UU_h^{n+1} \in \cS(\cT_x) \otimes \cS(\cT_y)$ satisfying
\begin{multline*}
  \frac{2\Delta t}{3d_\alpha}\int_{0}^{\cY}{y^{1-\alpha}\left(\nabla_{xy} \UU_h^{n+1}(t),\nabla_{xy} \VV_h\right)_{L^2(\Omega)}\,dy}
  +\left(\trace \UU_h^{n+1},\trace {\VV_h}\right)_{L^2(\Omega)} \\
  = \Big(\frac{4}{3}\trace \UU_h^{n} - \frac{1}{3}\trace \UU_h^{n-1},
    \trace {\VV_h}\Big)_{L^2(\Omega)}
    +\frac{2\Delta t}{3} \left(f(\trace{\UU}^{n}) - f(\trace{\UU}^{n-1}),\trace \VV_h\right)_{L^2(\Omega)}.
\end{multline*}
for all test functions $\VV_h \in \cS(\cT_x) \otimes \cS(\cT_y)$.
Note that the nonlinearity is treated explicitly.
Since we expect this contribution to be non-stiff this should not cause concern.
The fact that we replaced the infinite cylinder $\Omega \times \R_+$ with $\Omega \times (0,\cY)$ can be justified by the fact that
when expanding $\UU$ in an eigenbasis of the Laplacian, all contributions (except for the constant one) decay exponentially.
(See~\cite{tensor_fem,fractional_heat} for details in the Dirichlet case, the Neumann case behaves analogously).
The cutoff of the constant contribution does not impact the result since we used homogeneous Neumann conditions on the artificial boundary, which capture the constant mode exactly.
When referencing this algorithm, we will refer to it as the ``$hp$-method''.

\smallskip
In order to be able to compare the results obtained by our $hp$-FEM type approach, we also implemented and ran our experiments using a simpler spectral type method
proposed in \cite{spectral_code}.
It directly exploits the spectral definition of the fractional Laplacian by using a discrete cosine transform.
For time discretization we use an IMEX-type Euler method.
When referencing this algorithm, we will refer to it as the ``spectral--method''.

\subsection{Detecting an interface and measuring speeds}
In order to extract information on the behavior of the interfaces, we need to fix our methodology for
extracting the quantities of interest from the numerical solutions.
For simplicity, even when using a high order method for solving~\eqref{eq:AC} numerically, our measurements are always taken from a
piecewise linear approximation to the numerical solution.
This is done by interpolation on a finer grid.

Since the interfaces need a certain time to form, and we expect the interfaces to collapse after some time, with the scales of these timings depending on $\eps$ and $\alpha$, we use an (empirically determined) estimate for the time  until the collapse happens $t_{col}$, i.e., the solution becomes equal to the constant $1$ function.
When measuring we then assume that for $\frac{1}{4}t_{col}$ the interface has already formed and assume that afterwards, the interfaces remained stable for a duration of  $\frac{1}{10}t_{col}$ and we can do accurate measurements.

When computing the speed of the interface, we consider the speed of the zeros of the numerical solution.
In order to determine the width of the interface, we made the very simplistic assumption that an interface happens whenever the solution dips between $\max(u)-\delta$ and $\min(u)+\delta$.
For the images chosen here we used $\delta:=0.1$.
This gave good correspondence of the wave speed determined by the zeros of $u$ and the speed of the midpoint of the interface.
We also tested using $\delta=0.01$ and $\delta \approx \varepsilon^{\alpha}$, but $\delta=0.1$ seemed to give the most robust results.
Both width and speed are computed by taking the average over 100 measurements in the time-frame $\big(\frac{1}{4}t_{col},\frac{7}{20}t_{col}\big)$. Figure~\ref{fig:ex_interfaces} shows a typical situation.

\begin{figure}[h]
  \center
  \includegraphics[width=12cm]{figures/interfaces.png}
  \caption{Example for the detected interfaces and roots. }
\label{fig:ex_interfaces}
\end{figure}

In order to determine rates, we used MATLABs \textit{Curve Fitting Toolbox} for polynomial and power-law fits.

\subsection{Model problem}
For our computations we fixed the following parameters.
We considered the domain $\Omega:=(-L,L)$ with $L:=10$.
In order to derive a starting condition which is close to the theoretical
considerations of Section~\ref{sec:analysis}, we approximately computed
a layer solution $v$, as described in Section~\ref{subsec:layer_solution_on_R}.
Tacitly extending the function $v$ by a constant outside of $(-L,L)$,
the initial condition is then given by
\begin{align*}
  u_0(x):=v\Big(\frac{x-x_0}{\eps}\Big) + v\Big(\frac{-x-x_0}{\eps}\Big) +1.
\end{align*}
This function possesses two interfaces at $\pm x_0$, which was taken to be $x_0:=L/4$.
We will later on also discuss what happens for more general initial conditions.
We used a timestepping size of $\Delta t:=\min(0.05,10^{-4}t_{col})$, ensuring that
we make at least $10^4$ steps.

For the spectral method we used a grid with $N=262144$ points.
For the $hp$-method, in $\Omega$, we used a uniform grid of mesh size $\min(2L/\eps,8012)$ and polynomial degree $p_x=5$ in order to resolve the appearing steep flanks of the interfaces well. To discretize the artificial
variable, we used $p_y=8$ on a geometric grid with $8$ layers and a mesh grading factor of $\sigma=0.125$,
i.e., the grid consists of points of the form $8 (0.125)^{\ell}$ for $\ell=0, \dots 8$.

\subsection{Comparison of the two methods}
\label{subsec:comparision_methods}
In this section, we compare the different methods we used.
Since they are quite different (first order vs second order in time, different approaches to the fractional Laplacian), we expect that this
comparison gives a good indication of the real accuracy of our simulations.
Instead of comparing the discrete solutions directly, we compare the postprocessed quantities: interface speed, interface width and time--to--collapse.
The results are collected in Figure~\ref{fig:comparison_of_methods}.
As we can see, the error is reasonably small for all contributions, giving confidence that the behavior observed numerically matches the one of the exact solution. The only exception is for small values of $\alpha \approx 0.2$.  As we do not expect our discretization to be robust as $\alpha\to 0$ this has to be expected. At the same time, this behavior may explain the observed mismatch between our predicted behavior and the observations in this regime.

\begin{figure}[h]
  \begin{subfigure}[b]{0.5\textwidth}
    \includeTikzOrEps{figures/comparison_hp_bdf2_l10_new_initial_condition_spectral_l10_new_initial_condition_s_speed}
     \caption{Relative difference of the renormalized interface speed between the methods.}
  \end{subfigure}\qquad
  \begin{subfigure}[b]{0.5\textwidth}
    \includeTikzOrEps{figures/comparison_hp_bdf2_l10_new_initial_condition_spectral_l10_new_initial_condition_s_time}
    \caption{Relative difference of the time--to--collapse between the methods.}
  \end{subfigure} \\\vspace{2mm}
  \begin{subfigure}[b]{0.4\textwidth}
    \includeTikzOrEps{figures/comparison_hp_bdf2_l10_new_initial_condition_spectral_l10_new_initial_condition_s_width}
    \caption{Relative difference of the interface width between the methods.}
  \end{subfigure}
  \caption{Comparison of two numerical methods: \textit{hp}--method vs. spectral--method.
    Different lines represent different values of $\alpha \in (0.2,1.9)$. }
  \label{fig:comparison_of_methods}
\end{figure}

\section{Results on Metastability}
\label{sec:results}
In this section, we consider the evolution of two interfaces.
We present the results of our numerical experiments on intervals $[-L,L]$ and compare them with our analysis on the real line.
Therefore, it is natural to investigate the dependence of our quantities of interest on the size of the considered interval $[-L,L]$.
As can be seen in Figure~\ref{fig:dependence_on_domain}, neither the (modified) speed nor the width of the interface seem to depend strongly on the size of the domain.
Especially for larger domains, the results become indistinguishable.

\begin{figure}
  \centering
  \begin{subfigure}[b]{0.45\textwidth}
    \includeTikzOrEps{figures/old/domain_comparison_mod_speed}
    \caption{Modified interface speed $\widehat{s}$.}
    \label{fig:comp_speed_on_intervals}
  \end{subfigure}
  \begin{subfigure}[b]{0.45\textwidth}
    \includeTikzOrEps{figures/old/domain_comparison_width}
    \caption{Interface width $w$.}
    \label{fig:comp_width_on_intervals}
  \end{subfigure}
  \caption{Comparing results on intervals~$[-L,L]$ with $L=10,20,40$ and $\eps=0.1$.}
  \label{fig:dependence_on_domain}
\end{figure}

We illustrate the metastable behavior for solutions of the Allen-Cahn equation~\eqref{eq:AC} with well--prepared initial data~\eqref{eq:TL:simple} with $L=10$, $\eps=0.01$ and $\alpha=0.9$ in Figure~\ref{fig:metastability}. In this simulation, we can define a time--to--collision $\widetilde{T}_C$ as the first time that two interfaces collide, i.e. $\widetilde{T}_C$ is the time at which the solution becomes non-negative.
This happens at some time between $t=3671.53$ and $t=3694.24$, see Figure~\ref{fig:metastability}.
We studied the evolution of interfaces also via the ODE system~\eqref{ODE:centers} which governs the evolution of the centers of interfaces.
There the time--to--collision $T_C$ is defined as the first time that two centers collide, see~\eqref{cond:TC}.
Using the parameters of our simulation, $T_C$ as given in~\eqref{eq:TC} is approximately $3587.10$, where we used $d_0=2$ as in the simulation and computed $\gamma=1.28550$ as in Appendix~\ref{subsec:layer_solution_on_R}. Comparing $\widetilde{T}_C$ with $T_C$ shows that the collision of interfaces takes place at a little later than expected.

\begin{figure}[h]
  \begin{subfigure}[b]{0.5\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation2/out_0.pdf}
    \caption{solution at time $t=0$.}
  \end{subfigure}\qquad
  \begin{subfigure}[b]{0.5\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation2/out_1.pdf}
    \caption{solution at time $t=75.71$.}
  \end{subfigure}\qquad
  \begin{subfigure}[b]{0.5\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation2/out_30.pdf}
    \caption{solution at time $t=2271.05$.}
  \end{subfigure}
  \begin{subfigure}[b]{0.5\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation2/out_150.pdf}
    \caption{solution at time $t=3596.34$.}
  \end{subfigure}\\\vspace{2mm}
  \begin{subfigure}[b]{0.5\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation2/out_485.pdf}
    \caption{solution at time $t=3671.53$.}
  \end{subfigure}
  \begin{subfigure}[b]{0.5\textwidth}
    \includegraphics[width=\textwidth]{figures/simulation2/out_488.pdf}
    \caption{solution at time $t=3694.24$.}
  \end{subfigure}
  \caption{Evolution of solution for Allen-Cahn equation~\eqref{eq:AC} with well--prepared initial data~\eqref{eq:TL:simple} with $L=10$, $\eps=0.01$ and $\alpha=0.9$.}
  \label{fig:metastability}
\end{figure}

Next, we will first discuss the interface speed $s=s(\eps,\alpha)$ and, then, the interface width~$w=w(\eps,\alpha)$.

\subsection{Interface speed~$s=s(\eps,\alpha)$}
\label{sec:interface_speed}

Due to the analysis of the evolution of interfaces of~\eqref{eq:AC} on the real line, the (centers of) interfaces move according to the ODE system~\eqref{ODE:centers}.
In the case of two interfaces with centers $x_1$ and $x_2$, the interface speed satisfies approximately
\begin{equation} \label{speed_R}
 \begin{split}
 s_{\R} (\eps,\alpha)
  &\approx \Big| \frac{\eps^{1+\alpha}}{\alpha} 4 \frac{2^{\alpha} \Gamma((1+\alpha)/2)}{\sqrt{\pi} |\Gamma(-\alpha/2)|} \gamma \frac{x_1 -x_2}{|x_1 -x_2|^{1+\alpha}} \Big| \\
  &= \eps^{1+\alpha}\ \frac{4}{\alpha} \frac{2^{\alpha} \Gamma((1+\alpha)/2)}{\sqrt{\pi} |\Gamma(-\alpha/2)|}\ \gamma\ \frac{1}{|x_1 -x_2|^{\alpha}} \,.
 \end{split}
\end{equation}
where $\gamma := (\int_\R (v'(x))^2 \txtd x)^{-1}$ and $v$ is a basic layer solution of~\eqref{eq:BLS:1} satisfying $\lim_{x\pm \infty} v(x) =\pm 1$ and $v(0) =0$.

In order to cross-validate this analysis, we measured the speed numerically. Renormalizing according to our previous considerations, we compute a renormalized speed
\begin{align*}
  \widehat{s}(\eps,\alpha)
  &:=\frac{ C_{\alpha} |x_1-x_2|^{\alpha}}{\eps^{1+\alpha} } s(\eps,\alpha).
\end{align*}
with the normalization factor
\[
 C_\alpha
  :=\Bigg( \frac{4}{\alpha} \frac{2^{\alpha} \Gamma((1+\alpha)/2)}{\sqrt{\pi} |\Gamma(-\alpha/2)|}\Bigg)^{-1}
  =\frac{\Gamma(1-\alpha)}{2\sec(\alpha \frac{\pi}{2})} \,.
\]
If the relationship~\eqref{speed_R} holds true, we expect that the renormalized speed then behaves like the factor~$\gamma$.
In Figure~\ref{fig:normalized_speed}, we have plotted the relationship $\alpha \mapsto \widehat{s}(\eps,\alpha)$ for different values of~$\eps$.
We first note, that all the curves (roughly) correspond to each other which hints that the
$\eps$--dependence is properly captured by the asymptotics.
Comparing the resulting curve to the numerically computed value of~$\gamma$ then
confirms the overall relationship.

\begin{figure}
  \centering
  \includeTikzOrEps{figures/new_hp_bdf2_fixed_s_speed}
  \caption{Comparison of the renormalized speed~$\widehat{s}(\eps,\alpha)$ to $\gamma$.
  Blue: different lines represent different values of $\varepsilon \in (0.002,0.05)$.}
  \label{fig:normalized_speed}
\end{figure}

\subsection{Time--to--collapse $t_{col}$}
The previous computation of the interface speed still relied on locating the interface and measuring its speed using a number of sample points. This might cause some hard to account for inaccuracies, resulting in small discrepancies to the predicted behavior and the slightly noisy look of Figure~\ref{fig:normalized_speed}.
In order to arrive at a second, more robust result, we also looked at the
time it takes for the interfaces to collide and finally annihilate, what we call the \textit{time--to--collapse}~$t_{col}$ for solutions of PDE~\eqref{eq:AC}.
We compute it by comparing our numerical solution to the constant function
taking the value $1$ everywhere.
If the difference, as measured in the $L^2$-norm, drops below $10^{-6}$ we record the current time and plot the analog to Figure~\ref{fig:normalized_speed}.

The time--to--collapse~$t_{col}$ for solutions of the Allen-Cahn equation~\eqref{eq:AC} with well--prepared initial data~\eqref{eq:TL:simple} is approximately/greater than the time--to--collision~$T_C$ in the associated ODE system~\eqref{ODE:centers} (up to some additive terms, see Theorems~\ref{thm:TC} and~\ref{thm:Teps}):
\begin{align}
  t_{col}(\eps,\alpha)
  \geq \bigg(\frac{d_0}{\eps}\bigg)^{1+\alpha} \frac{C_\alpha}{2 (1+\alpha) \gamma}.
\end{align}
Computing the modified time--to--collapse as
\begin{equation*}
  \widehat{t}_{col}
  :=\bigg(\frac{\eps}{d_0}\bigg)^{1+\alpha} \frac{2 (1+\alpha)}{C_\alpha} t_{col} \,,
\end{equation*}
we expect that $\widehat{t}_{col} \geq \gamma^{-1}$ but to be of comparable size.
Figure~\ref{fig:normalized_time} confirms the estimate on~$\widehat{t}_{col}$.

\begin{figure}
  \centering
  \includeTikzOrEps{figures/new_hp_bdf2_fixed_s_time}
  \caption{ Comparison of the renormalized time--to--collapse~$\widehat{t}_{col}$ to $1/\gamma$.
    Blue: different lines represent different values of $\varepsilon  \in (0.002,0.05)$.}
  \label{fig:normalized_time}
\end{figure}

\subsection{Interface width $w=w(\eps,\alpha)$}
Next, we consider the width of individual interfaces.
For each fixed $\alpha$, we can plot the interface width as a function of $\eps$ in a log-log plot as shown in Figure~\ref{fig:interface_width}.

\begin{figure}
  \center
  \begin{subfigure}[b]{0.46\textwidth}
    \includeTikzOrEps{figures/hp_bdf2_l10_new_initial_condition_width}
  \end{subfigure}\quad
  \begin{subfigure}[b]{0.46\textwidth}
    \includeTikzOrEps{figures/hp_bdf2_l10_new_initial_condition_width_fit}
  \end{subfigure}

  \caption{Behavior of the interface width (hp-BDF2 method). Different lines represent different values of $\alpha$.
  Blue line: actual data; black dashed line: power law fit $w(\eps,\alpha)\approx \eps^{a(\alpha)} b(\alpha)$.}
  \label{fig:interface_width}
\end{figure}

We see that the behavior with respect to $\eps$ can be approximated well by a power law of the form $w(\eps,\alpha) \approx \eps^{a(\alpha)} b(\alpha)$. In the second subplot of Figure~\ref{fig:interface_width}, we plotted the behavior of the coefficient $a(\alpha)$ when varying $\alpha$ in $(0.2,2)$. Namely, we chose a uniform grid of 40 points between $0.2$ and $1.9$. We see that the behavior especially in the regime of larger $\alpha$, can be well modeled by $a(\alpha) \approx \kappa_1 \alpha^{-1}+\kappa_2$.
Empirically we determined the constants in the fit as $\kappa_1:=-0.168298$ and $\kappa_2:= 1.11709$.
Overall we see the empirical law
\begin{align*}
  w(\eps,\alpha)\approx b(\alpha) \eps^{-0.168298 \alpha^{-1} +  1.11709},
\end{align*}
which is included in Figure~\ref{fig:interface_width} as the dotted black line. In order to get an estimate on $b(\alpha)$, we again rescale $w$ to get
\begin{align}
  \label{eq:definition_what}
 \widehat{w}(\eps,\alpha):=w(\eps,\alpha) \eps^{0.168298 \alpha^{-1} - 1.11709} .
\end{align}
We plot the behavior of $\widehat{w}$ in Figure~\ref{fig:modified_width}.
We observe that, while the $\eps$-dependence is well captured for larger values
of $\alpha$, the behavior for small $\alpha$ is more erratic.

\begin{figure}
  \begin{subfigure}{0.5\textwidth}
    \includeTikzOrEps{figures/new_hp_bdf2_fixed_s_width_nocorrection}
    \caption{Behavior of $\widehat{w}(\eps,\alpha)$ in~\eqref{eq:definition_what}.}
    \label{fig:modified_width_a}
  \end{subfigure}

  \begin{subfigure}{0.5\textwidth}
    \includeTikzOrEps{figures/new_hp_bdf2_fixed_s_width}
    \caption{Same  as in~(a), but also correcting using $C_{\alpha}$.}
  \end{subfigure}
  \begin{subfigure}{0.5\textwidth}
    \includeTikzOrEps{figures/new_hp_bdf2_fixed_s_width_alternatecorrection}
    \caption{Same  as in~(a), but also correcting using $C_{\alpha}$ (alternative version).}
  \end{subfigure}

  \caption{Behavior of the modified interface width $\widehat{w}(\varepsilon,\alpha)$ with respect to $\alpha$.
    Different lines represent different values of $\varepsilon$.}
\label{fig:modified_width}
\end{figure}

\section{Summary and Outlook}
\label{sec:outlook}

Following the analysis of a fractional Allen-Cahn equation on~$\R$, we obtained formulas for the interface speed and time--to--collision and its dependence on $\eps$ and $\alpha$. In Section~\ref{sec:results} we checked the validity of these formulas for our model on a bounded interval with homogeneous Neumann boundary conditions. Moreover, we studied numerically the behavior of the interface width well before the two interfaces collide. Following the analysis by Patrizi and Valdinoci, we considered \textit{well--prepared initial data}, i.e., we started close to a metastable profile. Hence, we are now in the position to provide a good dynamical unterstanding of the fractional Allen-Cahn equation in its metastable regime. Yet, for the classical Allen-Cahn equation, it is well-known that there are four possible stages in the generation, propagation and annihilation of metastable patterns~\cite{Chen:2004}: Starting from a large class of oscillating initial data, there is phase separation, then the formation of a metastable pattern, next a slow motion of the metastable pattern and (a possible) annihilation of two interfaces. The latter stages may repeat until either a constant stable state or a single basic layer solution remains. Thus, we are still missing a more detailed analysis of the early stages of the generation of metastable patterns for the fractional Allen-Cahn equation, which remains a topic for future work.

Another natural question is to try to extend the approach in~\cite{CarrPego} to the case of a fractional Allen-Cahn equation, i.e., to re-consider the auxiliary BVP~\eqref{eq:BVP} as a problem in its own right. We may ask what happens if we replace the second-order derivative by the fractional Laplacian? Is this problem still somehow connected to metastability? Is the solution unique if we require positivity in the interior? If so, what is its shape? If it has a peak near zero, what is the scaling of the spike layer as a function of $\ell$, $\eps$ and $\alpha$? These questions all remain to be dealt with in future work.

Another important direction is to consider other potentials $F$. In this work we considered a balanced potential, i.e. $F(-1)=F(1)$. To fully characterize the dynamics of the fractional Allen-Cahn equation for general unbalanced potentials is an open question, but see~\cite[\S 5]{NecNepGol:2008} and \cite{AchKue:2015, CabSir:2015} as examples for encouraging results towards this goal.

\appendix

\section{Dependence on the size of the domain}

Since all our analytical insights are dependent on working on the full space $\R$, whereas our numerics is restricted to a finite interval, it is natural to investigate the dependence of our quantities of interest on the size of the considered interval $(-L,L)$.
As can be seen in Figure~\ref{fig:dependence_on_domain}, neither the (modified) speed nor the width of the interface seems to depend strongly on the size of the domain. Especially for larger domains, the results become indistinguishable.

\subsection{Basic layer solution on $\R$}
\label{subsec:layer_solution_on_R}

The layer solutions of~\eqref{eq:BLS:1}--\eqref{eq:BLS:2} are -- in fact -- the stationary traveling wave solutions of~\eqref{eq:AC}, which are unique up to translations (see~\cite{CabSir:2015}) and locally asymptotically stable, see~\cite{AchKue:2015, MaNiuWang:2019}.
Therefore we aim to approximate the layer solution as the stationary solution of~\eqref{eq:AC}.
We use again the numerical scheme on intervals $[-L,L]$ as presented in Section~\ref{sec:numerics}.
Starting our numerical scheme in Section~\ref{sec:numerics} with $\eps=1$ and an initial datum~$u_0$ given as the linear function connecting $-1$ and $1$, we let the simulation run until a stationary state is reached. Finally we compute
\[ \int_{-\infty}^{\infty}{|v'|^2 \txtd x}
 = \int_{-\infty}^{-L}{|v'|^2 \txtd x} +\int_{-L}^{L}{|v'|^2 \txtd x} +\int_{L}^{\infty}{|v'|^2 \txtd x}\,.
\]
The term $\int_{-L}^{L}{|v'|^2 \txtd x}$ is computed from our numerical solution.
Whereas for the first and third summand we use that the layer solution is an odd function (see \cite[Thm. 2.4]{CabSir:2015} and its tail behaves as
\[
 v(x) \approx 1 +p x^{-\alpha} \quad\text{for } x\gg L \quad\text{where } p = -\frac{\sec(\alpha\pi/2)}{2\Gamma(1-\alpha)}\,,
\]
see \cite[Eq. (23)]{NecNepGol:2008} and \cite[Thm. 2.7]{CabSir:2015}.
Thus,
\begin{align*}
 \int_{-\infty}^{-L}{|v'|^2 \txtd x}
  &=\int_{L}^{\infty}{|v'|^2 \txtd x} \\
  &=\int_{L}^{\infty}{p^2 \alpha^2 x^{-2\alpha -2} \txtd x} \\
  &=p^2 \alpha^2 \frac{-1}{2\alpha +1} x^{-2\alpha -1} \Big|_{L}^{\infty} \\
  &=p^2 \alpha^2 \frac{1}{2\alpha +1} L^{-2\alpha -1} \,.
\end{align*}
Summing up, we deduce
\[ \int_{-\infty}^{\infty}{|v'|^2 \txtd x}
 = \int_{-L}^{L}{|v'|^2 \txtd x} +2 p^2 \alpha^2 \frac{1}{2\alpha +1} L^{-2\alpha -1}\,.
\]
We again chose different values of $L$ in order to investigate the influence of the choice of bounded domain, but observe that for larger intervals, the influence of the cutoff becomes negligible.

\begin{figure}
  \begin{subfigure}{0.5\linewidth}
    \includeTikzOrEps{figures/layer_solution_norms}
  \end{subfigure}
  \begin{subfigure}{0.5\linewidth}
    \includeTikzOrEps{figures/layer_solution_norms_full}
  \end{subfigure}
  \caption{$H^1$-seminorm $\int_{-\infty}^{\infty}{|v'|^2 \txtd x}$ of a basic layer solution~$v$ of~\eqref{eq:BLS:1}--\eqref{eq:BLS:2}.}
\end{figure}

\textbf{Acknowledgments:} The authors (FA, JMM, AR) gratefully acknowledge financial support by the Austrian Science Fund (FWF) through
the research program ``Taming complexity in partial differential systems'' (grant SFB F65).
AR has in addition been funded through project P29197-N32. CK gratefully
acknowledges support via a Lichtenberg Professorship of the
VolkswagenStiftung as well as partial support of the SFB/TR 109
Discretization in Geometry and Dynamics funded by the German Science
Foundation (DFG).

\bibliographystyle{plain}
\bibliography{references}

config: arxiv


TITLE: In which ratio does the point $P$ divide the segment $\overline{AN}$?
QUESTION [0 upvotes]: In an arbitrary triangle $\triangle ABC$, let $M\in\overline{AC}$ s. t. $|AM|:|MC|=2:1$ and let $N\in\overline{BC}$ s. t. $|BN|:|NC|=1:2$. Let $P$ be the intersection point of the segments $\overline{AN}$ and $\overline{BM}$. In which ratio does the point $P$ divide the segment $\overline{AN}$?

My attempt:
I thought I could apply the intercept theorem to find the ratio in which the point $P$ divides the segment $\overline{BM}$ and then express $\overrightarrow{AP}$ as a linear combination of $\overrightarrow{BM}$ and some vector in $\triangle ABC$ linearly independent of $\overrightarrow{BM}$.
Let $S\in\overline{NC}$ s. t. $\overline{AN}\parallel\overline{MS}$. From the given ratios, it follows:
$|AM|=2\lambda,\ |MC|=\lambda,\ |BN|=\mu, |NC|=2\mu, \ \lambda,\mu\in\Bbb Q$.
By the intercept theorem,
$$\begin{aligned}&|SC|:|NS|=|MC|:|AM|=1:2\\\implies&|SC|=\nu,\ |NS|=2\nu,\ \nu\in\Bbb Q\\\implies&|NC|=|NS|+|SC|=3\nu=2\mu\implies\mu=\frac32\nu\end{aligned}$$
Then
$$\begin{aligned}&|BP|:|PM|=|BN|:|NS|=\frac{\mu}{2\nu}=\frac{\frac32\nu}{2\nu}=\frac34\\\implies&\overrightarrow{PM}=\frac47\overrightarrow{BP}=\frac47\left(\frac13\overrightarrow{AC}-\overrightarrow{BC}\right)\end{aligned}$$
but it doesn't seem I've accomplished anything by finding $\frac{|BP|}{|PM|}.$
It would be perfect if I could find $\frac{|AP|}{|PN|}$ the same way, but there isn't enough information to do that and compare that result with $\overrightarrow{AP}=\alpha\left(\overrightarrow{AC}-\overrightarrow{NC}\right),\ \alpha\in\Bbb Q$.

Another option was to consider a midpoint $T$ of the segment $\overline{NC}$, so $\overrightarrow{BT}=\frac23\overrightarrow{BC}$. Then
$$\overrightarrow{AP}=\alpha\left(\overrightarrow{AB}+\overrightarrow{AT}\right)$$
May I ask for advice on solving this task?
Thank you in advance!

REPLY [1 votes]: Let $\mathcal{H}_{X,k}$ denote a homothety with center at $X$ and a extension factor $k$.
We have a following fact (Theorem):
If $\mathcal{H}_{X,k_1}$ and $\mathcal{H}_{Y,k_2}$ are homotheties then their compostion $\mathcal{H}_{X,k_1}\circ \mathcal{H}_{Y,k_2}$ is again some homothety $\mathcal{H}_{S,k}$ with $k=k_1k_2$ (if $k\ne 1$) and it center $S$ lies on a line $XY$.
Since we have:
\begin{align}
\mathcal{H}_{M, {-2}}: & \;C \longmapsto A\\
\mathcal{H}_{B,{3}}: &\; N \longmapsto C\\
\end{align}
we see that $P$ is a center of homothety which takes $N$ to $A$ with ratio $-6$, so $AP:PN = 6:1$.


TITLE: Why are trace monoids cancellative?
QUESTION [3 upvotes]: Trace monoids are just partially commutative free monoids. In other words, this is the set of words where words which can be obtained by switching around certain pairs of letters are considered equivalent.
Formally, take a set of symbols $\Sigma$ and a relation (called the independence relation) $I \subseteq \Sigma \times \Sigma$ which is symmetric. Define a relation $\sim$ on the set of words $\Sigma^*$, so that $x \sim y$ if there exists $(a, b) \in I$ and $u_1, u_2 \in \Sigma^*$ such that $x = u_1 a b u_2$ and $y = u_1 b a u_2$. Let $\equiv$ be the transitive reflexive closure of $\sim$. We define the traces on $\Sigma$ induced by $I$ to be the quotient of $\Sigma^*$ under $\equiv$.
One can show that concatenation of $\Sigma^*$ is stable under $\equiv$, and thus $\Sigma^*/\equiv$ has a natural monoid structure.
I want to show that this monoid structure is cancellative. So, I want to show that $xy \equiv xz \implies y \equiv z$.
According to the discussion near (1.7) here, this is clear for $\sim$ and hence is clear for $\equiv$. I do not understand why the later part of their claim is straightforward.
Any help with this proof is appreciated.

REPLY [0 votes]: After some thought, I found an elementary proof of this myself.

Let us focus on left cancellativity. Right cancellativity should be symmetric.
It is enough to prove the following:
$$\forall a \in \Sigma, x, y \in \Sigma^*, a \cdot x \equiv a \cdot y \implies x \equiv y \qquad (1)$$
Left cancellativity can be proved from this lemma using a simple induction on the term on the left.

To prove (1), we will show the following:

Lemma (2) : As before, say $a \in \Sigma$, $y, x' \in \Sigma^*$ and that $ax' \equiv y$. Then, (a) $y$ can be decomposed (written) in the form $x_0 \cdot a \cdot x_1$ where $x_0, x_1 \in \Sigma^*$. So that, (b) $a$ does not occur in $x_0$ and (c) $x_0 \cdot x_1 \equiv x'$. Also, (d) $a$ commutes with every symbol in $x_0$.

(Note that (d) readily means that $a$ commutes with the word $x_0$ itself.)
Now, let us see why (2) implies (1). Say there is some word $a \cdot x = a \cdot y$ as in the hypothesis of (1). Using (2), we see that there is a way to decompose $a\cdot y$ into $x_0 \cdot a \cdot x_1$ following the conditions. Since $x_0$ cannot have $a$ in it, it must be empty, which means $x_1 = x$. By condition (c), we have $y \equiv x_0 \cdot x_1 = x_1 = x$, like we wanted.

Lemma (2) looks straight forward by visual inspection. Indeed, it can be proven directly by induction on the proof of $ax' \equiv y$.
I have, however found that choosing the correct induction principle here can be tricky. I recommend using the following induction principle.

Fix an element $x \in \Sigma^*$. Let $P$ be a property on elements of $\Sigma^*$. Now suppose conditions $(*)$ and $(**)$ below hold
$$\forall y, [x \sim y \implies P y] \qquad (*)$$
$$\forall y \; z, [x \equiv y \land P y \land y \sim z \implies P z] \qquad (**)$$
Then, the for every $y$ such that $x \equiv y$, it holds that $P y$.

Showing (**) in our case requires us to think of the situation where we have $x_0 \cdot a \cdot x_1 \sim z$. This involves a bit of case work, where we inspect where the transposed pair is located. This is slightly tedious, but can be done.

I have checked this proof with Coq, so I believe this should be correct. You may need to assume that the independence relation is irreflexive, or that for all $x, y \in \Sigma$, either $x = y$ or $x \neq y$.


TITLE: How much space do I need to sort my socks?
QUESTION [11 upvotes]: In my pile of finished laundry, there are $2n$ socks of $n$ types, comprising $n$ easily distinguishable pairs. I sort the socks into pairs by picking one sock at a time randomly from the pile and either (1) laying it in a row of distinct socks or (2) pairing it with a sock already in the row, if there is such a sock, and putting the pair away. What is the expected maximum number of socks in the row?
Let $N$ be the peak number of socks in the row. Clearly $1\leqslant N\leqslant n$. And obviously, if $n=1$, then $N=\mathrm E[N]=1$. For $n=2$, I get $\mathrm E[N]=\frac53$. As $n$ increases, the calculation gets more complicated, and I don't know whether $\mathrm E[N]$ has a reasonable explicit general formulation. If not, perhaps it still has a nice asymptotic expression for large $n$. My guess is that there are positive constants $\alpha$ and $\beta$ such that $\lim_{\,n\rightarrow\infty}\mathrm E[N]/n^\beta=\alpha$.

REPLY [7 votes]: Expanding my comments into an answer per request.
Let $N_n$ be the random variable for the peak number of socks for a given $n$. It is known that
$$\overline{N}_n := \mathbf{E}[N_n] = 1, \frac{5}{3}, \frac{7}{3}, \frac{311}{105}, \frac{3377}{945}, \frac{3943}{945}, \ldots
$$
or in general $\displaystyle\overline{N}_n = \frac{a_n}{b_n}$ where
$a_n$ = OEIS A225177 and $b_n = (2n-1)!! = \prod\limits_{k=1}^n (2k-1)$.
As of this moment, closed form of the sequence is an unknown.
For more infos, one can follow the reference in the OEIS A225177 or consult the M/J2 section of this article.
There are two things in above article I'll like to comment.

It described an algorithm to compute the probability distribution for $N_n$.
It mentioned as $n \to \infty$, the probability distribution of $N_n$ becomes normal with mean $\frac{n}{2}$ and variance $\frac{n}{4}$.

I have implemented an algorithm similar to what described in the article and performed some
numerical experiment of the probability distribution.
Let $\mathbf{Var}_n = \mathbf{E}\left[ \left(N_n - \overline{N}_n \right)^2 \right]$  be the variance for $N_n$. Following table summarizes what I get:
$$\begin{array}{|r|rr:rr:r|}
\hline
N_n & \overline{N}_n & \overline{N}_n - \frac{n}{2} &
\mathbf{Var}_n & \mathbf{Var}_n/n & \Delta_n(2)\\
\hline
 200 &   105.107856 &     5.107856 &    40.733642 &    20.366821\% &     2.175188\%\\
 400 &   206.615246 &     6.615246 &    84.774930 &    21.193732\% &     1.600965\%\\
 600 &   307.674159 &     7.674159 &   129.707716 &    21.617953\% &     1.291395\%\\
 800 &   408.517714 &     8.517714 &   175.150639 &    21.893830\% &     1.124633\%\\
1000 &   509.230343 &     9.230343 &   220.941243 &    22.094124\% &     1.009333\%\\
1200 &   609.853430 &     9.853430 &   266.990895 &    22.249241\% &     0.923493\%\\
1400 &   710.410713 &    10.410713 &   313.244447 &    22.374603\% &     0.856246\%\\
1600 &   810.917224 &    10.917224 &   359.664606 &    22.479038\% &     0.802341\%\\
1800 &   911.383157 &    11.383157 &   406.224651 &    22.568036\% &     0.757826\%\\
2000 &  1011.815776 &    11.815776 &   452.904598 &    22.645230\% &     0.719742\%\\
\hline
\end{array}$$

As one can see, as $n$ becomes large, the leading behavior of $\overline{N}_n$ is indeed $\frac{n}{2}$.
The difference between $\overline{N}_n$ and $\frac{n}{2}$ is a very slowly
growing function. The difference seems to grow like $O(n^\beta)$ with $\beta \sim 0.37$.
The variance $\mathbf{Var}_n$ does seem scale as $n$. However, the actual coefficient in front of $n$ is smaller than $\frac14$. It is not clear whether this is a finite size effect or the limiting coefficient is indeed different from $\frac14$.

About the claim the probability distribution of $N_n$ is normal, it actually works surprisingly well.
Let $\sigma_{n} = \sqrt{\mathbf{Var}_n}$ and $P_n(k) = \mathbf{Prob}[ N_n = k ]$ be the probability that $N_n$ equal to a particular value $k$. At the end is a picture comparing $P_n(k)$ versus a normal distribution for various $n$.

The $x$ axis is $\displaystyle\;x = \frac{k - \overline{N}_n}{\sigma_n}$.
The $y$ axis is the scaled probability $\sigma_n P_n(k)$.
The $\color{orange}{\verb/orange/}$ curve is the PDF $\displaystyle\;\rho(x) = \frac{e^{-x^2/2}}{\sqrt{2\pi}}$ for a normal distribution of unit variance.

As one can see, for all the cases tested, the scaled probability $P_n(k)$ nearly fall right on top of $\rho(x)$. To measure how close the computed probability differ from the PDF of a normal distribution. Let us define following quantity
$$\Delta_n(s) = \max\left\{ \left|\frac{\sigma_nP_n(k)}{\rho(x)} - 1\right| : |x| < s \right\}$$
This measures the biggest difference between the ratio $\frac{\text{scaled probability}}{\text{reference PDF}}$ and $1$ for $x$ within $s$ standard derivation.
As shown in last column of table above. As long as we are within 2 standard derivation,
the difference between the ratio and $1$ is at most a few percents and improves as $n$ increases.


TITLE: How to write a polynomial function that has the roots $-2$ and $\sqrt7$?
QUESTION [2 upvotes]: I need to write a polynomial function with integer coefficients that has the roots $-2$ and $\sqrt7$. I'm able to do this correctly when I'm given roots like $-3+i$ & $-3-i$, in which I set the roots equal to zero and then multiply them by one another. However, when I try this with $-2$ and $\sqrt7$ and multiply $x+2$ by $x-\sqrt7$, I get $x^2+2x-\sqrt7 x-2\sqrt7$. I don't know where to go from here, and I don't think that this is the correct next step. What do I do next?

REPLY [0 votes]: $x = \sqrt 7 \implies$
$x^2 = 7\implies$
$x^2 - 7 =0$.  So $x^2 -7$ has $\sqrt 7$ as a root.
And $x=-20\implies x+2 = 0$ so $x+2$ has $-2$ as a root.
If $x^2-7=0$ and $x+2=0$ then $(x^2-7)(x+2)=0$ and $(x^2-7)(x+2)=x^3+2x^2-7x-14$ has $\sqrt 7$ and $-2$ as roots.
(It also has $-\sqrt 7$ as a root.)
=====
By the way.  There is no 2nd degree polynomial with two roots $-2$ and $\sqrt 7$.
A second degree polynomial has at most two roots and if the two roots are $r_1$ and $r_2$ the polynomial is $(x-r_1)(x-r_2)= x^2 -(r_1+r_2)x + r_1r_2$ and if $r_1$ or $r_2$ are irrational we have no reason to assume $r_1+r_2$ or $r_1r_2$ are integers.
Is if $r_1,r_2 = -2,\sqrt 7$ they won't be.  But to counter $\sqrt 7$ as a root we can add a third root of $-\sqrt 7$ to get a third degree polynomial:
$(x-\sqrt 7)(x+\sqrt 7)(x+2) = x^3+2x^2-7x-14$.


TITLE: Justification for the nature of planet's orbit in gravitational field!
QUESTION [1 upvotes]: In kleppner Mechanics in the chapter central force he derived the polar form of orbit for gravitational force as illustrated below: (first two equations are derived from fundamentals of central force)

Often we are interested in the path of the particle, which means
knowing r as a function of θ rather than as a function of time.
We call r(θ) the orbit of the particle. (The term is used even if
the trajectory does not close on itself.) Dividing Eq. (9.14) by
Eq. (9.12) gives:

we  find  that  the  equation  for  the  orbit  Eq. (9.16)  becomes,using  indefinite  integrals,

Now my Questions:

why he takes dr/dt as positive(because all the terms in RHS in eq(9.12) is positive);? it is evidently can be negative and also from energy conservation it can be both positive and negative!
And as he has taken it(dr/dt) as positive that's why right side of eq**(9.16)** is also is positive! But for a conic section dθ/dr is not always positive! It is positive for some interval and negative for some interval(actually 50:50)! so eq(9.16) is only valid for a certain interval (according to the choice of my coordinate system); but he(the writter of the book) derived the general solution not only for a certain interval but the total interval! Where am I getting wrong or how this is happening?

According to me something is fishy in the integration of eq(9.16) : the value of the integral doesn not hold everywhere in the interval of θ.    But I am not able to understand the matter properly!

REPLY [0 votes]: Yes! it's not necessary that $dr/dt$ is positive, it can be negative or it might be positive for some interval of time, and then again positive. This can be seen from the trajectory of comets themselves. Any comet which is in a parabolic trajectory came to attractor (some planet) with $dr/dt<0$ and then swing back making $dr/dt$ to be positive.
So you can take
$$\frac{dr}{dt}=-\sqrt{\cdots}$$
It makes no difference. Try by yourself!

To show this you might  wanna check it with your final result
$$\frac{dr}{dt}=-\frac{r_0}{(1-\epsilon\cos\theta)^2}(\epsilon \sin\theta)\frac{d\theta}{dt}$$
As you  can see  the sign of $dr/dt$ determine by the $\dot{\theta}$ and $\theta$. Note that sign of $\dot{\theta}$  depend only on initial  condition (look at the relation $\dot{\theta}=l/\mu r^2).$
You will also find taking $\theta_0=-\pi /2$ is crucial!


TITLE: general solution to linear second order difference equation
QUESTION [1 upvotes]: Is there a general solution to difference equations of the form:
$$
u(n+2) + a(n)u(n+1) + b(n)u(n) = 0
$$
Thank you in advance

REPLY [0 votes]: Besides linear ODEs (http://en.wikipedia.org/wiki/Abel%27s_identity#Proof_using_Liouville.27s_formula)(Does this ODE question have closed form solution?), linear difference equations also have similar tricks:
In fact $u(n+2)+a(n)u(n+1)+b(n)u(n)=0$ can rewrite to the matrix form
$\begin{pmatrix}u(n+1)\\u(n+2)\end{pmatrix}=\begin{pmatrix}0&1\\-b(n)&-a(n)\end{pmatrix}\begin{pmatrix}u(n)\\u(n+1)\end{pmatrix}$
$\begin{pmatrix}u(n)\\u(n+1)\end{pmatrix}=\begin{pmatrix}0&1\\-b(n-1)&-a(n-1)\end{pmatrix}\begin{pmatrix}u(n-1)\\u(n)\end{pmatrix}$
This becomes to a matrix difference equation.
For the discrete version $\begin{pmatrix}u_n\\u_{n+1}\end{pmatrix}=\begin{pmatrix}0&1\\-b(n-1)&-a(n-1)\end{pmatrix}\begin{pmatrix}u_{n-1}\\u_n\end{pmatrix}$ ,
$\begin{pmatrix}u_n\\u_{n+1}\end{pmatrix}=\prod_{k=n-1}^c\begin{pmatrix}0&1\\-b(k)&-a(k)\end{pmatrix}\begin{pmatrix}u_c\\u_{c+1}\end{pmatrix}$
$\therefore u_n=C_1f(n)+C_2g(n)$ , where $C_1$ and $C_2$ are arbitrary constants, $f(n)$ and $g(n)$ are the top-left and top-right elements of the simplified matrix expression $\prod_{k=n-1}^c\begin{pmatrix}0&1\\-b(k)&-a(k)\end{pmatrix}$ respectively.
Now for the continuous version $\begin{pmatrix}u(n)\\u(n+1)\end{pmatrix}=\begin{pmatrix}0&1\\-b(n-1)&-a(n-1)\end{pmatrix}\begin{pmatrix}u(n-1)\\u(n)\end{pmatrix}$
$\therefore u(n)=T_1(n)f(n)+T_2(n)g(n)$ , where $T_1(n)$ and $T_2(n)$ are arbitrary periodic functions with period $1$, $f(n)$ and $g(n)$ are the top-left and top-right elements of the simplified and continuous extended matrix expression $\prod_{k=n-1}^c\begin{pmatrix}0&1\\-b(k)&-a(k)\end{pmatrix}$ respectively.
But the simplification process and the continuous extension process of $\prod_{k=n-1}^c\begin{pmatrix}0&1\\-b(k)&-a(k)\end{pmatrix}$ are super complicated.


TITLE: Show that the map is a diffeomorphism
QUESTION [2 upvotes]: Given a set $M = \{(x,y,z)\in \mathbb R^3 | x^2+y^2-z^2=1\}$ and a function $\alpha \colon M \to S^1 \times \mathbb R$,
$$\alpha(x, y, z) = \left(\frac{(x,y)}{\sqrt{1+z^2}},z\right)$$
Is there a smart method of showing that $\alpha$ is a diffeomorphism? I already showed that $M$ is an immersed submanifold of $\mathbb R^3$.
Is it enough to show that it is a bijection, the map is smooth, and the jacobian is nowhere 0?

REPLY [2 votes]: Writing $\alpha(x,y,z) = (u,v)$ with $u \in S^1 \subseteq \mathbb{R}^2$ and $v \in \mathbb{R}$, we have $\alpha^{-1}(u,v) = (\sqrt{1 + v^2}u,v)$. Now just verify that both $\alpha$ and $\alpha^{-1}$ are smooth.


TITLE: How to make correct system of equations to solve for the angles in this triangle?
QUESTION [0 upvotes]: I'm trying to solve this triangle for $X$. Thereby, I've tried to make correct system of equations. What would be the correct equations?
Here are the equations I can find

In $\triangle ABC$,  recalling that $\angle ACD = y$

$$48 + 24 + x  + 12 + y = 180$$

In $\triangle ADC$

$$84+x + y = 180$$
However, I'm getting the same equations.
Regards

REPLY [0 votes]: The key insight is to observe that $\angle(BDA)=\angle(BDC)=138$. From here in triangle $BDC$ we have $AB/\sin(138)=BD/\sin(24)$ and from triangle BCD we get
$BC/\sin(138)=BD/\sin(12)$.  Divingig both equations gives
$$BC=\frac{AB\sin(24)}{\sin(12)}.$$
Now use law of cosine in $ABC$ to get
$$AC^2=AB^2+BC^2-2AB\cdot BC\cos(48)
=BC^2\left(\frac{\sin^2(12)}{\sin^2(24)}+1
-2\frac{\sin(12)}{\sin(24)}{\cos(48)}\right).$$
Finally use law of sine in triangle $ABC$:
$$\frac{\sin(x+24)}{BC}=\frac{\sin(48)}{AC}$$
plug in $AB$ to get
$$\sin(x+48)=\frac{\sin(48)}{\sqrt{\dfrac{\sin^2(12)}{\sin^2(24)}+1
-2\dfrac{\sin(12)}{\sin(24)}{\cos(48)}}}$$