# Algorithms for MPS and MPO

## Summing MPS

`sum(MPS A1, MPS A2, Args args = Args::global()) -> MPS`

Return the sum of the MPS

`A1`

and`A2`

. The returned MPS will have an orthogonality center on site 1. Before being returned, the MPS representing the sum will be compressed using truncation parameters provided in the named arguments`args`

.The input MPS must have the same site indices. The link indices of the output MPS will have the same tags as the first input MPS.

Show Exampleauto A3 = sum(A1,A2,{"MaxDim",500,"Cutoff",1E-8});

`sum(vector<MPS> terms, Args args = Args::global()) -> MPS`

Returns the sum of all the MPS provided in the vector

`terms`

as a single MPS, using the truncation accuracy parameters (such as "Cutoff" or "MaxDim") provided in the named arguments`args`

to control the accuracy of the sum.This function uses a hierarchical, tree-like algorithm which first sums pairs of MPS, then pairs of pairs, etc. so that the largest bond dimensions are only reached toward the end of the process for maximum efficiency. Therefore using this algorithm can be much faster than calling the above two-argument

`sum`

function to sum the terms one at a time.Show Exampleauto terms = vector<MPS>(4); terms.at(0) = A0; terms.at(1) = A1; terms.at(2) = A2; terms.at(3) = A3; auto res = sum(terms,{"Cutoff",1E-8});

## Measuring Properties of MPS

`expect(MPS psi, SiteSet sites, string A) -> vector<Real>`

`expectC(MPS psi, SiteSet sites, string A) -> vector<Cplx>`

`expect(MPS psi, SiteSet sites, string A, vector<int> site_list) -> vector<Real>`

`expectC(MPS psi, SiteSet sites, string A, vector<int> site_list) -> vector<Cplx>`

Compute the expected value @@\langle\psi| \hat{A} |\psi\rangle@@ for an operator

`A`

on every site, returning a`vector`

of the results. The operator can be any operator name recognized by the`SiteSet`

, such as`"Sz"`

for the`SpinHalf`

SiteSet. The function`expectC`

returns a vector of complex values and should be used when complex results are expected.You can pass a vector

`site_list`

of integers as an optional last argument to obtain expectation values only for sites in that vector.Show Exampleauto exsz = expect(psi,sites,"Sz"); //Measure only sites 3,5,8 auto exsz = expect(psi,sites,"Sz",{{3,5,8}});

`expect(MPS psi, SiteSet sites, vector<string> ops) -> vector<vector<Real>>`

`expectC(MPS psi, SiteSet sites, vector<string> ops) -> vector<vector<Cplx>>`

`expect(MPS psi, SiteSet sites, vector<string> ops, vector<int> site_list) -> vector<vector<Real>>`

`expectC(MPS psi, SiteSet sites, vector<string> ops, vector<int> site_list) -> vector<vector<Cplx>>`

Compute the expected value of a set of operators

`ops`

for each site of the MPS`psi`

. Similar to calling the single-operator version of`expect`

above but is more efficient. The returned value is a vector of vectors, one for each operator, containing the expected values on each site.You can pass a vector

`site_list`

of integers as an optional last argument to obtain expectation values only for sites in that vector.Show Exampleauto ex = expect(psi,sites,{"Sz","Sx","Sy"}); //Measure only sites 3,5,8 auto ex = expect(psi,sites,{"Sz","Sx","Sy"},{{3,5,8}});

`correlationMatrix(MPS psi, SiteSet sites, string A, string B) -> vector<vector<Real>>`

`correlationMatrixC(MPS psi, SiteSet sites, string A, string B) -> vector<vector<Cplx>>`

Given an MPS

`psi`

, a`SiteSet`

, and two strings denoting operators computes the two-point correlation function matrix @@M_{ij} = \langle\psi| \hat{A}_i \hat{B}_j |\psi\rangle@@

using efficient MPS techniques. Returns the matrix M.Show Exampleauto czz = correlationMatrix(psi,sites,"Sz","Sz"); auto cpm = correlationMatrix(psi,sites,"S+","S-");

## Inner Products and Expectation Values

`inner(MPS y, MPS x) -> Real`

`innerC(MPS y, MPS x) -> Cplx`

Compute the exact inner product @@\langle y|x \rangle@@ of two MPS (the tensors of

`y`

will get conjugated). If the inner product is expected to be a complex number use`innerC`

.The algorithm used scales as @@m^3 d@@ where @@m@@ is a typical link dimension of the MPS and @@d@@ is the site dimension.

Note that if

`x`

and`y`

don't have the same site indices, this function will attempt to make them match.`inner(MPS y, MPO A, MPS x) -> Real`

`innerC(MPS y, MPO A, MPS x) -> Cplx`

Compute the exact inner product @@\langle y|A|x \rangle@@ of two MPS

`y`

and`x`

with respect to an MPO`A`

(the tensors of`y`

will get conjugated).The algorithm used scales as @@m^3\, k\,d + m^2\, k^2\, d^2@@ where @@m@@ is typical link dimension of the MPS, @@k@@ is the typical MPO dimension, and @@d@@ is the site dimension.

Note that

`A`

and`x`

must share a set of site indices. If the remaining site indices of`A`

are not shared with`y`

, this function will attempt to match them (i.e. it has the same behavior as`inner(y,Ax)`

if`Ax`

was the exact application of MPO`A`

to MPS`x`

).`inner(MPS y, MPO B, MPO A, MPS x) -> Real`

`innerC(MPS y, MPO B, MPO A, MPS x) -> Cplx`

Compute the exact inner product @@\langle y|BA|x \rangle@@ of two MPS

`y`

and`x`

with respect to two MPOs`B`

and`A`

(MPS`y`

will get conjugated).MPO

`A`

must share a set of site indices with MPS`x`

, and the other set of site indices with MPO`B`

. If the remaining set of site indices of`B`

are not shared with`y`

, with function will attempt to make them match.The algorithm used scales as @@m^3\, k^2\,d + m^2\, k^3\, d^2@@ where @@m@@ is typical bond dimension of the MPS, @@k@@ is the typical MPO dimension, and @@d@@ is the site dimension.

`inner(MPO B, MPS y, MPO A, MPS x) -> Real`

`innerC(MPO B, MPS y, MPO A, MPS x) -> Cplx`

Compute the exact inner product @@\langle By|A|x \rangle@@ (i.e. the inner product of of @@B|y \rangle@@ and @@A|x \rangle@@ ). MPO

`B`

and MPS`y`

will get conjugated.MPO

`A`

must share a set of site indices with MPS`x`

, and MPO`B`

must share a set of site indices with MPS`y`

. If the remaining site indices of`A`

and`B`

do not match with each other, the function will attempt to make them match.The algorithm used scales as @@m^3\, k^2\,d + m^2\, k^3\, d^2@@ where @@m@@ is typical bond dimension of the MPS, @@k@@ is the typical MPO dimension, and @@d@@ is the site dimension.

## Tracing an MPO

`trace(MPO A) -> Real`

`traceC(MPO A) -> Cplx`

Trace over the site indices of the MPO.

`trace(MPO A, MPO B) -> Real`

`traceC(MPO A, MPO B) -> Cplx`

Return the trace of the operator that would result from performing the exact contraction of MPO

`A`

with MPO`B`

. For each`j`

,`A(j)`

and`B(j)`

must share one or two site indices.Note that neither

`A`

or`B`

will get conjugated by this function.

## Multiplying MPOs

`nmultMPO(MPO A, MPO B, Args args = Args::global()) -> MPO`

Multiply MPOs

`A`

and`B`

, returning the results MPO. MPO tensors are multiplied one at a time from left to right and the resulting tensors are compressed using the truncation parameters (such as "Cutoff" and "MaxDim") provided through the named arguments`args`

.For each

`j`

, MPO tensors`A(j)`

and`B(j)`

must share a single site index. MPO`C`

will contain the site indices not shared by MPOs`A`

and`B`

. In addition, the link indices of MPO`C`

will have the same tags as the link indices of the MPO`A`

.Show Exampleauto sites = SiteSet(10,2); // Make trivial MPOs auto A = MPO(sites); auto B = MPO(sites); // Prime MPO A to ensure only one set of site indices are shared auto C = nmultMPO(prime(A),B,{"MaxDim",500,"Cutoff",1E-8}); auto s3 = sites(3); Print(hasInds(C(3),{s3,prime(s3,2)})); //print: true

## Applying MPO to MPS

`applyMPO(MPO A, MPS x, Args args = Args::global()) -> MPS`

Apply an MPO

`A`

to an MPS`x`

, resulting in an approximation to the MPS`y`

:

@@|y\rangle = A |x\rangle@@ .

The resulting MPS is returned. The algorithm used is chosen with the parameter "Method" in the named arguments`args`

.MPO

`A`

and MPS`x`

must share a set of site indices. The links of the output MPS will have the same tags as the links of the input MPS`x`

. The site indices of the output MPS will be the site indices of`A`

that are not shared with`x`

.The default algorithm used is the "density matrix" algorithm, chosen by setting the parameter "Method" to "DensityMatrix". If the input MPS has a typical bond dimension of @@m@@ and the MPO has typical bond dimension @@k@@ , this algorithm scales as @@m^3 k^2 + m^2 k^3@@ .

No approximation is made when applying the MPO, but after applying it the resulting MPS is compressed using the truncation parameters provided in the named arguments

`args`

.An alternative algorithm can be chosen by setting the parameter "Method" to "Fit". This is a sweeping algorithm that iteratively optimizes the resulting MPS @@|y\rangle@@ (analogous to DMRG). The scaling of the "Fit" method is @@m^3 k + m^2 k^2@@ where @@m@@ is the typical MPS bond dimension and @@k@@ is the typical MPO bond dimension. Thus this algorithm has better scaling in the MPO bond dimension @@k@@ compared to the "DensityMatrix" method, but is not guaranteed to converge (depending on the input MPO and MPS). To ensure convergence it is helpful to use an initial guess for the result MPS that is close to the actual one. Also it can be helpful (if appropriate) if the MPO is close to the identity, such as when time evolving using a small time step. The number of sweeps can be chosen with the parameter "Nsweep".

It is recommended to try the default "DensityMatrix" first because it is more reliable. Then, the "Fit" method can be tried if higher performance is required.

Named arguments recognized:

`"Method"`

— (default: "DensityMatrix") algorithm used for applying the MPO to the MPS. Currently available options are- "DensityMatrix"
- "Fit"

`"Cutoff"`

— (default: 1E-13) truncation error cutoff for compressing resulting MPS`"MaxDim"`

— maximum bond dimension of resulting compressed MPS`"Verbose"`

— (default: false) if true, prints extra output`"Normalize"`

— (default: false) choose whether or not to normalize the output wavefunction`"Nsweep"`

— (default: 1) sets the number of sweeps of the "Fit" algorithm

Show Example//Use the method "DensityMatrix" auto y1 = applyMPO(A,x,{"Method=","DensityMatrix","MaxDim=",100,"Cutoff=",1E-8}); //Use the method "Fit" with 5 sweeps auto y2 = applyMPO(A,x,{"Method=","Fit","MaxDim=",100,"Cutoff=",1E-8,"Nsweep=",5});

`applyMPO(MPO A, MPS x, MPS x0, Args args = Args::global()) -> MPS`

Similar to

`applyMPO`

above, but accepts a guess for the output wavefunction (the guess wavefunction`x0`

is not overwritten).MPO

`A`

and MPS`x`

must share a set of site indices. The site indices of`x0`

will be made to match the site indices of`A`

that are not shared by`x`

. The links of the output MPS will have the same tags as the links of the guess MPS`x0`

.Currently, this version of

`applyMPO`

only accepts "Fit" for the parameter "Method". Choosing a good guess state`x0`

can improve the convergence of the "Fit" method.Show Exampleauto sites = SiteSet(10,2); // Make trivial MPO and random MPS auto A = MPO(sites); auto x = randomMPS(sites); // Some other random starting state auto x0 = randomMPS(sites); //Use the method "Fit" with 5 sweeps and a guess state x0 auto y = applyMPO(A,x,x0,{"Method=","Fit","MaxDim=",100,"Cutoff=",1E-8,"Nsweep=",2});

`errorMPOProd(MPS y, MPO A, MPS x) -> Real`

Computes, without approximation, the difference @@||\, |y\rangle - A |x\rangle ||^2@@ , where

`A`

is an MPO that shares a set of site indices with MPS`x`

. This is especially useful for testing methods for applying an MPO to an MPS.`A`

and`x`

need to share a set of site indices. The function will attempt to match the remaining site indices of`A`

with the site indices of`y`

.Show Exampleauto sites = SiteSet(10,2); // Make trivial MPO and random MPS auto A = MPO(sites); auto x = randomMPS(sites); //Approximate A|x> auto y = applyMPO(A,x,{"MaxDim=",200,"Cutoff=",1E-12}); //Check Print(errorMPOProd(y,A,x)); //should be close to zero

*This page current as of version 3.0.0*

Back to Classes

Back to Main