MPO¶
full name: tenpy.networks.mpo.MPO
parent module:
tenpy.networks.mpo
type: class
Inheritance Diagram

Methods
|
|
|
Apply self to an MPS psi and compress psi in place. |
|
Applies an MPO to an MPS (in place) naively, without compression. |
|
Applies an MPO to an MPS (in place) with the zip-up method. |
|
Make a shallow copy of self. |
Return hermition conjugate copy of self. |
|
|
Repeat the unit cell for infinite MPS boundary conditions; in place. |
|
Calculate |
|
Calculate |
|
Calculate |
|
Calculate |
|
Extract a segment from the MPO. |
|
Initialize an MPO from grids. |
|
Load instance from a HDF5 file. |
|
Create a (finite) MPO wave packet representing |
|
Return index of IdL at bond to the left of site i. |
|
Return index of IdR at bond to the right of site i. |
|
Return W at site i. |
|
extract the full Hamiltonian as a |
|
group each blocklen subsequent tensors and return result as a new MPO. |
|
Modify self inplace to group sites. |
|
Check if self and other represent the same MPO to precision eps. |
|
Check if self is a hermitian MPO. |
|
Creates the U_I or U_II propagator. |
|
Creates the \(U_I\) propagator with W_I tensors. |
|
Creates the \(U_II\) propagator. |
|
Export self into a HDF5 file. |
|
Set W at site i. |
Sort virtual legs by charges. |
|
Sanity check, raises ValueErrors, if something is wrong. |
|
|
Calculate |
Class Attributes and Properties
Number of physical sites; for an iMPO the len of the MPO unit cell. |
|
Dimensions of the virtual bonds. |
|
List of local physical dimensions. |
|
Distinguish MPO vs iMPO. |
- class tenpy.networks.mpo.MPO(sites, Ws, bc='finite', IdL=None, IdR=None, max_range=None, explicit_plus_hc=False)[source]¶
Bases:
object
Matrix product operator, finite (MPO) or infinite (iMPO).
- Parameters
sites (list of
Site
) – Defines the local Hilbert space for each site.Ws (list of
Array
) – The matrices of the MPO. Should have labelswL, wR, p, p*
.bc ({'finite' | 'segment' | 'infinite'}) – Boundary conditions as described in
mps
.'finite'
requiresWs[0].get_leg('wL').ind_len = 1
.IdL ((iterable of) {int | None}) – Indices on the bonds, which correspond to ‘only identities to the left’. A single entry holds for all bonds.
IdR ((iterable of) {int | None}) – Indices on the bonds, which correspond to ‘only identities to the right’.
max_range (int | np.inf | None) – Maximum range of hopping/interactions (in unit of sites) of the MPO.
None
for unknown.explicit_plus_hc (bool) – If True, this flag indicates that the hermitian conjugate of the MPO should be computed and added at runtime, i.e., self is not (necessarily) hermitian.
- chinfo¶
The nature of the charge.
- Type
ChargeInfo
- sites¶
Defines the local Hilbert space for each site.
- Type
list of
Site
- bc¶
Boundary conditions as described in
mps
.'finite'
requiresWs[0].get_leg('wL').ind_len = 1
.- Type
{‘finite’ | ‘segment’ | ‘infinite’}
- IdL¶
Indices on the bonds (length L`+1), which correspond to ‘only identities to the left’. ``None` for bonds where it is not set. In standard form, this is 0 (except for unset bonds in finite case)
- Type
list of {int | None}
- IdR¶
Indices on the bonds (length L`+1), which correspond to ‘only identities to the right’. ``None` for bonds where it is not set. In standard form, this is the last index on the bond (except for unset bonds in finite case).
- Type
list of {int | None}
- max_range¶
Maximum range of hopping/interactions (in unit of sites) of the MPO.
None
for unknown.- Type
int | np.inf | None
- grouped¶
Number of sites grouped together, see
group_sites()
.- Type
- explicit_plus_hc¶
If True, this flag indicates that the hermitian conjugate of the MPO should be computed and added at runtime, i.e., self is not (necessarily) hermitian.
- Type
- save_hdf5(hdf5_saver, h5gr, subpath)[source]¶
Export self into a HDF5 file.
This method saves all the data it needs to reconstruct self with
from_hdf5()
.Specifically, it saves
sites
,chinfo
,max_range
(under these names),_W
as"tensors"
,IdL
as"index_identity_left"
,IdR
as"index_identity_right"
, andbc
as"boundary_condition"
. Moreover, it savesL
,explicit_plus_hc
andgrouped
as HDF5 attributes, as well as the maximum ofchi
under the namemax_bond_dimension
.
- classmethod from_hdf5(hdf5_loader, h5gr, subpath)[source]¶
Load instance from a HDF5 file.
This method reconstructs a class instance from the data saved with
save_hdf5()
.- Parameters
hdf5_loader (
Hdf5Loader
) – Instance of the loading engine.h5gr (
Group
) – HDF5 group which is represent the object to be constructed.subpath (str) – The name of h5gr with a
'/'
in the end.
- Returns
obj – Newly generated class instance containing the required data.
- Return type
cls
- classmethod from_grids(sites, grids, bc='finite', IdL=None, IdR=None, Ws_qtotal=None, legs=None, max_range=None, explicit_plus_hc=False)[source]¶
Initialize an MPO from grids.
- Parameters
sites (list of
Site
) – Defines the local Hilbert space for each site.grids (list of list of list of entries) – For each site (outer-most list) a matrix-grid (corresponding to
wL, wR
) with entries being or representing (seegrid_insert_ops()
) onsite-operators.bc ({'finite' | 'segment' | 'infinite'}) – Boundary conditions as described in
mps
.IdL ((iterable of) {int | None}) – Indices on the bonds, which correspond to ‘only identities to the left’. A single entry holds for all bonds.
IdR ((iterable of) {int | None}) – Indices on the bonds, which correspond to ‘only identities to the right’.
Ws_qtotal ((list of) total charge) – The qtotal to be used for each grid. Defaults to zero charges.
legs (list of
LegCharge
) – List of charges for ‘wL’ legs left of each W, L + 1 entries. The last entry should be the conjugate of the ‘wR’ leg, i.e. identical tolegs[0]
for ‘infinite’ bc. By default, determine the charges automatically. This is limited to cases where there are no “dangling open ends” in the MPO graph. (TheMPOGraph
can handle those cases, though.)max_range (int | np.inf | None) – Maximum range of hopping/interactions (in unit of sites) of the MPO.
None
for unknown.explicit_plus_hc (bool) – If True, the Hermitian conjugate of the MPO is computed at runtime, rather than saved in the MPO.
See also
grid_insert_ops
used to plug in entries of the grid.
tenpy.linalg.np_conserved.grid_outer
used for final conversion.
- classmethod from_wavepacket(sites, coeff, op, eps=1e-15)[source]¶
Create a (finite) MPO wave packet representing
sum_i coeff[i] op_i
.Note that we define it only for finite systems; a generalization to fininite systems is not straight forward due to normalization issues: the individual terms vanish in the thermodynamic limit!
- Parameters
Examples
Say you have fermions, so
op='Cd'
, and want to create a gaussian wave paket \(\sum_x \alpha_x c^\dagger_x\) with \(\alpha_x \propto e^{-0.5(x-x_0)^2/\sigma^2} e^{i k_0 x}\). Then you would useIndeed, we can apply this to a (vacuum) MPS and get the correct state:
- property L¶
Number of physical sites; for an iMPO the len of the MPO unit cell.
- property dim¶
List of local physical dimensions.
- property finite¶
Distinguish MPO vs iMPO.
True for an MPO (
bc='finite', 'segment'
), False for an iMPO (bc='infinite'
).
- property chi¶
Dimensions of the virtual bonds.
- enlarge_mps_unit_cell(factor=2)[source]¶
Repeat the unit cell for infinite MPS boundary conditions; in place.
- Parameters
factor (int) – The new number of sites in the unit cell will be increased from L to
factor*L
.
- group_sites(n=2, grouped_sites=None)[source]¶
Modify self inplace to group sites.
Group each n sites together using the
GroupedSite
. This might allow to do TEBD with a Trotter decomposition, or help the convergence of DMRG (in case of too long range interactions).- Parameters
n (int) – Number of sites to be grouped together.
grouped_sites (None | list of
GroupedSite
) – The sites grouped together.
- extract_segment(first, last)[source]¶
Extract a segment from the MPO.
- Parameters
- Returns
cp – A copy of self with “segment” boundary conditions.
- Return type
See also
tenpy.networks.mps.MPS.extract_segment
similar method for MPS.
- sort_legcharges()[source]¶
Sort virtual legs by charges. In place.
The MPO seen as matrix of the
wL, wR
legs is usually very sparse. This sparsity is captured by the LegCharges for these bonds not being sorted and bunched. This requires a tensordot to do more block-multiplications with smaller blocks. This is in general faster for large blocks, but might lead to a larger overhead for small blocks. Therefore, this function allows to sort the virtual legs by charges.
- make_U(dt, approximation='II')[source]¶
Creates the U_I or U_II propagator.
Approximations of MPO exponentials following [zaletel2015].
- Parameters
dt (float|complex) – The time step per application of the propagator. Should be imaginary for real time evolution!
approximation (
'I' | 'II'
) – Selects the approximation,make_U_I()
('I'
) ormake_U_II()
('II'
).
- Returns
U – The propagator, i.e. approximation \(U ~= exp(H*dt)\)
- Return type
MPO
- expectation_value(psi, tol=1e-10, max_range=100, init_env_data={})[source]¶
Calculate
<psi|self|psi>/<psi|psi>
(or density for infinite).For infinite MPS, it assumes that self is extensive, e.g. a Hamiltonian but not a unitary, and returns the expectation value density. For finite MPS, it just returns the total value.
This function is just a small wrapper around
expectation_value_finite()
,expectation_value_powermethod()
orexpectation_value_transfer_matrix()
.- Parameters
- Returns
exp_val – The expectation value of self with respect to the state psi. For an infinite MPS: the (energy) density per site.
- Return type
float/complex
- expectation_value_finite(psi, init_env_data={})[source]¶
Calculate
<psi|self|psi>/<psi|psi>
for finite MPS.
- expectation_value_TM(psi, tol=1e-10, init_env_data={})[source]¶
Calculate
<psi|self|psi>/<psi|psi> / L
from the MPOTransferMatrix.Only for infinite MPS, and assumes that the Hamiltonian is an extensive sum of (quasi)local terms, and that the MPO has all
IdL
andIdR
defined.Diagonalizing the
MPOTransferMatrix
allows to find energy densities for infinite systems even for hamiltonians with infinite (exponentially decaying) range.- Parameters
- Returns
exp_val – The expectation value density of self with respect to the state psi.
- Return type
float/complex
- expectation_value_power(psi, tol=1e-10, max_range=100)[source]¶
Calculate
<psi|self|psi>/<psi|psi>
with a power-method.Only for infinite MPS, and assumes that the Hamiltonian is an extensive sum of (quasi)local terms, and that the MPO has all
IdL
andIdR
defined. Only for infinite MPS.Instead of diagonalizing the MPOTransferMatrix like
expectation_value_TM()
, this method uses just considers terms of the MPO starting in the first unit cell and then continues to contract tensors until convergence. For infinite-range MPOs, this converges like a power-method (i.e. slower thanexpectation_value_TM()
), but for finite-range MPOs it’s likely faster, and conceptually cleaner.- Parameters
psi (
MPS
) – The state in which to calculate the expectation value.tol (float) – For infinite MPO containing exponentially decaying long-range terms, stop evaluating further terms if the terms in LP have norm < tol.
max_range (int) – Ignored for finite psi. Contract at most
self.L * max_range
sites, even if tol is not reached. In that case, issue a warning.
- Returns
exp_val – The expectation value of self with respect to the state psi. For an infinite MPS: the density per site.
- Return type
float/complex
- variance(psi, exp_val=None)[source]¶
Calculate
<psi|self^2|psi> - <psi|self|psi>^2
.Works only for finite systems. Ignores the
norm
of psi.Todo
This is a naive, expensive implementation contracting the full network. Try to follow arXiv:1711.01104 for a better estimate; would that even work in the infinite limit?
- Parameters
psi (
MPS
) – State for which the variance should be taken.exp_val (float/complex | None) – The result of
<psi|self|psi> = self.expectation_value(psi)
if known; otherwise obtained fromexpectation_value()
. (Set this to 0 to obtain only the part<psi|self^2|psi>
.)
- is_hermitian(eps=1e-10, max_range=None)[source]¶
Check if self is a hermitian MPO.
Shorthand for
self.is_equal(self.dagger(), eps, max_range)
.
- is_equal(other, eps=1e-10, max_range=None)[source]¶
Check if self and other represent the same MPO to precision eps.
To compare them efficiently we view self and other as MPS and compare the overlaps
abs(<self|self> + <other|other> - 2 Re(<self|other>)) < eps*(<self|self>+<other|other>)
- Parameters
other (
MPO
) – The MPO to compare to.eps (float) – Precision threshold what counts as zero.
max_range (None | int) – Ignored for finite MPS; for finite MPS we consider only the terms contained in the sites with indices
range(self.L + max_range)
. None defaults tomax_range
(orL
in case this is infinite or None).
- Returns
equal – Whether self equals other to the desired precision.
- Return type
- apply(psi, options)[source]¶
Apply self to an MPS psi and compress psi in place.
For infinite MPS, the assumed form of self is a product (e.g. a time evolution operator \(U= e^{-iH dt}\), not an (extensive) sum as a Hamiltonian would have. See A warning about infinite MPS for more details.
Options
- config ApplyMPO¶
option summary By default (``None``) this feature is disabled. [...]
chi_list_reactivates_mixer (from Sweep) in Sweep.sweep
If True, the mixer is reset/reactivated each time the bond dimension growth [...]
Whether to combine legs into pipes. This combines the virtual and [...]
Mandatory. [...]
init_env_data (from Sweep) in DMRGEngine.init_env
Dictionary as returned by ``self.env.get_initialization_data()`` from [...]
lanczos_params (from Sweep) in Sweep
Lanczos parameters as described in :cfg:config:`Lanczos`.
m_temp (from ZipUpApplyMPO) in MPO.apply_zipup
bond dimension will be truncated to `m_temp * chi_max`
max_sweeps (from VariationalCompression) in VariationalCompression
Minimum and maximum number of sweeps to perform for the compression.
min_sweeps (from VariationalCompression) in VariationalCompression
Minimum and maximum number of sweeps to perform for the compression.
Specifies which :class:`Mixer` to use, if any. [...]
mixer_params (from Sweep) in DMRGEngine.mixer_activate
Mixer parameters as described in :cfg:config:`Mixer`.
orthogonal_to (from Sweep) in DMRGEngine.init_env
Deprecated in favor of the `orthogonal_to` function argument (forwarded fro [...]
Number of sweeps to be performed without optimization to update the environment.
start_env_sites (from VariationalCompression) in VariationalCompression
Number of sites to contract for the inital LP/RP environment in case of inf [...]
tol_theta_diff (from VariationalCompression) in VariationalCompression
Stop after less than `max_sweeps` sweeps if the 1-site wave function change [...]
Truncation parameters as described in :cfg:config:`truncation`.
trunc_weight (from ZipUpApplyMPO) in MPO.apply_zipup
reduces cut for Schmidt values to `trunc_weight * svd_min`
-
option compression_method:
'SVD' | 'variational' | 'zip_up'
¶ Mandatory. Selects the method to be used for compression. For the SVD compression, trunc_params is the only other option used.
- option trunc_params: dict¶
Truncation parameters as described in
truncation
.
-
option compression_method:
- apply_naively(psi)[source]¶
Applies an MPO to an MPS (in place) naively, without compression.
This function simply contracts the W tensors of the MPO to the B tensors of the MPS, resulting in an MPS with bond dimension self.chi * psi.chi.
Warning
This function sets only a wild guess for the new singular values. You should either compress the MPS or at least call
canonical_form()
. If you useapply()
instead, this will be done automatically.- Parameters
psi (
MPS
) – The MPS to which self should be applied. Modified in place!
- apply_zipup(psi, options)[source]¶
Applies an MPO to an MPS (in place) with the zip-up method.
Described in Ref. [stoudenmire2010].
The ‘W’ tensors are contracted to the ‘B’ tensors with intermediate SVD compressions, truncated to bond dimensions chi_max * m_temp.
Warning
The MPS afterwards is only approximately in canonical form (under the assumption that self is close to unity). You should either compress the MPS or at least call
canonical_form()
. If you useapply()
instead, this will be done automatically.- Parameters
psi (
MPS
) – The MPS to which self should be applied. Modified in place!trunc_params (dict) – Truncation parameters as described in
truncation
.
Options
- config ZipUpApplyMPO¶
option summary bond dimension will be truncated to `m_temp * chi_max`
Truncation parameters as described in :cfg:config:`truncation`.
reduces cut for Schmidt values to `trunc_weight * svd_min`
- option trunc_params: dict¶
Truncation parameters as described in
truncation
.
- option m_temp: int¶
bond dimension will be truncated to m_temp * chi_max
- option trunc_weight: float¶
reduces cut for Schmidt values to trunc_weight * svd_min
- get_grouped_mpo(blocklen)[source]¶
group each blocklen subsequent tensors and return result as a new MPO.
Deprecated since version 0.5.0: Make a copy and use
group_sites()
instead.
- get_full_hamiltonian(maxsize=1000000.0)[source]¶
extract the full Hamiltonian as a
d**L``x``d**L
matrix.Deprecated since version 0.5.0: Use
tenpy.algorithms.exact_diag.ExactDiag.from_H_mpo()
instead.