quimb.tensor.tensor_1d#

Classes and algorithms related to 1D tensor networks.

Functions

expec_TN_1D(*tns[, compress, eps])

Compute the expectation of several 1D TNs, using transfer matrix compression if any are periodic.

gate_TN_1D(tn, G, where[, contract, tags, ...])

Act with the gate g on sites where, maintaining the outer indices of the 1D tensor netowork.

maybe_factor_gate_into_tensor(G, dp, ng, where)

set_default_compress_mode(opts[, cyclic])

superop_TN_1D(tn_super, tn_op[, ...])

Take a tensor network superoperator and act with it on a tensor network operator, maintaining the original upper and lower indices of the operator.

Classes

Dense1D(array[, phys_dim, tags, ...])

Mimics other 1D tensor network structures, but really just keeps the full state in a single tensor.

MatrixProductOperator(arrays[, shape, ...])

Initialise a matrix product operator, with auto labelling and tagging.

MatrixProductState(arrays, *[, shape, tags, ...])

Initialise a matrix product state, with auto labelling and tagging.

SuperOperator1D(arrays[, shape, ...])

A 1D tensor network super-operator class.

TNLinearOperator1D(*args, **kwargs)

A 1D tensor network linear operator like.

TensorNetwork1D(ts, *[, virtual, ...])

Base class for tensor networks with a one-dimensional structure.

TensorNetwork1DFlat(ts, *[, virtual, ...])

1D Tensor network which has a flat structure.

TensorNetwork1DOperator(ts, *[, virtual, ...])

TensorNetwork1DVector(ts, *[, virtual, ...])

1D Tensor network which overall is like a vector with a single type of site ind.

class quimb.tensor.tensor_1d.Dense1D(array, phys_dim=2, tags=None, site_ind_id='k{}', site_tag_id='I{}', **tn_opts)[source]#

Mimics other 1D tensor network structures, but really just keeps the full state in a single tensor. This allows e.g. applying gates in the same way for quantum circuit simulation as lazily represented hilbert spaces.

Parameters
  • array (array_like) – The full hilbert space vector - assumed to be made of equal hilbert spaces each of size phys_dim and will be reshaped as such.

  • phys_dim (int, optional) – The hilbert space size of each site, default: 2.

  • tags (sequence of str, optional) – Extra tags to add to the tensor network.

  • site_ind_id (str, optional) – String formatter describing how to label the site indices.

  • site_tag_id (str, optional) – String formatter describing how to label the site tags.

  • tn_opts – Supplied to TensorNetwork.

classmethod rand(n, phys_dim=2, dtype=<class 'float'>, **dense1d_opts)[source]#

Create a random dense vector ‘tensor network’.

class quimb.tensor.tensor_1d.MatrixProductOperator(arrays, shape='lrud', site_tag_id='I{}', tags=None, upper_ind_id='k{}', lower_ind_id='b{}', bond_name='', **tn_opts)[source]#

Initialise a matrix product operator, with auto labelling and tagging.

Parameters
  • arrays (sequence of arrays) – The tensor arrays to form into a MPO.

  • shape (str, optional) – String specifying layout of the tensors. E.g. ‘lrud’ (the default) indicates the shape corresponds left-bond, right-bond, ‘up’ physical index, ‘down’ physical index. End tensors have either ‘l’ or ‘r’ dropped from the string.

  • upper_ind_id (str) – A string specifiying how to label the upper physical site indices. Should contain a '{}' placeholder. It is used to generate the actual indices like: map(upper_ind_id.format, range(len(arrays))).

  • lower_ind_id (str) – A string specifiying how to label the lower physical site indices. Should contain a '{}' placeholder. It is used to generate the actual indices like: map(lower_ind_id.format, range(len(arrays))).

  • site_tag_id (str) – A string specifiying how to tag the tensors at each site. Should contain a '{}' placeholder. It is used to generate the actual tags like: map(site_tag_id.format, range(len(arrays))).

  • tags (str or sequence of str, optional) – Global tags to attach to all tensors.

  • bond_name (str, optional) – The base name of the bond indices, onto which uuids will be added.

add_MPO(other, inplace=False, compress=False, **compress_opts)[source]#

Add another MatrixProductState to this one.

apply(other, compress=False, **compress_opts)[source]#

Act with this MPO on another MPO or MPS, such that the resulting object has the same tensor network structure/indices as other.

For an MPS:

       | | | | | | | | | | | | | | | | | |
 self: A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A
       | | | | | | | | | | | | | | | | | |
other: x-x-x-x-x-x-x-x-x-x-x-x-x-x-x-x-x-x

                       -->

       | | | | | | | | | | | | | | | | | |   <- other.site_ind_id
  out: y=y=y=y=y=y=y=y=y=y=y=y=y=y=y=y=y=y

For an MPO:

       | | | | | | | | | | | | | | | | | |
 self: A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A
       | | | | | | | | | | | | | | | | | |
other: B-B-B-B-B-B-B-B-B-B-B-B-B-B-B-B-B-B
       | | | | | | | | | | | | | | | | | |

                       -->

       | | | | | | | | | | | | | | | | | |   <- other.upper_ind_id
  out: C=C=C=C=C=C=C=C=C=C=C=C=C=C=C=C=C=C
       | | | | | | | | | | | | | | | | | |   <- other.lower_ind_id

The resulting TN will have the same structure/indices as other, but probably with larger bonds (depending on compression).

Parameters
Return type

MatrixProductOperator or MatrixProductState

dot(other, compress=False, **compress_opts)#

Act with this MPO on another MPO or MPS, such that the resulting object has the same tensor network structure/indices as other.

For an MPS:

       | | | | | | | | | | | | | | | | | |
 self: A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A
       | | | | | | | | | | | | | | | | | |
other: x-x-x-x-x-x-x-x-x-x-x-x-x-x-x-x-x-x

                       -->

       | | | | | | | | | | | | | | | | | |   <- other.site_ind_id
  out: y=y=y=y=y=y=y=y=y=y=y=y=y=y=y=y=y=y

For an MPO:

       | | | | | | | | | | | | | | | | | |
 self: A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A-A
       | | | | | | | | | | | | | | | | | |
other: B-B-B-B-B-B-B-B-B-B-B-B-B-B-B-B-B-B
       | | | | | | | | | | | | | | | | | |

                       -->

       | | | | | | | | | | | | | | | | | |   <- other.upper_ind_id
  out: C=C=C=C=C=C=C=C=C=C=C=C=C=C=C=C=C=C
       | | | | | | | | | | | | | | | | | |   <- other.lower_ind_id

The resulting TN will have the same structure/indices as other, but probably with larger bonds (depending on compression).

Parameters
Return type

MatrixProductOperator or MatrixProductState

identity(**mpo_opts)[source]#

Get a identity matching this MPO.

property lower_inds#

An ordered tuple of the actual lower physical indices.

partial_transpose(sysa, inplace=False)[source]#

Perform the partial transpose on this MPO by swapping the bra and ket indices on sites in sysa.

Parameters
  • sysa (sequence of int or int) – The sites to transpose indices on.

  • inplace (bool, optional) – Whether to perform the partial transposition inplace.

Return type

MatrixProductOperator

permute_arrays(shape='lrud')[source]#

Permute the indices of each tensor in this MPO to match shape. This doesn’t change how the overall object interacts with other tensor networks but may be useful for extracting the underlying arrays consistently. This is an inplace operation.

Parameters

shape (str, optional) – A permutation of 'lrud' specifying the desired order of the left, right, upper and lower (down) indices respectively.

rand_state(bond_dim, **mps_opts)[source]#

Get a random vector matching this MPO.

trace(left_inds=None, right_inds=None)[source]#

Take the trace of this MPO.

class quimb.tensor.tensor_1d.MatrixProductState(arrays, *, shape='lrp', tags=None, bond_name='', site_ind_id='k{}', site_tag_id='I{}', **tn_opts)[source]#

Initialise a matrix product state, with auto labelling and tagging.

Parameters
  • arrays (sequence of arrays) – The tensor arrays to form into a MPS.

  • shape (str, optional) – String specifying layout of the tensors. E.g. ‘lrp’ (the default) indicates the shape corresponds left-bond, right-bond, physical index. End tensors have either ‘l’ or ‘r’ dropped from the string.

  • site_ind_id (str) – A string specifiying how to label the physical site indices. Should contain a '{}' placeholder. It is used to generate the actual indices like: map(site_ind_id.format, range(len(arrays))).

  • site_tag_id (str) – A string specifiying how to tag the tensors at each site. Should contain a '{}' placeholder. It is used to generate the actual tags like: map(site_tag_id.format, range(len(arrays))).

  • tags (str or sequence of str, optional) – Global tags to attach to all tensors.

  • bond_name (str, optional) – The base name of the bond indices, onto which uuids will be added.

add_MPS(other, inplace=False, compress=False, **compress_opts)[source]#

Add another MatrixProductState to this one.

bipartite_schmidt_state(sz_a, get='ket', cur_orthog=None)[source]#

Compute the reduced state for a bipartition of an OBC MPS, in terms of the minimal left/right schmidt basis:

     A            B
 .........     ...........
 >->->->->--s--<-<-<-<-<-<    ->   +-s-+
 | | | | |     | | | | | |         |   |
k0 k1...                          kA   kB
Parameters
  • sz_a (int) – The number of sites in subsystem A, must be 0 < sz_a < N.

  • get ({'ket', 'rho', 'ket-dense', 'rho-dense'}, optional) –

    Get the:

    • ’ket’: vector form as tensor.

    • ’rho’: density operator form, i.e. vector outer product

    • ’ket-dense’: like ‘ket’ but return qarray.

    • ’rho-dense’: like ‘rho’ but return qarray.

  • cur_orthog (int, optional) – If given, take as the current orthogonality center so as to efficienctly move it a minimal distance.

entropy(i, cur_orthog=None, method='svd')[source]#

The entropy of bipartition between the left block of i sites and the rest.

Parameters
  • i (int) – The number of sites in the left partition.

  • cur_orthog (int) – If given, the known current orthogonality center, to speed up the mixed canonization.

Return type

float

classmethod from_dense(psi, dims, site_ind_id='k{}', site_tag_id='I{}', **split_opts)[source]#

Create a MatrixProductState directly from a dense vector

Parameters
  • psi (array_like) – The dense state to convert to MPS from.

  • dims (sequence of int) – Physical subsystem dimensions of each site.

  • site_ind_id (str, optional) – How to index the physical sites, see MatrixProductState.

  • site_tag_id (str, optional) – How to tag the physical sites, see MatrixProductState.

  • split_opts – Supplied to tensor_split() to in order to partition the dense vector into tensors.

Return type

MatrixProductState

Examples

>>> dims = [2, 2, 2, 2, 2, 2]
>>> psi = rand_ket(prod(dims))
>>> mps = MatrixProductState.from_dense(psi, dims)
>>> mps.show()
 2 4 8 4 2
o-o-o-o-o-o
| | | | | |
gate_split(G, where, inplace=False, **compress_opts)[source]#

Apply a two-site gate and then split resulting tensor to retrieve a MPS form:

-o-o-A-B-o-o-
 | | | | | |            -o-o-GGG-o-o-           -o-o-X~Y-o-o-
 | | GGG | |     ==>     | | | | | |     ==>     | | | | | |
 | | | | | |                 i j                     i j
     i j

As might be found in TEBD.

Parameters
  • G (array) – The gate, with shape (d**2, d**2) for physical dimension d.

  • where ((int, int)) – Indices of the sites to apply the gate to.

  • compress_opts – Supplied to tensor_split().

See also

gate, gate_with_auto_swap

gate_with_auto_swap(G, where, inplace=False, cur_orthog=None, **compress_opts)[source]#

Perform a two site gate on this MPS by, if necessary, swapping and compressing the sites until they are adjacent, using gate_split, then unswapping the sites back to their original position.

Parameters
  • G (array) – The gate, with shape (d**2, d**2) for physical dimension d.

  • where ((int, int)) – Indices of the sites to apply the gate to.

  • cur_orthog (int, sequence of int, or 'calc') – If known, the current orthogonality center.

  • inplace (bond, optional) – Perform the swaps inplace.

  • compress_opts – Supplied to tensor_split().

See also

gate, gate_split

logneg_subsys(sysa, sysb, compress_opts=None, approx_spectral_opts=None, verbosity=0, approx_thresh=4096)[source]#

Compute the logarithmic negativity between subsytem blocks, e.g.:

                   sysa         sysb
                 .........       .....
... -o-o-o-o-o-o-A-A-A-A-A-o-o-o-B-B-B-o-o-o-o-o-o-o- ...
     | | | | | | | | | | | | | | | | | | | | | | | |
Parameters
  • sysa (sequence of int) – The sites, which should be contiguous, defining subsystem A.

  • sysb (sequence of int) – The sites, which should be contiguous, defining subsystem B.

  • eps (float, optional) – Tolerance to use when compressing the subsystem transfer matrices.

  • method (str or (str, str), optional) – Method(s) to use for laterally compressing the state then vertially compressing subsytems.

  • compress_opts (dict, optional) – If given, supplied to partial_trace_compress to govern how singular values are treated. See tensor_split.

  • approx_spectral_opts – Supplied to approx_spectral_function().

Returns

ln – The logarithmic negativity.

Return type

float

See also

MatrixProductState.partial_trace_compress, approx_spectral_function

magnetization(i, direction='Z', cur_orthog=None)[source]#

Compute the magnetization at site i.

measure(site, remove=False, outcome=None, renorm=True, cur_orthog=None, get=None, inplace=False)[source]#

Measure this MPS at site, including projecting the state. Optionally remove the site afterwards, yielding an MPS with one less site. In either case the orthogonality center of the returned MPS is min(site, new_L - 1).

Parameters
  • site (int) – The site to measure.

  • remove (bool, optional) –

    Whether to remove the site completely after projecting the measurement. If True, sites greater than site will be retagged and reindex one down, and the MPS will have one less site. E.g:

    0-1-2-3-4-5-6
           / / /  - measure and remove site 3
    0-1-2-4-5-6
                  - reindex sites (4, 5, 6) to (3, 4, 5)
    0-1-2-3-4-5
    

  • outcome (None or int, optional) – Specify the desired outcome of the measurement. If None, it will be randomly sampled according to the local density matrix.

  • renorm (bool, optional) – Whether to renormalize the state post measurement.

  • cur_orthog (None or int, optional) – If you already know the orthogonality center, you can supply it here for efficiencies sake.

  • get ({None, 'outcome'}, optional) – If 'outcome', simply return the outcome, and don’t perform any projection.

  • inplace (bool, optional) – Whether to perform the measurement in place or not.

Returns

  • outcome (int) – The measurement outcome, drawn from range(phys_dim).

  • psi (MatrixProductState) – The measured state, if get != 'outcome'.

normalize(bra=None, eps=1e-15, insert=None)[source]#

Normalize this MPS, optional with co-vector bra. For periodic MPS this uses transfer matrix SVD approximation with precision eps in order to be efficient. Inplace.

Parameters
  • bra (MatrixProductState, optional) – If given, normalize this MPS with the same factor.

  • eps (float, optional) – If cyclic, precision to approximation transfer matrix with. Default: 1e-14.

  • insert (int, optional) – Insert the corrective normalization on this site, random if not given.

Returns

old_norm – The old norm self.H @ self.

Return type

float

partial_trace(keep, upper_ind_id='b{}', rescale_sites=True)[source]#

Partially trace this matrix product state, producing a matrix product operator.

Parameters
  • keep (sequence of int or slice) – Indicies of the sites to keep.

  • upper_ind_id (str, optional) – The ind id of the (new) ‘upper’ inds, i.e. the ‘bra’ inds.

  • rescale_sites (bool, optional) – If True (the default), then the kept sites will be rescaled to (0, 1, 2, ...) etc. rather than keeping their original site numbers.

Returns

rho – The density operator in MPO form.

Return type

MatrixProductOperator

partial_trace_compress(sysa, sysb, eps=1e-08, method=('isvd', None), max_bond=(None, 1024), leave_short=True, renorm=True, lower_ind_id='b{}', verbosity=0, **compress_opts)[source]#

Perform a compressed partial trace using singular value lateral then vertical decompositions of transfer matrix products:

        .....sysa......     ...sysb....
o-o-o-o-A-A-A-A-A-A-A-A-o-o-B-B-B-B-B-B-o-o-o-o-o-o-o-o-o
| | | | | | | | | | | | | | | | | | | | | | | | | | | | |

                          ==> form inner product

        ...............     ...........
o-o-o-o-A-A-A-A-A-A-A-A-o-o-B-B-B-B-B-B-o-o-o-o-o-o-o-o-o
| | | | | | | | | | | | | | | | | | | | | | | | | | | | |
o-o-o-o-A-A-A-A-A-A-A-A-o-o-B-B-B-B-B-B-o-o-o-o-o-o-o-o-o

                          ==> lateral SVD on each section

          .....sysa......     ...sysb....
          /\             /\   /\         /\
  ... ~~~E  A~~~~~~~~~~~A  E~E  B~~~~~~~B  E~~~ ...
          \/             \/   \/         \/

                          ==> vertical SVD and unfold on A & B

                  |                 |
          /-------A-------\   /-----B-----\
  ... ~~~E                 E~E             E~~~ ...
          \-------A-------/   \-----B-----/
                  |                 |

With various special cases including OBC or end spins included in subsytems.

Parameters
  • sysa (sequence of int) – The sites, which should be contiguous, defining subsystem A.

  • sysb (sequence of int) – The sites, which should be contiguous, defining subsystem B.

  • eps (float or (float, float), optional) – Tolerance(s) to use when compressing the subsystem transfer matrices and vertically decomposing.

  • method (str or (str, str), optional) – Method(s) to use for laterally compressing the state then vertially compressing subsytems.

  • max_bond (int or (int, int), optional) – The maximum bond to keep for laterally compressing the state then vertially compressing subsytems.

  • leave_short (bool, optional) – If True (the default), don’t try to compress short sections.

  • renorm (bool, optional) – If True (the default), renomalize the state so that tr(rho)==1.

  • lower_ind_id (str, optional) – The index id to create for the new density matrix, the upper_ind_id is automatically taken as the current site_ind_id.

  • compress_opts (dict, optional) – If given, supplied to partial_trace_compress to govern how singular values are treated. See tensor_split.

  • verbosity ({0, 1}, optional) – How much information to print while performing the compressed partial trace.

Returns

rho_ab – Density matrix tensor network with outer_inds = ('k0', 'k1', 'b0', 'b1') for example.

Return type

TensorNetwork

permute_arrays(shape='lrp')[source]#

Permute the indices of each tensor in this MPS to match shape. This doesn’t change how the overall object interacts with other tensor networks but may be useful for extracting the underlying arrays consistently. This is an inplace operation.

Parameters

shape (str, optional) – A permutation of 'lrp' specifying the desired order of the left, right, and physical indices respectively.

ptr(keep, upper_ind_id='b{}', rescale_sites=True)[source]#

Alias of partial_trace().

schmidt_gap(i, cur_orthog=None, method='svd')[source]#

The schmidt gap of bipartition between the left block of i sites and the rest.

Parameters
  • i (int) – The number of sites in the left partition.

  • cur_orthog (int) – If given, the known current orthogonality center, to speed up the mixed canonization.

Return type

float

schmidt_values(i, cur_orthog=None, method='svd')[source]#

Find the schmidt values associated with the bipartition of this MPS between sites on either site of i. In other words, i is the number of sites in the left hand partition:

....L....   i
o-o-o-o-o-S-o-o-o-o-o-o-o-o-o-o-o
| | | | |   | | | | | | | | | | |
       i-1  ..........R..........

The schmidt values, S, are the singular values associated with the (i - 1, i) bond, squared, provided the MPS is mixed canonized at one of those sites.

Parameters
  • i (int) – The number of sites in the left partition.

  • cur_orthog (int) – If given, the known current orthogonality center, to speed up the mixed canonization.

Returns

S – The schmidt values.

Return type

1d-array

swap_site_to(i, f, cur_orthog=None, inplace=False, **compress_opts)[source]#

Swap site i to site f, compressing the bond after each swap:

      i       f
0 1 2 3 4 5 6 7 8 9      0 1 2 4 5 6 7 3 8 9
o-o-o-x-o-o-o-o-o-o      o-o-o-o-o-o-o-x-o-o
| | | | | | | | | |  ->  | | | | | | | | | |
Parameters
  • i (int) – The site to move.

  • f (int) – The new location for site i.

  • cur_orthog (int, sequence of int, or 'calc') – If known, the current orthogonality center.

  • inplace (bond, optional) – Perform the swaps inplace.

  • compress_opts – Supplied to tensor_split().

swap_sites_with_compress(i, j, cur_orthog=None, inplace=False, **compress_opts)[source]#

Swap sites i and j by contracting, then splitting with the physical indices swapped.

Parameters
  • i (int) – The first site to swap.

  • j (int) – The second site to swap.

  • cur_orthog (int, sequence of int, or 'calc') – If known, the current orthogonality center.

  • inplace (bond, optional) – Perform the swaps inplace.

  • compress_opts – Supplied to tensor_split().

class quimb.tensor.tensor_1d.SuperOperator1D(arrays, shape='lrkud', site_tag_id='I{}', outer_upper_ind_id='kn{}', inner_upper_ind_id='k{}', inner_lower_ind_id='b{}', outer_lower_ind_id='bn{}', tags=None, tags_upper=None, tags_lower=None, **tn_opts)[source]#

A 1D tensor network super-operator class:

0   1   2       n-1
|   |   |        |     <-- outer_upper_ind_id
O===O===O==     =O
|\  |\  |\       |\     <-- inner_upper_ind_id
  )   )   ) ...    )   <-- K (size of local Kraus sum)
|/  |/  |/       |/     <-- inner_lower_ind_id
O===O===O==     =O
|   | : |        |     <-- outer_lower_ind_id
      :
     chi (size of entangling bond dim)
Parameters

arrays (sequence of arrays) – The data arrays defining the superoperator, this should be a sequence of 2n arrays, such that the first two correspond to the upper and lower operators acting on site 0 etc. The arrays should be 5 dimensional unless OBC conditions are desired, in which case the first two and last two should be 4-dimensional. The dimensions of array can be should match the shape option.

class quimb.tensor.tensor_1d.TNLinearOperator1D(*args, **kwargs)[source]#

A 1D tensor network linear operator like:

         start                 stop - 1
           .                     .
         :-O-O-O-O-O-O-O-O-O-O-O-O-:                 --+
         : | | | | | | | | | | | | :                   |
         :-H-H-H-H-H-H-H-H-H-H-H-H-:    acting on    --V
         : | | | | | | | | | | | | :                   |
         :-O-O-O-O-O-O-O-O-O-O-O-O-:                 --+
left_inds^                         ^right_inds

Like TNLinearOperator, but performs a structured contract from one end to the other than can handle very long chains possibly more efficiently by contracting in blocks from one end.

Parameters
  • tn (TensorNetwork) – The tensor network to turn into a LinearOperator.

  • left_inds (sequence of str) – The left indicies.

  • right_inds (sequence of str) – The right indicies.

  • start (int) – Index of starting site.

  • stop (int) – Index of stopping site (does not include this site).

  • ldims (tuple of int, optional) – If known, the dimensions corresponding to left_inds.

  • rdims (tuple of int, optional) – If known, the dimensions corresponding to right_inds.

See also

TNLinearOperator

class quimb.tensor.tensor_1d.TensorNetwork1D(ts, *, virtual=False, check_collisions=True)[source]#

Base class for tensor networks with a one-dimensional structure.

property L#

The number of sites.

align(*, ind_ids=None, inplace=False)#

Align an arbitrary number of tensor networks in a stack-like geometry:

a-a-a-a-a-a-a-a-a-a-a-a-a-a-a-a-a-a
| | | | | | | | | | | | | | | | | | <- ind_ids[0] (defaults to 1st id)
b-b-b-b-b-b-b-b-b-b-b-b-b-b-b-b-b-b
| | | | | | | | | | | | | | | | | | <- ind_ids[1]
               ...
| | | | | | | | | | | | | | | | | | <- ind_ids[-2]
y-y-y-y-y-y-y-y-y-y-y-y-y-y-y-y-y-y
| | | | | | | | | | | | | | | | | | <- ind_ids[-1]
z-z-z-z-z-z-z-z-z-z-z-z-z-z-z-z-z-z
Parameters
  • tns (sequence of TensorNetwork) – The TNs to align, should be structured and either effective ‘vectors’ (have a site_ind_id) or ‘operators’ (have a up_ind_id and lower_ind_id).

  • ind_ids (None, or sequence of str) – String with format specifiers to id each level of sites with. Will be automatically generated like (tns[0].site_ind_id, "__ind_a{}__", "__ind_b{}__", ...) if not given.

  • inplace (bool) – Whether to modify the input tensor networks inplace.

Returns

tns_aligned

Return type

sequence of TensorNetwork

contract_structured(tag_slice, structure_bsz=5, inplace=False, **opts)[source]#

Perform a structured contraction, translating tag_slice from a slice or to a cumulative sequence of tags.

Parameters
  • tag_slice (slice or ...) – The range of sites, or for all.

  • inplace (bool, optional) – Whether to perform the contraction inplace.

Returns

The result of the contraction, still a TensorNetwork if the contraction was only partial.

Return type

TensorNetwork, Tensor or scalar

See also

contract, contract_tags, contract_cumulative

maybe_convert_coo(x)[source]#

Check if x is an integer and convert to the corresponding site tag if so.

property nsites#

The number of sites.

retag_sites(new_id, where=None, inplace=False)[source]#

Modify the site tags for all or some tensors in this 1D TN (without changing the site_tag_id).

site_tag(i)[source]#

The name of the tag specifiying the tensor at site i.

property site_tag_id#

The string specifier for tagging each site of this 1D TN.

property site_tags#

An ordered tuple of the actual site tags.

slice2sites(tag_slice)[source]#

Take a slice object, and work out its implied start, stop and step, taking into account cyclic boundary conditions.

Examples

Normal slicing:

>>> p = MPS_rand_state(10, bond_dim=7)
>>> p.slice2sites(slice(5))
(0, 1, 2, 3, 4)
>>> p.slice2sites(slice(4, 8))
(4, 5, 6, 7)

Slicing from end backwards:

>>> p.slice2sites(slice(..., -3, -1))
(9, 8)

Slicing round the end:

>>> p.slice2sites(slice(7, 12))
(7, 8, 9, 0, 1)
>>> p.slice2sites(slice(-3, 2))
(7, 8, 9, 0, 1)

If the start point is > end point (before modulo n), then step needs to be negative to return anything.

class quimb.tensor.tensor_1d.TensorNetwork1DFlat(ts, *, virtual=False, check_collisions=True)[source]#

1D Tensor network which has a flat structure.

amplitude(b)[source]#

Compute the amplitude of configuration b.

Parameters

b (sequence of int) – The configuration to compute the amplitude of.

Returns

c_b

Return type

scalar

as_cyclic(inplace=False)[source]#

Convert this flat, 1D, TN into cyclic form by adding a dummy bond between the first and last sites.

bond(i, j)[source]#

Get the name of the index defining the bond between sites i and j.

bond_size(i, j)[source]#

Return the size of the bond between site i and j.

calc_current_orthog_center()[source]#

Calculate the site(s) of the current orthogonality center.

Returns

The site, or min/max, around which this MPS is orthogonal.

Return type

int or (int, int)

canonize(where, cur_orthog='calc', bra=None)[source]#

Mixed canonize this TN. If this is a MPS, this implies that:

              i                      i
>->->->->- ->-o-<- -<-<-<-<-<      +-o-+
| | | | |...| | |...| | | | |  ->  | | |
>->->->->- ->-o-<- -<-<-<-<-<      +-o-+

You can also supply a set of indices to orthogonalize around, and a current location of the orthogonality center for efficiency:

      current                              where
      .......                              .....
>->->-c-c-c-c-<-<-<-<-<-<      >->->->->->-w-w-w-<-<-<-<
| | | | | | | | | | | | |  ->  | | | | | | | | | | | | |
>->->-c-c-c-c-<-<-<-<-<-<      >->->->->->-w-w-w-<-<-<-<
   cmin     cmax                           i   j

This would only move cmin to i and cmax to j if necessary.

Parameters
  • where (int or sequence of int) – Which site(s) to orthogonalize around. If a sequence of int then make sure that section from min(where) to max(where) is orthog.

  • cur_orthog (int, sequence of int, or 'calc') – If given, the current site(s), so as to shift the orthogonality ceneter as efficiently as possible. If ‘calc’, calculate the current orthogonality center.

  • bra (MatrixProductState, optional) – If supplied, simultaneously mixed canonize this MPS too, assuming it to be the conjugate state.

canonize_cyclic(i, bra=None, method='isvd', inv_tol=1e-10)[source]#

Bring this MatrixProductState into (possibly only approximate) canonical form at site(s) i.

Parameters
  • i (int or slice) – The site or range of sites to make canonical.

  • bra (MatrixProductState, optional) – Simultaneously canonize this state as well, assuming it to be the co-vector.

  • method ({'isvd', 'svds', ...}, optional) – How to perform the lateral compression.

  • inv_tol (float, optional) – Tolerance with which to invert the gauge.

compress(form=None, **compress_opts)[source]#

Compress this 1D Tensor Network, possibly into canonical form.

Parameters
  • form ({None, 'flat', 'left', 'right'} or int) – Output form of the TN. None left canonizes the state first for stability reasons, then right_compresses (default). 'flat' tries to distribute the singular values evenly – state will not be canonical. 'left' and 'right' put the state into left and right canonical form respectively with a prior opposite sweep, or an int will put the state into mixed canonical form at that site.

  • compress_opts – Supplied to Tensor.split().

compress_site(i, canonize=True, cur_orthog='calc', bra=None, **compress_opts)[source]#

Compress the bonds adjacent to site i, by default first setting the orthogonality center to that site:

     i                     i
-o-o-o-o-o-    -->    ->->~o~<-<-
 | | | | |             | | | | |
Parameters
  • i (int) – Which site to compress around

  • canonize (bool, optional) – Whether to first set the orthogonality center to site i.

  • cur_orthog (int, optional) – If given, the known current orthogonality center, to speed up the mixed canonization.

  • bra (MatrixProductState, optional) – The conjugate state to also apply the compression to.

  • compress_opts – Supplied to tensor_split().

expand_bond_dimension(new_bond_dim, rand_strength=0.0, bra=None, inplace=True)[source]#

Expand the bond dimensions of this 1D tensor network to at least new_bond_dim.

Parameters
  • new_bond_dim (int) – Minimum bond dimension to expand to.

  • inplace (bool, optional) – Whether to perform the expansion in place.

  • bra (MatrixProductState, optional) – Mirror the changes to bra inplace, treating it as the conjugate state.

  • rand_strength (float, optional) – If rand_strength > 0, fill the new tensor entries with gaussian noise of strength rand_strength.

Return type

MatrixProductState

left_canonize(stop=None, start=None, normalize=False, bra=None)[source]#

Left canonize all or a portion of this TN. If this is a MPS, this implies that:

              i              i
>->->->->->->-o-o-         +-o-o-
| | | | | | | | | ...  =>  | | | ...
>->->->->->->-o-o-         +-o-o-
Parameters
  • start (int, optional) – If given, the site to start left canonizing at.

  • stop (int, optional) – If given, the site to stop left canonizing at.

  • normalize (bool, optional) – Whether to normalize the state, only works for OBC.

  • bra (MatrixProductState, optional) – If supplied, simultaneously left canonize this MPS too, assuming it to be the conjugate state.

left_canonize_site(i, bra=None)[source]#

Left canonize this TN’s ith site, inplace:

    i                i
   -o-o-            ->-s-
... | | ...  ==> ... | | ...
Parameters
  • i (int) – Which site to canonize. The site at i + 1 also absorbs the non-isometric part of the decomposition of site i.

  • bra (None or matching TensorNetwork to self, optional) – If set, also update this TN’s data with the conjugate canonization.

left_compress(start=None, stop=None, bra=None, **compress_opts)[source]#

Compress this 1D TN, from left to right, such that it becomes left-canonical (unless absorb != 'right').

Parameters
  • start (int, optional) – Site to begin compressing on.

  • stop (int, optional) – Site to stop compressing at (won’t itself be an isometry).

  • bra (None or TensorNetwork like this one, optional) – If given, update this TN as well, assuming it to be the conjugate.

  • compress_opts – Supplied to Tensor.split().

left_compress_site(i, bra=None, **compress_opts)[source]#

Left compress this 1D TN’s ith site, such that the site is then left unitary with its right bond (possibly) reduced in dimension.

Parameters
  • i (int) – Which site to compress.

  • bra (None or matching TensorNetwork to self, optional) – If set, also update this TN’s data with the conjugate compression.

  • compress_opts – Supplied to Tensor.split().

right_canonize(stop=None, start=None, normalize=False, bra=None)[source]#

Right canonize all or a portion of this TN. If this is a MPS, this implies that:

      i                           i
   -o-o-<-<-<-<-<-<-<          -o-o-+
... | | | | | | | | |   ->  ... | | |
   -o-o-<-<-<-<-<-<-<          -o-o-+
Parameters
  • start (int, optional) – If given, the site to start right canonizing at.

  • stop (int, optional) – If given, the site to stop right canonizing at.

  • normalize (bool, optional) – Whether to normalize the state.

  • bra (MatrixProductState, optional) – If supplied, simultaneously right canonize this MPS too, assuming it to be the conjugate state.

right_canonize_site(i, bra=None)[source]#

Right canonize this TN’s ith site, inplace:

      i                i
   -o-o-            -s-<-
... | | ...  ==> ... | | ...
Parameters

i (int) –

Which site to canonize. The site at i - 1 also absorbs the

non-isometric part of the decomposition of site i.

braNone or matching TensorNetwork to self, optional

If set, also update this TN’s data with the conjugate canonization.

right_compress(start=None, stop=None, bra=None, **compress_opts)[source]#

Compress this 1D TN, from right to left, such that it becomes right-canonical (unless absorb != 'left').

Parameters
  • start (int, optional) – Site to begin compressing on.

  • stop (int, optional) – Site to stop compressing at (won’t itself be an isometry).

  • bra (None or TensorNetwork like this one, optional) – If given, update this TN as well, assuming it to be the conjugate.

  • compress_opts – Supplied to Tensor.split().

right_compress_site(i, bra=None, **compress_opts)[source]#

Right compress this 1D TN’s ith site, such that the site is then right unitary with its left bond (possibly) reduced in dimension.

Parameters
  • i (int) – Which site to compress.

  • bra (None or matching TensorNetwork to self, optional) – If set, update this TN’s data with the conjugate compression.

  • compress_opts – Supplied to Tensor.split().

shift_orthogonality_center(current, new, bra=None)[source]#

Move the orthogonality center of this MPS.

Parameters
  • current (int) – The current orthogonality center.

  • new (int) – The target orthogonality center.

  • bra (MatrixProductState, optional) – If supplied, simultaneously move the orthogonality center of this MPS too, assuming it to be the conjugate state.

singular_values(i, cur_orthog=None, method='svd')[source]#

Find the singular values associated with the ith bond:

....L....   i
o-o-o-o-o-l-o-o-o-o-o-o-o-o-o-o-o
| | | | |   | | | | | | | | | | |
       i-1  ..........R..........

Leaves the 1D TN in mixed canoncial form at bond i.

Parameters
  • i (int) – Which bond, or equivalently, the number of sites in the left partition.

  • cur_orthog (int) – If given, the known current orthogonality center, to speed up the mixed canonization, e.g. if sweeping this function from left to right would use i - 1.

Returns

svals – The singular values.

Return type

1d-array

class quimb.tensor.tensor_1d.TensorNetwork1DOperator(ts, *, virtual=False, check_collisions=True)[source]#
lower_ind(i)[source]#

The name of the lower (‘ket’) index at site i.

property lower_ind_id#

The string specifier for the lower phyiscal indices

property lower_inds#

An ordered tuple of the actual lower physical indices.

phys_dim(i=None, which='upper')[source]#

Get a physical index size of this 1D operator.

reindex_lower_sites(new_id, where=None, inplace=False)[source]#

Update the lower site index labels to a new string specifier.

Parameters
  • new_id (str) – A string with a format placeholder to accept an int, e.g. "ket{}".

  • where (None or slice) – Which sites to update the index labels on. If None (default) all sites.

  • inplace (bool) – Whether to reindex in place.

reindex_upper_sites(new_id, where=None, inplace=False)[source]#

Update the upper site index labels to a new string specifier.

Parameters
  • new_id (str) – A string with a format placeholder to accept an int, e.g. “ket{}”.

  • where (None or slice) – Which sites to update the index labels on. If None (default) all sites.

  • inplace (bool) – Whether to reindex in place.

to_dense(*inds_seq, **contract_opts)[source]#

Return the dense matrix version of this 1D operator, i.e. a qarray with shape (d, d).

upper_ind(i)[source]#

The name of the upper (‘bra’) index at site i.

property upper_ind_id#

The string specifier for the upper phyiscal indices

property upper_inds#

An ordered tuple of the actual upper physical indices.

class quimb.tensor.tensor_1d.TensorNetwork1DVector(ts, *, virtual=False, check_collisions=True)[source]#

1D Tensor network which overall is like a vector with a single type of site ind.

correlation(A, i, j, B=None, **expec_opts)[source]#

Correlation of operator A between i and j.

Parameters
  • A (array) – The operator to act with, can be multi site.

  • i (int or sequence of int) – The first site(s).

  • j (int or sequence of int) – The second site(s).

  • expec_opts – Supplied to expec_TN_1D().

Returns

C – The correlation <A(i)> + <A(j)> - <A(ij)>.

Return type

float

Examples

>>> ghz = (MPS_computational_state('0000') +
...        MPS_computational_state('1111')) / 2**0.5
>>> ghz.correlation(pauli('Z'), 0, 1)
1.0
>>> ghz.correlation(pauli('Z'), 0, 1, B=pauli('X'))
0.0
expec(*, compress=None, eps=1e-15)[source]#

Compute the expectation of several 1D TNs, using transfer matrix compression if any are periodic.

Parameters
  • tns (sequence of TensorNetwork1D) – The MPS and MPO to find expectation of. Should start and begin with an MPS e.g. (MPS, MPO, ...,  MPS).

  • compress ({None, False, True}, optional) – Whether to perform transfer matrix compression on cyclic systems. If set to None (the default), decide heuristically.

  • eps (float, optional) – The accuracy of the transfer matrix compression.

Returns

x – The expectation value.

Return type

float

gate(G, where, contract=False, tags=None, propagate_tags='sites', inplace=False, cur_orthog=None, **compress_opts)[source]#

Act with the gate g on sites where, maintaining the outer indices of the 1D tensor netowork:

contract=False       contract=True
    . .                    . .             <- where
o-o-o-o-o-o-o        o-o-o-GGG-o-o-o
| | | | | | |        | | | / \ | | |
    GGG
    | |


contract='split-gate'        contract='swap-split-gate'
      . .                          . .                      <- where
  o-o-o-o-o-o-o                o-o-o-o-o-o-o
  | | | | | | |                | | | | | | |
      G~G                          G~G
      | |                          \ /
                                    X
                                   / \

contract='swap+split'
        . .            <- where
  o-o-o-G=G-o-o-o
  | | | | | | | |

Note that the sites in where do not have to be contiguous. By default, site tags will be propagated to the gate tensors, identifying a ‘light cone’.

Parameters
  • tn (TensorNetwork1DVector) – The 1D vector-like tensor network, for example, and MPS.

  • G (array) – A square array to act with on sites where. It should have twice the number of dimensions as the number of sites. The second half of these will be contracted with the MPS, and the first half indexed with the correct site_ind_id. Sites are read left to right from the shape. A two-dimensional array is permissible if each dimension factorizes correctly.

  • where (int or sequence of int) – Where the gate should act.

  • contract ({False, 'split-gate', 'swap-split-gate',) –

    ‘auto-split-gate’, True, ‘swap+split’}, optional Whether to contract the gate into the 1D tensor network. If,

    • False: leave the gate uncontracted, the default

    • ’split-gate’: like False, but split the gate if it is two-site.

    • ’swap-split-gate’: like ‘split-gate’, but decompose the gate as if a swap had first been applied

    • ’auto-split-gate’: automatically select between the above three options, based on the rank of the gate.

    • True: contract the gate into the tensor network, if the gate acts on more than one site, this will produce an ever larger tensor.

    • ’swap+split’: Swap sites until they are adjacent, then contract the gate and split the resulting tensor, then swap the sites back to their original position. In this way an MPS structure can be explicitly maintained at the cost of rising bond-dimension.

  • tags (str or sequence of str, optional) – Tag the new gate tensor with these tags.

  • propagate_tags ({'sites', 'register', False, True}, optional) –

    Add any tags from the sites to the new gate tensor (only matters if contract=False else tags are merged anyway):

    • If 'sites', then only propagate tags matching e.g. ‘I{}’ and ignore all others. I.e. just propagate the lightcone.

    • If 'register', then only propagate tags matching the sites of where this gate was actually applied. I.e. ignore the lightcone, just keep track of which ‘registers’ the gate was applied to.

    • If False, propagate nothing.

    • If True, propagate all tags.

  • inplace – Perform the gate in place.

  • bool – Perform the gate in place.

  • optional – Perform the gate in place.

  • compress_opts – Supplied to split() if contract='swap+split' or gate_with_auto_swap() if contract='swap+split'.

Return type

TensorNetwork1DVector

Examples

>>> p = MPS_rand_state(3, 7)
>>> p.gate_(spin_operator('X'), where=1, tags=['GX'])
>>> p
<MatrixProductState(tensors=4, L=3, max_bond=7)>
>>> p.outer_inds()
('k0', 'k1', 'k2')
reindex_all(new_id, inplace=False)[source]#

Reindex all physical sites and change the site_ind_id.

reindex_sites(new_id, where=None, inplace=False)[source]#

Update the physical site index labels to a new string specifier. Note that this doesn’t change the stored id string with the TN.

Parameters
  • new_id (str) – A string with a format placeholder to accept an int, e.g. “ket{}”.

  • where (None or slice) – Which sites to update the index labels on. If None (default) all sites.

  • inplace (bool) – Whether to reindex in place.

site_ind(i)[source]#

Get the physical index name of site i.

property site_ind_id#

The string specifier for the physical indices

property site_inds#

An ordered tuple of the actual physical indices.

to_dense(*inds_seq, **contract_opts)[source]#

Return the dense ket version of this 1D vector, i.e. a qarray with shape (-1, 1).

quimb.tensor.tensor_1d.expec_TN_1D(*tns, compress=None, eps=1e-15)[source]#

Compute the expectation of several 1D TNs, using transfer matrix compression if any are periodic.

Parameters
  • tns (sequence of TensorNetwork1D) – The MPS and MPO to find expectation of. Should start and begin with an MPS e.g. (MPS, MPO, ...,  MPS).

  • compress ({None, False, True}, optional) – Whether to perform transfer matrix compression on cyclic systems. If set to None (the default), decide heuristically.

  • eps (float, optional) – The accuracy of the transfer matrix compression.

Returns

x – The expectation value.

Return type

float

quimb.tensor.tensor_1d.gate_TN_1D(tn, G, where, contract=False, tags=None, propagate_tags='sites', inplace=False, cur_orthog=None, **compress_opts)[source]#

Act with the gate g on sites where, maintaining the outer indices of the 1D tensor netowork:

contract=False       contract=True
    . .                    . .             <- where
o-o-o-o-o-o-o        o-o-o-GGG-o-o-o
| | | | | | |        | | | / \ | | |
    GGG
    | |


contract='split-gate'        contract='swap-split-gate'
      . .                          . .                      <- where
  o-o-o-o-o-o-o                o-o-o-o-o-o-o
  | | | | | | |                | | | | | | |
      G~G                          G~G
      | |                          \ /
                                    X
                                   / \

contract='swap+split'
        . .            <- where
  o-o-o-G=G-o-o-o
  | | | | | | | |

Note that the sites in where do not have to be contiguous. By default, site tags will be propagated to the gate tensors, identifying a ‘light cone’.

Parameters
  • tn (TensorNetwork1DVector) – The 1D vector-like tensor network, for example, and MPS.

  • G (array) – A square array to act with on sites where. It should have twice the number of dimensions as the number of sites. The second half of these will be contracted with the MPS, and the first half indexed with the correct site_ind_id. Sites are read left to right from the shape. A two-dimensional array is permissible if each dimension factorizes correctly.

  • where (int or sequence of int) – Where the gate should act.

  • contract ({False, 'split-gate', 'swap-split-gate',) –

    ‘auto-split-gate’, True, ‘swap+split’}, optional Whether to contract the gate into the 1D tensor network. If,

    • False: leave the gate uncontracted, the default

    • ’split-gate’: like False, but split the gate if it is two-site.

    • ’swap-split-gate’: like ‘split-gate’, but decompose the gate as if a swap had first been applied

    • ’auto-split-gate’: automatically select between the above three options, based on the rank of the gate.

    • True: contract the gate into the tensor network, if the gate acts on more than one site, this will produce an ever larger tensor.

    • ’swap+split’: Swap sites until they are adjacent, then contract the gate and split the resulting tensor, then swap the sites back to their original position. In this way an MPS structure can be explicitly maintained at the cost of rising bond-dimension.

  • tags (str or sequence of str, optional) – Tag the new gate tensor with these tags.

  • propagate_tags ({'sites', 'register', False, True}, optional) –

    Add any tags from the sites to the new gate tensor (only matters if contract=False else tags are merged anyway):

    • If 'sites', then only propagate tags matching e.g. ‘I{}’ and ignore all others. I.e. just propagate the lightcone.

    • If 'register', then only propagate tags matching the sites of where this gate was actually applied. I.e. ignore the lightcone, just keep track of which ‘registers’ the gate was applied to.

    • If False, propagate nothing.

    • If True, propagate all tags.

  • inplace – Perform the gate in place.

  • bool – Perform the gate in place.

  • optional – Perform the gate in place.

  • compress_opts – Supplied to split() if contract='swap+split' or gate_with_auto_swap() if contract='swap+split'.

Return type

TensorNetwork1DVector

Examples

>>> p = MPS_rand_state(3, 7)
>>> p.gate_(spin_operator('X'), where=1, tags=['GX'])
>>> p
<MatrixProductState(tensors=4, L=3, max_bond=7)>
>>> p.outer_inds()
('k0', 'k1', 'k2')
quimb.tensor.tensor_1d.superop_TN_1D(tn_super, tn_op, upper_ind_id='k{}', lower_ind_id='b{}', so_outer_upper_ind_id=None, so_inner_upper_ind_id=None, so_inner_lower_ind_id=None, so_outer_lower_ind_id=None)[source]#

Take a tensor network superoperator and act with it on a tensor network operator, maintaining the original upper and lower indices of the operator:

outer_upper_ind_id                           upper_ind_id
   | | | ... |                               | | | ... |
   +----------+                              +----------+
   | tn_super +---+                          | tn_super +---+
   +----------+   |     upper_ind_id         +----------+   |
   | | | ... |    |      | | | ... |         | | | ... |    |
inner_upper_ind_id|     +-----------+       +-----------+   |
                  |  +  |   tn_op   |   =   |   tn_op   |   |
inner_lower_ind_id|     +-----------+       +-----------+   |
   | | | ... |    |      | | | ... |         | | | ... |    |
   +----------+   |      lower_ind_id        +----------+   |
   | tn_super +---+                          | tn_super +---+
   +----------+                              +----------+
   | | | ... | <--                           | | | ... |
outer_lower_ind_id                           lower_ind_id
Parameters
  • tn_super (TensorNetwork) – The superoperator in the form of a 1D-like tensor network.

  • tn_op (TensorNetwork) – The operator to be acted on in the form of a 1D-like tensor network.

  • upper_ind_id (str, optional) – Current id of the upper operator indices, e.g. usually 'k{}'.

  • lower_ind_id (str, optional) – Current id of the lower operator indices, e.g. usually 'b{}'.

  • so_outer_upper_ind_id (str, optional) – Current id of the superoperator’s upper outer indices, these will be reindexed to form the new effective operators upper indices.

  • so_inner_upper_ind_id (str, optional) – Current id of the superoperator’s upper inner indices, these will be joined with those described by upper_ind_id.

  • so_inner_lower_ind_id (str, optional) – Current id of the superoperator’s lower inner indices, these will be joined with those described by lower_ind_id.

  • so_outer_lower_ind_id (str, optional) – Current id of the superoperator’s lower outer indices, these will be reindexed to form the new effective operators lower indices.

Returns

KAK – The tensornetwork of the superoperator acting on the operator.

Return type

TensorNetwork