cuquantum.cutensornet.experimental.contract_decompose¶
- cuquantum.cutensornet.experimental.contract_decompose(subscripts, *operands, algorithm=None, options=None, optimize=None, stream=None, return_info=False)[source]¶
Evaluate the compound expression for contraction and decomposition on the input operands.
The expression follows a combination of Einstein summation notation for contraction and the decomposition notation for decomposition (as in
decompose()
). The input represents a tensor network that will be contracted to form an intermediate tensor for subsequent decomposition operation, yielding two or three outputs depending on the final decomposition method. The expression requires explicit specification of modes for the input and output tensors (excludingS
for SVD method). The modes for the intermediate tensor are inferred based on the subscripts representing the output modes by using the implicit form of the Einstein summation expression (similar to the treatment innumpy.einsum
implicit mode).See the notes and examples for clarification.
- Parameters
subscripts – The mode labels (subscripts) defining the contraction and decomposition operation as a comma-separated sequence of characters. Unicode characters are allowed in the expression thereby expanding the size of the tensor network that can be specified.
algorithm – Specify the algorithm to perform the contraction and decomposition. Alternatively, a
dict
containing the parameters for theContractDecomposeAlgorithm
constructor can be provided. If not specified, the value will be set to the default-constructedContractDecomposeAlgorithm
object.operands – A sequence of tensors (ndarray-like objects). The currently supported types are
numpy.ndarray
,cupy.ndarray
, andtorch.Tensor
.options – Specify options for the tensor network as a
NetworkOptions
object. Alternatively, adict
containing the parameters for theNetworkOptions
constructor can also be provided. If not specified, the value will be set to the default-constructedNetworkOptions
object.optimize – This parameter specifies options for path optimization as an
OptimizerOptions
object. Alternatively, a dictionary containing the parameters for theOptimizerOptions
constructor can also be provided. If not specified, the value will be set to the default-constructedOptimizerOptions
object.stream – Provide the CUDA stream to use for the autotuning operation. Acceptable inputs include
cudaStream_t
(as Pythonint
),cupy.cuda.Stream
, andtorch.cuda.Stream
. If a stream is not provided, the current stream will be used.return_info – If true, information about the contraction and decomposition will also be returned as a :
ContractDecomposeInfo
object.
- Returns
For QR decomposition (default), if
return_info
isFalse
, the output tensors Q and R (ndarray-like objects) of the same type and on the same device as the input operand are returned as the result of the decomposition. Ifreturn_info
isTrue
, a 3-tuple of output tensors Q, R and aContractDecomposeInfo
object that contains information about the operations will be returned.For SVD decomposition, if
return_info
isFalse
, a 3-tuple of output tensors U, S and V (ndarray-like objects) of the same type as the input operand are returned as the result of the decomposition. Ifreturn_info
isTrue
, a 4-tuple of output tensors U, S, V and aContractDecomposeInfo
object that contains information about the operations will be returned. Note, depending on the choice ofpartition
, the returned S operand may beNone
. Also seepartition
.
- Return type
Depending on the decomposition setting specified in
algorithm
, the results returned may vary
The contract and decompose expression adopts a combination of Einstein summation notation for contraction and the decomposition notation introduced in
decompose()
. Thesubscripts
string is a list of subscript labels where each label refers to a mode of the corresponding operand. The subscript labels are separated by either comma or identifier->
. The subscript labels before the identifier->
are viewed as inputs, and the ones after are viewed as outputs, respectively. The requirements on the subscripts for SVD and QR decomposition are summarized below:For SVD and QR decomposition, the subscripts string is expected to contain more than one input and exactly two output labels (The modes for
S
is not needed in the case of SVD).One and only one identical mode is expected to exist in the two output mode labels.
The modes for the intermediate tensor which will be decomposed are inferred based on the subscripts representing the output modes by using the implicit form of the Einstein summation expression (similar to the treatment in
numpy.einsum
implicit mode). Therefore, assembling the input modes and the intermediate modes together should result in a validnumpy.einsum
expression (classical or generalized).
Examples
>>> # equivalent: >>> # q, r = numpy.linalg.qr(numpy.einsum('ij,jk->ik', a, b)) >>> q, r = contract_decompose('ij,jk->ix,xk', a, b)
>>> # equivalent: >>> # u, s, v = numpy.linalg.svd(numpy.einsum('ij,jk->ik', a, b), full_matrices=False) >>> u, s, v = tensor_decompose('ij,jk->ix,xk', a, algorithm={'qr_method':False, 'svd_method': {}})
For generalization to generic tensor network with multi-dimensional tensors (
a
,b
,c
are all rank-4 tensors). In this case, the intermediate modesijabe
is inferred from the output modesixeb
andjax
:>>> # equivalent: >>> # t = contract('ijc,cad,dbe->ijabe', a, b, c) >>> # u, s, v = tensor.decompose('ijabe->ixeb,jax', t, method=SVDMethod()) >>> u, s, v = contract_decompose('ijc,cad,dbe->ixeb,jax', a, b, c, algorithm={'qr_method': False, 'svd_method': True})
If the contract and decompose problem amounts to a ternary-operand gate split problem commonly seen in quantum circuit simulation (see Gate Split Algorithm for details), the user may be able to take advantage of optimized kernels from
cutensornetGateSplit()
by placing the gate operand as the last one in the input operands. In this case, QR decomposition can potentially be used to speed up the execution of contraction and SVD. This can be achieved by setting bothqr_method
and bothsvd_method
, as demonstrated below.Example
Applying a two-qubit gate to adjacent MPS tensors:
>>> a, _, b = contract_decompose('ipj,jqk,pqPQ->iPx,xQk', a, b, gate, algorithm={'qr_method':{}, 'svd_method':{}})
Broadcasting is supported for certain cases via ellipsis notation. One may add ellipses in the input modes to represent all the modes that are not explicitly specified in the labels. In such case, an ellipsis is allowed to appear in at most one of the output modes. If an ellipsis appears in one of the output modes, the implicit modes are partitioned onto the corresponding output. If no ellipsis is found in the output, the implicit modes will be summed over to construct the intermediate tensors.
Examples
Below are some examples based on two rank-4 tensors
a
andb
:>>> # equivalent: >>> # out = contract_decompose('ijab,abcd->ijx,xcd', a, b) >>> out = contract_decompose('ijab,ab...->ijx,x...', a, b) # intermediate modes being "ijcd"
>>> # equivalent: >>> # out = contract_decompose('ijab,abcd->ix,xj', a, b) >>> out = contract_decompose('ijab,ab...->ix,xj', a, b) # intermediate modes being "ij"
>>> # equivalent: >>> # out = contract_decompose('ijab,jkab->ix,xj', a, b) >>> out = contract_decompose('ij...,jk...->ix,xj', a, b) # intermediate modes being "ij"
>>> # equivalent: >>> # out = contract_decompose('ijab,jkab->ixab,xj', a, b) >>> out = contract_decompose('ij...,jk...->ix...,xj', a, b) # intermediate modes being "ijab"
Note that the number of modes that are implicitly represented by the ellipses must be the same for all occurrences.
Note
It is encouraged for users to maintain the library handle themselves so as to reduce the context initialization time:
from cuquantum import cutensornet as cutn from cuquantum.cutensornet.experimental import contract_decompose handle = cutn.create() q, r = contract_decompose(..., options={"handle": handle}, ...) # ... the same handle can be reused for further calls ... # when it's done, remember to destroy the handle cutn.destroy(handle)