Algebra with YASTN tensors#

Basic algebra operations with symmetric tensors#

See examples at Basic algebra operations.

Symmetric tensors can be added to and multiplied by a scalar through the usual operations +, -, *, and /. Element-wise raising to a power is done with the standard power operation **.

Simple element-wise operations#

Tensor.__abs__() → Tensor#

Return tensor with element-wise absolute values.

Can be called on tensor as abs(tensor).

Tensor.real() → Tensor#

Return tensor with imaginary part set to zero.

Note

Follows the behavior of the backend.real() when it comes to creating a new copy of the data or handling datatype dtype.

Tensor.imag() → Tensor#

Return tensor with real part set to zero.

Note

Follows the behavior of the backend.imag() when it comes to creating a new copy of the data or handling datatype dtype.

Tensor.sqrt() → Tensor#

Return tensor after applying element-wise square root for each tensor element.

Tensor.rsqrt(cutoff=0) → Tensor#

Return element-wise operation 1/sqrt(tensor).

The tensor elements with absolute value below the cutoff are set to zero.

Parameters:

cutoff (real scalar) – (element-wise) cutoff for inversion

Tensor.reciprocal(cutoff=0) → Tensor#

Return element-wise operation 1/tensor.

The tensor elements with absolute value below the cutoff are set to zero.

Parameters:

cutoff (real scalar) – (element-wise) cutoff for inversion

Tensor.exp(step=1.0) → Tensor#

Return element-wise exp(step * tensor).

Note

This applies only to non-empty blocks of tensor

Tensor.clip(a_min=None, a_max=None) → Tensor#

Return element-wise min(max(tensor, a_min), a_max), limiting tensor elements to the interval [a_min, a_max].

At least one of the bounds has to be provided; the other one is then not applied. For a_min > a_max, the upper bound takes precedence and all elements are set to a_max.

Not supported for complex and bool tensors.

Note

This applies only to non-empty blocks of tensor. In particular, a_min does not raise the elements of structurally absent blocks, which remain zero.

Parameters:

a_min, a_max (real scalar) – lower and upper bound for tensor elements

Tensor.__mul__(number) → Tensor#

Multiply tensor by a number, use: number * tensor.

Tensor.__pow__(exponent) → Tensor#

Element-wise exponent of tensor, use: tensor ** exponent.

Tensor.__truediv__(number) → Tensor#

Divide tensor by a scalar, use: tensor / number.

yastn.Tensor.__add__(a, b) → Tensor[source]#

Add two tensors, use: \(a + b\).

Signatures and total charges of two tensors should match.

yastn.Tensor.__sub__(a, b) → Tensor[source]#

Subtract two tensors, use: \(a - b\).

Signatures and total charges of two tensors should match.

yastn.add(*tensors, amplitudes=None, lazy_threshold=None, **kwargs) → Tensor[source]#

Linear combination of tensors with given amplitudes, \(\sum_i amplitudes[i] tensors[i]\).

Parameters:

  • tensors (Sequence[yastn.Tensor]) – Signatures and total charges of all tensors should match.

  • amplitudes (None | Sequence[Number]) – If None, all amplitudes are assumed to be one. Otherwise, the number of tensors and amplitudes should be the same. Individual amplitude can be None, which gives the same result as 1 but without an extra multiplication.

  • lazy_threshold (None | float) – As in yastn.tensordot(): if the blocks filled by the tensors are fewer than this fraction of the blocks allowed by the merged legs, the result stores only the filled ones. None takes the value from the tensors’ config.

yastn.Tensor.__lt__(a, number) → Tensor[bool][source]#

Logical tensor with elements less-than a number (if it makes sense for backend data tensors), use: mask = tensor < number

Intended for diagonal tensor to be applied as a truncation mask.

yastn.Tensor.__gt__(a, number) → Tensor[bool][source]#

Logical tensor with elements greater-than a number (if it makes sense for backend data tensors), use: mask = tensor > number

Intended for diagonal tensor to be applied as a truncation mask.

yastn.Tensor.__le__(a, number) → Tensor[bool][source]#

Logical tensor with elements less-than-or-equal-to a number (if it makes sense for backend data tensors), use: mask = tensor <= number

Intended for diagonal tensor to be applied as a truncation mask.

yastn.Tensor.__ge__(a, number) → Tensor[bool][source]#

Logical tensor with elements greater-than-or-equal-to a number (if it makes sense for backend data tensors), use: mask = tensor >= number

Intended for diagonal tensor to be applied as a truncation mask.

yastn.Tensor.bitwise_not(a) → Tensor[bool][source]#

Return tensor after applying bitwise not on each tensor element.

Note

Operation applies only to non-empty blocks of tensor with tensor data dtype that allows for bitwise operation, i.e. intended for masks used to truncate tensor legs.

yastn.allclose(a, b, rtol=1e-13, atol=1e-13) → bool[source]#

Check if a and b are identical within a desired tolerance. To be True, all tensors’ blocks and merge history have to be identical. If this condition is satisfied, execute backend.allclose function to compare tensors’ data.

Note that if two tensors differ by zero blocks, the function returns False. To resolve such differences, use (a - b).norm() < tol

Parameters:

  • a, b (yastn.Tensor) – Tensor for comparison.

  • rtol, atol (float) – Desired relative and absolute precision.

Tensor contractions#

See examples at Tensor contractions.

Tensor contractions are the main building blocks of tensor network algorithms. The functions below facilitate the computation of

  • Trace: \(B_{jl}= \sum_{i} T_{ijil}\) or, using Einstein’s summation convention, repeated indices as \(B_{jl} = T_{ijil}\).

  • Contractions: in the usual form \(C_{abc} = A_{aijb}{\times}B_{cij}\) and also outer products \(M_{abkl} = A_{ak}{\times}B_{bl}\).

or composition of such operations over several tensors.

Tensor.__matmul__(b) → Tensor#

Compute the tensor dot product using the @ operator.

The operation contracts the last leg of the first tensor with the first leg of the second tensor. It is equivalent to yastn.tensordot(a, b, axes=(a.ndim - 1, 0)).

yastn.tensordot(a, b, axes, conj=(0, 0), lazy_threshold=None) → Tensor[source]#

Compute tensor dot product of two tensors along specified axes.

Outgoing legs are ordered such that first ones are the remaining legs of the first tensor in the original order, followed by the remaining legs of the second tensor in the original order.

Parameters:

  • a, b (yastn.Tensor) – Tensors to contract.

  • axes (tuple[int, int] | tuple[Sequence[int], Sequence[int]]) – legs of both tensors to be contracted (for each, they are specified by int or tuple of ints) e.g. axes=(0, 3) to contract 0th leg of a with 3rd leg of b; axes=((0, 3), (1, 2)) to contract legs 0 and 3 of a with 1 and 2 of b, respectively.

  • conj (tuple[int, int]) – specify tensors to conjugate by: (0, 0), (0, 1), (1, 0), or (1, 1). The default is (0, 0), i.e., neither tensor is conjugated.

yastn.vdot(a, b, conj=(1, 0)) → Number[source]#

Compute scalar product \(\langle a|b \rangle\) between two tensors.

Parameters:

  • a, b (yastn.Tensor) – Tensors to contract.

  • conj (tuple[int, int]) – indicate which tensors to conjugate: (0, 0), (0, 1), (1, 0), or (1, 1). The default is (1, 0), i.e., tensor a is conjugated.

yastn.broadcast(a, *args, axes=0, lazy_threshold=None) → 'Tensor' | tuple['Tensor'][source]#

Compute tensordot product of diagonal tensor a with tensors in args.

Produce diagonal tensor if both are diagonal. Legs of the resulting tensors are ordered in the same way as those of tensors in args. It is used (in combination with yastn.transpose()) as a subroutine of yastn.tensordot() for contractions involving diagonal tensor.

Parameters:

  • a, args (yastn.Tensor) – a is diagonal tensor to be broadcasted.

  • axes (int | Sequence[int]) – legs of tensors in args to be multiplied by diagonal tensor a. Number of tensors provided in args should match the length of axes.

yastn.apply_mask(a, *args, axes=0) → 'Tensor' | tuple['Tensor'][source]#

Apply mask given by nonzero elements of diagonal tensor a on specified axes of tensors in args. Number of tensors in args is not restricted. The length of the list axes has to be matching with args.

Legs of resulting tensor are ordered in the same way as those of tensors in args. Bond dimensions of specified axes of args are truncated according to the mask a. Produce diagonal tensor if both are diagonal.

Parameters:

  • a, args (yastn.Tensor) – a is a diagonal tensor

  • axes (int | Sequence[int]) – leg of tensors in args where the mask is applied.

yastn.trace(a, axes=(0, 1), lazy_threshold=None) → Tensor[source]#

Compute trace of legs specified by axes.

Parameters:

axes (tuple[int, int] | tuple[Sequence[int], Sequence[int]]) – Legs to be traced out, e.g., axes=(0, 1); or axes=((2, 3, 4), (0, 1, 5)).

yastn.einsum(subscripts, *operands, order=None, swap=None, charge_swap=None) → Tensor[source]#

Execute a series of tensor contractions.

This covers trace, tensordot (including outer products), and transpose operations. It follows the notation of np.einsum() as closely as possible.

Parameters:

  • subscripts (str)

  • operands (Sequence[yastn.Tensor])

  • order (str) – Specify order in which repeated indices from subscipt are contracted. By default it follows alphabetic order.

  • swap (str) – Comma-separated pairs of subscript characters identifying pairs of legs where swap gate is applied, e.g., swap='ab,cd'.

  • charge_swap (Sequence[tuple[str, Sequence[int]]]) – Pairs of a subscript character and a charge, e.g., [('a', (1,))]; see charge_swap in yastn.ncon().

Example

yastn.einsum('*ij,jh->ih', t1, t2)

# matrix-matrix multiplication, where the first matrix is conjugated.
# Equivalent to

t1.conj() @ t2

yastn.einsum('ab,al,bm->lm', t1, t2, t3, order='ba')

# Contract along `b` first, and `a` second.
yastn.ncon(ts, inds, conjs=None, order=None, swap=None, release_cuda_cache=False, oom_retry=False, charge_swap=None) → Tensor[source]#

Execute a series of tensor contractions.

Parameters:

  • ts (Sequence[yastn.Tensor]) – list of tensors to be contracted.

  • oom_retry (bool) – If True (torch/CUDA backends), each contraction op is retried once after torch.cuda.empty_cache() when it raises a CUDA out-of-memory error. This reclaims reserved-but-unallocated (fragmented) cache back to the driver so a large contiguous allocation can succeed. The blocking empty_cache only runs on the rare OOM path; the common path is untouched. Default: False.

  • inds (Sequence[Sequence[int]]) – each inner tuple labels legs of respective tensor with integers. Positive values label legs to be contracted, with pairs of legs to be contracted denoted by the same integer label. Non-positive numbers label legs of the resulting tensor, in reversed order, i.e. -0 for the first outgoing leg, -1 for the second, -2 for the third, etc.

  • swap (Sequence[Sequence[int]]) – Sequence of two-element tuples identifying pairs of legs where swap gate is applied.

  • charge_swap (Sequence[tuple[int, Sequence[int]]]) – Sequence of pairs (ind, charge): a swap gate between the leg labelled ind and a one-dimensional leg of fixed charge, e.g., a fermionic string crossing that leg (see charge in yastn.swap_gate()). The charge is read in the symmetry of the tensors, and the fermionic flag of their configuration selects the components that enter the sign. A contracted leg is named by its label; the gate acts on one of its two ends, which is equivalent. Pairs naming the same leg multiply.

  • conjs (Sequence[int]) – For each tensor in ts contains either 0 or 1. If the value is 1, the tensor is conjugated.

  • order (Sequence[int]) – Order in which legs, marked by positive indices in inds, are contracted. If None, the legs are contracted following an ascending indices order. The default is None.

Note

yastn.ncon() and yastn.einsum() differ only by syntax.

Example

# matrix-matrix multiplication where the first matrix is conjugated

yastn.ncon([a, b], ((-0, 1), (1, -1)), conjs=(1, 0))

# outer product

yastn.ncon([a, b], ((-0, -2), (-1, -3)))
yastn.swap_gate(a, axes, charge=None) → Tensor[source]#

Return tensor after application of a swap gate.

The function’s action is controlled by the fermionic flag in the tensor config. Multiply blocks with odd charges on swapped legs by \(-1\). The fermionic flag selects which individual charges (in case of a direct product of a few symmetries) are tested for oddity, where the contributions from each selected charge get multiplied. See yastn.operators.SpinfulFermions for an example. For fermionic=True, all charges are considered. For fermionic=False, swap_gate returns a.

Parameters:

  • axes (Sequence[int | Sequence[int]]) – Tuple with groups of legs. Consecutive pairs of grouped legs that are to be swapped. For instance, axes = (0, 1) apply swap gate between 0th and 1st leg. axes = ((0, 1), (2, 3), 4, 5) swaps (0, 1) with (2, 3), and 4 with 5.

  • charge (Optional[Sequence[int] | Sequence[Sequence[int]]]) – If provided, the swap gate is applied between a virtual one-dimensional leg of specified charge, e.g., a fermionic string, and tensor legs specified in axes. In this case, there is no application of a swap gates between legs specified in axes. One can provide list of charges corresponding to each axes, of a single charge to be applied to all axes.

yastn.fkron(*operators, sites=None, application_order=None)[source]#

Return a Kronecker product of operators, including fermionic swap-gates.

Parameters:

  • operators (yastn.Tensor) – a sequence of rank-2 tensors

  • sites (Sequence[int] | None) – sites corresponding to the provided operators. Should be a permutation of 0, 1, …, len(operators) - 1. If None, assume 0, 1, …, len(operators) - 1. Site 0 is the first in the fermionic order.

  • application_order (Sequence[int] | None) – Order of applying operators, which might correspond to a sign change for fermionic operators. Should be a permutation of 0, 1, …, len(operators) - 1. If None, the last operator is applied first.

Results#

Order of outgoing legs, where sites 0, 1, … go from left to right, e.g., fkron(A, B, C, sites=(0, 1, 2)) gives

   1     3     5
   |     |     |
┌──┴─────┴─────┴──┐
|  A     B     C  |
└──┬─────┬─────┬──┘
   |     |     |
   0     2     4

Note

Contractions are metadata-heavy: matching blocks between operands and planning the output is pure Python work that caching reuses across calls (tensordot_f2m, tensordot_fc, tensordot_nf, tensordot_cutensor_*, ncon, …). In a loop that contracts the same block structure repeatedly, these caches are what keeps the Python cost flat.

Depending on tensordot_policy, yastn.tensordot() also fuses legs into matrices before calling into the backend. On a CUDA device that fusion runs through the path described in GPU execution: hybrid scatter/loop, so the knobs documented there apply to contractions as well as to explicit fuse_legs() calls.

Transposition#

See examples at Transposition.

yastn.transpose(a, axes=None)[source]#

Transpose tensor by permuting the order of its legs (spaces). Do not copy tensor data.

Parameters:

axes (Sequence[int]) – new order of legs. Has to be a valid permutation of (0, 1, ..., ndim-1) where ndim is tensor order (number of legs). By default is range(a.ndim)[::-1], which reverses the order of the axes.

property Tensor.T: Tensor#

Same as self.transpose().

yastn.moveaxis(a, source, destination) → Tensor[source]#

Change the position of an axis (or a group of axes) of the tensor. This is a convenience function for subset of possible permutations. It computes the corresponding permutation and calls yastn.transpose().

Makes a shallow copy of tensor data if the order is not changed.

Parameters:

source, destination (int | Sequence[int])

property Tensor.H: Tensor#

Same as self.T.conj(), i.e., transpose and conjugate.

Fusion of legs (reshaping)#

See examples at Fusion (reshaping).

Fusion of several vector spaces \(V_1,V_2,\ldots,V_n\) creates a new vector space as direct product \(W=V_1 \otimes V_2 \otimes \ldots \otimes V_n\), which is then indexed by a single index of dimension \(\prod_i {\rm dim}(V_i)\). Here multiplication depends on abelian symmetry, as the resulting total dimension is a sum of dimensions for effective charges. The inverse operation can split the fused space into its original constituents.

For dense tensors, the operation corresponds to reshaping.

Fusion can be used to vary compression between (unfused) symmetric tensors with many small non-zero blocks and tensors with several fused spaces having just few, but large non-zero blocks.

Tensor.fuse_legs(axes, mode=None) → Tensor#

Fuse groups of legs into effective legs, reducing the rank of the tensor.

Note

Fusion can be reverted back by yastn.Tensor.unfuse_legs()

First, the legs are permuted into desired order. Then, selected groups of consecutive legs are fused. The desired order of the legs is given by a tuple axes of leg indices where the groups of legs to be fused are denoted by inner tuples

axes=(0,1,(2,3,4),5)  keep leg order, fuse legs (2,3,4) into new leg
->   (0,1, 2,     3)
    __              __
0--|  |--3      0--|  |--3<-5
1--|  |--4  =>  1--|__|
2--|__|--5          |
                    2<-(2,3,4)


axes=(2,(3,1),4,(7,6),5,0)  permute indices, then fuse legs (3,1) into new leg
->   (0, 1,   2, 3,   4,5)  and legs (7,6) into another new leg
    __                   __
0--|  |--4        0->5--|  |--2<-4
1--|  |--5 =>           |  |--4<-5
2--|  |--6        2->0--|  |--3<-(7,6)
3--|__|--7    (3,1)->1--|__|

Two types of fusion are supported: meta and hard:

  • 'hard' changes both the structure and data by aggregating smaller blocks into larger ones. Such fusion allows to balance number of non-zero blocks and typical block size.

  • 'meta' performs the fusion only at the level of syntax, where it operates as a tensor with lower rank. Tensor structure and data (blocks) are not affected - apart from a transpose that may be needed for consistency.

It is possible to use both meta and hard fusion of legs on the same tensor. Applying hard fusion on tensor turns all previous meta fused legs into hard fused ones.

Parameters:

  • axes (Sequence[int | Sequence[int]]) – tuple of leg indices. Groups of legs to be fused together are accumulated within inner tuples.

  • mode (str) – can select 'hard' or 'meta' fusion. If None, uses default_fusion from tensor’s configuration. Configuration option force_fusion can be used to override mode (introduced for debugging purposes).

Tensor.unfuse_legs(axes) → Tensor#

Unfuse legs, reverting one layer of fusion.

If the tensor has been obtained by fusing some legs together, unfuse_legs can revert such fusion. The legs to be unfused are passed in axes as int or tuple[int] in case of more legs to be unfused. The unfused legs are inserted at the positions of the fused legs. The remaining legs are shifted accordingly

axes=2              unfuse leg 2 into legs 2,3,4
->   (0,1,2,3,4,5)
    __                 __
0--|  |--3         0--|  |--3
1--|__|        =>  1--|  |--4
    |              2--|__|--5<-3
    2=(2,3,4)


axes=(    2,      5  )  unfuse leg 2 into legs 2,3 and leg 5 into legs 6,7
->   (0,1,2,3,4,5,6,7)
          __                   __
      0--|  |--3           0--|  |--4
         |  |--4       =>  1--|  |--5
      1--|  |--5=(6,7)     2--|  |--6
(2,3)=2--|__|              3--|__|--7

Unfusing a leg obtained by fusing together other previously fused legs, unfuses only the last fusion.

axes=2              unfuse leg 2 into legs 2,3
->   (0,1,2,3,4)
    __                    __
0--|  |--3            0--|  |--3=(3,4)
1--|__|           =>  1--|  |
    |                 2--|__|--4<-3
    2=(2,3=(3,4))

fuse_legs may involve leg transposition, which is not undone by unfuse_legs.

Parameters:

axes (int | Sequence[int]) – leg(s) to unfuse.

Tensor.add_leg(axis=-1, s=-1, t=None, leg=None) → Tensor#

Creates a new tensor with an extra leg that carries the charge (or part of it) of the orignal tensor. This is achieved by the extra leg having a single charge sector of dimension D=1. The total charge of the tensor yastn.Tensor.n can be modified this way.

Makes a shallow copy of tensor data.

Parameters:

  • axis (int) – index of the new leg

  • s (int) – signature of the new leg, +1 or -1. The default is -1, where the leg charge is equal to the tensor charge for t=None.

  • t (int | Sequence[int]) – charge carried by the new leg. The default is None, which takes the total charge n of the original tensor resulting in a tensor with n=0.

  • leg (Optional[Leg]) – It is possible to provide a new leg directly. It has to be of dimension one but can contain information about the fusion of other dimension-one legs. If provided, it overrides information in s and t. The default is None.

Tensor.remove_leg(axis=-1) → Tensor#

Removes leg with a single charge sector of dimension one from tensor. The charge carried by that leg (if any) is added to the tensor’s total charge yastn.Tensor.n.

Makes a shallow copy of tensor data.

Parameters:

axis (int) – index of the leg to be removed.

Tensor.drop_leg_history(axes=None) → Tensor#

Drops information about original structure of fused or blocked legs that have been combined into a selected tensor leg(s).

Makes a shallow copy of tensor data.

Parameters:

axes (int | Sequence[int]) – index of the leg, or a group of legs. The default is None, which drops information from all legs.

GPU execution: hybrid scatter/loop#

Hard fusion (fuse_legs() with mode='hard') and its inverse unfuse_legs() are not just bookkeeping: they physically transpose and copy every block into (or out of) the fused buffer. On the torch backend this move has two implementations, and which one runs is chosen per call.

On CPU it is always a loop over blocks, one copy per block.

On GPU a per-block loop is a poor fit when a tensor has many small blocks: each block costs a kernel launch, and thousands of tiny launches dominate the actual data movement. The alternative — building an explicit index map and doing a single scatter/gather over the whole buffer — pays a one-off cost to construct that map, then moves everything in one kernel. Neither wins everywhere: the loop is bandwidth-optimal for large blocks (no index build at all), the scatter wins for many small ones.

YASTN therefore splits each call by block size rather than choosing globally. Blocks with at least YASTN_FUSE_SCATTER_THRESH elements go through the loop; smaller blocks are collected into one compact scatter (fusing) or gather (unfusing). Three regimes fall out:

  • all blocks large — pure per-block loop; no index map is built,

  • all blocks small — a single lean scatter/gather over the whole buffer,

  • mixed — the hybrid kernel: the large blocks loop while the small ones ride one scatter built over compact indices that cover only the real blocks.

Unfusing additionally falls back to the loop when the destination blocks do not tile the output buffer exactly.

The index maps are keyed on structure and cached, so a repeated fusion pattern builds them once — see Caching, entry pack_transpose_and_merge_params. That entry holds device-resident tensors; yastn.clear_cache() releases them.

Variable

Default

Effect

YASTN_FUSE_SCATTER_THRESH

65536 (2**16)

Per-block element count separating the loop (>= threshold) from the scatter/gather (< threshold). The default is roughly the point at which a single block saturates the memory interface of a datacentre GPU, which is also where the loop’s launch cost and the index build cost cross over. Lower it to push more blocks through the loop, raise it to push more through the scatter.

YASTN_FUSE_SCATTER_CHUNK

unset

Tile size for building the index map. Unset means a single tile of 2**27 elements. A positive value tiles the build, bounding its peak scratch memory. 0 is the escape hatch: force the per-block loop even on GPU, which is also the baseline to A/B against when checking whether the scatter path is helping.

Both variables are read on every call, so they can be changed at runtime without reimporting YASTN. They have no effect on CPU or on the NumPy backend. A negative or non-integer value is reported through warnings and the default is used instead — an invalid setting never raises.

All regimes produce bit-identical results; the choice is purely one of performance.

Conjugation of symmetric tensors#

See examples at Conjugation of symmetric tensors.

Tensor.conj() → Tensor#

Return conjugated tensor. In particular, change the sign of the signature s to -s, the total charge n to -n, and complex conjugate each block of the tensor.

Follows the behavior of the backend.conj() when it comes to creating a new copy of the data.

Tensor.conj_blocks() → Tensor#

Complex-conjugate all blocks leaving symmetry structure (signature, blocks charge, and total charge) unchanged.

Follows the behavior of the backend.conj() when it comes to creating a new copy of the data.

Tensor.flip_signature() → Tensor#

Change the signature of the tensor, s to -s or equivalently reverse the direction of in- and out-going legs, and also the total charge of the tensor n to -n. Does not complex-conjugate the elements of the tensor.

Creates a shallow copy of the data.

Tensor.flip_charges(axes=None) → Tensor#

Flip signs of charges and signatures on specified legs.

Flipping charges/signature of hard-fused legs is not supported.

Parameters:

axes (int | Sequence[int]) – index of the leg, or a group of legs. The default is None, which flips all legs.

Tensor norms#

yastn.linalg.norm(a, p='fro') → Number[source]#

Return the norm of the tensor.

Parameters:

p (str) – 'fro' for Frobenius norm; 'inf' for \(l^\infty\) (or supremum) norm.

Spectral decompositions and truncation#

See examples at Decompositions of symmetric tensors.

yastn.linalg.svd(a, axes=(0, 1), sU=1, nU=True, compute_uv=True, Uaxis=-1, Vaxis=0, policy='fullrank', fix_signs=False, svd_on_cpu=False, thresh=0.1, **kwargs) → tuple['Tensor', 'Tensor', 'Tensor'] | 'Tensor'[source]#

Split a tensor into \(a = U S V\) using an exact singular value decomposition (SVD), where the columns of U and the rows of V form orthonormal bases and S is a positive diagonal matrix.

Parameters:

  • axes (tuple[int, int] | tuple[Sequence[int], Sequence[int]]) – Specify two groups of legs between which to perform SVD, as well as their final order.

  • sU (int) – Signature of the new leg in U; equal to 1 or -1. The default is 1. V is going to have the opposite signature on the connecting leg.

  • nU (bool) – Whether or not to attach the charge of a to U. If False, it is attached to V. The default is True.

  • compute_uv (bool) – If True, compute and return U, S, V. If False, compute and return only S. The default is True.

  • Uaxis, Vaxis (int) – Specify which leg of U and V tensors are connecting with S. By default, it is the last leg of U and the first of V, in which case a = U @ S @ V.

  • policy (str) –

    Strategy for computing the SVD or a partial SVD.

    • (default) "fullrank" compute full SVD then truncate.

    • "lowrank" default policy for partial SVD.

      On NumPy backend uses block_arnoldi. On torch backend uses block_arnoldi.

    • "randomized" randomized SVD up to desired size in each block. Requires providing k_block in kwargs.

      Requires torch backend and uses torch.svd_lowrank.

    • "block_arnoldi" partial SVD using scipy’s svds arnoldi method. Requires providing k_block in kwargs.

    • "block_propack" partial SVD using scipy’s svds propack method. Requires providing k_block in kwargs.

    kwargs will be passed to those functions for non-default settings.

  • thresh (float) – For policy='block_arnoldi' or policy='block_propack', this sets the threshold on the minimal block size for applying a partial SVD solver instead of a full SVD. The default is thresh=0.1. If for a matrix of size \(N \times N\) N * thresh is smaller than the requested number of singular triples, a full SVD is applied.

  • fix_signs (bool) – Whether or not to fix phases in U and V, so that the largest element in each column of U is positive. Provide uniqueness of decomposition for non-degenerate cases. The default is False.

  • svd_on_cpu (bool) – GPU tends to be very slow when executing SVD. If True, the data will be copied to CPU for SVD, and the results will be copied back to the device. Nothing is done for data already residing on CPU. The default is False.

  • k_block (None (default) | int | dict) – When policy='lowrank', number of singular values to compute in each block. If D_block is provided, it is used instead to determine number of singular values to compute.

Return type:

U, S, V (when compute_uv=True) or S (when compute_uv=False)

yastn.linalg.svd_with_truncation(a, axes=(0, 1), sU=1, nU=True, Uaxis=-1, Vaxis=0, policy='fullrank', fix_signs=False, svd_on_cpu=False, tol=-inf, tol_block=-inf, D_total=inf, D_block=inf, largest_gap=False, eps_multiplet=None, hermitian=False, mask_f=None, **kwargs) → tuple['Tensor', 'Tensor', 'Tensor'][source]#

Split a tensor using an exact singular value decomposition (SVD) into \(a = U S V\), where the columns of U and the rows of V form orthonormal bases and S is a positive diagonal matrix.

The function allows optional truncation based on relative tolerance, block-wise bond dimension, and total bond dimension across all blocks (whichever gives the smaller total dimension).

Parameters:

  • axes (tuple[int, int] | tuple[Sequence[int], Sequence[int]]) – Specify two groups of legs between which to perform SVD, as well as their final order.

  • sU (int) – Signature of the new leg in U; equal to 1 or -1. The default is 1. V has the opposite signature on the connecting leg.

  • nU (bool) – Whether or not to attach the charge of a to U. If False, it is attached to V. The default is True.

  • Uaxis, Vaxis (int) – Specify which legs of U and V connect to S. By default, these are the last leg of U and the first leg of V.

  • policy (str) – "fullrank" or "lowrank" are allowed. For "fullrank" use a standard full (but reduced) SVD, while "lowrank" uses a randomized or truncated SVD and requires D_block or k_block in kwargs.

  • tol (float) – Relative tolerance with respect to the largest absolute value element of S.

  • tol_block (float) – Relative tolerance per block.

  • D_total (int) – Maximum number of elements kept across all blocks.

  • D_block (int | dict) – Maximum number of elements kept per block. It is also possible to provide a dictionary mapping charges to maximal number of elements in the charge sector.

  • largest_gap (bool) – If True, enlarge the truncation range specified by other arguments by shifting the cut to the largest gap between to-be-truncated singular values across all blocks. It provides a heuristic mechanism to avoid truncating part of a multiplet. If True, tol_block and D_block are ignored, as largest_gap is a global condition. The default is False.

  • eps_multiplet (float) – Relative tolerance on multiplet splitting. If relative difference between two consecutive elements of S is larger than eps_multiplet, these elements are not considered as part of the same multiplet. Partially truncated multiplets are truncated down. The default is None, when this scheme is not used. If True, tol_block and D_block are ignored, as eps_multiplet is a global condition. Cannot be used together with largest_gap scheme.

  • hermitian (bool) – If True, blocks related by hermitian conjugation are truncated equally, truncating down to the intersecting part. The default is False.

  • mask_f (None | function[yastn.Tensor] -> yastn.Tensor) – It is possible to provide a custom mask function, which provides a mechanism to pass such a function to many tensor network algorithms where the function truncation_mask is being called. If provided, it overrides the default function, and all other parameters are ignored. The default is None.

Return type:

U, S, V

yastn.linalg.qr(a, axes=(0, 1), sQ=1, Qaxis=-1, Raxis=0) → tuple['Tensor', 'Tensor'][source]#

Split a tensor using a reduced QR decomposition such that \(a = Q R\), with \(Q Q^\dagger = I\). The charge of R is zero, and the charge of a is carried by Q.

Parameters:

  • axes (tuple[int, int] | tuple[Sequence[int], Sequence[int]]) – Specify two groups of legs between which to perform QR, as well as their final order.

  • sQ (int) – signature of connecting leg in Q; equal 1 or -1. The default is 1. R is going to have opposite signature on connecting leg.

  • Qaxis, Raxis (int) – specify which leg of Q and R tensors are connecting to the other tensor. By default, it is the last leg of Q and the first leg of R.

Return type:

Q, R

yastn.linalg.eig(a, axes=(0, 1), sU=1, nU=True, compute_uv=True, Uaxis=-1, Vaxis=0, policy='fullrank', which='LM', **kwargs) → tuple['Tensor', 'Tensor', 'Tensor'] | 'Tensor'[source]#

Split a tensor into \(a = U S V\) using an exact eigenvalue decomposition (ED), where the columns of U and the rows of V satisfy biorthogonality, i.e. V @ U = I, and S is a diagonal matrix. Unlike in the symmetric/Hermitian case, U and V are not necessarily related and need not form orthonormal bases.

Parameters:

  • axes (tuple[int, int] | tuple[Sequence[int], Sequence[int]]) – Specify two groups of legs between which to perform ED, as well as their final order.

  • sU (int) – Signature of the new leg in U; equal to 1 or -1. The default is 1. V is going to have the opposite signature on the connecting leg.

  • nU (bool) – Whether or not to attach the charge of a to U. If False, it is attached to V. The default is True.

  • compute_uv (bool) – If True, compute and return U, S, V. If False, compute and return only S. The default is True.

  • Uaxis, Vaxis (int) – Specify which leg of U and V tensors are connecting with S. By default, it is the last leg of U and the first of V, in which case a = U @ S @ V.

  • policy (str) – "fullrank" or "lowrank" are allowed. Use standard ED for "fullrank". For "lowrank", uses dominant ED methods and requires providing D_block in kwargs. This employs scipy.sparse.linalg.eigs for numpy backend. kwargs will be passed to those functions for non-default settings.

  • which (str) – One of ['SR', 'LR'`, ``'SM', 'LM'] specifying how to order S: 'LM' : (default) sort by absolute value, largest first, 'SM' : sort by absolute value, smallest first, 'SR' : sort by real part, smallest first, 'LR' : sort by real part, largest first.

Return type:

U, S, V (when compute_uv=True) or S (when compute_uv=False)

yastn.linalg.eigh(a, axes, sU=1, Uaxis=-1, which='LR', policy='fullrank', **kwargs) → tuple['Tensor', 'Tensor'][source]#

Split a symmetric tensor using an exact eigenvalue decomposition, \(a = U S U^\dagger\).

The tensor is expected to be symmetric (Hermitian) with total charge 0.

Parameters:

  • axes (tuple[int, int] | tuple[Sequence[int], Sequence[int]]) – Specify two groups of legs between which to perform eigh, as well as their final order.

  • sU (int) – signature of connecting leg in U equal 1 or -1. The default is 1.

  • Uaxis (int) – specify which leg of U is the new connecting leg. By default, it is the last leg.

  • which (str) – One of ['SR', 'LR', 'SM', 'LM'] specifying how to order S: 'LM' : sort by absolute value, largest first, 'SM' : sort by absolute value, smallest first, 'SR' : (default) sort by real part, smallest first, 'LR' : sort by real part, largest first.

  • policy (str) – 'fullrank' : (default) use standard full eigenvalue decomposition. 'block_lanczos' : use partial eigenvalue decomposition via scipy.sparse.linalg.eigsh. Requires D_block or k_block in kwargs specifying the number of eigenvalues per block.

  • k_block (None (default) | int | dict) – When policy='block_lanczos', number of eigenvalues to compute in each block. If D_block is provided, it is used instead to determine number of eigenvalues to compute.

Return type:

S, U

yastn.linalg.eigh_with_truncation(a, axes, sU=1, Uaxis=-1, which='LR', policy='fullrank', tol=0, tol_block=0, D_block=inf, D_total=inf, largest_gap=False, mask_f=None, **kwargs) → tuple['Tensor', 'Tensor'][source]#

Split a symmetric tensor using an exact eigenvalue decomposition, \(a = U S U^\dagger\). Optionally truncate the resulting decomposition.

The tensor is expected to be symmetric (Hermitian) with total charge 0. Truncation can be based on relative tolerance, bond dimension of each block, and total bond dimension across all blocks (whichever gives smaller total dimension). Truncate based on tolerance only if some eigenvalues are positive – then all negative ones are discarded.

Parameters:

  • axes (tuple[int, int] | tuple[Sequence[int], Sequence[int]]) – Specify two groups of legs between which to perform eigh, as well as their final order.

  • sU (int) – signature of connecting leg in U equal 1 or -1. The default is 1.

  • Uaxis (int) – specify which leg of U is the new connecting leg. By default, it is the last leg.

  • which (str) – One of ['SR', 'LR'`, ``'SM', 'LM'] specifying how to order S: 'LM' : sort by absolute value, largest first, 'SM' : sort by absolute value, smallest first, 'SR' : (default) sort by real part, smallest first, 'LR' : sort by real part, largest first.

  • policy (str) – "fullrank" : Use standard full ED for "fullrank" and then truncate. kwargs will be passed to those functions for non-default settings.

  • tol (float) – relative tolerance of eigen-values below which to truncate across all blocks.

  • tol_block (float) – relative tolerance of eigen-values below which to truncate within individual blocks.

  • D_block (int) – largest number of eigen-values to keep in a single block.

  • D_total (int) – largest total number of eigen-values to keep.

  • mask_f (function[yastn.Tensor] -> yastn.Tensor) – custom truncation-mask function. If provided, it overrides all other truncation-related arguments.

Return type:

S, U

yastn.linalg.truncation_mask(S, which='LR', tol=-inf, tol_block=-inf, D_total=inf, D_block=inf, largest_gap=False, eps_multiplet=None, hermitian=False, mask_f=None, **kwargs) → Tensor[bool][source]#

Generate a mask tensor from a diagonal tensor S. The mask can then be used for truncation.

Parameters:

  • S (yastn.Tensor) – Diagonal tensor with spectrum.

  • which (str) – Which values to keep from ['LM', 'LR', 'SR', 'SM']: 'LR' : largest real part (the default), 'LM' : largest magnitude, 'SR' : smallest real part, 'SM' : smallest magnitude.

  • tol (float) – Relative tolerance with respect to the largest absolut value element of S.

  • tol_block (float) – Relative tolerance per block.

  • D_total (int) – Maximum number of elements kept across all blocks.

  • D_block (int | dict) – Maximum number of elements kept per block. It is also possible to provide a dictionary mapping charges to maximal number of elements in the charge sector.

  • largest_gap (bool) – If True, enlarge the truncation range specified by other arguments by shifting the cut to the largest gap between to-be-truncated singular values across all blocks. It provides a heuristic mechanism to avoid truncating part of a multiplet. If True, tol_block and D_block are ignored, as largest_gap is a global condition. The default is False.

  • eps_multiplet (float) – Relative tolerance on multiplet splitting. If relative difference between two consecutive elements of S is larger than eps_multiplet, these elements are not considered as part of the same multiplet. Partially truncated multiplets are truncated down. The default is None, when this scheme is not used. If True, tol_block and D_block are ignored, as eps_multiplet is a global condition. Cannot be used together with largest_gap scheme.

  • hermitian (bool) – If True, blocks related by hermitian conjugation are truncated equally, truncating down to the intersecting part. The default is False.

  • mask_f (None | function[yastn.Tensor] -> yastn.Tensor) – It is possible to provide a custom mask function, which provides a mechanism to pass such a function to many tensor network algorithms where the function truncation_mask is being called. If provided, it overrides the default function, and all other parameters are ignored. The default is None.

yastn.linalg.entropy(a, alpha=1, tol=1e-12) → Number[source]#

Compute the entropy from probabilities encoded in the diagonal tensor a.

The sum of a is normalized to 1, but its correctness is not checked otherwise. Base-2 logarithms are used. For an empty or zero tensor, 0 is returned.

Parameters:

  • alpha (float) – Order of Renyi entropy. alpha=1 (the default) is von Neumann entropy: \(-{\rm Tr}(a \cdot {\rm log2}(a))\) otherwise: \(\frac{1}{1-alpha} {\rm log2}({\rm Tr}(a^{alpha}))\)

  • tol (float) – Discard all probabilities smaller than tol during calculation.

Auxiliary#

Eliminating individual blocks

Tensor.remove_zero_blocks(rtol=1e-12, atol=0) → Tensor#

Remove blocks whose entries are below a cutoff.

The cutoff combines an absolute tolerance and a relative tolerance with respect to the largest element of the tensor.

Tensor.remove_random_blocks(number, keep_legs=True) → Tensor#

Randomly remove a number of blocks from a tensor.

The function attempts to remove number blocks from the tensor structure, selecting them randomly. If keep_legs is True, any attempted removal that would change the tensor legs is rejected. The number of actually removed blocks can be smaller than requested. If no blocks are removed, the original tensor is returned. Otherwise, a new tensor with updated structure and copied data is returned.

The function is useful mostly for testing.

Parameters:

  • number (int) – Number of attempts to remove random block.

  • keep_legs (bool) – If True, ensures that tensor legs are not changed, i.e. all leg charges appear in some blocks. The default is True.

Methods called by Krylov-based algorithms.

Tensor.expand_krylov_space(f, tol, ncv, hermitian, V, H=None, **kwargs)#

Expand the Krylov base up to ncv states or until reaching desired tolerance tol. Implementation for yastn.Tensor.