quimb.experimental.operatorbuilder ================================== .. py:module:: quimb.experimental.operatorbuilder .. autoapi-nested-parse:: Tools for defining and constructing sparse operators with: * arbitrary geometries, * numba acceleration, * support for symmetic sectors, * efficient parallelization, and optionally producing: * sparse matrix form * matrix product operators, * dict of local gates form * VMC 'coupled configs' form Currently only supports composing operators which are sums of products of diagonal or anti-diagonal real dimension 2 operators. TODO:: - [ ] support for non-diagonal and qudit operators (lower priority) - [ ] product of operators generator (e.g. for PEPS DMRG) DONE:: - [x] fix sparse matrix being built in opposite direction - [x] complex and single precision support (lower priority) - [x] use compact bitbasis - [x] design interface for HilbertSpace / OperatorBuilder interaction - [x] automatic symbolic jordan wigner transformation - [x] numba accelerated coupled config - [x] general definition and automatic 'symbolic' jordan wigner - [x] multithreaded sparse matrix construction - [x] LocalHam generator (e.g. for simple update, normal PEPS algs) - [x] automatic MPO generator Submodules ---------- .. toctree:: :maxdepth: 1 /autoapi/quimb/experimental/operatorbuilder/builder/index /autoapi/quimb/experimental/operatorbuilder/configcore/index /autoapi/quimb/experimental/operatorbuilder/hilbertspace/index /autoapi/quimb/experimental/operatorbuilder/models/index /autoapi/quimb/experimental/operatorbuilder/pepobuilder/index Classes ------- .. autoapisummary:: quimb.experimental.operatorbuilder.HilbertSpace quimb.experimental.operatorbuilder.SparseOperatorBuilder Functions --------- .. autoapisummary:: quimb.experimental.operatorbuilder.fermi_hubbard_from_edges quimb.experimental.operatorbuilder.fermi_hubbard_spinless_from_edges quimb.experimental.operatorbuilder.heisenberg_from_edges quimb.experimental.operatorbuilder.rand_operator quimb.experimental.operatorbuilder.get_mat Package Contents ---------------- .. py:class:: HilbertSpace(sites, order=None, sector=None, symmetry=None) Take a set of 'sites' (any sequence of sortable, hashable objects), and map this into a 'register' or linearly indexed range, optionally using a particular ordering. A symmetry and sector can also be specified, which will change the size of the Hilbert space and how the valid configurations are enumerated. :param sites: The sites to map into a linear register. If an integer, simply use ``range(sites)``. :type sites: int or sequence of hashable objects :param order: If provided, use this to order the sites. If a callable, it should be a sorting key. If a sequence, it should be a permutation of the sites, and ``key=order.index`` will be used. :type order: callable or sequence of hashable objects, optional :param sector: The sector of the Hilbert space. If None, no sector is assumed. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space if any. If `None` and a `sector` is provided, the symmetry will be inferred from the sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional .. py:attribute:: _order :value: None .. py:attribute:: _sites .. py:attribute:: _mapping_inv .. py:attribute:: _mapping .. py:attribute:: _size :value: None .. py:attribute:: _pt :value: None .. py:method:: set_ordering(order=None) .. py:method:: from_edges(edges, order=None) :classmethod: Construct a HilbertSpace from a set of edges, which are pairs of sites. .. py:property:: sites The ordered tuple of all sites in the Hilbert space. .. py:property:: sector The sector of the Hilbert space. .. py:method:: get_sector_numba(sector=None, symmetry=None) The sector of the Hilbert space, in 'numba form'. A non-default symmetry and sector can be provided. :param sector: The sector of the Hilbert space. If None, the default sector is used. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space. If None, the default symmetry is used, or inferred from the supplied sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :returns: * **sector** (*ndarray[int64]*) -- The sector of the Hilbert space, in 'numba form'. This is a 1D array of length 1, 2 or 4, depending on the symmetry. * **symmetry** (*str*) -- The symmetry of the Hilbert space. This is one of "None", "Z2", "U1", or "U1U1". .. py:property:: symmetry The symmetry of the Hilbert space. .. py:property:: nsites The total number of sites in the Hilbert space. .. py:method:: get_pascal_table() Get a sufficiently large pascal table for this Hilbert space. .. py:method:: get_size(sector=None, symmetry=None) Get the size of the Hilbert space, optionally given a non-default symmetry and sector. :param sector: The sector of the Hilbert space. If None, the default sector is used. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space. If None, the default symmetry is used, or inferred from the supplied sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional .. py:property:: size Get the size of this Hilbert space, taking into account the symmetry and sector. .. py:method:: site_to_reg(site) Convert a site to a linear register index. .. py:method:: reg_to_site(reg) Convert a linear register index back to a site. .. py:method:: has_site(site) Check if this HilbertSpace contains a given site. .. py:method:: rank_to_flatconfig(rank) Convert a rank (linear index) into a flat configuration. :param rank: The rank (linear index) to convert. :type rank: int :returns: **flatconfig** -- A flat configuration, with the occupation number or spin state of each site in the order given by this ``HilbertSpace``. :rtype: ndarray[uint8] .. py:method:: flatconfig_to_rank(flatconfig) Convert a flat configuration into a rank (linear index). :param flatconfig: A flat configuration, with the occupation number or spin state of each site in the order given by this ``HilbertSpace``. :type flatconfig: ndarray[uint8] :returns: **rank** -- The rank (linear index) of the flat configuration in the Hilbert space. :rtype: int .. py:method:: config_to_flatconfig(config) Turn a configuration into a flat configuration, assuming the order given by this ``HilbertSpace``. :param config: A dictionary mapping sites to their occupation number / spin. :type config: dict[hashable, int] :returns: **flatconfig** -- A flat configuration, with the occupation number or spin state of each site in the order given by this ``HilbertSpace``. :rtype: ndarray[uint8] .. py:method:: flatconfig_to_config(flatconfig) Turn a flat configuration into a configuration, assuming the order given by this ``HilbertSpace``. :param flatconfig: A flat configuration, with the occupation number or spin state of each site in the order given by this ``HilbertSpace``. :type flatconfig: ndarray[uint8] :returns: **config** -- A dictionary mapping sites to their occupation number / spin state. :rtype: dict[hashable, int] .. py:method:: rank_to_config(rank) Convert a rank (linear index) into a configuration. :param rank: The rank (linear index) to convert. :type rank: int :returns: **config** -- A dictionary mapping sites to their occupation number / spin state. :rtype: dict[hashable, int] .. py:method:: config_to_rank(config) Convert a configuration into a rank (linear index). :param config: A dictionary mapping sites to their occupation number / spin state. :type config: dict[hashable, int] :returns: **rank** -- The rank (linear index) of the configuration in the Hilbert space. :rtype: int .. py:method:: rand_rank(seed=None) Get a random rank (linear index) in the Hilbert space. :param seed: The random seed or generator to use. If None, a new generator will be created with default settings. :type seed: None, int or numpy.random.Generator, optional :returns: **rank** -- A random rank in the Hilbert space. :rtype: int64 .. py:method:: rand_flatconfig(seed=None) Get a random flat configuration. :param seed: The random seed or generator to use. If None, a new generator will be created with default settings. :type seed: None, int or numpy.random.Generator, optional :returns: **flatconfig** -- A flat configuration, with the occupation number or spin state of each site in the order given by this ``HilbertSpace``. :rtype: ndarray[uint8] .. py:method:: rand_config(seed=None) Get a random configuration. :param seed: The random seed or generator to use. If None, a new generator will be created with default settings. :type seed: None, int or numpy.random.Generator, optional :returns: **config** -- A dictionary mapping sites to their occupation number / spin state. :rtype: dict[hashable, np.uint8] .. py:method:: __repr__() .. py:function:: fermi_hubbard_from_edges(edges, t=1.0, U=1.0, mu=0.0, order=None, sector=None, symmetry=None, hilbert_space=None, dtype=None) Create a Fermi-Hubbard Hamiltonian on the graph defined by ``edges``. The Hamiltonian is given by: .. math:: H = -t \sum_{\{i,j\}}^{|E|} \sum_{\sigma \in \uparrow, \downarrow} \left( c_{\sigma,i}^\dagger c_{\sigma,j} + c_{\sigma,j}^\dagger c_{\sigma,i} \right) + U \sum_{i}^{|V|} n_{\uparrow,i} n_{\downarrow,i} - \mu \sum_{i}^{|V|} \left( n_{\uparrow,i} + n_{\downarrow,i} \right) where :math:`\{i,j\}` are the edges of the graph, and :math:`c_{\sigma,i}` is the fermionic annihilation operator acting on site :math:`i` with spin :math:`\sigma`. The Jordan-Wigner transformation is used to implement fermionic statistics. :param edges: The edges, as pairs of hashable 'sites', that define the graph. Multiple edges are allowed, but will be treated as a single edge. :type edges: Iterable[tuple[hashable, hashable]] :param t: The hopping amplitude. Default is 1.0. :type t: float, optional :param U: The on-site interaction strength. Default is 1.0. :type U: float, optional :param mu: The chemical potential. Default is 0.0. :type mu: float, optional :param order: If provided, use this to order the sites. If a callable, it should be a sorting key. If a sequence, it should be a permutation of the sites, and ``key=order.index`` will be used. :type order: callable or sequence of hashable objects, optional :param sector: The sector of the Hilbert space. If None, no sector is assumed. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space if any. If `None` and a `sector` is provided, the symmetry will be inferred from the sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :param hilbert_space: The Hilbert space to use. If not given, one will be constructed automatically from the edges. This overrides the ``order``, ``symmetry``, and ``sector`` parameters. :type hilbert_space: HilbertSpace, optional :param dtype: The data type of the Hamiltonian. If None, a default dtype will be used, np.float64 for real and np.complex128 for complex. :type dtype: {None, str, type}, optional :returns: **H** -- The Hamiltonian as a SparseOperatorBuilder object. :rtype: SparseOperatorBuilder .. py:function:: fermi_hubbard_spinless_from_edges(edges, t=1.0, V=0.0, mu=0.0, order=None, sector=None, symmetry=None, hilbert_space=None, dtype=None) Create a spinless Fermi-Hubbard Hamiltonian on the graph defined by ``edges``. The Hamiltonian is given by: .. math:: H = -t \sum_{\{i,j\}}^{|E|} \left( c_i^\dagger c_j + c_j^\dagger c_i \right) + V \sum_{\{i,j\}}^{|E|} n_i n_j - \mu \sum_{i}^{|V|} n_i where :math:`\{i,j\}` are the edges of the graph, and :math:`c_i` is the fermionic annihilation operator acting on site :math:`i`. The Jordan-Wigner transformation is used to implement fermionic statistics. :param edges: The edges, as pairs of hashable 'sites', that define the graph. Multiple edges are allowed, but will be treated as a single edge. :type edges: Iterable[tuple[hashable, hashable]] :param t: The hopping amplitude. Default is 1.0. :type t: float, optional :param V: The nearest neighbor interaction strength. Default is 0.0. :type V: float, optional :param mu: The chemical potential. Default is 0.0. :type mu: float, optional :param order: If provided, use this to order the sites. If a callable, it should be a sorting key. If a sequence, it should be a permutation of the sites, and ``key=order.index`` will be used. :type order: callable or sequence of hashable objects, optional :param sector: The sector of the Hilbert space. If None, no sector is assumed. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space if any. If `None` and a `sector` is provided, the symmetry will be inferred from the sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :param hilbert_space: The Hilbert space to use. If not given, one will be constructed automatically from the edges. This overrides the ``order``, ``symmetry``, and ``sector`` parameters. :type hilbert_space: HilbertSpace, optional :param dtype: The data type of the Hamiltonian. If None, a default dtype will be used, np.float64 for real and np.complex128 for complex. :type dtype: {None, str, type}, optional :returns: **H** -- The Hamiltonian as a SparseOperatorBuilder object. :rtype: SparseOperatorBuilder .. py:function:: heisenberg_from_edges(edges, j=1.0, b=0.0, order=None, sector=None, symmetry=None, hilbert_space=None, dtype=None) Create a Heisenberg Hamiltonian on the graph defined by ``edges``. .. math:: H = \sum_{\{i,j\}}^{|E|} \left( J_x S^x_i S^x_j + J_y S^y_i S^y_j + J_z S^z_i S^z_j \right) + \sum_{i}^{|V|} \left( B_x S^x_i + B_y S^y_i + B_z S^z_i \right) where :math:`\{i,j\}` are the edges of the graph, and :math:`S^x_i` is the spin-1/2 operator acting on site :math:`i` in the x-direction, etc. Note positive values of :math:`J` correspond to antiferromagnetic coupling here, and the magnetic field is in the z-direction by default. :param edges: The edges, as pairs of hashable 'sites', that define the graph. Multiple edges are allowed, and will be treated as a single edge. :type edges: Iterable[tuple[hashable, hashable]] :param j: The Heisenberg exchange coupling constant(s). If a single float is given, it is used for all three terms. If a tuple of three floats is given, they are used for the xx, yy, and zz terms respectively. Note that positive values of ``j`` correspond to antiferromagnetic coupling. :type j: float or tuple[float, float, float], optional :param b: The magnetic field strength(s). If a single float is given, it is used taken as a z-field. If a tuple of three floats is given, they are used for the x, y, and z fields respectively. :type b: float or tuple[float, float, float], optional :param order: If provided, use this to order the sites. If a callable, it should be a sorting key. If a sequence, it should be a permutation of the sites, and ``key=order.index`` will be used. :type order: callable or sequence of hashable objects, optional :param sector: The sector of the Hilbert space. If None, no sector is assumed. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space if any. If `None` and a `sector` is provided, the symmetry will be inferred from the sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :param hilbert_space: The Hilbert space to use. If not given, one will be constructed automatically from the edges. This overrides the ``order``, ``symmetry``, and ``sector`` parameters. :type hilbert_space: HilbertSpace, optional :param dtype: The data type of the Hamiltonian. If None, a default dtype will be used, np.float64 for real and np.complex128 for complex. :type dtype: {None, str, type}, optional :returns: **H** -- The Hamiltonian as a SparseOperatorBuilder object. :rtype: SparseOperatorBuilder .. py:function:: rand_operator(n, m, k, kmin=None, seed=None, ops='XYZ') Generate a random operator with n qubits and m terms. Each term is a sum of k operators acting on different qubits. The operators are chosen randomly from the set {X, Y, Z, +, -, n}. The coefficients are drawn from a normal distribution. :param n: The number of qubits. :type n: int :param m: The number of terms in the operator. :type m: int :param k: The number of operators in each term. :type k: int :param kmin: The minimum number of operators in each term. If not given, kmin = k. :type kmin: int, optional :param seed: The random seed for reproducibility. :type seed: int, optional :param ops: The set of operators to choose from. :type ops: str, optional :returns: The random operator as a SparseOperatorBuilder object. :rtype: SparseOperatorBuilder .. py:function:: get_mat(op, dtype=None) .. py:class:: SparseOperatorBuilder(terms=(), hilbert_space: quimb.experimental.operatorbuilder.hilbertspace.HilbertSpace = None, dtype=None) Object for building operators with sparse structure. Specifically, a sum of terms, where each term is a product of operators, where each of these local operators acts on a single site and has at most one entry per row. :param terms: The terms to initialize the builder with. ``add_term`` is simply called on each of these. :type terms: sequence, optional :param hilbert_space: The Hilbert space to build the operator in. If this is not supplied then a minimal Hilbert space will be constructed from the sites used, when required. The default symmetry and sector to build operators in is inherited from this Hilbert space, but can be overridden. :type hilbert_space: HilbertSpace :param dtype: A default data type for the coefficients of the operator. If not provided, will be automatically determined at building time based on the terms in the operator. If the operator is complex, will be set to ``np.complex128``. If the operator is real, will be set to ``np.float64``. Individual building methods can override this. :type dtype: numpy.dtype or str, optional .. py:attribute:: _sites_used .. py:attribute:: _hilbert_space :value: None .. py:attribute:: _terms_raw .. py:attribute:: _terms_transformed :value: None .. py:attribute:: _terms :value: None .. py:attribute:: _dtype :value: None .. py:attribute:: _coupling_maps .. py:attribute:: _cache .. py:property:: sites_used A tuple of the sorted coordinates/sites seen so far. .. py:property:: hilbert_space :type: quimb.experimental.operatorbuilder.hilbertspace.HilbertSpace The Hilbert space of the operator. Created from the sites seen so far if not supplied at construction. .. py:property:: nsites The total number of coordinates/sites seen so far. .. py:method:: site_to_reg(site) Get the register / linear index of coordinate ``site``. .. py:method:: reg_to_site(reg) Get the site of register / linear index ``reg``. .. py:property:: terms_raw A tuple of the raw terms seen so far, as a mapping from operator strings to coefficients. .. py:method:: _get_terms_current() .. py:method:: _reset_caches() .. py:method:: _finalize() .. py:method:: _get_terms_final() .. py:property:: terms A tuple of the, possibly transformed, terms seen so far. .. py:property:: nterms The total number of terms seen so far. .. py:property:: locality The locality of the operator, the maximum support of any term. .. py:property:: iscomplex Whether the operator has complex terms. .. py:method:: add_term(*coeff_ops, atol=1e-12) Add a term to the operator. :param coeff: The overall coefficient of the term. :type coeff: float, optional :param ops: The operators of the term, together with the sites they act on. Each term should be a pair of ``(operator, site)``, where ``operator`` can be: - ``'x'``, ``'y'``, ``'z'``: Pauli matrices - ``'sx'``, ``'sy'``, ``'sz'``: spin operators (i.e. scaled Pauli matrices) - ``'+'``, ``'-'``: creation/annihilation operators - ``'n'``, ``'sn'``, or ``'h'``: number, symmetric number (n - 1/2) and hole (1 - n) operators And ``site`` is a hashable object that represents the site that the operator acts on. :type ops: sequence of tuple[str, hashable] .. py:method:: __iadd__(term) .. py:method:: __isub__(term) .. py:method:: jordan_wigner_transform() Transform the terms in this operator by pre-prending pauli Z strings to all creation and annihilation operators. This is always applied directly to the raw terms, so that the any fermionic ordering is respected. Any further transformations (e.g. simplification or pauli decomposition) should thus be applied after this transformation. Note this doesn't decompose +, - into (X + iY) / 2 and (X - iY) / 2, it just prepends Z strings. Call `pauli_decompose` after this to get the full decomposition. The placement of the Z strings is defined by the ordering of the hilbert space, by default, the sorted order of the site labels. Calling this multiple times will reset the transformed terms and any simplifications or decompositions that have been applied to the transformed terms will be lost. .. py:method:: simplify(atol=1e-12) Simplify terms so that each term is a string of single operators acting on different sites. This assumes that operators on different sites commute, and thus is always performed after transformations such as Jordan-Wigner. .. py:method:: pauli_decompose(atol=1e-12, use_zx=False) Transform the terms in this operator by decomposing them into Pauli strings. .. py:method:: show(filler='.') Print an ascii representation of the terms in this operator. .. py:method:: get_dtype(dtype=None) Calculate the numpy data type of the operator to use. :param dtype: The data type of the coefficients. If not provided, will be automatically determined based on the terms in the operator. :type dtype: numpy.dtype or str, optional :returns: **dtype** -- The data type of the coefficients. :rtype: numpy.dtype .. py:method:: get_coupling_map(dtype=None) Build and cache the coupling map for the specified dtype. :param dtype: The data type of the coefficients. If not provided, will be automatically determined based on the terms in the operator. :type dtype: numpy.dtype, optional :returns: **coupling_map** -- The operator defined as tuple of flat arrays. :rtype: tuple[ndarray] .. py:method:: flatconfig_coupling(flatconfig, dtype=None) Get an array of other configurations coupled to the given individual ``flatconfig`` by this operator, and the corresponding coupling coefficients. This is for use with VMC for example. :param flatconfig: The linear array of the configuration to get the coupling for. :type flatconfig: array[np.uint8] :param dtype: The data type of the matrix. If not provided, will be automatically determined based on the terms in the operator. :type dtype: numpy.dtype, optional :returns: * **coupled_flatconfigs** (*ndarray[np.uint8]*) -- Each distinct flatconfig coupled to ``flatconfig``. * **coeffs** (*ndarray[dtype]*) -- The corresponding coupling coefficients. .. py:method:: config_coupling(config, dtype=None) Get a list of other configurations coupled to ``config`` by this operator, and the corresponding coupling coefficients. This is for use with VMC for example. :param config: The configuration to get the coupling for. :type config: dict[site, int] :returns: * **coupled_configs** (*list[dict[site, np.uint8]]*) -- Each distinct configuration coupled to ``config``. * **coeffs** (*list[dtype]*) -- The corresponding coupling coefficients. .. py:method:: build_coo_data(sector=None, symmetry=None, dtype=None, parallel=False) Build the raw data for a sparse matrix in COO format, optionally in parallel. :param sector: The sector of the Hilbert space. If None, the default sector is used. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space. If None, the default symmetry is used, or inferred from the supplied sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :param dtype: The data type of the matrix. If not provided, will be automatically determined based on the terms in the operator. :type dtype: numpy.dtype, optional :param parallel: Whether to build the matrix in parallel (multi-threaded). If True, will use number of threads equal to the number of available CPU cores. If False, will use a single thread. If an integer is provided, it will be used as the number of threads to use. :type parallel: bool or int, optional :returns: * **data** (*array*) -- The data entries for the sparse matrix in COO format. * **rows** (*array*) -- The row indices for the sparse matrix in COO format. * **cols** (*array*) -- The column indices for the sparse matrix in COO format. * **d** (*int*) -- The total number of basis states. .. py:method:: build_sparse_matrix(sector=None, symmetry=None, dtype=None, stype='csr', parallel=False) Build a sparse matrix in the given format. Optionally in parallel. :param sector: The sector of the Hilbert space. If None, the default sector is used. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space. If None, the default symmetry is used, or inferred from the supplied sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :param dtype: The data type of the matrix. If not provided, will be automatically determined based on the terms in the operator. :type dtype: numpy.dtype, optional :param stype: The sparse matrix format to use. Can be one of 'coo', 'csr', 'csc', 'bsr', 'lil', 'dok', or 'dia'. Default is 'csr'. :type stype: str, optional :param parallel: Whether to build the matrix in parallel (multi-threaded). :type parallel: bool, optional :rtype: scipy.sparse matrix .. py:method:: build_dense(sector=None, symmetry=None, dtype=None, parallel=False) Get the dense (`numpy.ndarray`) matrix representation of this operator. :param sector: The sector of the Hilbert space. If None, the default sector is used. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space. If None, the default symmetry is used, or inferred from the supplied sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :param dtype: The data type of the matrix. If not provided, will be automatically determined based on the terms in the operator. :type dtype: numpy.dtype, optional :param parallel: Whether to build the matrix in parallel (multi-threaded). If True, will use number of threads equal to the number of available CPU cores. If False, will use a single thread. If an integer is provided, it will be used as the number of threads to use. :type parallel: bool or int, optional :returns: **A** -- The dense matrix representation of this operator. :rtype: numpy.ndarray .. py:method:: matvec(x, out=None, sector=None, symmetry=None, dtype=None, parallel=False) Apply this operator lazily (i.e. without constructing a sparse matrix) to a vector. This uses less memory but is much slower. :param x: The vector to apply the operator to. :type x: array :param out: An array to store the result in. If not provided, a new array will be created. :type out: array, optional :param sector: The sector of the Hilbert space. If None, the default sector is used. The implicit size should match that of `x`. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space. If None, the default symmetry is used, or inferred from the supplied sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :param dtype: The data type of the matrix. If not provided, will be automatically set as the same as the input vector. :type dtype: numpy.dtype, optional :param parallel: Whether to apply the operator in parallel (multi-threaded). If True, will use number of threads equal to the number of available CPU cores. If False, will use a single thread. If an integer is provided, it will be used as the number of threads to use. Uses `num_threads` more memory but is faster. :type parallel: bool or int, optional :returns: **out** -- The result of applying the operator to the vector. :rtype: array .. py:method:: aslinearoperator(sector=None, symmetry=None, dtype=None, parallel=False) Get a `scipy.sparse.linalg.LinearOperator` for this operator. This is a lazy representation of the operator, which uses `matvec` to apply the operator to a vector. Less memory is required than constructing the full sparse matrix, but it is significantly slower. Note currently the operator is assumed to be hermitian. :param sector: The sector of the Hilbert space. If None, the default sector is used. :type sector: {None, str, int, ((int, int), (int, int))}, optional :param symmetry: The symmetry of the Hilbert space. If None, the default symmetry is used, or inferred from the supplied sector if possible. :type symmetry: {None, "Z2", "U1", "U1U1"}, optional :param dtype: The data type of the matrix. If not provided, will be automatically determined based on the terms in the operator. :type dtype: numpy.dtype, optional :param parallel: Whether to apply the operator in parallel (multi-threaded). If True, will use number of threads equal to the number of available CPU cores. If False, will use a single thread. If an integer is provided, it will be used as the number of threads to use. Uses `num_threads` more memory but is faster. :type parallel: bool or int, optional :returns: **Alo** -- The linear operator representation of this operator. :rtype: scipy.sparse.linalg.LinearOperator .. py:method:: build_local_terms(dtype=None) Get a dictionary of local terms, where each key is a sorted tuple of sites, and each value is the local matrix representation of the operator on those sites. For use with e.g. tensor network algorithms. Note terms acting on the same sites are summed together and the size of each local matrix is exponential in the locality of that term. :returns: **Hk** -- The local terms. :rtype: dict[tuple[hashable], numpy.ndarray] .. py:method:: build_state_machine_greedy(atol=1e-12) .. py:method:: draw_state_machine(method='greedy', figsize='auto', G=None) Draw the fintie state machine for this operator, as if buildling the MPO. .. py:method:: build_mpo(method='greedy', dtype=None, **mpo_opts) Build a matrix product operator (MPO) representation of this operator. :param method: The method to use for building the MPO. Currently only "greedy" is supported. :type method: str, optional :param dtype: The data type of the MPO. If not supplied, will be chosen automatically based on the terms in the operator. :type dtype: type, optional :param mpo_opts: Additional options to pass to the MPO constructor. See `MatrixProductOperator` for details. :type mpo_opts: keyword arguments :returns: **mpo** -- The MPO representation of this operator. :rtype: MatrixProductOperator .. py:method:: __repr__() .. py:method:: build_matrix_ikron(**ikron_opts) Build either the dense or sparse matrix of this operator via explicit calls to `ikron`. This is a slower but useful alternative testing method.