#  Anomalous Conductivity of Two-dimensional Dirac Electrons in Organic
Conductor under Pressure

Yoshikazu Suzumura 1 and Masao Ogata 2 E-mail:
<EMAIL_ADDRESS>of Physics1Department of Physics
Nagoya University Nagoya University Nagoya 464-8602 Nagoya 464-8602 Japan
2Department of Physics Japan
2Department of Physics University of Tokyo University of Tokyo Bunkyo Bunkyo
Tokyo 113-0033 Tokyo 113-0033 Japan
Japan

###### Abstract

The electric conductivity of Dirac electrons in the organic conductor
$\alpha$-(BEDT-TTF)2I3 [BEDT-TTF = bis(ethylenedithio)tetrathiafulvalene]
under pressure has been examined using a two-dimensional tight-binding (TB)
model with both impurity and electron–phonon (e–p) scatterings. We study an
anomalous temperature dependence of the conductivity, which shows a crossover
from $\sigma_{x}<\sigma_{y}$ at low temperatures [region (I)] to
$\sigma_{x}>\sigma_{y}$ at high temperatures [region (II)]. $\sigma_{y}$ and
$\sigma_{x}$ are diagonal conductivities parallel and perpendicular to a
stacking axis of molecules, respectively. The effect of Dirac cone tilting is
dominant in region (I), whereas the anisotropy of the velocity of the Dirac
cone is dominant in region (II). Such behavior is further examined by
calculating the deviation of principal axes due to the off-diagonal
conductivity $\sigma_{xy}$ and a nearly constant conductivity at high
temperatures due to the e–p scattering, which is the extension of the previous
result of the simple two-band model [Phys. Rev. B 98,161205 (2018)]. The
relevance to experiments of organic conductors is discussed.

##  1 Introduction

Since the discovery of two-dimensional massless Dirac fermions, [1] extensive
studies have been performed on various materials. Among them, noticeable
phenomena of Dirac electrons in molecular crystals[2] have been studied in the
organic conductor $\alpha$-(BEDT-TTF)2I3[3] [BEDT-
TTF=bis(ethylenedithio)tetrathiafulvalene]. After noting that the density of
states (DOS) vanishes linearly at the Fermi energy, [4] the two-dimensional
Dirac cone with a zero-gap state (ZGS) [5] was found using a tight-binding
(TB) model, where transfer energies are estimated from the extended Hückel
method.[6] The existence of such Dirac cone was verified by first-principles
DFT calculation, [7] which has been used for studying further $\alpha$-(BEDT-
TTF)2I3 under hydrostatic pressure.[8]

There are common features among organic conductors with isostructure salts,
[9, 10] $\alpha$-D2I3 (D = ET, STF, and BETS), where ET = BEDT-TTF, STF =
bis(ethylenedithio)diselenadithiafuluvalene), and BETS =
bis(ethylenedithio)tetraselenafulvalene. These salts display an energy band
with a Dirac cone, [5, 11, 12, 13] and the resistivity at high temperatures
shows a nearly constant behavior, [9, 10, 14, 15, 16, 17, 18] whereas the
conventional metal shows the linearly increasing one. Such unconventional
behavior was also observed in Dirac electrons with nodal line semimetals of
single-component molecular conductors [19, 20, 21, 22, 23] Thus, it has been
considered that the nearly constant behavior of resistivity at high
temperatures is attributable to the intrinsic property of Dirac electrons.

Thus, the conductivity of Dirac electrons has been studied in simple models. A
two-band model shows that the static conductivity at absolute zero temperature
remains finite with a universal value, i.e., independent of the magnitude of
impurity scattering owing to a quantum effect.[24] The conductivity increases
with increasing doping concentration. The effect of Dirac cone tilting at
absolute zero was studied previously.[25] It provides the anisotropic
conductivity and the deviation of the current from the applied electric field.
[25] At finite temperatures, on the other hand, the conductivity depends on
the magnitude of the impurity scattering, $\Gamma$, which is proportional to
the inverse of the life-time by the disorder. With increasing temperature
($T$), the conductivity remains unchanged for $T\ll\Gamma$, whereas it
increases for $\Gamma\ll T$.[26] Noting that $\Gamma\sim$ 0.0003 eV for
organic conductors,[2] a monotonic increase in the conductivity at finite
temperature $T>0.0005$ eV is expected. However the measured conductivity (or
resistivity) on the above organic conductor shows an almost constant behavior
at high temperatures. To comprehend such an exotic phenomenon, the acoustic
phonon scatterings has been proposed as a possible mechanism, which was
studied using a simple two-band model of the Dirac cone without tilting. [27]
Although the scattering shows a reasonable suppression of the conductivity at
high temperatures, it is unclear if such a model sufficiently explains the
conductivity of the actual organic conductor, where the Dirac cone shows
deviation from the linear spectrum and then the calculation must be performed
beyond the model.

In this paper, we study the anomalous temperature dependence of the
anisotropic conductivity of $\alpha$-(BEDT-TTF)2I3 using the TB model [28, 29,
30, 8] and acoustic phonons. Thus, the present paper is a natural extension of
the previous study [27] in which a simplified model was used without the
details of the lattice structure of $\alpha$-(BEDT-TTF)2I3. By using the TB
model, we can automatically take into account the effects of tilting and
anisotropy discussed at absolute zero [25] to explore the temperature
dependence. It will be shown that the presence of acoustic phonons gives rise
to conductivity being nearly constant at high temperatures. We will also
demonstrate a crossover of the conductivity from a quantum regime at low
temperatures to a classical regime at high temperatures by focusing on a
competition between the tilted Dirac cone and the anisotropic velocity of the
cone. We use the TB model with transfer energies of $\alpha$-(BEDT-TTF)2I3
under pressure, which are well known compared with those of other
isostructural salts, BETS and STF. Since the tilted Dirac cone is also
obtained for BETS [12] and STF [13], it is expected that the present result
provides a common feature for these salts, i.e., the nearly constant
conductivity at high temperatures.

The paper is organized as follows. In Sect. 2, the model and formulation are
given for both the uniaxial and hydrostatic pressures, where the latter case
is examined by adding site potential. In Sect. 3, after examining the chemical
potential and density of states (DOS), the conductivity is calculated for
$\alpha$-(BEDT-TTF)2I3 under both uniaxial and hydrostatic pressures. Section
4 is devoted to summary and discussion.

##  2 Model and Formulation

Figure 1: (Color online) Crystal structure, where there are four molecules A,
A’, B, and C in the unit cell (dot-dashed line), which forms a square lattice.
Transfer energies are shown by $a_{1},\cdots,b_{4}$ for the nearest neighbor
(NN) sites and $a_{1d}$, $a_{3d}$, and $a_{4d}$ for the next-nearest neighbor
(NNN) sites. The cross denotes an inversion center between A and A’.

We consider a two-dimensional Dirac electron system per spin, which is given
by

$H_{\rm total}=H_{0}+H_{1}+H_{\rm p}+H_{\rm e-p}+H_{imp}\;.$ (1)

$H_{0}$ describes a TB model of the organic conductor $\alpha$-(BEDT-TTF)2I3
consisting of four molecules per unit cell (Fig. 1), where $H_{1}$ describes a
site potential obtained from a mean field of short-range repulsive
interactions.[8] $H_{\rm p}$ and $H_{\rm e-p}$ denote an acoustic phonon and
an electron-phonon (e–p) interaction, respectively. $H_{\rm imp}$ is the
impurity potential. The terms $H_{0}+H_{\rm p}+H_{\rm e-p}$ provide the
Fröhlich Hamiltonian [31] applied to the present Dirac electron system. The
unit of the energy is taken as eV.

###  2.1 Energy band

First, we derive the energy band for $H=H_{0}+H_{1}$ and the associated
quantities. A TB model, $H_{0}$, is expressed as

$\displaystyle H_{0}$ $\displaystyle=$
$\displaystyle\sum_{i,j=1}^{N}\sum_{\alpha,\beta=1}^{4}t_{i,j;\alpha,\beta}a^{\dagger}_{i,\alpha}a_{j,\beta}$
(2) $\displaystyle=$
$\displaystyle\sum_{\bm{k}}\sum_{\alpha,\beta=1}^{4}t_{\alpha,\beta}(\bm{k})a^{\dagger}_{\alpha}(\bm{k})a_{\beta}(\bm{k})\;,$

where $a^{\dagger}_{i,\alpha}$ denotes a creation operator of an electron of
molecule $\alpha$ [= A(1), A’(2), B(3), and C(4)] in the unit cell at the
$i$-th lattice site. $N$ is the total number of square lattice sites and
$t_{i,j;\alpha,\beta}$ are the transfer energies for the nearest neighbor (NN)
and next-nearest neighbor (NNN) sites, which are shown in Fig. 1. A Fourier
transform for the operator $a_{j,\alpha}$ is given by
$a_{j,\alpha}=1/N^{1/2}\sum_{\bm{k}}a_{\alpha}(\bm{k})\exp[i\bm{k}\cdot\bm{r}_{j}]$,
where $\bm{k}=(k_{x},k_{y})$ and the lattice constant is taken as unity. The
quantity $H_{1}$ corresponds to a site potential, $V_{\alpha}$, acting on the
$\alpha$ site, where $V_{\rm A}=V_{\rm A^{\prime}}$ due to an inversion
symmetry around the cross in Fig. 1. The Hamiltonian $H_{1}$ is obtained as
(Appendix A)

$\displaystyle H_{1}$ $\displaystyle=$
$\displaystyle\sum_{\alpha}(V_{\alpha}-V_{\rm A})\hat{n}_{\alpha}$ (3)
$\displaystyle=$ $\displaystyle\tilde{V}_{\rm B}\hat{n}_{\rm B}+\tilde{V}_{\rm
C}\hat{n}_{\rm C}\;,$

where $\tilde{V}_{\alpha}$ denotes a potential measured from that of the A
site and $\hat{n}_{\alpha}=a^{\dagger}_{\alpha}(\bm{k})a_{\beta}(\bm{k})$.
From Eqs. (2) and (3), $H$ is written as[8]

$\displaystyle
H=\sum_{\bm{k}}\sum_{\alpha,\beta}a^{\dagger}_{\alpha}(\bm{k})h_{\alpha,\beta}a_{\beta}(\bm{k})\;,$
(4)

with $h_{\alpha,\beta}$ being the matrix element (Appendix A) and is
diagonalized as

$\displaystyle
H=\sum_{\bm{k}}\sum_{\gamma}c_{\gamma}^{\dagger}(\bm{k})E_{\gamma}(\bm{k})c_{\gamma}(\bm{k})\;,$
(5a) where $E_{1}(\bm{k})>E_{2}(\bm{k})>E_{3}(\bm{k})>E_{4}(\bm{k})$ and
$\displaystyle\sum_{\beta}h_{\alpha,\beta}(\bm{k})d_{\beta\gamma(\bm{k})}$
$\displaystyle=$ $\displaystyle E_{\gamma}(\bm{k})d_{\alpha\gamma}(\bm{k})\;,$
(5b) $\displaystyle c_{\gamma}(\bm{k})$ $\displaystyle=$
$\displaystyle\sum_{\alpha}d_{\alpha\gamma}(\bm{k})a_{\alpha}(\bm{k})\;.$ (5c)

The Dirac point ($\bm{k}_{\rm D}$) is calculated from

$\displaystyle E_{1}(\bm{k}_{\rm D})=E_{2}(\bm{k}_{\rm D})=\epsilon_{\rm
D}\;.$ (6)

The ZGS is obtained when $\epsilon_{\rm D}$ becomes equal to the chemical
potential at $T=0$.

From $E_{\gamma}$, the local density $n_{\alpha}$ including both spin
$\uparrow$ and $\downarrow$ is calculated as

$\displaystyle n_{\alpha}$ $\displaystyle=$
$\displaystyle\frac{2}{N}\sum_{\bm{k}}\left<\hat{n}_{\alpha}(\bm{k})\right>_{H}$
(7) $\displaystyle=$
$\displaystyle\frac{2}{N}\sum_{\bm{k}}\sum_{\gamma}d^{*}_{\alpha\gamma}(\bm{k})d_{\alpha\gamma}(\bm{k})f(E_{\gamma}(\bm{k})-\mu)\;,$

which is determined self-consistently. $n_{\rm A}=n_{\rm A^{\prime}}$ owing to
transfer energies being symmetric with respect to the inversion center between
A and A’ in Fig. 1. In Eq. (7), $f(\epsilon)=1/(\exp[\epsilon/T]+1)$ with $T$
being temperature in the unit of eV and $k_{\rm B}=1$. The chemical potential
$\mu$ is determined from the three-quarter-filled condition, which is given by

$\displaystyle\frac{1}{N}\sum_{\bm{k}}\sum_{\gamma}f(E_{\gamma}(\bm{k})-\mu)=\int_{-\infty}^{\infty}{\rm
d}\omega D(\omega)f(\omega)=3\;,$ (8)

where

$\displaystyle D(\omega)$ $\displaystyle=$
$\displaystyle\frac{1}{N}\sum_{\bm{k}}\sum_{\gamma}\delta(\omega-
E_{\gamma}(\bm{k}))\;.$ (9)

$D(\omega)$ denotes DOS per spin and per unit cell, which satisfies $\int{\rm
d}\omega D(\omega)=4$. Note that $n_{\rm A}+n_{\rm A^{\prime}}+n_{\rm
B}+n_{\rm C}=6$ from Eq. (8). We use $\mu(T)$ at finite $T$ and $\mu$ =
$\mu(0)$ at $T$=0.

In Eq. (1), the third term denotes the harmonic phonon given by $H_{\rm
p}=\sum_{\bm{q}}\omega_{\bm{q}}b_{\bm{q}}^{\dagger}b_{\bm{q}}$ with
$\omega_{\bm{q}}=v_{s}|\bm{q}|$ and $\hbar$ =1, whereas the fourth term is the
e–p interaction expressed as [31]

$H_{\rm
e-p}=\sum_{\bm{k},\gamma}\sum_{\bm{q}}\alpha_{\bm{q}}c_{\gamma}(\bm{k}+\bm{q})^{\dagger}c_{\gamma}(\bm{k})\phi_{\bm{q}}\;,$
(10)

with $\phi_{\bm{q}}=b_{\bm{q}}+b_{-\bm{q}}^{\dagger}$. We introduce a coupling
constant $\lambda=|\alpha_{\bm{q}}|^{2}/\omega_{\bm{q}}$, which becomes
independent of $|\bm{q}|$ for small $|\bm{q}|$. The e–p scattering is
considered within the same band (i.e., intraband) owing to the energy
conservation with $v\gg v_{s}$, where $v\simeq 0.05$ [8] denotes the averaged
velocity of the Dirac cone. The last term of Eq. (1), $H_{\rm imp}$, denotes a
normal impurity scattering, which is introduced to obtain the finite
conductivity and to avoid the infinite conductivity in the presence of only
the e–p interaction [32].

###  2.2 Conductivity

By using the component of the wave function $d_{\alpha\gamma}$ in Eq. (5c), we
calculate the conductivity per spin and per site as[33]

$\displaystyle\sigma_{\nu\nu^{\prime}}(T)$ $\displaystyle=$
$\displaystyle\frac{e^{2}}{\pi\hbar
N}\sum_{\bm{k}}\sum_{\gamma,\gamma^{\prime}}v^{\nu}_{\gamma\gamma^{\prime}}(\bm{k})^{*}v^{\nu^{\prime}}_{\gamma^{\prime}\gamma}(\bm{k})$
$\displaystyle\int_{-\infty}^{\infty}d\epsilon\left(-\frac{\partial
f(\epsilon)}{\partial\epsilon}\right)$ $\displaystyle\times$
$\displaystyle\frac{\Gamma_{\gamma}}{(\epsilon-\xi_{\bm{k}\gamma})^{2}+\Gamma_{\gamma}^{2}}\times\frac{\Gamma_{\gamma^{\prime}}}{(\epsilon-\xi_{\bm{k}\gamma^{\prime}})^{2}+\Gamma_{\gamma^{\prime}}^{2}}\;,$
$\displaystyle v^{\nu}_{\gamma\gamma^{\prime}}(\bm{k})$ $\displaystyle=$
$\displaystyle\sum_{\alpha\beta}d_{\alpha\gamma}(\bm{k})^{*}\frac{\partial
h_{\alpha\beta}}{\partial k_{\nu}}d_{\beta\gamma^{\prime}}(\bm{k})\;,$ (12)

where $\nu=x$ and $y$. $h=2\pi\hbar$ and $e$ denote Planck’s constant and
electric charge, respectively. The quantity $\Gamma_{\gamma}$ denotes the
damping of the electron of the $\gamma$ band given by

$\displaystyle\Gamma_{\gamma}=\Gamma+\Gamma_{\rm ph}^{\gamma}\;,$ (13)

where the first term comes from the impurity scattering and the second term
corresponding to the phonon scattering is given by [27] (Appendix B)

$\displaystyle\Gamma_{\rm ph}^{\gamma}$ $\displaystyle=$ $\displaystyle
C_{0}R\times T|\xi_{\gamma,\bm{k}}|\;,$ (14a) $\displaystyle R$
$\displaystyle=$ $\displaystyle\frac{\lambda}{\lambda_{0}}\;,$ (14b)

with $C_{0}=6.25\lambda_{0}/(2\pi v^{2})$. For $v\simeq 0.05$ and
$\lambda_{0}/2\pi v=0.1$, we obtain $C_{0}\simeq$ 12.5 (eV)-1. $R$ denotes a
normalized e–p coupling constant.

In the following, we denote $\sigma_{x}$, $\sigma_{y}$, and $\sigma_{xy}$
instead of $\sigma_{xx}(T)$, $\sigma_{yy}(T)$, and $\sigma_{xy}(T)$ for
simplicity. In terms of $\sigma_{x}$, $\sigma_{y}$, and $\sigma_{xy}$, the
current $(j_{x},j_{y})$ obtained from a response to an external electric field
$(E_{x},E_{y})$ is written as

$\displaystyle\begin{pmatrix}j_{x}\\\
j_{y}\end{pmatrix}=\begin{pmatrix}\sigma_{x}&\sigma_{xy}\\\
\sigma_{xy}&\sigma_{y}\end{pmatrix}\begin{pmatrix}E_{x}\\\
E_{y}\end{pmatrix}\;.$ (15)

The principal axis of the Dirac cone has an angle $\phi$ measured from the
$k_{y}$ axis, where $-\pi/2<\phi<\pi/2$. When we denote the current and
electric field in this axis direction as $j_{x}^{\prime}$ and
$E_{x}^{\prime}$, we obtain

$\displaystyle\begin{pmatrix}j_{x}^{\prime}\\\
j_{y}^{\prime}\end{pmatrix}=\begin{pmatrix}\sigma_{-}&0\\\
0&\sigma_{+}\end{pmatrix}\begin{pmatrix}E_{x}^{\prime}\\\
E_{y}^{\prime}\end{pmatrix}\;,$ (16a) where
$\displaystyle\begin{pmatrix}j_{x}^{\prime}\\\
j_{y}^{\prime}\end{pmatrix}=\begin{pmatrix}\cos\phi&\sin\phi\\\
-\sin\phi&\cos\phi\end{pmatrix}\begin{pmatrix}j_{x}\\\ j_{y}\end{pmatrix}\;,$
(16b) $\displaystyle\begin{pmatrix}E_{x}^{\prime}\\\
E_{y}^{\prime}\end{pmatrix}=\begin{pmatrix}\cos\phi&\sin\phi\\\
-\sin\phi&\cos\phi\end{pmatrix}\begin{pmatrix}E_{x}\\\ E_{y}\end{pmatrix}\;.$
(16c)

The relations between ($\sigma_{-},\sigma_{+}$, $\phi$) and
$(\sigma_{x},\sigma_{y},\sigma_{xy})$ are

$\displaystyle\sigma_{x}$ $\displaystyle=$
$\displaystyle(\cos^{2}\phi)\;\sigma_{-}+(\sin^{2}\phi)\;\sigma_{+}\;,$ (17a)
$\displaystyle\sigma_{y}$ $\displaystyle=$
$\displaystyle(\sin^{2}\phi)\;\sigma_{-}+(\cos^{2}\phi)\;\sigma_{+}\;,$ (17b)
$\displaystyle\sigma_{xy}$ $\displaystyle=$
$\displaystyle\sin\phi\cos\phi\;(\sigma_{-}-\sigma_{+})\;,$ (17c)

or inversely,

$\displaystyle\tan 2\phi$ $\displaystyle=$
$\displaystyle\frac{2\sigma_{xy}}{\sigma_{x}-\sigma_{y}}\;,$ (18a)
$\displaystyle\sigma_{-}$ $\displaystyle=$
$\displaystyle\frac{1}{2}[\sigma_{x}+\sigma_{y}-\sqrt{(\sigma_{x}-\sigma_{y})^{2}+4\sigma_{xy}^{2}}]\;,$
(18b) $\displaystyle\sigma_{+}$ $\displaystyle=$
$\displaystyle\frac{1}{2}[\sigma_{x}+\sigma_{y}+\sqrt{(\sigma_{x}-\sigma_{y})^{2}+4\sigma_{xy}^{2}}]\;.$

Note that $\sigma_{xy}$ is not a Hall conductivity and
$\sigma_{xy}=\sigma_{yx}$ holds as in Eq. (15) in the case of zero magnetic
field.[34] $\sigma_{xy}$ is finite when $\sigma_{-}\neq\sigma_{+}$. The sign
of $\phi$ is chosen such that $\phi<0$ for $\sigma_{xy}>0$ and $\phi>0$ for
$\sigma_{xy}<0$, where $0<|\phi|<\pi/4$ for $\sigma_{y}>\sigma_{x}$ and
$\pi/4<|\phi|<\pi/2$ for $\sigma_{x}>\sigma_{y}$.

##  3 Results

We calculate the conductivity for the TB model with transfer energies shown in
Fig. 1. The direction of molecular stacking is given by the $y$ axis and that
perpendicular to the stacking is given by the $x$ axis. The NN transfer
energies are examined for both uniaxial pressure and hydrostatic pressure,
while those of NNN sites are added only for hydrostatic pressure.[7] In the
following calculations, the conductivity is normalized by $e^{2}/\hbar$ and
the energy is scaled by eV. First, we examine the case of uniaxial pressure
($P=6$ kbar) without site potential and then the case of hydrostatic pressure
with site potential to comprehend the similarity and dissimilarity.

###  3.1  $\alpha$-(BEDT-TTF)2I3 under uniaxial pressure

For the uniaxial pressure $P$ (kbar), we take only transfer energies with NN
neighbor sites, $t=a_{1},\cdots,b_{4}$ (eV), which are estimated by the
extended Hückel method and interpolation, [6, 5]

$\displaystyle t(P)=t(0)(1+K_{t}P)\;,$ (19)

where $t(0)$ is given by $a_{1}(0),\cdots,b_{4}(0)$ = $-0.028$, $-0.048$,
0.020, 0.123, 0.140, 0.062, and 0.025, and the corresponding $K_{t}$ is
written as 0.089, 0.167, $-0.025$, 0, 0.011, and 0.032, respectively. Since
the ZGS is obtained for $P>3$ kbar, we choose $P$ = 6 kbar as a typical
pressure.

Figure 2(a) shows two bands of $E_{1}(\bm{k})$ and $E_{2}(\bm{k})$, which
touch at Dirac points $\pm\bm{k}_{\rm D}=\pm(0.57,0.30)\pi$ with an energy
$\epsilon_{\rm D}$ = $\mu=0.178$ corresponding to the three-quarter-filled
band. The ranges of the energy of the conduction and valence bands
$E_{1}(\bm{k})$ and $E_{2}(\bm{k})$ are given by
$0<E_{1}(\bm{k})-\epsilon_{\rm D}<0.19$ and $-0.11<E_{2}(\bm{k})-\epsilon_{\rm
D}<0$, respectively. Such ZGS shows the relation $E_{2}(Y)<\epsilon_{\rm
D}<E_{1}(X)<E_{1}(M)$, where $\Gamma$, X, Y, and M are TRIMs given by
$\Gamma=(0,0)\pi$, $X=(1,0)\pi$, $Y=(0,1)\pi$, and $M=(1,1)\pi$. Figure 2(b)
shows contour plots of $E_{1}(\bm{k})-E_{2}(\bm{k})$ as the function of
$\delta\bm{k}=\bm{k}-\bm{k}_{\rm D}$ in a small region around $\bm{k}_{\rm
D}$. The contour lines form almost isotropic circles, suggesting that the
velocity of the Dirac cone is isotropic. Figure 2(c) shows
$E_{1}(\bm{k})-\epsilon_{\rm D}$. The Dirac point is located at $(\delta
k_{x},\delta k_{y})$ = (0,0). This contour suggests a tilted Dirac cone and
shows a slight deviation from the ellipse. It is also found that the cone
shows a slight rotation clockwise from the $k_{x}$ axis, which plays a crucial
role in the transport property, as shown later. Figure 2(d) shows
$E_{2}(\bm{k})-\epsilon_{\rm D}$. The Dirac point is located at (0,0). The
contour of $E_{2}(\bm{k})-\epsilon_{\rm D}$ also shows a tilted Dirac cone and
a slight deviation from the ellipse. We define a phase $\phi_{1}(<0)$ (
$\phi_{2}$) as a tilting angle of $E_{1}(\bm{k})$ ($E_{2}(\bm{k})$) measured
from the $k_{x}$ axis. Since $E_{1}(\bm{k})$ and $E_{2}(\bm{k})$ form a pair
of Dirac cones, $\phi_{2}-\phi_{1}=\pi$ for $\bm{k}$ in the limit of the Dirac
point. The deviation from the limiting value increases with increasing
$|\delta\bm{k}_{\rm D}|$. The tilting parameter, which is a ratio of the
tilting velocity to that of the Dirac cone, is given by $\eta\sim 0.8$. Figure
2(e) shows contour plots of $E_{1}(\bm{k})+E_{2}(\bm{k})-2\epsilon_{\rm D}$.
The bright line, $E_{1}(\bm{k})+E_{2}(\bm{k})-2\epsilon_{\rm D}=0$, is
perpendicular to the tilting direction of the cone.

Figure 2: (Color online) (a) Conduction and valence bands given by
$E_{1}(\bm{k})$ (upper band) and $E_{2}(\bm{k})$ (lower band) at $P$= 6 kbar.
Two bands contact at the Dirac points $\pm\bm{k}_{\rm D}=\pm(0.57,0.30)\pi$
with an energy $\epsilon_{\rm D}=\mu=0.178$. (b) Contour plots of
$E_{1}(\bm{k})-E_{2}(\bm{k})$ with the energy range [0, 0.106], where
$\delta\bm{k}(=\bm{k}-\bm{k}_{\rm D})$. The yellow line (outermost bright
line) corresponds to $E_{1}(\bm{k})-E_{2}(\bm{k})$ = 0.03. (c) Contour plots
of $E_{1}(\bm{k})-\epsilon_{\rm D}$ with the range [0, 0.096], where the Dirac
point exists at (0,0). The yellow line corresponds to
$E_{1}(\bm{k})-\epsilon_{\rm D}$ = 0.005. (d) Contour plots of
$E_{2}(\bm{k})-\epsilon_{\rm D}$ with the range [-0.076, 0]. The Dirac point
exists at (0,0). The outermost bright line corresponds to
$E_{2}(\bm{k})-\epsilon_{\rm D}$ = -0.005. (e) Contour plots of
$E_{1}(\bm{k})+E_{2}(\bm{k})-2\epsilon_{\rm D}$ with the range [-0.062,
0.086]. The bright color line denotes
$E_{1}(\bm{k})+E_{2}(\bm{k})-2\epsilon_{\rm D}=0$.

Figure 3 shows the $T$ dependence of the chemical potential $\mu(T)$ for $P$ =
6 (solid line) and a hydrostatic pressure $P_{\rm hydro}$ (dashed line)
examined in the next subsection, where the corresponding DOS as a function of
$\omega-\mu$ is shown in the inset. With increasing $T$, $\mu$ varies slowly,
suggesting that the $T$ dependence of $\mu$ on $\sigma_{\nu}$ is negligibly
small. In fact, the $T$ dependence of $\mu$ in Eq. (14a) may be ignored for
$\sigma_{\nu}(T)$ at low temperatures. We verified that $\mu(T)$ can be
replaced by $\mu(0)$ in $\xi_{\gamma,\bm{k}}$ for $0<T<0.015$ ($0<T<0.013$),
which is the range of temperatures of the following numerical calculation of
$\sigma_{\nu}$ under the pressure of $P$= 6 kbar ($P_{\rm hydro}$). The DOS
close to the chemical potential shows a clear deviation from the linear
dependence with respect to $\omega-\mu$, where $\mu$ is the chemical potential
at $T$=0. Such asymmetry of DOS is ascribed to the nonlinear spectrum of the
Dirac cone, where $E_{1}(\bm{k})$ ($E_{2}(\bm{k})$) is convex downward
(upward) on the tangent plane at the Dirac point. The slight decrease in the
chemical potential at low temperatures comes from such asymmetry. The increase
in $\mu$ above the minimum occurs since the van Hove singularity at $E_{2}(Y)$
below the chemical potential has a large peak compared with that of
$E_{1}(X)$.

Figure 3: (Color online) Chemical potential ($\mu$) as a function of
temperature ($T$) for $P$ =6 (solid line) and $P=P_{\rm hydro}$ (dashed line).
It is found that $\mu(T)$ for $P$=6 ($P=P_{\rm hydro}$) takes a minimum at
$T\sim 0.004$ ($T\sim 0.006$). The inset denotes the corresponding DOS as a
function of $\omega-\mu=\omega^{\prime}$ for $P$ = 6, where D($-\omega$’)(red
dotted line) is compared with D($\omega$’) showing the deviation from the
linear dependence.

Figure 4 shows the temperature dependence of conductivity of $\sigma_{\nu}$
($\nu=x,y,xy,+,-$) without the e–p interaction, where $P$ = 6 and $\Gamma$ =
0.0005. Note that the solid line was already obtained in Ref.
Katayama2006_cond but $\sigma_{xy}$ was absent since the energy band with the
Dirac cone was found later.[5] We discuss $\sigma_{\nu}(T)$ by dividing
temperatures into two regions: (I) low temperatures with $\sigma_{y}$ being
larger than $\sigma_{x}$ and (II) high temperatures with $\sigma_{x}$ being
larger than $\sigma_{y}$.

First, we examine region (I), where the quantum region is found, leading to
the finite conductivity even at $T$ = 0. In fact, the conductivity in the zero
limit of $T$ is given by $\sigma_{x}(0)\simeq 0.06$ and $\sigma_{y}(0)\simeq
0.09$, which are compared with that of the universal value of
$1/2\pi^{2}=0.051$ obtained for the simple Dirac cone. The slightly larger
value in the present case comes from the tilting of the Dirac cone. [25] An
equality $\sigma_{y}>\sigma_{x}$ at low temperatures is understood as follows.
Figure 2(b) shows a nearly isotropic velocity of the Dirac cone, whereas Figs.
2(c) and 2(d) present a large tilting of the Dirac cone almost along the
$k_{x}$ direction. Our previous calculation of a tilted Dirac cone at $T$ = 0
shows that the conductivity becomes maximum for the direction perpendicular to
a tilting axis. [25] This finding explains $\sigma_{y}>\sigma_{x}$ in the
present case. Such inequality still holds at finite temperatures above
$T>\Gamma$. As a reference, the conductivity of the isotropic Dirac cone
without tilting, $\sigma_{\rm iso}$, is shown by the dot-dashed line for
$0<T<0.005$. With increasing $T$, $\sigma_{\rm iso}$, which is slightly lower
than $\sigma_{x}$ at $T=0$, becomes larger than $\sigma_{x}$. Thus, the effect
of tilting is found in region (I), where $\sigma_{y}>\sigma_{\rm
iso}>\sigma_{x}$. By comparing the solid line ($\Gamma$ = 0.0005) with the
dotted line ($\Gamma$ = 0.001), we found that the effect of tilting is reduced
for large $\Gamma$.

Next, we examine region (II), where $\sigma_{x}>\sigma_{y}$ is found for
$T>0.012$. In region (II) the effect of the tilting is neglected but the
difference between $\sigma_{x}$ and $\sigma_{y}$ still exists owing to the
anisotropy of the velocity of the cone. Since the DOS is averaged owing to a
factor $(-\partial f(\epsilon)/\partial\epsilon)$ in Eq. (LABEL:eq:sigma)
having a broad distribution with a width $\sim T$, this resembles a classical
conductivity with a finite chemical potential being $\sim T$. Thus, we find
$\sigma_{x}>\sigma_{y}$, since the average velocity is larger for the
$x$-direction than for the $y$-direction because the transfer energies are
larger for the $x$-direction than for the $y$-direction. Such crossover occurs
at lower $T$ for larger $\Gamma$. We also note that such crossover occurs at
higher $T$ for $P$ = 8 (not shown here).

Figure 4: (Color online) $T$ dependence of conductivity in the absence of the
e–p interaction at $P$ = 6 with fixed $\Gamma$ = 0.0005, where the solid lines
denote $\sigma_{x}$, $\sigma_{y}$, and $\sigma_{xy}$ and the dashed lines
denote $\sigma_{\pm}$. The dotted lines correspond to $\Gamma=0.001$. As a
reference, the conductivity of the isotropic Dirac cone without tilting,
$\sigma_{\rm iso}$, is shown by the (red) dot-dashed line for $0<T<0.005$.
Principal values of $\sigma_{-}$ and $\sigma_{+}$ are given by Eqs. (18b) and
(LABEL:eq:eq21c), respectively, whereas $\phi$ is given by Eq. (18a). The
inset shows the phase $\phi$, which is an angle of the principal axis of
$\sigma_{-}$ measured from the $k_{x}$ axis.

The dashed lines in Fig. 4 show the $T$ dependence of the principal value
$\sigma_{\pm}$. In the limit of $T=0$, the principal conductivity
$\sigma_{+}$, which gives the maximum of the conductivity, becomes
$\simeq\sigma_{y}$ and is perpendicular to the tilted Dirac cone. Note that
this direction is almost parallel to the line of
$E_{1}(\bm{k})+E_{2}(\bm{k})=0$ (see Fig. 2(e)). For $T\rightarrow 0$, we
obtain $\sigma_{+}(0)\simeq 0.092$ and $\sigma_{-}(0)\simeq 0.0601$, which
give $\sigma_{+}(0)/\sigma_{0}\simeq 1.82$ and $\sigma_{-}(0)/\sigma_{0}\simeq
1.32$ with $\sigma_{0}=1/(2\pi^{2})$. These results are consistent with the
analytical results of the tilted Dirac cone, where $\sigma_{+}/\sigma_{0}$ =
1.84 and $\sigma_{-}/\sigma_{0}$ =1.19 for $\eta$ =084. [25] The behavior of
$\sigma_{\pm}$ being linear in $T$ resembles that of DOS around $\omega=\mu$
in the inset of Fig. 3, which is found for $|\omega-\mu|<0.015$. The inset
shows the $T$ dependence of $\phi$, where $\phi<0$ owing to $\sigma_{xy}>0$
and $\sigma_{y}>\sigma_{x}$ at low temperatures, as seen from Eq. (18a). With
increasing $T$, $\phi$ decreases and becomes smaller than $-\pi/4$ at a
temperature corresponding to $\sigma_{x}=\sigma_{y}$, where the axis close to
$\sigma_{-}$ changes from the $k_{x}$ axis to the $k_{y}$ axis. Thus, the
principal axes rotate clockwise with an angle $\phi(<0)$ (see Eq. (16c)).

Figure 5: (Color online) Conductivity $\sigma_{\nu}$ ($\nu=x,y$, and $xy$) in
the presence of the e–p interaction with $P$ = 6 and $\Gamma$ = 0.0005. The
solid, dashed, and dot-dashed lines denote $\sigma_{x}$, $\sigma_{y}$, and
$\sigma_{xy}$, respectively, where $R$ denotes a normalized e–p coupling
constant defined by the ratio $\lambda/\lambda_{0}$ (Eq. (14b)). $R$ = 1
corresponds to $\lambda=0.1$. At high $T$, $\sigma_{x}$ shows almost constant
behavior and $\sigma_{y}$ shows a broad maximum.  Figure 6: (Color online)
Conductivity $\sigma_{\pm}$ with $P$ = 6 and $\Gamma$ = 0.0005 in the presence
of the e–p interaction. The solid and dashed lines correspond to $\sigma_{+}$
and $\sigma_{-}$, respectively, which show a broad maximum. The inset denotes
the corresponding $\phi$.

As shown in Fig. 4 ( $P$=6), $\sigma_{\nu}$ increases monotonically as a
function of $T$, which is different from the experiment showing nearly
constant behavior at high temperatures. [2] Such an exotic $T$ dependence of
$\sigma_{\nu}$ is examined next by taking account of the e–p interaction,
which is expected to reduce $\sigma_{\nu}$. By using Eq. (23b) and (14a), we
calculate $\sigma_{\nu}$ of Eq. (LABEL:eq:sigma), where $\Gamma$ in the
absence of the e–p interaction is replaced by
$\Gamma_{\gamma}(=\Gamma+\Gamma_{\rm ph}^{\gamma})$. Owing to the $T$
dependence of $\Gamma_{\rm ph}^{\gamma}$, $\Gamma$ is dominant at low $T$
whereas $\Gamma_{\rm ph}^{\gamma}$ is dominant at high $T$. Note that such
crossover with increasing $T$ depends on $R$.

Figure 5 shows the $T$ dependence of $\sigma_{\nu}$ ($\nu=x,y,$ and $xy$) in
the presence of the e–p interaction with some choices of $R$. The effect of
$R\not=0$ appears when $\sigma_{\nu}$ deviates from the $T$-linear behavior.
Compared with $\sigma_{\nu}$ with $R$ = 0 (Fig. 4), $\sigma_{\nu}$ is reduced
noticeably. At temperatures around $T\sim 0.015$, $\sigma_{x}$ becomes nearly
constant, whereas $\sigma_{y}$ takes a broad maximum at lower temperatures.
The quantity $\sigma_{y}-\sigma_{x}$ decreases noticeably at high
temperatures, suggesting that the effect of the tilting of the Dirac cone
decreases with increasing temperature. The $R$ dependence on $\sigma_{xy}$ is
small. The crossover temperature corresponding to $\sigma_{x}=\sigma_{y}$
decreases with increasing $R$, since the reduction of $\sigma_{y}$ is larger
than that of $\sigma_{x}$. The nearly constant behavior of $\sigma_{x}$ is
understood as follows. With increasing $T$, $\sigma_{\nu}$ without the e–p
interaction ($R=0$) increases linearly owing to the DOS obtained from the
Dirac cone. Such a linear increase is suppressed for $R\not=0$, since the
noticeable effect of the acoustic phonon emerges at finite temperatures. In
fact, the electron is scattered by both normal impurity ($\Gamma$) and the e–p
interaction ($\Gamma_{\rm ph}^{\gamma}$), and the latter becomes dominant at
high temperatures as seen from Eq. (14a). However, compared with a case of the
conventional metal with a Fermi surface, the effect of the e–p scattering in
the case of the Dirac cone close to the three-quarter-filled band is strongly
reduced owing to a constraint by the energy-momentum conservation.[27] This
could possibly give rise to a nearly constant behavior or a broad maximum in
$\sigma_{\nu}$, owing to a competition between the enhancement by DOS of the
Dirac cone and the suppression by the e–p interaction.

Figure 6 shows the $T$ dependence of the principal values $\sigma_{\pm}$
corresponding to Fig. 5. Since $\sigma_{+}$ and $\sigma_{-}$ give the upper
and lower bounds of the conductivity, these quantities are convenient for
comprehending the results of experiments even when the precise relationship
between the crystal axis and the direction of the applied electric field is
unknown. In the inset, the angle between the $k_{x}$ axis and the tilting axis
is shown, where $\phi<0$ for $\sigma_{xy}>0$ and $\sigma_{y}>\sigma_{x}$. The
angle $|\phi|$ increases and exceeds $\pi/4$ at a certain temperature with
$\sigma_{x}=\sigma_{y}$, implying that the axis of $\sigma_{-}$ becomes closer
to the $k_{y}$ axis at high temperatures. Such a behavior can be understood on
the basis of a simplified model, [27] where
$\sigma_{+}\sim\sigma_{-}\sim\sigma$. Note that $\Gamma_{\rm ph}^{\gamma}$ is
obtained in Eq. (14a) and $\sigma\simeq a_{\nu}^{\prime}10^{3}T/\Gamma$ with
$a_{\nu}^{\prime}=o(0.1)$ without e–p interaction. By taking $\Gamma$ replaced
by $\Gamma+\Gamma_{\rm ph}^{\gamma}$ and employing an idea
$<|\xi_{\gamma,\bm{k}}|>\sim T$ with $<>$ being an average value in the
summation of Eq. (24), we obtain

$\displaystyle\sigma\simeq\frac{a_{\nu}^{\prime}10^{3}T}{1+C_{0}RT^{2}/\Gamma}\;,$
(20)

with $C_{0}$ = 12.5 and $\Gamma$ = 0.0005. From Eq. (20), it is found that a
maximum of $\sigma$ as a function of $T$ is obtained by a competition of the
DOS (the numerator) and the e–p interaction (the denominator). Equation (20)
suggests that $\sigma$ decreases with increasing $R$.

###  3.2  $\alpha$-(BEDT-TTF)2I3 under hydrostatic pressures

We examine Dirac electrons under hydrostatic pressure using the TB model with
NN ($a_{1},\dots,b_{4}$) and NNN ($a_{1d},\cdots,a_{4d}$) transfer energies
given by[7] $a_{1}=-0.0267$, $a_{2}=-0.0511$, $a_{3}=0.0323$, $b_{1}=0.1241$,
$b_{2}=0.1296$, $b_{3}=0.0513$, $b_{4}=0.0152$, $a_{1d}=0.0119$,
$a_{3d}=0.0046$, and $a_{4d}=0.0060$. The site potentials given by Eqs. (22)
and (22) [8] are estimated as $\tilde{V}_{B}=0.0511$, $\tilde{V}_{C}=0.0032$
(Appendix A).

Band energies $E_{j}(\bm{k})$ under hydrostatic pressure, which have Dirac
cones similar to Fig. 2, are as follows. Figure 7(a) shows the conduction and
valence bands ($0<E_{1}(\bm{k})-\epsilon_{\rm D}<0.15$ and
$-0.09<E_{2}(\bm{k})-\epsilon_{\rm D}<0$) in the first Brillouin zone, which
touch at the Dirac points $\pm\bm{k}_{\rm D}=\pm(0.69,0.44)\pi$ with
$\epsilon_{\rm D}$ =0.172. Compared with Fig. 2(a), the band width of both
$E_{1}(\bm{k})$ and $E_{2}(\bm{k})$ is slightly large and the Dirac points
move away from the $\Gamma$ point. In Fig. 7(b), the energy difference between
$E_{1}$ and $E_{2}$ is shown. The contour lines form ellipsoids, suggesting
that the anisotropy of the velocity of the Dirac cone is larger than that of
Fig. 2(b), e.g., an ellipse with the ratio of the major to minor axes being
$\simeq 1.2$. Figure 7(c) shows that the tilting angle $\phi(>0)$ of
$E_{1}(\bm{k})-\epsilon_{\rm D}$ measured from the $k_{x}$ axis has the sign
opposite to that of Fig. 2(c). As shown later, such difference in the sign
results in the difference in the current direction, i.e., the rotation with
respect to the applied electric field of $E_{x}$ or $E_{y}$. Figure 7(d) shows
$E_{2}(\bm{k})-\epsilon_{\rm D}$ forming a pair of Dirac cones with that of
Fig. 7(c), where the deviation of the tilting axis from the $k_{x}$ axis is
opposite to that of Fig. 2(d). The tilting direction is understood from Fig.
7(e), where the energy of the Dirac point is given by $E_{1}(\bm{k}_{\rm
D})-\epsilon_{\rm D}=E_{2}(\bm{k}_{\rm D})-\epsilon_{\rm D}=0$ and the tilting
direction of the cone is almost perpendicular to the zero line
[$E_{1}(\bm{k})+E_{2}(\bm{k})-2\epsilon_{\rm D}=0$]. By comparing Figs. 2 with
Figs. 7, we find the opposite rotation of the tilted cone.

Figure 7: (Color online) Energy bands of $E_{1}(\bm{k})$ and $E_{2}(\bm{k})$
under hydrostatic pressure, where the Dirac points are given by
$\pm\bm{k}_{\rm D}=\pm(0.69,0.44)\pi$ and $\mu=\epsilon_{\rm D}$ =0.172. (a)
Conduction and valence bands given by $E_{1}(\bm{k})$ (upper band) and
$E_{2}(\bm{k})$ (lower band). (b) Contour plots of
$E_{1}(\bm{k})-E_{2}(\bm{k})$ as a function of
$\delta\bm{k}=\bm{k}-\bm{k}_{\rm D}$ with the energy range [0, 0.09]. The
yellow line (outermost bright line) corresponds to
$E_{1}(\bm{k})-E_{2}(\bm{k})$ = 0.03. (c) Contour plots of
$E_{1}(\bm{k})-\epsilon_{\rm D}$ as a function of $\delta\bm{k}$ with the
range [0, 0.066]. The yellow line corresponds to $E_{1}(\bm{k})-\epsilon_{\rm
D}$ = 0.005. (d) Contour plots of $E_{2}(\bm{k})-\epsilon_{\rm D}$ as a
function of $\delta\bm{k}$ with the range [-0.056, 0]. The outermost colored
region corresponds to $E_{1}(\bm{k})-\epsilon_{\rm D}=-0.005$. (e) Contour
plots of $E_{1}(\bm{k})+E_{2}(\bm{k})-2\epsilon_{\rm D}$ with the range
[-0.045, 0.045]. The bright color line denotes
$E_{1}(\bm{k})+E_{2}(\bm{k})-2\epsilon_{\rm D}$ = 0.

Now we examine the electric conductivity. The parameters for $P_{\rm hydro}$
are taken as $\Gamma$ = 0.0005 and $\tilde{V}_{B}$= 0.0511, $\tilde{V}_{C}$ =
0.0032 for the site potentials, which are obtained from Eqs. (22) and (22).
Figure 8 shows the conductivity $\sigma_{\nu}$ without the e–p interaction,
where the solid line (dashed line) denotes $\sigma_{x}$, $\sigma_{y}$, and
$\sigma_{xy}$ ($\sigma_{\pm}$). Compared with the case of the uniaxial
pressure (Fig. 4), the sign with $\sigma_{xy}<0$ is opposite, as also seen
from a comparison of Figs. 7(e) and 2(e), showing that the angle of the zero
line $E_{1}+E_{2}-2\epsilon_{\rm D}=0$ from the $k_{y}$ axis is opposite. This
difference originates from the transfer energies. A noticeable difference
between hydrostatic pressure and uniaxial pressure exists in $\sigma_{\nu}$ at
low temperatures, where $\sigma_{y}>\sigma_{x}$ in Fig. 4 and
$\sigma_{y}\simeq\sigma_{x}$ in Fig. 8. The origin of
$\sigma_{y}\simeq\sigma_{x}$ is understood as follows from Figs. 7(b)–7(d).
Figure 7(b) shows an anisotropy of the contour implying that the ratio of the
maximum to minimum velocities of the Dirac cone is $\simeq 1.2$, which results
in the ratio of the maximum and minimum conductivities being $\simeq 1.4$.
[25] On the other hand, the tilted cone gives $\sigma_{y}/\sigma_{x}\simeq
1.4$ for the tilting parameter $\eta\sim 0.9$. [25] The competition of these
two ingredients provides an almost isotropic conductivity at low temperatures.
The conductivity at high temperatures is given by the classical region, where
the effect of the tilting can be neglected and the larger transfer energies
along the $x$ direction gives a larger velocity for the $k_{x}$ axis as seen
from Fig. 7(b).

Figure 8 also shows the principal value $\sigma_{\pm}$ of the conductivity,
where $\phi(>0)$ in the inset denotes an angle between the $\sigma_{-}$ axis
and the $k_{x}$ axis. With increasing $T$, the difference between $\sigma_{+}$
and $\sigma_{-}$ increases but is smaller than that shown in Fig. 4. This
comes from a difference at low temperatures, where $\sigma_{y}-\sigma_{x}$ for
the uniaxial pressure is larger than that for the hydrostatic pressure.

We discuss $\sigma_{\nu}$ with transfer energies of hydrostatic pressure but
without site potentials, which is shown by the dot-dashed line in Fig. 8. In
this case, $\sigma_{y}\simeq\sigma_{x}$ but $\sigma_{xy}>0$. This comes from
the fact that the ellipse corresponding to Fig. 7(b) is slightly rotated
anticlockwise (not shown here). Thus, it turns out that
$\sigma_{y}\simeq\sigma_{x}$ comes from transfer energies of hydrostatic
pressure and that $\sigma_{xy}>0$ is attributable to the absence of the
interaction.

Figure 8: (Color online) Conductivity $\sigma_{\nu}$ with $\nu=x,y$, and $xy$
(solid line), and $\nu=\pm$ (dashed line) under hydrostatic pressure. Note
that $\sigma_{xy}<0$ and $\phi/\pi(>0)$ (inset) from Eq. (18a). The dot-dashed
line denotes $\sigma_{\nu}$ without the site potential, where $\mu(0)$ = 0.168

Now, we examine the effect of the e–p scattering on the $T$ dependence of
conductivity with some choices of $R$, Eq. (14b), where the effect of the site
potential is included in the energy band. Figure 9 shows $\sigma_{\nu}$ with
the fixed $R$ = 0.5 and 1, in which sufficient suppression of $\sigma_{\nu}$
by $R$ is seen compared with that of Fig. 8 ($R$ = 0). The crossover
temperature, where $\sigma_{y}>\sigma_{x}$ at low temperatures and
$\sigma_{x}>\sigma_{y}$ at high temperatures, decreases with increasing $R$.
The increase toward a constant behavior at high $T$ is seen for $\sigma_{x}$,
whereas a broad maximum is found for $\sigma_{y}$. $\sigma_{xy}$ as a function
of $T$ takes a minimum followed by the change of the sign for large R.
Compared with Fig. 4, the crossover temperature of $\sigma_{y}=\sigma_{x}$ is
lower for the hydrostatic pressure, where $\sigma_{y}$ is sufficiently reduced
and $\sigma_{x}$ is less reduced.

Figure 9: (Color online) Conductivity $\sigma_{\nu}$ ($\nu$ = $x,y,$ and $xy$)
under hydrostatic pressure in the presence of the e–p interaction with choices
of $R$ = 0.5 and 1.

##  4 Summary and Discussion

We calculated the $T$ dependence of conductivity $\sigma_{\nu}$ ($\nu=x$ and
$y$) of Dirac electrons in the organic conductor $\alpha$-(BEDT-TTF)2I3 under
both uniaxial and hydrostatic pressures to find that the anisotropic
conductivity exhibits a crossover from region (I) including the quantum regime
at low temperatures to region (II) showing the classical regime at high
temperatures. The former with $\sigma_{y}>\sigma_{x}$ comes from the tilted
Dirac cone, whereas the latter with $\sigma_{x}>\sigma_{y}$ originates from
the anisotropy of the velocity of the cone, and a nearly constant conductivity
comes from the phonon scattering on the Dirac cone. The presence of the off-
diagonal component ($\sigma_{xy}$), which is associated with the deviation of
the tilting axis of the Dirac cone from the $k_{x}$ or $k_{y}$ axis, results
in the principal axis with clockwise or anticlockwise rotation depending on
the sign of $\sigma_{xy}$.

Finally, we compare our result with that of the experiment. The temperature
dependence of resistance (corresponding to the inverse of the conductivity)
under hydrostatic pressure shows a nearly constant behavior at high
temperatures and a minimum at low temperatures, whereas the minimum is
invisible under uniaxial pressure.[2] Our results show the broad maximum for
$\sigma_{y}$ and monotonic variation for $\sigma_{x}$. Although we obtain the
qualitative coincidence, the details of the correspondence between them, e.g.,
the direction of measurement, are needed for the quantitative comparison. It
also remains a future problem to clarify if the conductivity under pressure
suggests the presence of site potentials. The validity of our present
calculation may be examined by the measurement of the deviation angle of the
principal axis, i.e., clockwise (anticlockwise).

###### Acknowledgements.

We thank R. Kato for useful discussions on the effect of the e–p interaction.
This work was supported by Grants-in-Aid for Scientific Research from the
Japan Society for the Promotion of Science (Grant No. JP18H01162), and by JST-
Mirai Program (Grant No. JPMJMI19A1).

##  Appendix A Matrix element of H

The quantity $h_{i,j}$, which denotes the matrix element of $H$ and is the
same as that of Ref.Katayama_EPJ, is given as

$\displaystyle h_{12}(\bm{k})$ $\displaystyle=$ $\displaystyle
a_{3}+a_{2}Y\;,$ (21a) $\displaystyle h_{13}(\bm{k})$ $\displaystyle=$
$\displaystyle b_{3}+b_{2}X\;,$ (21b) $\displaystyle h_{14}(\bm{k})$
$\displaystyle=$ $\displaystyle b_{4}Y+b_{1}XY\;,$ (21c) $\displaystyle
h_{23}(\bm{k})$ $\displaystyle=$ $\displaystyle b_{2}+b_{3}X\;,$ (21d)
$\displaystyle h_{24}(\bm{k})$ $\displaystyle=$ $\displaystyle
b_{1}+b_{4}X\;,$ (21e) $\displaystyle h_{34}(\bm{k})$ $\displaystyle=$
$\displaystyle 2a_{1}\;,$ (21f) $\displaystyle h_{11}(\bm{k})$
$\displaystyle=$ $\displaystyle t_{22}(\bm{k})=a_{1d}(Y+\bar{Y})\;,$ (21g)
$\displaystyle h_{33}(\bm{k})$ $\displaystyle=$ $\displaystyle
a_{3d}(Y+\bar{Y})+\tilde{V}_{\rm B}\;,$ (21h) $\displaystyle h_{44}(\bm{k})$
$\displaystyle=$ $\displaystyle a_{4d}(Y+\bar{Y})+\tilde{V}_{\rm C}\;,$ (21i)

and $h_{ij}(\bm{k})=h_{ji}^{*}(\bm{k})$, where $X=\exp[ikx]=\bar{X}^{*}$ and
$Y=\exp[iky]=\bar{Y}^{*}$,

The matrix elements of $H_{1}$ are $\tilde{V}_{B}$ for $\alpha=\beta=3$,
$\tilde{V}_{C}$ for $\alpha=\beta=4$, and zero otherwise. The site potentials
$\tilde{V}_{\rm B}$ and $\tilde{V}_{\rm C}$ in $H_{1}$ are given by the
mean–field of short-range repulsive interactions, [8]

$\displaystyle\tilde{V}_{\rm B}$ $\displaystyle=$
$\displaystyle(n_{B}-n_{A})U/2$ $\displaystyle+$ $\displaystyle
2V_{a}(n_{C}-n_{A})+2V_{b}(2n_{A}-n_{B}-n_{C})\;,$
$\displaystyle\tilde{V}_{\rm C}$ $\displaystyle=$
$\displaystyle(n_{C}-n_{A})U/2$ $\displaystyle+$ $\displaystyle
2V_{a}(n_{B}-n_{A})+2V_{b}(2n_{A}-n_{B}-n_{C})\;,$

where $U$ is the on-site repulsive interaction and $V_{a}$ ($V_{b}$) denotes
the nearest neighbor interaction along the $y$ ($x$) axis. $n_{\alpha}$
denotes a local density corresponding to an electron number per unit cell at
the $\alpha$ site and is determined self-consistently. In the present case,
the site potentials given by Eqs. (22) and (22) are estimated as
$\tilde{V}_{B}=0.0511$ and $\tilde{V}_{C}=0.0032$, respectively, for $U$= 0.4,
$V_{a}$=0.17, and $V_{b}$=0.05, where $n_{A}=n_{A^{\prime}}$ = 1.46, $n_{B}$ =
1.37, and $n_{C}$ = 1.71, and $\mu$ = 0.172 at $T$=0 .

##  Appendix B Damping by phonon scattering

The damping of electrons of the $\gamma$ band, which is defined by
$\Gamma_{\gamma}$, is obtained from the electron Green function[35] expressed
as

$\displaystyle G_{\gamma}(\bm{k},i\omega_{n})^{-1}$ $\displaystyle=$
$\displaystyle i\omega_{n}-E_{\gamma,\bm{k}}+\mu+i\Gamma_{\gamma}\;,$ (23a)
$\displaystyle\Gamma_{\gamma}$ $\displaystyle=$
$\displaystyle\Gamma+\Gamma_{\rm ph}^{\gamma}\;,$ (23b)

where $\Gamma_{\rm ph}^{\gamma}=-{\rm
Im}\Sigma_{\gamma}(\bm{k},E_{\gamma,\bm{k}}-\mu)$ with
$\Sigma_{\gamma}(\bm{k},E_{\gamma,\bm{k}}-\mu)$ being a self-energy given by
the e–p interaction. The real part of the self-energy can be neglected for
doping at low concentrations. [27] The quantity $\Gamma$ comes from another
self-energy by the impurity scattering. Note that $\Gamma_{\rm ph}^{\gamma}$
does not depend on $\Gamma$, and that the ratio $\Gamma_{\rm
ph}^{\gamma}/\Gamma$ is crucial to the determination of the $T$ dependence of
the conductivity. The quantity
$\Sigma_{\gamma}(\bm{k},\omega)=\Sigma_{\gamma}(\bm{k},i\omega_{n})$ with
$i\omega_{n}\rightarrow\omega+0$ is estimated as [35]

$\displaystyle\Sigma_{\gamma}(\bm{k},i\omega_{n})=T\sum_{m}\sum_{\bm{q}}\;|\alpha_{q}|^{2}$
$\displaystyle\times\frac{1}{i\omega_{n+m}-\xi_{\gamma,\bm{k}+\bm{q}}}\times\frac{2\omega_{\bm{q}}}{\omega_{m}^{2}+\omega_{\bm{q}}^{2}}\;,$
(24)

which is a product of electron and phonon Green functions.
$\omega_{n}=(2n+1)\pi T$, $\omega_{m}=2\pi mT$ with $n$ and $m$ being
integers. $\xi_{\gamma,\bm{k}}=E_{\gamma,\bm{k}}-\mu$. Applying the previous
result,[27] we obtain Eqs. (14a) and (14b).

## References

  * [1] K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, M. I. Katsnelson, I. V. Grigorieva, S. V. Dubonos, and A. A. Firsov, Nature 438, 197 (2005).
  * [2] K. Kajita, Y. Nishio, N. Tajima, Y. Suzumura, and A. Kobayashi, J. Phys. Soc. Jpn. 83, 072002 (2014).
  * [3] T. Mori, A. Kobayashi, Y. Sasaki, H. Kobayashi, G. Saito, and H. Inokuchi, Chem. Lett. 13, 957 (1984).
  * [4] A. Kobayashi, S. Katayama, K. Noguchi, and Y. Suzumura, J. Phys. Soc. Jpn. 73, 3135 (2004).
  * [5] S. Katayama, A. Kobayashi, and Y. Suzumura, J. Phys. Soc. Jpn. 75, 054705 (2006).
  * [6] R. Kondo, S. Kagoshima, and J. Harada, Rev. Sci. Instrum. 76, 093902 (2005).
  * [7] H. Kino and T. Miyazaki, J. Phys. Soc. Jpn. 75, 034704 (2006).
  * [8] S. Katayama, A. Kobayashi, and Y. Suzumura, Eur. Phys. J. B 67, 139 (2009).
  * [9] M. Inokuchi, H. Tajima, A. Kobayashi, H. Kuroda, R. Kato, T. Naito, and H. Kobayashi, Synth. Met. 56, 2495 (1993).
  * [10] M. Inokuchi, H. Tajima, A. Kobayashi, T. Ohta, H. Kuroda, R. Kato, T. Naito, and H. Kobayashi, Bull. Chem. Soc. Jpn. 68, 547 (1995)
  * [11] R. Kondo, S. Kagoshima, N. Tajima, and R. Kato, J. Phys. Soc. Jpn. 78, 114714 (2009).
  * [12] T. Morinari and Y. Suzumura, J. Phys. Soc. Jpn. 83, 094701 (2014).
  * [13] T. Naito, R. Doi, and Y. Suzumura, J. Phys. Soc. Jpn. 89, 023701 (2020).
  * [14] K. Kajita, T. Ojiro, H. Fujii, Y. Nishio, H. Kobayashi, A. Kobayashi, and R. Kato, J. Phys. Soc. Jpn. 61, 23 (1992).
  * [15] N. Tajima, M. Tamura, Y. Nishio, K. Kajita, and Y. Iye, J. Phys. Soc. Jpn. 69, 543 (2000).
  * [16] N. Tajima, A. Ebina-Tajima, M. Tamura, Y. Nishio, and K. Kajita, J. Phys. Soc. Jpn. 71, 1832 (2002).
  * [17] N. Tajima, S. Sugawara, M. Tamura, R. Kato, Y. Nishio, and K. Kajita, EPL 80, 47002 (2007).
  * [18] D. Liu, K. Ishikawa, R. Takehara, K. Miyagawa, M. Tanuma, and K. Kanoda, Phys. Rev. Lett. 116, 226401 (2016).
  * [19] R. Kato, H. Cui, T. Tsumuraya, T. Miyazaki, and Y. Suzumura, J. Am. Chem. Soc. 139, 1770 (2017).
  * [20] R. Kato and Y. Suzumura, J. Phys. Soc. Jpn. 86, 064705 (2017).
  * [21] Y. Suzumura, H. B. Cui, and R. Kato, J. Phys. Soc. Jpn. 87, 084702 (2018).
  * [22] B. Zhou, S. Ishibashi, T. Ishii, T. Sekine, R. Takehara, K. Miyagawa, K. Kanoda, E. Nishibori, and A. Kobayashi, Chem. Commun. 55, 3327 (2019).
  * [23] R. Kato and Y. Suzumura, J. Phys. Soc. Jpn. 89, 044713 (2020).
  * [24] N. H. Shon and T. Ando, J. Phys. Soc. Jpn. 67, 2421 (1998).
  * [25] Y. Suzumura, I. Proskurin, and M. Ogata J. Phys. Soc. Jpn. 83, 023701 (2014).
  * [26] N. M. R. Peres, F. Guinea, and A. H. Castro Neto, Phys. Rev. B 83,125411 (2006).
  * [27] Y. Suzumura and M. Ogata, Phys. Rev. B 98, 161205 (2018).
  * [28] A. Kobayashi, S. Katayama, Y. Suzumura, and H. Fukuyama, J. Phys. Soc. Jpn. 76, 034711 (2007).
  * [29] M. O. Goerbig, J.-N. Fuchs, G. Montambaux, and F. Pi$\acute{\rm e}$chon, Phys. Rev. B 78, 045415 (2008).
  * [30] A. Kobayashi, Y. Suzumura, and H. Fukuyama, J. Phys. Soc. Jpn. 77, 064718 (2008).
  * [31] H. Fröhlich, Proc. Phys. Soc. A 223, 296 (1954).
  * [32] T. Holstein, Ann. Phys. 29, 410 (1964).
  * [33] S. Katayama, A. Kobayashi, and Y. Suzumura, J. Phys. Soc. Jpn. 75, 023708 (2006).
  * [34] For example, see Eq. (6.9) in R. Kubo, J. Phys. Soc. Jpn. 12, 570 (1957).
  * [35] A. A. Abrikosov, L. P. Gorkov, and I. E. Dzyaloshinskii, Methods of Quantum Field Theory in Statistical Physics (Prentice-Hall, Englewood Cliffs, NJ, 1963).


# All-passive multiple-place optical phase noise cancellation

Liang Hu  State Key Laboratory of Advanced Optical Communication Systems and
Networks, Department of Electronic Engineering, Shanghai Jiao Tong University,
Shanghai 200240, China  Corresponding author<EMAIL_ADDRESS>Ruimin Xue
State Key Laboratory of Advanced Optical Communication Systems and Networks,
Department of Electronic Engineering, Shanghai Jiao Tong University, Shanghai
200240, China  Xueyang Tian  State Key Laboratory of Advanced Optical
Communication Systems and Networks, Department of Electronic Engineering,
Shanghai Jiao Tong University, Shanghai 200240, China  Guiling Wu  State Key
Laboratory of Advanced Optical Communication Systems and Networks, Department
of Electronic Engineering, Shanghai Jiao Tong University, Shanghai 200240,
China  Shanghai Key Laboratory of Navigation and Location-Based Services,
Shanghai 200240, China  Jianping Chen  State Key Laboratory of Advanced
Optical Communication Systems and Networks, Department of Electronic
Engineering, Shanghai Jiao Tong University, Shanghai 200240, China  Shanghai
Key Laboratory of Navigation and Location-Based Services, Shanghai 200240,
China

###### Abstract

We report on the realization of delivering coherent optical frequency to
multiple places based on passive phase noise cancellation over a bus topology
fiber network. This technique mitigates any active servo controller on the
main fiber link and at arbitrary access places as opposed to the conventional
technique, in which an active phase compensation circuit has to be adopted to
stabilize the main fiber link. Although the residual fiber phase noise power
spectral density (PSD) in the proposed technique turns out to be a factor of 7
higher than that of in the conventional multiple-access technique when the
access place is close to the end of the fiber link, it could largely suppress
the phase noise introduced by the servo bumps, improve the response speed and
phase recovery time, and minimize hardware overhead in systems with many
stations and connections without the need of the active servo circuits
including phase discriminators and active compensators. The proposed technique
could considerably simplify future efforts to make precise optical frequency
signals available to many users, as required by some large-scale science
experiments.

††journal: ol

Atomic optical clocks have rapidly grown in the last decade and are proving to
be a powerful tool for investigation of fundamental and applied physics [1, 2,
3, 4, 5]. Precision clock networks are of particular interest in precision
measurements and fundamental physics tests, such as general relativity,
temporal variation of the fundamental constant [6], searching for dark matter,
gravitational waves and physics beyond the Standard Model [7, 8, 9, 10], as
well as providing innovative quantum technologies for other branches of
science [11]. Most of the applications mentioned above require high-precision
clock networks with multiple places over fiber links. To achieve this aim,
active compensation schemes as first demonstrated in 1994 by Ma et al. have
been proposed to cancel the fiber-induced phase drift and implement highly
stable optical frequency distribution [12], which generally utilizes the phase
error from a round-trip probe signal to achieve the feedback control of a
voltage-controlled oscillator (VCO) via a servo controller [12]. One
intriguing question is how to distribute a reference optical frequency to many
users simultaneously in a cost-effective and robust way.

Figure 1: Schematic diagram of our multi-access optical frequency
dissemination with passive phase correction over a fiber link. We evaluate the
system performance by measuring the beatnote between the local light reference
and the output of the access place. AOM: acousto-optic modulator, FRM: Faraday
mirror, FM: frequency mixer, EDFA: erbium-doped fiber amplifier, OC: optical
coupler, FD: frequency divider, PD: photo-detector, PC: power combiner, DDS:
direct digital synthesizer. The black solid and dashed arrows represent the
directions of the forward and backward lights, respectively. The blue solid
arrows are the direction of the RF signals.

In the past years, many works have demonstrated that coherence optical phase
can be delivered to multiple users with the help of servo controllers. Grosche
et al. have first proposed and demonstrated to extract the ultrastable signal
for multiple users along the active phase stabilized main link [13], enabling
to tap the fiber anywhere with the same precision level as that achieved at
the main link output. Alternatively, a branching optical fiber network with
phase noise correction at each output end has been proposed and experimentally
demonstrated [14, 15]. All these existing demonstrators for the multiple-
access applications have to adopt at least one active servo controller for the
main fiber link [16, 13, 17, 18] or for each branch fiber link [14, 15, 19].
Consequently, the overall performance of the above mentioned techniques is
mainly dependent on whether the servo controller’s parameters set properly
[20]. For example, there will be occasional interruptions where a phase lock
temporarily loses lock, corresponding to unpredictable jumps in the optical
signals. These interruptions will reduce the effective averaging time to bias
the measurement [12]. Moreover, the effect of the servo bumps on the main
fiber link will pass to each access place and corrupt the spectral purity
received by each place.

In our previous work, we have demonstrated optical frequency transfer with
passive phase stabilization over a bus or a ring fiber network [21, 19]. The
technique possesses the advantages of an unlimited compensation precision and
a fast compensation speed and free from the effect of servo bumps on the
spectral purity [21, 19]. However, the feasibility and adaptability of the
passive phase noise cancellation technique for the multiple-access application
needs to be theoretically and experimentally studied. Therefore, the
investigation of the novel multiple-access optical frequency transfer scheme
which implement alternative phase correction techniques other than the more
commonly used conventional, active, phase correction technique is seeing
increasing demand. Improving the robustness and sensitivity of multiple-access
optical frequency transfer devices, as well as understanding and
characterizing the limitations of the novel multiple-access optical frequency
transfer scheme with passive phase correction is an important step towards the
goal of heralding a new generation of viable precision optical frequency
transfer devices. Theses devices could be employed in the above mentioned
applications [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11].

In this article, we extend the study of the optical frequency transfer
technique based on our previous passive phase correction technique [21, 19],
demonstrating its application as a multiple-access optical frequency transfer
technique. We experimentally study the optical frequency transfer stabilities,
phase noise power spectral densities (PSDs) and accuracies for the two access
places at the most symmetric 50/50 km ($L_{a}/L_{b}$) one and the relative
asymmetric 70/30 km one, over a total fiber link of $L=100$ km. In comparison
with our previous work [21, 19, 22], in which we have demonstrated multiple-
place optical frequency transfer over the star and ring fiber networks with
the passive phase noise cancellation technique, here we demonstrated that a
high performance multiple-place optical frequency transfer over the widely
adopted bus fiber topology with all-passive phase noise cancellation at the
access places and the main fiber link. The proposed technique could increase
the adaptability to incorporate the optical frequency transfer technique into
any existing communication networks with different topologies.

Figure 1 shows the schematic diagram of multiple-place optical frequency
transfer with a simple extraction along the passively stabilized main link.
The main optical link aims at regenerating a coherent optical phase at the
output end of the fiber [21] and at any places along the fiber as the input
end of the fiber. The principle of passively stabilizing the main fiber link
can be found in [21, 19]. In brief, an ultrastable reference $\nu$ at the
local site is upshifted by frequency of $\omega_{l}$ with an acousto-optic
modulator (AOM) denoted as AOM1 and then injects into the fiber. At the
output, the light is downshifted by frequency of $\omega_{r}$ with the AOM2
(here $\omega_{l}>\omega_{r}$), and part of light is sent back with a Faraday
mirror (FRM). The round-trip signal is mixed with the input ultrastable laser
using an interferometer consisting of an optical coupler (OC) and another FRM.
The beatnote radio frequency (RF) signal is twice the sum of the AOMs’
frequencies and exhibits the round-trip phase noise, $2\phi_{p}$, with
$\phi_{p}=\phi_{a}+\phi_{b}$, where $\phi_{a}$ and $\phi_{b}$ are the noise of
the fiber spans of length of $L_{a}$ and $L_{b}$, respectively. Afterwards, we
divide the beatnote with a factor of 2 and then mix the divided signal with an
assistant frequency of $\omega_{s}$. The lower sideband of the mixed signal,
$\omega_{s}-\omega_{l}+\omega_{r}$, is used. Afterwards, the mixed signal
together with $\omega_{l}$ is fed into the electrical port of the AOM1.
Finally at the output end, the phase fluctuations of the optical frequency
$\nu+\omega_{s}-\omega_{l}$ are automatically cancelled.

At each access place, the configuration is similar with the one presented in
[13], at a distance $L_{a}$ from the input end and $L_{b}$ from the output
end, a $2\times 2$ optical coupler is adopted to extract both the forward and
backward signals from the main fiber. The forward signal can be expressed as,

$\begin{split}S_{F}(\omega)\propto&\cos((\nu+\omega_{l})t+\phi_{a})\\\
&+\cos((\nu+\omega_{s}-\omega_{l}+\omega_{r})t-\phi_{p}+\phi_{a}).\end{split}$
(1)

Similarly, the backward signal has a form of

$\begin{split}S_{B}(\omega)\propto&\cos((\nu-2\omega_{r}+\omega_{l})t+\phi_{p}+\phi_{b})\\\
&+\cos((\nu+\omega_{s}-\omega_{l}-\omega_{r})t+\phi_{b}).\end{split}$ (2)

The beatnote frequencies between $\nu+\omega_{l}$ and
$\nu-2\omega_{r}+\omega_{l}$ and between
$\nu+\omega_{s}-\omega_{l}+\omega_{r}$ and
$\nu+\omega_{s}-\omega_{l}-\omega_{r}$, result in the same frequency of
$2\omega_{r}$ with the phase fluctuation $2\phi_{b}$. The signal is divided by
a factor of 2, filtered, and drives an AOM (AOM3) in order to correct for the
frequency and phase fluctuations of the forward signal. The forward signal
after passing through the AOM3, is thus downshifted to,

$\begin{split}S_{F^{\prime}}(\omega)\propto&\cos((\nu-\omega_{r}+\omega_{l})t+\phi_{p})\\\
&+\cos((\nu+\omega_{s}-\omega_{l})t).\end{split}$ (3)

A similar phase noise cancelled signal can be obtained on the backward
extracted signal with an opposite frequency shifter. Thus, the phase noise of
the remote site and the access sites are automatically cancelled.

As demonstrated by Williams et al., the phase noise rejection capability is
limited by the propagation delay [23]. The residual phase noise PSD,
$S_{\text{access}}(f)$, at each access place in terms of the single-pass free-
running phase noise PSD, $S_{\text{fiber}}(f)$, and the single-pass fiber link
$L$ propagation delay, $\tau_{0}$, by compensating the forward optical signal
can be expressed as [23, 22],

$S_{\text{access}}(f)\simeq{(2\pi
f\tau_{0})^{2}}\left[1+2\frac{L_{a}}{L}-\frac{2}{3}\left(\frac{L_{a}}{L}\right)^{2}\right]S_{\text{fiber}}(f).$
(4)

Note that, in comparison with the conventional multiple-access technique, the
phase noise rejection capability of both the conventional and proposed
techniques is proportional to $\tau^{2}_{0}$ and the uncompensated single-pass
fiber phase noise PSD, $S_{\text{fiber}}(\omega)$. Although, the residual
fiber phase noise PSD in the proposed technique turns out to be a factor of 7
worse than that of in the conventional multiple-access scheme at the output of
the fiber link ($L_{a}=L$) [16, 19], it could largely suppress the phase noise
introduced by the servo bumps (see Fig. 2), which could significantly increase
as the increase of the fiber link [24], and simplify the hardware overhead in
systems with many stations and connections [21].

Figure 2: Measured phase noise PSDs of a the 70/30 km free-running access
place (blue curve), b the 50/50 km free-running access place (orange curve), c
the compensated 70/30 km access place (black curve), d the compensated 50/50
km access place (magenta curve), e the 70/30 km theoretical delay-limited
value (gray curve), and f the 50/50 km theoretical delay-limited value (cyan
curve). Strong servo bumps can be effectively removed in the proposed scheme.

We implement the proposed extraction scheme in Fig. 1 on a 100 km fiber link.
The signal source that has been adopted in this experiment is a narrow-
linewidth optical source (NKT X15) at a frequency of near 193 THz with a
typical linewidth of 100 Hz. At an arbitrary access place, a $2\times 2$
optical coupler (coupling ratio 50%: 50%) is used to tap both forward and
backward lights. To demonstrate the phase noise rejection capability at any
places along the fiber link, we choose two representative places, which are at
the most symmetric one, 50/50 km ($L_{a}/L_{b}$), and the relative asymmetric
one, 70/30 km. To compensate the total loss of 23 dB of the 100 km fiber link
system, we implement a home-made bidirectional Erbium-doped fiber amplifier
(EDFA) at the middle of the fiber link. The angular frequencies at the local
and remote sites are $\omega_{l}=2\pi\times 75$ MHz, $\omega_{r}=2\pi\times
45$ MHz and $\omega_{s}=2\pi\times 115$ MHz, respectively. The AOM1 is
simultaneously supplied by two angular frequencies of $2\pi\times 85$ MHz and
$2\pi\times 75$ MHz, and the AOM3 is supplied by the angular frequency of
$2\pi\times 45$ MHz. Consequently, the beatnotes between the input and
extracted sites have a frequency of $2\pi\times 40$ MHz. This beatnote is
recorded simultaneously with a dead-time free counter with a gate time of 1 s
and non-averaging $\Pi$-type operation. At the same time, we implement the
phase noise measurement by feeding the heterodyne beatnote signal together
with a stable frequency reference to a phase detector. The phase fluctuations
$S_{\phi}(\omega)$ at the phase detector output are measured with a fast
Fourier transform (FFT) analyzer [22, 25].

Figure 3: Measured fractional frequency instability, calculated from
$\Pi$-type data with the overlapping Allan deviation (ADEV), of the 100-km
free-running main link (red triangles), the stabilized 70/30 km access place
output (black squares), the compensated 50/50 km access place output (blue
cicles), and the noise floor, in which the fiber links are replaced by the
short fibers, of the access place output (gray diamonds)

Figure 2 shows the phase noise power spectral densities (PSDs) of the
stabilized 70/30 km access place (c, black curve) and the stabilized 50/50 km
access place (d, magenta curve), respectively. The phase noise PSDs of the
free-running 70/30 km access place (a, blue curve) and the 50/50 km access
place (b, orange curve) are also shown. The free-running curve is typical for
optical fibre links, with a noise of approximately $100$ rad${}^{2}/$Hz at 1
Hz for the 70/30 km access place, scaling down with a slope of about $f^{-2}$,
and reaching around $10^{-3}$ rad${}^{2}/$Hz after 100 Hz. Their noises of the
free-running 70/30 km (a, blue curve) and 50/50 km access places (b, orange
curve) slightly differ because the two measurements are performed at different
times. At the same time, both passive phase stabilized cases are also very
similar with the phase noise PSDs about $10^{-3}$ rad${}^{2}/$Hz between $1$
and $100$ Hz. We can clearly see that the noise correction is limited by the
main fiber propagation delay approximately $\sim 150$ Hz, which is compatible
with the theoretical bandwidth of 189 Hz given by
$1/(4\sqrt{\text{7}}\tau_{0})$. This limit is the same for both access places
because the bandwidth of the extraction is limited by the longer delay, namely
the main fiber link delay [16]. More importantly, we can clearly see that
compared with the conventional multiple-access optical frequency transfer
technique [13, 14, 15], the strong servo bumps are effectively suppressed in
our passive phase noise cancellation technique as observed in our previous
work [19, 21, 22].

Figure 3 displays the fractional frequency stability of the free-running and
passively stabilized configurations, calculated from $\Pi$-type data with
overlapping Allan deviation (ADEV). The stability of the 70/30 km (50/50 km)
access place is $3.1\times 10^{-15}$ ($2.0\times 10^{-15}$) at 1 s averaging
time, decreases as a slope of approximately $\tau^{-1}$ and reaches a floor of
approximately $4.1\times 10^{-19}$ ($3.3\times 10^{-19}$) at $10^{4}$ s. As
calculated by Eq. 4, the ratio of the stability of the $50/50$ km and $30/70$
km places should be $\mathcal{R}=0.64$. In our experiment, we also obtain the
ratio of $\mathcal{R}\simeq 0.64$ ($2.0\times 10^{-15}/3.1\times 10^{-15}$).
As a comparison, we also measured the noise floor for the access output as
displayed in Fig. 3. The long-term stability is mainly attributed to the
thermal noise on uncompensated fibre paths, due to imperfect length adjustment
and thermal stabilisation in the extraction optical set-up, or in the main
fiber link [22, 25]. Compared with our previous work in which we characterized
the system performance at the output of the main fiber link [21], the
ultrastable laser can be thus transferred through the secondary and main links
without significant degradation.

Figure 4: Frequency comparison between input and access place frequencies
after 70 km over the 100 km fiber. (b) 84,014 data points were taken with
dead-time free $\Pi$-type frequency counters with a 1 s gate time (green
points, left axis). We calculated unweighted mean ($\Pi$-type) values for all
cycle-slip free 1,000 s long segments, resulting in 84 data points (black
dots, right frequency axis, enlarged scale). Histograms (brown bars) and
Gaussian fits (red curves) for (a) frequency values taken with one second gate
time and (c) 84 phase coherent 1,000-second frequency averages.

Complementary to the characterization of the stability and phase noise, the
accuracy has to be throughout-fully examined by calculating the mean value of
the end-to-end beatnote frequency offset. Figure 4(b) shows the frequency
deviation of the beatnote’s data for the 70/30 km access place, recorded with
a 1 s gate time and $\Pi$-type counters, over successive 84,014 s (green
point, left axis) and the arithmetic mean of all cycle-slip free 1,000 s
intervals (black dots, right axis). Histograms (brown bars) and Gaussian fits
(red curves) of the frequency deviation for the access place after 70 km are
also illustrated in Fig. 4(a) and (c). According to the Gaussian fit in Fig.
4(c), the calculated results demonstrate that the mean frequency is shifted by
-28.9 $\mu$Hz ($-1.5\times 10^{-19}$). The standard deviation of the 1,000 s
data points is 493 $\mu$Hz ($2.5\times 10^{-18}$), which is a factor of 1,000
smaller than the ADEV at 1 s as expected for this $\Pi$-type evaluation.
Considering the long-term stability of frequency transfer as illustrated in
Fig. 3, we conservatively estimate the accuracy of the transmitted optical
signal as shown in the last data point of the ADEV, resulting in a relative
frequency accuracy of $4.1\times 10^{-19}$.

Adopting the same procedure, the mean frequency offset and the standard
deviation for the 50/50 km place were calculated using the 1000 s point with
the total 89,690 $\Pi$-type counter data to be 52.2 $\mu$Hz ($2.7\times
10^{-19}$) and 454 $\mu$Hz ($2.3\times 10^{-18}$), respectively. Taking into
account the long-term ADEV at 10,000 s of the data set for the 50/50 km access
place of $3.3\times 10^{-19}$, we conservatively estimate that the mean
frequency offset is $2.7\times 10^{-19}$ with a statistical uncertainty of
$3.3\times 10^{-19}$ for the 50/50 km access place. We can conclude that there
is no systematic frequency shift arising in the all-passive multiple-place
phase noise cancellation setup at a level of a few $10^{-19}$.

In summary, we have presented a new method, making a stable optical frequency
available at any arbitrary access places along the fiber link with passive
phase noise cancellation. In comparison with previous work, here we
demonstrated that a high performance multiple-place optical frequency transfer
over the widely adopted bus fiber topology with the open-loop design,
mitigating some technical difficulties in conventional active multiple-access
phase noise cancellation. We experimentally demonstrate transferring of
optical frequency to two different access places. After being compensated,
delivering an optical frequency to different places with the relative
frequency instability in terms of ADEV measured by the $\Pi$-mode frequency
counter can be as low as $3.1\times 10^{-15}$ at 1 s and $4.1\times 10^{-19}$
at 10,000 s. The frequency uncertainty of the light after transferring through
the fiber relative to that of the input light is a few $10^{-19}$ for the
access places over the 100 km fiber link.

Funding. This research was supported by the National Natural Science
Foundation of China (NSFC) (61905143, 61627817).

Disclosures. The authors declare no conflicts of interest.

## References

  * [1] B. J. Bloom, T. L. Nicholson _et al._ , Nature 506, 71 (2014).
  * [2] E. Oelker, R. B. Hutson _et al._ , Nature Photonics 13, 714 (2019).
  * [3] W. F. McGrew, X. Zhang, R. J. Fasano _et al._ , Nature 564, 87 (2018).
  * [4] N. Huntemann, C. Sanner _et al._ , Physical Review Letters 116, 063001 (2016).
  * [5] S. M. Brewer, J.-S. Chen _et al._ , Physical Review Letters 123, 033201 (2019).
  * [6] R. H. Parker, C. Yu, W. Zhong _et al._ , Science 360, 191 (2018).
  * [7] A. Derevianko and M. Pospelov, Nat. Phys. 10, 933 (2014).
  * [8] K. Van Tilburg, N. Leefer _et al._ , Phys. Rev. Lett. 115, 011802 (2015).
  * [9] A. Arvanitaki, J. Huang _et al._ , Phys. Rev. D 91, 015015 (2015).
  * [10] L. Hu, N. Poli _et al._ , Phys. Rev. Lett. 119, 263601 (2017).
  * [11] M. S. Safronova, D. Budker _et al._ , Rev. Mod. Phys. 90, 025008 (2018).
  * [12] L.-S. Ma, P. Jungner, J. Ye, and J. L. Hall, Opt. Lett. 19, 1777 (1994).
  * [13] G. Grosche, Opt. Lett. 39, 2545 (2014).
  * [14] S. W. Schediwy _et al._ , Optics Letters 38, 2893 (2013).
  * [15] L. Wu, Y. Jiang, C. Ma _et al._ , Opt. Lett. 41, 4368 (2016).
  * [16] A. Bercy, S. Guellati-Khelifa _et al._ , J. Opt. Soc. Am. B 31, 678 (2014).
  * [17] P. Krehlik, L. Sliwczynski, L. Buczek, and M. Lipinski, IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control 60, 1804 (2013).
  * [18] G. Grosche, German patent application DE 10 (2010).
  * [19] L. Hu, X. Tian, G. Wu _et al._ , Journal of Lightwave Technology 38, 3644 (2020).
  * [20] K. H. Ang _et al._ , IEEE Transactions on Control Systems Technology 13, 559 (2005).
  * [21] L. Hu, X. Tian _et al._ , Opt. Lett. 45, 4308 (2020).
  * [22] L. Hu, X. Tian _et al._ , Journal of Lightwave Technology 38, 5916 (2020).
  * [23] P. A. Williams _et al._ , J. Opt. Soc. Am. B 25, 1284 (2008).
  * [24] S. M. Raupach, A. Koczwara _et al._ , Phys. Rev. A 92, 021801 (2015).
  * [25] X. Tian, L. Hu, G. Wu, and J. Chen, Journal of Lightwave Technology 38, 4270 (2020).