# Cross-Over between universality classes in a magnetically disordered
metallic wire.

Guillaume Paulin  Institut für Theoretische Physik, Universität zu K öln, Z
ülpicher Str. 77, 50937 K öln, Deutschland  David Carpentier  CNRS -
Laboratoire de Physique de l’Ecole Normale Supérieure de Lyon,
46, Allée d’Italie, 69007 Lyon, France

###### Abstract

In this article we present numerical results of conduction in a disordered
quasi-1D wire in the possible presence of magnetic impurities. Our analysis
leads us to the study of universal properties in different conduction regimes
such as the localized and metallic ones. In particular, we analyse the cross-
over between universality classes occuring when the strength of magnetic
disorder is increased. For this purpose, we use a numerical Landauer approach,
and derive the scattering matrix of the wire from electron’s Green’s function

##  1 Introduction

Interplay between disorder and quantum interferences leads to one of the most
remarkable phenomenon in condensed matter : the Anderson localization of
waves. The possibility to probe directly the properties of this localization
with cold atoms[1, 2] have greatly renewed the interest on this fascinating
physics. In this paper, we focus on the particular situation where electrons
encounter two kinds of disorder: a usual scalar potential at the origin of
diffusion, and a magnetic potential, arising from a collection of frozen
random magnetic moments. This situation is naturally realized experimentally
in the study of transport properties of metallic spin glass wires [3, 4, 5,
6]. In these wires, the spins freeze at low temperatures when entering the
spin glass phase due to the frustrating magnetic couplings. In this glassy
phase, and neglecting any residual Kondo effect in this regime, the impurities
act effectively as a (weak) magnetic potential. We study numerically the
effect of both types of disorder on the statistical properties of the wire
conductance. In particular, we will focus on the experimentally relevant
crossover of (weak) localization properties of the wire as a function of the
magnetic disorder strength.

One dimensional disordered electronic systems are always localized. Following
the scaling theory [7] this implies that by increasing the length $L_{x}$ of
the wire for a fixed amplitude of disorder, its typical conductance ultimately
reaches vanishingly small values. The localization length $\xi$ separates
metallic regime for small length $L_{x}\ll\xi$ from the asymptotic insulating
regime. In the present paper, we focus on several universal properties of both
metallic and insulating regime of these wires in the simultaneous presence of
two kinds of disorder. The first type corresponds to scalar potentials induced
by the impurities, for which the system has time reversal symmetry (TRS) and
spin rotation degeneracy. In this class the Hamiltonian belongs to the so-
called Gaussian Orthogonal Ensemble (GOE) of the Random Matrix Theory
classification [8] (RMT), corresponding to the class AI in the modern
classification of Anderson universality classes (see e.g. [9]). If impurities
do have a spin, the TRS is broken as well as spin rotation invariance. The
Hamiltonian is then a unitary matrix, which corresponds in RMT to the Gaussian
Unitary Ensemble (GUE) with the breaking of Kramers degeneracy [10], and to
the Anderson class A [9]. However, for the experimentally relevant case of a
magnetic potential weaker than the scalar potential, the system is neither
described by the GUE class, nor by the GOE class, but extrapolates in between.
This intermediate regime, of particular relevance experimentally, is the main
object of study of the present paper. Moreover the present work paved the way
towards a numerical study of the correlation of conductances in the cross over
regime [11].
This paper is organized as follows : in section 2, the model and the numerical
method used will be described in details. In section 3 we identify the
localized and metallic regime of transport of the system. The localization
length will be determined by two different methods and the cross over between
both universality classes (GOE, GUE) will be highlighted. In section 4, the
insulating regime is studied with a particular focus on the statistical
distribution of the conductance, which allows us to highlight universal
behavior. In section 5, we turn to the study of the metallic regime, and
perform a careful analysis of each of the first three cumulants of the
statistical distribution of the conductance. We focus on the universal
properties of conductance fluctuations, and the non-analyticity of the
complete distribution is discussed. Finally section 6 is devoted to the
conclusion.

##  2 The model and the method

###  2.1 The model

In this paper, we study numerically the scaling of transport properties of
wires in the presence of magnetic and scalar disorders. We will focus on the
regime of phase coherent transport, reached experimentally at low temperature
(see in particular [6]). In this regime, the phase coherence length $L_{\phi}$
which phenomenologically accounts for inelastic scattering of electrons on
impurities [12] is larger than (or comparable with) the wire’s length $L_{x}$,
so that phase coherence for the propagating electrons is preserved in the
whole sample. Note that this phase coherence length $L_{\phi}$ includes in
particular a contribution from inelastic scattering on the non frozen magnetic
impurities through a Kondo dephasing, which is strongly reduced when entering
the magnetic glass phase [6].

We describe the behavior of electrons inside the disordered wire using a
tight-binding Anderson lattice model with two kinds of disorder potentials:

$\displaystyle\mathcal{H}=\sum_{<i,j>,s}t_{ij}c_{j,s}^{\dagger}c_{i,s}$
$\displaystyle+$
$\displaystyle\sum_{i,s}v_{i}c_{i,s}^{\dagger}c_{i,s}+J\sum_{i,s,s^{\prime}}{\vec{S}}_{i}.{\vec{\sigma}}_{s,s^{\prime}}c_{i,s}^{\dagger}c_{i,s^{\prime}}.$
(1)

$t_{ij}$ is the hopping term from site $i$ to $j$. In the following, $t_{ij}$
will take two different values: $t_{ij}=t_{//}$ in the longitudinal $x$
direction and $t_{ij}=t_{\perp}$ in the transverse $y$ direction. The scalar
disorder potential $V=\\{v_{i}\\}_{i}$ is diagonal in electron-spin space. We
choose the $v_{i}$ to be random scalars uniformly distributed in the interval
$[-W/2,W/2]$. In this work, we have chosen without loss of generality to fix
$t_{//}=1$ so that all energy scales are relative to the bandwidth $t_{//}=1$,
and the amplitude of disorder $W=0.6$. In eq. (1) $s,s^{\prime}$ label the
$SU(2)$ spin of electrons and the $\vec{S}_{i}$ account for spins of the
frozen magnetic impurities. A realistic choice for these frozen spins in a
spin glass phase consist in considering classical spins with random
orientations, thereby neglecting any small spatial correlation in a spin glass
configuration [13]. The coupling $J$ between the electron spins and the
magnetic impurities fixes the amplitude of the magnetic disorder. In this
work, we will monitor the behavior of the transport properties of the sample
as a function of this amplitude $J$, from $J=0$ to $J=0.4$. Indeed, variation
of the amplitude of magnetic disorder $J$ allows to extrapolate from GOE /
class AI for $J=0$ to GUE / class A for $J\neq 0$. We will also show in
section 5 that, as a bonus, the presence of the magnetic disorder allows also
an numerically easier settlement of the universal metallic regime of the wire.

For a given realization of both scalar $\\{v_{i}\\}_{i}$ and magnetic disorder
$\\{\vec{S}_{i}\\}_{i}$, the Landauer conductance of this model on a 2D square
lattice of size $L_{x}\times L_{y}$ is evaluated numerically using a recursive
Landaueur method described in details in the next paragraph.

###  2.2 The method

We have chosen to use a numerical method based on the lattice model (1) as
opposed to e.g. random matrix or Dorokhov-Mello-Pereyra-Kumar (DMPK) [14, 15]
method so as to provide a numerical study allowing for an easy comparison with
the experimental situation[6]. Moreover, this method allows for further
natural developments such as the study of the conductance change upon magnetic
impurities spin flipping, which would be difficult to reach by alternative
method. Starting from a lattice model such as (1), the natural method
providing the conductance of a finite size sample is based on the Landauer
formalism [16].

We consider a two-terminal setup, with electrodes connected to the wire at
$x=0$ and at $x=L_{x}$. These electrodes are described as semi-infinite
ribbons with the same transverse geometry as the sample, and described with
(1) but without randomness. Electrons are then confined in the transverse
$y$-direction via a potential that has the form:

$\displaystyle V(y)$ $\displaystyle=$ $\displaystyle 0\hskip
14.22636pt\mathrm{if}\hskip 14.22636pt0\leq y\leq L_{y}\ ;\ V(y)=\infty\hskip
14.22636pt\mathrm{otherwise}$ (2)

This confining potential in the $y$ direction leads to the appearance in the
electrodes of $N$ modes propagating in the $x$-direction. The complete wave
function of an electron in this tight binding lattice model reads then:

$\psi(x,y)=\phi_{n}(y)e^{\imath k_{x}x},$ (3)

where $k_{x}$ is the momentum of electrons in the longitudinal direction and

$\phi_{n}(y)=\sqrt{\frac{2}{N_{y}+1}}\sin\left(\frac{n\pi y}{N_{y}+1}\right).$
(4)

We used $N_{y}=L_{y}$ in units of lattice spacing. The group velocity of this
mode reads:

$v_{n}=2\frac{t_{//}}{\hbar}\sin\left(k_{x}\right).$ (5)

where $t_{//}$ is the longitudinal hopping amplitude. This velocity depends on
the momentum $k_{x}$ which is determined for a constant energy by the
dispersion relation

$E_{n}=\mu-2t_{\perp}\cos\left(\frac{n\pi}{N_{y}+1}\right)-2t_{//}\cos(k_{x}).$
(6)

At given energy $E-\mu$, we end up with the following relation for the
longitudinal part of the momentum of the electron:

$k_{x}=\arccos\left(\frac{\mu-E}{2t_{//}}-\frac{t_{\perp}}{t_{//}}\cos\left(\frac{n\pi}{N_{y}+1}\right)\right).$
(7)

To optimize the efficiency of the numerical study, we fix $t_{//}=2t_{\perp}$
and stay near the band center, avoiding in particular the presence of
fluctuating states studied in [17].

To compute the conductance of such a wire, the Landauer-Büttiker formalism of
coherent transport is used [18]. It allows to relate the dimensionless
conductance $g$ of a diffusive wire with the scattering matrix $T$:

$g=\sum_{\mathrm{modes}\,m,n}T_{mn},$ (8)

where $T_{mn}$ is the transmission coefficient between modes $m$ and $n$. The
dimensionless conductance $g$ is defined from the conductance $G$ of the
system as $g=G/(2e^{2}/h)$ when spin degeneracy is present ($J=0$), and
$g=G/(e^{2}/h)$ otherwise ($J\neq 0$).

Following Fisher and Lee[19] we relate the scattering matrix to the electronic
retarded Green’s functions of the system through:

$\displaystyle t_{mn}$ $\displaystyle=$
$\displaystyle\imath\hbar\sqrt{v_{n}v_{m}}\sum_{y,y^{\prime}=1}^{N_{y}}\phi_{n}(y)\mathcal{G}^{R}(y,x=0|y^{\prime},x=L_{x})\phi_{m}(y^{\prime}).$
(9)

with $T=Tr(t^{\dagger}t)$, $v_{n}$ and $\phi_{n}$ (resp. $v_{m}$ and
$\phi_{m}$) are the group velocity and the eigen wave function of propagating
mode $n$ (resp. $m$). Mode $n$ belongs to the left lead whereas mode $m$
belongs to the right one. In (9) $\mathcal{G}^{R}(y,x=0|y^{\prime},x=L_{x})$
represents the retarded Green’s function of an electron between the point
$(x=0,y)$ in the left contact between the system and the electrode and the
point $(x=L_{x},y^{\prime})$ in the right contact. Consequently, by
calculating the retarded Green’s functions of the system only between its both
sides, we are able to determine the dimensionless conductance of the wire. An
efficient method to calculate $\mathcal{G}^{R}(y,x=0|y^{\prime},x=L_{x})$,
which takes advantage of the quasi-one dimensional nature of the ribbon
consists in obtaining it recursively[20]. In figure 1 the principle of the
method is sketched: using a Dyson equation we deduce the boundary Green’s
function of a system of size $n+1$ from the corresponding Green’s function of
a subsystem of size $n$, and the exact Green’s function of the $n+1$ row. This
allows to perform matrix inversion only of the simple row system. At both
initial and final steps, we reconnect the system to the semi-infinite
electrodes (see fig. 2) described by a standard self energy:

$\mathcal{G}^{R}_{\mathrm{bound}}(y_{1},y_{2})=-\frac{1}{t_{//}}\sum_{n=1}^{N_{y}}\phi_{n}(y_{1})e^{\imath
k_{x}}\phi_{n}(y_{2}).$ (10)

The last step consists in combining the Green’s functions of the sample of
desired longitudinal size with the one of the right lead, which is also given
by equation (10).

Figure 1:  Principle of recursive calculation of retarded Green’s functions of
the wire. Use of a Dyson: $G^{R}_{n+1}=G^{R}_{n}+G^{R}_{n}VG^{R}_{n+1}$. At
each step the longitudinal length $L_{x}$ is increased by one lattice spacing.

Figure 2:  Boundary conditions: the wire is connected to two leads represented
by two semi infinite metallic wires.

The unitarity of the corresponding scattering matrix $\\{T_{mn}\\}_{m,n}$ is
used to monitor the accuracy of the numerical method. Such test yielded
typical relative error of order $10^{-4}$ for a system size $L_{x}=1600$ and
$L_{y}=40$. This method allows one to compute the conductance of a wire of
length $L_{x}$ and of width $L_{y}$ for any given configuration of scalar
disorder $V$ and for any configuration of frozen classical spins
$\\{\vec{S}_{i}\\}_{i}$.

In the next sections, we study the properties of this conductance for one
given configuration of magnetic disorder but for many different configurations
of scalar disorder. Universal properties are identified by varying the
transverse length $L_{y}$ from $10$ to $80$, with the aspect ratio
$L_{x}/L_{y}$ taken from $1$ to $6000$. Typical number of configurations of
scalar disorder $V$ we used were $N_{d}=5000$, with exceptions for the study
of the localization properties where for $L_{y}=10$ we sampled the conductance
distribution for $50000$ different configurations of disorder.

##  3 Localization length

###  3.1 Determination of $\xi$

We start our analysis by a determination of the parameters corresponding to
the metallic and localized regimes, through a careful determination of the
localization length of the system. While experimentally the only accessible
regime of a phase coherent wire is the metallic regime, numerically this
regime is difficult to reach and describe quantitatively, as opposed to the
localized regime. This is related to the extreme reduction in the number of
propagating modes in the numerical system which is associated with a
corresponding reduction of the localization length. Hence in order to clearly
identify the conditions to access universal properties of the transport in the
metallic regime, we start by a detailed determination of this localization
length in the experimentally relevant crossover situation. Afterwards we will
take the opportunity of the present study to describe other characteristics of
the localized regimes in the crossover situation, before turning to our main
interest : the universal metallic regime.

The localization length separates short wires of metallic behavior from a
insulating long wires. A first method to access to the localization length
$\xi$ from the conductance consists in considering the scaling behavior of the
typical conductance $g_{typ}$ defined as:

$g_{typ}=e^{\langle\log g\rangle},$ (11)

where $\langle\cdot\rangle$ represents the average over the different
configurations of scalar disorder $V$. This typical conductance decays
exponentially with the longitudinal length of the wire[21, 22] :

$g_{typ}\sim e^{-\frac{2L_{x}}{\xi}}$ (12)

in the regime of long wires $L_{x}\gg\xi$.

Figure 3:  Evolution of $\langle\log g\rangle$ as a function of longitudinal
size for different transverse length ($L_{y}=10,20,30,40$) and $J=0$. The
linear part of the curve allows one to get the localization length $\xi(J=0)$
from the scaling form in the insulating regime $\langle\log
g\rangle=-\frac{2L_{x}}{\xi}$.

Figure 4:  Evolution of $\langle\log g\rangle$ as a function of longitudinal
size for different transverse length ($L_{y}=10,20,30,40$) and $J=0.05$. The
linear part of the curve allows one to get the localization length
$\xi(J=0.05)$ from the scaling form in the insulating regime $\langle\log
g\rangle=-\frac{2L_{x}}{\xi}$.

Figure 5:  Evolution of $\langle\log g\rangle$ as a function of longitudinal
size for different transverse length ($L_{y}=10,20,30,40$) and $J=0.2$. The
linear part of the curve allows one to get the localization length
$\xi(J=0.2)$ from the scaling form in the insulating regime $\langle\log
g\rangle=-\frac{2L_{x}}{\xi}$.

Figure 6:  Evolution of $\langle\log g\rangle$ as a function of longitudinal
size for different transverse length ($L_{y}=10,20,30,40$) and $J=0.4$. The
linear part of the curve allows one to get the localization length
$\xi(J=0.4)$ from the scaling form in the insulating regime $\langle\log
g\rangle=-\frac{2L_{x}}{\xi}$.

Figures 3, 4, 5 and 6 show the behavior of $\langle\log g\rangle$ as a
function of longitudinal length $L_{x}$ for different widths $L_{y}$.
Different curves correspond to different values of magnetic disorder. The
linear fit of the large length part of the curves allow for a precise
determination of the corresponding localization length for each value of
$L_{y}$ and $J$.
A second method to determine this localization length consists in considering
the Lyapunov exponent $\gamma$ of the transfer matrix of the system, following
a standard random matrix theory approach[23, 21, 24]. This exponent can be
deduced from the conductance as

$\gamma(L_{x})=\frac{1}{2L_{x}}\log\left(1+\frac{1}{g(L_{x})}\right),$ (13)

and the localization length follows from its asymptotic behavior :

$\xi^{-1}=\lim_{L_{x}\to\infty}\gamma(L_{x}).$ (14)

On figure 7 and 8 we have plotted the Lyapunov exponent versus the inverse of
the longitudinal length for different values of magnetic disorder. Different
curves correspond to different widths of the wire.

Figure 7:  Evolution of Lyapunov exponent with the inverse of longitudinal
length in semi-log plot. Circles correspond to $J=0$, squares to $J=0.05$,
diamonds to $J=0.2$ and triangles to $J=0.4$. The value of the transverse
length is $L_{y}=10$. The localization length can be extrapolated from the
value of $\gamma$ for $L_{x}\to\infty$.

Figure 8:  Evolution of Lyapunov exponent with the inverse of longitudinal
length in semi-log plot. Circles correspond to $J=0$, squares to $J=0.05$,
diamonds to $J=0.2$ and triangles to $J=0.4$. The value of the transverse
length is $L_{y}=20$. The localization length can be extrapolated from the
value of $\gamma$ for $L_{x}\to\infty$.

With this method, a simple extrapolation of the curve is necessary to obtain
$\xi$, without any fit. the value of the inverse of the localization length
for an infinite wire. Nevertheless, we have found that this method shows less
accuracy than the preceding one: as shown on figures 7 or 8, the Lyapunov
exponent is still varying for the longest longitudinal length. We find that
both methods give fully compatible results while the Lyapunov exponent method
requires much larger system sizes than the typical conductance method for a
given required accuracy.

A first manifestation of the universality of the Anderson localization classes
appears through the dependance of $\xi$ on the transverse length $L_{y}$ (or
the number of propagating modes). It is expected to follow[21, 25]:

$\xi=(\beta L_{y}+2-\beta)l_{e},$ (15)

with $l_{e}$ the elastic mean free path and $\beta$ encodes the universal
class of the model : $\beta=1$ corresponds to the orthogonal universality
class GOE while $\beta=2$ for GUE. Note that this change in $\beta$ is
accompanied by an artificial doubling of the number of transverse modes
$N_{y}\equiv L_{y}\to 2N_{y}$ due to the breaking of Kramers degeneracy [21].
This effective factor 4 when breaking the spin rotation symmetry has been
discussed in ref. [26, 27] in details, when discussing the magnetic field
dependance of this localization length, in comparison with random matrix and
Non-linear sigma models. Comparison of numerical localization lengths for
different $J$ with (15) is shown in fig 9. Excellent agreement is found for
$J=0$ (GOE class, $\beta=1$) and with the GUE class for $J\geq 0.2$. For
intermediate values of $J\neq 0$ we observe a crossover between the two
extreme GOE and GUE laws, which cannot be described by eq. (15), and for which
no analytical work exists to our knowledge.

From these results, we also notice that the localization regime is reached for
much longer wires in the presence of magnetic impurities (GUE case) than
without (GOE): localization is hampered by the presence of these magnetic
impurities. As we will discuss below, this property helps in observing
numerically the universal weak localization regime and the associated
universal conductance fluctuations.

Figure 9:  Evolution of localization length as a function of transverse
length. $l_{e}$ is the mean free path of the diffusive sample. Different
behavior of the localization length if $J=0$ or $J\neq 0$. Inset : Scaling of
the typical conductance $\langle\log g\rangle=-\frac{2L_{x}}{\xi}$.

###  3.2 The Insulating and Metallic regimes

The localization length discriminates between both insulating and metallic
regimes: the ribbon behaves indeed as a metal ($g\gg 1$) for lengths
$L_{x}\ll\xi$ and as an insulator ($g\ll 1$) if $L_{x}\gg\xi$. In both
asymptotic regimes the shape of the Probability Density Function (PDF) of the
conductance which is known: it is Log-normal for insulating wires[21] and
gaussian for metallic ones. By varying the longitudinal length we can study
the evolution of this PDF from a gaussian to a log-normal distribution, as
shown on figure 10. This plot is done for a given value of $L_{y}$ and $J$.
One can notice that in the metallic regime, the PDF is very well approximated
by a gaussian for relatively small wires : the gaussian regimes is easily
reached, whereas it takes length much larger than the localization length for
the distribution to become log-normal in the insulating regime. This point
will be discussed more precisely below on the cumulants of this PDF.

Figure 10:  Evolution of the statistical distribution of the conductance for
different longitudinal sizes for $L_{y}=10$ and $J=0$. $N$ is the number of
disorder configurations used. Plain lines are gaussian (if $\langle
g\rangle>1$) or log-normal (if $\langle g\rangle<1$) interpolations of
numerical data.

The insulating regime is then characterized by $\langle g\rangle<1$ and the
metallic one by $\langle g\rangle>1$.

In order to characterize samples by the average $\langle g\rangle$, and in
particular plots higher cumulants as a function of $\langle g\rangle$, we now
turn to a short study of the behavior of this first cumulant as a function of
the system size. On figure 11 we have plotted $\langle g\rangle(L_{x})$ for
different values of magnetic disorder $J$ (hence different localization
lengths). These curves approximately collapse when plotted against the scaling
variable $L_{x}/\xi$, as show on fig. 12. we remind the reader that for $J=0$,
the average conductance is defined by $g=G/2G_{0}$, which explains why $J=0$
and $J=0.4$ curves coincide on fig. 11: according to figure 9, only for these
values of magnetic disorder universality classes are reached.

Figure 11:  Evolution of average conductance versus longitudinal length for
different values of magnetic disorder for $L_{y}=10$ and $J=0,0.05,0.1,0.2$
and $0.4$. The line $\langle g\rangle=1$ is plotted as a frontier between
insulating and metallic regimes.

Figure 12:  Evolution of average conductance versus the scaling variable
$L_{x}/\xi$, for $L_{y}=10$ and $J=0,0.05,0.1,0.2$ and $0.4$. All curves
nearly collapse in one single curve. The frontier between insulating and
metallic is drawn again.

This study of $\langle g\rangle(L_{x})$ allows to proceed in the study of
higher cumulants of the PDF of $g$ and test prediction of the single parameter
scaling of distributions[28].

##  4 Universal Insulating Regime

###  4.1 Probability Density Functions

In the insulating regime $L_{x}\geq\xi$, we expect a Log-normal conductance
statistical distribution, as seen previously. However in the weakly insulating
regime $\langle g\rangle\lesssim 1$ this log-normal asymptotic form is not
reached. Instead, as shown in fig. 10 we find a non-analytical behavior of
$P(g)$ in agreement with [29, 30, 31, 32].

Figure 13: Comparison of Probability density functions (PDF) of conductance
for $J=0$ (plain curves) and $J=0.2$ (dashed curves). Plots are performed for
different values of average conductance. (a): $\langle g\rangle(J=0)=0.84$ and
$\langle g\rangle(J=0.2)=0.79$. (b): $\langle g\rangle(J=0)=0.67$ and $\langle
g\rangle(J=0.2)=0.62$. (c): $\langle g\rangle(J=0)=0.45$ and $\langle
g\rangle(J=0.2)=0.42$. (d): $\langle g\rangle(J=0)=0.21$ and $\langle
g\rangle(J=0.2)=0.18$.

Figure 14: PDF of conductance for $\langle g\rangle<1(J=0.2)$ and $\langle
g\rangle>1(J=0)$ and Gaussian interpolations. $L_{y}=10$.

In order to study the dependance of this non-analyticity on the universality
class, we have plotted on figure 13 the distribution $P(g)$ for similar values
of $\langle g\rangle$ but different magnetic strengths $J=0$ (GOE) and $J=0.2$
(GUE). The shapes of these distributions are highly similar if $\langle
g\rangle\ll 1$, showing that distributions for $J=0$ and $J\neq 0$ reaches the
same Log-normal distribution at large system sizes, in agreement with the
super-universality scenario[33]. In the intermediate regime ($\langle
g\rangle\approx 1$), shapes are symmetry dependent. Moreover we find that the
non-analyticity appears for different values of conductance (close to $1$) and
the rate of the exponential decay [29] in the metallic regime seems to differ
from one class to the other (see for instance curves (a) or (b)). Finally, the
plot on figure 14 represents the distribution of conductance of mean value
just above and below the threshold $\langle g\rangle=1$. Plain lines represent
gaussian interpolations with a mean and a variance given by the first and the
second cumulant of each numerical conductance distribution. For $\langle
g\rangle>1$, the gaussian interpolation approximates very well the full
distribution while as soon as $\langle g\rangle<1$, the gaussian approximation
only applies in the tail $g\geq 1$ of the distribution of conductance[32]. The
shape of the distribution for $g<1$ (Figure 13) agrees qualitatively with the
numerical results of [34] (see in particular their fig. 4), Unfortunately, a
more accurate comparison proves to be difficult due to the lack of analytical
description of the distribution. To quantize further these results on the
whole distribution of $\log g$ we now turn to a quantitative study of the
second and third cumulant of this distribution.

###  4.2 Study of cumulants

This conductance distribution converges to the Log-normal only deeply in the
insulating regime, the convergence being very slow (much slower than in the
metallic regime). This qualitative result is confirmed by the study of
moments: in the insulating regime the second cumulant is expected to
follow[35, 25]:

$\langle\left(\log g-\langle\log g\rangle\right)^{2}\rangle=\langle(\log
g)^{2}\rangle_{c}=-2\langle\log g\rangle,$ (16)

Our numerical results are in agreement with this scaling with however very
slow convergence towards this law: corrections are measurable even for the
largest system size where the system is deeply in the localized state, as
shown on fig. 15. More precisely, we find that for the deep insulating regime
$\langle(\log g)^{2}\rangle_{c}=-1.88\langle\log g\rangle$ slope $-1.88$, with
a slight discrepancy with (16).

Figure 15:  Plot of the variance of $\log g$ as a function of the mean for the
orthogonal ($J=0$) and unitary ($J=0.2$) case, for $L_{y}=10$ and $20$. The
slope of the linear fit is equal to $-1.88$. This plot shows also super-
universality as the behavior of the second cumulant does not depend neither on
geometry of the wire nor on the universality class.

This plot also shows that in the deep insulating regime the behavior of the
second cumulant as a function of the first one does not depend on the value of
magnetic disorder: both curves for $J=0$ and $J=0.2$ follow the same law. This
is in agreement with our previous result on statistical distributions: there
is a super-universal behavior in the deep insulating regime.

Finally we have studied the third cumulant of $\log g$ scaled in figure 16 as
a function of the first cumulant. The linear behavior for each value of
magnetic disorder in the deep insulating regime ($\langle\log g\rangle<-4$) is
in agreement with the single parameter scaling. We find that contrary to the
second cumulant the coefficient of proportionality between the skewness and
the average depends on the symmetry of disorder, which denotes a lack of
super-universality concerning this cumulant. For instance dots and diamonds
(which correspond to the case $J=0$) have the same behavior, as opposed to the
case $J=0.2$ (squares or triangles). A systematic study of this point with
even larger system sizes and other numerical methods would be of high interest
but is definitely beyond the scope of the present paper.

Figure 16:  Plot of the skewness of $\log g$ as a function of the mean for the
orthogonal ($J=0$) and unitary ($J=0.2$) case, for $L_{y}=10$ and $20$.

The study of cumulants of the distribution of $\log g$ confirm the single
parameter scaling of the distribution, with a slight discrepancy concerning
the value of the coefficient of proportionality between second and first
cumulant. Moreover, super-universality has been highlighted concerning the
second cumulant but is lacking concerning the third one.

###  4.3 Comparison with exact results in the cross-over regime

For localization in wires connected to ideal contacts, exact formula for the
average conductance [36] and conductance fluctuations [37] have been derived
for the two universal orthogonal and unitary classes. These formula are of
particular interest in the intermediate regime between the metallic
($L_{x}\ll\xi$) and the deeply localized ($L_{x}\gg\xi$) regime. We have
numerically evaluated the formula (3.105) and (3.106) of ref. [37] and we
compare them with our numerical data in Fig. 17 and Fig. 18. We find an
excellent agreement between the $J=0$ data and the orthogonal exact formulae
on one hand, and the $J\neq 0$ data and the unitary formulae on the other
hand. This comparison naturally breaks down for small sizes where the quasi-1d
assumption for diffusion breaks down and a non universal regime takes place.
As shown on Fig. 18, the exact unitary behavior is recovered for $J=0.05$
beyond the cross-over length $L_{m}$ (see next section). Whether the scale
dependance of $\langle g\rangle$ and $\langle\delta g^{2}\rangle$ in the
cross-over regime ($J$ small such that $L_{m}\simeq\xi$) is amenable to an
exact formula along the lines of ref. [37] is an open question of high
interest.

Figure 17:  Comparison of the numerical average conductance $\langle g\rangle$
as a function of the reduced length $L_{x}/\xi$ with the exact expressions of
ref. [37] for the orthogonal ($J=0$) and unitary ($J=0.2$) case, for various
transverse sizes $L_{y}$.

Figure 18:  Comparison of the numerical average conductance fluctuations
$\langle\delta g^{2}\rangle$ as a function of the reduced length $L_{x}/\xi$
with the exact expressions of ref. [37] for the orthogonal ($J=0$) and unitary
($J=0.4$) case, for various transverse sizes $L_{y}$.

##  5 Universal Metallic Regime

We now focus on the universal metallic regime described by weak localization.
By definition weak localization corresponds to metallic diffusion, expected
for lengths of wire $l_{e}\ll L_{x}\ll\xi$. For this regime to be reached, we
thus need to increase $\xi$ through an increase of the number of transverse
modes $L_{y}$ with all other parameters fixed (see eq. (15)). Moreover, for a
fixed geometry, this regime will be easier to reach with magnetic impurities
than without. As we saw on fig. 13, the shape of the PDF of conductance is a
truncated gaussian in this regime. In the following we study its three first
cumulants quantitatively, starting with the variance.

###  5.1 Conductance Fluctuations and Universal Crossover

In the weak localization regime, the conductance fluctuations $\langle
g^{2}\rangle_{c}$ are expected to be independent on the size of the system,
and only depend on the universal localization class of the model. The figure
19 shows that for a suitable value of transverse length, the system reaches
the expect plateau in conductance fluctuations. The value of the plateau
identifies with the expected values ($1/15$ and $4/15$ for GUE and GOE
respectively) with a high accuracy. However, the presence of this plateau
depends strongly on the value of transverse length $L_{y}$ and on the magnetic
disorder. Figures 20 and 21 illustrates this point with further details: these
plots show the conductance fluctuations as a function of longitudinal size for
$J=0$ and $J=0.2$ and for two values of transverse length $L_{y}$. On the
first plot, for both values of magnetic disorder the universal plateau arises,
whereas it appears only for $J=0.2$ if $L_{y}=40$.

The evolution of this variance of the PDF of $g$ depends on the longitudinal
length through [12]:

$\displaystyle\langle\delta g^{2}\rangle=\langle g^{2}\rangle_{c}$
$\displaystyle=$
$\displaystyle\frac{1}{4}F\left(0\right)+\frac{3}{4}F\left(x\sqrt{\frac{4}{3}}\right)+\frac{1}{4}F\left(x\sqrt{2}\right)+\frac{1}{4}F\left(x\sqrt{\frac{2}{3}}\right),$
(17)

where $x=L_{x}/L_{m}$ and the scaling function $F(x)$ depends only on
dimension[12, 11]. The universality occurs in this equation in the two limit
$x=0$ corresponding to $J=0$ or $x\gg 1$ where $F(x)\to 0$. In both cases the
variance becomes geometry independent. Moreover this expression shows that the
whole crossover between the two classes is described by a universal crossover
function parametrized solely by the length $L_{m}$, called the magnetic
dephasing length scale [12]. In the figure 22, these conductance fluctuations
are plotted as a function of longitudinal length $L_{x}$ for different values
of $J$. A single parameter fit by (17) provides the determination of the
magnetic dephasing length $L_{m}$ as a function of the magnetic disorder $J$.
The Universal behavior is also highlighted for strong enough magnetic
disorder. The determination of this scattering length allows a precise study
of average conductance, and in particular the study of the classical part as
described in the next sub section below. Fig. 23 shows the scaling form of
these fluctuations (as a function of $L_{x}/L_{m}(J)$) in excellent agreement
with the theory (17). Moreover, for long wires (and large values of $J$)
conductance fluctuations are no longer $L_{x}$ dependent and equal to $1/15$,
in units of $G_{0}^{2}$. This is the so-called Universal Conductance
Fluctuations (UCF) regime which is precisely identified numerically in the
present work.

Figure 19:  Second cumulant of $g$ as a function of first cumulant of $g$
showing the universal behavior in the metallic regime. Different curves
correspond to different transverse lengths ($L_{y}=10,40,80$) or magnetic
disorder ($J=0,0.2$).

Figure 20:  Second cumulant of $g$ as a function of longitudinal length, for
$L_{y}=80$. Dots represent data for $J=0$ and squares for $J=0.2$. The value
of UCF is shown in each symmetry class (with or without magnetic disorder).
UCF regime is reached in both cases.

Figure 21:  Second cumulant of $g$ as a function of longitudinal length, for
$L_{y}=40$. Dots represent data for $J=0$ and squares for $J=0.2$. The value
of UCF is shown in each symmetry class (with or without magnetic disorder. UCF
regime is reached in the case of magnetic impurities but not for scalar
impurities.

Figure 22:  Variance of $g$ as a function of longitudinal size. UCF are shown.
Different curves correspond to different values of magnetic disorder $J$
($J=0.025\to 0.4$). Transverse length $L_{y}=40$. The only free parameter in
analytical fits is the magnetic length $L_{m}$.

Figure 23:  Variance of $g$ as a function of $L_{x}/L_{m}$. Transverse length
$L_{y}=40$.

The plot of figure 19 confirms analytical results from [31] both qualitatively
in the shape of the curves and quantitatively in the values of fluctuations in
both universality classes. In our study, values of UCF are reached with a
maximal error of $1\%$ for GOE and $3\%$ for GUE with respect to the
analytical value of the UCF in the regime independent of $\langle g\rangle$
(i.e with much higher precision than e.g [38] and [30]).

Let us finally comment the work of Z. Qiao et al. [33] who have performed a
similar numerical Landuer study of 1D transport in various universality
classes, focusing mostly on the metallic regime. While both our study agree on
the UCF (although we have a higher accuracy for $\beta=1$), we did not find
evidence for a second universal conductance plateau. This result would
definitely deserves further study.

###  5.2 The average conductance

The main contribution to the average conductance is of classical origin.
However a weak localization correction must be added when the quantum behavior
of electrons is taken into account [12]. This quantum part manifests itself in
the magneto-conductance behavior of long wires (larger than the the phase
coherence length) where a weak magnetic field is sufficient to destroy this
quantum correction by dephasing the various diffusing path with respect to
each other[12]. We can write this average conductance as (without any magnetic
field):

$\langle g\rangle(J,L_{x},L_{m})=g_{class}(J,L_{x})+\delta
g_{WL}(L_{x},L_{m}),$ (18)

where $g_{class}$ is the classical part of the conductance and $\delta g_{WL}$
is the weak localization correction. For a quasi one dimensional system the
quantum correction reads the simple form[12]

$\delta
g_{WL}=\sum_{n=1}^{\infty}\left(\frac{-1/\pi^{2}}{n^{2}+2\left(\frac{L_{x}}{L_{m}}\right)^{2}}-\frac{3/\pi^{2}}{n^{2}+\frac{2}{3}\left(\frac{L_{x}}{L_{m}}\right)^{2}}\right).$
(19)

The knowledge of the magnetic length $L_{m}$ we gained in the previous study
of conductance fluctuations can now be used to completely characterize this
weak localization contribution to the conductance. By subtracting the
corresponding contribution to the average conductance, we obtain the classical
conductance, plotted in figure 24 as a function of longitudinal length.

Figure 24:  Evolution of classical conductance $g_{class}=\langle
g\rangle-\delta g$ with the longitudinal length for $L_{y}=40$ and
$J=0.05,0.1,0.2,0.4$. Numerical data are well fitted using the conductivity as
the only free parameter. different data correspond to different values of
magnetic disorder. In inset is shown the set of numerical data for $J=0.2$.
The dotted line is the analytical fit without taking into account the contact
conductance. Plain line is the same analytical fit as in the main plot.

Taking into account the contact resistance[39] in the two terminal setup, the
expected expression for this classical conductivity reads

$g_{class}(J,L_{x})=\frac{1}{\frac{1}{L_{y}}+\frac{L_{x}}{\sigma_{0}(J)}}.$
(20)

Figure 24 shows this expression plotted for different values of magnetic
disorder. The corresponding fitting parameter $\sigma_{0}$ is plotted as a
function of $J$ in figure 25.

Figure 25:  Evolution of conductivity $\sigma_{0}$ with magnetic disorder $J$
for $L_{y}=40$. Dots are numerical data. Plain line is the theory given by
Einstein relation and Matthiesen rule for the conductivity. Agreement is good
especially at low magnetic disorder.

To compare these results with theory, consider the Einstein relation which
links the (Einstein) conductivity to the diffusion constant:

$\sigma_{0}=se^{2}\rho_{0}(\varepsilon_{F})D,$ (21)

where $s$ is the spin degeneracy and $\rho_{0}(\varepsilon_{F})$ the
electronic density of states at the Fermi level. By definition, the diffusion
coefficient reads, for non magnetic impurities:

$D=v_{F}^{2}\tau_{e},$ (22)

with $v_{F}$ the Fermi velocity and $\tau_{e}$ the elastic scattering time. It
can be related to the scalar disorder by:

$\tau_{e}=\frac{1}{2\pi\rho_{0}n_{i}v_{0}^{2}}$ (23)

where $n_{i}$ is the impurity density and $v_{0}^{2}=W^{2}/12$. For more than
one diffusive process, it is compulsory to use the Matthiesen rule that
modifies scattering time $\tau_{e}$ in the following way:

$\frac{1}{\tau_{e}}\to\frac{1}{\tau_{e}}+\frac{1}{\tau_{m}},$ (24)

where $\tau_{m}=L_{m}^{2}/D$ and is related to the magnetic disorder:

$\tau_{m}=\frac{1}{2\pi\rho_{0}n_{i}J^{2}\langle S^{2}\rangle}.$ (25)

Using this allows one to get the $J$ dependance of the Einstein conductivity:

$\sigma_{0}(J)=\frac{\sigma_{0}(J=0)}{1+\frac{3}{W^{2}}J^{2}}.$ (26)

In figure 25 we compare this expression with numerical evaluation of the
conductivity. The good agreement between both curves provides an additional
check of the correct determination of the magnetic dephasing length $L_{m}$.
Including the magnetic disorder dependance of the diffusion coefficient
(through Matthiesen rule), we obtain a perturbative expression to second order
in $J$ for this magnetic dephasing length:

$L_{m}(J)=\sqrt{D(J)\tau_{m}(J)}\propto\frac{1}{J\sqrt{\frac{W^{2}}{12}+\frac{J^{2}}{4}}}.$
(27)

Figure 26:  Evolution of magnetic dephasing length $L_{m}$ with magnetic
disorder $J$ for $L_{y}=40$. Dots are numerical data extracted from the study
of UCF and diamonds are numerical data extracted from the study of conductance
correlations for two different spin configurations. Plain and dotted lines are
analytical fits from perturbation theory at second order in $J$. Error bars
are smaller than dots and diamonds sizes.

On figure 26 we have plotted the numerical evaluation of the magnetic length
as a function of magnetic disorder and the corresponding fit with eq.(27). We
notice that it is also possible to obtain the magnetic length via the study of
correlations of conductance, i.e via the study of $\langle
g(V,\\{\vec{S}_{i}^{(1)}\\}_{i})g(V,\\{\vec{S}_{i}^{(2)}\\}_{i})\rangle_{c}$,
which goes beyond the scope of the present paper[11].

###  5.3 The third cumulant

Finally we consider the third cumulant of the distribution of conductance.
According to the analytical study of [31], this cumulant decays to zero in a
universal way as $\langle g\rangle$ increases. Here in figure 28 we find a
dependance of this decrease on the symmetry class: for GOE $\langle
g^{3}\rangle_{c}$ goes to zero in a monotonous way whereas it decreases,
changes its sign and then goes to zero in GUE case. For $\langle g\rangle>4$
numerical errors are dominant, then this part of the curve is irrelevant.
Moreover, for GUE this decrease seems to be universal whereas it depends on
the transverse length for GOE.

Figure 27:  Plot of $\langle g^{3}\rangle_{c}$ as a function of $\langle
g\rangle$ in the metallic regime, averages are performed with various a number
of configurations $N_{des}$. Convergence curves are shown for $L_{y}=10$ and
$J=0$ or $J=0.2$.

On figure 27, is represented the convergence of the skewness when increasing
the number of configurations used to perform averages $N_{des}$ for both GOE
and GUE. Plots show a good enough convergence of averages to conclude that the
third cumulant of conductance is not zero for all values of $\langle
g\rangle$. Notice that the maximal number of averages is 50000.

Figure 28:  Plot of $\langle g^{3}\rangle_{c}$ as a function of $\langle
g\rangle$ in the metallic regime for $L_{y}=10,80$ and $J=0,0.2$.

Moreover this fast vanishing of the third cumulant confirms the faster
convergence of the whole distribution towards the gaussian, compared to what
happens in the insulating regime. Based on our numerical results, we cannot
confirm nor refute the expected law $\langle g^{3}\rangle_{c}\propto 1/\langle
g\rangle^{n}$, with $n=2$ in GOE and $n=3$ in GUE [40, 41].

##  6 Conclusion

To conclude we have conducted extensive numerical studies of electronic
transport in the presence of random frozen magnetic moments. Comparing and
extending previous analytical and numerical studies, we have identified the
insulating and metallic regimes described by the universality classes GOE and
GUE. We have paid special attention to the dependance on this symmetry of
cumulants of the distribution of conductance in both metallic and insulating
universal regimes. In particular, we have identified with high accuracy the
domain of universal conductance fluctuations, and determined its extension in
the present model. We have also determined precisely the so-called magnetic
length $L_{m}$ which represents the elastic scattering length of the spin on
magnetic impurities. This length is of primary importance in experiments as it
controls the crossover between universality classes. This work paves the way
for further studies of transport in metals with frozen magnetic impurities as
we have clearly identified the range of the parameters to access the
experimentally relevant metallic diffusive regime. One possible extension
consists in considering evolution of the statistics of conductance as the
magnetic disorder is varied, e.g. by rotating or flipping the spins of
impurities. Comparing the conductance obtained in both spin configurations
mimics the experimental measurement of the conductance of a low temperature
canonical spin glass after two successive quenches [42, 11], without the
necessary restrictions of analytical approaches[11]. Experimentally, this
approach could give access to fundamental properties of a spin glass, that
have never been measured.

We thank X. Waintal for useful discussions. This work was supported by the ANR
grants QuSpins and Mesoglass. All numerical calculations were performed on the
computing facilities of the ENS-Lyon calculation center (PSMN).

## References

  * [1] Juliette Billy, Vincent Josse, Zhanchun Zuo, Alain Bernard, Ben Hambrecht, Pierre Lugan, David Clément, Laurent Sanchez-Palencia, Philippe Bouyer, and Alain Aspect.  Direct observation of anderson localization of matter waves in a controlled disorder.  Nature, 453:891–894, 2008.
  * [2] Giacomo Roati, Chiara D’Errico, Leonardo Fallani, Marco Fattori, Chiara Fort, Matteo Zaccanti, Giovanni Modugno, Michele Modugno, and Massimo Inguscio.  Anderson localization of a non-interacting bose–einstein condensate.  Nature, 453:895–898, 2008.
  * [3] P.G.N. de Vegvar, L.P. Lévy, and T.A. Fulton.  Conductance fluctuations of mesoscopic spin glasses.  Phys. Rev. Lett., 66:2380, 1991.
  * [4] J. Jaroszynski, J. Wrobel, G. Karczewski, T. Wojtowicz, and T. Dietl.  Magnetoconductance noise and irreversibilities in submicron wires of spin-glass n+-cdmnte.  Phys. Rev. Lett., 80:5635, 1998.
  * [5] G. Neuttiens, J. Eom, C. Strunk, H. Pattyn, C. Van Haesendonck, Y. Bruynseraede, and V. Chandrasekhar.  Thermoelectric effects in mesoscopic aufe spin-glass wire.  Europhys. Lett., 42:185, 1998.
  * [6] T. Capron, A. Perrat-Mabilon, C. Peaucelle, T. Meunier, D. Carpentier, L. P. Levy, C. Baeuerle, and L. Saminadayar.  Remanence effects in the electrical resistivity of spin glasses.  EPL, 93(2), JAN 2011.
  * [7] E. Abrahams, P.W. Anderson, D.C. Licciardello, and T.V. Ramakrishnan.  Scaling theory of localization: Absence of quantum diffusion in two dimensions.  Phys. Rev. Lett., 42:673, 1979.
  * [8] F. Evers and A. D. Mirlin.  Anderson transitions.  Rev. Mod. Phys., 80:1355, 2008.
  * [9] Shinsei Ryu, Andreas P Schnyder, Akira Furusaki, and Andreas W W Ludwig.  Topological insulators and superconductors: tenfold way and dimensional hierarchy.  New Journal of Physics, 12(065010), 2010.
  * [10] A. D. Mirlin.  Statistics of enrgy levels and eigenfunctions in disordered systems.  Phys. Rep., 326:259, 2000.
  * [11] G. Paulin, T. Capron, D. Carpentier, T. Meunier, C. Bäuerle, L.P. Lévy, and L. Saminadayar.  Spin flip induced conductance fluctuations.  preprint, 2011.
  * [12] E. Akkermans and G. Montambaux.  Mesoscopic Physics of electrons and photons.  Cambridge University Press, 2007.
  * [13] M. Mézard, G. Parisi, and M. Virasoro.  Spin Glass Theory and Beyond.  World Scientific, 1987.
  * [14] O. N. Dorokhov.  Zh. Eksp. Teor. Fiz., 85:1040, 1983.
  * [15] P. A. Mello, P. Pereyra, and N. Kumar.  Macroscopic approach to multichannel disordered conductors.  Ann. Phys., 181:290–317, 1988.
  * [16] R. Landauer.  Spatial variation of currents and fields due to localized scatterers in metallic conduction.  IBM J. Res. Dev., 1:223, 1957.
  * [17] L. I. Deych, M. V. Erementchouk, and A. A. Lisyansky.  Scaling in the one-dimensional anderson localization problem in the region of fluctuation states.  Phys. Rev. Lett., 90:126601, 2003.
  * [18] M. Buttiker, Y. Imry, R. Landauer, and S. Pinhas.  Generalized many-channel conductance formula with application to small rings.  Phys. Rev. B, 31:6207–6215, 1985.
  * [19] D. S. Fisher and P.A. Lee.  Relation between conductivity and transmission matrix.  Phys. Rev. B, 23:6851–6854, 1981.
  * [20] A. Croy, R. A. Roemer, and M. Schreiber.  Localization of electronic states in amorphous materials: recursive green’s function method and the metal-insulator transition at e¡¿0.  In K. Hoffman and A. Meyer, editors, Lecture notes in Computational Science and Engineering, volume Lecture notes in Computational Science and Engineering. Springer (Berlin), 2006.
  * [21] C. W. J. Beenakker.  Random-matrix theory of quantum transport.  Rev. Mod. Phys., 69:731, 1997.
  * [22] K. Slevin, P. Markos, and T. Ohtsuki.  Reconciling conductance fluctuations and the scaling theory of localization.  Phys. Rev. Lett., 86(16):3594, 2001.
  * [23] P. Markos.  Weak disorder expansion of lyapunov exponents of products of random matrices: A degenerate theory.  Journal of Statistical Physics, 70(3/4):899, 1993.
  * [24] J-L. Pichard.  Random matrix theory of scattering in chaotic and disordered media.  In P. Sebbah, editor, Waves and Imaging through Complex Media. Kluwer Academin Publishers, 2001.
  * [25] J-L. Pichard.  In NATO ASI SERIES B254 B. KRAMER, editor, Quantum Coherence in Mesoscopic systems, page 369. PLENUM, NEW YORK, 1991.
  * [26] I. Larkin.  Zh. Eksp. Teor. Fiz., 85:764, 1983.
  * [27] I.V. Lerner and Y. Imry.  Magnetic-field dependence of the localization length in anderson insulators.  Europhys. Lett., 29:49–54, 1995.
  * [28] A. MacKinnon and B. Kramer.  One-parameter scaling of localization length and conductance in disordered systems.  Phys. Rev. Lett., 47(21):1546, 1981.
  * [29] K. A. Muttalib and P. Wölfle.  ”one-sided” log-normal distribution of conductances for a disordered quantum wire.  Phys. Rev. Lett., 83(15):3013, 1999.
  * [30] P. Markos.  Dimension dependance of the conductance distribution in the nonmetallic regime.  Phys. Rev. B, 65:104207, 2002.
  * [31] L.S. Froufe-Pérez, P. Garcia-Mochales, P.A. Serena, P.A. Mello, and J.J. Saenz.  Conductance distributions in quasi-one-dimensional disordered wires.  Phys. Rev. Lett., 89(24):246403, 2002.
  * [32] K. A. Muttalib, P. Wölfle, A. García-Martín, and V. A. Gopar.  Nonanalyticity in the distribution of conductances in quasi-one-dimensional wires.  Europhys. Lett., 61:95, 2003.
  * [33] Z. Qiao, Y. Xing, and J. Wang.  New universal conductance fluctuation of mesoscopic systems in the crossover regime from metal to insulator.  ArXiv:0910.3475v1, 2009.
  * [34] V. A. Gopar, K. A. Muttalib, and P. Wölfle.  Conductance distribution in disordered quantum wires: Crossover between the metallic and insulating regimes.  Phys. Rev. B, 66:174204, 2002.
  * [35] C. W. J. Beenakker and H. van Houten.  Quantum transport in semiconductor nanostructures.  Solid State Physics, 44:1, 1991.
  * [36] M.R. Zirnbauer.  Super fourier analysis and localization in disordered wires.  Phys. Rev. Lett., 69:1584, 1992.
  * [37] A. D. Mirlin, A. Mueller-Groeling, and M. R. Zirnbauer.  Conductance fluctuations of disordered wires: Fourier analysis on supersymmetric spaces.  Ann. Phys., 236:325, 1994.
  * [38] M. Cieplak, B.R. Bulka, and T. Dietl.  Universal conductance fluctuations in spin glasses.  Phys. Rev. B, 44:12337, 1991.
  * [39] H. L. Engquist and P.W. Anderson.  Definition and measurement of the electrical and thermal resistances.  Phys. Rev. B, 24:1151, 1981.
  * [40] A. M. S. Macêdo.  Random-matrix approach to the quantum-transport theory of disordered conductors.  Phys. Rev. B, 49:1858, 1994.
  * [41] M.C.W. van Rossum, Igor V. Lerner, Boris L. Altshuler, and Th. M. Nieuwenhuizen.  Deviations from the gaussian distribution of mesoscopic conductance fluctuations.  Phys. Rev. B, 55:4710, 1997.
  * [42] D. Carpentier and E. Orignac.  Measuring overlaps in mesoscopic spin glasses via conductance fluctuations.  Phys. Rev. Lett., 100:057207, 2008.


# Quintessential Axions from a New Confining Force

Jihn E. Kim, Younggeun Kim and Soonkeon Nam  Department of Physics, Kyung Hee
University, 26 Gyungheedaero, Dongdaemun-Gu, Seoul 02447, Republic of Korea

###### Abstract

Almost massless mesons created at the condensation scale of extra nonabelian
gauge group can be candidates of dark energy source. This can be another
motivation for introducing an additional confining force.

Dark energy, Quintessential axion, Supersymmetry, Chiral symmetry

###### pacs:

11.25.Mj, 11.30.Er, 11.25.Wx, 12.60.Jv

##  I Introduction

The origin of dark energy (DE) in the Universe is most mysterious in the
current cosmology. At present, the DE density $(0.003\,\,\mathrm{eV})^{4}$ is
synonymous to the cosmological constant (CC). In the multiverse distinguished
by CC, the anthropic choice allows some universes to evolve to the ones
similar to our currently nearly flat and empty Universe to have creatures to
quest about it now Weinberg88 . On the other hand, if the current DE is a
transient one, parametrized by a varying $\theta$, the classically evolving
scalar field $\theta$ is interpreted as the source for the current domination
of the energy density of the Universe FriemanWaga98 ; Carroll98 . The name
‘quintessential axion’ KimNilles03 ; Hill02 ; KimRMP10  is by identifying
$\theta$ as the pseudo-Goldstone boson of spontaneously broken global symmetry
U(1)DE.

Light pseudoscalar particles are the lampposts to the road leading to physics
scales much above their masses. The well-known example is the pion triplet
looked at low energy strong interaction scale, which guided toward the
SU(2)${}_{L}\times$SU(2)R chiral symmetry at a higher energy scale Nambu61 ;
GellMann68 . Another is the very light axion which hints at the intermediate
scale KimPRL79 ; KNS18 . In this pseudoscalar scenario, the scale where the
symmetry is explicit is at the scale $E\gtrsim f\sim\Lambda^{2}/m$, for the
pseudoscalar mass $m$, where $\Lambda^{4}$ is a typical energy density
contributed by such a pseudoscalar. Before obtaining mass of quintessential
axion by the energy density perturbation of order $\Lambda^{4}$, the
pseudoscalar degree is the phase degree $\theta(x)$ of some unitary operators
such that they have kinetic energy terms in quantum field theory below the
defining scale $f$ of the pseudoscalar $a$. Below the scale $f$, the global
transformation is expressed as $\theta\to a/f$ where $f$ is the decay
constant. Above the scale $f$, the phase in $e^{i\,a/f}$ is not a dynamical
field, i.e. not depending on $x$, but it still represents a phase direction of
the global symmetry. On the other hand, if a pseudoscalar is originally
present in the theory, then $f$ is the defining scale of that theory. In this
sense, string scale $M_{s}$ itself is the decay constant for axions from
string theory.

Quintessential axion requires Carroll98 :

  * •

The decay constant is near the Planck scale $M_{\rm P}$, and

  * •

Mass is about $10^{-32\,}\,\mathrm{eV}$.

Since the Planck mass is the mass defining scale, the first condition looks
easy to be implemented. But, the singlets beyond the standard model may be
required to interact with light fields, in which case a judicious care is
needed to allow the VEV of singlets still remaining near the Planck scale.
Implementing the second condition needs some symmetries such that the breaking
scale is sufficiently small.

If the pseudoscalar mass $m$ is less than 1 MeV, it can rarely decay to two
photons and two neutrinos. Using the $\pi^{0}$ decay rate to two photons
$\Gamma(\pi^{0}\to 2\gamma)\simeq(\alpha_{\rm
em}^{2}/64\pi^{3})(M_{\pi^{0}}/f_{\pi})^{2}M_{\pi^{0}}$ (at leading order
Holstein02 ) and to two neutrinos $1.2\times 10^{-5}\Gamma(\pi^{0}\to
2\gamma)$ Shrock79 , we can get an idea on the lifetime of the pseudoscalar.
For it to be still present in our Universe, therefore we require
$\Gamma^{-1}>4.3\times 10^{17}{\rm s}=1/(1.53\times 10^{-33}\,\mathrm{eV})$.
Using the above $f-m$ relation, we get an idea on $m$ from the condition that
it survives until now, using the QCD scale $\Lambda$,

$\begin{split}\Gamma^{-1}\simeq\frac{64\pi^{3}\,\Lambda_{\rm
QCD}^{4}}{\alpha_{\rm em}^{2}m^{4}}\frac{1}{m}>\frac{1}{1.53\times
10^{-33}\,\mathrm{eV}}\end{split}$ (1)

or $m<65\,\,\mathrm{eV}$ for $\Lambda_{\rm QCD}=380\,\,\mathrm{MeV}$.111For
the very light axion, the specific coupling reduces it further to 24 eV
KimRMP10 . If we take $\Lambda^{4}$ as the current energy density of the
Universe, then pseudoscalars with mass less than $3.3\times
10^{-4}\,\,\mathrm{eV}$ survives until now.222Ultra-light axions (ULAs)
KimPRD16 ; Witten17  was suggested for the galactic scale structures.

We will be using two words for describing breaking the global symmetry. One is
the spontaneous breaking scale $f$ mentioned above. This is the vacuum
expectation value (VEV) of the scalar field which carries the non-vanishing
quantum number of the global symmetry and is determined from the interaction
terms respecting the global symmetry. The other is the explicit breaking terms
which are interaction terms explicitly breaking the global symmetry in
question and provide the curvature in the direction of $a$.

We will consider the QCD axion also because it is present in most models. Not
only providing a solution of the strong CP problem, the QCD axion can work as
a source of dark matter (DM) in the Universe PWW83 ; AS83 ; DF83 ; Bae08 ,
which depends on the confining scale $\Lambda_{\rm QCD}$ and the axion decay
constant $f_{a}$. So, let us start with building a field theoretic model
housing both the QCD and quintessential axions. Working out an explicit
example, we employ supersymmetry because supersymmetry reduces the number of
couplings and hence the discussion for generating symmetry breaking terms is
much simpler.

Let the gauge group for a new confining force be $\cal G$ and the flavor
symmetry of ex-quarks under $\cal G$ be U$(N)\times$U$(N)$ which reduces to
SU$(N)\times$SU$(N)$ by removing the heavy singlet meson obtaining mass by the
anomaly contribution. The light mesons
$\Pi^{i}_{j}(i=1,\cdots,N;j=1,\cdots,N)$ with Tr $\Pi=0$, belonging to the
adjoint representation ($\bf N^{2}-1$) of SU$(N)_{A}\subset$
SU$(N)\times$SU$(N)$, obtain mass by explicit breaking of the flavor symmetry
of ex-quarks.

One obvious scale breaking the global symmetry is the condensation scale
$\Lambda$ of the ex-quarks,

$\begin{split}\langle\overline{Q}_{L}T^{a\,i}_{j}Q_{L}\rangle=\Lambda^{3}e^{i\Pi^{ai}_{j}/f}\end{split}$
(2)

where $\overline{Q}_{L}T^{a\,i}_{j}Q_{L}$ is interpreted as the matrix for
adjoint representation $T^{a}$ sandwiched between the relevant fermion fields
$\overline{\psi}_{R}T^{a}\psi_{L}$. To construct a model based on discrete
symmetry, supersymmetric extension is simpler and we will discuss
supersymmetric extension below.

In the ${\cal N}=1$ supersymmetric extension, we include scalar partners of
ex-quarks also. These scalar quarks are defined to have the chirality of the
corresponding ex-quarks denoted by tilded fields. The
SU$(N)_{L}\times$SU$(N)_{R}$ representation of scalar ex-quarks are

$\begin{split}Q_{L}=\begin{pmatrix}a_{1}\\\ a_{2}\\\ \cdot\\\ \cdot\\\
\cdot\\\
a_{N}\end{pmatrix}_{L},~{}~{}~{}~{}\overline{Q}_{L}=\begin{pmatrix}\bar{b}_{1}\\\
\bar{b}_{2}\\\ \cdot\\\ \cdot\\\ \cdot\\\
\bar{b}_{N}\end{pmatrix}_{L}.\end{split}$ (3)

Usually the same notations $Q_{L}$ and $Q_{R}$ are used for the corresponding
superfields also,
$Q_{L}+\sqrt{2}\tilde{Q}_{L}\,\vartheta+F_{Q_{L}}\vartheta^{2}$, etc., and
interaction is given by $\int d^{4}xd^{2}\vartheta W$. In supersymmetric
extension, the R-charges of chiral superfields determine possible terms in the
superpotential $W$. For dynamical SUSY breaking via the new confining gauge
group ${\cal G}$, the gaugino condensation scale $\Lambda$ is expected around
$10^{13\,}\,\textrm{GeV}$ Nilles82 . With U(1)DE, we also consider the
breaking scale of U(1)PQ for the QCD axion KNS18 . For simplicity, let us
assume the VEV of the singlet scalar field $\langle\sigma\rangle(\equiv V)$ is
the QCD axion scale $f_{a}\approx 10^{10\,}\,\textrm{GeV}$. In addition to
singlet fields, below we will consider Higgs fields also.

In supersymmetric models, condensation is possible for bilinears of fermions
and also bilinears of scalar ex-quarks (bi-scalar). Condensation of bi-scalars
does not break supersymmetry. This is a merit introducing quintessential axion
in the phases of bi-scalar condensation. So, $f$ for the bi-scalar
condensation can be different from the gaugino condensation scale $\Lambda$.

Condensation of bi-scalars is expected when the gauge coupling of the new non-
abelian force becomes strong at $\Lambda$. Let the bi-scalar singlet be $X$,

$\begin{split}\overline{Q}_{L}Q_{L}\equiv X.\end{split}$ (4)

As mentioned above, nonzero $X$ does not break supersymmetry. If we consider a
superpotential in terms of $X$,

$\begin{split}&W=\Lambda X-\frac{1}{2M}X^{2}+\cdots.\\\ \end{split}$ (5)

where $\Lambda$ and $1/M$ can result from the VEVs of some singlet fields and
$M$ in the following is used for the Planck mass $M_{\rm P}=2.43\times
10^{18\,}\,\textrm{GeV}$. So, the first minimum of the potential appears from
$V=\left(\Lambda-\frac{1}{M_{\rm P}}X\right)^{2}$ and the next minimum can be
far above the Planck scale. Then, the decay constant $f=\sqrt{X}$ is expected
at a median of $\Lambda$ and $M_{\rm P}$.

Now let us consider the explicit breaking terms. Anomalies of non-abelian
gauge groups are inevitable breaking terms, and the QCD anomaly is used for
the explicit breaking term of the U(1)PQ global symmetry. If we consider a
phase of the SM singlet scalar field, the weak gauge group in the SM may work
for an explicit breaking term of the U(1)DE global symmetry Kim21PLB . But our
quintessence fields are mesons in the phase of bi-scalars, or scalar ex-
quarks, and hence the SU(2)W anomaly cannot work for our quintessence fields.
The only method is from the superpotential made of scalar fields. To have a
sufficiently small numerical value for the vacuum energy, we consider the VEV
scale of the SM Higgs fields. A possible superpotential is

$\begin{split}\Delta
W=\frac{1}{M^{n+3}}\overline{Q}_{L}Q_{L}\left(H_{u}H_{d}\right)^{2}\sigma^{n},\end{split}$
(6)

where we started from $\left(H_{u}H_{d}\right)^{2}$. If we started with one
power of $H_{u}H_{d}$, too large a value of $n$ would be required. Considering
condensation of the hidden sector quark $Q$ in Eq. (6), the vacuum energy
density is

$\begin{split}\frac{4}{M^{n+3}}\Lambda^{3}(v_{u}v_{d})^{2}V^{n}=\frac{4}{M^{n+3}}\Lambda^{3}\frac{v_{d}^{4}}{\cos\beta^{4}}V^{n}\simeq(0.003\,\mathrm{eV})^{4}\end{split}$
(7)

where $\tan\beta=v_{u}/v_{d}$.

Figure 1:  Potential generated by Yukawa terms breaking U(1)DE. At the
intersection of the blue curve and the $f_{u}=1$ line, $v_{d}$ is
$25.6\,\textrm{GeV}$.

The top quark mass $m_{t}$ is $f_{u}v_{u}/\sqrt{2}$ and the bottom quark mass
$m_{b}$ is $f_{d}v_{d}/\sqrt{2}$. For the perturbative calculation to be
possible, we limit the study in the region $f_{u}\leq 1$ and $f_{d}\leq 1$.
Equation (7) has one parameter $\tan\beta$. For a given $\tan\beta$, relative
values of $f_{u}$ and $f_{d}$ are given, as shown in Fig. 2. The blue curve in
Fig. 2 is interpreted as the current DE, $(0.003\,\mathrm{eV})^{4}$, for
$f_{d}=1$.

As an example, taking $v_{d}=20\,\textrm{GeV}$ and
$m_{t}=173\,\,\textrm{GeV}$, Eq. (7) is estimated as,

$\begin{split}\frac{4}{M^{n+3}}\Lambda^{3}\frac{v_{d}^{4}}{\cos^{4}\beta}V^{n}\simeq\tan^{4}\beta\left(\frac{V}{M}\right)^{n}\left(\frac{20}{2.43\times
10^{18}}\right)^{4}\Lambda^{3}\cdot 0.97\times
10^{28}\,\mathrm{eV}=10^{-10}\,\mathrm{eV}^{4},\end{split}$ (8)

leading to

$\begin{split}\Lambda=\left\\{\begin{array}[]{l}2.7\times
10^{+14}\,\textrm{GeV},~{}\textrm{for }n=1\\\ 1.7\times
10^{+17}\,\textrm{GeV},~{}\textrm{for }n=2\\\ \end{array}\right.\end{split}$
(9)

This example is shown as the red bullet in Fig. 2.

| Representation under ${\cal G}\equiv$SU(${\cal N}$) | SU$(2)_{W}\times$U(1)Y | ${{\bf Z}}_{6R}$ |
---|---|---|---|---
$Q_{L}$ | ${\cal N}$ | ${\bf 1}$ | $+1$ |
$\overline{Q}_{L}$ | $\overline{\cal N}$ | ${\bf 1}$ | $-1$ |
$H_{u}$ | ${\bf 1}$ | ${\bf 2}_{+1/2}$ | $+3$ |
$H_{d}$ | ${\bf 1}$ | ${\bf 2}_{-1/2}$ | $+2$ |
$\sigma$ | ${\bf 1}$ | ${\bf 1}$ | $+4$ |
$S$ | ${\bf 1}$ | ${\bf 1}$ | $+5$ |
Table 1: ${{\bf Z}}_{6R}$ quantum numbers of relevant chiral superfileds
appearing in Eq. (6).

In Table 1, an example for the ${{\bf Z}}_{6R}$ discrete symmetry is shown,
from which $n=1$ is chosen in Eq. (9). The ${{\bf Z}}_{6R}$ is a discrete
subgroup of R-symmetry U(1)R. Because the antisymmetric parameter $\vartheta$
carries –1 unit of U(1)R charge, superpotential $W$ with +2 units of U(1)R
charge survives in the integration $\int d^{2}\vartheta\,W$. To assign tree
level Yukawa couplings of the SM quarks, we have the following quantum numbers

$\begin{split}q_{L}:0,~{}u_{L}^{c}:-1,~{}d_{L}^{c}:0,\end{split}$ (10)

such that both $q_{L}u_{L}^{c}H_{u}$ and $q_{L}d_{L}^{c}H_{d}$ carry two units
(modulo 6) of U(1)R charge. To obtain a TeV scale $\mu$ term KimNilles84 , one
can introduce the following superpotential

$\begin{split}W_{\mu}=\frac{(10^{10\,}\,\textrm{GeV})^{2}}{M}\,H_{u}H_{d},\end{split}$
(11)

but should forbid the dimensions 2 and 3 superpotential terms
$H_{u}H_{d},H_{u}H_{d}\sigma$ and $H_{u}H_{d}S$. Equation (11) can be
satisfied with $\langle S\rangle\approx\langle\sigma\rangle\approx
10^{10\,}\,\textrm{GeV}$. Then, singlets $S$ and $\sigma$ lead to the
following superpotential

$\begin{split}W=-\alpha\sigma
S^{2}+\frac{\varepsilon}{M}S^{4}-\frac{x}{M^{2}}\sigma
S^{2}Q_{L}\overline{Q}_{L}+\cdots,\end{split}$ (12)

where the neglected terms $\cdots$ denote other terms including
$Q_{L}\overline{Q}_{L}$. The relations between the VEVs of $\sigma$ and $S$
are

$\begin{split}&\frac{\partial W}{\partial\sigma}:\to
Q_{L}\overline{Q}_{L}=-\frac{\alpha M^{2}}{x}\\\ &\frac{\partial W}{\partial
S}:\to(x\frac{Q_{L}\overline{Q}_{L}}{M^{2}}+\alpha)\sigma=\frac{2\varepsilon}{M}S^{2}.\end{split}$
(13)

Relations in (13) give exact SUSY. But SUSY is broken dynamically by the
condensation of ex-quarks at $\Lambda\simeq 10^{13\,}\,\textrm{GeV}$. This
SUSY breaking effect is added in the first SUSY relation,

$\begin{split}-\alpha
S^{2}-\frac{x}{M^{2}}S^{2}Q_{L}\overline{Q}_{L}+\delta_{1}\Lambda^{2}=0,\\\
\frac{2\varepsilon S^{3}}{M}-\alpha
S\sigma-\frac{x}{M^{2}}SQ_{L}\overline{Q}_{L}\sigma+\frac{\delta_{2}\Lambda^{2}}{2}=0,\\\
2\frac{\varepsilon}{M}S^{4}+(\frac{\delta_{2}S-2\delta_{1}\sigma}{2})\Lambda^{2}=0.\end{split}$
(14)

Solutions of $\sigma$ and $S$ satisfying Eq. (14) are shown in Fig. 2. So, the
needed $\mu$ term is around the point $(1,1)$ in Fig. 2.

The above U(1)R symmetry is a kind of the PQ symmetry because the U(1)R
quantum number of $H_{u}H_{d}$ is non-zero. At the scale $V_{\rm
PQ}\equiv\sqrt{\sigma^{2}+S^{2}}$, therefore, the PQ symmetry is spontaneously
broken and there results the QCD axion discussed in the literature KimRMP10 .
Because the quantum numbers of $\sigma$ and $S$ are relatively prime, all the
axionic vacua are connected and the domain wall number is 1 KNS18 ; ChoiKim85
. This QCD axion is responsible for the DM in the Universe. On the other hand,
the mesons resulting from the new confining force are responsible for the DE
of the Universe.

Figure 2: Solutions of $\sigma$ and $S$ satisfying Eq. (14).

Gravity effects of breaking global and discrete symmetries can be given by the
terms explicitly breaking them. Since supergravity with superpotential $W$ has
a fine-tuning problem of adjusting the cosmological constant Nilles04 , we
lack any understanding on the magnitude of the cosmological constant in
supergravity. So, no attempt is made to introduce explicit DE terms. However,
we can discuss the gravity effects to the parameters in the gauge part,
especially on the $\bar{\theta}$ of QCD, as done in Barr92 . Here, we consider
only the potential $V$ satisfying the ${\bf Z}_{6R}$ symmetry, but not the
superpotential $W$. Possible ${\bf Z}_{6R}$ invariant terms are

$\begin{split}\sigma^{3},~{}{\sigma^{*}}^{3},~{}S^{2}{\sigma^{*}},~{}S^{6},~{}{S^{*}}^{6},~{}(\sigma\sigma^{*})^{n}~{}{\rm
for}~{}n=1,2,\cdots,~{}(SS^{*})^{n}~{}{\rm for}~{}n=1,2,\cdots,~{}{\rm
etc.}\end{split}$ (15)

Note that the $(\sigma\sigma^{*})^{n}$ and $(SS^{*})^{n}$ terms respect the PQ
symmetry and do not shift $\bar{\theta}$. The leading terms in shifting
$\bar{\theta}$ in Eq. (15) are the cubic terms. But these cubic terms do not
arise from the superpotential and hence appear with supersymmetry breaking
coefficients, i.e. the gaugino condensation $\langle{\cal G}{\cal G}\rangle$
of the new confining force. Thus the $\bar{\theta}$ shift is estimated from

$\begin{split}\frac{1}{M^{2}}\langle{\cal G}{\cal
G}\rangle\sigma^{3},\cdots\end{split}$ (16)

which can be compared to the axion (or Pontryagin number) shift by $\langle
G_{\mu\nu}\tilde{G}^{\mu\nu}\rangle\approx\Lambda^{4}$,

$\begin{split}\Delta\bar{\theta}\approx\frac{\sigma^{3}}{M^{2}\Lambda}\approx
1.7\times 10^{-20},\end{split}$ (17)

for $\Lambda=10^{13}\,\textrm{GeV}$ and $\sigma=10^{10}\,\textrm{GeV}$. So,
the $\bar{\theta}$ shift by gravity effects falls in the allowed region,
$|\bar{\theta}|\lesssim 10^{-11}$ KimRMP10 . In addition, the higher order
terms will give a potential of the form given in Fig. 3. Restrcting the vaccum
for $a$ below $M_{\rm P}$ as shown in the left-hand side of the dash line, the
nonzero VEV of $a$ is near the Planck scale. In addition, note that the
effects of pseudo-Goldstone bosons may induce terms expected in general
relativityBludman77 . We note, however, that all such effects with
supersymmetry are included in Eq. (12).

Figure 3: Shape of the potential by Eq. (15).

In conclusion, we propsed a mechanism creating quintessential axions that are
almost massless mesons created at the condensation scale of a new nonabelian
gauge group. Since an extra nonabelian gauge group has been proposed for
dymnamical supersymmetry breaking, our mechanism is better suited in
supersymmetric models.

###### Acknowledgements.

This work is supported in part by the National Research Foundation (NRF) grant
NRF-2018R1A2A3074631.

## References

  * (1)
  * (2) S. Weinberg, _Anthropic Bound on the Cosmological Constant_ , Phys. Rev. Lett. 59 (1987) 2607 [doi: 10.1103/PhysRevLett.59.2607].
  * (3) J. A. Frieman and I. Waga, _Constraints from high redshift supernovae upon scalar field cosmologies_ , Phys. Rev. D 57 (1998) 4642 [doi: 10.1103/PhysRevD.57.4642].
  * (4) S. M. Carroll, _Quintessence and the rest of the world_ , Phys. Rev. Lett. 81 (1998) 3067 [e-Print: astro-ph/9806099 [astro-ph]].
  * (5) J. E. Kim and H. P. Nilles, _A quintessential axion_ , Phys. Lett. B 553 (2003) 1 [e-print: hep-ph/0210402 [hep-ph]].
  * (6) C. T. Hill, _Natural theories of ultralow mass pseudo Nambu-Goldstone bosons: Axions and quintessence_ , Phys. Rev. D 66 (2002) 075010 [e-print:hep-ph/0205237].
  * (7) J. E. Kim and G. Carosi, _Axions and the strong CP problem_ , Rev. Mod. Phys. 82 (2010) 557 [e-print: 0807.3125 [hep-ph]].
  * (8) Y. Nambu and G. Jona-Lasinio, _Dynamical Model of Elementary Particles Based on an Analogy with Superconductivity. 1._ , Phys. Rev. D 122 (1961) 345 [doi:10.1103/PhysRev.122.345].
  * (9) M. Gell-Mann, R. J. Oakes, and B. Renner, _Behavior of current divergences under SU(3) $\times$SU(3)_, Phys. Rev. D 175 (1968) 2195 [doi:10.1103/PhysRev.175.2195].
  * (10) J. E. Kim, _Weak Interaction Singlet and Strong CP Invariance_ , Phys. Rev. Lett. 43 (1979) 103 [doi:10.1103/ PhysRevLett.43.103].
  * (11) For a review, see, J. E. Kim, S. Nam, and Y. Semertzidis, _Fate of global symmetries in the Universe: QCD axion, quintessential axion and trans-Planckian inflaton decay-constant_ , Int. J. Mod. Phys. A 33 (2018) 1830002 [e-print: 1712.08648 [hep-ph]].
  * (12) J. L. Goity, A. M. Bernstein, and B. R. Holstein, _Decay $\pi^{0}\to\gamma\gamma$ to next to leading order in chiral perturbation theory_, Phys. Rev. D 66 (2002) 076014 [e-print:hep-ph/0206007].
  * (13) R. E. Shrock and M. Voloshin, _Bounds on quark mixing angles from the decay $K_{L}\to\mu\bar{\mu}$_, Phys. Lett. B 87 (1979) 375 [doi:10.1016/0370-2693(79)90557-4].
  * (14) J. E. Kim and D. J. E. Marsh, _An ultralight pseudoscalar boson_ , Phys. Rev. D 93 (2016) 025027 [e-Print: 1510.01701 [hep-ph]].
  * (15) L. Hui, J. P. Ostriker, S. Tremaine, and E. Witten, _Ultralight scalars as cosmological dark matter_ , Phys. Rev. D 95 (2017) 043541 [e-print:1610.08297 [astro-ph.CO]].
  * (16) J. Preskill, M. B. Wise, and F. Wilczek, _Cosmology of the Invisible Axion_ , Phys. Lett. B 120 (1983) 127 [doi:10.1016/0370-2693(83)90637-8].
  * (17) L. F. Abbott and P. Sikivie, _A Cosmological Bound on the Invisible Axion_ , Phys. Lett. B 120 (1983) 133 [doi:10.1016/0370-2693(83)90638-X].
  * (18) M. Dine and W. Fischler, _The Not So Harmless Axion_ , Phys. Lett. B 120 (1983) 137 [doi:10.1016/0370-2693(83)90639-1].
  * (19) See a fit, K. J. Bae, J-H. Huh, and J. E. Kim, _Update of axion CDM energy density_ , JCAP 09 (2008) 005 [e-print: 0806.0497 [hep-ph]].
  * (20) H. P. Nilles, _Dynamically Broken Supergravity and the Hierarchy Problem_ , Phys. Lett. B 115 (1982) 193 [doi: 10.1016/0370-2693(82)90642-6].
  * (21) J. E. Kim, _Anomalies and parities for quintessential and ultra-light axions_ , Phys. Lett. B 817 (2021) 136248 [e-Print: 2102.01795 [hep-ph]].
  * (22) M. Bronstein, Phys. Z. Sowjetunion 3 (1933) 73;
M. $\ddot{\rm O}$zer, M. O. Taha, Nucl. Phys. B 287 (1987) 797 ;
B. Ratra and P. J. E. Peebles, _Cosmological consequences of a rolling
homogeneous scalar field_ , Phys. Rev. D 37 (1988) 3406
[doi:10.1103/PhysRevD.37.3406];
C. Wetterich, _Cosmologies with variable Newton’s ‘constant’_ , Nucl. Phys. B
302 (1988) 645 [doi:10.1016/0550-3213(88)90192-7];
H. Gies and C. Wetterich, _Renormalization flow from UV to IR degrees of
freedom_ , hep-ph/0205226;
J. A. Frieman, C. T. Hill, and R. Watkins, _Late-time cosmological phase
transitions: Particle-physics models and cosmic evolution_ , Phys. Rev. D 46
(1992) 1226 [doi:10.1103/PhysRevD.46.1226];
R. Caldwell, R. Dave, and P. J. Steinhardt, _Cosmological imprint of an energy
component with general equation of state_ , Phys. Rev. Lett. 80 (1998) 1582
[doi: 10.1103/PhysRevLett.80.1582];
P. Binetruy, _Models of dynamical supersymmetry breaking and quintessence_ ,
Phys. Rev. D 60 (1999) 063502 [e-Print: hep-ph/9810553 [hep-ph]];
C. Kolda and D. H. Lyth, _Quintessential difficulties_ , Phys. Lett. B 458
(1999) 197 [e-Print: hep-ph/9811375 [hep-ph]];
T. Chiba, _Quintessence, the gravitational constant, and gravity_ , Phys. Rev.
D 60 (1999) 083508 [doi: 10.1103/PhysRevD.60.083508];
P. Brax and J. Martin, Phys. Lett. B 468 (1999) 40 [];
A. Masiero, M. Pietroni, F. Rosati, _SUSY QCD and quintessence_ , Phys. Rev. D
61 (2000) 023504 [e-Print: hep-ph/9905346 [hep-ph]];
J. E. Kim, _Model dependent axion as quintessence with almost massless hidden
sector quarks_ , JHEP 06 (2000) 016 [e-Print: hep-ph/9907528 [hep-ph]];
M. C. Bento, O. Bertolami, Gen. Relativ. Gravit. 31 (1999) 1461 [e-Print: gr-
qc/9905075 [gr-qc]];
F. Perrotta, C. Baccigalupi, S. Matarrase, _Extended quintessence_ , Phys.
Rev. D 61 (2000) 023507 [e-Print: astro-ph/9906066 [astro-ph]];
J. E. Kim, _Model dependent axion as quintessence with almost massless hidden
sector quarks_ , JHEP 06 (2000) 016 [e-Print: hep-ph/9907528 [hep-ph]];
A. Arbey, J. Lesgourgues, P. Salati, _Cosmological constraints on
quintessential halos_ , Phys. Rev. D 65 (2002) 083514 [e-Print: astro-
ph/0112324 [astro-ph]].

  * (23) K. Choi, S. H. Im, and C. S. Shin, _Recent Progress in the Physics of Axions and Axion-Like Particles_ , Ann. Rev. Nucl. Part. Sci. 71 (2021) 225 [e-Print: 2012.05029 [hep-ph]].
  * (24) J. E. Kim, H. P. Nilles, and M. Peloso, _Completing natural inflation_ , JCAP 01 (2005) 005 [e-Print: hep-ph/0409138 [hep-ph]].
  * (25) For a recent survay, see, A. Lagu$\ddot{\rm e}$, J. R. Bond, R. Hlozek, K. K. Rogers, D. J. E. Marsh, and D. Grin, _Constraining Ultralight Axions with Galaxy Survaeys_ , eprint:2104.07802[astro-ph.CO].
  * (26) T. W. B.Kibble, _Topology of cosmic domains and strings_ , Jour. of Physics A: Mathematical and General 9 (1976) 1387.
  * (27) P. Sikivie, _Of Axions, Domain Walls and the Early Universe_ , Phys. Rev. Lett. 48 (1982) 1156 [doi:10.1103/ PhysRevLett.48.1156].
  * (28) S. M. Barr, K. Choi and J. E. Kim, _Axion Cosmology in Superstring Models_ , Nucl. Phys. B 283 (1987) 591 [doi:10.1016/0550-3213(87)90288-4].
  * (29) J. E. Kim and H. P. Nilles, _The $\mu$ problem and the strong CP problem_, Phys. Lett. B 138 (1984) 150 [doi:10.1016/0370-2693(84)91890-2].
  * (30) K. Choi and J. E. Kim, _Domain Walls in Superstring Models_ , Phys. Rev. Lett. 55 (1985) 2637 [doi:10.1103/PhysRevLett.55.2637].
  * (31) H P. Nilles, _Supersymmetry, Supergravity and Particle Physics_ , Phys. Rep. 110 (1984) 1 [doi:10.1016/0370-1573(84)90008-5].
  * (32) S. M. Barr and D. Seckel, _Planck scale corrections to axion models_ , Phys. Rev. D 46 (1992) 539 [doi: 10.1103/PhysRevD.46.539].
  * (33) S. A. Bludman and M. A. Ruderman, _Induced Cosmological Constant Expected above the Phase Transition Restoring the Broken Symmetry_ , Phys. Rev. Lett. 38 (1977) 255 [doi:10.1103/PhysRevLett.38.255].