Tensor network influence functionals in the continuous-time limit: connections to
quantum embedding, bath discretization, and higher-order time propagation
Gunhee Park,
1
Nathan Ng,
2
David R. Reichman,
2
and Garnet Kin-Lic Chan
3
1
Division of Engineering and Applied Science, California Institute of Technology, Pasadena, CA 91125, USA
2
Department of Chemistry, Columbia University, New York, New York 10027, United States
3
Division of Chemistry and Chemical Engineering,
California Institute of Technology, Pasadena, California 91125, USA
We describe two developments of tensor network influence functionals (in particular, influence
functional matrix product states (IF-MPS)) for quantum impurity dynamics within the fermionic
setting of the Anderson impurity model. The first provides the correct extension of the IF-MPS to
continuous time by introducing a related mathematical object, the boundary influence functional
MPS. The second connects the dynamics described by a compressed IF-MPS to that of a quantum
embedding method with a time-dependent effective bath undergoing non-unitary dynamics. Using
these concepts, we implement higher-order time propagators for the quench dynamics of the Ander-
son impurity model within the boundary IF-MPS formalism. The calculations illustrate the ability
of the current formulation to efficiently remove the time-step error in standard discrete-time IF-MPS
implementations as well as to interface with state-vector propagation techniques. They also show
the advantages of IF-MPS dynamics, with its associated highly compact effective bath dynamics,
over state-vector propagation with a static bath discretization.
I. INTRODUCTION
Quantum impurity models, such as the Anderson im-
purity model (AIM), consist of an interacting impu-
rity coupled to a (possibly non-interacting) bath. They
provide simple settings in which to study quantum
many-body physics, and the challenge of describing non-
equilibrium dynamics in such models has spurred the de-
velopment of many computational techniques [1–12].
One way to describe quantum impurity dynamics is
to use the equation of motion of the impurity reduced
density operator obtained by tracing out the bath. The
influence of the bath is clearly expressed in a path integral
language via the Feynman-Vernon influence functional
(IF) [13], which reweights the paths in the path integral.
Although the influence functional for a non-interacting
bath with linear coupling can be expressed in a compact
mathematical form, its effect on the impurity dynamics
must still be numerically approximated in practice [14–
25].
Recently, tensor network methods have been explored
within the IF language to overcome some limitations of
other numerical IF techniques, and in particular, to cap-
ture longer time memory effects [26–36]. By treating the
IF as a temporal wavefunction expressed as a tempo-
ral matrix product state (MPS), one can exploit low en-
tanglement in time, leading to a compact representation
of certain temporal correlations. Such tensor network-
based IF methods have been successfully demonstrated
in several contexts, including the spin-boson model [28–
32], hard-core bosons [31], one-dimensional spin chains
[26, 27, 33], and more recently, in the context of the AIM
[34–36].
The current work is concerned with two developments
of the tensor network IF approach. Although the ideas
generalize to IFs of any linearly-coupled non-interacting
bath, for concreteness we work specifically with the AIM.
The first development removes the time-step error in the
standard IF treatment obtained from a Trotterized rep-
resentation of the path integral. Formulating the prob-
lem in the continuous-time limit yields a continuous MPS
version of the IF-MPS. We show that the standard IF-
MPS in fact does not have a useful continuous-time limit,
and instead one must consider a closely related object,
the boundary IF-MPS, to define this limit. We provide
an explicit construction of the continuous-time boundary
IF-MPS for the AIM.
The second development is concerned with the rela-
tionship between the (boundary) IF-MPS and discrete
bath dynamics. Standard state-vector methods for the
dynamics of the AIM use a discrete bath with a fixed
set of bath energies and couplings [1–5]. We show that
impurity dynamics defined by an IF-MPS of fixed bond
dimension is equivalent to propagating in a space of ef-
fective bath orbitals in Liouville space, where the bath
energies and couplings dynamically vary with time. This
naturally connects tensor network influence functionals
to dynamical quantum embedding theories, such as real-
time density matrix embedding, which also utilize a dy-
namical set of bath energies and couplings [37]. An im-
portant difference, however, is that the IF-MPS defines
a non-unitary dynamics in the bath.
We implement the boundary IF-MPS propagation for
the AIM in both the discrete-time and continuous-time
formulations. Directly comparing to standard static bath
discretizations, we show that the IF-MPS method con-
verges extremely rapidly with respect to the effective
bath size, and shows none of the discretization artifacts
of standard bath discretizations. Further, the combina-
tion of the above two insights suggests that the bound-
ary IF-MPS dynamics in the continuum limit can be ef-
ficiently implemented using standard state-vector time-
propagation techniques. We use this to implement a
arXiv:2401.12460v1 [cond-mat.str-el] 23 Jan 2024
2
high-order Runge-Kutta time-propagator for boundary
IF-MPS dynamics and demonstrate the high-order error
with time-step. In contrast, we show that higher than
first-order Trotter methods in the standard discrete-time
IF-MPS still suffer from first-order time-step errors, due
to the IF-MPS bond truncation. By Trotterizing the cor-
rect continuous-time dynamics, we derive a corrected ver-
sion of the Trotterized error with the correct time-order
scaling, significantly improving on the second-order Trot-
ter formulation currently used in tensor network IF ap-
proaches.
The paper is organized as follows.
In Sec. II, we
first recapitulate the formulation of the IF-MPS for non-
interacting fermionic baths, describing in detail the for-
malism we use here in terms of number-conserving Slater
determinants. Through this picture, we establish the
connection between the IF-MPS dynamics and Liouville
state-vector propagation of an impurity coupled to a set
of effective bath orbitals, that is, the dynamics of a quan-
tum embedding of the impurity. In Sec. III, we analyze
the continuous-time limit of the IF-MPS and show that
the correct object to consider is the boundary IF-MPS
and we provide its continuous-time limit. Using the state-
vector formalism, we rewrite the continuous-time propa-
gation in terms of a differential equation of motion for the
Liouville state-vector in the quantum embedding space.
In Sec. IV, we provide numerical results using the bound-
ary IF-MPS for the single impurity Anderson model and
compare the discretization errors associated with a static
set of bath orbitals with the dynamic set of bath orbitals
defined by the IF-MPS. We further analyze the time-step
errors from both the discrete-time and continuous-time
formulations. In Sec. V, we conclude with discussions of
some implications of our results.
II. INFLUENCE FUNCTIONAL THEORY IN
DISCRETE TIME
A. Influence Functional for Single Impurity
Anderson Model
We consider a single impurity Anderson model,
ˆ
H
=
ˆ
H
S
+
ˆ
H
SB
+
ˆ
H
B
,
ˆ
H
S
=
U
ˆ
n
↑
ˆ
n
↓
+
X
σ
ε
σ
ˆ
n
σ
,
ˆ
H
SB
=
X
i,σ
t
i
ˆ
c
†
i,σ
ˆ
d
σ
+ h.c.
,
ˆ
H
B
=
X
i,σ
E
i
ˆ
c
†
i,σ
ˆ
c
i,σ
,
(1)
where
ˆ
d
†
σ
(ˆ
c
†
i,σ
) creates a fermion of spin
σ
∈{↑
,
↓}
in the
impurity (bath) orbitals and ˆ
n
σ
=
ˆ
d
†
σ
ˆ
d
σ
is the number
density operator of the impurity orbital of spin
σ
. Note
that
ˆ
H
SB
and
ˆ
H
B
are of noninteracting (quadratic) form.
Hereafter, we assume a discrete and finite set of bath or-
bitals and also assume that the impurity is initially de-
coupled from the bath, ˆ
ρ
(0) = ˆ
ρ
S
(0)
⊗
ˆ
ρ
B
. We further
assume a Gaussian initial bath, in particular, the ther-
mal state ˆ
ρ
B
∝
e
−
β
ˆ
H
B
with inverse temperature
β
. For
quadratic operators, we omit the hats when expressing
their matrix elements in a single-particle basis.
While the time evolution of the density operator can
be described by the von Neumann equation,
i
d
dt
ˆ
ρ
= [
ˆ
H,
ˆ
ρ
]
,
(2)
we instead will adopt the super-fermion representation of
Liouville space [17, 38–41]. This uses a super-Fock space
with twice the number of orbitals as the original Hilbert
space, obtained by applying a particle-hole transforma-
tion to the bra of the density operator, and assuming
the resulting operator acts on a vacuum to produce a
state. We denote states in the super-Fock space using a
double bra/ket notation. For example, the initial den-
sity operator is written as
|
ρ
(0)
⟩⟩
=
|
ρ
S
(0)
⟩⟩⊗|
ρ
B
⟩⟩
. The
von Neumann equation now becomes a Hamiltonian time
evolution with the Liouville operator,
ˆ
L
,
i
d
dt
|
ρ
⟩⟩
=
ˆ
L
|
ρ
⟩⟩
.
(3)
The Liouville operator for the Anderson impurity model
takes the form (using tildes on the operators from the
particle-hole transformed Fock space),
ˆ
L
=
ˆ
L
S
+
ˆ
L
SB
+
ˆ
L
B
,
ˆ
L
S
=
U
ˆ
n
↑
ˆ
n
↓
−
U
(1
−
ˆ
̃
n
↑
)(1
−
ˆ
̃
n
↓
)
+
X
σ
ε
σ
(ˆ
n
σ
+
ˆ
̃
n
σ
)
,
ˆ
L
SB
=
X
i,σ
t
i
ˆ
c
†
i,σ
ˆ
d
σ
+
t
i
ˆ
̃
c
†
i,σ
ˆ
̃
d
σ
+ h.c.
,
ˆ
L
B
=
X
i,σ
E
i
(ˆ
c
†
i,σ
ˆ
c
i,σ
+
ˆ
̃
c
†
i,σ
ˆ
̃
c
i,σ
)
,
(4)
where we have omitted all constant terms. We will collec-
tively refer to the tilde and non-tilde creation operators
as ˆ
a
†
i,σ
,
ˆ
a
†
i,σ
=
(
ˆ
c
†
i,σ
1
≤
i
≤
N
b
ˆ
̃
c
†
i
−
N
b
,σ
N
b
+ 1
≤
i
≤
2
N
b
,
(5)
for
N
b
bath orbitals.
The initial bath state in the super-Fock space,
|
ρ
B
⟩⟩
,
is given by a Slater determinant
|
ρ
B
⟩⟩
=
Y
i,σ
f
+
(
E
i,σ
,β
)ˆ
c
†
i,σ
+
f
−
(
E
i,σ
,β
)
ˆ
̃
c
†
i,σ
|
0
⟩⟩
(6)
where
|
0
⟩⟩
is a vacuum state in the super-Fock space,
f
+
(
E,β
) = (1 +
e
βE
)
−
1
, and
f
−
(
E,β
) = 1
−
f
+
(
E,β
).
3
We are often interested in the reduced density operator of
the impurity, ˆ
ρ
S
= Tr
B
( ˆ
ρ
). In the super-Fock space, the
trace over the bath is equivalent to taking the overlap
with a trace vector that can be expressed as a Slater
determinant,
|
Tr
B
⟩⟩
=
Y
i,σ
ˆ
c
†
i,σ
+
ˆ
̃
c
†
i,σ
|
0
⟩⟩
,
(7)
and the amplitude of a configuration
s
in
ρ
S
can be ex-
pressed as,
⟨⟨
s
|
ρ
S
⟩⟩
= (
⟨⟨
s
|⊗⟨⟨
Tr
B
|
)
|
ρ
⟩⟩
.
(8)
The discretized time evolution of the density operator
with timestep ∆
t
can be expressed via a second-order
Trotter decomposition
|
ρ
S
(
t
N
)
⟩⟩
= Tr
B
e
−
i
2
ˆ
L
S
∆
t
e
−
i
ˆ
L
SB
∆
t
e
−
i
2
ˆ
L
S
∆
t
N
t
|
ρ
(0)
⟩⟩
= Tr
B
h
ˆ
U
S
ˆ
U
SB
ˆ
U
S
···
ˆ
U
SB
ˆ
U
S
|
ρ
(0)
⟩⟩
i
(9)
where
t
N
t
=
N
t
∆
t
and
ˆ
U
S
(
ˆ
U
SB
) is the time evolution
operator for
ˆ
L
S
(
ˆ
L
SB
and
ˆ
L
B
). The tensor network di-
agram for this time evolution is shown in Fig. 1, where
we see that the
ˆ
U
S
tensor is applied only within the im-
purity
S
, whereas the tensor for
ˆ
U
SB
is applied to both
S
and
B
.
The influence functional (IF) tensor is defined as the
tensor arising from the contraction of all the bath degrees
of freedom in
ˆ
U
SB
, the initial bath density operator, and
the trace vector (Fig. 1). This IF tensor is indexed by
the state of the impurity orbitals before and after each
ˆ
U
SB
, for
N
t
time-steps (denoted
s
i
m
and
s
f
m
, respectively,
for the
m
-th time-step) and so is indexed by the configu-
ration of 2
N
t
impurity orbitals. In addition, we mention
that the IF tensors for different spins are constructed
separately thanks to the absence of spin-mixing terms in
ˆ
L
SB
, and hence, the spin indices are omitted for brevity.
Due to the one-dimensional structure of the tempo-
ral axis, the IF can be rewritten in terms of
N
t
tensors
that are directly obtained from
ˆ
U
SB
by inserting the
bath configurations between each time-step. Denoting
the impurity configurations,
s
= (
s
i
1
,s
f
1
,
···
,s
i
N
t
,s
f
N
t
),
and the corresponding IF tensor element as
I
(
s
), the ma-
trix product state (MPS) representation of the IF can be
written as,
I
(
s
) =
l
T
·
A
s
i
N
t
,s
f
N
t
N
t
···
A
s
i
1
,s
f
1
1
·
r,
(10)
where the matrix elements for
A
m
are given by
A
s
i
m
,s
f
m
m
b
m
,b
m
−
1
=
⟨⟨
s
f
m
,b
m
|
ˆ
U
SB
|
s
i
m
,b
m
−
1
⟩⟩
,
(11)
and
b
m
is the bath configuration after
m
applications of
ˆ
U
SB
,
r
b
0
=
⟨⟨
b
0
|
ρ
B
(0)
⟩⟩
, and
l
b
N
t
=
⟨⟨
Tr
B
|
b
N
t
⟩⟩
.
In the above MPS representation, the bond dimension
is given by the dimension of the super-Fock bath space,
ˆ
ρ
S
ˆ
ρ
B
Tr
B
IF
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
FIG. 1. A schematic diagram for the real-time evolution
of the Anderson Impurity model in Liouville space after
N
t
= 5 time-steps. The initial density operator is given by
ˆ
ρ
(0) = ˆ
ρ
S
⊗
ˆ
ρ
B
, described by a vectorized state in a Liouville
space. The time evolution of the density operator is described
by time evolution operators following a second-order Trotter
decomposition, alternating between
ˆ
U
S
(squares) and
ˆ
U
SB
(rectangles). After the time evolution, the bath degrees of
freedom are traced out, which is equivalent to applying the
trace tensor, Tr
B
, to the bath. The influence functional (IF)
tensor corresponds to the tensor after the contraction of the
bath modes in all
ˆ
U
SB
, ˆ
ρ
B
, and Tr
B
tensors.
which, in many cases, is too large to deal with directly.
Previous studies [35, 36] have made use of the fermionic
Gaussian properties of the IF (arising from the linear cou-
pling and quadratic bath) to find a compressed form of
the MPS representation [42, 43]. In this work, we will use
a slightly different language to formulate the MPS com-
pression in terms of finding Schmidt vectors in the bath.
This algorithm (described in Sec. II E) is closely related
to that in Ref. [42] and improves on the computational
scaling in Ref. [43].
B. From Influence Functionals to State-Vector
Propagation
In this work, we will often switch between two equiv-
alent pictures: dynamics encoded by a compressed IF-
MPS, and a state-vector propagation corresponding to a
quantum embedding. To understand this mapping, we
first describe the conventional MPS compression scheme
(i.e. without using any Gaussian properties of the bath)
and explain how it can be interpreted as a type of pro-
jected bath dynamics. (We recall that the IF-MPS can,
in general, be expressed through Eq. 10 and 11, for ar-
bitrary system-bath quantum dynamics, i.e. even for an
interacting bath).
4
b
1
b
2
(
s
i
1
, s
f
1
)
(
s
i
2
, s
f
2
)
(
s
i
3
, s
f
3
)
b
3
(
s
i
4
, s
f
4
)
b
4
(
s
i
5
, s
f
5
)
A
1
A
2
A
3
A
4
A
5
B
1
A
3
A
4
A
5
G
1
A
1
G
−
1
1
A
2
G
2
L
2
G
−
1
2
G
3
L
3
G
−
1
3
G
4
G
−
1
4
L
4
L
5
̃
b
2
̃
b
3
̃
A
1
A
3
A
4
A
5
G
1
A
1
G
−
1
1
A
2
G
2
̃
A
2
G
−
1
2
G
3
̃
A
3
G
−
1
3
G
4
G
−
1
4
̃
A
4
̃
A
5
̃
b
1
̃
b
4
FIG. 2. General matrix product state (MPS) compression scheme for an influence functional MPS (IF-MPS). The bond
dimension of the initial uncompressed MPS is given by the dimension of the bath Liouville space, indexed by bath configurations,
b
m
(left). The MPS needs to be converted into canonical form for an optimal truncation. The gauge matrices,
G
m
, are inserted
in the bath Liouville space to convert the MPS into the canonical form (middle). Gauge matrices transform the original matrix
elements of the MPS into left-normalized matrices, denoted as
L
m
=
G
m
A
m
G
−
1
m
−
1
. Afterward, projectors are inserted into
the bath subspace that contains the bath states with the largest singular values (right). The projected bath configurations are
denoted as
̃
b
m
.
In conventional MPS compression, the MPS is first
transformed to a canonical form to enable an optimal
truncation of the bond dimension. The canonical form
is defined using the gauge degrees of freedom in the
MPS [44],
A
s
i
m
,s
f
m
m
→
L
s
i
m
,s
f
m
m
=
G
m
A
s
i
m
,s
f
m
m
G
−
1
m
−
1
,
(12)
X
s
i
m
,s
f
m
L
s
i
m
,s
f
m
†
m
L
s
i
m
,s
f
m
m
=
I,
(13)
where the matrices,
L
s
i
m
,s
f
m
m
, satisfying Eq. 13, are called
left-normalized matrices. The gauge matrices,
G
m
, are
inserted in the bond space, or the Liouville bath space
of the IF-MPS (Fig. 2). The MPS that is composed of
the left-normalized matrices is called left-canonical. The
left-canonical MPS is then
I
(
s
) =
L
s
i
N
t
,s
f
N
t
N
t
···
L
s
i
2
,s
f
2
2
B
s
i
1
,s
f
1
1
,
(14)
where
B
s
i
1
,s
f
1
1
=
G
1
A
s
i
1
,s
f
1
1
. After choosing this gauge,
the MPS is compressed by iteratively applying truncated
singular value decompositions (SVD) from right to left.
We can group together the effect of the gauging and
compression together with the system-bath evolution to
define new matrices of the IF-MPS,
̃
A
m
,
̃
A
s
i
m
,s
f
m
m
̃
b
m
,
̃
b
m
−
1
=
⟨⟨
s
f
m
,
̃
b
m
|P
m
ˆ
G
m
ˆ
U
SB
ˆ
G
−
1
m
−
1
P
m
−
1
|
s
i
m
,
̃
b
m
−
1
⟩⟩
,
(15)
where
P
m
denotes the projectors onto the bath states as-
sociated with the largest singular values (or equivalently,
the largest eigenvalues of the bath density matrix) after
the gauge transformation (Fig. 2), and
̃
b
m
denotes the
projected bath configurations.
̃
A
m
can be viewed as defining a (non-unitary) evolution
in the Liouville space of the system and the bath. Alter-
natively, the MPS gauging and compression procedure
can be seen as a pure projected bath dynamics in the Li-
ouville space,
ˆ
G
−
1
m
P
m
ˆ
G
m
, inserted between the system-
bath evolution
ˆ
U
SB
. The coarse-graining (i.e. projec-
tion) of the bath degrees of freedom is referred to as
a quantum embedding, and consequently, the truncated
IF-MPS dynamics is a quantum embedding scheme with
a dynamically evolved bath, similar to (real-time) den-
5
sity matrix embedding theory (DMET) [37, 45–47], but
with the important difference that the dynamics of the
projected bath is non-unitary, as
ˆ
G
m
is non-unitary.
As the time-step ∆
t
→
0, we are led to the continuous-
time limit of the IF-MPS. We consider the subtleties of
continuous-time construction in Sec. III A. However, the
above shows that we can also view time evolution as be-
ing performed on an embedded state-vector (“wavefunc-
tion”) in the projected Liouville space (embedding space)
{|
s
i
m
̃
b
m
⟩⟩}
. We can then formulate the continuous-time
dynamics in terms of the equations of motion for the em-
bedded wavefunction and the bath projectors (or equiv-
alently, the bath density matrix). The dynamics of such
a bath density matrix has previously been considered in
the MPS language in Ref. [27] for a 1D spin chain, which
derived the dissipative contribution of the system-bath
coupling to the density matrix dynamics.
We will be interested in the above constructions for the
case of a fermionic Gaussian bath where we can replace
the discussion of many-body states and density matrices
with orbitals and 1-particle reduced density matrices (1-
RDM). We now turn to the formulation of the IF-MPS
operations in terms of these quantities.
C. Schmidt Decomposition of Slater determinants
In this section, we review the Schmidt decomposition
of Slater determinants [46, 48]. A Slater determinant is
given by
|
ψ
⟩
=
N
occ
Y
p
=1
ˆ
c
†
p
|
0
⟩
,
ˆ
c
†
p
=
X
i
C
ip
ˆ
a
†
i
(16)
where ˆ
c
†
p
is a creation operator of the
N
occ
occupied or-
bitals, ˆ
a
†
i
is a creation operator of orthonormal orbitals
in the basis of
n
sites, and
C
ip
is the orbital coefficient
matrix. Given a bipartite Hilbert space,
H
=
H
A
⊗H
B
,
where the first
n
A
orbitals belong to subsystem A and the
other
n
B
=
n
−
n
A
orbitals belong to subsystem B, the
Schmidt decomposition of the Slater determinant can be
obtained by diagonalizing the one-particle reduced den-
sity matrices (1-RDM), Γ
ij
=
⟨
ψ
|
ˆ
a
†
j
ˆ
a
i
|
ψ
⟩
, of the subsys-
tems. Assuming
n
A
< n
B
and
N
occ
> n
A
, the Schmidt
decomposition can be written as,
|
ψ
⟩
=
n
A
Y
k
=1
√
ν
k
ˆ
c
†
A,k
+
√
1
−
ν
k
ˆ
c
†
B,k
N
occ
Y
l
=
n
A
+1
ˆ
c
†
B,l
|
0
⟩
(17)
where
ν
k
(1
−
ν
k
) denote the eigenvalues of the 1-RDM of
A
(
B
) with values between 0 and 1, ˆ
c
†
A,k
(ˆ
c
†
B,k
) create the
corresponding eigenmodes, and ˆ
c
†
B,l
create eigenmodes
with eigenvalue 1 of the 1-RDM of
B
.
Based on the above, the orbitals in
B
can be classi-
fied into three different categories: (1)
entangled
orbitals,
ˆ
c
†
B,k
, 1
≤
k
≤
n
A
, which are entangled with
A
, (2)
core
orbitals, ˆ
c
†
B,l
,
n
A
+ 1
≤
l
≤
N
occ
, which are fully occu-
pied in
B
and not entangled with
A
, (3)
virtual
orbitals,
which are unoccupied and so do not appear in Eq. 17 and
also are not entangled with
A
.
Note that the Slater determinant with all
ν
k
=
1
2
cor-
responds to a maximally entangled fermionic state,
|
φ
⟩
,
where the reduced density operator of subsystem
A
is
proportional to the identity. It is possible to write the
1-RDM of
A
from
|
ψ
⟩
as that of a maximally entangled
fermionic state
|
φ
⟩
after a gauge transformation within
the subsystem
A
, i.e.
ˆ
G
A
|
φ
⟩
, where
ˆ
G
A
= exp
X
k
log
g
k
ˆ
c
†
A,k
ˆ
c
A,k
!
, g
k
=
r
ν
k
1
−
ν
k
,
(18)
up to a normalization constant factor, assuming 0
< ν
k
<
1 for all
k
.
ˆ
G
A
satisfies the following
ˆ
G
A
ˆ
c
†
A,k
ˆ
G
−
1
A
=
g
k
ˆ
c
†
A,k
,
ˆ
G
A
ˆ
c
A,k
ˆ
G
−
1
A
=
g
−
1
k
ˆ
c
A,k
(19)
This gauge transformation is related to the gauge trans-
formation introduced in Sec. II B because the state
|
φ
⟩
is ‘left-normalized’ with respect to the subsystem
A
. We
will use this gauge transformation to convert Slater de-
terminants into left-normalized forms in the next section.
The Schmidt decomposition can be truncated by treat-
ing the entangled orbitals with
ν
k
≈
1 (
ν
k
≈
0) as
core (virtual) orbitals, retaining the orbitals with larger
p
ν
k
(1
−
ν
k
). In other words, a projection operator on
the entangled orbital space that keeps only
ν
k
close to
1
2
can be applied to truncate the Schmidt decomposition.
This truncation scheme is sometimes called a “mode”
truncation [42, 49, 50] and has been utilized in the con-
text of tensor network truncations of fermionic Gaussian
states.
After a truncation to
n
ent
entangled orbitals, the Slater
determinant can be written as,
|
ψ
⟩
=
n
ent
Y
k
=1
√
ν
k
ˆ
c
†
A,k
+
√
1
−
ν
k
ˆ
c
†
B,k
×
n
ent
+
n
A,c
Y
l
=
n
ent
+1
ˆ
c
†
A,l
n
ent
+
n
B,c
Y
l
=
n
ent
+1
ˆ
c
†
B,l
|
0
⟩
,
(20)
where
n
A,c
(
n
B,c
) is the number of core orbitals in
A
(
B
).
D. Influence Functional tensors as Slater
determinants
We now describe how the bipartitions of the exact in-
fluence functional tensor can be expressed as Slater deter-
minants, which are obtained by propagating a finite num-
ber of steps forward in time from
|
ρ
B
(0)
⟩⟩
or backward in
time from
⟨⟨
Tr
B
|
. Thanks to the noninteracting nature
of the bath, both the IF and its partitions correspond to
fermionic Gaussian states, specifically, Bardeen-Cooper-
Schrieffer (BCS) states [34–36]. Despite this fact, we will
6
ρ
B
(0)
ρ
B
(0)
ρ
B
(0)
(a)
(b)
Tr
B
(c)
A
S
B
B
B
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
ˆ
U
SB
|
Ψ
R,
1
〉〉
|
Ψ
R,
3
〉〉
|
Ψ
L,N
t
−
3
〉〉
φ
s
i
1
s
f
1
b
1
FIG. 3. (a) A right IF state after one time-step,
|
Ψ
R,
1
⟩⟩
. It
represents a partial contraction of the bath degrees of freedom
between the time evolution operator,
ˆ
U
SB
, and the initial
bath density operator,
ρ
B
(0). It can be expressed as a Slater
determinant by introducing a maximally entangled fermionic
state,
φ
, which has impurity and auxiliary fermionic orbitals
denoted by
S
and
A
, respectively. Coupling the maximally
entangled state to the time evolution operator on the input
impurity orbitals yields the Slater determinant state. (b) A
right IF state after the
m
-th time-step,
|
Ψ
R,m
⟩⟩
, (
m
= 3 in
the figure) can also be represented as a Slater determinant
by introducing maximally entangled fermionic states at each
time-step. (c) A left IF state,
|
Ψ
L,m
⟩⟩
, (
m
=
N
t
−
3 in the
figure) can be constructed by contracting the time evolution
operator from the top downwards, and can also be repre-
sented as a Slater determinant by inserting maximally en-
tangled fermionic states starting from the top.
prefer to work with Slater determinants and convert the
BCS states to Slater determinants through a particle-
hole transformation. This is because when building the
matrix elements of the IF-MPS in Eq. 15, the basis of
Slater determinants will allow us to use number sym-
metry, which significantly reduces the prefactors in the
numerical computations.
We define the right and left bipartitions of IF (right
and left IF for short)
Q
m
ˆ
U
SB
m
· |
ρ
B
(0)
⟩⟩
and
⟨⟨
Tr
B
| ·
Q
m
ˆ
U
SB
m
, respectively, where the products follow the or-
der in the MPS representation and here,
·
indicates the
partial contraction of the intermediate bath configura-
tions. We start with the first right IF,
ˆ
U
SB
· |
ρ
B
(0)
⟩⟩
for simplicity. Its tensor elements are determined from
⟨⟨
s
f
1
,b
1
|
ˆ
U
SB
(
|
s
i
1
⟩⟩⊗|
ρ
B
(0)
⟩⟩
), and its external degrees of
freedom are
s
i
1
,
s
f
1
, and
b
1
(see Fig. 3(a)). Its occupation
numbers satisfy the relationship,
n
(
s
f
1
) +
n
(
b
1
)
−
n
(
s
i
1
) =
n
(
ρ
B
(0))
,
(21)
which is not consistent with a number-conserving state.
After applying the particle-hole transformation to the
input configuration,
s
i
1
→
̄
s
i
1
, (for example, 00
→
11,
01
→
10, 10
→
01, and 11
→
00), we obtain,
n
(
s
f
1
) +
n
(
b
1
) +
n
( ̄
s
i
1
) =
n
(
ρ
B
(0)) + 2
,
(22)
where
n
( ̄
s
i
1
) = 2
−
n
(
s
i
1
). Therefore, after the transforma-
tion, the right partition IF is a Slater determinant with
n
(
ρ
B
(0)) + 2 occupied orbitals.
Formally, the particle-hole transformations can be ex-
pressed in terms of fermionic tensor network contractions
by inserting maximally entangled fermionic states,
|
φ
⟩⟩
=
Y
s
1
√
2
( ˆ
α
†
s
+
ˆ
d
†
s
)
|
0
⟩⟩
,
(23)
where the index
s
denotes the input impurity orbitals
from the super-Fock Liouville space and ˆ
α
†
s
indicates a
creation operator of auxiliary fermionic orbitals. Note
that the occupation numbers of the auxiliary fermionic
orbitals have particle-hole transformed occupation num-
bers compared to the original input impurity orbitals.
The right IF state representation of
ˆ
U
SB
·|
ρ
B
(0)
⟩⟩
, which
we will denote
|
Ψ
R
⟩⟩≡|
Ψ
R,
1
⟩⟩
, is then written as follows
(Fig. 3a),
|
Ψ
R,
1
⟩⟩
=
ˆ
I
A
⊗
ˆ
U
SB
|
φ
⟩⟩⊗|
ρ
B
(0)
⟩⟩
,
(24)
where
ˆ
I
A
is the identity operator on the auxiliary
fermionic orbitals.
Because the initial state,
|
φ
⟩⟩ ⊗
|
ρ
B
(0)
⟩⟩
is given by a Slater determinant and
ˆ
I
A
⊗
ˆ
U
SB
is number-conserving, we see again that the right IF
state is also a Slater determinant. The right IF state for
Q
m
ˆ
U
SB
m
·|
ρ
B
(0)
⟩⟩
can be similarly expressed by insert-
ing maximally entangled fermionic states at each time-
step and coupling one of the modes to
ˆ
U
SB
, as shown
schematically in Fig. 3b.
To construct the compression of the IF-MPS later, we
require the bath 1-RDMs. Given the right IF state at
the
m
-th time-step,
|
Ψ
R,m
⟩⟩
, we define the bath 1-RDM
of the right IF state (Fig. 4a) as,
Γ
R,m
ij
=
⟨⟨
Ψ
R,m
|
ˆ
a
†
j
ˆ
a
i
|
Ψ
R,m
⟩⟩
⟨⟨
Ψ
R,m
|
Ψ
R,m
⟩⟩
,
(25)
where the indices
i
and
j
refer to the super-Fock bath
orbitals. The 1-RDM at the next time-step can be com-
puted from the evolution of
|
φ
⟩⟩⊗|
Ψ
R,m
⟩⟩
under
ˆ
I
A
⊗
ˆ
U
SB
.
7
The bath 1-RDM after the evolution can be written as,
Γ
R,m
+1
ij
=
U
SB
Γ
R,m
U
SB
†
ij
+
U
ASB
Γ
φ
U
ASB
†
ij
,
(26)
where Γ
φ
denotes the 1-RDM of
|
φ
⟩⟩
and
U
ASB
denotes
a single particle-basis representation of
ˆ
I
A
⊗
ˆ
U
SB
. The
diagrammatic representation of the evolution of the bath
1-RDM is drawn in Fig. 4b.
We can construct the left IF state in an analogous way,
but where the state propagates in the inverse (negative)
time direction using
ˆ
U
SB
†
(Fig. 4c). Denoting the left IF
state and its bath 1-RDM,
|
Ψ
L
⟩⟩
and Γ
L
, respectively, the
bath 1-RDMs at successive time-steps (from top down-
wards) are related by
Γ
L,m
−
1
ij
=
U
SB
†
Γ
L,m
U
SB
ij
+
U
ASB
†
Γ
φ
U
ASB
ij
,
(27)
Note that the initial state is given by the trace vector,
|
Ψ
L,N
t
⟩⟩
=
|
Tr
B
⟩⟩
, which is also given by a Slater deter-
minant as in Eq. 7. With both the left and right IF state,
the IF-MPS can be written as
⟨⟨
Ψ
L,m
|·|
Ψ
R,m
⟩⟩
for any
m
.
E. Influence Functional Matrix Product State
Compression
We now revisit the MPS compression of the IF-
MPS, described in Sec. II B, but now utilizing the non-
interacting nature of the bath, which allows all the steps
to be expressed at the level of orbitals and single-particle
quantities.
In Sec II B, the optimal compression of the MPS re-
quired the IF-MPS to be in a canonical form. This was
achieved by transforming matrices into left-normalized
matrices by inserting gauge matrices into the MPS. As
discussed in Sec. II D the partitions of the IF-MPS for
a non-interacting bath are Slater determinants, and the
gauge transformation (Eq. 18) to convert Slater determi-
nants to maximally entangled fermionic pairs, which are
left-normalized, was introduced in Sec. II C.
We therefore have all the ingredients to convert the IF-
MPS to canonical form. We start with the left IF state,
|
Ψ
L
⟩⟩
, and determine the gauge transformation,
ˆ
G
, from
the eigenvalues and eigenvectors of the bath 1-RDM, Γ
L
ij
.
We denote the eigenvalues of Γ
L
ij
as
ν
k
, where
k
indexes
the eigenvalues, and the rotation matrix as
R
L
ik
, whose
columns are the eigenvectors of Γ
L
in a single-particle
basis. With this gauge transformation, we can represent
the IF-MPS as,
⟨⟨
Ψ
L
|·|
Ψ
R
⟩⟩
=
⟨⟨
φ
|·
ˆ
G
·|
Ψ
R
⟩⟩
. The gauge
matrix in the single-particle basis,
G
, can be written as,
G
ki
=
g
k
R
L
∗
ik
,
(28)
where
g
k
=
p
ν
k
/
(1
−
ν
k
) (Eq. 18). Note that
g
k
diverges
when
ν
k
→
1, so in practice, we regularize
ν
k
with a small
threshold
ε
so that
ν
k
=
ε
when
ν
k
< ε
and
ν
k
= 1
−
ε
when 1
−
ν
k
< ε
. A similar regularization scheme for
the 1-RDM has been used in multiconfiguration time-
dependent Hartree theory [51, 52] and real-time density
matrix embedding theory [37, 47].
The gauge matrices need to be absorbed into the right
IF state to further canonicalize the IF-MPS. We call the
state,
|
Ψ
G
⟩⟩
=
ˆ
G
·|
Ψ
R
⟩⟩
, a gauge-transformed right IF
state and its bath 1-RDM, Γ
G
. The gauge transformation
is a non-unitary transformation, so we always normalize
the state when computing Γ
G
.
Γ
G
ij
=
⟨⟨
Ψ
G
|
ˆ
a
†
j
ˆ
a
i
|
Ψ
G
⟩⟩
⟨⟨
Ψ
G
|
Ψ
G
⟩⟩
=
h
G
√
Γ
R
(
√
Γ
R
G
†
G
√
Γ
R
+
I
−
Γ
R
)
−
1
√
Γ
R
G
†
i
ij
(29)
where Γ
R
is the bath 1-RDM of the
|
Ψ
R
⟩⟩
. A detailed
derivation of this expression is included in the Supple-
mental Material.
Subsequently, we decompose the gauge-transformed
right IF state into a right-normalized maximally entan-
gled fermionic state and another gauge transformation,
and this gauge transformation contains the singular val-
ues of the IF-MPS at this bipartition. Taking the largest
singular values in the truncated SVD in the conventional
MPS compression is equivalent to taking the eigenvectors
of Γ
G
with eigenvalues closest to
1
2
(or large
ν
G
k
(1
−
ν
G
k
))
as done in the mode truncation approach introduced in
Sec. II C. We will call the bath orbitals from the selected
eigenvectors of Γ
G
the effective bath orbitals.
Hence, the procedure for IF-MPS compression at the
orbital level is as follows: (1) Diagonalize the bath 1-
RDM of the gauge-transformed right IF state, Γ
G
, and
obtain its eigenvalues,
ν
G
k
, and rotation matrix,
R
G
. (2)
Take the eigenvectors with the
N
eff
largest
ν
G
k
(1
−
ν
G
k
),
where
N
eff
is a number of the effective bath orbitals we
select, defining the truncated rotation matrix
R
eff
. (3)
Construct the projection operator,
P
, from
R
eff
, i.e. the
eigenvectors corresponding to the effective bath orbitals.
The configurations from the other bath orbitals in the
projected subspace are fixed to be either fully occupied
(core orbitals,
ν
G
k
≈
1) or unoccupied (virtual orbitals,
ν
G
k
≈
0).
For the Liouville time evolution of the Anderson impu-
rity model with the initial Gaussian thermal bath, it is
possible to prove that if
ν
G
k
is an eigenvalue of Γ
G
, so is
1
−
ν
G
k
(Supplemental Material). In this case, we obtain
the same number of core and virtual orbitals and choose
the effective bath orbitals symmetrically by occupancy.
The tensor elements of the IF-MPS can then be com-
puted from Eq. 15, after applying the gauge transforma-
tion and the projectors,
P
m
, for each
m
-th time-step,
⟨⟨
s
f
m
,
̃
b
m
|P
m
ˆ
G
m
ˆ
U
SB
ˆ
G
−
1
m
−
1
P
m
−
1
|
s
i
m
,
̃
b
m
−
1
⟩⟩
,
(30)
with the configurations of the effective bath orbitals,
̃
b
m
.
This tensor element defines an effective time evolution
operator for the embedded wavefunction defined on the
impurity and effective bath orbitals. The diagrammatic
8
=
Γ
R
=
Γ
R,m
+1
Γ
R,m
Γ
R
Γ
L
Γ
R
φ
ˆ
G
ˆ
G
Γ
G
φ
=
=
ˆ
U
SB
→
ˆ
U
SB
ˆ
G
m
+1
ˆ
G
−
1
m
(a)
(b)
(c)
(d)
|
1
〉
⊗
n
=
|
0
〉
⊗
n
=
FIG. 4.
(a)
The bath 1-RDM of the right IF state, Γ
R
, is obtained by tracing out the impurity orbitals, i.e. by contracting
the state with its complex conjugate, within the impurity space.
(b)
The bath 1-RDM of the right IF state after the
m
-th
time-step, Γ
R,m
can be updated to Γ
B,m
+1
by applying the time evolution operator and tracing out the impurity orbitals.
(c)
A gauge transformation,
ˆ
G
, is extracted from the bath 1-RDM of the left IF state, Γ
L
, by converting it to a maximally
entangled fermionic state,
φ
. The gauge-transformed right IF state is obtained by applying the gauge transformation,
ˆ
G
, to
the bath, and its bath 1-RDM is transformed to Γ
G
.
(d)
The impurity-bath time evolution operator,
ˆ
U
SB
, is approximated
by an effective time evolution operator by projecting the bath to the effective bath orbitals after the gauge transformation,
P
m
+1
ˆ
G
m
+1
ˆ
U
SB
ˆ
G
−
1
m
P
m
. The filled and unfilled rectangles represent the core and virtual orbital spaces and these spaces are
projected into fully occupied and unoccupied states, respectively.
representation for the effective time evolution operator is
illustrated in Fig. 4d.
The above tensor elements can be efficiently computed
from determinant formulae. Since the configurations in
the core and virtual orbitals are fixed, we can compute
the determinant of block matrices, keeping the core and
virtual orbital block matrices fixed. Therefore, the com-
putational complexity to compute determinants for all
configurations is
O
(
N
3
b
+ 2
2
N
eff
N
3
eff
), where the first term
corresponds to the computation of the determinant and
inverse of the core and virtual orbital block matrices and
the second term corresponds to the computation of the
determinant of the effective bath orbital block matrices.
The cost for computing the full set of MPS tensor ele-
ments at time-step
N
t
is
O
(
N
3
b
N
t
+ 2
2
N
eff
N
3
eff
N
t
), which
is
linear
in the number of time-steps
N
t
.
III. CONTINUOUS-TIME FORMULATION OF
IF-MPS
A. Boundary Influence functional tensor
Defining the continuous-time limit of the influence
functional tensor network is in principle one way to elim-
inate the time-step error from the standard second-order
Trotter decomposition, and in numerical applications al-
lows for the introduction of a wide variety of higher-order
differential integrators. However, as shown in [33], the
IF-MPS shows a nonphysical entanglement entropy scal-
ing in the limit of ∆
t
→
0, as the entanglement entropy
always scales to zero. This suggests that the continuous-
time limit requires a more careful treatment.
In particular, the formalism of continuous matrix prod-
uct states (cMPS) [53, 54] describes a quantum wave-
function of continuous variables that (in general) sup-
ports a finite entanglement entropy as the discretiza-
tion approaches the continuum limit. In this section, we
show that the usual IF-MPS does not support a standard
cMPS representation in the continuous-time limit, and
instead a closely related object, the
boundary
influence
functional MPS should be used. The boundary influence
functional MPS is implicitly used in transverse contrac-
tion [26, 55] (e.g., Ref. [27] states that this becomes a
continuous MPS in the continuum limit, without pro-
viding an explicit construction) and has previously been
used in influence functional calculations with interacting
baths [31].
Consider the continuous-time limit ∆
t
→
0, where
ˆ
U
SB
can be expressed as,
ˆ
U
SB
=
ˆ
I
−
i
(
ˆ
L
SB
+
ˆ
L
B
)∆
t,
(31)
and the corresponding matrix elements of the IF-MPS,
A
s
i
,s
f
, are,
A
0
,
0
=
I
−
i
ˆ
L
B
∆
t
A
0
,
1
=
A
1
,
0
=
−
i
ˆ
L
SB
∆
t
A
1
,
1
=
I
−
i
ˆ
L
B
∆
t.
(32)