View
Online
Export
Citation
CrossMark
RESEARCH ARTICLE
|
DECEMBER 18 2023
B 2
: A comprehensive open source framework to develop
and apply state-of-the-art DMRG algorithms in electronic
structure and beyond
Huanchen Zhai
;
Henrik R. Larsson
;
Seunghoon Lee
;
Zhi-Hao Cui
;
Tianyu Zhu
;
Chong Sun
;
Linqing Peng
;
Ruojing Peng
;
Ke Liao
;
Johannes Tölle
;
Junjie Y
ang
;
Shuoxue Li
;
Garnet Kin-Lic Chan
J. Chem. Phys.
159, 234801 (2023)
https://doi.org/10.1063/5.0180424
19 December 2023 20:42:18
The Journal
of Chemical Physics
ARTICLE
pubs.aip.org/aip/jcp
B
LOCK
2: A comprehensive open source framework
to develop and apply state-of-the-art DMRG
algorithms in electronic structure and beyond
Cite as: J. Chem. Phys.
159
, 234801 (2023); doi: 10.1063/5.0180424
Submitted: 9 October 2023
•
Accepted: 16 November 2023
•
Published Online: 18 December 2023
Huanchen Zhai,
a)
Henrik R. Larsson,
b)
Seunghoon Lee,
c)
Zhi-Hao Cui,
d)
Tianyu Zhu,
e)
Chong Sun,
f)
Linqing Peng,
Ruojing Peng,
Ke Liao,
g)
Johannes Tölle,
Junjie Yang,
Shuoxue Li,
and Garnet Kin-Lic Chan
h)
AFFILIATIONS
Division of Chemistry and Chemical Engineering, California Institute of Technology, Pasadena, California 91125, USA
Note:
This paper is part of the JCP Special Topic on Modular and Interoperable Software for Chemical Physics.
a)
Electronic mail:
hczhai.ok@gmail.com
b)
Present address:
Department of Chemistry and Biochemistry, University of California, Merced, California 95343, USA.
c)
Present address:
Department of Chemistry, Seoul National University, Seoul 151-747, South Korea.
d)
Present address:
Department of Chemistry, Columbia University, New York, New York 10027, USA.
e)
Present address:
Department of Chemistry, Yale University, New Haven, Connecticut 06520, USA.
f)
Present address:
Department of Chemistry, Rice University, Houston, Texas 77005, USA.
g)
Present address:
Department of Physics and Arnold Sommerfeld Center for Theoretical Physics, Ludwig-Maximilians-
Universität München, D-80333 Munich, Germany.
h)
Author to whom correspondence should be addressed:
gkc1000@gmail.com
ABSTRACT
BLOCK
2 is an open source framework to implement and perform density matrix renormalization group and matrix product state algo-
rithms. Out-of-the-box it supports the eigenstate, time-dependent, response, and finite-temperature algorithms. In addition, it carries special
optimizations for
ab initio
electronic structure Hamiltonians and implements many quantum chemistry extensions to the density matrix
renormalization group, such as dynamical correlation theories. The code is designed with an emphasis on flexibility, extensibility, and effi-
ciency and to support integration with external numerical packages. Here, we explain the design principles and currently supported features
and present numerical examples in a range of applications.
Published under an exclusive license by AIP Publishing.
https://doi.org/10.1063/5.0180424
I. INTRODUCTION
The
density
matrix
renormalization
group
(DMRG)
algorithm
1–3
is now a widely used numerical approach to
obtain the low-energy eigenstates and other properties of a wide
range of systems of interest in electronic structure, such as strongly
correlated low-dimensional models,
4–8
large active space problems
in quantum chemistry,
9–18
electronic quantum dynamics,
19,20
and
even models to benchmark quantum computing protocols.
21,22
Over the last several decades, there have been many developments
and extensions of the original algorithm for these applications.
23–25
For example, dynamical electron correlation,
26,27
the effects
of finite temperature,
28,29
relativistic effects,
30
time-dependent
phenomena,
31
and frequency-dependent
32
response properties can
now all be treated within the DMRG framework.
In parallel with these theoretical advances, we have also wit-
nessed the success of many different DMRG implementations for
electronic structure,
11,13,33–44
which together have helped demon-
strate the power of DMRG in different settings. However, as more
extensions of the DMRG algorithm have appeared and as DMRG
techniques have grown to encompass more applications, so has the
complexity grown in many implementations.
The new open source
BLOCK
2 code
45
introduced here rep-
resents an effort to provide users with a comprehensive set of
state-of-the-art DMRG algorithms used today in electronic structure
and other applications while enabling access to the development and
J. Chem. Phys.
159
, 234801 (2023); doi: 10.1063/5.0180424
159
, 234801-1
Published under an exclusive license by AIP Publishing
19 December 2023 20:42:18
The Journal
of Chemical Physics
ARTICLE
pubs.aip.org/aip/jcp
extensions of these methods. To this end, the core design principles
in
BLOCK
2 aim to balance efficiency, flexibility, and simplicity. For
example, we have the following:
(i)
For efficiency and flexibility, we provide a special highly
optimized DMRG implementation for the standard quantum
chemistry Hamiltonian, but many of the optimizations also
accelerate computations with arbitrary sparse Hamiltonians,
thus enabling production level applications involving custom
Hamiltonians.
(ii)
For efficiency and simplicity, we provide both optimized and
unoptimized implementations of most features in the code,
so a developer can always start with the simplest version to
explore new ideas or compare different implementations to
understand our optimization techniques. The framework sup-
ports both the matrix product state (MPS)/matrix product
operator (MPO) and renormalized operator views of DMRG
algorithms, balancing the ease of expressing algorithms with
optimized performance.
(iii)
For flexibility and simplicity, while the algorithms in
BLOCK
2
support general Hamiltonian and operator expressions, e.g.,
arbitrary-order interactions, arbitrary-order density matrix
computations, and arbitrary-order excitations, the compu-
tations can be specified through a simple symbolic Python
syntax.
We note that to achieve the above design goals, many imple-
mentations in
BLOCK
2 use techniques that to our knowledge have
not been published in the literature. The following does not attempt
to provide a full description of the implementation, which after all
comprises
∼
200 000 lines of code. Instead, we introduce a sample of
the existing features of
BLOCK
2, some ideas behind the implementa-
tion, and how the capabilities are used in applications and algorithm
development.
II. PROGRAM FEATURES
A. Standard DMRG
BLOCK
2 provides a highly optimized and flexible implemen-
tation of the DMRG algorithm for solving various aspects of the
electronic structure problem using
ab initio
and model Hamiltoni-
ans. This section presents an overview of the most common features
of the code used in electronic structure calculations.
1. Ab initio DMRG and symmetries
The most common application of the
ab initio
DMRG algo-
rithm is to compute a many-body ground state wavefunction with
some desired symmetries.
2,35,41,42,46–48
One first requires the Hamil-
tonian expressed in a basis of orthogonal orbitals, which form the
sites in DMRG. Such orbitals, and the corresponding Hamiltonian
integrals, can be computed by many quantum chemistry packages,
for example, P
Y
SCF.
49,50
The Hamiltonian expressed in the chosen orbital basis typically
has a number of symmetries.
BLOCK
2 has built-in efficient support
for many different symmetries [including particle-number U(1),
34
total spin SU(2),
35
projected spin U(1),
34
K
-point Z
n
,
51
SO(4) in
the Hubbard model,
52,53
and the Abelian point group (D
2h
and
its subgroups
24
)] and is designed to be easily extensible to other
symmetries.
FIG. 1.
The scaling of DMRG wall time [using 24 central processing unit (CPU)
cores] per sweep with respect to the MPS bond dimension
M
for the (36o, 48e)
active space of a protonated iron sulfur dimer
[
Fe
2
S
(
CH
3
)(
SCH
3
)
4
]
3
−
. The
active space model is from Ref. 55. Spin-adapted uses particle-number and SU(2)
symmetry, while non-spin-adapted uses particle-number and projected-spin U(1)
symmetry.
Figure 1 illustrates the computational scaling with bond dimen-
sion of the
ab initio
DMRG algorithm implemented in
BLOCK
2 using
different spin symmetries. We note that the actual scaling of
ab initio
DMRG is lower than the formal
O
(
M
3
)
scaling (where
M
is the
MPS bond dimension), mainly due to the use of symmetries and
block-sparse matrix operations.
54
2. Sweep algorithms
Standard DMRG optimizations (sweeps) update the wave-
function in one of two ways, termed the two-site and one-site
algorithms.
2,56,57
The local Hilbert space in DMRG refers to the Fock
space of the site, which is of, e.g., dimension
d
=
4 for a spatial
orbital and
d
=
2 for a spin orbital. Assuming each DMRG site is
a spatial orbital, the two-site algorithm uses the product space of the
larger two-site local Hilbert space
d
=
16, and the left and right block
renormalized states to represent the DMRG renormalized wavefunc-
tion (which helps escape from local minima), while the one-site
algorithm substitutes the
d
=
4 local Hilbert space, lowering mem-
ory and computational costs. An effective “half-site” calculation with
d
=
2 can also be performed by using the one-site algorithm with
spin orbital sites. One can also fuse sites in the sweep [by contract-
ing adjacent tensors in the Hamiltonian matrix product operator
(MPO)] before performing DMRG, allowing
BLOCK
2 to perform
larger site sweep optimizations.
We note that the dimension of the local Hilbert space on
different sites need not be the same. For example, DMRG algorithms
with one or two large sites
58
use a single site to represent all the
virtual and/or core orbitals to efficiently treat dynamic correlation.
In such cases,
BLOCK
2 supports mixing the one-site and two-site
sweep algorithms to balance efficiency and accuracy.
58
3. Excited states
The ground-state DMRG algorithm can be extended to opti-
mize excited states, and such excited-state algorithms are available
in
BLOCK
2. The excited-state algorithms implemented in
BLOCK
2
take the following forms:
J. Chem. Phys.
159
, 234801 (2023); doi: 10.1063/5.0180424
159
, 234801-2
Published under an exclusive license by AIP Publishing
19 December 2023 20:42:18
The Journal
of Chemical Physics
ARTICLE
pubs.aip.org/aip/jcp
(i)
The state-averaged DMRG approach.
59
In this formalism,
there is a common set of rotation matrices for all states, but
the center site carries a unique tensor (renormalized wave-
function) for each state. This algorithm is normally a robust
and convenient choice for tens of roots, but for a given bond
dimension, the accuracy decreases as the number of roots is
increased.
(ii)
The projected orthogonalization approach.
40,60,61
To improve
the quality of the states obtained from the state-averaged
DMRG, one can further refine each excited state (as an inde-
pendent MPS) by projecting out lower-lying states when
optimizing the central site.
(iii)
The level-shift approach.
33,39
Here, we compute the ground
and excited states one-by-one with a modified level-shifted
Hamiltonian defined as
ˆ
H
′
k
=
ˆ
H
+
k
−
1
∑
i
=
1
w
i
∣
Ψ
i
⟩⟨
Ψ
i
∣
,
(1)
where
∣
Ψ
i
⟩
are the converged states below the targeted excited state
∣
Ψ
k
⟩
and
w
i
are the energy level shifts. This can be combined with
the state-averaged ansatz to determine batches of states at a time.
(iv)
The harmonic Davidson approach.
59
By changing the stan-
dard Davidson solver to the harmonic Davidson solver, one
can directly target interior excited states close to a given
energy without computing all the states below. In the current
implementation of
BLOCK
2, this approach is not as stable as
the above approaches for dense spectra.
4. Relativistic effects
Relativistic effects can be classified into scalar relativistic and
spin-dependent relativistic effects,
62
where the mean-field level
scalar relativistic effects can be described inexpensively.
63
B
LOCK
2
supports several strategies to include the scalar and spin-dependent
relativistic effects.
(i)
Relativistic DMRG.
64–66
We can directly solve for eigenstates
of four-component Dirac–Coulomb, Dirac–Coulomb–Gaunt,
or Dirac–Coulomb–Breit Hamiltonians.
67
Here, DMRG is
performed in complex number mode where both the MPS and
MPO are complex.
(ii)
One-step SOC-DMRG.
68
Within the two-component for-
malism of relativistic quantum chemistry, we can include
spin–orbit coupling (SOC) by adding a SOC term to the
scalar-relativistic Hamiltonian and then perform DMRG on
this two-component Hamiltonian with a complex MPO and
MPS. Alternatively, we can simplify the SOC term within
the spin–orbit mean-field (SOMF) approximation,
69
giving
rise to complex one-electron and real two-electron integrals.
One can then represent the Hamiltonian as the sum of two
MPOs where the expensive two-electron interaction is con-
tained only in the second MPO, which is completely real.
68
One then runs
BLOCK
2 in the hybrid complex and real arith-
metic DMRG mode, where the real and imaginary parts of the
effective wavefunction share the same real rotation matrices
in the MPS.
(iii)
Two-step SOC-DMRG.
30,70,71
When the SOC effect is weak,
we can first compute eigenstates of the scalar-relativistic
Hamiltonian. Then, SOC effects can be added by carrying out
the state-interaction of a small number of eigenstates.
FIG. 2.
The energy difference between representative excited states computed
from the one-step and two-step SOC-DMRG approaches for the (7o, 9e) active
space of the
[
DyCl
6
]
3
−
cluster. Adapted from Fig. 4(a) in Ref. 68.
Figure 2 shows a numerical comparison of the accuracy of the
1-step and 2-step SOC-DMRG treatments.
5. Non-standard Hamiltonians and operators
In
addition
to
computations
involving
the
standard
Hermitian quantum chemistry Hamiltonian with up to two-body
interaction terms,
BLOCK
2 supports the use of other types of
non-standard Hamiltonians and operators in DMRG algorithms
and/or MPO/MPS operations.
(i)
General Hamiltonians with high-order interactions. Using the
techniques in Sec. II C,
BLOCK
2 supports DMRG computa-
tions with arbitrary parity-preserving (i.e., even number of
creation and annihilation operator) Hamiltonians. For exam-
ple, DMRG calculations can be carried out for Hamiltoni-
ans with three-body or up to
n
-body terms. We note that
higher than two-body terms appear commonly in quantum
chemistry in effective Hamiltonians obtained from similar-
ity transformations
72
or canonical transformations
73–76
of the
ab initio
Hamiltonians.
(ii)
Normal-ordered Hamiltonians.
77
The DMRG algorithm can
work instead with the normal-ordered version of a Hamilto-
nian where the vacuum expectation value is separated out.
This can make the DMRG algorithm with reduced floating
point precision numerically more stable.
(iii)
Non-Hermitian Hamiltonians.
72,78–80
This allows DMRG to
be used, for example, with Hamiltonians that arise from
many-body similarity transformations, such as in coupled
cluster theory or the transcorrelated method.
(iv)
Anti-Hermitian operators.
BLOCK
2 has special support for
anti-Hermitian operators. Exponentials of anti-Hermitian
operators appear in various applications, most notably, when
performing orbital rotations of the MPS, or in canoni-
cal transformation theory.
74
The explicit support for anti-
Hermitian operators allows such unitaries to be expressed
as an exponential of a real MPO, thus avoiding complex
arithmetic.
(v)
General MPO and MPS operations with arbitrary operators.
BLOCK
2 supports MPO and MPS operations with arbitrary
J. Chem. Phys.
159
, 234801 (2023); doi: 10.1063/5.0180424
159
, 234801-3
Published under an exclusive license by AIP Publishing
19 December 2023 20:42:18
The Journal
of Chemical Physics
ARTICLE
pubs.aip.org/aip/jcp
operators, including non-number-conserving operators, as
described in Sec. II D. For example, to support the MPS
compression approach for strongly-correlated NEVTP2,
81
one needs to apply a partial Hamiltonian to the reference
wavefunction to get the perturber wavefunction,
∣
Ψ
′
i
⟩
=
(
∑
abc
v
iabc
a
†
a
a
†
b
a
c
) ∣
Ψ
0
⟩
.
(2)
In
BLOCK
2, this is done by building an MPO for the operator
part (as indicated by the parenthesis) in the above expression.
6. Initial guesses
BLOCK
2 provides several ways to generate initial MPS guesses
for various algorithms. It also provides an implementation of pertur-
bative noise,
56,57
which helps DMRG algorithms escape from poor
initial guesses.
(i)
A random MPS with the desired bond dimension can be
generated.
24
(ii)
One can supply an estimate of the occupancy information at
each site from a cheaper method [e.g., coupled cluster sin-
gles and doubles (CCSD)].
BLOCK
2 then uses this information
to construct the initial quantum number distribution in the
MPS (by assuming that the probability distribution of the
occupancy at each site is independent from each other).
(iii)
We can generate an initial guess from another (more approx-
imate) MPS using fitting. For example, we can use an MPS
parameterized and optimized in the singles and doubles
excitation space to initialize the MPS in the full Hilbert space.
(iv)
We can construct an initial guess MPS from a given lin-
ear combination of determinants or Configuration State
Functions (CSFs).
7. Orbital reordering algorithms
For a given bond dimension, the accuracy of DMRG is affected
by the orbital ordering.
82
BLOCK
2 provides several options to choose
good orbital orderings in DMRG.
(i)
For molecules with non-trivial point group symmetries, it can
group the orbitals by irreducible representation.
83,84
(ii)
For less symmetric systems, it can choose an ordering using
the Fiedler vector
85
of an appropriate metric matrix. This
ordering method is very inexpensive.
86,87
(iii)
The Fiedler ordering can be further optimized using a more
expensive genetic algorithm (GA),
87
as used in the S
TACK
-
B
LOCK
code.
60
In
BLOCK
2, this GA is reimplemented with
greatly increased efficiency.
The Fiedler and GA orderings require one to define a metric
to measure the correlation between pairs of orbitals.
BLOCK
2 can
use the absolute value of the exchange matrix
87
or the two-orbital
mutual information [computed from the one- and two-orbital
density matrices (1 and 2ODMs)]
88
for this purpose.
8. DMRG-CASSCF
For typical quantum chemistry problems, it is often not pos-
sible to treat all orbitals in DMRG. One can combine DMRG
with the Complete Active Space Self-Consistent Field (CASSCF)
method
89,90
to use DMRG to approximate the FCI problem within
the active space while optimizing the active space orbitals.
91–95
Using
BLOCK
2, we can perform DMRG-CASSCF through interfaces pro-
vided by external quantum chemistry packages, such as P
Y
SCF.
49,50
With P
Y
SCF, we currently support DMRG-CASSCF energy opti-
mization
92
for spin restricted and spin unrestricted orbitals and
DMRG-CASSCF gradients and geometry optimization
94
for spin
restricted orbitals.
B. High performance DMRG
Parallelization in
ab initio
DMRG is essential to realistic chem-
ical problems. In
BLOCK
2, we have implemented a hierarchy of
different parallel and related DMRG strategies.
1. Standard parallelization strategies
In a recent report,
96
we summarized five standard parallel
DMRG strategies implemented in
BLOCK
2: (i) parallelism over
the sum of sub-Hamiltonians,
97
(ii) parallelism over sites,
98,99
(iii)
parallelism over normal and complementary operators,
24,100
(iv)
parallelism over symmetry sectors,
11,101
and (v) parallelism within
dense matrix multiplications.
101,102
For large scale calculations, it is
essential to combine these strategies to achieve high performance.
96
We have implemented most of these strategies in
BLOCK
2 not
only for the ground state optimization problem but also for other
types of calculations, such as time evolution and the calculation of
response properties. The user can change the parallelization strategy
at runtime to maximize the efficiency of custom workflows.
Figure 3 shows a numerical example of the performance of
parallel DMRG.
2. Integral distribution strategies
In large
ab initio
calculations, the memory to store the two-
electron integrals of the Hamiltonian can be a bottleneck. In
BLOCK
2, we provide several integral/MPO distribution strategies
to balance the CPU, disk, and memory consumption in such
calculations.
FIG. 3.
Speed-up of average wall time per sweep relative to the
N
core
=
448 case
for different MPS bond dimensions using parallel DMRG algorithms for solving the
ground state in a (108o, 30e) active space of the benzene molecule
103
using a
cc-pVDZ basis.
104
Adapted from Fig. 4 in Ref. 96.
J. Chem. Phys.
159
, 234801 (2023); doi: 10.1063/5.0180424
159
, 234801-4
Published under an exclusive license by AIP Publishing
19 December 2023 20:42:18
The Journal
of Chemical Physics
ARTICLE
pubs.aip.org/aip/jcp
(i)
In the default strategy for small- to medium-sized problems,
integrals are replicated across different nodes, which pro-
vides the most flexibility for the DMRG algorithm. In this
case, we have the freedom to switch between the normal/
complementary (NC) parallelization strategy and the comple-
mentary/normal (CN) parallelization strategy near the mid-
dle site of a sweep,
11
which leads to higher computational
efficiency.
(ii)
For larger calculations where integral replication is problem-
atic, we provide an alternative strategy where each node only
stores the portion of integrals required by that node. This
strategy reduces the communication and memory costs at the
price of slightly reduced CPU efficiency.
96
(iii)
When the number of orbitals is large and localization has been
performed, there are many small or close-to-zero integral ele-
ments.
BLOCK
2 can store the integrals in a compressed format
to further reduce memory consumption using a floating point
number compression algorithm introduced in Ref. 96. (Note
that this only affects the memory consumption before the
Hamiltonian MPO is built.)
(iv)
For even larger scale applications, the MPOs constructed from
the sub-Hamiltonians can also consume a large amount of
memory on each node (since each MPO must contain at least
the same amount of information as contained in the integrals).
To alleviate this problem, we can store the MPO primarily on
disk and only the MPO tensors required at the given one or
two sites in an iteration of the sweep algorithm are loaded into
memory. Note that MPO storage is not a problem in tradi-
tional DMRG implementations, such as S
TACK
B
LOCK
, where
the action of the MPO is performed on the fly.
(v)
At the cost of introducing some error, the bond dimension
of the sub-Hamiltonian MPOs can be compressed during the
construction, as described in Sec. II C.
3. Parallelization strategies for specific tasks
The aforementioned parallelization strategies are common to
a wide range of DMRG algorithms. We have additionally imple-
mented more specialized parallelization strategies for some specific
computational tasks.
(i)
The generation of the Hamiltonian matrix elements in
large-site DMRG calculations (used, e.g., to treat multirefer-
ence dynamical correlation
58
) is implemented with shared-
memory parallelization over the determinants or CSFs.
(ii)
Green’s function calculations in dynamical DMRG
105
are par-
allelized in a distributed manner over the site index of the
one-particle Green’s function matrix element.
(iii)
Stochastic perturbative DMRG
106
is implemented with a
hybrid shared-memory and distributed parallelization to
sample determinants or CSFs.
(iv)
The GA algorithm for orbital reordering
87
is implemented
with distributed parallelization over the different initial
populations in a multi-start style GA algorithm.
(v)
The MPS compression step in the DMRG version of second
order strongly contracted
N
-Electron Valence States for Mul-
tireference Perturbation Theory (SC-NEVPT2)
81
is imple-
mented with distributed parallelization over virtual orbital
indices.
C. General operators and MPO construction
algorithms
BLOCK
2 provides full support for algorithms involving general
Hamiltonians (and other operators) expressed as MPOs. The imple-
mentation in
BLOCK
2 is carefully designed to support such general
operations while not (severely) sacrificing efficiency. As discussed in
Ref. 97, the original
ab initio
DMRG algorithm
107
implicitly provides
an efficient implementation of the expectation value of a sparse MPO
but does not explicitly construct the MPO. The advantages of an
explicit MPO object for relevant operators and Hamiltonians is that
a wider variety of algorithms can be easily implemented. However,
naive implementations of the MPO for complicated operators, such
as the
ab initio
Hamiltonian, lead to the wrong cost scaling.
97
Con-
sequently,
BLOCK
2 uses the techniques described below to achieve
high efficiency.
1. Symbolic MPO approach for ab initio Hamiltonians
To obtain performance competitive with traditional imple-
mentations of
ab initio
DMRG (such as the one found in S
TACK
-
B
LOCK
),
34,35
B
LOCK
2 provides a symbolic MPO class for
ab initio
Hamiltonians. Here, the multiplication of the MPO matrices is inter-
preted as a set of formulas that reproduce the DMRG equations in
traditional
ab initio
implementations. Symbolic MPOs are available
for different spin symmetries:
S
2
(SU2),
S
z
, and no spin symmetry,
and support the Normal/Complementary (NC) partition, the Com-
plementary/Normal (CN) partition, and the more efficient mixed
partition with a NC to CN transformation near the middle site.
11
In
addition, simplification rules are applied during the symbolic mul-
tiplication of the MPO matrices to use the index permutation and
transposition symmetry of the normal and complementary oper-
ators
11
(such as
ˆ
A
i j
=
−
ˆ
A
ji
, where
ˆ
A
i j
is defined as
a
†
i
a
†
j
). This
reduces the number of non-redundant renormalized operators that
have to be stored and computed during DMRG sweeps. With this
optimization, the DMRG algorithm implemented via the symbolic
MPO is more memory efficient than with MPOs constructed using
other approaches. In addition, because the symbols correspond to
the labels of renormalized operators in a traditional DMRG algo-
rithm, we can also implement conventional DMRG parallelization
strategies in the literature within this MPO approach.
97,100
2. General MPO construction for Hamiltonians
BLOCK
2 also implements general operator and Hamiltonian
MPOs efficiently. At the interface level, the user writes the definition
of the operator as a symbolic expression using characters in a string
to represent elementary operators, such as
ˆ
a
†
and
ˆ
a
. If non-standard
elementary operators are required, the user can define the matrix
representation of each elementary operator in a Python script.
Bosonic operators
108,109
and general Pauli strings
84
are supported.
Starting from such an arbitrary symbolic expression of a many-body
operator (such as the Hamiltonian) and the matrix representation
of the involved elementary operators (in any local orthonormal
basis), the code generates the MPO automatically. One can then
use the MPO to optimize or evaluate expectation values and carry
out other DMRG computations. Both spin-adapted and non-spin-
adapted MPOs can be generated. Listing 1 is an example Python
script for performing DMRG of the Hubbard–Holstein model using
BLOCK
2.
J. Chem. Phys.
159
, 234801 (2023); doi: 10.1063/5.0180424
159
, 234801-5
Published under an exclusive license by AIP Publishing
19 December 2023 20:42:18
The Journal
of Chemical Physics
ARTICLE
pubs.aip.org/aip/jcp
LISTING 1.
A Python script to setup the Hubbard–Holstein model Hamiltonian and perform DMRG to find the ground state using
BLOCK
2.
In the script, we use the characters
c, d, C, D, E,
and
F
(as an example) to represent the fermionic elementary operators
a
†
α
,
a
α
,
a
†
β
, and
a
β
and the bosonic elementary operators
b
†
and
b
, respectively. The expected output of this script (the ground state energy) is
−
6.956 893.
from pyblock2.driver.core import DMRGDriver, SymmetryTypes, MPOAlgorithmTypes
import numpy as np
N_SITES_ELEC, N_SITES_PH, N_ELEC
=
4, 4, 4
N_PH, U, OMEGA, G
=
11, 2, 0.25, 0.5
L
=
N_SITES_ELEC
+
N_SITES_PH
driver
=
DMRGDriver(scratch
=
"./tmp", symm_type
=
SymmetryTypes.SZ, n_threads
=
4)
driver.initialize_system(n_sites
=
L, n_elec
=
N_ELEC, spin
=
0)
# [Part A] Set states and matrix representation of operators in local Hilbert space
site_basis, site_ops
=
[], []
Q
=
driver.bw.SX # quantum number wrapper (n_elec, 2
∗
spin, point group irrep)
for k in range(L):
if k
<
N_SITES_ELEC:
# electron part
basis
=
[(Q(0, 0, 0), 1), (Q(1, 1, 0), 1), (Q(1, -1, 0), 1), (Q(2, 0, 0), 1)] # [0ab2]
ops
=
{
"": np.array([[1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 1, 0], [0, 0, 0, 1]]), # identity
"c": np.array([[0, 0, 0, 0], [1, 0, 0, 0], [0, 0, 0, 0], [0, 0, 1, 0]]), # alpha
+
"d": np.array([[0, 1, 0, 0], [0, 0, 0, 0], [0, 0, 0, 1], [0, 0, 0, 0]]), # alpha
"C": np.array([[0, 0, 0, 0], [0, 0, 0, 0], [1, 0, 0, 0], [0, -1, 0, 0]]), # beta
+
"D": np.array([[0, 0, 1, 0], [0, 0, 0, -1], [0, 0, 0, 0], [0, 0, 0, 0]]), # beta
}
else:
# phonon part
basis
=
[(Q(0, 0, 0), N_PH)]
ops
=
{
"": np.identity(N_PH), # identity
"E": np.diag(np.sqrt(np.arange(1, N_PH)), k
=
-1), # ph
+
"F": np.diag(np.sqrt(np.arange(1, N_PH)), k
=
1), # ph
}
site_basis.append(basis)
site_ops.append(ops)
# [Part B] Set Hamiltonian terms in Hubbard-Holstein model
driver.ghamil
=
driver.get_custom_hamiltonian(site_basis, site_ops)
b
=
driver.expr_builder()
# electron part
b.add_term("cd", np.array([[i, i
+
1, i
+
1, i] for i in range(N_SITES_ELEC - 1)]).ravel(), -1)
b.add_term("CD", np.array([[i, i
+
1, i
+
1, i] for i in range(N_SITES_ELEC - 1)]).ravel(), -1)
b.add_term("cdCD", np.array([[i, i, i, i] for i in range(N_SITES_ELEC)]).ravel(), U)
# phonon part
b.add_term("EF", np.array([[i
+
N_SITES_ELEC, ]
∗
2 for i in range(N_SITES_PH)]).ravel(), OMEGA)
# interaction part
b.add_term("cdE", np.array([[i, i, i
+
N_SITES_ELEC] for i in range(N_SITES_ELEC)]).ravel(), G)
b.add_term("cdF", np.array([[i, i, i
+
N_SITES_ELEC] for i in range(N_SITES_ELEC)]).ravel(), G)
b.add_term("CDE", np.array([[i, i, i
+
N_SITES_ELEC] for i in range(N_SITES_ELEC)]).ravel(), G)
b.add_term("CDF", np.array([[i, i, i
+
N_SITES_ELEC] for i in range(N_SITES_ELEC)]).ravel(), G)
# [Part C] Perform DMRG
mpo
=
driver.get_mpo(b.finalize(adjust_order
=
True), algo_type
=
MPOAlgorithmTypes.FastBipartite)
mps
=
driver.get_random_mps(tag
=
"KET", bond_dim
=
250, nroots
=
1)
energy
=
driver.dmrg(mpo, mps, n_sweeps
=
10, bond_dims
=
[250]
∗
4
+
[500]
∗
4,
noises
=
[1e-4]
∗
4
+
[1e-5]
∗
4
+
[0], thrds
=
[1e-10]
∗
8, dav_max_iter
=
30, iprint
=
2)
print("DMRG energy
=
%20.15f" % energy)
J. Chem. Phys.
159
, 234801 (2023); doi: 10.1063/5.0180424
159
, 234801-6
Published under an exclusive license by AIP Publishing
19 December 2023 20:42:18