quimb.tensor.tensor_1d ====================== .. py:module:: quimb.tensor.tensor_1d .. autoapi-nested-parse:: Classes and algorithms related to 1D tensor networks. Attributes ---------- .. autoapisummary:: quimb.tensor.tensor_1d.align_TN_1D Classes ------- .. autoapisummary:: quimb.tensor.tensor_1d.TensorNetwork1D quimb.tensor.tensor_1d.TensorNetwork1DVector quimb.tensor.tensor_1d.TensorNetwork1DOperator quimb.tensor.tensor_1d.TensorNetwork1DFlat quimb.tensor.tensor_1d.MatrixProductState quimb.tensor.tensor_1d.MatrixProductOperator quimb.tensor.tensor_1d.Dense1D quimb.tensor.tensor_1d.SuperOperator1D quimb.tensor.tensor_1d.TNLinearOperator1D Functions --------- .. autoapisummary:: quimb.tensor.tensor_1d.expec_TN_1D quimb.tensor.tensor_1d.maybe_factor_gate_into_tensor quimb.tensor.tensor_1d.gate_TN_1D quimb.tensor.tensor_1d.superop_TN_1D quimb.tensor.tensor_1d.parse_cur_orthog quimb.tensor.tensor_1d.convert_cur_orthog quimb.tensor.tensor_1d.set_default_compress_mode Module Contents --------------- .. py:data:: align_TN_1D .. py:function:: expec_TN_1D(*tns, compress=None, eps=1e-15) Compute the expectation of several 1D TNs, using transfer matrix compression if any are periodic. :param tns: The MPS and MPO to find expectation of. Should start and begin with an MPS e.g. ``(MPS, MPO, ..., MPS)``. :type tns: sequence of TensorNetwork1D :param compress: Whether to perform transfer matrix compression on cyclic systems. If set to ``None`` (the default), decide heuristically. :type compress: {None, False, True}, optional :param eps: The accuracy of the transfer matrix compression. :type eps: float, optional :returns: **x** -- The expectation value. :rtype: float .. py:function:: maybe_factor_gate_into_tensor(G, phys_dim, nsites, where) .. py:function:: gate_TN_1D(tn, G, where, contract=False, tags=None, propagate_tags='sites', info=None, inplace=False, cur_orthog=None, **compress_opts) Act with the gate ``G`` on sites ``where``, maintaining the outer indices of the 1D tensor network:: 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'. :param tn: The 1D vector-like tensor network, for example, and MPS. :type tn: TensorNetwork1DVector :param G: 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. :type G: array :param where: Where the gate should act. :type where: int or sequence of int :param contract: '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. :type contract: {False, 'split-gate', 'swap-split-gate', :param tags: Tag the new gate tensor with these tags. :type tags: str or sequence of str, optional :param propagate_tags: 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. :type propagate_tags: {'sites', 'register', False, True}, optional :param inplace: Perform the gate in place. :param bool: Perform the gate in place. :param optional: Perform the gate in place. :param compress_opts: Supplied to :meth:`~quimb.tensor.tensor_core.Tensor.split` if ``contract='swap+split'`` or :meth:`~quimb.tensor.tensor_1d.MatrixProductState.gate_with_auto_swap` if ``contract='swap+split'``. :rtype: TensorNetwork1DVector .. seealso:: :py:obj:`MatrixProductState.gate_split` .. rubric:: Examples >>> p = MPS_rand_state(3, 7) >>> p.gate_(spin_operator('X'), where=1, tags=['GX']) >>> p >>> p.outer_inds() ('k0', 'k1', 'k2') .. py:function:: 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) 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 :param tn_super: The superoperator in the form of a 1D-like tensor network. :type tn_super: TensorNetwork :param tn_op: The operator to be acted on in the form of a 1D-like tensor network. :type tn_op: TensorNetwork :param upper_ind_id: Current id of the upper operator indices, e.g. usually ``'k{}'``. :type upper_ind_id: str, optional :param lower_ind_id: Current id of the lower operator indices, e.g. usually ``'b{}'``. :type lower_ind_id: str, optional :param so_outer_upper_ind_id: Current id of the superoperator's upper outer indices, these will be reindexed to form the new effective operators upper indices. :type so_outer_upper_ind_id: str, optional :param so_inner_upper_ind_id: Current id of the superoperator's upper inner indices, these will be joined with those described by ``upper_ind_id``. :type so_inner_upper_ind_id: str, optional :param so_inner_lower_ind_id: Current id of the superoperator's lower inner indices, these will be joined with those described by ``lower_ind_id``. :type so_inner_lower_ind_id: str, optional :param so_outer_lower_ind_id: Current id of the superoperator's lower outer indices, these will be reindexed to form the new effective operators lower indices. :type so_outer_lower_ind_id: str, optional :returns: **KAK** -- The tensornetwork of the superoperator acting on the operator. :rtype: TensorNetwork .. py:function:: parse_cur_orthog(cur_orthog='calc', info=None) .. py:function:: convert_cur_orthog(fn) .. py:class:: TensorNetwork1D(ts=(), *, virtual=False, check_collisions=True) Bases: :py:obj:`quimb.tensor.tensor_arbgeom.TensorNetworkGen` Base class for tensor networks with a one-dimensional structure. .. py:attribute:: _NDIMS :value: 1 .. py:attribute:: _EXTRA_PROPS :value: ('_site_tag_id', '_L') .. py:attribute:: _CONTRACT_STRUCTURED :value: True .. py:method:: _compatible_1d(other) Check whether ``self`` and ``other`` are compatible 2D tensor networks such that they can remain a 2D tensor network when combined. .. py:method:: combine(other, *, virtual=False, check_collisions=True) Combine this tensor network with another, returning a new tensor network. If the two are compatible, cast the resulting tensor network to a :class:`TensorNetwork1D` instance. :param other: The other tensor network to combine with. :type other: TensorNetwork1D or TensorNetwork :param virtual: Whether the new tensor network should copy all the incoming tensors (``False``, the default), or view them as virtual (``True``). :type virtual: bool, optional :param check_collisions: Whether to check for index collisions between the two tensor networks before combining them. If ``True`` (the default), any inner indices that clash will be mangled. :type check_collisions: bool, optional :rtype: TensorNetwork1D or TensorNetwork .. py:property:: L The number of sites, i.e. length. .. py:property:: nsites The number of sites. .. py:method:: gen_site_coos() Generate the coordinates of all possible sites. .. py:method:: site_tag(i) The name of the tag specifiying the tensor at site ``i``. .. py:method:: slice2sites(tag_slice) Take a slice object, and work out its implied start, stop and step, taking into account cyclic boundary conditions. .. rubric:: 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. .. py:method:: maybe_convert_coo(x) Check if ``x`` is an integer and convert to the corresponding site tag if so. .. py:method:: contract_structured(tag_slice, structure_bsz=5, inplace=False, **opts) Perform a structured contraction, translating ``tag_slice`` from a ``slice`` or `...` to a cumulative sequence of tags. :param tag_slice: The range of sites, or `...` for all. :type tag_slice: slice or ... :param inplace: Whether to perform the contraction inplace. :type inplace: bool, optional :returns: The result of the contraction, still a ``TensorNetwork`` if the contraction was only partial. :rtype: TensorNetwork, Tensor or scalar .. seealso:: :py:obj:`contract`, :py:obj:`contract_tags`, :py:obj:`contract_cumulative` .. py:method:: compute_left_environments(**contract_opts) Compute the left environments of this 1D tensor network. :param contract_opts: Supplied to :meth:`~quimb.tensor.tensor_core.TensorNetwork.contract`. :returns: Environments indexed by the site they are to the left of, so keys run from (1, ... L - 1). :rtype: dict[int, Tensor] .. py:method:: compute_right_environments(**contract_opts) Compute the right environments of this 1D tensor network. :param contract_opts: Supplied to :meth:`~quimb.tensor.tensor_core.TensorNetwork.contract`. :returns: Environments indexed by the site they are to the right of, so keys run from (0, ... L - 2). :rtype: dict[int, Tensor] .. py:method:: _repr_info() General info to show in various reprs. Sublasses can add more relevant info to this dict. .. py:class:: TensorNetwork1DVector(ts=(), *, virtual=False, check_collisions=True) Bases: :py:obj:`TensorNetwork1D`, :py:obj:`quimb.tensor.tensor_arbgeom.TensorNetworkGenVector` 1D Tensor network which overall is like a vector with a single type of site ind. .. py:attribute:: _EXTRA_PROPS :value: ('_site_tag_id', '_site_ind_id', '_L') .. py:method:: reindex_sites(new_id, where=None, inplace=False) Update the physical site index labels to a new string specifier. Note that this doesn't change the stored id string with the TN. :param new_id: A string with a format placeholder to accept an int, e.g. "ket{}". :type new_id: str :param where: Which sites to update the index labels on. If ``None`` (default) all sites. :type where: None or slice :param inplace: Whether to reindex in place. :type inplace: bool .. py:attribute:: reindex_sites_ .. py:method:: site_ind(i) Get the physical index name of site ``i``. .. py:method:: gate(*args, inplace=False, **kwargs) .. py:attribute:: gate_ .. py:method:: expec(*args, **kwargs) .. py:method:: correlation(A, i, j, B=None, **expec_opts) Correlation of operator ``A`` between ``i`` and ``j``. :param A: The operator to act with, can be multi site. :type A: array :param i: The first site(s). :type i: int or sequence of int :param j: The second site(s). :type j: int or sequence of int :param expec_opts: Supplied to :func:`~quimb.tensor.tensor_1d.expec_TN_1D`. :returns: **C** -- The correlation `` + - ``. :rtype: float .. rubric:: 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 .. py:class:: TensorNetwork1DOperator(ts=(), *, virtual=False, check_collisions=True) Bases: :py:obj:`TensorNetwork1D`, :py:obj:`quimb.tensor.tensor_arbgeom.TensorNetworkGenOperator` Base class for tensor networks with a one-dimensional structure. .. py:attribute:: _EXTRA_PROPS :value: ('_site_tag_id', '_upper_ind_id', '_lower_ind_id', '_L') .. py:method:: reindex_lower_sites(new_id, where=None, inplace=False) Update the lower site index labels to a new string specifier. :param new_id: A string with a format placeholder to accept an int, e.g. ``"ket{}"``. :type new_id: str :param where: Which sites to update the index labels on. If ``None`` (default) all sites. :type where: None or slice :param inplace: Whether to reindex in place. :type inplace: bool .. py:attribute:: reindex_lower_sites_ .. py:method:: reindex_upper_sites(new_id, where=None, inplace=False) Update the upper site index labels to a new string specifier. :param new_id: A string with a format placeholder to accept an int, e.g. "ket{}". :type new_id: str :param where: Which sites to update the index labels on. If ``None`` (default) all sites. :type where: None or slice :param inplace: Whether to reindex in place. :type inplace: bool .. py:attribute:: reindex_upper_sites_ .. py:function:: set_default_compress_mode(opts, cyclic=False) .. py:class:: TensorNetwork1DFlat(ts=(), *, virtual=False, check_collisions=True) Bases: :py:obj:`TensorNetwork1D` 1D Tensor network which has a flat structure. .. py:attribute:: _EXTRA_PROPS :value: ('_site_tag_id', '_L') .. py:method:: left_canonize_site(i, bra=None) Left canonize this TN's ith site, inplace:: i i -o-o- ->-s- ... | | ... ==> ... | | ... :param i: Which site to canonize. The site at i + 1 also absorbs the non-isometric part of the decomposition of site i. :type i: int :param bra: If set, also update this TN's data with the conjugate canonization. :type bra: None or matching TensorNetwork to self, optional .. py:method:: right_canonize_site(i, bra=None) Right canonize this TN's ith site, inplace:: i i -o-o- -s-<- ... | | ... ==> ... | | ... :param i: 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. :type i: int .. py:method:: left_canonicalize(stop=None, start=None, normalize=False, bra=None, inplace=False) Left canonicalize all or a portion of this TN (i.e. sweep the orthogonality center to the right). If this is a MPS, this implies that:: i i >->->->->->->-o-o- +-o-o- | | | | | | | | | ... => | | | ... >->->->->->->-o-o- +-o-o- :param start: If given, the site to start left canonicalizing at. :type start: int, optional :param stop: If given, the site to stop left canonicalizing at. :type stop: int, optional :param normalize: Whether to normalize the state, only works for OBC. :type normalize: bool, optional :param bra: If supplied, simultaneously left canonicalize this MPS too, assuming it to be the conjugate state. :type bra: MatrixProductState, optional :param inplace: Whether to perform the operation inplace. If ``bra`` is supplied then it is always modifed inplace. :type inplace: bool, optional :rtype: TensorNetwork1DFlat .. py:attribute:: left_canonicalize_ .. py:attribute:: left_canonize .. py:method:: right_canonicalize(stop=None, start=None, normalize=False, bra=None, inplace=False) Right canonicalize all or a portion of this TN (i.e. sweep the orthogonality center to the left). If this is a MPS, i i -o-o-<-<-<-<-<-<-< -o-o-+ ... | | | | | | | | | -> ... | | | -o-o-<-<-<-<-<-<-< -o-o-+ :param start: If given, the site to start right canonizing at. :type start: int, optional :param stop: If given, the site to stop right canonizing at. :type stop: int, optional :param normalize: Whether to normalize the state. :type normalize: bool, optional :param bra: If supplied, simultaneously right canonicalize this MPS too, assuming it to be the conjugate state. :type bra: MatrixProductState, optional :param inplace: Whether to perform the operation inplace. If ``bra`` is supplied then it is always modifed inplace. :type inplace: bool, optional :rtype: TensorNetwork1DFlat .. py:attribute:: right_canonicalize_ .. py:attribute:: right_canonize .. py:method:: canonize_cyclic(i, bra=None, method='isvd', inv_tol=1e-10) Bring this MatrixProductState into (possibly only approximate) canonical form at site(s) ``i``. :param i: The site or range of sites to make canonical. :type i: int or slice :param bra: Simultaneously canonize this state as well, assuming it to be the co-vector. :type bra: MatrixProductState, optional :param method: How to perform the lateral compression. :type method: {'isvd', 'svds', ...}, optional :param inv_tol: Tolerance with which to invert the gauge. :type inv_tol: float, optional .. py:method:: shift_orthogonality_center(current, new, bra=None) Move the orthogonality center of this MPS. :param current: The current orthogonality center. :type current: int :param new: The target orthogonality center. :type new: int :param bra: If supplied, simultaneously move the orthogonality center of this MPS too, assuming it to be the conjugate state. :type bra: MatrixProductState, optional .. py:method:: canonicalize(where, cur_orthog='calc', info=None, bra=None, inplace=False) Gauge this MPS into mixed canonical form, implying:: i i >->->->->- ->-o-<- -<-<-<-<-< +-o-+ | | | | |...| | |...| | | | | -> | | | >->->->->- ->-o-<- -<-<-<-<-< +-o-+ You can also supply a min/max of sites to orthogonalize around, and a current location of the orthogonality center for efficiency:: current where ....... ..... >->->-c-c-c-c-<-<-<-<-<-< >->->->->->-w-<-<-<-<-<-< | | | | | | | | | | | | | -> | | | | | | | | | | | | | >->->-c-c-c-c-<-<-<-<-<-< >->->->->->-w-<-<-<-<-<-< cmin cmax i j This would only move ``cmin`` to ``i`` and ``cmax`` to ``j`` if necessary. :param where: Which site(s) to orthogonalize around. If a sequence of int then make sure that section from min(where) to max(where) is orthog. :type where: int or sequence of int :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param bra: If supplied, simultaneously mixed canonicalize this MPS too, assuming it to be the conjugate state. :type bra: MatrixProductState, optional :param inplace: Whether to perform the operation inplace. If ``bra`` is supplied then it is always modifed inplace. :type inplace: bool, optional :returns: The mixed canonical form MPS. :rtype: MatrixProductState .. py:attribute:: canonicalize_ .. py:attribute:: canonize .. py:method:: left_compress_site(i, bra=None, **compress_opts) Left compress this 1D TN's ith site, such that the site is then left unitary with its right bond (possibly) reduced in dimension. :param i: Which site to compress. :type i: int :param bra: If set, also update this TN's data with the conjugate compression. :type bra: None or matching TensorNetwork to self, optional :param compress_opts: Supplied to :meth:`Tensor.split`. .. py:method:: right_compress_site(i, bra=None, **compress_opts) Right compress this 1D TN's ith site, such that the site is then right unitary with its left bond (possibly) reduced in dimension. :param i: Which site to compress. :type i: int :param bra: If set, update this TN's data with the conjugate compression. :type bra: None or matching TensorNetwork to self, optional :param compress_opts: Supplied to :meth:`Tensor.split`. .. py:method:: left_compress(start=None, stop=None, bra=None, **compress_opts) Compress this 1D TN, from left to right, such that it becomes left-canonical (unless ``absorb != 'right'``). :param start: Site to begin compressing on. :type start: int, optional :param stop: Site to stop compressing at (won't itself be an isometry). :type stop: int, optional :param bra: If given, update this TN as well, assuming it to be the conjugate. :type bra: None or TensorNetwork like this one, optional :param compress_opts: Supplied to :meth:`Tensor.split`. .. py:method:: right_compress(start=None, stop=None, bra=None, **compress_opts) Compress this 1D TN, from right to left, such that it becomes right-canonical (unless ``absorb != 'left'``). :param start: Site to begin compressing on. :type start: int, optional :param stop: Site to stop compressing at (won't itself be an isometry). :type stop: int, optional :param bra: If given, update this TN as well, assuming it to be the conjugate. :type bra: None or TensorNetwork like this one, optional :param compress_opts: Supplied to :meth:`Tensor.split`. .. py:method:: compress(form=None, **compress_opts) Compress this 1D Tensor Network, possibly into canonical form. :param form: 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. :type form: {None, 'flat', 'left', 'right'} or int :param compress_opts: Supplied to :meth:`Tensor.split`. .. py:method:: compress_site(i, canonize=True, info=None, bra=None, **compress_opts) 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~<-<- | | | | | | | | | | :param i: Which site to compress around :type i: int :param canonize: Whether to first set the orthogonality center to site ``i``. :type canonize: bool, optional :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param bra: The conjugate state to also apply the compression to. :type bra: MatrixProductState, optional :param compress_opts: Supplied to :func:`~quimb.tensor.tensor_core.tensor_split`. .. py:method:: bond(i, j) Get the name of the index defining the bond between sites i and j. .. py:method:: bond_size(i, j) Return the size of the bond between site ``i`` and ``j``. .. py:method:: bond_sizes() .. py:method:: amplitude(b) Compute the amplitude of configuration ``b``. :param b: The configuration to compute the amplitude of. :type b: sequence of int :returns: **c_b** :rtype: scalar .. py:method:: singular_values(i, info=None, method='svd') 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``. :param i: Which bond, or equivalently, the number of sites in the left partition. :type i: int :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :returns: **svals** -- The singular values. :rtype: 1d-array .. py:method:: expand_bond_dimension(new_bond_dim, rand_strength=0.0, bra=None, inplace=True) Expand the bond dimensions of this 1D tensor network to at least ``new_bond_dim``. :param new_bond_dim: Minimum bond dimension to expand to. :type new_bond_dim: int :param inplace: Whether to perform the expansion in place. :type inplace: bool, optional :param bra: Mirror the changes to ``bra`` inplace, treating it as the conjugate state. :type bra: MatrixProductState, optional :param rand_strength: If ``rand_strength > 0``, fill the new tensor entries with gaussian noise of strength ``rand_strength``. :type rand_strength: float, optional :rtype: MatrixProductState .. py:method:: count_canonized() .. py:method:: calc_current_orthog_center() Calculate the site(s) of the current orthogonality center. :returns: The min/max of sites around which the TN is currently orthogonal. :rtype: (int, int) .. py:method:: as_cyclic(inplace=False) Convert this flat, 1D, TN into cyclic form by adding a dummy bond between the first and last sites. .. py:method:: show(max_width=None) .. py:class:: MatrixProductState(arrays, *, sites=None, L=None, shape='lrp', tags=None, site_ind_id='k{}', site_tag_id='I{}', **tn_opts) Bases: :py:obj:`TensorNetwork1DVector`, :py:obj:`TensorNetwork1DFlat` Initialise a matrix product state, with auto labelling and tagging. :param arrays: The tensor arrays to form into a MPS. :type arrays: sequence of arrays :param sites: Construct the MPO on these sites only. If not given, enumerate from zero. Should be monotonically increasing and match ``arrays``. :type sites: sequence of int, optional :param L: The number of sites the MPO should be defined on. If not given, this is taken as the max ``sites`` value plus one (i.e.g the number of arrays if ``sites`` is not given). :type L: int, optional :param shape: 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. :type shape: str, optional :param tags: Global tags to attach to all tensors. :type tags: str or sequence of str, optional :param site_ind_id: 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)))``. :type site_ind_id: str :param site_tag_id: 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)))``. :type site_tag_id: str .. py:attribute:: _EXTRA_PROPS :value: ('_site_tag_id', '_site_ind_id', 'cyclic', '_L') .. py:attribute:: _L .. py:attribute:: _site_ind_id :value: 'k{}' .. py:attribute:: _site_tag_id :value: 'I{}' .. py:attribute:: cyclic .. py:method:: from_fill_fn(fill_fn, L, bond_dim, phys_dim=2, sites=None, cyclic=False, shape='lrp', site_ind_id='k{}', site_tag_id='I{}', tags=None) :classmethod: Create an MPS by supplying a 'filling' function to generate the data for each site. :param fill_fn: A function with signature ``fill_fn(shape : tuple[int]) -> array_like``. :type fill_fn: callable :param L: The number of sites. :type L: int :param bond_dim: The bond dimension. :type bond_dim: int :param phys_dim: The physical dimension(s) of each site, if a sequence it will be cycled over. :type phys_dim: int or Sequence[int], optional :param sites: Construct the MPS on these sites only. If not given, enumerate from zero. :type sites: None or sequence of int, optional :param cyclic: Whether the MPS should be cyclic (periodic). :type cyclic: bool, optional :param shape: What specific order to layout the indices in, should be a sequence of ``'l'``, ``'r'``, and ``'p'``, corresponding to left, right, and physical indices respectively. :type shape: str, optional :param site_ind_id: How to label the physical site indices. :type site_ind_id: str, optional :param site_tag_id: How to tag the physical sites. :type site_tag_id: str, optional :param tags: Global tags to attach to all tensors. :type tags: str or sequence of str, optional :rtype: MatrixProductState .. py:method:: from_dense(psi, dims=2, tags=None, site_ind_id='k{}', site_tag_id='I{}', **split_opts) :classmethod: Create a ``MatrixProductState`` directly from a dense vector :param psi: The dense state to convert to MPS from. :type psi: array_like :param dims: Physical subsystem dimensions of each site. If a single int, all sites have this same dimension, by default, 2. :type dims: int or sequence of int :param tags: Global tags to attach to all tensors. :type tags: str or sequence of str, optional :param site_ind_id: How to index the physical sites, see :class:`~quimb.tensor.tensor_1d.MatrixProductState`. :type site_ind_id: str, optional :param site_tag_id: How to tag the physical sites, see :class:`~quimb.tensor.tensor_1d.MatrixProductState`. :type site_tag_id: str, optional :param split_opts: Supplied to :func:`~quimb.tensor.tensor_core.tensor_split` to in order to partition the dense vector into tensors. ``absorb='left'`` is set by default, to ensure the compression is canonical / optimal. :rtype: MatrixProductState .. rubric:: 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 | | | | | | .. py:method:: add_MPS(other, inplace=False, **kwargs) Add another MatrixProductState to this one. .. py:attribute:: add_MPS_ .. py:method:: permute_arrays(shape='lrp') 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. :param shape: A permutation of ``'lrp'`` specifying the desired order of the left, right, and physical indices respectively. :type shape: str, optional .. py:method:: normalize(bra=None, eps=1e-15, insert=None) 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. :param bra: If given, normalize this MPS with the same factor. :type bra: MatrixProductState, optional :param eps: If cyclic, precision to approximation transfer matrix with. Default: 1e-14. :type eps: float, optional :param insert: Insert the corrective normalization on this site, random if not given. :type insert: int, optional :returns: **old_norm** -- The old norm ``self.H @ self``. :rtype: float .. py:method:: gate_split(G, where, inplace=False, **compress_opts) 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. :param G: The gate, with shape ``(d**2, d**2)`` for physical dimension ``d``. :type G: array :param where: Indices of the sites to apply the gate to. :type where: (int, int) :param compress_opts: Supplied to :func:`~quimb.tensor.tensor_split`. .. seealso:: :py:obj:`gate`, :py:obj:`gate_with_auto_swap` .. py:attribute:: gate_split_ .. py:method:: swap_sites_with_compress(i, j, info=None, inplace=False, **compress_opts) Swap sites ``i`` and ``j`` by contracting, then splitting with the physical indices swapped. If the sites are not adjacent, this will happen multiple times. :param i: The first site to swap. :type i: int :param j: The second site to swap. :type j: int :param cur_orthog: If known, the current orthogonality center. :type cur_orthog: int, sequence of int, or 'calc' :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param inplace: Perform the swaps inplace. :type inplace: bond, optional :param compress_opts: Supplied to :func:`~quimb.tensor.tensor_core.tensor_split`. .. py:attribute:: swap_sites_with_compress_ .. py:method:: swap_site_to(i, f, info=None, inplace=False, **compress_opts) 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 >->->->->->->-x-<-< | | | | | | | | | | -> | | | | | | | | | | :param i: The site to move. :type i: int :param f: The new location for site ``i``. :type f: int :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param inplace: Perform the swaps inplace. :type inplace: bond, optional :param compress_opts: Supplied to :func:`~quimb.tensor.tensor_core.tensor_split`. .. py:attribute:: swap_site_to_ .. py:method:: gate_with_auto_swap(G, where, info=None, swap_back=True, inplace=False, **compress_opts) 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. :param G: The gate, with shape ``(d**2, d**2)`` for physical dimension ``d``. :type G: array :param where: Indices of the sites to apply the gate to. :type where: (int, int) :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param swap_back: Whether to swap the sites back to their original position after applying the gate. If not, for sites ``i < j``, the site ``j`` will remain swapped to ``i + 1``, and sites between ``i + 1`` and ``j`` will be shifted one place up. :type swap_back: bool, optional :param inplace: Perform the swaps inplace. :type inplace: bond, optional :param compress_opts: Supplied to :func:`~quimb.tensor.tensor_core.tensor_split`. .. seealso:: :py:obj:`gate`, :py:obj:`gate_split` .. py:attribute:: gate_with_auto_swap_ .. py:method:: gate_with_submpo(submpo, where=None, method='direct', transpose=False, info=None, inplace=False, inplace_mpo=False, **compress_opts) Apply an MPO, which only acts on a subset of sites, to this MPS, compressing the MPS with the MPO only on the minimal set of sites covering `where`, keeping the MPS form:: │ │ │ A───A─A │ │ │ -> │ │ │ │ │ │ │ │ >─>─O━O━O━O─<─< │ │ │ │ │ │ │ │ o─o─o─o─o─o─o─o :param submpo: The MPO to apply. :type submpo: MatrixProductOperator :param where: The range of sites the MPO acts on, will be inferred from the support of the MPO if not given. :type where: sequence of int, optional :param method: The compression method to use. :type method: {'direct", 'dm', 'zipup', 'zipup-first', 'fit'}, optional :param transpose: Whether to transpose the MPO before applying it. By default the lower inds of the MPO are contracted with the MPS, if transposed the upper inds are contracted. :type transpose: bool, optional :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param inplace: Whether to perform the application and compression inplace. :type inplace: bool, optional :param compress_opts: Supplied to :func:`~quimb.tensor.tensor_1d_compress.tensor_network_1d_compress`. :rtype: MatrixProductState .. py:attribute:: gate_with_submpo_ .. py:method:: gate_with_mpo(mpo, method='direct', transpose=False, inplace=False, inplace_mpo=False, **compress_opts) Gate this MPS with an MPO and compress the result with one of various methods back to MPS form:: │ │ │ │ │ │ │ │ A─A─A─A─A─A─A─A │ │ │ │ │ │ │ │ -> │ │ │ │ │ │ │ │ O━O━O━O━O━O━O━O │ │ │ │ │ │ │ │ o─o─o─o─o─o─o─o :param mpo: The MPO to apply. :type mpo: MatrixProductOperator :param max_bond: A maximum bond dimension to keep when compressing. :type max_bond: int, optional :param cutoff: A singular value cutoff to use when compressing. :type cutoff: float, optional :param method: The compression method to use. :type method: {'direct", 'dm', 'zipup', 'zipup-first', 'fit', ...}, optional :param transpose: Whether to transpose the MPO before applying it. By default the lower inds of the MPO are contracted with the MPS, if transposed the upper inds are contracted. :type transpose: bool, optional :param inplace: Whether to perform the compression inplace. :type inplace: bool, optional :param inplace_mpo: Whether the modify the MPO in place, a minor performance gain. :type inplace_mpo: bool, optional :param compress_opts: Other options supplied to :func:`~quimb.tensor.tensor_1d_compress.tensor_network_1d_compress`. :rtype: MatrixProductState .. py:attribute:: gate_with_mpo_ .. py:method:: gate_nonlocal(G, where, dims=None, method='direct', info=None, inplace=False, **compress_opts) Apply a potentially non-local gate to this MPS by first decomposing it into an MPO, then compressing the MPS with MPO only on the minimal set of sites covering `where`. :param G: The gate to apply. :type G: array_like :param where: The sites to apply the gate to. :type where: sequence of int :param max_bond: A maximum bond dimension to keep when compressing. :type max_bond: int, optional :param cutoff: A singular value cutoff to use when compressing. :type cutoff: float, optional :param dims: The factorized dimensions of the gate ``G``, which should match the physical dimensions of the sites it acts on. Calculated if not supplied. If a single int, all sites are assumed to have this same dimension. :type dims: sequence of int, optional :param method: The compression method to use. :type method: {'direct", 'dm', 'zipup', 'zipup-first', 'fit', ...}, optional :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param inplace: Whether to perform the compression inplace. :type inplace: bool, optional :param compress_opts: Supplied to :func:`~quimb.tensor.tensor_1d_compress.tensor_network_1d_compress`. :rtype: MatrixProductState .. py:attribute:: gate_nonlocal_ .. py:method:: flip(inplace=False) Reverse the order of the sites in the MPS, such that site ``i`` is now at site ``L - i - 1``. .. py:method:: magnetization(i, direction='Z', info=None) Compute the magnetization at site ``i``. .. py:method:: schmidt_values(i, info=None, method='svd') 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. :param i: The number of sites in the left partition. :type i: int :param info: If given, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. :type info: dict, optional :returns: **S** -- The schmidt values. :rtype: 1d-array .. py:method:: entropy(i, info=None, method='svd') The entropy of bipartition between the left block of ``i`` sites and the rest. :param i: The number of sites in the left partition. :type i: int :param info: If given, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. :type info: dict, optional :rtype: float .. py:method:: schmidt_gap(i, info=None, method='svd') The schmidt gap of bipartition between the left block of ``i`` sites and the rest. :param i: The number of sites in the left partition. :type i: int :param info: If given, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. :type info: dict, optional :rtype: float .. py:method:: partial_trace_to_mpo(keep, upper_ind_id='b{}', rescale_sites=True) Partially trace this matrix product state, producing a matrix product operator. :param keep: Indicies of the sites to keep. :type keep: sequence of int or slice :param upper_ind_id: The ind id of the (new) 'upper' inds, i.e. the 'bra' inds. :type upper_ind_id: str, optional :param rescale_sites: If ``True`` (the default), then the kept sites will be rescaled to ``(0, 1, 2, ...)`` etc. rather than keeping their original site numbers. :type rescale_sites: bool, optional :returns: **rho** -- The density operator in MPO form. :rtype: MatrixProductOperator .. py:method:: partial_trace(*_, **__) Partially trace this tensor network state, keeping only the sites in ``keep``, using compressed contraction. :param keep: The sites to keep. :type keep: iterable of hashable :param max_bond: The maximum bond dimensions to use while compressed contracting. :type max_bond: int :param optimize: The contraction path optimizer to use, should specifically generate contractions paths designed for compressed contraction. :type optimize: str or PathOptimizer, optional :param flatten: Whether to force 'flattening' (contracting all physical indices) of the tensor network before contraction, whilst this makes the TN generally more complex to contract, the accuracy is usually improved. If ``'all'`` also flatten the tensors in ``keep``. :type flatten: {False, True, 'all'}, optional :param reduce: Whether to first 'pull' the physical indices off their respective tensors using QR reduction. Experimental. :type reduce: bool, optional :param normalized: Whether to normalize the reduced density matrix at the end. :type normalized: bool, optional :param symmetrized: Whether to symmetrize the reduced density matrix at the end. This should be unecessary if ``flatten`` is set to ``True``. :type symmetrized: {'auto', True, False}, optional :param rehearse: Whether to perform the computation or not:: - False: perform the computation. - 'tn': return the tensor network without running the path optimizer. - 'tree': run the path optimizer and return the ``cotengra.ContractonTree``.. - True: run the path optimizer and return the ``PathInfo``. :type rehearse: {False, 'tn', 'tree', True}, optional :param contract_compressed_opts: Additional keyword arguments to pass to :meth:`~quimb.tensor.tensor_core.TensorNetwork.contract_compressed`. :type contract_compressed_opts: dict, optional :returns: **rho** -- The reduce density matrix of sites in ``keep``. :rtype: array_like .. py:method:: ptr(*_, **__) .. py:method:: partial_trace_to_dense_canonical(where, normalized=True, info=None, **contract_opts) Compute the dense local reduced density matrix by canonicalizing around the target sites and then contracting the local tensors. Note this moves the orthogonality around inplace, and records it in `info`. :param where: The site or sites to compute the reduced density matrix for. :type where: int or tuple[int] :param normalized: Explicitly normalize the local reduced density matrix. :type normalized: bool, optional :param info: If supplied, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param contract_opts: Passed to `tensor_contract` when computing the reduced local density matrix. :rtype: array_like .. py:method:: local_expectation_canonical(G, where, normalized=True, info=None, **contract_opts) Compute a local expectation value (via forming the reduced density matrix). Note this moves the orthogonality around inplace, and records it in `info`. :param G: The local operator to compute the expectation of. :type G: array_like :param where: The site or sites to compute the expectation at. :type where: int or tuple[int] :param normalized: Explicitly normalize the local reduced density matrix. :type normalized: bool, optional :param info: If supplied, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param contract_opts: Passed to `tensor_contract` when computing the reduced local density matrix. :rtype: float .. py:method:: compute_local_expectation_canonical(terms, normalized=True, return_all=False, info=None, inplace=False, **contract_opts) Compute many local expectations at once, via forming the relevant reduced density matrices via canonicalization. This moves the orthogonality around inplace, and records it in `info`. :param terms: The local terms to compute values for. :type terms: dict[int or tuple[int], array_like] :param normalized: Explicitly normalize each local reduced density matrix. :type normalized: bool, optional :param return_all: Whether to return each expectation in `terms` separately or sum them all together (the default). :type return_all: bool, optional :param info: If supplied, will be used to infer and store various extra information. Currently, the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param inplace: Whether to perform the required canonicalizations inplace. :type inplace: bool, optional :param contract_opts: Supplied to :meth:`~quimb.tensor.tensor_core.TensorNetwork.contract` when contracting the local density matrices. :returns: The expecetation value(s), either summed or for each term if `return_all=True`. :rtype: float or dict[in or tuple[int], float] .. seealso:: :py:obj:`compute_local_expectation_via_envs`, :py:obj:`local_expectation_canonical`, :py:obj:`partial_trace_to_dense_canonical` .. py:method:: compute_local_expectation_via_envs(terms, normalized=True, return_all=False, **contract_opts) Compute many local expectations at once, via forming the relevant local overlaps using left and right environments formed via contraction. This does not require any canonicalization and can be quicker if the canonical center is not already aligned. :param terms: The local terms to compute values for. :type terms: dict[int or tuple[int], array_like] :param normalized: Explicitly normalize each local reduced density matrix. :type normalized: bool, optional :param return_all: Whether to return each expectation in `terms` separately or sum them all together (the default). :type return_all: bool, optional :param contract_opts: Supplied to :meth:`~quimb.tensor.tensor_core.TensorNetwork.contract` when contracting the local overlaps. :returns: The expecetation value(s), either summed or for each term if `return_all=True`. :rtype: float or dict[int or tuple[int], float] .. seealso:: :py:obj:`compute_local_expectation_canonical`, :py:obj:`compute_left_environments`, :py:obj:`compute_right_environments` .. py:method:: compute_local_expectation(terms, normalized=True, return_all=False, method='canonical', info=None, inplace=False, **contract_opts) Compute many local expectations at once. :param terms: The local terms to compute values for. :type terms: dict[int or tuple[int], array_like] :param normalized: Explicitly normalize each local term. :type normalized: bool, optional :param return_all: Whether to return each expectation in `terms` separately or sum them all together (the default). :type return_all: bool, optional :param method: The method to use to compute the local expectations. - 'canonical': canonicalize around the sites of interest and contract the local reduced density matrices, moving the canonical center around as needed. - 'envs': form the local overlaps using left and right environments and contract these directly. This can be quicker if the canonical center is not already aligned. :type method: {'canonical', 'envs'}, optional :param info: If supplied, and `method=="canonical"`, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. Its input value can be ``"calc"``, a single site, or a pair of sites representing the min/max range, inclusive. It will be updated to the actual range after. :type info: dict, optional :param inplace: If `method=="canonical"`, whether to perform the required canonicalizations inplace or on a copy of the state. :type inplace: bool, optional :param contract_opts: Supplied to :meth:`~quimb.tensor.tensor_core.TensorNetwork.contract` when contracting the local overlaps or density matrices. :returns: The expecetation value(s), either summed or for each term if `return_all=True`. :rtype: float or dict[int or tuple[int], float] .. seealso:: :py:obj:`compute_local_expectation_canonical`, :py:obj:`compute_local_expectation_via_envs` .. py:method:: bipartite_schmidt_state(sz_a, get='ket', info=None) 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 :param sz_a: The number of sites in subsystem A, must be ``0 < sz_a < N``. :type sz_a: int :param get: 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``. :type get: {'ket', 'rho', 'ket-dense', 'rho-dense'}, optional :param info: If given, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. :type info: dict, optional .. py:method:: _do_lateral_compress(mps, kb, section, leave_short, ul, ll, heps, hmethod, hmax_bond, verbosity, compressed, **compress_opts) :staticmethod: .. py:method:: _do_vertical_decomp(mps, kb, section, sysa, sysb, compressed, ul, ur, ll, lr, vmethod, vmax_bond, veps, verbosity, **compress_opts) :staticmethod: .. py:method:: 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) 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. :param sysa: The sites, which should be contiguous, defining subsystem A. :type sysa: sequence of int :param sysb: The sites, which should be contiguous, defining subsystem B. :type sysb: sequence of int :param eps: Tolerance(s) to use when compressing the subsystem transfer matrices and vertically decomposing. :type eps: float or (float, float), optional :param method: Method(s) to use for laterally compressing the state then vertially compressing subsytems. :type method: str or (str, str), optional :param max_bond: The maximum bond to keep for laterally compressing the state then vertially compressing subsytems. :type max_bond: int or (int, int), optional :param leave_short: If True (the default), don't try to compress short sections. :type leave_short: bool, optional :param renorm: If True (the default), renomalize the state so that ``tr(rho)==1``. :type renorm: bool, optional :param lower_ind_id: The index id to create for the new density matrix, the upper_ind_id is automatically taken as the current site_ind_id. :type lower_ind_id: str, optional :param compress_opts: If given, supplied to ``partial_trace_compress`` to govern how singular values are treated. See ``tensor_split``. :type compress_opts: dict, optional :param verbosity: How much information to print while performing the compressed partial trace. :type verbosity: {0, 1}, optional :returns: **rho_ab** -- Density matrix tensor network with ``outer_inds = ('k0', 'k1', 'b0', 'b1')`` for example. :rtype: TensorNetwork .. py:method:: logneg_subsys(sysa, sysb, compress_opts=None, approx_spectral_opts=None, verbosity=0, approx_thresh=2**12) 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- ... | | | | | | | | | | | | | | | | | | | | | | | | :param sysa: The sites, which should be contiguous, defining subsystem A. :type sysa: sequence of int :param sysb: The sites, which should be contiguous, defining subsystem B. :type sysb: sequence of int :param eps: Tolerance to use when compressing the subsystem transfer matrices. :type eps: float, optional :param method: Method(s) to use for laterally compressing the state then vertially compressing subsytems. :type method: str or (str, str), optional :param compress_opts: If given, supplied to ``partial_trace_compress`` to govern how singular values are treated. See ``tensor_split``. :type compress_opts: dict, optional :param approx_spectral_opts: Supplied to :func:`~quimb.approx_spectral_function`. :returns: **ln** -- The logarithmic negativity. :rtype: float .. seealso:: :py:obj:`MatrixProductState.partial_trace_compress`, :py:obj:`approx_spectral_function` .. py:method:: measure(site, remove=False, outcome=None, renorm=True, info=None, get=None, seed=None, inplace=False) 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)``. :param site: The site to measure. :type site: int :param remove: 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 :type remove: bool, optional :param outcome: Specify the desired outcome of the measurement. If ``None``, it will be randomly sampled according to the local density matrix. :type outcome: None or int, optional :param renorm: Whether to renormalize the state post measurement. :type renorm: bool, optional :param info: If given, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. :type info: dict, optional :param get: If ``'outcome'``, simply return the outcome, and don't perform any projection. :type get: {None, 'outcome'}, optional :param seed: A random seed or generator to use. :type seed: None, int, or np.random.Generator, optional :param inplace: Whether to perform the measurement in place or not. :type inplace: bool, optional :returns: * **outcome** (*int*) -- The measurement outcome, drawn from ``range(phys_dim)``. * **psi** (*MatrixProductState*) -- The measured state, if ``get != 'outcome'``. .. py:attribute:: measure_ .. py:method:: sample_configuration(seed=None, info=None) Sample a configuration from this MPS. :param seed: A random seed or generator to use. :type seed: None, int, or np.random.Generator, optional :param info: If given, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. :type info: dict, optional .. py:method:: sample(C, seed=None, info=None) Generate ``C`` samples rom this MPS, along with their probabilities. :param C: The number of samples to generate. :type C: int :param seed: A random seed or generator to use. :type seed: None, int, or np.random.Generator, optional :param info: If given, will be used to infer and store various extra information. Currently the key "cur_orthog" is used to store the current orthogonality center. :type info: dict, optional :Yields: * **config** (*sequence of int*) -- The sample configuration. * **omega** (*float*) -- The probability of this configuration. .. py:class:: MatrixProductOperator(arrays, *, sites=None, L=None, shape='lrud', tags=None, upper_ind_id='k{}', lower_ind_id='b{}', site_tag_id='I{}', **tn_opts) Bases: :py:obj:`TensorNetwork1DOperator`, :py:obj:`TensorNetwork1DFlat` Initialise a matrix product operator, with auto labelling and tagging. :param arrays: The tensor arrays to form into a MPO. :type arrays: sequence of arrays :param sites: Construct the MPO on these sites only. If not given, enumerate from zero. Should be monotonically increasing and match ``arrays``. :type sites: sequence of int, optional :param L: The number of sites the MPO should be defined on. If not given, this is taken as the max ``sites`` value plus one (i.e.g the number of arrays if ``sites`` is not given). :type L: int, optional :param shape: 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. :type shape: str, optional :param tags: Global tags to attach to all tensors. :type tags: str or sequence of str, optional :param upper_ind_id: 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)))``. :type upper_ind_id: str :param lower_ind_id: 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)))``. :type lower_ind_id: str :param site_tag_id: 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)))``. :type site_tag_id: str .. py:attribute:: _EXTRA_PROPS :value: ('_site_tag_id', '_upper_ind_id', '_lower_ind_id', 'cyclic', '_L') .. py:attribute:: _L :value: None .. py:attribute:: _upper_ind_id :value: 'k{}' .. py:attribute:: _lower_ind_id :value: 'b{}' .. py:attribute:: _site_tag_id :value: 'I{}' .. py:attribute:: cyclic .. py:method:: from_fill_fn(fill_fn, L, bond_dim, phys_dim=2, sites=None, cyclic=False, shape='lrud', tags=None, upper_ind_id='k{}', lower_ind_id='b{}', site_tag_id='I{}') :classmethod: Create an MPO by supplying a 'filling' function to generate the data for each site. :param fill_fn: A function with signature ``fill_fn(shape : tuple[int]) -> array_like``. :type fill_fn: callable :param L: The number of sites. :type L: int :param bond_dim: The bond dimension. :type bond_dim: int :param phys_dim: The physical dimension(s) of each site, if a sequence it will be cycled over. :type phys_dim: int or Sequence[int], optional :param sites: Construct the MPO on these sites only. If not given, enumerate from zero. :type sites: None or sequence of int, optional :param cyclic: Whether the MPO should be cyclic (periodic). :type cyclic: bool, optional :param shape: 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. :type shape: str, optional :param tags: Global tags to attach to all tensors. :type tags: str or sequence of str, optional :param upper_ind_id: A string specifiying how to label the upper physical site indices. Should contain a ``'{}'`` placeholder. :type upper_ind_id: str :param lower_ind_id: A string specifiying how to label the lower physical site indices. Should contain a ``'{}'`` placeholder. :type lower_ind_id: str :param site_tag_id: How to tag the physical sites. Should contain a ``'{}'`` placeholder. :type site_tag_id: str, optional :rtype: MatrixProductState .. py:method:: from_dense(A, dims=2, sites=None, L=None, tags=None, site_tag_id='I{}', upper_ind_id='k{}', lower_ind_id='b{}', **split_opts) :classmethod: Build an MPO from a raw dense matrix. :param A: The dense operator, it should be reshapeable to ``(*dims, *dims)``. :type A: array :param dims: The physical subdimensions of the operator. If any integer, assume all sites have the same dimension. If a sequence, the dimension of each site. Default is 2. :type dims: int, sequence of int, optional :param sites: The sites to place the operator on. If None, will place it on first `len(dims)` sites. :type sites: sequence of int, optional :param L: The total number of sites in the MPO, if the operator represents only a subset. :type L: int, optional :param tags: Global tags to attach to all tensors. :type tags: str or sequence of str, optional :param site_tag_id: The string to use to label the site tags. :type site_tag_id: str, optional :param upper_ind_id: The string to use to label the upper physical indices. :type upper_ind_id: str, optional :param lower_ind_id: The string to use to label the lower physical indices. :type lower_ind_id: str, optional :param split_opts: Supplied to :func:`~quimb.tensor.tensor_core.tensor_split`. :rtype: MatrixProductOperator .. py:method:: fill_empty_sites(mode='full', phys_dim=None, fill_array=None, inplace=False) Fill any empty sites of this MPO with identity tensors, adding size 1 bonds or draping existing bonds where necessary such that the resulting tensor has nearest neighbor bonds only. :param mode: Whether to fill in all sites, including at either end, or simply the minimal range covering the min to max current sites present. :type mode: {'full', 'minimal'}, optional :param phys_dim: The physical dimension of the identity tensors to add. If not specified, will use the upper physical dimension of the first present site. :type phys_dim: int, optional :param fill_array: The array to use for the identity tensors. If not specified, will use the identity array of the same dtype as the first present site. :type fill_array: array, optional :param inplace: Whether to perform the operation inplace. :type inplace: bool, optional :returns: The modified MPO. :rtype: MatrixProductOperator .. py:attribute:: fill_empty_sites_ .. py:method:: add_MPO(other, inplace=False, **kwargs) .. py:attribute:: add_MPO_ .. py:method:: _apply_mps(other, compress=False, contract=True, **compress_opts) .. py:method:: _apply_mpo(other, compress=False, contract=True, **compress_opts) .. py:method:: apply(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). :param other: The object to act on. :type other: MatrixProductOperator or MatrixProductState :param compress: Whether to compress the resulting object. :type compress: bool, optional :param compress_opts: Supplied to :meth:`TensorNetwork1DFlat.compress`. :rtype: MatrixProductOperator or MatrixProductState .. py:attribute:: dot .. py:method:: permute_arrays(shape='lrud') 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. :param shape: A permutation of ``'lrud'`` specifying the desired order of the left, right, upper and lower (down) indices respectively. :type shape: str, optional .. py:method:: trace(left_inds=None, right_inds=None) Take the trace of this MPO. .. py:method:: partial_transpose(sysa, inplace=False) Perform the partial transpose on this MPO by swapping the bra and ket indices on sites in ``sysa``. :param sysa: The sites to transpose indices on. :type sysa: sequence of int or int :param inplace: Whether to perform the partial transposition inplace. :type inplace: bool, optional :rtype: MatrixProductOperator .. py:method:: rand_state(bond_dim, **mps_opts) Get a random vector matching this MPO. .. py:method:: identity(**mpo_opts) Get a identity matching this MPO. .. py:method:: show(max_width=None) .. py:class:: Dense1D(array, phys_dim=2, tags=None, site_ind_id='k{}', site_tag_id='I{}', **tn_opts) Bases: :py:obj:`TensorNetwork1DVector` 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. :param array: The full hilbert space vector - assumed to be made of equal hilbert spaces each of size ``phys_dim`` and will be reshaped as such. :type array: array_like :param phys_dim: The hilbert space size of each site, default: 2. :type phys_dim: int, optional :param tags: Extra tags to add to the tensor network. :type tags: sequence of str, optional :param site_ind_id: String formatter describing how to label the site indices. :type site_ind_id: str, optional :param site_tag_id: String formatter describing how to label the site tags. :type site_tag_id: str, optional :param tn_opts: Supplied to :class:`~quimb.tensor.tensor_core.TensorNetwork`. .. py:attribute:: _EXTRA_PROPS :value: ('_site_ind_id', '_site_tag_id', '_L') .. py:attribute:: _L .. py:attribute:: _site_ind_id :value: 'k{}' .. py:attribute:: _site_tag_id :value: 'I{}' .. py:method:: rand(n, phys_dim=2, dtype=float, **dense1d_opts) :classmethod: Create a random dense vector 'tensor network'. .. py:class:: 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) Bases: :py:obj:`TensorNetwork1D` 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) :param 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. :type arrays: sequence of arrays .. py:attribute:: _EXTRA_PROPS :value: ('_site_tag_id', '_outer_upper_ind_id', '_inner_upper_ind_id', '_inner_lower_ind_id',... .. py:attribute:: _L .. py:attribute:: _outer_upper_ind_id :value: 'kn{}' .. py:attribute:: _inner_upper_ind_id :value: 'k{}' .. py:attribute:: _inner_lower_ind_id :value: 'b{}' .. py:attribute:: _outer_lower_ind_id :value: 'bn{}' .. py:attribute:: _site_tag_id :value: 'I{}' .. py:attribute:: cyclic .. py:method:: rand(n, K, chi, phys_dim=2, herm=True, cyclic=False, dtype=complex, **superop_opts) :classmethod: .. py:property:: outer_upper_ind_id .. py:property:: inner_upper_ind_id .. py:property:: inner_lower_ind_id .. py:property:: outer_lower_ind_id .. py:class:: TNLinearOperator1D(tn, left_inds, right_inds, start, stop, ldims=None, rdims=None, is_conj=False, is_trans=False) Bases: :py:obj:`scipy.sparse.linalg.LinearOperator` 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 :class:`~quimb.tensor.tensor_core.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. :param tn: The tensor network to turn into a ``LinearOperator``. :type tn: TensorNetwork :param left_inds: The left indicies. :type left_inds: sequence of str :param right_inds: The right indicies. :type right_inds: sequence of str :param start: Index of starting site. :type start: int :param stop: Index of stopping site (does not include this site). :type stop: int :param ldims: If known, the dimensions corresponding to ``left_inds``. :type ldims: tuple of int, optional :param rdims: If known, the dimensions corresponding to ``right_inds``. :type rdims: tuple of int, optional .. seealso:: :py:obj:`TNLinearOperator` .. py:attribute:: tn .. py:attribute:: tags .. py:attribute:: is_conj :value: False .. py:attribute:: is_trans :value: False .. py:attribute:: _conj_linop :value: None .. py:attribute:: _adjoint_linop :value: None .. py:attribute:: _transpose_linop :value: None .. py:method:: _matvec(vec) Default matrix-vector multiplication handler. If self is a linear operator of shape (M, N), then this method will be called on a shape (N,) or (N, 1) ndarray, and should return a shape (M,) or (M, 1) ndarray. This default implementation falls back on _matmat, so defining that will define matrix-vector multiplication as well. .. py:method:: _matmat(mat) Default matrix-matrix multiplication handler. Falls back on the user-defined _matvec method, so defining that will define matrix multiplication (though in a very suboptimal way). .. py:method:: copy(conj=False, transpose=False) .. py:method:: conj() .. py:method:: _transpose() Default implementation of _transpose; defers to rmatvec + conj .. py:method:: _adjoint() Hermitian conjugate of this TNLO. .. py:method:: to_dense() .. py:method:: toarray() .. py:property:: A