quimb.tensor.tensor_dmrg ======================== .. py:module:: quimb.tensor.tensor_dmrg .. autoapi-nested-parse:: DMRG-like variational algorithms, but in tensor network language. Exceptions ---------- .. autoapisummary:: quimb.tensor.tensor_dmrg.DMRGError Classes ------- .. autoapisummary:: quimb.tensor.tensor_dmrg.MovingEnvironment quimb.tensor.tensor_dmrg.DMRG quimb.tensor.tensor_dmrg.DMRG1 quimb.tensor.tensor_dmrg.DMRG2 quimb.tensor.tensor_dmrg.DMRGX Functions --------- .. autoapisummary:: quimb.tensor.tensor_dmrg.get_default_opts quimb.tensor.tensor_dmrg.get_cyclic_canonizer quimb.tensor.tensor_dmrg.parse_2site_inds_dims Module Contents --------------- .. py:function:: get_default_opts(cyclic=False) Get the default advanced settings for DMRG. :returns: * **default_sweep_sequence** (*str*) -- How to sweep. Will be repeated, e.g. "RRL" -> RRLRRLRRL..., default: R. * **bond_compress_method** (*{'svd', 'eig', ...}*) -- Method used to compress sites after update. * **bond_compress_cutoff_mode** (*{'sum2', 'abs', 'rel'}*) -- How to perform compression truncation. * **bond_expand_rand_strength** (*float*) -- In DMRG1, strength of randomness to expand bonds with. Needed to avoid singular matrices after expansion. * **local_eig_tol** (*float*) -- Relative tolerance to solve inner eigenproblem to, larger = quicker but more unstable, default: 1e-3. Note this can be much looser than the overall tolerance, the starting point for each local solve is the previous state, and the overall accuracy comes from multiple sweeps. * **local_eig_ncv** (*int*) -- Number of inner eigenproblem lanczos vectors. Smaller can mean quicker. * **local_eig_backend** (*{None, 'AUTO', 'SCIPY', 'SLEPC'}*) -- Which to backend to use for the inner eigenproblem. None or 'AUTO' to choose best. Generally ``'SLEPC'`` best if available for large problems, but it can't currently handle ``LinearOperator`` Neff as well as ``'lobpcg'``. * **local_eig_maxiter** (*int*) -- Maximum number of inner eigenproblem iterations. * **local_eig_ham_dense** (*bool*) -- Force dense representation of the effective hamiltonian. * **local_eig_EPSType** (*{'krylovschur', 'gd', 'jd', ...}*) -- Eigensovler tpye if ``local_eig_backend='slepc'``. * **local_eig_norm_dense** (*bool*) -- Force dense representation of the effective norm. * **periodic_segment_size** (*float or int*) -- How large (as a proportion if float) to make the 'segments' in periodic DMRG. During a sweep everything outside this (the 'long way round') is compressed so the effective energy and norm can be efficiently formed. Tradeoff: longer segments means having to compress less, but also having a shorter 'long way round', meaning that it needs a larger bond to represent it and can be 'pseudo-orthogonalized' less effectively. 0.5 is the largest fraction that makes sense. Set to >= 1.0 to not use segmentation at all, which is better for small systems. * **periodic_compress_method** (*{'isvd', 'svds'}*) -- Which method to perform the transfer matrix compression with. * **periodic_compress_norm_eps** (*float*) -- Precision to compress the norm transfer matrix in periodic systems. * **periodic_compress_ham_eps** (*float*) -- Precision to compress the energy transfer matrix in periodic systems. * **periodic_compress_max_bond** (*int*) -- The maximum bond to use when compressing transfer matrices. * **periodic_nullspace_fudge_factor** (*float*) -- Factor to add to ``Heff`` and ``Neff`` to remove nullspace. * **periodic_canonize_inv_tol** (*float*) -- When psuedo-orthogonalizing, an inverse gauge is generated that can be very ill-conditioned. This factor controls cutting off the small singular values of the gauge to stop this. * **periodic_orthog_tol** (*float*) -- When psuedo-orthogonalizing, if the local norm is within this distance to 1 (pseudo-orthogonoalized), then the generalized eigen decomposition is *not* used, which is much more efficient. If set too large the total normalization can become unstable. .. py:class:: MovingEnvironment(tn, begin, bsz, *, cyclic=False, segment_callbacks=None, ssz=0.5, eps=1e-08, method='isvd', max_bond=-1, norm=False) Helper class for efficiently moving the effective 'environment' of a few sites in a 1D tensor network. E.g. for ``begin='left', bsz=2``, this initializes the right environments like so:: n - 1: ●─●─●─ ─●─●─● │ │ │ │ │ │ H─H─H─ ... ─H─H─H │ │ │ │ │ │ ●─●─●─ ─●─●─● n - 2: ●─●─●─ ─●─●─╮ │ │ │ │ │ ● H─H─H─ ... ─H─H─H │ │ │ │ │ ● ●─●─●─ ─●─●─╯ n - 3: ●─●─●─ ─●─╮ │ │ │ │ ●● H─H─H─ ... ─H─HH │ │ │ │ ●● ●─●─●─ ─●─╯ ... 0 : ●─●─╮ │ │ ●● ●●● H─H─HH...HHH │ │ ●● ●●● ●─●─╯ which can then be used to efficiently generate the left environments as each site is updated. For example if ``bsz=2`` and the environements have been shifted many sites into the middle, then ``MovingEnvironment()`` returns something like:: <---> bsz sites ╭─●─●─╮ ●●●●● │ │ ●●●●●●● HHHHH─H─H─HHHHHHH ●●●●● │ │ ●●●●●●● ╰─●─●─╯ 0 ... i i+1 ... n-1 For periodic systems ``MovingEnvironment`` approximates the 'long way round' transfer matrices. E.g consider replacing segment B (to arbitrary precision) with an SVD:: ╭───────────────────────────────────────────────╮ ╰─A─A─A─A─A─A─A─A─A─A─A─A─B─B─B─B─B─B─B─B─B─B─B─╯ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ ==> ╭─A─A─A─A─A─A─A─A─A─A─A─A─B─B─B─B─B─B─B─B─B─B─B─╮ ╰───────────────────────────────────────────────╯ ╭┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄╮ ┊ ╭─A─A─A─A─A─A─A─A─A─A─A─A─╮ ┊ ==> ╰┄┄╯ ╰─A─A─A─A─A─A─A─A─A─A─A─A─╯ ^ ^ segment_start segment_stop - 1 ╭┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄╮ ┊ ╭─A─A─╮ ┊ ==> ╰┄┄╯ ╰─A─A─╯ ... <-bsz-> ╭┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄╮ ┊ ╭─A─A─╮ ┊ ==> ╰~~╯ ╰─A─A─╯ i i+1 -----sweep---------> Can then contract and store left and right environments for efficient sweeping just as in non-periodic case. If the segment is long enough (50+) sites, often only 1 singular value is needed, and thus the efficiency is the same as for OBC. :param tn: The 1D tensor network, should be closed, i.e. an overlap of some sort. :type tn: TensorNetwork :param begin: Which side to start at and sweep from. :type begin: {'left', 'right'} :param bsz: The number of sites that form the the 'non-environment', e.g. 2 for DMRG2. :type bsz: int :param ssz: The size of the segment to use, if float, the proportion. Default: 1/2. :type ssz: float or int, optional :param eps: The tolerance to approximate the transfer matrix with. See :meth:`~quimb.tensor.TensorNetwork.replace_with_svd`. :type eps: float, optional :param cyclic: Whether this is a periodic ``MovingEnvironment``. :type cyclic: bool, optional :param segment_callbacks: Functions with signature ``callback(start, stop, self.begin)``, to be called every time a new segment is initialized. :type segment_callbacks: sequence of callable, optional :param method: How to perform the transfer matrix compression if PBC. See :meth:`~quimb.tensor.TensorNetwork.replace_with_svd`. :type method: {'isvd', 'svds', ...}, optional :param max_bond: If > 0, the maximum bond of the compressed transfer matrix. :type max_bond: , optional :param norm: If True, treat this ``MovingEnvironment`` as the state overlap, which enables a few extra checks. :type norm: bool, optional .. rubric:: Notes Does not necessarily need to be an operator overlap tensor network. Useful for any kind of sweep where only local tensor updates are being made. Note that *only* the current site is completely up-to-date and can be modified with changes meant to propagate. .. py:attribute:: tn .. py:attribute:: begin .. py:attribute:: bsz .. py:attribute:: cyclic :value: False .. py:attribute:: L .. py:attribute:: _site_tag_id .. py:method:: site_tag(i) .. py:method:: init_segment(begin, start, stop) Initialize the environments in ``range(start, stop)`` so that one can start sweeping from the side defined by ``begin``. .. py:method:: init_non_segment(start, stop) Compress and label the effective env not in ``range(start, stop)`` if cyclic, else just add some dummy left and right end pieces. .. py:method:: move_right() .. py:method:: move_left() .. py:method:: move_to(i) Move this effective environment to site ``i``. .. py:method:: __call__() Get the current environment. .. py:function:: get_cyclic_canonizer(k, b, inv_tol=1e-10) Get a function to use as a callback for ``MovingEnvironment`` that approximately orthogonalizes the segments of periodic MPS. .. py:function:: parse_2site_inds_dims(k, b, i) Sort out the dims and inds of:: ---O---O--- | | For use in 2 site algorithms. .. py:exception:: DMRGError Bases: :py:obj:`Exception` Common base class for all non-exit exceptions. .. py:class:: DMRG(ham, bond_dims, cutoffs=1e-09, bsz=2, which='SA', p0=None) Density Matrix Renormalization Group variational groundstate search. Some initialising arguments act as defaults, but can be overidden with each solve or sweep. See :func:`~quimb.tensor.tensor_dmrg.get_default_opts` for the list of advanced options initialized in the ``opts`` attribute. :param ham: The hamiltonian in MPO form. :type ham: MatrixProductOperator :param bond_dims: The bond-dimension of the MPS to optimize. If ``bsz > 1``, then this corresponds to the maximum bond dimension when splitting the effective local groundstate. If a sequence is supplied then successive sweeps iterate through, then repeate the final value. E.g. ``[16, 32, 64] -> (16, 32, 64, 64, 64, ...)``. :type bond_dims: int or sequence of ints. :param cutoffs: The cutoff threshold(s) to use when compressing. If a sequence is supplied then successive sweeps iterate through, then repeate the final value. E.g. ``[1e-5, 1e-7, 1e-9] -> (1e-5, 1e-7, 1e-9, 1e-9, ...)``. :type cutoffs: dict-like :param bsz: Number of sites to optimize for locally i.e. DMRG1 or DMRG2. :type bsz: {1, 2} :param which: Whether to search for smallest or largest real part eigenvectors. :type which: {'SA', 'LA'}, optional :param p0: If given, use as the initial state. :type p0: MatrixProductState, optional .. attribute:: state The current, optimized state. :type: MatrixProductState .. attribute:: energy The current most optimized energy. :type: float .. attribute:: energies The total energy after each sweep. :type: list of float .. attribute:: local_energies The local energies per sweep: ``local_energies[i, j]`` contains the local energy found at the jth step of the (i+1)th sweep. :type: list of list of float .. attribute:: total_energies The total energies per sweep: ``local_energies[i, j]`` contains the total energy after the jth step of the (i+1)th sweep. :type: list of list of float .. attribute:: opts Advanced options e.g. relating to the inner eigensolve or compression, see :func:`~quimb.tensor.tensor_dmrg.get_default_opts`. :type: dict .. attribute:: (bond_sizes_ham) If cyclic, the sizes of the energy environement transfer matrix bonds, per segment, per sweep. :type: list[list[int]] .. attribute:: (bond_sizes_norm) If cyclic, the sizes of the norm environement transfer matrix bonds, per segment, per sweep. :type: list[list[int]] .. py:attribute:: L .. py:attribute:: phys_dim .. py:attribute:: bsz :value: 2 .. py:attribute:: which :value: 'SA' .. py:attribute:: cyclic .. py:attribute:: _b .. py:attribute:: ham .. py:attribute:: TN_energy .. py:attribute:: energies :value: [] .. py:attribute:: local_energies :value: [] .. py:attribute:: total_energies :value: [] .. py:attribute:: opts .. py:method:: _set_bond_dim_seq(bond_dims) .. py:method:: _set_cutoff_seq(cutoffs) .. py:property:: energy .. py:property:: state .. py:method:: _canonize_after_1site_update(direction, i) Compress a site having updated it. Also serves to move the orthogonality center along. .. py:method:: _eigs(A, B=None, v0=None) Find single eigenpair, using all the internal settings. .. py:method:: print_energy_info(Heff=None, loc_gs=None) .. py:method:: print_norm_info(i=None) .. py:method:: form_local_ops(i, dims, lix, uix) Construct the effective Hamiltonian, and if needed, norm. .. py:method:: post_check(i, Neff, loc_gs, loc_en, loc_gs_old) Perform some checks on the output of the local eigensolve. .. py:method:: _update_local_state_1site(i, direction, **compress_opts) Find the single site effective tensor groundstate of:: >->->->->-/|\-<-<-<-<-<-<-<-< /|\ <-- uix | | | | | | | | | | | | | | / | \ H-H-H-H-H--H--H-H-H-H-H-H-H-H = L--H--R | | | | | i| | | | | | | | | \i| / >->->->->-\|/-<-<-<-<-<-<-<-< \|/ <-- lix And insert it back into the states ``k`` and ``b``, and thus ``TN_energy``. .. py:method:: _update_local_state_2site(i, direction, **compress_opts) Find the 2-site effective tensor groundstate of:: >->->->->-/| |\-<-<-<-<-<-<-<-< /| |\ | | | | | | | | | | | | | | | / | | \ H-H-H-H-H--H-H--H-H-H-H-H-H-H-H = L--H-H--R | | | | | i i+1| | | | | | | | \ | | / >->->->->-\| |/-<-<-<-<-<-<-<-< \| |/ i i+1 And insert it back into the states ``k`` and ``b``, and thus ``TN_energy``. .. py:method:: _update_local_state(i, **update_opts) Move envs to site ``i`` and dispatch to the correct local updater. .. py:method:: sweep(direction, canonize=True, verbosity=0, **update_opts) Perform a sweep of optimizations, either rightwards:: optimize --> ... >->-o-<-<-<-<-<-<-<-<-<-<-<-<-< | | | | | | | | | | | | | | | | H-H-H-H-H-H-H-H-H-H-H-H-H-H-H-H | | | | | | | | | | | | | | | | >->-o-<-<-<-<-<-<-<-<-<-<-<-<-< or leftwards (`direction='L'`):: <-- optimize ... >->->->->->->->->->->->->-o-<-< | | | | | | | | | | | | | | | | H-H-H-H-H-H-H-H-H-H-H-H-H-H-H-H | | | | | | | | | | | | | | | | >->->->->->->->->->->->->-o-<-< After the sweep the state is left or right canonized respectively. :param direction: Sweep from left to right (->) or right to left (<-) respectively. :type direction: {'R', 'L'} :param canonize: Canonize the state first, not needed if doing alternate sweeps. :type canonize: bool, optional :param verbosity: Show a progress bar for the sweep. :type verbosity: {0, 1, 2}, optional :param update_opts: Supplied to ``self._update_local_state``. .. py:method:: sweep_right(canonize=True, verbosity=0, **update_opts) .. py:method:: sweep_left(canonize=True, verbosity=0, **update_opts) .. py:method:: _print_pre_sweep(i, direction, max_bond, cutoff, verbosity=0) Print this before each sweep. .. py:method:: _compute_post_sweep() Compute this after each sweep. .. py:method:: _print_post_sweep(converged, verbosity=0) Print this after each sweep. .. py:method:: _check_convergence(tol) By default check the absolute change in energy. .. py:method:: solve(tol=0.0001, bond_dims=None, cutoffs=None, sweep_sequence=None, max_sweeps=10, verbosity=0, suppress_warnings=True) Solve the system with a sequence of sweeps, up to a certain absolute tolerance in the energy or maximum number of sweeps. :param tol: The absolute tolerance to converge energy to. :type tol: float, optional :param bond_dims: Overide the initial/current bond_dim sequence. :type bond_dims: int or sequence of int :param cutoffs: Overide the initial/current cutoff sequence. :type cutoffs: float of sequence of float :param sweep_sequence: String made of 'L' and 'R' defining the sweep sequence, e.g 'RRL'. The sequence will be repeated until ``max_sweeps`` is reached. :type sweep_sequence: str, optional :param max_sweeps: The maximum number of sweeps to perform. :type max_sweeps: int, optional :param verbosity: How much information to print about progress. :type verbosity: {0, 1, 2}, optional :param suppress_warnings: Whether to suppress warnings about non-convergence, usually due to the intentional low accuracy of the inner eigensolve. :type suppress_warnings: bool, optional :returns: **converged** -- Whether the algorithm has converged to ``tol`` yet. :rtype: bool .. py:class:: DMRG1(ham, which='SA', bond_dims=None, cutoffs=1e-08, p0=None) Bases: :py:obj:`DMRG` Simple alias of one site ``DMRG``. .. py:class:: DMRG2(ham, which='SA', bond_dims=None, cutoffs=1e-08, p0=None) Bases: :py:obj:`DMRG` Simple alias of two site ``DMRG``. .. py:class:: DMRGX(ham, p0, bond_dims, cutoffs=1e-08, bsz=1) Bases: :py:obj:`DMRG` Class implmenting DMRG-X [1], whereby local effective energy eigenstates are chosen to maximise overlap with the previous step's state, leading to convergence on an mid-spectrum eigenstate of the full hamiltonian, as long as it is perturbatively close to the original state. [1] Khemani, V., Pollmann, F. & Sondhi, S. L. Obtaining Highly Excited Eigenstates of Many-Body Localized Hamiltonians by the Density Matrix Renormalization Group Approach. Phys. Rev. Lett. 116, 247204 (2016). :param ham: The hamiltonian in MPO form, should have ~area-law eigenstates. :type ham: MatrixProductOperator :param p0: The initial MPS guess, e.g. a computation basis state. :type p0: MatrixProductState :param bond_dims: See :class:`DMRG`. :type bond_dims: int or sequence of int :param cutoffs: See :class:`DMRG`. :type cutoffs: float or sequence of float .. attribute:: k The current, optimized state. :type: MatrixProductState .. attribute:: energies The list of energies after each sweep. :type: list of float .. py:attribute:: TN_energy2 .. py:attribute:: variances .. py:attribute:: _target_energy .. py:attribute:: opts .. py:property:: variance .. py:method:: form_local_ops(i, dims, lix, uix) Construct the effective Hamiltonian, and if needed, norm. .. py:method:: _update_local_state_1site_dmrgx(i, direction, **compress_opts) Like ``_update_local_state``, but re-insert all eigenvectors, then choose the one with best overlap with ``eff_ovlp``. .. py:method:: _update_local_state(i, **update_opts) Move envs to site ``i`` and dispatch to the correct local updater. .. py:method:: sweep(direction, canonize=True, verbosity=0, **update_opts) Perform a sweep of the algorithm. :param direction: Sweep from left to right (->) or right to left (<-) respectively. :type direction: {'R', 'L'} :param canonize: Canonize the state first, not needed if doing alternate sweeps. :type canonize: bool, optional :param verbosity: Show a progress bar for the sweep. :type verbosity: {0, 1, 2}, optional :param update_opts: Supplied to ``self._update_local_state``. .. py:method:: _compute_post_sweep() Compute this after each sweep. .. py:method:: _print_post_sweep(converged, verbosity=0) Print this after each sweep. .. py:method:: _check_convergence(tol) By default check the absolute change in energy.