# On the quantization of AB phase in nonlinear systems

Xi Liu  NUS Graduate School - Integrative Sciences and Engineering Programme
(ISEP), National University of Singapore, Singapore 119077, Singapore  Qing-
hai Wang  Department of Physics, National University of Singapore, Singapore
117551, Singapore  Jiangbin Gong<EMAIL_ADDRESS>Department of Physics,
National University of Singapore, Singapore 117551, Singapore  Centre for
Quantum Technologies, National University of Singapore, Singapore 117543,
Singapore

###### Abstract

Self-intersecting energy band structures in momentum space can be induced by
nonlinearity at the mean-field level, with the so-called nonlinear Dirac cones
as one intriguing consequence. Using the Qi-Wu-Zhang model plus power law
nonlinearity, we systematically study in this paper the Aharonov-Bohm (AB)
phase associated with an adiabatic process in the momentum space, with two
adiabatic paths circling around one nonlinear Dirac cone. Interestingly, for
and only for Kerr nonlinearity, the AB phase experiences a jump of $\pi$ at
the critical nonlinearity at which the Dirac cone appears or disappears,
whereas for all other powers of nonlinearity the AB phase always changes
continuously with the nonlinear strength. Our results may be useful for
experimental measurement of power-law nonlinearity and shall motivate further
fundamental interest in aspects of geometric phase and adiabatic following in
nonlinear systems.

##  I Introduction

The dynamics depicted by a nonlinear discretized Schödinger equation (NDSE)
can be extremely rich, including the emergence of many-dimensional chaos,
solitons, and breathers etc. The problem can be much reduced by assuming the
translational invariance of a wave under consideration. With this assumption
the main physics is about the features of Bloch waves, the associated energy
bands, and how they respond to changes in the parameters of a nonlinear
system. Interestingly, the nonlinear Bloch bands of NDSE can induce gapless
band structures absent in linear systems, such as 2-dimensional (2D) noninear
Dirac cones [1] induced by Kerr nonlinearity [2]. Even more peculiar, such
nonlinear Dirac cones are formed by exotic nonlinear energy bands in a
subregime of the Brillioun zone [3, 1, 4, 5, 6, 7, 8].

As a close analog to a setting in real space to measure the Aharonov-Bohm (AB)
phase around a singularity point with magnetic flux, let us now imagine two
adiabatic paths, in the momentum space, circling around a band-crossing point.
If we adiabatically change the Bloch momentum, so as to guide the Bloch wave
to evolve along the two adiabatic paths, the final phase difference thus
generated between the two adiabatic paths is termed the nonlinear AB phase
[1]. One may naïvely think of the following: provided that the dynamical
phases between the two adiabatic paths are identical and hence have zero
contribution to the phase difference of interest, the obtained AB phase would
be just the Berry phase associated with the band degeneracy point. The actual
physics turns out to be more interesting than just a Berry phase. Because of
nonlinearity, any small deviation of the adiabatically following state from
the instantaneous Bloch wave causes a tiny correction to the dynamical phase,
and accumulation of such tiny corrections over the entire adiabatic protocol
yields an unfamiliar geometrical phase on top of the expected Berry phase.
Remarkably, as a possible means of topological charaterization of nonlinear
Dirac cones, it is found in Ref. [1] that the nonlinear AB phase around
nonlinear Dirac cones induced by Kerr nonlinearity added to the so-called Qi-
Wu-Zhang (QWZ) model [9] is quantized in $\pi$, whereas the Berry phase is not
quantized (thus in sharp contrast to a variery of linear systems, where the
Berry phase around a Dirac cone is quantized in $\pi$ [10, 11, 12, 13]).
Echoing with the finding in [1], Ref. [4] found $\pi$-quantization of a
nonlinear Zak phase and Ref. [14] further confirmed the $\pi$-quantization of
the nonlinear AB phase around a nodal line induced by Kerr nonlinearity.

The goal of this work is entirely focused on aspects of the nonlinear AB phase
around Dirac cones induced by general power law nonlinearity [15, 16, 17, 18,
19, 20, 21, 22, 23]. In this way, it becomes possible to answer whether the
previously obtained AB phase quantization is unique to Kerr nonlinearity and
if so, why there is such uniqueness. Using the QWZ model [9] as the linear
limit, we are able to analytically show that Kerr nonlinearity happens to be a
critical case among all kinds of power law nonlinearity. Specifically, for any
nonlinearity other than the cubic order, the $\pi$-quantization of nonlinear
AB phase does not exist. Our analytical results are further confirmed by
direct numerical simulations.

##  II Hamiltonian and energy spectrum

The momentum-space Hamiltonian is composed of a QWZ model with power law
nonlinearity characterized by a parameter $p$:

$\displaystyle\widehat{H}(\psi)$ $\displaystyle=J_{1}\sin
k_{1}\sigma_{1}+J_{2}\sin
k_{2}\sigma_{2}+\beta(k_{1},k_{2})\sigma_{3}+g\begin{bmatrix}|\psi_{1}|^{2p}&0\\\
0&|\psi_{2}|^{2p}\end{bmatrix},$ (1)

where $\sigma_{i}$ are Pauli matrices and $\psi_{a}$ are two components of the
wavefunction, $\psi=\begin{bmatrix}\psi_{1}\\\ \psi_{2}\end{bmatrix}$. The
normalization of the wavefunction means that
$|\psi_{1}|^{2}+|\psi_{2}|^{2}=1$. The nonlinearity parameter $p$ is a
nonnegative real number. The Kerr nonlinearity corresponds to $p=1$. The
parameters $k_{1}$ and $k_{2}$ are two quasimomenta, whose values will be
adiabatically tuned in order to implement an actual adiabatic protocol to
generate the nonlinear AB phase.

To solve the nonliner eigenvalue problem,

$\displaystyle\widehat{H}(\psi)\ket{\psi}$ $\displaystyle=E\ket{\psi},$ (2)

we introduce a real parameter $x$ as

$\displaystyle\psi_{1}$
$\displaystyle=\sqrt{\frac{1+x}{2}},\qquad\psi_{2}=\sqrt{\frac{1-x}{2}}e^{i\varphi}.$
(3)

We will see later that the angular variable $\varphi$ is the same as in the
Fig. 2. It turns out that $x$ is the central quantity for expressing energy,
dynamical phase, Berry phase, and nonlinear AB phase. It can be shown that the
instantaneous eigenenergy is

$\displaystyle E$
$\displaystyle=\frac{\beta}{x}+\frac{g}{x}\left[\left(\frac{1+x}{2}\right)^{p+1}-\left(\frac{1-x}{2}\right)^{p+1}\right],$
(4)

where $x$ satisfies the following algebraic equation,

$\displaystyle\frac{1-x^{2}}{x^{2}}\left\\{\beta+\frac{g}{2}\left[\left(\frac{1+x}{2}\right)^{p}-\left(\frac{1-x}{2}\right)^{p}\right]\right\\}^{2}$
$\displaystyle=|\gamma|^{2},$ (5)

with $\gamma:=J_{1}\sin k_{1}-iJ_{2}\sin k_{2}$.

In order to have a Dirac point in the energy spectrum, the energy must be
doubly degenerate at $k_{1}=k_{2}=0$. Since $\gamma=0$ at the this point, $x$
must satisfy

$\displaystyle\beta(0,0)+\frac{g}{2}\left[\left(\frac{1+x}{2}\right)^{p}-\left(\frac{1-x}{2}\right)^{p}\right]$
$\displaystyle=0.$ (6)

For simplicity, we choose

$\displaystyle J_{1}$ $\displaystyle=J_{2}:=B,$ (7)
$\displaystyle\beta(k_{1},k_{2})$ $\displaystyle=B(-1+\cos k_{1}+\cos k_{2}).$
(8)

Hence $\beta(0,0)=B$. It is clear that the nonlinearity strength $g$ and
energy $E$ can be scaled in terms of $B$. Energy spectra with $p=1,1.5,2$ and
$g=2.5B$ are shown in Fig. 1, where Dirac cone is clearly visible around the
origin. A perturbative analysis of energy spectrum near the Dirac cone can be
found in Appendix A.

Figure 1: Nonlinear band structure for small momenta in the vicinity of the
origin, i.e., for small values of $|k_{1}|$ and $|k_{2}|$, with the power law
nonlinearity parameter $p=1.5,2,2.5$ and the nonlinear strength parameter
$g=2.5B$. See the main text for details of the system parameters. The Dirac
cone emerges from the lower energy band.

##  III Dynamics of Adiabatic Following

To obtain the nonlinear AB phase, let us consider two adiabatic paths along a
small circle around the origin $k_{1}=k_{2}=0$. As shown in Fig. 2, starting
at the same point S, along each path the system is guided to move along one
half of the perimeter of the circle using the same amount of time. The two
adiabatic paths are “recombined" at the end of the evolution at point N. As
introduced in Sec. I, the phase difference acquired by the system between two
adiabatic paths is called the nonlinear AB phase. Clearly, the nonlinear AB
phase here is the sum of the dynamical phase difference and the Berry phase
associated with the closed loop around the band-degeneracy point. We shall
study below the possible AB phase quantization for a varying nonlinearity
strength $g$ and for different nonlinear parameters $p$. The quasimomenta
$k_{1}$ and $k_{2}$ associated with two spatial dimensions are parameterized
by $\varphi$ and will be made to adiabatically change.

At the starting point S, the system is assumed to be prepared in the Bloch
eigenstate at momentum space location S. As the system adiabatically evolves
along the path SEN or SWN, the time-evolving state deviates from the
instantaneous eigenstate along the path, with the tiny deviation at the order
of the adiabatic parameter $\varepsilon$. The slower the rate of adiabatic
change is, the smaller $\varepsilon$ is, and the less the deviation. Here
nonlinearity plays a key role. That is, the dynamical phase also obtains a
correction at the order of $\varepsilon$. Since the total evolution time is of
order $O(\varepsilon^{-1})$, the $O(\varepsilon)$ term in this phase
correction will contribute an $\varepsilon$-independent term through
accumulation, yielding a geometric phase term out of the dynamical phase. This
will not occur in linear terms because such correction accumulated over the
entire adiabatic process is at most of the order of $\varepsilon$, which
vanishes for sufficiently slow adiabatic protocols.

Figure 2: Dynamical paths in the momentum space. Path SWN and path SEN are
symmetric halves of the perimeter of the circle. The system starts its
adiabatic following at point S, and ends at point N. The two paths are
parameterized by $\varphi$ in the main text.

The dynamics of the states is governed by the time-dependent Schrödinger
equation,

$\displaystyle i\ket{\dot{\Psi}}$ $\displaystyle=\widehat{H}(\Psi)\ket{\Psi},$
(9)

where the Hamiltonian is given by Eq. (1) with $\psi$ being replaced by
$\Psi$. Here the overhead dot denotes the time derivative. We will solve this
equation up to the order of $\varepsilon$ as described above. Through the
lengthy computation as illustrated in Appendix B, we obtain the instantaneous
change rate of the overall phase of a time-evolving state as

$\displaystyle\dot{\theta}$
$\displaystyle\sim-E-\frac{1-x}{2}\dot{\varphi}+gp\frac{x(1-x^{2})}{4\Delta}\left[\left(\frac{1+x}{2}\right)^{p}-\left(\frac{1-x}{2}\right)^{p}\right]\dot{\varphi},$
(10)

with

$\displaystyle\Delta$
$\displaystyle:=\beta+\frac{g}{2}\left[(1-px+px^{2})\left(\frac{1+x}{2}\right)^{p}-(1+px+px^{2})\left(\frac{1-x}{2}\right)^{p}\right].$
(11)

We recognize that the circular integration of the second term in Eq. (10) is
nothing but the Berry phase $\theta_{B}$, because it assumes the same form as
in the linear limit. The rest of the phase is from the dynamical phase
$\theta_{D}$, which contains two parts: the first part comes from the
instantaneous eiegnenergy $E$ and the second part from the third term in Eq.
(10) as a new contribution from the nonlinearity. Specifically,

$\displaystyle\theta_{B}$ $\displaystyle:=-\oint\frac{1-x}{2}d\varphi,$ (12)
$\displaystyle\theta_{D}$ $\displaystyle:=-\int
Edt+gp\int\frac{x(1-x^{2})}{4\Delta}\left[\left(\frac{1+x}{2}\right)^{p}-\left(\frac{1-x}{2}\right)^{p}\right]d\varphi.$
(13)

In the event that the Dirac cone does exist at the point $k_{1}=k_{2}=0$, the
obtained phase difference between the two adiabatic paths described in Fig. 2
then becomes the nonlinear AB phase $\theta_{AB}$. Since the two adiabatic
paths are symmetric by construction and that they take the same amount of
time, the leading term in Eq. (13) contributes the same in each of the two
paths. Thus, the difference of the dynamical phases between two paths comes
from the second term of Eq. (13) only. Thus, the total nonlinear AB phase is

$\displaystyle\theta_{AB}$ $\displaystyle:=\theta_{B}+\delta\theta_{D}$
$\displaystyle\sim-\pi(1-x)+\pi
gp\frac{x(1-x^{2})}{2\Delta}\left[\left(\frac{1+x}{2}\right)^{p}-\left(\frac{1-x}{2}\right)^{p}\right].$
(14)

Note that we take into account that the paths are chosen to be close to the
Dirac cone (so that the cones indeed have linear dispersion relations),
namely, $|k_{1}|$ and $|k_{2}|$ are small at all times. The leading behavior
of the dynamical phase difference term is then found to be

$\displaystyle\delta\theta_{D}$ $\displaystyle\sim\pi
gp\frac{x_{0}(1-x_{0}^{2})}{2\Delta_{0}}\left[\left(\frac{1+x_{0}}{2}\right)^{p}-\left(\frac{1-x_{0}}{2}\right)^{p}\right],$
(15)

where $\Delta_{0}$ is $\Delta$ evaluated at $x=x_{0}$ and $k_{1}=k_{2}=0$,
$x_{0}$ is the solution of Eq. (6), and

$\displaystyle\Delta_{0}$
$\displaystyle=-\frac{gpx_{0}}{2}\left[(1-x_{0})\left(\frac{1+x_{0}}{2}\right)^{p}+(1+x_{0})\left(\frac{1-x_{0}}{2}\right)^{p}\right].$
(16)

For the Berry phase, the leading behavior is

$\displaystyle\theta_{B}$ $\displaystyle\sim-\pi(1-x)\sim-\pi(1-x_{0}).$ (17)

Figure 3: Dynamical, Berry and the nonlinear AB phase plotted against
nonlinearity strength $g$ for different values of power law nonlinearity
parameter $p$. The solid lines are theoretical values, and the dots are
numerical verification. Only for Kerr nonlinearity $p=1$, the AB phase has a
quantized jump of $\pi$ at the critical value $g=\pm 2B$ and stays at $\pi$
for $|g|>2B$.

As detailed in Appendix A, For $|g|>2B$, a nonlinear Dirac cone is located at
the origin. For $|g|<2B$, the only possible solutions to Eq. (5) are $x=\pm 1$
and there is no Dirac cone. For $g\in(0,2B)$, we can hence assign $x_{0}=-1$,
and for $g\in(-2B,0)$, we may assign $x_{0}=1$. With this convention, it is
clear to see that $\theta_{B}$ is constantly $0$ (mod $2\pi$) for
$g\in(-2B,2B)$. The Berry phase $\theta_{B}$ becomes nonzero and changes
continuously for $|g|>2B$. For each $p$, as we continuously tune $g$, $x_{0}$
can be easily solved numerically using Eq. (6), thus obtaining the theoretical
values of the leading terms of the dynamical phase, Berry phase and AB phase
around the origin. We also numerically solve the evolution using the
Schrödinger equation Eq. (9) along the two paths, and compute the dynamical
phase, AB phase and Berry phase using numerical solutions of the evolution.
The evolution is computed using an operator-splitting algorithm. The results
are presented in Fig. 3.

In each plot, solid lines are theoretical values, while dots on the solid
lines are computed from numerical evolutions. In Fig. 3(a), for any
$p=0.5,1,1.5,2,2.5,3$, the dynamical phase around the origin is $0$ for
$g\in(-2B,2B)$. At the critical value $g=\pm 2B$ where the Dirac cone appears,
for $p=0.5$, the Dirac cone changes continuously with respect to $g$. For
$p=1$, there is a quantized jump of $\pm\pi$ at $g=\pm 2B$. For
$p=1.5,2,2.5,3$, there is a quantized jump of $\pm 2\pi$ at the critical value
$g=\pm 2B$ (so this is equivalent to no change). In Fig. 3(b), the Berry phase
(modulo $2\pi$) is identically $0$ for $g\in(-2B,2B)$, and changes
continuously with respect to $g$. In Fig 3(c), the AB phase (modulo $2\pi$) is
the sum of the dynamical phase in Fig. 3(a) and the Berry phase in Fig. 3(b).
Only for $p=1$, the AB phase has a quantized jump of $\pi$ at the critical
value $g=\pm 2B$ and stays at $\pi$ for $|g|>2B$, as discovered by Ref. [1].
For all other values of $p$, the AB phase changes continuously with respect to
$g$. The special behavior of $p=1$ is because of the fact that $p=1$ is a
critical value for the limit $\lim_{g\to\pm 2B^{\pm}}\delta\theta_{D}$, as
will be explained in the next section.

##  IV Mechanism of the jump of AB phase at $g=\pm 2B$ for Kerr nonlinearity

For $p>1$, we can factor out a factor $(1-x_{0}^{2})$ from $\Delta_{0}$ which
cancels the same factor in the numerator of $\delta\theta_{D}$,

$\displaystyle\delta\theta_{D}(p>1)$
$\displaystyle\sim-\pi\frac{(1+x_{0})^{p}-(1-x_{0})^{p}}{(1+x_{0})^{p-1}+(1-x_{0})^{p-1}},$
(18)

which equals $\mp 2\pi$ or equivalently zero since $x_{0}=\pm 1$, for $|g|=2B$
or when the Dirac cone starts to appear.

likewise, for $p=1$, we have

$\displaystyle\delta\theta_{D}(p=1)$ $\displaystyle\sim-x_{0}\pi,$ (19)

which equals $\mp\pi$ since $x_{0}=\pm 1$, for $|g|=2B$.

Finally, for $0<p<1$,

$\displaystyle\delta\theta_{D}(0<p<1)$
$\displaystyle\sim-\pi(1-x_{0}^{2})^{1-p}\frac{(1+x_{0})^{p}-(1-x_{0})^{p}}{(1+x_{0})^{1-p}+(1-x_{0})^{1-p}},$
(20)

which vanishes for $x_{0}=\pm 1$, for $|g|=2B$.

Calculations above make it clear that the nonlinear AB phase associated with
Kerr nonlinearity ($p=1$) is most special as the extra nonlinearity-induced
correction to dynamical phase experiences a $\pi$ jump when the Dirac cone
appears. What is intriguing for Kerr nonlinearity is that the nonlinear AB
phase stays quantized at $\pi$ for $|g|>2B$, as $\theta_{B}$ and
$\delta\theta_{D}$ happen to be complementary to each other, as shown in Eqs.
(17) and (19). For all other forms of power-law nonlinearity, there is no such
jump, $\pi$-quantization is thus absent, and consequently, the nonlinear AB
phase only changes continuously with respect to $g$. This finally explains why
in Fig. 3 only the nonlinear AB phase for Kerr nonlinearity ($p=1$) displays a
quantization plateau for $|g|>2B$.

##  V Conclusion

In this paper, we analytically and computationally examined the so-called
nonlinear AB phase around Dirac cones induced by power-law nonlinearity added
to the QWZ model often used for studies of topological band structures. With
our analytical results, we are able to explain why the nonlinear AB phase has
a quantized jump of $\pi$ when Dirac cone starts to appear or disappear, for
and only for Kerr nonlinearity. In the context of nonlinear AB phase that can
be in principle measured in experiments, Kerr nonlinearity is thus identified
as a critical form of nonlinearity.

###### Acknowledgements.

J.G. is grateful to Prof. Giulio Casati for his many years of guidance,
interaction and collaboration. J.G. would also like to thank Prof. Chushun
Tian for very useful discussions. The computational work for this article was
fully performed on resources of the National Supercomputing Centre, Singapore
(https://www.nscc.sg).

##  Appendix A Eigenvalue problem

The instantaneous eigenenergy and eigenstate satisfy the Schrödinger equation,

$\displaystyle\begin{bmatrix}\beta(k_{1},k_{2})+g|\psi_{1}|^{2p}&\gamma(k_{1},k_{2})\\\
\gamma^{*}(k_{1},k_{2})&-\beta(k_{1},k_{2})+g|\psi_{2}|^{2p}\end{bmatrix}\begin{bmatrix}\psi_{1}\\\
\psi_{2}\end{bmatrix}=E\begin{bmatrix}\psi_{1}\\\ \psi_{2}\end{bmatrix}.$ (21)

In terms of the two components of the state, we have

$\displaystyle\gamma\psi_{2}$
$\displaystyle=(E-\beta-g|\psi_{1}|^{2p})\psi_{1},$ (22)
$\displaystyle\gamma^{*}\psi_{1}$
$\displaystyle=(E+\beta-g|\psi_{2}|^{2p})\psi_{2}.$ (23)

Plugging the expressions in Eq. (3), we see that

$\displaystyle\gamma e^{i\varphi}$
$\displaystyle=\left[E-\beta-g\left(\frac{1+x}{2}\right)^{p}\right]\sqrt{\frac{1+x}{1-x}}.$
(24)

Since the right hand side of the above equation is real, we recognize that the
phase variable $\varphi$ is simply the opposite of the phase of $\gamma$,

$\displaystyle\varphi=-\textrm{arg}(\gamma).$ (25)

Recall that $\gamma=B(\sin k_{1}-i\sin k_{2})$ in our choice, this means that
$\varphi$ is the same angle illustrated in Fig. 2 for sufficiently small
$|k_{1}|$ and $|k_{2}|$.

Multiplying $\psi_{1}^{*}$ on both sides of Eq. (22), multiplying
$\psi_{2}^{*}$ on Eq. (23) and taking complex conjugate, then subtracting the
two equations, one ontains

$\displaystyle\beta+g(|\psi_{1}|^{2p+2}-|\psi_{2}|^{2p+2})$
$\displaystyle=E(|\psi_{1}|^{2}-|\psi_{2}|^{2}).$ (26)

In terms of the parameter $x$ defined in Eq. (3), we get the instantaneous
eigenenergy as in Eq. (4).

One can then multiply the two equations in (22) and (23) together. Eliminating
the common factor $\psi_{1}\psi_{2}$, we arrive at

$\displaystyle|\gamma|^{2}$
$\displaystyle=E^{2}-gE(|\psi_{1}|^{2p}+|\psi_{2}|^{2p})-\beta^{2}-\beta
g(|\psi_{1}|^{2p}-|\psi_{2}|^{2p})+g^{2}|\psi_{1}|^{2p}|\psi_{2}|^{2p}.$ (27)

Further using Eq. (3), we obtain the equation satisfied by the variable $x$ in
Eq. (5).

Apparently, if the Dirac cone exists, at the Dirac point $k_{1}=k_{2}=0$, the
energy is doubly degenerate. As a result, $x$ is also doubly degenerate.
Namely, it must satisfy Eq. (6) with $\beta(0,0)=B$. That is

$\displaystyle 2^{p+1}B+g\left[(1+x)^{p}-(1-x)^{p}\right]$ $\displaystyle=0.$
(28)

Denote the solution of the above equation as $x_{0}$, i.e.,

$\displaystyle\left(\frac{1+x_{0}}{2}\right)^{p}-\left(\frac{1-x_{0}}{2}\right)^{p}$
$\displaystyle=-\frac{2B}{g}.$ (29)

Note that the left hand side of the above equation is a monotonically
increasing function of $x_{0}$ as $x_{0}\in[-1,1]$, with a minimum of $-1$ and
a maximum of $+1$. Therefore,

$\displaystyle-1\leq-\frac{2B}{g}\leq 1.$ (30)

This means that

$\displaystyle|g|\geq 2B.$ (31)

This is the necessary condition for a Dirac cone to exist.

It is also of interest to use the perturbation theory to solve the
eigenenergies near the Dirac cone. For sufficiently small $|k_{1}|$ and
$|k_{2}|$, we let

$\displaystyle x$ $\displaystyle\sim x_{0}+\chi,$ (32) $\displaystyle\beta$
$\displaystyle\sim B+\rho,$ (33) $\displaystyle|\gamma|^{2}$
$\displaystyle\sim 0+\eta,$ (34)

where $\chi$ is at least in the first order in $k_{1}$ and $k_{2}$, and $\rho$
and $\eta$ are at least in the second order in $k_{1}$ and $k_{2}$. Plugging
Eq. (32), Eq. (33) and Eq. (34) into Eq. (5), we have

$\displaystyle\frac{g^{2}p^{2}}{2^{2p+2}}\frac{1-x_{0}^{2}}{x_{0}^{2}}\left[(1+x_{0})^{p-1}+(1-x_{0})^{p-1}\right]^{2}\chi^{2}$
$\displaystyle\sim\eta.$ (35)

To this order, we get the correction to the parameter $x$,

$\displaystyle\chi$
$\displaystyle\sim\pm\frac{2^{p+1}}{gp}\frac{x_{0}}{(1+x_{0})^{p-1}+(1-x_{0})^{p-1}}\frac{\sqrt{J_{1}^{2}k_{1}^{2}+J_{2}^{2}k_{2}^{2}}}{\sqrt{1-x_{0}^{2}}}.$
(36)

Plugging this into Eq. (4), we find the expression for the eigenenergy,

$\displaystyle E$ $\displaystyle\sim E_{0}\left(1+p\frac{\chi}{x_{0}}\right),$
(37)

where the nonperturbed eigenenergy is

$\displaystyle E_{0}$
$\displaystyle=\frac{g}{2}\left[\left(\frac{1+x_{0}}{2}\right)^{p}+\left(\frac{1-x_{0}}{2}\right)^{p}\right].$
(38)

We can see clearly from the expansion of $E$ that there is a Dirac cone
structure at the origin, provided $|x_{0}|<1$, which corresponds to $|g|>2B$.
For $|g|<2B$, the system contains two smooth energy bands. At the critical
value $g=2B$ ($g=-2B$), a kink will develop on the lower (upper) band at
$k_{1}=k_{2}=0$. Once $g>2B$ ($g<-2B$), a 2D self-intersection structure,
i.e., a nonlinear Dirac cone, will appear from the lower (upper) band, whose
vertex is at $k_{1}=k_{2}=0$. This is true for any $p>0$. We show 5 plots with
different values of nonlinearity in Fig. 4, along section $k_{1}=0$ and with
$-0.1\pi\leq k_{2}\leq 0.1\pi$. In each plot, red dots are perturbative
eigenenergies around the Dirac point (or at the origin for $|g|\leq 2B$),
while blue lines are numerical solutions by solving Eq. (5) exactly. We can
see that the perturbative solutions perfectly match the numerical solutions
for sufficiently small $|k_{2}|$.

Figure 4: Numerical and perturbative solutions of $E$ for $p=2$ with
$g/B=-2.5,-2,1,2,2.5$, $B=2$ along section $k_{1}=0$. In these plots, the red
dots are perturbative eigenenergies near the Dirac cone and the blue lines are
numerical solutions.

##  Appendix B Dynamics around the Dirac cone

We solve the following Schrödinger equation perturbatively in terms of
adiabatic parameter $\varepsilon$ [1],

$\displaystyle i\partial_{t}\begin{bmatrix}\Psi_{1}\\\
\Psi_{2}\end{bmatrix}=\begin{bmatrix}\beta&\gamma\\\
\gamma^{*}&-\beta\end{bmatrix}\begin{bmatrix}\Psi_{1}\\\
\Psi_{2}\end{bmatrix}+g\begin{bmatrix}|\Psi_{1}|^{2p}\Psi_{1}\\\
|\Psi_{2}|^{2p}\Psi_{2}\end{bmatrix}.$ (39)

During the adiabatic following process, the quasimomenta are tuned
adiabatically,

$\displaystyle k_{1}$ $\displaystyle=k_{1}(\varepsilon t),\qquad
k_{2}=k_{2}(\varepsilon t),$ (40)

with $0<\varepsilon\ll 1$. Let

$\displaystyle\ket{\Psi}$ $\displaystyle=\begin{bmatrix}\Psi_{1}\\\
\Psi_{2}\end{bmatrix}=e^{i\theta}(\ket{\psi}+\varepsilon\ket{\phi})=e^{i\theta}\begin{bmatrix}\psi_{1}\\\
\psi_{2}\end{bmatrix}+\varepsilon e^{i\theta}\begin{bmatrix}\phi_{1}\\\
\phi_{2}\end{bmatrix},$ (41)

with

$\displaystyle\dot{\theta}$ $\displaystyle\sim-E+\varepsilon\alpha,$ (42)

where $\psi_{a}$ are the solutions to the eigenvalue problem in Eq. (2), and
$\varepsilon\phi_{a}$ are the first order corrections. As we shall see, the
dynamical phase comes from both $E$ and $\varepsilon\alpha$, and geometric
phase comes from only $\varepsilon\alpha$.

The solution to the adiabatic process deviates from the instantaneous
eigensolution by the order of $\varepsilon$, multiplied by a phase term. Note
that the increment of phase $\theta$ over a small duration of time $dt$ may
deviate from the contribution of dynamical phase by $O(\varepsilon)$, but the
accumulation of the $O(\varepsilon)$ term over the total time
$O(\varepsilon^{-1})$ has a contribution of $O(1)$.

The Hamiltonian can be expand in the power of $\varepsilon$ accordingly,

$\displaystyle\widehat{H}(\Psi)$
$\displaystyle\sim\widehat{H}(\psi)+\varepsilon\hat{h}(\psi,\phi),$ (43)

where $\widehat{H}(\psi)$ is given in Eq. (1) and $\hat{h}$ depends on both
$\psi$ and $\phi$ with a diagonal form,

$\displaystyle\hat{h}(\psi,\phi)$
$\displaystyle=gp\begin{bmatrix}|\psi_{1}|^{2p-2}(\psi_{1}^{*}\phi_{1}+\psi_{1}\phi_{1}^{*})&0\\\
0&|\psi_{2}|^{2p-2}(\psi_{2}^{*}\phi_{2}+\psi_{2}\phi_{2}^{*})\end{bmatrix}.$
(44)

To compute $\hat{h}$, we need to expand $|\Psi|^{2p}$. To the first order in
$\varepsilon$, we have

$\displaystyle|\psi_{a}+\varepsilon\phi_{a}|^{2p}$
$\displaystyle\sim|\psi_{a}|^{2p}\left[1+2\varepsilon
p\operatorname{Re}\left(\frac{\phi_{a}}{\psi_{a}}\right)\right]$
$\displaystyle=|\psi_{a}|^{2p}+\varepsilon
p|\psi_{a}|^{2p-2}(\psi_{a}^{*}\phi_{a}+\psi_{a}\phi_{a}^{*}).$ (45)

Plugging Eqs. (41), (42) and (45) into the time-dependent Schrödinger equation
in Eq. (39), up to the first order in $\varepsilon$, we obtain

$\displaystyle E\ket{\psi}-\varepsilon\alpha\ket{\psi}+\varepsilon
E\ket{\phi}+i\ket{\dot{\psi}}+i\varepsilon\ket{\dot{\phi}}$
$\displaystyle\sim\widehat{H}(\psi)\ket{\psi}+\varepsilon\hat{h}\ket{\psi}+\varepsilon\widehat{H}(\psi)\ket{\phi}.$
(46)

Note that the time derivative brings a factor of $\varepsilon$ because we are
in the adiabatic regime, thus the term $\varepsilon\ket{\dot{\phi}}$ is
actually in the order of $\varepsilon^{2}$ and it can be discarded. Apply the
instantaneous eigenvalue equation in Eq. (2), we get the equation for
$\ket{\phi}$,

$\displaystyle\left[E-\widehat{H}(\psi)\right]\ket{\phi}$
$\displaystyle\sim\alpha\ket{\psi}-\frac{i}{\varepsilon}\ket{\dot{\psi}}+\hat{h}\ket{\psi}.$
(47)

Multiply $\bra{\psi}$ from the left to Eq. (47), we get

$\displaystyle\varepsilon\alpha$
$\displaystyle=i\braket{\psi}{\dot{\psi}}-\varepsilon\bra{\psi}\hat{h}\ket{\psi}.$
(48)

After a lengthy calculation, we find the solution to the correction of
wavefunction as

$\displaystyle\varepsilon\phi_{1}$
$\displaystyle=-\frac{x(1-x)\sqrt{1+x}}{4\sqrt{2}}\frac{\dot{\varphi}}{\Delta}-i\frac{x}{4\sqrt{2(1+x)}}\frac{\dot{x}}{\Delta^{\prime}},$
(49) $\displaystyle\varepsilon\phi_{2}$
$\displaystyle=\frac{x(1+x)\sqrt{1-x}}{4\sqrt{2}}\frac{\dot{\varphi}}{\Delta}e^{i\varphi}+i\frac{x}{4\sqrt{2(1-x)}}\frac{\dot{x}}{\Delta^{\prime}}e^{i\varphi},$
(50)

where we introduce two quantities, $\Delta$ in Eq. (11) and $\Delta^{\prime}$
is defined as

$\displaystyle\Delta^{\prime}$
$\displaystyle:=\beta+\frac{g}{2}\left[\left(\frac{1+x}{2}\right)^{p}-\left(\frac{1-x}{2}\right)^{p}\right].$
(51)

It turns out that the $\dot{x}/\Delta^{\prime}$ terms do not contribute to
$\hat{h}$,

$\displaystyle\hat{h}=gp\frac{x(1-x^{2})}{4\Delta}\frac{\dot{\varphi}}{\varepsilon}\begin{bmatrix}-\left(\frac{1+x}{2}\right)^{p-1}&0\\\
0&\left(\frac{1-x}{2}\right)^{p-1}\end{bmatrix}.$ (52)

Putting all together, the change rate of the overall phase is given in Eq.
(10).

## References

  * Bomantara  _et al._ [2017] R. W. Bomantara, W. Zhao, L. Zhou, and J. Gong, Nonlinear dirac cones, Phys. Rev. B 96, 121406 (2017).
  * New [2011] G. New,  _Introduction to Nonlinear Optics_ (Cambridge University Press, 2011).
  * Wu and Niu [2000] B. Wu and Q. Niu, Nonlinear landau-zener tunneling, Phys. Rev. A 61, 023402 (2000).
  * Tuloup  _et al._ [2020] T. Tuloup, R. W. Bomantara, C. H. Lee, and J. Gong, Nonlinearity induced topological physics in momentum space and real space, Phys. Rev. B 102, 115411 (2020).
  * Liu  _et al._ [2002] J. Liu, L. Fu, B.-Y. Ou, S.-G. Chen, D.-I. Choi, B. Wu, and Q. Niu, Theory of nonlinear landau-zener tunneling, Phys. Rev. A 66, 023404 (2002).
  * Witthaut  _et al._ [2006] D. Witthaut, E. M. Graefe, and H. J. Korsch, Towards a generalized landau-zener formula for an interacting bose-einstein condensate in a two-level system, Phys. Rev. A 73, 063609 (2006).
  * Zhang  _et al._ [2008a] Q. Zhang, P. Hänggi, and J. Gong, Two-mode bose-einstein condensate in a high-frequency driving field that directly couples the two modes, Phys. Rev. A 77, 053607 (2008a).
  * Zhang  _et al._ [2008b] Q. Zhang, P. Hänggi, and J. Gong, Nonlinear landau–zener processes in a periodic driving field, New Journal of Physics 10, 073008 (2008b).
  * Qi  _et al._ [2006] X.-L. Qi, Y.-S. Wu, and S.-C. Zhang, Topological quantization of the spin hall effect in two-dimensional paramagnetic semiconductors, Phys. Rev. B 74, 085308 (2006).
  * Novoselov  _et al._ [2005] K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, M. I. Katsnelson, I. V. Grigorieva, S. V. Dubonos, and A. A. Firsov, Two-dimensional gas of massless dirac fermions in graphene, Nature 438, 197 (2005).
  * Zhang  _et al._ [2005] Y. Zhang, Y.-W. Tan, H. L. Stormer, and P. Kim, Experimental observation of the quantum hall effect and berry's phase in graphene, Nature 438, 201 (2005).
  * Ando  _et al._ [1998] T. Ando, T. Nakanishi, and R. Saito, Berry’s phase and absence of back scattering in carbon nanotubes, Journal of the Physical Society of Japan 67, 2857 (1998), https://doi.org/10.1143/JPSJ.67.2857 .
  * Mikitik and Sharlai [1999] G. P. Mikitik and Y. V. Sharlai, Manifestation of berry’s phase in metal physics, Phys. Rev. Lett. 82, 2147 (1999).
  * Tuloup  _et al._ [2022] T. Tuloup, R. W. Bomantara, and J. Gong, Topological characteristics of gap closing points in nonlinear weyl semimetals, Phys. Rev. B 106, 195411 (2022).
  * Milovanov  _et al._ [2021] A. V. Milovanov, J. J. Rasmussen, and G. Dif-Pradalier, Self-consistent model of the plasma staircase and nonlinear schrödinger equation with subquadratic power nonlinearity, Phys. Rev. E 103, 052218 (2021).
  * Wazwaz [2006] A.-M. Wazwaz, Exact solutions for the fourth order nonlinear schrodinger equations with cubic and power law nonlinearities, Mathematical and Computer Modelling 43, 802 (2006).
  * Kilic and Inc [2017] B. Kilic and M. Inc, Optical solitons for the schrödinger-hirota equation with power law nonlinearity by the bäcklund transformation, Optik 138, 64 (2017).
  * Sulem and Sulem [2004] C. Sulem and P. L. Sulem, The nonlinear schrödinger equation : self-focusing and wave collapse (2004).
  * Osman  _et al._ [2019] M. Osman, D. Lu, and M. M. Khater, A study of optical wave propagation in the nonautonomous schrödinger-hirota equation with power-law nonlinearity, Results in Physics 13, 102157 (2019).
  * Mirzazadeh  _et al._ [2017] M. Mirzazadeh, M. Ekici, Q. Zhou, and A. Biswas, Exact solitons to generalized resonant dispersive nonlinear schrödinger’s equation with power law nonlinearity, Optik 130, 178 (2017).
  * Dai  _et al._ [2017] C.-Q. Dai, X.-F. Zhang, Y. Fan, and L. Chen, Localized modes of the (n+1)-dimensional schrödinger equation with power-law nonlinearities in pt-symmetric potentials, Communications in Nonlinear Science and Numerical Simulation 43, 239 (2017).
  * Mirzazadeh  _et al._ [2014] M. Mirzazadeh, M. Eslami, B. F. Vajargah, and A. Biswas, Optical solitons and optical rogons of generalized resonant dispersive nonlinear schrödinger’s equation with power law nonlinearity, Optik 125, 4246 (2014).
  * Biswas and Konar [2006] A. Biswas and S. Konar, Introduction to non-kerr law optical solitons (2006).


n}\left\\{h(Z_{i},Z_{j})-h(Z_{i},Z^{\prime}_{j})-h(Z^{\prime}_{i},Z_{j})+h(Z^{\prime}_{i},Z^{\prime}_{j})\right\\}$
$\displaystyle=$ $\displaystyle
EE_{\epsilon}\sup\limits_{h\in\mathcal{H}}\frac{1}{n(n-1)}\sum_{1\leq i\neq
j\leq
n}\epsilon_{i}\epsilon_{j}\left\\{h(Z_{i},Z_{j})-h(Z_{i},Z^{\prime}_{j})-h(Z^{\prime}_{i},Z_{j})+h(Z^{\prime}_{i},Z^{\prime}_{j})\right\\},$

where the inequality follows from Jensen’s inequality and the last equality is
due to Fubini’s Theorem and a direct corollary of the symmetrization
inequality in Sherman (1994).

We prove the second inequality in (4.3) through the chaining technique. Let
$\alpha_{0}=b$ and for any $t\in\mathbb{N}_{+}$, let $\alpha_{t}=2^{-t}b$. For
each $t$, let $\mathcal{C}_{t}$ be a $\alpha_{t}$-cover of $\mathcal{H}$
w.r.t. $d_{U,n}(\cdot,\cdot)$ such that
$|\mathcal{C}_{t}|=N(\alpha_{t},\mathcal{H},d_{U,n})$. For each
$h\in\mathcal{H}$, there exists a function $\hat{h}_{t}\in\mathcal{C}_{t}$
such that $d_{U,n}(h,\hat{h}_{t})\leq\alpha_{t}$. Let $\hat{h}_{0}\equiv 0$
and for any $T\in\mathbb{N}_{+}$, we have the following chaining expression
for $h$:

$h=h-\hat{h}_{T}+\sum_{t=1}^{T}(\hat{h}_{t}-\hat{h}_{t-1}).$

Let

$G(h,Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j})=h(Z_{i},Z_{j})-h(Z_{i},Z^{\prime}_{j})-h(Z^{\prime}_{i},Z_{j})+h(Z^{\prime}_{i},Z^{\prime}_{j}).$

Hence for any $T\in\mathbb{N}_{+}$, we have

$\displaystyle\mathcal{S}_{n}^{\textup{Deg}}(\mathcal{H})$ $\displaystyle=$
$\displaystyle
E_{\epsilon}\sup\limits_{h\in\mathcal{H}}\frac{1}{n(n-1)}\sum_{1\leq i\neq
j\leq
n}\epsilon_{i}\epsilon_{j}G(h,Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j})$
(S.4) $\displaystyle\leq$ $\displaystyle
E_{\epsilon}\sup\limits_{h\in\mathcal{H}}\frac{1}{n(n-1)}\sum_{1\leq i\neq
j\leq
n}\epsilon_{i}\epsilon_{j}G(h-\hat{h}_{T},Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j})$
$\displaystyle+$
$\displaystyle\sum_{t=1}^{T}E_{\epsilon}\sup\limits_{h\in\mathcal{H}}\frac{1}{n(n-1)}\sum_{1\leq
i\neq j\leq
n}\epsilon_{i}\epsilon_{j}G(\hat{h}_{t}-\hat{h}_{t-1},Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j})$
$\displaystyle\leq$
$\displaystyle\alpha_{T}+\sum_{t=1}^{T}E_{\epsilon}\sup\limits_{h\in\mathcal{H}}\frac{1}{n(n-1)}\sum_{1\leq
i\neq j\leq
n}\epsilon_{i}\epsilon_{j}G(\hat{h}_{t}-\hat{h}_{t-1},Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j}),$

where the last inequality is due to

$\displaystyle\frac{1}{n(n-1)}\sum_{1\leq i\neq j\leq
n}\epsilon_{i}\epsilon_{j}G(h-\hat{h}_{T},Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j})$
$\displaystyle\leq$ $\displaystyle\left\\{\frac{1}{n(n-1)}\sum_{1\leq i\neq
j\leq
n}(\epsilon_{i}\epsilon_{j})^{2}\right\\}^{1/2}\left\\{\frac{1}{n(n-1)}\sum_{1\leq
i\neq j\leq
n}G^{2}(h-\hat{h}_{T},Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j})\right\\}^{1/2}$
$\displaystyle=$ $\displaystyle d_{U,n}(h,\hat{h}_{T})\leq\alpha_{T}.$

Now the second term in (S.4) is the summation of empirical Rademacher chaos
w.r.t. the function classes
$\mathcal{C}_{t}-\mathcal{C}_{t-1}=\\{h_{t}-h_{t-1}:h_{t}\in\mathcal{C}_{t},h_{t-1}\in\mathcal{C}_{t-1}\\},t=1,2,\ldots,T$.
For any measurable function $h\in\mathcal{H}$, let

$M(h)=(a^{h}_{ij})\in\mathbb{R}^{n\times
n}~{}\text{with}~{}a^{h}_{ij}=0,~{}\text{if}~{}i=j~{}\text{and
otherwise}~{}a^{h}_{ij}=G(\hat{h}_{t}-\hat{h}_{t-1},Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j}).$

Then for any $h_{t}\in\mathcal{C}_{t},h_{t-1}\in\mathcal{C}_{t-1}$,

$\displaystyle 1/\sqrt{n(n-1)}\|M(\hat{h}_{t}-\hat{h}_{t-1})\|_{F}$
$\displaystyle=$ $\displaystyle\left\\{\frac{1}{n(n-1)}\sum_{1\leq i\neq j\leq
n}G^{2}(\hat{h}_{t}-\hat{h}_{t-1},Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j})\right\\}^{1/2}$
$\displaystyle=$ $\displaystyle d_{U,n}(\hat{h}_{t},\hat{h}_{t-1})$
$\displaystyle\leq$ $\displaystyle
d_{U,n}(h,\hat{h}_{t})+d_{U,n}(h,\hat{h}_{t-1})$ $\displaystyle\leq$
$\displaystyle\alpha_{t}+\alpha_{t-1}$ $\displaystyle=$ $\displaystyle
3\alpha_{t}.$

Applying Lemma C.4 to $\mathcal{C}_{t}-\mathcal{C}_{t-1},t=1,2,\ldots,T$, we
have there exists a universal constant $C$ such that

$\displaystyle\mathcal{S}_{n}^{\textup{Deg}}(\mathcal{H})$ $\displaystyle\leq$
$\displaystyle\alpha_{T}+\sum_{t=1}^{T}E_{\epsilon}\sup\limits_{h\in\mathcal{H}}\frac{1}{n(n-1)}\sum_{1\leq
i\neq j\leq
n}\epsilon_{i}\epsilon_{j}G(\hat{h}_{t}-\hat{h}_{t-1},Z_{i},Z_{j},Z^{\prime}_{i},Z^{\prime}_{j})$
$\displaystyle\leq$
$\displaystyle\alpha_{T}+\frac{3C\sum_{t=1}^{T}\alpha_{t}\log(|\mathcal{C}_{t}|\cdot|\mathcal{C}_{t-1}|)}{\sqrt{n(n-1)}}$
$\displaystyle\leq$
$\displaystyle\alpha_{T}+\frac{6C\sum_{t=1}^{T}\alpha_{t}\log(|\mathcal{C}_{t}|)}{\sqrt{n(n-1)}}$
$\displaystyle=$
$\displaystyle\alpha_{T}+\frac{12C\sum_{t=1}^{T}(\alpha_{t}-\alpha_{t+1})\log(N(\alpha_{t},\mathcal{H},d_{U,n}))}{\sqrt{n(n-1)}}$
$\displaystyle\leq$
$\displaystyle\alpha_{T}+\frac{12C\int_{\alpha_{T+1}}^{\alpha_{0}}\log(N(u,\mathcal{H},d_{U,n}))du}{\sqrt{n(n-1)}}$
$\displaystyle\leq$
$\displaystyle\alpha_{T}+\frac{24C\int_{\alpha_{T+1}}^{\alpha_{0}}\log(N(u,\mathcal{H},d_{U,n}))du}{n}.$

For any $0\textless\delta\textless b$, let
$T=\min\\{K\in\mathbb{N}_{+}:2^{-K}b\leq\delta\\}$, then

$\displaystyle\mathcal{S}_{n}^{\textup{Deg}}(\mathcal{H})$ $\displaystyle\leq$
$\displaystyle\delta+\frac{24C\int_{\delta/2}^{b}\log(N(u,\mathcal{H},d_{U,n}))du}{n}$
$\displaystyle\leq$
$\displaystyle\delta+\frac{24C\int_{\delta/4}^{b}\log(\mathcal{N}_{2n(n-1)}(u,\|\cdot\|_{\infty},\mathcal{H}))du}{n}$

where the last inequality follows from (S.1). This completes the proof of
Lemma 4.3. $\Box$

######  Lemma C.9.

(a) Let
$\mathcal{D}\otimes\mathcal{R}=\\{h(z)=D(y,R(x)):D:\mathbb{R}^{d_{Y}+d_{0}}\to\mathbb{R}\in\mathcal{D},R:\mathbb{R}^{d_{X}}\to\mathbb{R}^{d_{0}}\in\mathcal{R}\\}$.
If $\mathcal{D}$ is a class of functions implemented by a ReLU network with
width $\mathcal{W}_{\mathcal{D}}$ and depth $\mathcal{L}_{\mathcal{D}}$, and
$\mathcal{R}$ is a class of functions implemented by a ReLU network with width
$\mathcal{W}_{\mathcal{R}}$ and depth $\mathcal{L}_{\mathcal{R}}$, then each
element of $\mathcal{D}\otimes\mathcal{R}$ can be implemented by a ReLU
network with width
$\max\\{\mathcal{W}_{\mathcal{D}},\mathcal{W}_{\mathcal{R}}+2d_{Y},2(d_{Y}+d_{0})\\}$
and depth $\mathcal{L}_{\mathcal{D}}+\mathcal{L}_{\mathcal{R}}+1$.

(b) Let
$\mathcal{Q}\circ\mathcal{R}=\\{h(z)=Q(R(x)):Q:\mathbb{R}^{d_{0}}\to\mathbb{R}^{d_{\textup{out}}}\in\mathcal{Q},R:\mathbb{R}^{d_{X}}\to\mathbb{R}^{d_{0}}\in\mathcal{R}\\}$.
If $\mathcal{Q}$ is a class of functions implemented by a ReLU network with
width $\mathcal{W}_{\mathcal{Q}}$ and depth $\mathcal{L}_{\mathcal{Q}}$, and
$\mathcal{R}$ is a class of functions implemented by a ReLU network with width
$\mathcal{W}_{\mathcal{R}}$ and depth $\mathcal{L}_{\mathcal{R}}$, then each
element of $\mathcal{Q}\circ\mathcal{R}$ can be implemented by a ReLU network
with width
$\max\\{\mathcal{W}_{\mathcal{Q}},\mathcal{W}_{\mathcal{R}},2d_{0}\\}$ and
depth $\mathcal{L}_{\mathcal{Q}}+\mathcal{L}_{\mathcal{R}}+1$.

###### Proof of Lemma C.9.

The proof follows from the method in Subsection B.1.1 in Nakada and Imaizumi
(2020), which shows how to construct a neural network from sub-neural networks
through concatenation, parallelization, and realization of identity function.
$\Box$

##  Appendix D Additional numerical studies and experiment details

This Appendix contains Table S.1 which includes some additional simulation
results.

Table S.1: Distance correlation (DC), average prediction errors (APE), and
their standard errors (based on 6-fold validation and best hyperparameters
selected by validation method) for $p=30$

|  | Model I | Model II | Model III | Model IV
---|---|---|---|---|---
| Method | DC | APE | DC | APE | DC | APE | DC | APE
S(i) | MSRL | .95(.00) | 0.30(.00) | .86(.01) |  0.30(.01) | .92(0.01) |  0.54(.05) | .94(.00) |  0.31(.01)
SIR | .97(.00) |  0.25(.01) | .00(.00) | 0.79(.01) | .20(0.00) | 1.70(.02) | .09(.02) | 1.22(.01)
SAVE | .96(.00) |  0.25(.01) | .00(.00) | 0.79(.01) | .20(0.01) | 1.71(.04) | .08(.02) | 1.23(.01)
GSIR | .94(.00) | 0.32(.01) | .84(.01) | 0.32(.02) | .74(0.01) | 0.72(.01) | .65(.04) | 0.34(.01)
GSAVE | .96(.00) | 0.27(.01) | .35(.10) | 0.63(.05) | .42(0.04) | 1.29(.06) | .32(.04) | 0.86(.09)
S(ii) | MSRL | .93(.00) | 0.29(.00) | .88(.01) |  0.31(.01) | .87(0.03) |  0.60(.08) | .92(.01) | 0.43(.05)
SIR | .95(.00) |  0.25(.00) | .00(.00) | 0.90(.04) | .16(0.02) | 1.51(.06) | .05(.02) | 1.44(.10)
SAVE | .95(.00) |  0.25(.00) | .00(.00) | 0.90(.04) | .15(0.02) | 1.54(.06) | .06(.02) | 1.44(.10)
GSIR | .93(.01) | 0.30(.01) | .87(.01) | 0.32(.01) | .72(0.03) | 0.68(.03) | .80(.03) |  0.38(.04)
GSAVE | .95(.00) | 0.26(.00) | .42(.02) | 0.68(.04) | .33(0.04) | 1.32(.04) | .30(.08) | 0.99(.13)
S(iii) | MSRL | .99(.00) | 0.30(.01) | .97(.00) |  0.38(.01) | .94(0.01) |  2.01(.21) | .97(.01) |  0.66(.08)
SIR | .99(.00) |  0.25(.00) | .00(.00) | 2.03(.04) | .33(0.02) | 6.36(.23) | .02(.01) | 3.65(.07)
SAVE | .99(.00) | 0.27(.00) | .00(.00) | 2.04(.03) | .32(0.02) | 6.41(.25) | .03(.01) | 3.63(.06)
GSIR | .84(.04) | 0.97(.14) | .93(.02) | 0.52(.06) | .71(0.04) | 2.35(.31) | .72(.01) | 0.96(.21)
GSAVE | .99(.00) | 0.34(.03) | .66(.13) | 1.16(.24) | .64(0.01) | 3.85(.42) | .36(.02) | 2.47(.09)
S(iv) | MSRL | .96(.00) | 0.28(.01) | .89(.01) |  0.33(.01) | .90(0.02) |  0.87(.08) | .93(.01) | 0.41(.02)
SIR | .97(.00) |  0.25(.00) | .00(.00) | 0.99(.05) | .19(0.02) | 2.38(.20) | .09(.01) | 1.44(.09)
SAVE | .97(.00) |  0.25(.00) | .00(.00) | 0.99(.06) | .18(0.03) | 2.39(.21) | .09(.01) | 1.44(.09)
GSIR | .95(.00) | 0.29(.01) | .89(.01) |  0.33(.01) | .68(0.06) | 0.93(.20) | .75(.04) |  0.36(.05)
GSAVE | .95(.00) | 0.29(.01) | .70(.02) | 0.54(.03) | .53(0.06) | 1.52(.22) | .50(.06) | 0.81(.06)

###  D.1 Experimental Details

Python codes will be made publicly available on a GitHub repository, we still
Here we provide the details for the numerical experimental details. In all the
experiments, we set $\lambda=2$ and use the Adam optimizer (Kingma and Ba,
2014) in Pytorch with learning rate lr and weight decay parameter wd specified
in each experiment. The activation function for the neural networks is the
LeakyReLU. The batch size in all experiments is $512$.

###  D.2 Simulation for Model (2.10)

The training data sample size is $4000$ and no ES criterion is applied. The
maximum number of epoches is $3000$. When the reference distribution is
$\gamma_{U}$, $\mathcal{R}$ has $2$ hidden layer with width $(64,64)$. When
the reference distribution is $\gamma_{S}$, $\mathcal{R}$ has $3$ hidden layer
with width $(64,64,64)$. Both $\mathcal{D}$ and $\mathcal{Q}$ has $2$ hidden
layers with widths $(32,32)$. lr is $0.003$ and wd is $0$. For GSIR and GSAVE,
we use some default tuning parameters. And the kernel estimates are based on
$10000$ new samples.

###  D.3 Simulation Study in Subsection 7.1

In this experiment, the total $6000$ data points are first divided into $5000$
training-validation data points and $1000$ testing data points and then $5000$
training-validation data points are divided into $4000$ training data points
and $1000$ validation data points.

In MSRL training, ES criterion is applied and it is the empirical distance
correlation (DC) calculated based on the validation data. The maximum number
of epochs is $1000$ and wd is $0.0001$. The patience is $200$, where patience
is the number of epochs until termination if no improvement is made on the
validation dataset. The learning procedure is replicated $10$ times with
different initial weights for the parameters of neural networks and we choose
the one with the largest empirical DC on the validation dataset among the ten
representation estimate as the ultimate representation estimate. $\mathcal{R}$
has $3$ hidden layers with widths $(32,16,8)$. Both $\mathcal{D}$ and
$\mathcal{Q}$ has $2$ hidden layers with widths $(16,8)$. For Model I & II,
$\textit{lr}=0.001$, $d_{0}=1$. For Model III $\textit{lr}=0.0003$, $d_{0}=2$.
For Model IV, $\textit{lr}=0.001$, $d_{0}=2$.

For SIR and SAVE, the tuning parameter is the slice number sn. For each
$\textit{sn}\in\\{5,10,15,20,25,30\\}$, we use the validation data to
calculate the empirical DC between the response and the representation learned
from training data. We use the slice number with the largest empirical DC as
the ultimate choice for the tuning parameter sn. For GSIR and GSAVE, we apply
the tuning parameter selection approach provided in Li (2018) to the
validation data.

###  D.4 Real Datasets

For each data split, $20\%$ of data is used as the testing data and then from
the other $80\%$ of the data, we randomly select $3000$ data points as the
validation data and the rest part is the training data.

For GSIR and GSAVE, limited by our computational capability, we only subsample
$4000$, $2000$, $3000$ data points from the above training, validation, and
testing data as the new training, validation, and testing data, and conduct
the estimation procedure we have stated in subsection D.3

In the superconductivty dataset, for our MSRL, $\mathcal{R}$ has $3$ hidden
layers with widths $(128,128,128)$. Both $\mathcal{D}$ and $\mathcal{Q}$ has
$2$ hidden layers with widths $(64,64)$. The maximum number of epochs is
$2000$, $\textit{lr}=0.0001,\textit{wd}=0.001$ and the patience is $400$. The
number of different initial weights is still $10$. For the goal of data
visualization, we randomly sample $800$ data points from the testing data for
one fold to conduct plotting.

In the pole-telecommunication dataset, for MSRL, $\mathcal{R}$ has $2$ hidden
layers with widths $(30,25)$. Both $\mathcal{D}$ and $\mathcal{Q}$ has $2$
hidden layers with widths $(16,8)$. The maximum number of epochs is $1000$,
$\textit{lr}=0.001,\textit{wd}=0.0001$ and the patience is $200$. The number
of different initial weights is still $10$. For the goal of data
visualization, we randomly sample $400$ data points from the testing data for
one fold to conduct plotting.

###  D.5 The extra plots for the superconductivity and pole-telecommunication
datasets

In Figure S.1, we plot the critical temperature against each component of the
learned representation with $d_{0}=5$ for SIR and SAVE. We fit a quadratic
model to check how the components impact the response.

SIR

SAVE

Figure S.1: SIR and SAVE: The plots of the critical temperature against each
component of the learned representation by SIR and SAVE for the
superconductivity dataset. The red line is a fitting quadratic model.

In Figure S.2, we plot the transformed data based on each pair of the
components of the learned representation with $d_{0}=5$ for SIR and SAVE.

SIR

SAVE

Figure S.2: The plots of the transformed data based on each pair of the
components of the learned representation by SIR and SAVE for the Pole-
Telecommunication dataset. Bule ’O’ and red ’$\Delta$’ means that $Y\geq 50$
and $Y\textless 50$, respectively.

###  D.6 Intrinsic Dimensionality Determination of Model (7.1)

The number of total data points is $8000$. $8000$ data points are first
divided into $6000$ training-validation data points and $2000$ testing sdata
points and then $6000$ training-validation data points are divided into $4000$
training data points and $2000$ validation data points. $\mathcal{R}$ has $1$
hidden layer with width $128$. Both $\mathcal{D}$ and $\mathcal{Q}$ also has
$1$ hidden layer with width $64$. The maximum number of epochs is $2000$,
$\textit{lr}=0.001,\textit{wd}=0.0001$ and the patience is $200$.

## References

  * Ali and Silvey (1966) [author] Ali, Syed MumtazS. M. and Silvey, Samuel DS. D. (1966). A general class of coefficients of divergence of one distribution from another. Journal of the Royal Statistical Society: Series B (Methodological) 28 131–142.
  * Anthony and Bartlett (1999) [author] Anthony, MartinM. and Bartlett, Peter L.P. L. (1999). Neural Network Learning: Theoretical Foundations. Cambridge University Press.
  * Bartlett et al. (2019) [author] Bartlett, Peter L.P. L., Harvey, NickN., Liaw, ChristopherC. and Mehrabian, AbbasA. (2019). Nearly-tight VC-dimension and pseudodimension bounds for piecewise linear neural networks. Journal of Machine Learning Research 20 1 – 17.
  * Berrett and Samworth (2019) [author] Berrett, T BT. B. and Samworth, R JR. J. (2019). Nonparametric independence testing via mutual information. Biometrika 106 547-566.
  * Cook (1998) [author] Cook, R. D.R. D. (1998). Regression Graphics: Ideas for Studying Regressions through Graphics. New York: Wiley.
  * Cook (2007) [author] Cook, R. DennisR. D. (2007). Fisher lecture: dimension reduction in regression. Statistical Science 22 1 – 26.
  * Cook (2018) [author] Cook, R. DennisR. D. (2018). Principal components, sufficient dimension reduction, and envelopes. Annual Review of Statistics and Its Application 5 533-559.
  * Cook and Weisberg (1991) [author] Cook, R. D.R. D. and Weisberg, S.S. (1991). Discussion of ”Sliced inverse regression for dimension reduction” by K.-C. Li. Journal of the American Statistical Association 86 328–332.
  * Cover and Thomas (2006) [author] Cover, TMT. and Thomas, Joy AJ. A. (2006). Elements of Information Theory, 2 ed. Wiley-Interscience.
  * Esposito, Gastpar and Issa (2020) [author] Esposito, Amedeo RobertoA. R., Gastpar, MichaelM. and Issa, IbrahimI. (2020). Robust generalization via f-mutual information. 2020 IEEE International Symposium on Information Theory (ISIT) 2723-2728.
  * Farrell, Liang and Misra (2021) [author] Farrell, Max H.M. H., Liang, TengyuanT. and Misra, SanjogS. (2021). Deep neural networks for estimation and inference. Econometrica 89 181 – 213.
  * Fukumizu et al. (2009) [author] Fukumizu, KenjiK., Gretton, ArthurA., Lanckriet, GertG., Schölkopf, BernhardB. and Sriperumbudur, Bharath K.B. K. (2009). Kernel choice and classifiability for RKHS embeddings of probability distributions. In Advances in Neural Information Processing Systems 22\.
  * Goodfellow et al. (2014) [author] Goodfellow, IanI., Pouget-Abadie, JeanJ., Mirza, MehdiM., Xu, BingB., Warde-Farley, DavidD., Ozair, SherjilS., Courville, AaronA. and Bengio, YoshuaY. (2014). Generative adversarial nets. Advances in Neural Information Processing Systems 27\.
  * Hamidieh (2018) [author] Hamidieh, KamK. (2018). A data-driven statistical model for predicting the critical temperature of a superconductor. Computational Materials Science 154 346–354.
  * Jiao et al. (2021) [author] Jiao, YulingY., Shen, GuohaoG., Lin, YuanyuanY. and Huang, JianJ. (2021). Deep nonparametric regression on approximately low-dimensional manifolds. arXiv:2104.06708.
  * Kingma and Ba (2014) [author] Kingma, Diederik PD. P. and Ba, JimmyJ. (2014). Adam: A method for stochastic optimization. arXiv:1412.6980.
  * Lee, Li and Chiaromonte (2013) [author] Lee, Kuang-YaoK.-Y., Li, BingB. and Chiaromonte, FrancescaF. (2013). A general theory for nonlinear sufficient dimension reduction: Formulation and estimation. The Annals of Statistics 41 221 – 249.
  * Li (1991) [author] Li, K. C.K. C. (1991). Sliced inverse regression for dimension reduction. Journal of the American Statistical Association 86 316–327.
  * Li (1992) [author] Li, K. C.K. C. (1992). On principal Hessian directions for data visualization and dimension reduction: another application of Stein’s lemma. Journal of the American Statistical Association 87 1025–1039.
  * Li (2018) [author] Li, BingB. (2018). Sufficient Dimension Reduction: Methods and Applications with R. Chapman and Hall/CRC.
  * Li, Artemiou and Li (2011) [author] Li, BingB., Artemiou, AndreasA. and Li, LexinL. (2011). Principal support vector machines for linear and nonlinear sufficient dimension reduction. The Annals of Statistics 39 3182 – 3210.
  * Li and Song (2017) [author] Li, BingB. and Song, JunJ. (2017). Nonlinear sufficient dimension reduction for functional data. The Annals of Statistics 45 1059 – 1095.
  * Liese and Vajda (2006) [author] Liese, FriedrichF. and Vajda, IgorI. (2006). On divergences and informations in statistics and information theory. Information Theory, IEEE Transactions on 52 4394 - 4412.
  * Ma and Zhu (2013) [author] Ma, Y. Y.Y. Y. and Zhu, L. P.L. P. (2013). A review on dimension reduction. International Statistical Review 81 134–150.
  * Nakada and Imaizumi (2020) [author] Nakada, RyumeiR. and Imaizumi, MasaakiM. (2020). Adaptive approximation and estimation of deep neural network with intrinsic dimensionality. Journal of Machine Learning Research 21 1 – 38.
  * Nguyen, Wainwright and Jordan (2010) [author] Nguyen, XuanLongX., Wainwright, MartinM. and Jordan, MichaelM. (2010). Estimating divergence functionals and the likelihood ratio by convex risk minimization. IEEE Transactions on Information Theory 56 5847 – 5861.
  * Nowozin, Cseke and Tomioka (2016) [author] Nowozin, SebastianS., Cseke, BotondB. and Tomioka, RyotaR. (2016). f-gan: Training generative neural samplers using variational divergence minimization. Advances in Neural Information Processing Systems 29\.
  * Polyanskiy and Wu (2017) [author] Polyanskiy, Y.Y. and Wu, Y.Y. (2017). Lecture notes on Information Theory. MIT (6.441), UIUC (ECE 563), Yale (STAT 664).
  * Sason (2019) [author] Sason, IgalI. (2019). On data-processing and majorization inequalities for f-divergences with applications. Entropy 21 1022\.
  * Shen, Yang and Zhang (2020) [author] Shen, ZuoweiZ., Yang, HaizhaoH. and Zhang, ShijunS. (2020). Deep network approximation characterized by number of neurons. Communications in Computational Physics 28 1768–1811.
  * Sherman (1994) [author] Sherman, Robert P.R. P. (1994). Maximal inequalities for degenerate $U$-processes with applications to optimization estimators. The Annals of Statistics 22 439 – 459.
  * Shin et al. (2017) [author] Shin, Seung JunS. J., Wu, YichaoY., Zhang, Hao HelenH. H. and Liu, YufengY. (2017). Principal weighted support vector machines for sufficient dimension reduction in binary classification. Biometrika 104 67-81.
  * Székely, Rizzo and Bakirov (2007) [author] Székely, Gábor JG. J., Rizzo, Maria LM. L. and Bakirov, Nail KN. K. (2007). Measuring and testing dependence by correlation of distances. The Annals of Statistics 35 2769–2794.
  * Tishby, Pereira and Bialek (1999) [author] Tishby, NaftaliN., Pereira, Fernando C.F. C. and Bialek, WilliamW. (1999). The information bottleneck method. In Proc. of the 37-th Annual Allerton Conference on Communication, Control and Computing 368-377.
  * Tishby and Zaslavsky (2015) [author] Tishby, NaftaliN. and Zaslavsky, NogaN. (2015). Deep learning and the information bottleneck principle. In 2015 IEEE Information Theory Workshop (ITW).
  * Tsybakov (2008) [author] Tsybakov, AlexandreA. (2008). Introduction to Nonparametric Estimation. Springer Science & Business Media.
  * Vershynin (2018) [author] Vershynin, RomanR. (2018). High-Dimensional Probability: An Introduction with Applications in Data Science. Cambridge University Press.
  * Villani (2009) [author] Villani, CédricC. (2009). Optimal transport: old and new 338\. Springer.
  * Weiss and Indurkhya (1995) [author] Weiss, SholomS. and Indurkhya, NitinN. (1995). Rule-based machine learning methods for functional prediction. Journal of Artificial Intelligence Research 3 383–403.
  * Xia et al. (2002) [author] Xia, Y. C.Y. C., Tong, H.H., Li, W. K.W. K. and Zhu, L. X.L. X. (2002). An adaptive estimation of dimension reduction space. Journal of the Royal Statistical Society, Series B 64 363–410.
  * Ying and Campbell (2010) [author] Ying, YimingY. and Campbell, ColinC. (2010). Rademacher chaos complexities for learning the kernel problem. Neural computation 22 2858-86.
  * Zhu, Miao and Peng (2006) [author] Zhu, L. X.L. X., Miao, B. Q.B. Q. and Peng, H.H. (2006). On sliced inverse regression with high-dimensional covariates. Journal of the American Statistical Association 101 630-643.
  * Zhu, Zhu and Feng (2010) [author] Zhu, L. P.L. P., Zhu, L. X.L. X. and Feng, Z. H.Z. H. (2010). Dimension reduction in regressions through cumulative slicing estimation. Journal of the American Statistical Association 105 1455–1466.


# Evaluation of flux expulsion and flux trapping sensitivity of SRF cavities
fabricated from cold work Nb sheet with successive heat treatment††thanks:
Work supported by the U.S. Department of Energy, Office of Science, Office of
Nuclear Physics under contract DE-AC05-06OR23177 and U.S.Department of Energy,
Office of Science, Office of High Energy Physics under Awards No. DE-SC
0009960.

B. D. Khanal1 and P. Dhakal2
1Center for Accelerator Science<EMAIL_ADDRESS>Department of Physics
Old Dominion University  Norfolk  VA 23529  USA
2Thomas Jefferson National Accelerator Facility  Newport News  VA 23606  USA

###### Abstract

The main source of RF losses leading to lower quality factor of
superconducting radio-frequency cavities is due to the residual magnetic flux
trapped during cool-down. The loss due to flux trapping is more pronounced for
cavities subjected to impurities doping. The flux trapping and its sensitivity
to rf losses are related to several intrinsic and extrinsic phenomena. To
elucidate the effect of re-crystallization by high temperature heat treatment
on the flux trapping sensitivity, we have fabricated two 1.3 GHz single cell
cavities from cold-worked Nb sheets and compared with cavities made from
standard fine-grain Nb. Flux expulsion ratio and flux trapping sensitivity
were measured after successive high temperature heat treatments. The cavity
made from cold worked Nb showed better flux expulsion after 800 ∘C/3h heat
treatment and similar behavior when heat treated with additional 900 ∘C/3h and
1000 ∘C/3h. In this contribution, we present the summary of flux expulsion,
trapping sensitivity, and RF results.

##  1 Introduction

Niobium has been the material of choice for superconducting radio-frequency
(SRF) cavities not only because of low power loss at the inner surface of the
cavities’ inner wall but also its high ductility which makes easier to
fabricate the complex structures [1].The niobium is elemental superconductor
with highest critical temperature, $T_{c}\sim$ 9.25 K and highest critical
field, $H_{c}\sim$ 200 mT. The performance is measured in terms of the quality
factor ($Q_{0}$) which defined as the ratio of stored energy inside the
cavities to the power dissipation on the inner wall of the cavities per radio
frequency (RF) cycle as a function of accelerating gradient ($E_{acc}$). The
ambient magnetic flux trapping during the cooldown is one of the prominent
factors causing the degradation of quality factor in cavities. The trapped
flux in the form of vortices oscillates in the presence of RF field and
dissipate energy. The field dependence of RF losses due to trapped vortices is
much stronger than the ohmic-type loss [2, 3]. The ambient flux trapping and
the flux trapping sensitivity to rf losses are related to several extrinsic
and intrinsic phenomena. The primary host sites of flux trapping are the
materials defects, dislocations, impurities, normal conducting precipitates.
For instance, the RF loss due to flux trapping can be minimized by maximizing
the flux expulsion when the cavity transition to the superconducting state
from normal conducting during the cavity cool-down by creating a large thermal
gradient across the cavity surface. Intrinsically, we can minimize the flux
trapping by minimizing the defects, dislocations, impurities with different
temperature heat treatments followed by chemical and mechanical polishing and
by high pressure rinse with de-ionized water. It has been demonstrated that
several different pinning mechanism plays a role to the rf losses due to
vortices [4]. Studies showed that doped cavities are more vulnerable to the
vortex dissipation loss due to the presence of the dopant on cavities rf
surface [5, 6, 7]. The flux expulsion can be maximized by increasing the
annealing temperature [8]. The increase in annealing temperature minimizes the
pinning centers by removing the clusters of dislocations and impurities. In
addition, the metallurgical state with larger grain size is expected as the
annealing temperature is increased. Fine-grain recrystallized microstructure
with an average grain size of 10–50 $\mu$m leads to flux trapping even with a
lack of dislocation structures in grain interiors [9]. Thus, it is important
to consider the crystallize structure of the niobium before the fabrications
and during the cavity processing [10]. In this contribution, we have
fabricated two single cell cavity, one from cavity grade SRF Nb with grain
size specified to ASTM 4-6 and other from the cold worked sheet with no
specified grain size. The cavity were processed together for chemical
polishing with electropolishing and successive annealing at 800, 900 and 1000
∘C/3h heat treatment. The flux expulsion ratio, flux trapping sensitivity and
$Q_{0}(B_{p})$ at 2.0 K were measured.

##  2 Fabrication and Surface Preparation

The SRF grade and cold-worked Nb sheets with residual resistivity ratio $\geq$
300 were purchased from Ningxia OTIC, China and 1.3 GHz TESLA shaped cavities
were fabricated at Zanon Research & Innovation Srl, Italy using standard
practice of deep drawing to half cells, trimming, machining of the iris and
equator of the half-cells, electron beam welding of the beam tubes (made from
low purity niobium). The cavity labeled TE1-05 was made from SRF grade Nb and
TE1-06 was made from cold-worked sheet. The cavities’ inner surface of $\sim$
150 $\mu m$ was removed in horizontal electropolishing setup using a mixture
of electronic grade $HF:H_{2}SO_{4}$ = 1:9 at a constant voltage of $\sim 14$
V, temperature of 15 - 20 ∘C and a speed of 1 rpm at Jefferson Lab.

Additional $\sim$25 $\mu$m inner surface was removed as a final surface
preparation using electropolishing after each successive heat treatments at
800 ∘C, 900 ∘C and 1000 ∘C for 3 hours.

##  3 Experimental setup

After the clean room assembly and leak check done, the cavity was transferred
to the vertical staging area. Three flux-gate magnetometers (FGM) were
fastened at the equator 120∘ apart each and parallel to the axis of the cavity
and six temperature sensors: two sensors at top iris, two sensors at equator,
and two sensors at bottom iris were attached. All temperature sensors were
attached 180 degree apart as shown in Fig. 1.

Figure 1: Experimental set up with FGMs and temperature sensors.

The cavity was moved to the Dewar and the temperature and magnetic field was
monitored with Cryocon 18i and flux gate magnetometer, respectively during the
thermal cycle (cool down and warm up above 10 K). The thermal cycles were done
at constant residual magnetic field of $\sim$ 10 mG with the help of
compensation coil as shown in Fig. 2. Several thermal cycles were performed in
order to measure the flux expulsion ratio.

Figure 2: The temperature and magnetic field during the cavity cooldown.
Depending on the temperature difference between iris-iris, the jump in
magnetic field is changed.

After the flux expulsion measurement, the cavity is filled with liquid helium
with residual magnetic field 0 and 20 mG in Dewar with cooldown condition of
full flux trapping and minimizing flux expulsion by changing the thermal
gradient across the cavity surface (iris-iris). The RF measurements consists
of $Q_{0}(T)$ at low peak RF field $B_{p}\sim 15mT$ at $Q_{0}(B_{p})$ at 2.0
K.

##  4 Results and Discussions

###  4.1 Flux expulsion ratio

The flux expulsion ratio is defined as the ratio of magnetic field measured in
superconducting state ($B_{sc}$) to the magnetic field measured in normal
conducting state $B_{n}$ as;

Flux expulsion ratio =$\frac{B_{sc}}{B_{n}}$

The flux expulsion ratio of 1 corresponds to the full flux trapping and ratio
of $\sim 1.7$ measured at outer surface of cavity’s equator corresponds to the
full flux expulsion. The ratio 1.7 calculated by COMSOL simulation assuming a
complete Meissner state [4].

Figure 3 shows the flux expulsion ratio ($B_{sc}/B_{n}$) as a function of the
temperature gradient (dT/ds) along the cavity surface. The cavity made from
cold work Nb showed better flux expulsion compared to cavity made from SRF
grade Nb after 800 ∘C/3h heat treatment. After additional annealing with 900
∘C/3h (not shown here) and 1000 ∘C/3h, the flux expulsion ratio for both
cavities is similar. The flux expulsion ratio quickly approaches to full
expulsion limit for dT/ds > 0.05, which corresponds the temperature difference
between the cavity irises to be $\sim$ 1 K.

Figure 3: Flux expulsion ratio as a function of temperature gradient across
the cavity surface (iris-iris) for cavity under study after 800 ∘C/3h and
additional 900 ∘C/3h and 1000 ∘C/3h annealing.

###  4.2 RF Results

Figure 4 shows the $R_{s}$ as a function of (1/T) for cavity TE1-05 and TE1-06
for two different cooldown condition with different residual magnetic field in
Dewar. The cavity was cooldown through $T_{c}$ with temperature gradient dT/ds
> 0.2 when the residual field in Dewar is $\sim$ 0 mG, where as the cooldown
through $T_{C}$ is kept dT/ds $\sim$ 0 when the residual magnetic field is
$\sim$ 20 mG. The $R_{s}(T)$ data were fitted using the model developed in
Ref. [11] to extract the temperature dependent $R_{BCS}$ and temperature
independent residual resistance $R_{i}$. The flux trapping sensitivity is
calculated using,

$S=\frac{R_{i,B_{2}}-R_{i,B_{1}}}{B_{2}-B_{1}}.$ (1)

where $B_{1}\approx$ 0 mG and $B_{2}\approx$ 20 mG.

Figure 4: $R_{s}$ vs. 1/T for cavity (a) TE1-05 and (b) TE1-06 after different
annealing and cooldown conditions. The solid lines are fits of the $R_{s}(T)$
using the model described in Ref. [11].

Figure 5 shows the $Q_{0}$ as a function of $E_{acc}$ at 2.0 K for both
cavities with different annealing and cooldown conditions. As expected,
$Q_{0}(E_{acc})$ shows the high field Q-slope starting $\sim$ 26 MV/m. All
test RF were limited by quench. Both cavities showed similar performance. The
cavity TE1-06 showed some multipacting $\sim$ 20-25 MV/m after 800 and 900 ∘
heat treatments. However, no multipacting were observed after 1000 ∘ heat
treatments. It is to be noted that the final surface preparation for all RF
tests were $\sim 25\leavevmode\nobreak\ \mu m$ electropolishing.

Figure 5: $Q_{0}$ vs. $E_{acc}$ at 2.0 K for cavity (a) TE1-05 and (b) TE1-06
after different annealing and cooldown conditions. All RF tests were limited
by quench. The error in $Q_{0}$ and $E_{acc}$ is < 10 % and < 5 %
respectively.

Figure 6 shows the flux trapping sensitivity as a function of cumulative
annealing temperature. The flux trapping sensitivity decreases with the
increase in annealing temperature. The largest change in sensitivity was
observed for cavity TE1-06. The flux trapping sensitivity for TE1-06 is higher
than TE1-05 after initial 800 $\circ$C/3h heat treatment. After additional 900
and 1000 $\circ$C/3h heat treatment, the sensitivity is very similar.

Figure 6: Flux trapping sensitivity as a function of annealing conditions. The
final surface preparation is electropolishing for all tests.

##  5 Discussion

As shown in Fig. 3, the flux expulsion for cavity made from cold work Nb
showed a better expulsion when the cavity was heat treated at 800 ∘C/3h
compared to SRF grade Nb. After additional 900 ∘C/3h and 1000 ∘C/3h heat
treatment, the initial microstructure state of the cavity sheet no longer
affects the flux expulsion performance [12]. On the other hand, the flux
trapping sensitivity depends on the final surface preparation of the cavities
before the rf test [13]. The poor flux expulsion will result in higher flux
sensitivity if the final surface preparation is other than electropolishing.
When the RF surface is modified with impurity doping, the flux trapping
sensitivity increases [4]. Thus, efficient flux expulsion is required to
minimize the additional rf loss related to the trapped flux. It is possible to
achieve the full flux expulsion limit on cavities fabricated with cold work Nb
after 800 ∘C/3h heat treatment if one can maintain the temperature gradient
between the cavities irises > 0.25 K/cm, which is not the case with cavities
made from SRF grade Nb. However, the cavity treated at 800 ∘C/3h requires
larger temperature gradient (dT/ds > 0.2 K/cm) to achieve the flux expulsion
ratio of $B_{sc}/B_{n}$ > 1.55, whereas with an additional 900 ∘C/3h and 1000
∘C/3h heat treatment the same expulsion ratio can be achieve for dT/ds > 0.5
K/cm. The reason that larger temperature gradient required to achieve the
maximum flux expulsion in SRF cavities showed that the recrystallization
process could be different on deformed structures vs. the sample coupons.
Additionally, significant number of magnetic flux pinning centers that are
present in the micro-structure pinned the flux weakly and need higher thermal
gradient to expel them away from the pinning centers. In fact, recent
recrystallization experiments on 30 % cold rolled Nb sheets showed a strain
path dependence on recrystallization [14].

To achieve the high flux expulsion ratio and hence the high quality factor in
nitrogen doped SRF cavities, the cavities used in LCLS-II project were heat
treated as high as 975 ∘C [15] for some Nb sheet from Ningxia. The initial
grain size specified with ASTM $\sim$ 5 showed better flux expulsion when
treated at 900 ∘C/3h compared to ASTM $\sim$ 7 requiring 975 ∘C/3h to achieve
the same level of flux expulsion. The electron back scattered imaging on
sample coupons showed that the full recrystallization and grain growth occur
after 800 ∘C/3h with average grain size of 141$\pm 6\mu$m. However, some
smaller grains of the order of 30 $\mu$m were observed along the cross-section
of the sample within 3 mm sheet thickness. The sample subjected 900 ∘C/3h has
average grain size of 211$\pm 13\leavevmode\nobreak\ \mu$m with an absence of
smaller grains after the 900 ∘C/3h [16]. It was previously demonstrated that
flux exit is harder in the regions of fine grins ($<50\leavevmode\nobreak\
\mu$m) than the regions of larger grains ($\geq$ 100 $\mu$m) [9].

The higher flux trapping sensitivity in TE1-06 cavity compared to TE1-05 after
initial 800 ∘C/3h, even with better flux expulsion ratio suggest the presence
of higher pinning center within the RF penetration depth. This results is in
consistent with the presence of smaller grains size on sample coupons after
800 ∘C/3h [16]. Furthermore, it confirms that the flux expulsion depends on
the bulk micro-structure, whereas flux trapping sensitivity depends on the
details of surface preparations within the RF penetration depth.

##  6 Summary

There is a correlation of flux expulsion ratio and flux trapping sensitivity
with the grain size, recrystallization as a result of high temperature heat
treatment. The recrystallization study on the cut out of cavity half cell to
understand the recrystallization with respect to the strain applied during
cavity fabrications. New cavities are being fabricated with know percentage of
cold work to further understand the controlled recrystalization to the flux
expulsion and flux trapping sensitivity.

##  7 ACKNOWLEDGEMENTS

We would like to acknowledge Dr. Gianluigi Ciovati and Dr. Shreyas
Balachandran for several discussion. We would like to acknowledge Jefferson
Lab technical staff members for the cavity fabrication, processing, and
cryogenic support during RF test.

## References

  * [1] H. Padamsee, " 50 years of success for SRF accelerators-a review", Supercond. Sci. Technol., 30, 053003 (2017).
  * [2] G. Ciovati, and A. Gurevich, "Evidence of high-field radio-frequency hot spots due to trapped vortices in niobium cavities", Phys. Rev. ST Accel. Beams 11, 122001 (2008).
  * [3] B. D. Khanal, P. Dhakal, and G. Ciovati, "Quench Detection in a Superconducting Radio Frequency Cavity with Combined Temperature and Magnetic Field Mapping", this workshop, paper: SUSPB016.
  * [4] P. Dhakal, G. Ciovati, and A. Gurevich, "Flux expulsion in niobium superconducting radiofrequency cavities of different purity and essential contributions to the flux sensitivity", Phys. Rev. ST Accel. Beams, 23, 023102 (2020).
  * [5] S. Posen et al., "Role of magnetic flux expulsion to reach $Q_{0}>3\times 10^{10}$in superconducting rf cryomodule", Phys. Rev. Accel.Beams, 22, 032001 (2019).
  * [6] D. Gonnella, J. Kaufman, and M. Liepe, "Impact of nitrogen doping of niobium superconducting cavities on the sensitivity of surface resistance to trapped magnetic flux", J.Appl. Phys. 119, 073904 (2016).
  * [7] M. Martinello et al., "Effect of interstitial impurities on the field dependent microwave surface resistance of niobium", Appl.Phys. Lett. 109, 062601 (2016).
  * [8] S. Posen et al., "Efficient expulsion of magnetic flux in superconducting radiofrequency cavities for high $Q_{0}$ applications", J. Appl. Phys. 119, 213903 (2016).
  * [9] S. Balachndran et al., "Direct evidence of microstructure dependence of magnetic flux trapping in niobium", Sci. Rep. 11, 5364 (2021).
  * [10] C. Z. Antoine, "Influence of crystalline structure on rf dissipation in superconducting niobium" Phys. Rev. Accel. Beams 22, 034801 (2019).
  * [11] G. Ciovati, P. Dhakal, and A. Gurevich, "Decrease of the surface resistance in superconducting niobium resonator cavities by the microwave field", Appl. Phys. Lett. 104 (9), 092601 (2014).
  * [12] B. D. Khanal, S. Balachandran, P. J. Lee, S Chetri, and P. Dhakal,"Magnetic flux expulsion in superconducting radio-frequency niobium cavities made from cold worked niobium",  in the Proc. of NAPAC22, Albaquarque, New MExico, USA.
  * [13] B. D. Khanal and P. Dhakal, "Insight to the Duration of 120 ∘C Baking on the Performance of SRF Niobium Cavities," IEEE Transactions on Applied Superconductivity 33 3, 1-6 (2023).
  * [14] Z. L. Thune, C. McKinney, N. G. Fleming, and T. R. Bieler, "The influence of strain path and heat treatment variations on recrystallization in cold-rolled high-purity niobium polycrystals" IEEE Transactions on Applied Superconductivity, 33(5), 1-4 (2023).
  * [15] A. D. Palczewski, D. Gonnella, O. S. Melnychuk, and D. A. Sergatskov, “Study of Flux Trapping Variability between Batches of Tokyo Denkai Niobium used for the LCLS-II Project and Subsequent 9-cell RF Loss Distribution between the Batches”, in Proc. SRF’19, Dresden, Germany, Jun.-Jul. 2019, pp. 570-575.
  * [16] S. Balachandran et al., "Microstructure development in a cold worked SRF niobium sheet after heat treatments", this workshop, paper ID: MOPMB041.