††thanks: These authors contributed equally.††thanks: These authors
contributed equally.

# Efficient simulation of open quantum systems coupled to a reservoir through
multiple channels

Hanggai Nuomin  Department of Chemistry, Duke University, Durham, North
Carolina 27708, United States  Jiaxi Wu  Kuang Yaming Honors School, Nanjing
University, Nanjing, Jiangsu 210023, China  Department of Physics, Duke
University, Durham, North Carolina 27708, United States  Peng Zhang
<EMAIL_ADDRESS>Department of Chemistry, Duke University, Durham, North
Carolina 27708, United States  David N. Beratan<EMAIL_ADDRESS>Department of Chemistry, Duke University, Durham, North Carolina 27708, United
States  Department of Physics, Duke University, Durham, North Carolina 27708,
United States  Department of Biochemistry, Duke University, Durham, North
Carolina 27710, United States

###### Abstract

The simulation of open quantum systems coupled to a reservoir through multiple
channels remains a substantial challenge. This kind of open quantum system
arises when considering the radiationless decay of excited states that are
coupled to molecular vibrations, for example. We use the chain mapping
strategy in the interaction picture to study systems linearly coupled to a
harmonic bath through multiple interaction channels. In the interaction
picture, the bare bath Hamiltonian is removed by a unitary transformation (the
system-bath interactions remain), and a chain mapping transforms the bath
modes to a new basis. The transformed Hamiltonian contains time-dependent
local system-bath couplings. The open quantum system is coupled to a limited
number of (transformed) bath modes in the new basis. As such, the entanglement
generated by the system-bath interactions is local, making efficient dynamical
simulations possible with matrix product states. We use this approach to
simulate singlet fission, using a generalized spin-boson Hamiltonian. The
electronic states are coupled to a vibrational bath both diagonally and off-
diagonally. This approach generalizes the chain mapping scheme to the case of
multi-channel system-bath couplings, enabling the efficient simulation of this
class of open quantum systems using matrix product states.

##  I Introduction

Dynamical simulations of open quantum systems have been significantly advanced
by the development of chain-mapped time-dependent density-matrix
renormalization group (DMRG) methods. These methods enhance the efficiency of
numerical simulations [1, 2, 3, 4, 5, 6, 7, 8]. A bath continuum is mapped to
a semi-infinite chain of modes (e.g., bosons or fermions) with nearest-
neighbor couplings using an orthogonal transformation in the chain mapping
method. As a result, the open quantum system in the transformed basis is
coupled only to the first mode of the chain. It was found that the chain
Hamiltonian has significantly slower entanglement growth during its time
evolution than the untransformed Hamiltonian [9]. In the transformed basis,
the quantum state of the system and bath lives in a small Hilbert space
spanned by a few dozen important states (adaptively generated by the DMRG
algorithm). As a consequence, the chain representation is particularly well
suited for efficient DMRG simulations. However, chain mapping is only possible
when a single spectral density describes the open quantum system. That is, for
a system interacting with a bath through many kinds of system-bath couplings
(i.e., multi-channel couplings [10]), the mapping of the bath to a nearest-
neighbor chain is not feasible. Efficiently simulating open quantum systems of
this class using DMRG, or equivalently, using matrix product states (MPS),
remains an open challenge [11]; however, open quantum systems with multi-
channel couplings are ubiquitous in chemical dynamics. Bridge-mediated
electron transfer in donor-bridge-acceptor complexes [12, 13, 14, 15] provides
a familiar example of such open quantum systems, as the donor and acceptor
electronic energies and couplings are influenced by nuclear motions. For this
reason, developing efficient algorithms to simulate the time evolution of
multi-channel open quantum systems is particularly timely.

In earlier studies [9], the mapped chain Hamiltonian was transformed into the
interaction picture with respect to the reservoir (bath) Hamiltonian, referred
to as the interaction picture chain (IC) Hamiltonian. In the IC Hamiltonian,
the system-bath coupling coefficients are time-dependent, in contrast to the
constant coefficients before transformation (i.e., the Schrödinger picture).
The IC time-dependent coupling coefficients are found to be localized, and the
system-bath interactions reside in a narrow window of bath modes and have a
travelling wave pattern [16] that resembles a chain structure. That is, the
system interacts with a small number of bath modes at all times, enabling
efficient simulations.

Here, we apply the interaction picture chain mapping to multi-channel open
quantum systems and study the pattern of time-dependent system-bath couplings.
As an example, we simulate the population dynamics and entanglement growth of
a singlet fission process [17] mediated by two types (diagonal and off-
diagonal) of vibronic couplings (two distinct spectral densities and one
bath), as an example. The singlet fission system was studied earlier by Zhao
et al. using the multiple Davydov ansatz [18, 19]. Our simulations show
slightly different dynamics from theirs, while both methods are exact in
principle, given sufficiently tight convergence thresholds. The disagreement
in simulation results requires a careful reexamination of the system studied
in Ref [18].

Our method extends the chain-mapped time dependent DMRG in the interaction
picture to the multi-channel open quantum systems. The efficiency of the
method is ensured by the localized time dependent system-bath couplings. Two
variants of chain mappings (Lanczos and block Lanczos) are presented in Sec.
II. In Sec. III, we apply both chain mappings to the singlet fission system,
and calculate its time evolution.

##  II THEORY AND METHODS

A general spin-Boson Hamiltonian is

$\displaystyle\begin{split}\hat{H}=&\Delta_{x}\hat{\sigma}_{x}+\hat{\sigma}_{x}\otimes\int{d\omega\,h_{x}(\omega)(\hat{a}_{\omega}^{\dagger}+\hat{a}_{\omega})}+\Delta_{z}\hat{\sigma}_{z}\\\
&+\hat{\sigma}_{z}\otimes\int{d\omega\,h_{z}(\omega)(\hat{a}_{\omega}^{\dagger}+\hat{a}_{\omega})}+\int{d\omega\,\omega\hat{a}_{\omega}^{\dagger}\hat{a}_{\omega}}.\end{split}$
(1)

$\hat{\sigma}_{x}$ and $\hat{\sigma}_{z}$ are the Pauli matrices, $\Delta_{x}$
and $\Delta_{z}$ are prefactors, $\hat{a}_{\omega}$
($\hat{a}^{\dagger}_{\omega}$) is the annihilation (creation) operator for
bosons with frequency $\omega$; $h_{z}(\omega)$ and $h_{x}(\omega)$ are the
frequency-dependent couplings between the spin and the bath (through
$\hat{\sigma}_{x}$ and $\hat{\sigma}_{z}$). We term them as $x$-interactions
and $z$-interactions, respectively. The Hamiltonian is discretized as

$\displaystyle\begin{split}\hat{H}=&\Delta_{x}\hat{\sigma}_{x}+\underbrace{\hat{\sigma}_{x}\otimes\sum_{i=0}^{N}{\xi_{i}(\hat{a}_{i}^{\dagger}+\hat{a}_{i})}}_{x\text{-interactions}}\\\
&+\Delta_{z}\hat{\sigma}_{z}+\underbrace{\hat{\sigma}_{z}\otimes\sum_{i=0}^{N}{\zeta_{i}(\hat{a}_{i}^{\dagger}+\hat{a}_{i})}}_{{}_{z\text{-interactions}}}+\underbrace{\sum_{i=0}^{N}\omega_{i}\hat{a}_{i}^{\dagger}\hat{a}_{i}}_{\mathrm{bath}}.\end{split}$
(2)

In the discretization, we use the Legendre polynomial discretization [20] and
choose the same discretization nodes (i.e., discrete frequencies) for the two
spectral densities $J_{x}=\pi h^{2}_{x}(\omega)$ and $J_{z}(\omega)=\pi
h^{2}_{z}(\omega)$. Any one of the two sets of discrete system-bath coupling
coefficients $\\{\xi_{i}\\}$ and $\\{\zeta_{i}\\}$ can be used to generate a
chain mapping, which allows the system-bath interaction terms characterized by
this set of couplings to be transformed to a nearest-neighbor chain form [20,
7, 21]. We refer to the chain mapping generated by $\\{\xi_{i}\\}$ as mapping
$\mathbf{X}$, and the mapping generated by $\\{\zeta_{i}\\}$ as mapping
$\mathbf{Z}$. It is generally impossible to transform both the $x$ and
$z$-type system-bath interaction terms to the chain form simultaneously,
unless the vectors $(\xi_{1},\ldots,\xi_{N})$ and
$(\zeta_{0},\ldots,\zeta_{N})$ are parallel, so that the $\mathbf{X}$ and
$\mathbf{Z}$ mappings are identical [22]. The optimal nearest-neighbor chain
form of the Hamiltonian in Eq. (2) therefore does not exist. Although one can
simulate the open quantum system without transformation (i.e., the star form
as in Eq. (2), see details in 111Since the spin interacts with each boson, the
couplings connect the spin to the bosons, forming a star pattern.), the
numerical cost is large [9]. A strategy to circumvent this difficulty, and to
enable efficient MPS simulations, is to transform the Hamiltonian of Eq. (2)
into the interaction picture with respect to the bath,

$\displaystyle\begin{split}\hat{H}_{I}=&\Delta_{x}\hat{\sigma}_{x}+\hat{\sigma}_{x}\otimes\sum_{j=0}^{N}{\xi_{j}(\hat{a}_{j}^{\dagger}e^{i\omega_{j}t}+\hat{a}_{j}e^{-i\omega_{j}t})}\\\
&+\Delta_{z}\hat{\sigma}_{z}+\hat{\sigma}_{z}\otimes\sum_{j=0}^{N}{\zeta_{j}(\hat{a}_{j}^{\dagger}e^{i\omega_{j}t}+\hat{a}_{j}e^{-i\omega_{j}t})},\end{split}$
(3)

and to use either of the two chain mappings ($\mathbf{X}$ or $\mathbf{Z}$) to
transform the interaction-picture Hamiltonian.

Mapping $\mathbf{X}$, when applied to the $x$-coupling terms, can produce a
set of time-dependent couplings that resemble the desired chain form, as shown
in Ref [9]. For the $z$-couplings, the chain mapping $\mathbf{X}$ cannot
transform the $z$-interaction terms into a nearest-neighbor interacting chain
form. However, in the interaction picture, mapping $\mathbf{X}$ produces time
dependent $z$-type couplings with the localized traveling-wave pattern (see
Fig. 2), which is essential for efficiently simulating the time-dependent
Hamiltonian Eq. (2) (see Ref [9]).

Chain mapping $\mathbf{X}$ can be generated as follows. We transform
$\mathbf{w}=\operatorname{diag}(\omega_{0},\cdots,\omega_{N})$ into a
tridiagonal chain form by applying an orthogonal transformation to
$\mathbf{w}$: $\mathbf{P}^{T}\mathbf{w}\mathbf{P}=\mathbf{W}$. The orthogonal
matrix $\mathbf{P}$ is then used to generate a new set of bath modes:
$\hat{b}^{\dagger}_{j}=\sum_{i}P_{ij}\hat{a}^{\dagger}_{i}$ (and
$\hat{b}_{j}=\sum_{i}P_{ij}\hat{a}_{i}$). The bath Hamiltonian can be written
in terms of the new operators $\hat{b}^{\dagger}_{i}$ and $\hat{b}_{i}$ as

$\displaystyle
H_{b}=\sum_{i}{\omega_{i}\hat{a}_{i}^{\dagger}\hat{a}_{i}}=\mathbf{b^{\dagger}}\mathbf{W}\mathbf{b}$
(4)

where
$\mathbf{\hat{b}^{\dagger}}=(\hat{b}_{0}^{\dagger},\ldots,\hat{b}_{N}^{\dagger})$
and $\mathbf{b}=(\hat{b}_{0},\ldots,\hat{b}_{N})^{T}$ represent the new bath
modes, and the tridiagonal matrix $\mathbf{W}$ is

$\displaystyle\mathbf{W}=\begin{pmatrix}\alpha_{0}&\beta_{1}&&&\\\
\beta_{1}&\alpha_{1}&\beta_{2}&&\\\ &\ddots&\ddots&\ddots&\\\ \
&&\beta_{N-1}&\alpha_{N-1}&\beta_{N}\\\ &&&\beta_{N}&\alpha_{N}\end{pmatrix}.$
(5)

By further requiring the $x$-interaction terms to become the interactions
between only the system and the zeroth mode:
$\hat{\sigma}_{x}\otimes\sum_{i}\xi_{i}(\hat{a}^{\dagger}_{i}+\hat{a}_{i})\to\hat{\sigma}_{x}\otimes\tilde{\xi}_{0}(\hat{b}^{\dagger}_{0}+\hat{b}_{0})$,
with $\tilde{\xi}_{0}=\sqrt{\sum_{i}\xi_{i}^{2}}$, we obtain the first column
of $\mathbf{P}$: $P_{0i}=\xi_{i}/\tilde{\xi}_{0}$ [24, 25, 7, 8]. The
orthogonal matrix $\mathbf{P}$ is then determined by solving

$\displaystyle\mathbf{wP}=\mathbf{PW}.$ (6)

for the other columns of $\mathbf{P}$.

The orthogonal matrix $\mathbf{P}$ produced by the above procedure links the
old operators $(\hat{a}^{\dagger}_{i},\hat{a}_{i})$ to the new ones
$(\hat{b}^{\dagger}_{j},\hat{b}_{j})$. One can replace the old operators in
Eq. (3) with the new operators and obtain a set of time-dependent couplings:

$\displaystyle\begin{split}\hat{H}_{I}=&\Delta_{x}\hat{\sigma}_{x}+\hat{\sigma}_{x}\otimes\sum_{jk=0}^{N}P_{jk}{\xi_{j}(\hat{b}_{k}^{\dagger}e^{i\omega_{j}t}+\hat{b}_{k}e^{-i\omega_{j}t})}\\\
&+\Delta_{z}\hat{\sigma}_{z}+\hat{\sigma}_{z}\otimes\sum_{jk=0}^{N}P_{jk}{\zeta_{j}(\hat{b}_{k}^{\dagger}e^{i\omega_{j}t}+\hat{b}_{k}e^{j\omega_{i}t})}\\\
=&\Delta_{x}\hat{\sigma}_{x}+\hat{\sigma}_{x}\otimes\sum_{k=0}^{N}[\xi_{k}(t)\hat{b}_{k}^{\dagger}+\xi^{*}_{k}(t)\hat{b}_{k}]\\\
&+\Delta_{z}\hat{\sigma}_{z}+\hat{\sigma}_{z}\otimes\sum_{k=0}^{N}[\zeta_{k}(t)\hat{b}_{k}^{\dagger}+\zeta^{*}_{k}(t)\hat{b}_{k}]\end{split}$
(7)

where $\xi_{k}(t)=\sum_{j}P_{jk}\xi_{j}e^{-i\omega_{j}t}$ and
$\zeta_{k}(t)=\sum_{j}P_{jk}\zeta_{j}e^{-i\omega_{j}t}$. $\mathbf{P}$ is not
the only possible orthogonal transformation that can be used to transform the
bath Hamiltonian of Eq. (4); rather, any orthogonal matrix can be used to
transform the diagonal harmonic bath Hamiltonian and may be used in Eq. (2) to
obtain a new interaction-picture Hamiltonian. However, different orthogonal
transformations produce different numerical efficiencies. As shown in [9], the
numerical efficiency of simulating the time-dependent Hamiltonians of Eq. (7)
depends on the properties of the orthogonal matrix ($\mathbf{P}$ or other
possible orthogonal matrices).

The Hamiltonian [Eq. (7)] has a natural ordering of sites for the MPS
representation of wave functions. The first site of the MPS is the spin; the
others represent bath oscillators, and the oscillators with smaller indices
are closer to the spin site, as shown in Fig. 1.

Figure 1: The MPS of the spin and bath in Eq. (7):
$\ket{\psi_{sb}}=\sum_{\sigma_{s}n_{0}\cdots
n_{N}}\mathbf{A}^{\sigma_{s}}\mathbf{A}^{n_{0}}\cdots\mathbf{A}^{n_{N}}$.
$\sigma_{s}$ is the spin basis index. $n_{0},\ldots,n_{N}$ are indices of bath
mode eigenbasis. The left-most matrix $\textbf{A}^{\sigma_{s}}$ corresponds to
the spin, while $\textbf{A}^{n_{0}},\ldots,\textbf{A}^{n_{N}}$ correspond to
bath modes 0 to $N$, respectively.

Any orthogonal matrix can replace the orthogonal matrix $\mathbf{P}$ without
changing the dynamics of the spin, as noted above. However, the numerical
efficiency of a MPS simulation depends crucially on the choice of the
orthogonal matrix. Here, we provide an alternative to $\mathbf{P}$. We
construct an orthogonal matrix with its first two columns orthogonal to the
vector space spanned by $\\{\zeta_{i}\\}$ and $\\{\xi_{i}\\}$ in Eq. (2) [26,
24]. One can require the orthogonal matrix (denoted $\mathbf{Q}$) to transform
the diagonal bath Hamiltonian into a chain of oscillators with nearest and
next-nearest-neighbor couplings so that the desired traveling wave structure
still appears, thus reducing the entanglement growth over the course of time.
The orthogonal matrix $\mathbf{Q}$ is then uniquely determined by solving the
following linear equation recursively [27, 28]:

$\displaystyle\begin{split}&\mathbf{Q}^{T}\mathbf{w}\mathbf{Q}=\\\
&\begin{pmatrix}\alpha_{0}&\beta_{1}&\kappa_{2}&0&&&\\\
\beta_{1}&\alpha_{1}&\beta_{2}&\kappa_{3}&0&&\\\
\kappa_{2}&\beta_{2}&\alpha_{2}&\beta_{3}&\kappa_{4}&0&&\\\
0&\kappa_{3}&\beta_{3}&\alpha_{3}&\beta_{4}&\kappa_{5}&0\\\
&\ddots&\ddots&\ddots&\ddots&\ddots&\ddots&\ddots\\\
&&0&\kappa_{N-3}&\beta_{N-3}&\alpha_{N-3}&\beta_{N-2}&\kappa_{N-1}&0\\\
&&&0&\kappa_{N-2}&\beta_{N-2}&\alpha_{N-2}&\beta_{N-1}&\kappa_{N}\\\
&&&&0&\kappa_{N-1}&\beta_{N-1}&\alpha_{N-1}&\beta_{N}\\\
&&&&&0&\kappa_{N}&\beta_{N}&\alpha_{N}\end{pmatrix}.\end{split}$ (8)

The orthogonal matrix $\mathbf{Q}$ can replace $\mathbf{P}$ (the Lanczos
matrix) in Eq. (7). The matrix $\mathbf{Q}$ is known as the block Lanczos
matrix [28].

To show that a single transformation ($\mathbf{P}$ or $\mathbf{Q}$) can
generate two sets of couplings with a traveling-wave pattern [9], we use two
distinct spectral densities for $x$ and $z$-couplings to obtain the
discretized Hamiltonian [Eq. (2)]. The discretized Hamiltonian is then
transformed by the orthogonal transformation generated by the Lanczos method
regarding either one of the two couplings [Fig. 2, panel (b)], or, instead, by
the transformation $\mathbf{Q}$ generated by the block Lanczos method using
two spectral densities simultaneously [Fig. 2, panel (c)]. The parameters of
the spectral densities are noted in the captions of Fig. 2. The traveling-wave
pattern appears in the time-dependent couplings ($\xi_{i}(t)$ and
$\zeta_{i}(t)$) of both the $x$ and $z$-interactions. This pattern ensures the
efficiency of the interaction-picture MPS simulations for the chain-mapped
Hamiltonian with two distinct spectral densities, since the interactions
between the system and the bath are localized spatially to a relatively narrow
window of bath oscillators [9].

Figure 2: (a): The spectral densities $J_{z}$ and $J_{x}$ characterizing the
$z$ and $x$ system-bath couplings.
$J_{z}(\omega)=L(\omega,\Omega_{1},\eta,\lambda)+L(\omega,\Omega_{2},\eta,\lambda)+L(\omega,\Omega_{3},\eta,\lambda)$
with $\Omega_{1}=2$, $\Omega_{2}=5$, $\Omega_{3}=10$, $\eta=1.5$, and
$\lambda=1$. The Lorentzian function is
$L(\omega,\Omega,\eta,\lambda)=\frac{2\lambda\Omega^{2}\eta\omega}{(\omega^{2}-\Omega^{2})^{2}+\eta^{2}\omega^{2}}$.
$J_{x}(\omega)=\lambda\omega e^{-\omega/\omega_{c}}$ with $\lambda=2$,
$\omega_{c}=5$. (b): The absolute values of the Lanczos-transformed
$z$-couplings $|\zeta_{j}(t)|$ and $x$-couplings $|\xi_{j}(t)|$ for bath
oscillators (indexed by $n$, see Eq. (7) with the definitions of $\xi_{j}(t)$
and $\zeta_{j}(t)$). The Lanczos transformation is generated by the
$z$-couplings $\\{\zeta_{i}\\}$ [Eq. (2)]. (c): The absolute values of the
block Lanczos transformed $z$ couplings $|\zeta_{j}(t)|$ and $x$-couplings
$|\xi_{j}(t)|$ for bath oscillators (see Eq. (7) and Eq. (8)).

##  III Application

The spin-boson model is widely used to describe transitions between molecular
electronic states coupled to nuclear vibrations, environmental phonons, and
photon modes of the radiation field [29]. This approach is frequently used to
simulate non-adiabatic dynamics, including bridge-mediated electron transfer
[30, 31, 32].

To explore the validity of the method developed in Sec.II, we apply the method
to a spin-boson system with both diagonal and non-diagonal couplings. In
organic crystals where intramolecular and intermolecular vibrations induce
diagonal and off-diagonal exciton-phonon couplings [33, 34], two spectral
densities are required to describe the coupling. We model singlet fission,
where a singlet exciton state $\ket{S_{1}}$ evolves into two triplet excitons
$\ket{TT}$ [35]. In some organic materials, singlet fission is facilitated by
an intermediate charge transfer state that couples the singlet and triplets
states [17]. The rate of singlet fission is influenced by both diagonal and
off-diagonal vibronic interactions. The Hamiltonian is represented by Eq. (2),
including the Hamiltonian of the exciton states, harmonic nuclear vibrations,
and the interactions between them. Here, we use the model singlet fission
Hamiltonian of Huang et al. [17]:

$\displaystyle\begin{split}\hat{H}=&\Delta_{x}(\ket{S_{1}}\bra{TT}+h.c.)+\Delta_{z}\ket{S_{1}}\bra{S_{1}}+\sum_{i}\hbar\omega_{i}\hat{a}^{\dagger}_{i}\hat{a}_{i}+\hat{H}_{sb},\\\
\hat{H}_{sb}=&\big{(}\sqrt{\lambda_{\smash{S_{1}}}}\ket{S_{1}}\bra{S_{1}}+\sqrt{\lambda_{\smash{TT}}}\ket{TT}\bra{TT}\big{)}\otimes\sum_{i}\zeta_{i}(\hat{a}^{\dagger}_{i}+\hat{a}_{i})\\\
&+\sqrt{\lambda_{\smash{o.d.}}}\big{(}\ket{S_{1}}\bra{TT}+h.c.\big{)}\otimes\sum_{i}\xi_{i}(\hat{a}^{\dagger}_{i}+\hat{a}_{i}).\end{split}$

$\lambda$ values characterize the magnitude of system-bath couplings. Two
Lorentzian spectral densities $J_{z},\,J_{x}$ define the diagonal and off-
diagonal couplings, $\\{\xi_{i}\\}$ and $\\{\zeta_{i}\\}$, respectively:

$\displaystyle J_{z}(\omega)$
$\displaystyle=\pi\sum_{i}\zeta_{i}^{2}\delta(\omega-\omega_{i})=\frac{4\gamma_{diag}\omega_{diag}^{2}\omega}{(\omega^{2}-\omega_{diag}^{2})^{2}+4\gamma_{diag}^{2}\omega^{2}},$
(9) $\displaystyle J_{x}(\omega)$
$\displaystyle=\pi\sum_{i}\xi_{i}^{2}\delta(\omega-\omega_{i})=\frac{4\gamma_{o.d.}\omega_{o.d.}^{2}\omega}{(\omega^{2}-\omega_{o.d.}^{2})^{2}+4\gamma_{o.d.}^{2}\omega^{2}}.$
(10)

$\omega_{diag}$ and $\omega_{o.d.}$ are vibrational frequency centers of the
Lorentzian spectral densities, and
$\gamma_{diag}=\gamma_{o.d.}=1\,\mathrm{ps}^{-1}$ are vibrational relaxation
rates. Values of the other parameters appear in Table 1.

Quantity | Value
---|---
Energy difference $\Delta_{z}$ | $100\,\mathrm{meV}$
Electronic coupling $\Delta_{x}$ | $20\,\mathrm{meV}$
Coupling strength $\lambda_{S_{1}}$ | $0.7\hbar\omega_{diag}$
Coupling strength $\lambda_{S_{TT}}$ | $1.4\hbar\omega_{diag}$
Coupling strength $\lambda_{o.d.}$ | $0.1\hbar\omega_{o.d.}$
Frequency $\omega_{diag}$, $\omega_{o.d.}$ | See Fig. 3
Spectral density cutoff | [0, $800\,\mathrm{cm}^{-1}$]
Bath temperature | $0\,\mathrm{K}$
Number of bath oscillators | $300$
SVD cutoff for MPS | $10^{-4}$
Energy levels of bath modes | 160
Table 1: Parameters used in the singlet fission simulations. Figure 3: Time
evolution of the singlet state population. The $x$-axis shows the time in
picosecond, and the $y$-axis shows the population of the singlet state. The
left column shows simulation results for the Lanczos transformation, while the
right column presents population evolution of block Lanczos transformation.
Each panel shows the simulation results for diagonal vibrational frequencies
of $20,\,40,\,60,\,80\,\mathrm{meV}$. For the block Lanczos method,
$\omega_{diag}$ and $\omega_{o.d.}$ must have distinct values, otherwise the
matrix $\mathbf{Q}$ in Eq. (8) would not become orthogonal.

We ran the singlet fission simulation with an initial state
$\ket{\varphi(0)}=\ket{S_{1}}\ket{0}$, where $\ket{0}$ indicates that all bath
harmonic oscillators are in their ground state at $T=0\,\mathrm{K}$ and $t=0$.
Varying the parameters $\omega_{diag},\,\omega_{o.d.}$ in $J_{z}$ and $J_{x}$,
we compute the population evolution of the singlet state for different central
mode frequencies ($\omega_{diag}$ and $\omega_{o.d.}$) of intermolecular and
intramolecular vibrations, as shown in Fig. 3. Our Lanczos and block Lanczos
simulation results are consistent in the time evolution of singlet state
population, and the time evolution behavior is similar to previous findings
for the singlet fission system [17]. However, we find different singlet state
populations after $1\,\mathrm{ps}$, which continues to drop in our simulation
but oscillates around $0.4$ in the earlier simulations. The discrepancy may be
due to the inadequate multiplicity in their multi-$D_{2}$ ansatz [17].

In the upper panels of Fig.3, where the off-diagonal coupling is weak compared
to the diagonal coupling, the singlet fission dynamics is dominated by
diagonal vibration modes and produces a resonant transition when
$\hbar\omega_{diag}=80\,\mathrm{meV}$. However, when $\hbar\omega_{o.d.}$ is
increased to $60\,\mathrm{meV}$, the off-diagonal coupling begins to break the
strict condition of resonant transitions in the singlet fission dynamics and
slows the transition rate accordingly.

Figure 4: Entanglement growth of MPS bonds (i.e, bipartite divisions) with
time. The evolution is computed using the Lanczos and block Lanczos
transformed Hamiltonians [Eq. (7)]. The y-axes of the panels are entanglement.
The time evolution direction is labeled by the arrow in the first panel. A
bond index of $0$ represents the bond between the spin and first bath mode,
while other indices represent bonds between bath modes. The block Lanczos
method shows higher values of entanglement (see Sec.III for the definition),
especially for the
$\hbar\omega_{diag}=80\,\mathrm{meV},\,\hbar\omega_{o.d.}=60\,\mathrm{meV}$
panels. These results demonstrate that MPS computations using the Lanczos
transformation have a faster simulation speed than the block Lanczos
transformation.

In a MPS, a quantum system is divided into two parts, with entanglement
between them. The division is referred to as a bond, which introduces the
concept of entanglement between the two parts. The entanglement for the
division (bond) measures how many states in one of the two parts are needed to
represent the wave function faithfully. As such, the larger the entanglement,
the slower the simulation becomes. Fig. 4 shows the entanglement growth of MPS
bonds with time evolution, providing comparisons between Lanczos and block
Lanczos algorithms. The growth of the MPS entanglement behaves slightly
differently for the two methods. However, the overall time-evolution patterns
of entanglement are similar for the Lanczos and block Lanczos methods, showing
a propagation of entanglement away from the spin. We conclude from Fig. 4 that
the Lanczos method provides a faster simulation speed than is found using the
block Lanczos algorithm. This is understood because the entanglement is
generated not only by nearest-neighbor interactions but also by next-nearest-
neighbor interactions.

##  IV Conclusions

We developed a DMRG-based approach to simulate multi-channel open quantum
system with two distinct spectral densities. We discretized the spectral
densities of the system-bath couplings using the Legendre polynomial
discretization method, producing two sets of coupling strengths
$\\{\xi_{i}\\}$ and $\\{\zeta_{i}\\}$. Then, we transformed the Hamiltonian to
the interaction picture, and applied an orthogonal matrix to perform chain
mapping and to generate new sets of modes. We developed two kinds of
orthogonal matrices, with one transforming the diagonal vibrational
frequencies into a tridiagonal matrix (generated by the Lanczos algorithm),
and the other resulting in a block tridiagonal matrix (generated by the block
Lanczos algorithm). The Lanczos algorithm, producing the chain mapping from
either $\\{\xi_{i}\\}$ or $\\{\zeta_{i}\\}$, generates only nearest-neighbor
interactions between chain modes. The block Lanczos method, in contrast,
generates additional next-nearest-neighbor interactions, and produces the
orthogonal matrix taking both $\\{\xi_{i}\\}$ and $\\{\zeta_{i}\\}$ into
account. In the interaction-picture chain mapping scenario, the interactions
between the system and bath are limited to a localized space within the MPS
chain, substantially accelerating the simulation.

We applied the two chain mapping techniques to simulate singlet fission, where
the singlet and triplet states are coupled non-adiabatically to a bath through
two Lorenztian spectral densities. By varying the central vibronic frequencies
of the spectral densities, we found different singlet fission rates that
characterize the decay rate of the singlet state population. The singlet
fission findings are in qualitative agreement with previous results but differ
at longer times. We also analyzed the efficiencies of the Lanczos and block
Lanczos methods by examining the entanglement evolution of the MPS bonds. We
found that the block Lanczos method generally has a higher entanglement value
than the Lanczos approach, leading to the higher efficiency for the Lanczos
method. The lower efficiency of the block Lanczos approach is attributed to
the introduction of next-nearest-neighbor interactions, since next-nearest-
neighbor interactions entangle a bath mode with many other bath modes.

###### Acknowledgements.

This material is based upon work supported by the National Science Foundation
under Grant No. 1955138.

## References

  * White [1993] S. R. White, Phys. Rev. B 48, 10345 (1993).
  * White [1992] S. R. White, Phys. Rev. Lett. 69, 2863 (1992).
  * White and Feiguin [2004] S. R. White and A. E. Feiguin, Phys. Rev. Lett. 93, 076401 (2004).
  * Xu  _et al._ [2022] Y. Xu, Z. Xie, X. Xie, U. Schollwöck, and H. Ma, JACS Au 2, 335 (2022).
  * Liu  _et al._ [2022a] K. T. Liu, F.-F. Song, D. N. Beratan, and P. Zhang, Phys. Rev. B 106, 104306 (2022a).
  * Ren  _et al._ [2022] J. Ren, W. Li, T. Jiang, Y. Wang, and Z. Shuai, Wiley Interdiscip. Rev. Comput. , e1614 (2022).
  * Chin  _et al._ [2010] A. W. Chin, Á. Rivas, S. F. Huelga, and M. B. Plenio, J. Math. Phys. 51, 092109 (2010).
  * Prior  _et al._ [2010] J. Prior, A. W. Chin, S. F. Huelga, and M. B. Plenio, Phys. Rev. Lett. 105, 050404 (2010).
  * Liu  _et al._ [2022b] K. T. Liu, D. N. Beratan, and P. Zhang, Phys. Rev. A 105, 032406 (2022b).
  * Broadbent  _et al._ [2012] C. J. Broadbent, J. Jing, T. Yu, and J. H. Eberly, Ann. Phys. 327, 1962 (2012).
  * Nüßeler  _et al._ [2020] A. Nüßeler, I. Dhand, S. F. Huelga, and M. B. Plenio, Physical Review B 101, 155134 (2020).
  * Li  _et al._ [2020] X. Li, J. Valdiviezo, S. D. Banziger, P. Zhang, T. Ren, D. N. Beratan, and I. V. Rubtsov, Phys. Chem. Chem. 22, 9664 (2020).
  * Lin  _et al._ [2009] Z. Lin, C. M. Lawrence, D. Xiao, V. V. Kireev, S. S. Skourtis, J. L. Sessler, D. N. Beratan, and I. V. Rubtsov, J. Am. Chem. Soc. 131, 18060 (2009).
  * Skourtis  _et al._ [1993] S. S. Skourtis, D. N. Beratan, and J. N. Onuchic, Chem. Phys. 176, 501 (1993).
  * Beratan  _et al._ [2015] D. N. Beratan, C. Liu, A. Migliore, N. F. Polizzi, S. S. Skourtis, P. Zhang, and Y. Zhang, Acc. Chem. Res. 48, 474 (2015).
  * Polyakov [2022] E. A. Polyakov, Phys. Rev. B 105, 054306 (2022).
  * Huang  _et al._ [2017] Z. Huang, Y. Fujihashi, and Y. Zhao, J. Phys. Chem. Lett. 8, 3306 (2017), pMID: 28673087, https://doi.org/10.1021/acs.jpclett.7b01247 .
  * Zhao  _et al._ [2022] Y. Zhao, K. Sun, L. Chen, and M. Gelin, Wiley Interdiscip. Rev. Comput. 12, e1589 (2022).
  * Shen  _et al._ [2022] K. Shen, K. Sun, and Y. Zhao, J. Phys. Chem. A 126, 2706 (2022).
  * de Vega  _et al._ [2015] I. de Vega, U. Schollwöck, and F. A. Wolf, Phys. Rev. B 92, 155126 (2015).
  * Lacroix  _et al._ [2021] T. Lacroix, A. Dunnett, D. Gribben, B. W. Lovett, and A. Chin, Phys. Rev. A 104, 052204 (2021).
  * Dunnett  _et al._ [2021] A. J. Dunnett, D. Gowland, C. M. Isborn, A. W. Chin, and T. J. Zuehlsdorff, J. Chem. Phys. 155, 144112 (2021).
  * Note [1] Since the spin interacts with each boson, the couplings connect the spin to the bosons, forming a star pattern.
  * Cederbaum  _et al._ [2005] L. S. Cederbaum, E. Gindensperger, and I. Burghardt, Phys. Rev. Lett. 94, 113003 (2005).
  * Tamascelli  _et al._ [2019] D. Tamascelli, A. Smirne, J. Lim, S. F. Huelga, and M. B. Plenio, Phys. Rev. Lett. 123, 090402 (2019).
  * Gindensperger  _et al._ [2006] E. Gindensperger, I. Burghardt, and L. S. Cederbaum, J. Chem. Phys. 124, 144103 (2006).
  * Noachtar  _et al._ [2022] D. D. Noachtar, J. Knörzer, and R. H. Jonsson, Phys. Rev. A 106, 013702 (2022).
  * Qiao  _et al._ [2005] S. Qiao, G. Liu, and W. Xu, in  _Advanced Signal Processing Algorithms, Architectures, and Implementations XV_, Vol. 5910, edited by F. T. Luk, International Society for Optics and Photonics (SPIE, 2005) p. 591010.
  * Nitzan [2006] A. Nitzan,  _Chemical Dynamics in Condensed Phases: Relaxation, Transfer and Reactions in Condensed Molecular Systems_ (Oxford University Press, 2006) Chap. 12.
  * Sun  _et al._ [2016] K.-W. Sun, Y. Fujihashi, A. Ishizaki, and Y. Zhao, J. Chem. Phys. 144, 204106 (2016), https://doi.org/10.1063/1.4950888 .
  * Hait  _et al._ [2017] D. Hait, M. G. Mavros, and T. Van Voorhis, J. Chem. Phys. 147, 014108 (2017), https://doi.org/10.1063/1.4990739 .
  * Acharyya  _et al._ [2020] N. Acharyya, R. Ovcharenko, and B. P. Fingerhut, J. Chem. Phys. 153, 185101 (2020), https://doi.org/10.1063/5.0027976 .
  * Berkelbach  _et al._ [2013a] T. C. Berkelbach, M. S. Hybertsen, and D. R. Reichman, J. Chem. Phys. 138, 114102 (2013a), https://doi.org/10.1063/1.4794425 .
  * Berkelbach  _et al._ [2013b] T. C. Berkelbach, M. S. Hybertsen, and D. R. Reichman, J. Chem. Phys. 138, 114103 (2013b), https://doi.org/10.1063/1.4794427 .
  * Smith and Michl [2010] M. B. Smith and J. Michl, Chem. Rev. 110, 6891 (2010), pMID: 21053979, https://doi.org/10.1021/cr1002613 .


# Phase structure of the Born-Infeld-anti-de Sitter black holes probed by non-
local observables

Xiao-Xiong Zeng1,2, Xian-Ming Liu3,4, Li-Fang Li5
1 School of Material Science and Engineering, Chongqing Jiaotong University,
Chongqing
400074, China
2 Institute of Theoretical Physics, Chinese Academy of Sciences, Beijing
100190, China
3 Center for Theoretical Physics, Massachusetts Institute of Technology,
Cambridge,
MA 02139,USA
4 Center for Theoretical Physics, School of Sciences, Hubei University for
Nationalities, Enshi, Hubei 445000, China
5 State Key Laboratory of Space Weather, Center for Space Science and Applied
Research,
Chinese Academy of Sciences, Beijing 100190, China
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract:

With the non-local observables such as two point correlation function and
holographic entanglement entropy, we probe the phase structure of the Born-
Infeld-anti-de Sitter black holes. We find for the case $bQ>0.5$, the phase
structure is similar to that of the Reissner-Nordström-AdS black hole, namely
the black hole undergoes a Hawking-Page phase transition, a first order phase
transition, and a second order phase transition. While for the case $bQ<0.5$,
we find there is a new branch for the infinitesimally small black hole so that
a pseudo phase transition emerges besides the original first order phase
transition. For the first order phase transition and the pseudo phase
transition, the equal area law is checked, and for the second order phase
transition, the critical exponent of the analogous heat capacity is obtained
in the neighborhood of the critical points. All the results show that the
phase structure of the non-local observables is the same as that of the
thermal entropy regardless of the size of the boundary region in the field
theory.

##  1 Introduction

Investigation on phase transition of an AdS space time has attracted attention
of many theoretical physicists recently. The main motivation maybe steps from
the existence of the AdS/CFT correspondence [1, 2, 3], which relates a black
hole in the AdS space time to a thermal system without gravity. For in this
case, some interesting but intractable phenomena in strongly coupled system
become to be tractable easily in the bulk. To probe these fascinating
phenomena in field theory, one should employ some non-local observables such
as two point correlation function, Wilson loop, and holographic entanglement
entropy, which are dual to the geodesic length, minimal area surface, and
minimal surface area in the bulk individually. It has been shown that these
observables can probe the non-equilibrium thermalization behavior [4, 5, 6, 7,
8, 9, 10, 11, 12, 13, 14], superconducting phase transition [15, 16, 17, 18,
19, 20, 21, 22], and cosmological singularity [23, 24].

In this paper, we intend to use the non-local observables to probe the phase
structure of the Born-Infeld-anti-de Sitter black holes. Usually, phase
structure of a black hole is understood from the viewpoint of thermodynamics.
For an uncharged AdS black hole, it has been found that there is a phase
transition between the thermal gas in AdS space and Schwarzschild AdS black
holes [25], which was interpreted as the confinement/deconfinement phase
transition in the dual gauge field theory later [26]. As the charge is endowed
with, the AdS black hole will undergo a Van der Waals-like phase transition
before it reaches the stable state in the entropy-temperature plane [27].
Specifically speaking, there exists a critical charge, and for the case that
the charge of the black hole is smaller than the critical charge, the black
hole undertakes a first order phase transition. As the charge increases to the
critical charge, the phase transition is second order, while the charge
exceeds the critical charge, there is not phase transition and the black hole
is always stable. The Van der Waals-like phase transition can be observed in
many circumstances. In [28], it was found that a 5-dimensional neutral Gauss-
Bonnet black hole demonstrates the Van der Waals-like phase transition in the
$T-\alpha$ plane, where $T$ is the Hawking temperature and $\alpha$ is the
Gauss-Bonnet coupling parameter. In [29], the Van der Waals-like phase
transition was also observed in the $Q-\Phi$ plane, where $Q$ is electric
charge and $\Phi$ is the chemical potential. Treating the negative
cosmological constant as the pressure $P$ and its conjugate as the
thermodynamical volume $V$, the Van der Waals-like phase transition has also
been observed in the $P-V$ plane recently[30, 31, 32, 33, 34, 35, 36, 37, 38].

Very recently, with the entanglement entropy as a probe, [40] investigated the
phase structure of the Reissner-Nordström-AdS black hole and found that there
was also a Van der Waals-like phase transition in the entanglement entropy-
temperature plane 111 In fact, the Van der Waals phase transition has also
been observed for a zero temperature system in the framework of holography
[39].. They also obtained the critical exponent of the heat capacity for the
second order phase transition in the neighborhood of the critical points. To
further confirm the similarity of the phase structure between the thermal
entropy and entanglement entropy, [41] checked the equal area law later and
found that it held for the first order phase transition in the entanglement
entropy-temperature plane. Now [40] has been generalized to the extended space
time [42], massive gravity [43] , as well as Weyl gravity [44], and all the
results showed that the entanglement entropy exhibited the same phase
structure as that of the thermal entropy.

In this paper, besides the entanglement entropy, we will employ the two point
correlation function to probe the phase structure of the black holes. We
choose the Born-Infeld-anti-de Sitter black holes as the gravity model, which
is a solution of the Einstein-Born-Infeld action. There have been many works
to study phase structure of the Born-Infeld-anti-de Sitter black holes [45,
46, 47, 48]. The results showed that both the Born-Infeld parameter $b$ and
charge $Q$ affect the phase structure, and to assure the existence of a non-
extremal black hole, one should impose the condition $bQ>0.5$. In this paper,
we find for the case $bQ>0.5$, the phase structure of the black hole is
similar to that of the Reissner-Nordström-AdS black hole in the entropy-
temperature plane. While for the case $bQ<0.5$, there is a novel phase
structure, which has not been observed previously. Specially speaking, a new
branch emerges compared with that of the Reissner-Nordström-AdS black hole so
that there are two unstable regions and two phase transition temperature
correspondingly. All these phase structure are probed by the nonlocal
observables such as two point correlation function as well as holographic
entanglement entropy, and phase structure of the nonlocal observables are
found to be the same as that of the thermal entropy.

Our paper is outlined as follows. In the next section, we will review the
thermodynamic properties of the Born-Infeld-anti-de Sitter black hole firstly,
and then study its phase structure in $T-S$ plane in a fixed charge ensemble.
In Section 3, we will employ the two point correlation function and
holographic entanglement entropy to probe the phase structure of the Born-
Infeld-anti-de Sitter black hole. In each subsection, the equal area law is
checked and the critical exponent of the heat capacity is obtained. The last
section is devoted to discussions and conclusions.

##  2 Thermodynamic phase transition of the Born-Infeld-anti-de Sitter black
hole

###  2.1 Review of the Born-Infeld-anti-de Sitter black hole

The 4-dimensional Born-Infeld AdS black hole is a solution of the following
action [49]

$S=\int d^{4}x\sqrt{-g}[\frac{R-3\Lambda}{16\pi G}+\frac{b^{2}}{4\pi
G}(1-\sqrt{1+\frac{2F}{b^{2}}})],$ (1)

in which $F=\frac{1}{4}F_{\mu\nu}F^{\mu\nu}$, $R$ is scalar curvature, $G$ is
the gravitational constant, $\Lambda$ is the cosmological constant that
relates to the AdS radius as $\Lambda=-3/l^{2}$, and $b$ is the Born-Infeld
parameter which relates to the string tension $\alpha^{\prime}$ as
$b=1/(2\pi\alpha^{\prime})$. Explicitly, the solutions of the Born-Infeld AdS
black hole can be written as [50, 51]

$ds^{2}=-f(r)dt^{2}+f^{-1}(r)dr^{2}+r^{2}(d\theta^{2}+\sin^{2}\theta
d\phi^{2}),$ (2)

where

$f(r)=\frac{4Q^{2}\,_{2}F_{1}\left(\frac{1}{4},\frac{1}{2};\frac{5}{4};-\frac{Q^{2}}{b^{2}r^{4}}\right)}{3r^{2}}+\frac{2b^{2}r^{2}}{3}\left(1-\sqrt{\frac{Q^{2}}{b^{2}r^{4}}+1}\right)+\frac{r^{2}}{l^{2}}-\frac{2M}{r}+1,$
(3)

in which ${}_{2}F_{1}$ is the hypergeometric function. From Eq.(3), we know
that in the limit $b\rightarrow\infty$, $Q\neq 0$, the solution reduces to the
the Reissner-Nordström-AdS black hole, and in the limit $Q\rightarrow 0$, it
reduces to the Schwarzschild AdS black hole.

The ADM mass of the black hole, defined by $f(r_{+})=0$, is given by

$M=\frac{4l^{2}Q^{2}\,_{2}F_{1}\left(\frac{1}{4},\frac{1}{2};\frac{5}{4};-\frac{Q^{2}}{b^{2}r_{+}^{4}}\right)-2b^{2}l^{2}r_{+}^{4}\sqrt{\frac{Q^{2}}{b^{2}r_{+}^{4}}+1}+2b^{2}l^{2}r_{+}^{4}+3l^{2}r_{+}^{2}+3r_{+}^{4}}{6l^{2}r_{+}},$
(4)

in which $r_{+}$ is the event horizon of the black hole. The Hawking
temperature defined by $T=\frac{f^{\prime}(r)}{4\pi}\mid_{r_{+}}$ can be
written as

$T=\frac{l^{2}\left(1-2b^{2}r_{+}^{2}\left(\sqrt{\frac{Q^{2}}{b^{2}r_{+}^{4}}+1}-1\right)\right)+3r_{+}^{2}}{4\pi
l^{2}r_{+}},$ (5)

where we have used Eq.(4). In addition, according to the Bekenstein-Hawking
entropy area relation, we can get the black hole entropy

$S=\pi r_{+}^{2}.$ (6)

Inserting Eq.(6) into Eq.(5), we can get the entropy temperature relation,
namely

$T=\frac{-2b^{2}l^{2}S\sqrt{\frac{\pi^{2}Q^{2}}{b^{2}S^{2}}+1}+2b^{2}l^{2}S+\pi
l^{2}+3S}{4\pi^{3/2}l^{2}\sqrt{S}}.$ (7)

It is obvious that besides the charge $Q$, the Born-Infeld parameter $b$ also
affects the phase structure of this space time in the $T-S$ plane.

###  2.2 Phase transition of thermal entropy

Based on Eq.(7), we will discuss the phase structure of the Born-Infeld-anti-
de Sitter black hole. We are interested in how $b$ or $Q$ affects it. For a
fixed charge $Q$, the effect of $b$ on the phase structure is shown in Figure
1. It is obvious that for the small $b$, the phase structure resembles as that
of the Schwarzschild AdS black hole, and for large $b$, it resembles as that
of the Reissner-Nordström-AdS black hole as expected. From (b) of Figure 1, we
can see that as $b$ increases, the minimum temperature of the space time
decreases, and the entropy, also the event horizon, becomes smaller. That is,
the Born-Infeld parameter $b$ promotes the formation of an AdS space time.

(a) $Q=0.15$ (b) $Q=0.2$

Figure 1: Relations between entropy and temperature for different $b$ at a
fixed $Q$. The curves from top to down correspond to the case $b$ varies from
1.5 to 3.5 with step 0.25.

More interestingly, for the smaller charge, we find there is a novel phase
transition, which is labeled by the red solid line in (a) of Figure 1. From
this curve, we know that there are two unstable regions. Namely the black hole
will undergo the following transition: unstable-stable-unstable-stable.
Compared with that of the Reissner-Nordström-AdS black hole, a new branch
emerges at the onset of the phase transition. The influence of the charge on
the phase structure for a fixed Born-Infeld parameter is plotted in Figure 2.

(a) $b=2$ (b) $b=4$

Figure 2: Relations between entropy and temperature for different $Q$ at a
fixed $b$. The curves from top to down in (a) correspond to that $Q$ varies
from 0.05 to 0.3 with step 0.02, and in (b) they correspond to the case $Q$
varies from 0.065 to 0.195 with step 0.01.

Apparently, the phase structure is similar to that of the Schwarzschild AdS
black hole for the small charges, and Reissner-Nordström-AdS black hole for
large charges. From (a) of Figure 2, we observe that the larger the charge is,
the lower the minimum temperature will be. In other words, the charge will
promote the formation of an AdS space time, which has the same effect as that
of the Born-Infeld parameter on the phase structure of the black hole.

(a) $b=5$ (b) $b=5$

Figure 3: Relations between entropy and temperature for different $Q$ at a
fixed $b$. The curves from top to down in (a) correspond to that $Q$ varies
from 0.02 to 0.2 with step 0.02, and in (b) they correspond to that $Q$ varies
from 0.015 to 0.19 with step 0.002.

For the case $b=4$ in (b) of Figure 2, we see that for some intermediate
values of charge $Q$, there is also a similar phase structure as that in (b)
in Figure 1, which is labeled by the red solid line in (b) of Figure 2. Note
that however, for the large $b$, the novel phase structure disappears, which
is shown in Figure 3. For in this case, the Born-Infeld parameter is large
enough, so that the phase structure resembles completely as that of the
Reissner-Nordström-AdS black hole.

Next, we will study detailedly the phase structure of the Born-Infeld AdS
black hole for the case $b=4$ and $b=5$ respectively. To finish it, we should
first find the critical charge for a fixed Born-Infeld parameter by the
condition

$(\frac{\partial T}{\partial S})_{Q}=(\frac{\partial^{2}T}{\partial
S^{2}})_{Q}=0.$ (8)

However, from Eq.(5), we find it is hard to get the analytical result
directly. Taking the case $b=5$ as an example, we will show how to get it
numerically next. We first plot a series of curves for different charges in
the $T-S$ plane, which is shown in (a) of Figure 3, and read off the charge
which satisfies likely the condition $(\frac{\partial T}{\partial S})_{Q}=0$.
With that rough value, we plot a bunch of curves in the $T-S$ plane once with
smaller step so that we can get the most likely value of charge. From (b) of
Figure 3, we find the critical charge should be $0.168<Q_{c}<0.17$, which are
labeled by the red dashed lines in (b) of Figure 3. Lastly, we adjust the
value of $Q_{c}$ by hand to find the only solution $S_{c}$ which satisfies
Eq.(8), which produces $Q_{c}=0.168678344129$, $S_{c}=0.510691$. Substituting
these critical values into Eq.(5), we further get the critical temperature
$T_{c}=0.259444$.

(a) $b=5$ (b) $b=5$

Figure 4: Relations between entropy and temperature for different $Q$ at
$b=5$. The red solid line and dashed line in (a) correspond to the minimum
temperature and Hawking-Page phase transition temperature. The red dashed
lines from top to down in (b) correspond to the first order phase transition
temperature and second order phase transition temperature.

Having obtained the critical charge, we can plot the isocharge curves for the
case $b=5$ in the $T-S$ plane, which is shown in Figure 4. We observe that the
black hole undergoes a Hawking-Page phase transition and a Van der Waals-like
phase transition. Specifically speaking, for the case $Q=0$, there is a
minimum temperature $T_{0}=\frac{\sqrt{3}}{2\pi}$ [52], which is indicated by
the red dashed line in (a) of Figure 4. When the temperature is higher than
$T_{0}$, there are two additional black hole branches. The small branch is
unstable while the large branch is stable. The Hawking-Page phase transition
occurs at the temperature given by $T_{1}=\frac{1}{\pi}$ [52], which is
indicated by the red dotted line. The Hawking-Page phase transition also can
be observed in the $F-T$ plane, in which $F$ is the Helmholtz free energy
defined by222Note that here we have implicitly chosen the pure AdS as the
reference spacetime because the free energy vanishes for pure AdS by this
formula.

$F=\frac{8l^{2}Q^{2}\,_{2}F_{1}\left(\frac{1}{4},\frac{1}{2};\frac{5}{4};-\frac{Q^{2}}{b^{2}r_{+}^{4}}\right)+l^{2}r_{+}^{2}\left(2b^{2}r_{+}^{2}\left(\sqrt{\frac{Q^{2}}{b^{2}r_{+}^{4}}+1}-1\right)+3\right)-3r_{+}^{4}}{12l^{2}r_{+}}.$
(9)

From (a) of Figure 5, we know that there is a minimum temperature $T_{0}$, and
above this temperature, there are two branches. The lower branch is stable
always. The Hawking-Page phase transition occurs at $T_{1}$ for in this case
the free energy vanishes.

For the case $Q\neq 0$, the phase structure is similar to that of the Van der
Waals phase transition, which is shown in (b) of Figure 4. The solid blue
lines from top to down correspond to the isocharges for the case
$Q=0.11,0.168678344129,0.21$. We can see that black holes endowed with
different charges have different phase structures. For the small charge, there
is an unstable black hole interpolating between the stable small hole and
stable large hole.

(a) $Q=0,b=5$ (b) $Q=1.1,b=5$ (c) $Q=0.168678344,b=5$

Figure 5: Relations between the free energy and temperature for different
charges at $b=5$.

The small stable hole will jump to the large stable hole at the critical
temperature $T_{a}$. As the charge increases to the critical charge, the small
hole and the large hole merge into one and squeeze out the unstable phase so
that an inflection point emerges. The divergence of the heat capacity in this
case implies that the phase transition is second order. As the charge exceeds
the critical charge, we simply have one stable black hole at each temperature.
The Van der Waals-like phase transition can also be observed from the $F-T$
relation. From (b) of Figure 5, we observe that there is a swallowtail
structure, which corresponds to the unstable phase in the top curve in (b)
Figure 4. The critical temperature $T_{a}=0.2825$ is apparently the value of
the horizontal coordinate of the junction between the small black hole and the
large black hole. As the temperature is lower than the critical temperature
$T_{a}$, the free energy of the small black hole is lowest, so the small hole
is stable. As the temperature is higher than $T_{a}$, the free energy of the
large black hole is lowest, so that the large hole dominates thereafter. The
non-smoothness of the junction indicates that the phase transition is first
order. From (c) of Figure 5, we know that there is an inflection point, which
corresponds to the inflection point in the middle curve in (b) of Figure 4.
The horizontal coordinate of the inflection point corresponds to the second
order phase transition temperature $T_{c}$.

(a) $b=4$ (b) $b=4$

Figure 6: Relations between entropy and temperature for different $Q$ at
$b=4$. The red solid line and dashed line in (a) correspond to pseudo phase
transition temperature and first order phase transition temperature. The red
dashed lines from top to down in (b) correspond to the first order phase
transition temperature and second order phase transition temperature.

Similarly, we also can study the phase structure for the case $b=4$ in the
$T-S$ plane. To finish it, we should first find the critical charge. Adopting
the same strategy as that of the case $b=5$, we find $Q_{c}=0.16987452395$,
$S_{c}=0.502683752119928$, $T_{c}=0.259166$. For the case $Q=0$, the phase
structure is the same as that in (a) of Figure 4. We will not repeat it here.
For the case $Q\neq 0$, we find besides the Van der Waals-like phase
transition, there is a novel phase transition for $Q=0.115$ as observed in
Figure 1 and Figure 2, which is plotted in (a) of Figure 6. That is, at the
beginning of the phase transition, an infinitesimally small black hole branch
emerges. It is unstable and its life is very short so that it becomes to a
small stable black hole quickly. In addition, we find there are two unstable
region, and correspondingly in the $F-T$ plane, we observe two swallowtail
structures, which is shown in (a) of Figure 7. The horizontal coordinate of
the junction of the swallowtail corresponds to the phase transition
temperature, thus besides the first order phase transition temperature between
the small black hole and large black hole, labeled by $T_{2}$, there is a new
phase transition temperature, labeled by $T_{3}$. It should be stressed that
the new phase transition can not take place actually, for above the critical
temperature $T_{2}$, the free energy of the large black hole is lowest so that
the space time is dominated by the large black hole, which can be seen in (a)
of Figure 7. Thus, we call thus phase transition pseudo phase transition.

(a) $Q=0.115,b=4$ (b) $Q=0.13,b=4$ (c) $Q=0.169874523952,b=4$

Figure 7: Relations between the free energy and temperature for different
charges at $b=4$.

We observe that there is also a Van der Waals-like phase transition, which is
plotted in (b) of Figure 6. The solid blue lines from top to down correspond
to the isocharges for the case $Q=0.13,0.16987452395,0.21$. The first order
phase transition temperature $T_{a}$ can be read off from the horizontal
coordinate of the junction of the swallowtail structure in (b) of Figure 7,
and the second order phase transition temperature $T_{c}$ can be read off from
the horizontal coordinate of the inflection point in (c) of Figure 7.

As done in [41], we will also check numerically whether Maxwell’s equal area
law holds for the first order phase transition and pseudo phase transition,
which states

$A_{L}\equiv\int_{S_{min}}^{S_{max}}T(S,Q)dS=T_{x}(S_{max}-S_{min})\equiv
A_{R},$ (10)

in which $T(S,Q)$ is defined in Eq.(7), $T_{x}$ is the phase transition
temperature, and $S_{min}$, $S_{max}$ are the smallest and largest values of
unstable region which satisfies the equation $T(S,Q)=T_{x}$. Usually for this
equation, there are three roots $S_{1},S_{2},S_{3}$, so $S_{min}=S_{1}$,
$S_{max}=S_{3}$, $T_{x}=T_{a}$. But for the case $b=4,Q=0.115$, there are four
roots. One can see from (a) of Figure 6 that $S_{min}=S_{1}$, $S_{max}=S_{3}$
for the pseudo phase transition with phase transition temperature
$T_{x}=T_{3}$, and $S_{min}=S_{2}$, $S_{max}=S_{4}$ for the first order phase
transition with phase transition temperature $T_{x}=T_{2}$. The calculated
results for the first order phase transition and pseudo phase transition are
listed in Table 1.

| $S_{min}$ | $S_{max}$ | $A_{L}$ | $A_{R}$
---|---|---|---|---
$b=5,Q=0.11,T_{a}=0.2825$ | 0.039159 | 1.7506571 | 0.4837 | 0.4835
$b=4,Q=0.13,T_{a}=0.2757$ | 0.0685377 | 1.46304 | 0.3840 | 0.3840
$b=4,Q=0.115,T_{2}=0.2817$ | 0.0210019 | 1.72464 | 0.480 | 0.480
$b=4,Q=0.115,T_{3}=0.3166$ | 0.00212107 | 0.182767 | 0.05717 | 0.05719
Table 1: Check of the equal area law in the $T-S$ plane.

From this table, we can see that $A_{L}$ equals $A_{R}$ for different $b$ in
our numeric accuracy. The equal area law therefore holds.

For the second order phase transition, we know that the heat capacity is
divergent near the critical point and the critical exponent is $-2/3$. As
stated in [40], near the critical point, there is always a linear relation

$\displaystyle\log\mid T-T_{c}\mid=3\log\mid S-S_{c}\mid+constant,$ (11)

with 3 the slope. It is not difficult to show that for the case $b=4,5$ in our
gravity model, the temperature and entropy also satisfy this linear relation
near the critical point. Next, we will take Eq.(11) as a reference to check
whether there is a similar relation for the second order phase transition in
the two point correlation function-temperature plane as well as entanglement
entropy-temperature plane.

##  3  Phase transition in the framework of holography

Having obtained the phase structure of thermal entropy of the Born-Infeld AdS
black hole in the $T-S$ plane, we will study the phase structure of two point
correlation function and entanglement entropy in the fled theory to see
whether they have the similar phase structure and critical behavior.

###  3.1 Phase structure probed by two point correlation function

According to the AdS/CFT correspondence, in the large $\Delta$ limit, the
equal time two point correlation function can be written as [53]

$\langle{\cal{O}}(t_{0},x_{i}){\cal{O}}(t_{0},x_{j})\rangle\approx
e^{-\Delta{L}},$ (12)

where $\Delta$ is the conformal dimension of scalar operator $\cal{O}$ in the
dual field theory, $L$ is the length of the bulk geodesic between the points
$(t_{0},x_{i})$ and $(t_{0},x_{j})$ on the AdS boundary. In our gravity model,
we can simply choose $(\phi=\frac{\pi}{2},\theta=0)$ and
$(\phi=\frac{\pi}{2},\theta=\theta_{0})$ as the two boundary points. Then with
$\theta$ to parameterize the trajectory, the proper length is given by

$\displaystyle
L=\int_{0}^{\theta_{0}}\mathcal{L}(r(\theta),\theta)d\theta,~{}~{}\mathcal{L}=\sqrt{\frac{\dot{r}^{2}}{f(r)}+r^{2}},$
(13)

in which $\dot{r}=dr/d\theta$. Imagining $\theta$ as time, and treating
$\mathcal{L}$ as the Lagrangian, one can get the equation of motion for
$r(\theta)$ by making use of the Euler-Lagrange equation. Then with the
following boundary conditions

$\displaystyle\dot{r}(0)=0,r(0)=r_{0},$ (14)

we can get the numeric result of $r(\theta)$. To explore whether the size of
the boundary region affects the phase structure, we will choose
$\theta_{0}=0.2,0.3$ as two examples. Note that for a fixed $\theta_{0}$, the
geodesic length is divergent, so it should be regularized by subtracting off
the geodesic length in pure AdS with the same boundary region, denoted by
$L_{0}$333$L_{0}$ can be obtained easily for there is an analytical result for
$r(\theta)$ in the pre AdS, namely
$r_{AdS}(\theta)=l[(\frac{\cos\theta}{\cos\theta_{0}})^{2}-1]^{-1/2}$ [56,
57]. . To achieve this, we are required to set a UV cutoff for each case,
which is chosen to be $r(0.199)$ and $r(0.299)$, respectively for our two
examples. The regularized geodesic length is labeled as $\delta L\equiv
L-L_{0}$. In addition, during the numerics, we will set the AdS radius $l$ to
be 1.

(a) $b=5,\theta_{0}=0.2$ (b) $b=5,\theta_{0}=0.2$

Figure 8: Isocharges in the $\delta L-T$ plane for the case
$b=5,\theta_{0}=0.2$. The red solid line and dashed line in (a) correspond to
the minimum temperature and Hawking-Page phase transition temperature. The red
solid lines from top to down in (b) correspond to the first order phase
transition temperature and second order phase transition temperature.

(a) $b=5,\theta_{0}=0.3$ (b) $b=5,\theta_{0}=0.3$

Figure 9: Isocharges in the $\delta L-T$ plane for the case
$b=5,\theta_{0}=0.3$. The red solid line and dashed line in (a) correspond to
the minimum temperature and Hawking-Page phase transition temperature. The red
solid lines from top to down in (b) correspond to the first order phase
transition temperature and second order phase transition temperature.

(a) $b=4,\theta_{0}=0.2$ (b) $b=4,\theta_{0}=0.2$

Figure 10: Isocharges in the $\delta L-T$ plane for the case
$b=4,\theta_{0}=0.2$. The red solid line and dashed line in (a) correspond to
the minimum temperature and Hawking-Page phase transition temperature. The red
solid lines from above to down in (b) correspond to the first order phase
transition temperature and second order phase transition temperature.

(a) $b=4$ (b) $b=4$

Figure 11: Isocharges in the $\delta L-T$ plane for the case
$b=4,\theta_{0}=0.3$. The red solid line and dashed line in (a) correspond to
the minimum temperature and Hawking-Page phase transition temperature. The red
solid lines from top to down in (b) correspond to the first order phase
transition temperature and second order phase transition temperature.

We plot the relation between $T$ and $\delta L$ for the case $b=5$ in Figure 8
and Figure 9 . The solid blue line in (a) corresponds to the isocharge for the
case $Q=0$, and the dashed blue lines from top to down in (b) correspond to
that $Q=0.11,0.168678344129,0.21$. From these figures, one can see that
$\delta L$ demonstrates a similar phase structure as that of the thermal
entropy in Figure 4. Namely the two point correlation function can probe both
the Hawking-Page phase transition and Van der Waals-like phase transition.
Precisely speaking, we find that the minimum temperature $T_{0}$ as well as
Hawking-Page phase transition temperature $T_{1}$ in $(a)$, the first order
phase transition temperature $T_{a}$, and second order phase transition
temperature $T_{c}$ in $(b)$ are exactly the same as those in $T-S$ plane.
This conclusion will not be affected as $\theta_{0}$ varies in a reasonable
region, which can be seen from Figure 8 and Figure 9.

With the two point correlation function, we also can probe the phase structure
of the case $b=4$, which is shown in Figure 10 and Figure 11. The solid blue
line in (a) corresponds to the isocharge for the case $Q=0.115$, and the solid
blue lines from top to down in (b) correspond to that
$Q=0.13,0.16987452395,0.21$. As the case as the thermal entropy in the $T-S$
plane, we find there ia also a novel phase structure besides the Van der
Waals-like phase transition in the $\delta L-T$ plane. There are also two
unstable region, correspondingly two phase transition temperature, labeled by
$T_{2}$, $T_{3}$. For the phase transition temperature mentioned above, it is
easy to check $T_{0}$ by locating the position of local minimum. But in order
to locate $T_{2}$, $T_{3}$, and $T_{a}$ precisely, we are required to examine
the equal area law for the first order phase transition as well as pseudo
phase transition. And to locate $T_{c}$, we should obtain the critical
exponent $-2/3$ for the second order phase transition, which are documented as
follows.

In the $\delta L-T$ plane, we define the equal area law as

$\displaystyle A_{L}\equiv\int_{\delta L_{min}}^{\delta L_{max}}T(\delta
L)d\delta L=T_{x}(\delta L_{max}-\delta L_{min})\equiv A_{R},$ (15)

in which $T(\delta L)$ is an Interpolating Function obtained from the numeric
result, $T_{x}$ is the phase transition temperature, and $\delta L_{min}$,
$\delta L_{max}$ are the smallest and largest values of the unstable region
which satisfy the equation $T(\delta L)=T_{x}$. Similar to the equal area law
in the $S-T$ plane, we know that for $b=5,Q=0.11$ and $b=4,Q=0.13$,
$T_{x}=T_{a}$, while for $b=4,Q=0.115$, $T_{x}=T_{2}$ for the first order
phase transition and $T_{x}=T_{3}$ for the pseudo phase transition. For
different $b$ and $\theta_{0}$, the calculated results are listed in Table 2.

| $\delta L_{min}$ | $\delta L_{max}$ | $A_{L}$ | $A_{R}$
---|---|---|---|---
$b=5,Q=0.11,T_{a}=0.2825,\theta_{0}=0.2$ | 0.00005901 | 0.0003329 | 0.0000777 | 0.0000774
$b=5,Q=0.11,T_{a}=0.2825,\theta_{0}=0.3$ | 0.0001982 | 0.001124 | 0.000263 | 0.000261
$b=4,Q=0.13,T_{a}=0.2757,\theta_{0}=0.2$ | 0.00007227 | 0.0002890 | 0.0000598 | 0.0000596
$b=4,Q=0.13,T_{a}=0.2757,\theta_{0}=0.3$ | 0.0002453 | 0.0009755 | 0.000201 | 0.000201
$b=4,Q=0.115,T_{2}=0.2817,\theta_{0}=0.2$ | 0.00005509 | 0.0003290 | 0.0000774 | 0.0000772
$b=4,Q=0.115,T_{2}=0.2817,\theta_{0}=0.3$ | 0.0001939 | 0.001110 | 0.000259 | 0.000258
$b=4,Q=0.115,T_{3}=0.3166,\theta_{0}=0.2$ | 0.00005485 | 0.00008727 | 0.0000103 | 0.0000103
$b=4,Q=0.115,T_{3}=0.3166,\theta_{0}=0.3$ | 0.0001842 | 0.00003476 | 0.0000349 | 0.0000348
Table 2: Check of the equal area law in the $T-\delta L$ plane.

It is obvious that $A_{L}$ equals $A_{R}$ for different $b$ and $\theta_{0}$
in a reasonable accuracy. That is, in the $T-\delta L$ plane, the equal area
law holds, and it is independent of the Born-Infeld parameter as well as the
size of the boundary region.

In order to get the critical exponent for the second order phase transition in
the $T-\delta L$ plane, we should first define an analogous heat capacity

$\displaystyle C=T\frac{\partial\delta L}{\partial T}.$ (16)

Provided a similar relation as showed in Eq.(11) is satisfied, one can get the
critical exponent immediately. So next, we are interested in the logarithm of
the quantities $T-T_{c}$, $\delta L-\delta L_{c}$, in which $T_{c}$ is the
critical temperature mentioned previously, and $L_{c}$ is obtained numerically
by the equation $T(\delta L)=T_{c}$. We plot the relation between $\log\mid
T-T_{c}\mid$ and $\log\mid\delta L-\delta L_{c}\mid$ for different
$\theta_{0}$ and $b$ in Figure 12, where these straight lines can be fitted as

$\log\mid T-T_{c}\mid=\begin{cases}25.0567+3.0113\log\mid\delta L-\delta
L_{c}\mid,&$for$~{}b=5,\theta_{0}=0.2,\\\ 21.207+3.0402\log\mid\delta L-\delta
L_{c}\mid,&$for$~{}b=5,\theta_{0}=0.3,\\\ 23.8282+3.01531\log\mid\delta
L-\delta L_{c}\mid,&$for$~{}b=4,\theta_{0}=0.2,\\\
20.3704+3.03817\log\mid\delta L-\delta
L_{c}\mid,&$for$~{}b=4,\theta_{0}=0.3.\\\ \end{cases}$ (17)

It is obvious that the slope is about 3, which indicates that the critical
exponent is $-2/3$ nearly for the analogous heat capacity, and the phase
transition is also second order with the phase transition temperature $T_{c}$.

(a) $b=5,\theta_{0}=0.2$ (b) $b=5,\theta_{0}=0.3$ (c) $b=4,\theta_{0}=0.2$ (d)
$b=4,\theta_{0}=0.3$

Figure 12: Relations between $\log\mid T-T_{c}\mid$ and $\log\mid\delta
L-\delta L_{c}\mid$ for different $b$ and $\theta_{0}$.

###  3.2 Phase structure probed by holographic entanglement entropy

According to the formula in [54, 55], holographic entanglement entropy can be
given by the area $A_{\Sigma}$ of a minimal surface $\Sigma$ anchored on the
boundary entangling surface $\partial\Sigma$, namely

$S=\frac{A_{\Sigma}(t)}{4G},$ (18)

where $G$ is the Newton’s constant. We will take the region $\Sigma$ to be a
spherical cap on the boundary delimited by $\theta\leq\theta_{0}$444To avoid
contamination of the entanglement entropy by the thermal entropy, we will
choose a small entangling region as pointed out in [40, 58]. . Then based on
the definition of area and (2), (18) can be rewritten as

$\displaystyle
S=\frac{\pi}{2}\int_{0}^{\theta_{0}}r\sin\theta\sqrt{\frac{(r^{\prime})^{2}}{f(r)}+r^{2}},$
(19)

in which $r^{\prime}=dr/d\theta$. Making use of the Euler-Lagrange equation,
one can get the equation of motion of $r(\theta)$. With the boundary
conditions in Eq.(14), we can obtain the numeric result of $r(\theta)$
immediately. Note that the entanglement entropy is divergent at the boundary,
so it should be regularized by subtracting off the entanglement entropy in
pure AdS with the same entangling surface and boundary values. We label the
regularized entanglement entropy as $\delta S$. We choose the size of the
boundary region as $\theta_{0}=0.2$ and set the UV cutoff in the dual field
theory to be $r(0.199)$.

The isocharges for the cases $b=5$ and $b=4$ on the entanglement entropy-
temperature plane are plotted in Figure 13 and Figure 14. As the same as that
of the two point correlation function, for $b=5$, we are interested in that
the black holes have charges $Q=0,0.11,0.168678344129,0.21$, and for $b=4$, we
choose $Q=0.115,0.13,0.16987452395,0.21$. From Figure 13, we know that the
entanglement entropy exhibits a Hawking-Page phase transition, and a Van der
Waals-like phase transition as the charge increases form 0 to 0.21. The
entanglement entropy also exhibits the novel phase structure as that of the
thermal entropy as well as two point correlation function, which is shown in
(a) of Figure 14.

(a) $b=5,\theta_{0}=0.2$ (b) $b=5,\theta_{0}=0.2$

Figure 13: Isocharges in the $\delta S-T$ plane for the case
$b=5,\theta_{0}=0.2$. The red solid line and dashed line in (a) correspond to
the minimum temperature and Hawking-Page phase transition temperature. The red
solid lines from top to down in (b) correspond to the first order phase
transition temperature and second order phase transition temperature.

(a) $b=4,\theta_{0}=0.2$ (b) $b=4,\theta_{0}=0.2$

Figure 14: Isocharges in the $\delta S-T$ plane for the case
$b=4,\theta_{0}=0.2$. The red solid line and dashed line in (a) correspond to
the minimum temperature and Hawking-Page phase transition temperature. The red
solid lines from above to down in (b) correspond to the first order phase
transition temperature and second order phase transition temperature.

We will also employ Maxwell’s equal area law to locate the first order phase
transition temperature, namely $T_{a}$ in (b) of Figure 13 and Figure 14,
$T_{2}$ in (a) of Figure 14, and pseudo phase transition temperature $T_{3}$.
In the $\delta S-T$ plane, we define the equal area law as

$\displaystyle A_{L}\equiv\int_{\delta S_{min}}^{\delta S_{max}}T(\delta
S)d\delta S=T_{x}(\delta S_{max}-\delta S_{min})\equiv A_{R},$ (20)

in which $T(\delta S)$ is an Interpolating Function obtained from the numeric
result, and $\delta S_{min}$, $\delta S_{max}$ are the smallest and largest
roots of the unstable regions which satisfy the equation $T(\delta S)=T_{x}$.
Surely the values of the phase transition temperature $T_{x}$ depends on the
choice of $Q$ and $b$. For different $Q$ and $b$, the results of $\delta
S_{min}$, $\delta S_{max}$ and $A_{L}$, $A_{R}$ are listed in Table (3).

| $\delta S_{min}$ | $\delta S_{max}$ | $A_{L}$ | $A_{R}$
---|---|---|---|---
$b=5,Q=0.11,T_{a}=0.2825,\theta_{0}=0.2$ | 0.00064474 | 0.0036781 | 0.000860 | 0.000857
$b=4,Q=0.13,T_{a}=0.2757,\theta_{0}=0.2$ | 0.00076078 | 0.0032308 | 0.000681 | 0.000681
$b=4,Q=0.115,T_{2}=0.2817,\theta_{0}=0.2$ | 0.00060499 | 0.0036363 | 0.000857 | 0.000853
$b=4,Q=0.115,T_{3}=0.3166,\theta_{0}=0.2$ | 0.00060237 | 0.00095971 | 0.000113 | 0.000112
Table 3: Check of the equal area law in the $T-\delta S$ plane.

It is obvious that $A_{L}$ equals nearly $A_{R}$ for a fixed $Q$ and $b$. That
is, the equal area law is also valid for the first order phase transition and
pseudo phase transition in the entanglement entropy-temperature plane. This
result is the same as that of the thermal entropy as well as two point
correlation function.

To confirm that $T_{c}$ is the second order phase transition temperature in
the $\delta S-T$ plane, we should explore whether it satisfies a similar
relation as Eq.(11). For different $b$ and $Q$, the relation between $\log\mid
T-T_{c}\mid$ and $\log\mid\delta S-\delta S_{c}\mid$ are plotted in Figure 14,
in which $S_{c}$ is the critical entropy obtained numerically by the equation
$T(\delta S)=T_{c}$. The analytical result of these curves can be fitted as

$\log\mid T-T_{c}\mid=\begin{cases}17.3605+2.91954\log\mid\delta S-\delta
S_{c}\mid,&$for$~{}b=5,\theta_{0}=0.2,\\\ 17.2131+2.93553\log\mid\delta
S-\delta S_{c}\mid,&$for$~{}b=4,\theta_{0}=0.2.\\\ \end{cases}$ (21)

We find for a fixed $b$ and $Q$, the slope is always about 3, which is
consistent with that of the thermal entropy. That is, the entanglement entropy
also exhibits a second order phase transition, with phase transition
temperature $T_{c}$, as that of the thermal entropy and two point correlation
function.

(a) $b=5,\theta_{0}=0.2$ (b) $b=4,\theta_{0}=0.2$

Figure 15: Relations between $\log\mid T-T_{c}\mid$ and $\log\mid\delta
S-\delta S_{c}\mid$ for different $b$ and $\theta_{0}$.

##  4 Concluding remarks

Investigation on the phase structure of a back hole is of great importance for
it reveals whether and how a black hole emerges. Usually, it is realized by
studying the relation between some thermodynamic quantities in a fixed
ensemble. In this paper, we found that the phase structure of a back hole also
can be probed by the two point correlation function and holographic
entanglement entropy, which provides a new strategy to understand the phase
structure of the black holes from the viewpoint of holography.

Specially speaking, we first investigated the phase structure of the Born-
Infeld-anti-de Sitter black holes in the $T-S$ plane in the fixed charge
ensemble. We found that the phase structure of the black hole depends on not
only the value of $Q$ or $b$, but also the combination of $bQ$. For the case
$b=5$, the black hole resembles as the Reissner-Nordström-AdS black hole, and
it undergoes the Van der Waals-like phase transition as its charge satisfies
the non-extremal condition $bQ>0.5$. For the case $b=4$, besides the Van der
Waals-like phase transition, we also observed a novel phase structure for the
case $bQ<0.5$. With the two point correlation function and holographic
entanglement entropy, we further probed the phase structure of the Born-
Infeld-anti-de Sitter black holes and found that both the probes exhibited the
same phase structure as that of the thermal entropy for cases $bQ>0.5$ and
$bQ<0.5$ regardless of the size of the boundary region. This conclusion was
reinforced by checking the equal area law for the first order phase transition
as well as pseudo phase transition, and calculating the critical exponent of
the analogous heat capacity for the second order phase transition.

In previous investigation on the phase structure of the Born-Infeld-anti-de
Sitter black holes, authors concentrated mainly on that of $bQ>0.5$ [45, 46,
47, 48]. In this paper we found that the phase structure for the case $bQ<0.5$
is also interesting. As $bQ\ll 0.5$, the black hole resembles as the
Schwarzschild AdS black hole as expected, and there is a minimum temperature.
We found the larger the charge or the Born-Infeld parameter is, the smaller
the minimum temperature will be. That is, both the charge and the Born-Infeld
parameter promote the formation of an AdS black hole. As $bQ$ approaches to
0.5, a new extremal small black hole branch emerges compared with that of the
Reissner-Nordström-AdS black hole, so that there are two unstable regions and
correspondingly two phase transition temperature. The high temperature phase
transition was found to be pseudo for in this case the space time was
dominated by the large black hole, which was observed from the $F-T$ relation.

## Acknowledgements

We would like to thank Professor Rong-Gen Cai for his helpful suggestions.
This work is supported by the National Natural Science Foundation of China
(Grant No. 11365008, Grant No. 11205226, Grant No. 11575270), and Natural
Science Foundation of Education Committee of Chongqing (Grant No. KJ1500530).

## References

  * [1] J. M. Maldacena, _Large N limit of superconformal field theories and supergravity_ , Int. J. Theor. Phys. 38, 1113 (1999)
  * [2] E. Witten,_Anti-de Sitter space and holography_ , Adv. Theor. Math. Phys. 2, 253 (1998) [hep-th/9802150]
  * [3] S. S. Gubser, I. R. Klebanov, A. M. Polyakov, _Gauge theory correlators from noncritical string theory_ , Phys. Lett. B 428, 105 (1998) [hep-th/9802109]
  * [4] V. Balasubramanian et al.,_Thermalization of Strongly Coupled Field Theories_ , Phys. Rev. Lett. 106, 191601 (2011) [arXiv:1012.4753 [hep-th]]
  * [5] V. Balasubramanian et al., _Holographic Thermalization_ , Phys. Rev. D 84, 026010 (2011) [arXiv:1103.2683 [hep-th]]
  * [6] D. Galante and M. Schvellinger, _Thermalization with a chemical potential from AdS spaces_ , JHEP 1207, 096 (2012) [arXiv:1205.1548 [hep-th]]
  * [7] E. Caceres and A. Kundu, _Holographic Thermalization with Chemical Potential_ , JHEP 1209, 055 (2012) [arXiv:1205.2354 [hep-th]]
  * [8] X. X. Zeng, B. W. Liu, _Holographic thermalization in Gauss-Bonnet gravity_ , Phys. Lett. B 726, 481 (2013) [arXiv:1305.4841 [hep-th]]
  * [9] X. X. Zeng, X. M. Liu, B. W. Liu, _Holographic thermalization with a chemical potential in Gauss-Bonnet gravity_ , JHEP 03, 031 (2014) [arXiv:1311.0718 [hep-th]]
  * [10] X. X. Zeng, D. Y. Chen, L. F. Li,_Holographic thermalization and gravitational collapse in the spacetime dominated by quintessence dark energy_ , Phys. Rev. D 91, 046005 (2015) [arXiv:1408.6632 [hep-th]]
  * [11] H. Liu and S. J. Suh,_Entanglement Tsunami: Universal Scaling in Holographic Thermalization_ , Phys. Rev. Lett. 112, 011601 (2014) [arXiv:1305.7244 [hep-th]].
  * [12] S. J. Zhang, E. Abdalla, _Holographic Thermalization in Charged Dilaton Anti-de Sitter Spacetime_ , Nuclear Physics B 896, 569 (2015) [arXiv:1503.07700 [hep-th]]
  * [13] A. Buchel, R. C. Myers, A. v. Niekerk,_Nonlocal probes of thermalization in holographic quenches with spectral methods_ , JHEP 02,017 (2015) [arXiv:1410.6201 [hep-th]]
  * [14] B. Craps et al., _Gravitational collapse and thermalization in the hard wall model_ , JHEP 02, 120 (2014) [arXiv:1311.7560 [hep-th]]
  * [15] T. Albash, C. V. Johnson, _Holographic Studies of Entanglement Entropy in Superconductors_ , JHEP 1205, 079 (2012) [arXiv:1202.2605 [hep-th]]
  * [16] R. G. Cai, S. He, L. Li, Y. L. Zhang,_Holographic Entanglement Entropy in Insulator/Superconductor Transition_ , JHEP 1207, 088 (2012) [arXiv:1203.6620 [hep-th]]
  * [17] R. G. Cai, L. Li, L. F. Li, R. K. Su, _Entanglement Entropy in Holographic P-Wave Superconductor/Insulator Model_ , JHEP 1306, 063 (2013) [arXiv:1303.4828 [hep-th]]
  * [18] L. F. Li, R. G. Cai, L. Li, C. Shen,_Entanglement entropy in a holographic p-wave superconductor model_ , Nucl. Phys. B 894, 15-28 (2015) [arXiv:1310.6239 [hep-th]]
  * [19] R. G. Cai, S. He, L. Li, Y. L. Zhang,_Holographic Entanglement Entropy in Insulator/Superconductor Transition_ , JHEP 1207, 088 (2012) [arXiv:1203.6620 [hep-th]]
  * [20] X. Bai, B. H. Lee, L. Li, J. R. Sun, H. Q. Zhang, _Time Evolution of Entanglement Entropy in Quenched Holographic Superconductors_ , JHEP 040, 66 (2015) [arXiv:1412.5500 [hep-th]]
  * [21] R. G. Cai, L. Li, L. F. Li, R. Q. Yang, _Introduction to Holographic Superconductor Models_ , Sci China-Phys Mech Astron, 58, 060401 (2015) [arXiv:1502.00437 [hep-th]]
  * [22] Y. Ling, P. Liu, C. Niu, J. P. Wu, Z. Y. Xian, _Holographic Entanglement Entropy Close to Quantum Phase Transitions_ [arXiv:1502.03661 [hep-th]]
  * [23] N. Engelhardt, T. Hertog, G. T. Horowitz, _Holographic Signatures of Cosmological Singularities_ , Phys. Rev. Lett. 113, 121602 (2014) [arXiv:1404.2309 [hep-th]]
  * [24] N. Engelhardt, T. Hertog, G. T. Horowitz, Further Holographic Investigations of Big Bang Singularities, JHEP 1507, 044 (2015) [arXiv:1503.08838 [hep-th]]
  * [25] S. W. Hawking and D. N. Page, _Thermodynamics Of Black Holes In Anti-De Sitter Space_ , Commun. Math. Phys. 87, 577 (1983)
  * [26] E. Witten, _Anti-de Sitter space, thermal phase transition, and confinement in gauge theories_ , Adv. Theor. Math. Phys. 2, 505 (1998) [arXiv:hep-th/9803131]
  * [27] A. Chamblin, R. Emparan, C. V. Johnson and R. C. Myers, _Charged AdS black holes and catastrophic holography_ , Phys. Rev. D 60, 064018 (1999) [hep-th/9902170]
  * [28] R. G. Cai,_Gauss-Bonnet Black Holes in AdS Spaces_ , Phys. Rev. D 65, 084014 (2002) [arXiv:hep-th/0109133]
  * [29] C. Niu, Y. Tian, and X. N. Wu, _Critical Phenomena and Thermodynamic Geometry of RN-AdS Black Holes_ , Phys. Rev. D 85, 024017 (2012) [arXiv:1104.3066 [hep-th]]
  * [30] D. Kastor, S. Ray, J. Traschen, _Enthalpy and the Mechanics of AdS Black Holes_ , Class. Quant. Grav. 26, 195011 (2009) [arXiv:0904.2765 [hep-th]]
  * [31] D. Kubiznak and R. B. Mann,_P-V criticality of charged AdS black holes_ , JHEP 1207, 033 (2012) [arXiv:1205.0559 [hep-th]]
  * [32] J. Xu, L. M. Cao, Y. P. Hu, _P-V criticality in the extended phase space of black holes in massive gravity_ , Phys. Rev. D 91, 124033 (2015) [arXiv:1506.03578 [gr-qc]]
  * [33] R. G. Cai, L. M. Cao, L. Li, R. Q. Yang, _P-V criticality in the extended phase space of Gauss-Bonnet black holes in AdS space_ , JHEP 09, 005 (2013) [arXiv:1306.6233 [gr-qc]]
  * [34] S. H. Hendi, A. Sheykhi, S. Panahiyan, B. Eslam Panah, _Phase transition and thermodynamic geometry of Einstein-Maxwell-dilaton black holes_ , Phys. Rev. D 92, 064028 (2015) [arXiv:1509.08593 [hep-th]]
  * [35] R. A. Hennigar, W. G. Brenna, R. B. Mann, _P-V Criticality in Quasitopological Gravity_ , JHEP 1507, 077 (2015) [arXiv:1505.05517 [hep-th]]
  * [36] S. W. Wei, Y. X. Liu, _Clapeyron equations and fitting formula of the coexistence curve in the extended phase space of charged AdS black holes_ , Phys. Rev. D 91, 044018 (2015) [arXiv:1411.5749 [hep-th]]
  * [37] J. X. Mo, W. B. Liu, _P-V Criticality of Topological Black Holes in Lovelock-Born-Infeld Gravity_ , Eur. Phys. J. C 74, 2836 (2014) [arXiv:1401.0785 [hep-th]]
  * [38] E. Spallucci and A. Smailagic, _Maxwell’s equal area law for charged Anti-deSitter black holes_ , Phys. Lett. B 723, 436 (2013) [arXiv:1305.3379 [hep-th]]
  * [39] F. Bigazzi, A. L. Cotrone, C. Nunez, A.Paredes, _Heavy quark potential with dynamical flavors: a first order transition_ , Phys. Rev. D 78, 114012 (2008) [arXiv:0806.1741 [hep-th]]
  * [40] C. V. Johnson, _Large N Phase Transitions, Finite Volume, and Entanglement Entropy_ , JHEP 1403, 047 (2014) [arXiv:1306.4955 [hep-th]]
  * [41] P. H. Nguyen, _An equal area law for the van der Waals transition of holographic entanglement entropy_ , JHEP 12, 139 (2015) [arXiv:1508.01955 [hep-th]]
  * [42] E. Caceres, P. H. Nguyen, J. F. Pedraza, _Holographic entanglement entropy and the extended phase structure of STU black holes_ , JHEP 1509, 184, (2015) [arXiv:1507.06069 [hep-th]]
  * [43] X. X. Zeng, H. Zhang, L. F. Li, _Phase transition of holographic entanglement entropy in massive gravity_ , arXiv:1511.00383 [gr-qc]
  * [44] A. Dey, S. Mahapatra, T. Sarkar,_Thermodynamics and Entanglement Entropy with Weyl Corrections_ , arXiv:1512.07117 [hep-th]
  * [45] Y. S. Myung, Y. W. Kim, Y. J. Park,_Thermodynamics and phase transitions in the Born-Infeld-anti-de Sitter black holes_ , Phys. Rev. D 78, 084002 (2008) [arXiv:0805.0187 [gr-qc]]
  * [46] A. Lala, D. Roychowdhur,_Ehrenfest s scheme and thermodynamic geometry in Born-Infeld AdS black holes_ , Phys. Rev. D 86, 084027 (2012) [arXiv:1111.5991 [gr-qc]]
  * [47] S. H. Hendi, S. Panahiyan, _Thermodynamic instability of topological black holes in Gauss-Bonnet gravity with a generalized electrodynamics_ , Phys. Rev. D 90, 124008 (2014) [arXiv:1501.05481 [gr-qc]]
  * [48] D. C. Zou, S. J. Zhang, B. Wang, _Critical behavior of Born-Infeld AdS black holes in the extended phase space thermodynamics_ , Phys. Rev. D 89, 044002 (2014) [arXiv:1311.7299 [hep-th]]
  * [49] D. A. Rasheed,_Non-Linear Electrodynamics: Zeroth and First Laws of Black Hole Mechanics_ , arXiv:hep-th/9702087
  * [50] R. G. Cai, D. W. Pang, A. Wang, _Born-Infeld Black Holes in (A)dS Spaces_ , Phys. Rev. D 70, 124034 (2004) [arXiv:hep-th/0410158]
  * [51] T. K. Dey, _Born-Infeld black holes in the presence of a cosmological constant,_ Phys. Lett. B 595, 484 (2004) [arXiv:hep-th/0406169].
  * [52] R. Banerjee, D. Roychowdhury,_Thermodynamics of phase transition in higher dimensional AdS black holes_ , JHEP 11, 004 (2011) [arXiv:1109.2433 [gr-qc]]
  * [53] V. Balasubramanian and S. F. Ross,_Holographic particle detection_ , Phys. Rev. D 61, 044007 (2000) [arXiv:hep-th/9906226]
  * [54] S. Ryu and T. Takayanagi, _Holographic derivation of entanglement entropy from AdS/CFT_ , Phys. Rev. Lett. 96, 181602 (2006) [hep-th/0603001]
  * [55] S. Ryu, T. Takayanagi,_Aspects of Holographic Entanglement Entropy_ , JHEP 0608, 045 (2006) [hep-th/0605073]
  * [56] D. D. Blanco, H. Casini, L.-Y. Hung, and R. C. Myers, _Entropy and Holography_ , JHEP 08, 060 (2013) [arXiv:1305.3182 [hep-th]]
  * [57] H. Casini, M. Huerta, and R. C. Myers, _a derivation of holographic entanglement entropy_ , JHEP 1105, 036 (2011) arXiv:1102.0440 [hep-th]
  * [58] V. E. Hubeny, H. Maxfield, M. Rangamani, and E. Tonni,_Holographic entanglement plateaux_ , JHEP 08, 092 (2013) [arXiv:1306.4004 [hep-th]]