Decompositions#

vec_to_mps#

tensorkrowch.decompositions.vec_to_mps(vec, n_batches=0, rank=None, cum_percentage=None, cutoff=None, renormalize=False)[source]#

Splits a vector into a sequence of MPS tensors via consecutive SVD decompositions. The resultant tensors can be used to instantiate a MPS with boundary = "obc".

The number of resultant tensors and their respective physical dimensions depend on the shape of the input vector. That is, if one expects to recover a MPS with physical dimensions

\[d_1 \times \cdots \times d_n\]

the input vector will have to be provided with that shape. This can be done with reshape.

If the input vector has batch dimensions, having as shape

\[b_1 \times \cdots \times b_m \times d_1 \times \cdots \times d_n\]

the number of batch dimensions \(m\) can be specified in n_batches. In this case, the resultant tensors will all have the extra batch dimensions. These tensors can be used to instantiate a MPSData with boundary = "obc".

To specify the bond dimension of each cut done via SVD, one can use the arguments rank, cum_percentage and cutoff. If more than one is specified, the resulting rank will be the one that satisfies all conditions.

Parameters
  • vec (torch.Tensor) – Input vector to decompose.

  • n_batches (int) – Number of batch dimensions of the input vector. Each resultant tensor will have also the corresponding batch dimensions. It should be between 0 and the rank of vec.

  • rank (int, optional) – Number of singular values to keep.

  • cum_percentage (float, optional) –

    Proportion that should be satisfied between the sum of all singular values kept and the total sum of all singular values.

    \[\frac{\sum_{i \in \{kept\}}{s_i}}{\sum_{i \in \{all\}}{s_i}} \ge cum\_percentage\]

  • cutoff (float, optional) – Quantity that lower bounds singular values in order to be kept.

  • renormalize (bool) – Indicates whether nodes should be renormalized after SVD/QR decompositions. If not, it may happen that the norm explodes as it is being accumulated from all nodes. Renormalization aims to avoid this undesired behavior by extracting the norm of each node on a logarithmic scale after SVD/QR decompositions are computed. Finally, the normalization factor is evenly distributed among all nodes of the MPS.

Return type

List[torch.Tensor]

mat_to_mpo#

tensorkrowch.decompositions.mat_to_mpo(mat, rank=None, cum_percentage=None, cutoff=None, renormalize=False)[source]#

Splits a matrix into a sequence of MPO tensors via consecutive SVD decompositions. The resultant tensors can be used to instantiate a MPO with boundary = "obc".

The number of resultant tensors and their respective input/output dimensions depend on the shape of the input matrix. That is, if one expects to recover a MPO with input/output dimensions

\[in_1 \times out_1 \times \cdots \times in_n \times out_n\]

the input matrix will have to be provided with that shape. Thus it must have an even number of dimensions. To accomplish this, it may happen that some input/output dimensions are 1. This can be done with reshape.

To specify the bond dimension of each cut done via SVD, one can use the arguments rank, cum_percentage and cutoff. If more than one is specified, the resulting rank will be the one that satisfies all conditions.

Parameters
  • mat (torch.Tensor) – Input matrix to decompose. It must have an even number of dimensions.

  • rank (int, optional) – Number of singular values to keep.

  • cum_percentage (float, optional) –

    Proportion that should be satisfied between the sum of all singular values kept and the total sum of all singular values.

    \[\frac{\sum_{i \in \{kept\}}{s_i}}{\sum_{i \in \{all\}}{s_i}} \ge cum\_percentage\]

  • cutoff (float, optional) – Quantity that lower bounds singular values in order to be kept.

  • renormalize (bool) – Indicates whether nodes should be renormalized after SVD/QR decompositions. If not, it may happen that the norm explodes as it is being accumulated from all nodes. Renormalization aims to avoid this undesired behavior by extracting the norm of each node on a logarithmic scale after SVD/QR decompositions are computed. Finally, the normalization factor is evenly distributed among all nodes of the MPS.

Return type

List[torch.Tensor]