quimb.tensor.tensor_dmrg#
DMRG-like variational algorithms, but in tensor network language.
Functions
|
Get a function to use as a callback for |
|
Get the default advanced settings for DMRG. |
|
Sort out the dims and inds of. |
Classes
|
Density Matrix Renormalization Group variational groundstate search. |
|
Simple alias of one site |
|
Simple alias of two site |
|
Class implmenting DMRG-X [1], whereby local effective energy eigenstates are chosen to maximise overlap with the previous step's state, leading to convergence on an mid-spectrum eigenstate of the full hamiltonian, as long as it is perturbatively close to the original state. |
|
Helper class for efficiently moving the effective 'environment' of a few sites in a 1D tensor network. E.g. for |
Exceptions
- class quimb.tensor.tensor_dmrg.DMRG(ham, bond_dims, cutoffs=1e-09, bsz=2, which='SA', p0=None)[source]#
Density Matrix Renormalization Group variational groundstate search. Some initialising arguments act as defaults, but can be overidden with each solve or sweep. See
get_default_opts()
for the list of advanced options initialized in theopts
attribute.- Parameters
ham (MatrixProductOperator) – The hamiltonian in MPO form.
bond_dims (int or sequence of ints.) – The bond-dimension of the MPS to optimize. If
bsz > 1
, then this corresponds to the maximum bond dimension when splitting the effective local groundstate. If a sequence is supplied then successive sweeps iterate through, then repeate the final value. E.g.[16, 32, 64] -> (16, 32, 64, 64, 64, ...)
.cutoffs (dict-like) – The cutoff threshold(s) to use when compressing. If a sequence is supplied then successive sweeps iterate through, then repeate the final value. E.g.
[1e-5, 1e-7, 1e-9] -> (1e-5, 1e-7, 1e-9, 1e-9, ...)
.bsz ({1, 2}) – Number of sites to optimize for locally i.e. DMRG1 or DMRG2.
which ({'SA', 'LA'}, optional) – Whether to search for smallest or largest real part eigenvectors.
p0 (MatrixProductState, optional) – If given, use as the initial state.
- state#
The current, optimized state.
- Type
- energies#
The total energy after each sweep.
- Type
list of float
- local_energies#
The local energies per sweep:
local_energies[i, j]
contains the local energy found at the jth step of the (i+1)th sweep.- Type
list of list of float
- total_energies#
The total energies per sweep:
local_energies[i, j]
contains the total energy after the jth step of the (i+1)th sweep.- Type
list of list of float
- opts#
Advanced options e.g. relating to the inner eigensolve or compression, see
get_default_opts()
.- Type
- (bond_sizes_ham)
If cyclic, the sizes of the energy environement transfer matrix bonds, per segment, per sweep.
- (bond_sizes_norm)
If cyclic, the sizes of the norm environement transfer matrix bonds, per segment, per sweep.
- form_local_ops(i, dims, lix, uix)[source]#
Construct the effective Hamiltonian, and if needed, norm.
- post_check(i, Neff, loc_gs, loc_en, loc_gs_old)[source]#
Perform some checks on the output of the local eigensolve.
- solve(tol=0.0001, bond_dims=None, cutoffs=None, sweep_sequence=None, max_sweeps=10, verbosity=0)[source]#
Solve the system with a sequence of sweeps, up to a certain absolute tolerance in the energy or maximum number of sweeps.
- Parameters
tol (float, optional) – The absolute tolerance to converge energy to.
bond_dims (int or sequence of int) – Overide the initial/current bond_dim sequence.
cutoffs (float of sequence of float) – Overide the initial/current cutoff sequence.
sweep_sequence (str, optional) – String made of ‘L’ and ‘R’ defining the sweep sequence, e.g ‘RRL’. The sequence will be repeated until
max_sweeps
is reached.max_sweeps (int, optional) – The maximum number of sweeps to perform.
verbosity ({0, 1, 2}, optional) – How much information to print about progress.
- Returns
converged – Whether the algorithm has converged to
tol
yet.- Return type
- sweep(direction, canonize=True, verbosity=0, **update_opts)[source]#
Perform a sweep of optimizations, either rightwards:
optimize --> ... >->-o-<-<-<-<-<-<-<-<-<-<-<-<-< | | | | | | | | | | | | | | | | H-H-H-H-H-H-H-H-H-H-H-H-H-H-H-H | | | | | | | | | | | | | | | | >->-o-<-<-<-<-<-<-<-<-<-<-<-<-<
or leftwards (direction=’L’):
<-- optimize ... >->->->->->->->->->->->->-o-<-< | | | | | | | | | | | | | | | | H-H-H-H-H-H-H-H-H-H-H-H-H-H-H-H | | | | | | | | | | | | | | | | >->->->->->->->->->->->->-o-<-<
After the sweep the state is left or right canonized respectively.
- Parameters
direction ({'R', 'L'}) – Sweep from left to right (->) or right to left (<-) respectively.
canonize (bool, optional) – Canonize the state first, not needed if doing alternate sweeps.
verbosity ({0, 1, 2}, optional) – Show a progress bar for the sweep.
update_opts – Supplied to
self._update_local_state
.
- class quimb.tensor.tensor_dmrg.DMRG1(ham, which='SA', bond_dims=None, cutoffs=1e-08, p0=None)[source]#
Simple alias of one site
DMRG
. Density Matrix Renormalization Group variational groundstate search. Some initialising arguments act as defaults, but can be overidden with each solve or sweep. Seeget_default_opts()
for the list of advanced options initialized in theopts
attribute.- Parameters
ham (MatrixProductOperator) – The hamiltonian in MPO form.
bond_dims (int or sequence of ints.) – The bond-dimension of the MPS to optimize. If
bsz > 1
, then this corresponds to the maximum bond dimension when splitting the effective local groundstate. If a sequence is supplied then successive sweeps iterate through, then repeate the final value. E.g.[16, 32, 64] -> (16, 32, 64, 64, 64, ...)
.cutoffs (dict-like) – The cutoff threshold(s) to use when compressing. If a sequence is supplied then successive sweeps iterate through, then repeate the final value. E.g.
[1e-5, 1e-7, 1e-9] -> (1e-5, 1e-7, 1e-9, 1e-9, ...)
.bsz ({1, 2}) – Number of sites to optimize for locally i.e. DMRG1 or DMRG2.
which ({'SA', 'LA'}, optional) – Whether to search for smallest or largest real part eigenvectors.
p0 (MatrixProductState, optional) – If given, use as the initial state.
- state#
The current, optimized state.
- Type
- energies#
The total energy after each sweep.
- Type
list of float
- local_energies#
The local energies per sweep:
local_energies[i, j]
contains the local energy found at the jth step of the (i+1)th sweep.- Type
list of list of float
- total_energies#
The total energies per sweep:
local_energies[i, j]
contains the total energy after the jth step of the (i+1)th sweep.- Type
list of list of float
- opts#
Advanced options e.g. relating to the inner eigensolve or compression, see
get_default_opts()
.- Type
- (bond_sizes_ham)
If cyclic, the sizes of the energy environement transfer matrix bonds, per segment, per sweep.
- class quimb.tensor.tensor_dmrg.DMRG2(ham, which='SA', bond_dims=None, cutoffs=1e-08, p0=None)[source]#
Simple alias of two site
DMRG
. Density Matrix Renormalization Group variational groundstate search. Some initialising arguments act as defaults, but can be overidden with each solve or sweep. Seeget_default_opts()
for the list of advanced options initialized in theopts
attribute.- Parameters
ham (MatrixProductOperator) – The hamiltonian in MPO form.
bond_dims (int or sequence of ints.) – The bond-dimension of the MPS to optimize. If
bsz > 1
, then this corresponds to the maximum bond dimension when splitting the effective local groundstate. If a sequence is supplied then successive sweeps iterate through, then repeate the final value. E.g.[16, 32, 64] -> (16, 32, 64, 64, 64, ...)
.cutoffs (dict-like) – The cutoff threshold(s) to use when compressing. If a sequence is supplied then successive sweeps iterate through, then repeate the final value. E.g.
[1e-5, 1e-7, 1e-9] -> (1e-5, 1e-7, 1e-9, 1e-9, ...)
.bsz ({1, 2}) – Number of sites to optimize for locally i.e. DMRG1 or DMRG2.
which ({'SA', 'LA'}, optional) – Whether to search for smallest or largest real part eigenvectors.
p0 (MatrixProductState, optional) – If given, use as the initial state.
- state#
The current, optimized state.
- Type
- energies#
The total energy after each sweep.
- Type
list of float
- local_energies#
The local energies per sweep:
local_energies[i, j]
contains the local energy found at the jth step of the (i+1)th sweep.- Type
list of list of float
- total_energies#
The total energies per sweep:
local_energies[i, j]
contains the total energy after the jth step of the (i+1)th sweep.- Type
list of list of float
- opts#
Advanced options e.g. relating to the inner eigensolve or compression, see
get_default_opts()
.- Type
- (bond_sizes_ham)
If cyclic, the sizes of the energy environement transfer matrix bonds, per segment, per sweep.
- class quimb.tensor.tensor_dmrg.DMRGX(ham, p0, bond_dims, cutoffs=1e-08, bsz=1)[source]#
Class implmenting DMRG-X [1], whereby local effective energy eigenstates are chosen to maximise overlap with the previous step’s state, leading to convergence on an mid-spectrum eigenstate of the full hamiltonian, as long as it is perturbatively close to the original state.
[1] Khemani, V., Pollmann, F. & Sondhi, S. L. Obtaining Highly Excited Eigenstates of Many-Body Localized Hamiltonians by the Density Matrix Renormalization Group Approach. Phys. Rev. Lett. 116, 247204 (2016).
- Parameters
ham (MatrixProductOperator) – The hamiltonian in MPO form, should have ~area-law eigenstates.
p0 (MatrixProductState) – The initial MPS guess, e.g. a computation basis state.
- k#
The current, optimized state.
- Type
- energies#
The list of energies after each sweep.
- Type
list of float
- form_local_ops(i, dims, lix, uix)[source]#
Construct the effective Hamiltonian, and if needed, norm.
- sweep(direction, canonize=True, verbosity=0, **update_opts)[source]#
Perform a sweep of the algorithm.
- Parameters
direction ({'R', 'L'}) – Sweep from left to right (->) or right to left (<-) respectively.
canonize (bool, optional) – Canonize the state first, not needed if doing alternate sweeps.
verbosity ({0, 1, 2}, optional) – Show a progress bar for the sweep.
update_opts – Supplied to
self._update_local_state
.
- class quimb.tensor.tensor_dmrg.MovingEnvironment(tn, begin, bsz, *, cyclic=False, segment_callbacks=None, ssz=0.5, eps=1e-08, method='isvd', max_bond=- 1, norm=False)[source]#
Helper class for efficiently moving the effective ‘environment’ of a few sites in a 1D tensor network. E.g. for
begin='left', bsz=2
, this initializes the right environments like so:n - 1: ●─●─●─ ─●─●─● │ │ │ │ │ │ H─H─H─ ... ─H─H─H │ │ │ │ │ │ ●─●─●─ ─●─●─● n - 2: ●─●─●─ ─●─●─╮ │ │ │ │ │ ● H─H─H─ ... ─H─H─H │ │ │ │ │ ● ●─●─●─ ─●─●─╯ n - 3: ●─●─●─ ─●─╮ │ │ │ │ ●● H─H─H─ ... ─H─HH │ │ │ │ ●● ●─●─●─ ─●─╯ ... 0 : ●─●─╮ │ │ ●● ●●● H─H─HH...HHH │ │ ●● ●●● ●─●─╯
which can then be used to efficiently generate the left environments as each site is updated. For example if
bsz=2
and the environements have been shifted many sites into the middle, thenMovingEnvironment()
returns something like:<---> bsz sites ╭─●─●─╮ ●●●●● │ │ ●●●●●●● HHHHH─H─H─HHHHHHH ●●●●● │ │ ●●●●●●● ╰─●─●─╯ 0 ... i i+1 ... n-1
For periodic systems
MovingEnvironment
approximates the ‘long way round’ transfer matrices. E.g consider replacing segment B (to arbitrary precision) with an SVD:╭───────────────────────────────────────────────╮ ╰─A─A─A─A─A─A─A─A─A─A─A─A─B─B─B─B─B─B─B─B─B─B─B─╯ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ ==> ╭─A─A─A─A─A─A─A─A─A─A─A─A─B─B─B─B─B─B─B─B─B─B─B─╮ ╰───────────────────────────────────────────────╯ ╭┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄╮ ┊ ╭─A─A─A─A─A─A─A─A─A─A─A─A─╮ ┊ ==> ╰┄<BL │ │ │ │ │ │ │ │ │ │ │ │ BR>┄╯ ╰─A─A─A─A─A─A─A─A─A─A─A─A─╯ ^ ^ segment_start segment_stop - 1 ╭┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄╮ ┊ ╭─A─A─╮ ┊ ==> ╰┄<BL │ │ AAAAAAAAAAAAAAAAAAAAABR>┄╯ ╰─A─A─╯ ... <-bsz-> ╭┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄╮ ┊ ╭─A─A─╮ ┊ ==> ╰~<BLAAAAAAAAAAA │ │ AAAAAAAABR>~╯ ╰─A─A─╯ i i+1 -----sweep--------->
Can then contract and store left and right environments for efficient sweeping just as in non-periodic case. If the segment is long enough (50+) sites, often only 1 singular value is needed, and thus the efficiency is the same as for OBC.
- Parameters
tn (TensorNetwork) – The 1D tensor network, should be closed, i.e. an overlap of some sort.
begin ({'left', 'right'}) – Which side to start at and sweep from.
bsz (int) – The number of sites that form the the ‘non-environment’, e.g. 2 for DMRG2.
ssz (float or int, optional) – The size of the segment to use, if float, the proportion. Default: 1/2.
eps (float, optional) – The tolerance to approximate the transfer matrix with. See
replace_with_svd()
.cyclic (bool, optional) – Whether this is a periodic
MovingEnvironment
.segment_callbacks (sequence of callable, optional) – Functions with signature
callback(start, stop, self.begin)
, to be called every time a new segment is initialized.method ({'isvd', 'svds', ...}, optional) – How to performt the transfer matrix compression. See
replace_with_svd()
.max_bond (, optional) – If > 0, the maximum bond of the compressed transfer matrix.
norm (bool, optional) – If True, treat this
MovingEnvironment
as the state overlap, which enables a few extra checks.
Notes
Does not necessarily need to be an operator overlap tensor network. Useful for any kind of sweep where only local tensor updates are being made. Note that only the current site is completely up-to-date and can be modified with changes meant to propagate.
- init_non_segment(start, stop)[source]#
Compress and label the effective env not in
range(start, stop)
if cyclic, else just add some dummy left and right end pieces.
- quimb.tensor.tensor_dmrg.get_cyclic_canonizer(k, b, inv_tol=1e-10)[source]#
Get a function to use as a callback for
MovingEnvironment
that approximately orthogonalizes the segments of periodic MPS.
- quimb.tensor.tensor_dmrg.get_default_opts(cyclic=False)[source]#
Get the default advanced settings for DMRG.
- Returns
default_sweep_sequence (str) – How to sweep. Will be repeated, e.g. “RRL” -> RRLRRLRRL…, default: R.
bond_compress_method ({‘svd’, ‘eig’, …}) – Method used to compress sites after update.
bond_compress_cutoff_mode ({‘sum2’, ‘abs’, ‘rel’}) – How to perform compression truncation.
bond_expand_rand_strength (float) – In DMRG1, strength of randomness to expand bonds with. Needed to avoid singular matrices after expansion.
local_eig_tol (float) – Relative tolerance to solve inner eigenproblem to, larger = quicker but more unstable, default: 1e-3. Note this can be much looser than the overall tolerance, the starting point for each local solve is the previous state, and the overall accuracy comes from multiple sweeps.
local_eig_ncv (int) – Number of inner eigenproblem lanczos vectors. Smaller can mean quicker.
local_eig_backend ({None, ‘AUTO’, ‘SCIPY’, ‘SLEPC’}) – Which to backend to use for the inner eigenproblem. None or ‘AUTO’ to choose best. Generally
'SLEPC'
best if available for large problems, but it can’t currently handleLinearOperator
Neff as well as'lobpcg'
.local_eig_maxiter (int) – Maximum number of inner eigenproblem iterations.
local_eig_ham_dense (bool) – Force dense representation of the effective hamiltonian.
local_eig_EPSType ({‘krylovschur’, ‘gd’, ‘jd’, …}) – Eigensovler tpye if
local_eig_backend='slepc'
.local_eig_norm_dense (bool) – Force dense representation of the effective norm.
periodic_segment_size (float or int) – How large (as a proportion if float) to make the ‘segments’ in periodic DMRG. During a sweep everything outside this (the ‘long way round’) is compressed so the effective energy and norm can be efficiently formed. Tradeoff: longer segments means having to compress less, but also having a shorter ‘long way round’, meaning that it needs a larger bond to represent it and can be ‘pseudo-orthogonalized’ less effectively. 0.5 is the largest fraction that makes sense. Set to >= 1.0 to not use segmentation at all, which is better for small systems.
periodic_compress_method ({‘isvd’, ‘svds’}) – Which method to perform the transfer matrix compression with.
periodic_compress_norm_eps (float) – Precision to compress the norm transfer matrix in periodic systems.
periodic_compress_ham_eps (float) – Precision to compress the energy transfer matrix in periodic systems.
periodic_compress_max_bond (int) – The maximum bond to use when compressing transfer matrices.
periodic_nullspace_fudge_factor (float) – Factor to add to
Heff
andNeff
to remove nullspace.periodic_canonize_inv_tol (float) – When psuedo-orthogonalizing, an inverse gauge is generated that can be very ill-conditioned. This factor controls cutting off the small singular values of the gauge to stop this.
periodic_orthog_tol (float) – When psuedo-orthogonalizing, if the local norm is within this distance to 1 (pseudo-orthogonoalized), then the generalized eigen decomposition is not used, which is much more efficient. If set too large the total normalization can become unstable.