quimb.core

Core functions for manipulating quantum objects.

Module Attributes

njit(*args, **kws)

Numba no-python jit, but obeying cache setting.

pnjit(*args, **kws)

Numba no-python jit, but obeying cache setting, with optional parallel target, depending on environment variable ‘QUIMB_NUMBA_PARALLEL’.

vectorize([ftylist_or_function, target, …])

Numba vectorize, but obeying cache setting.

pvectorize([ftylist_or_function])

Numba vectorize, but obeying cache setting, with optional parallel target, depending on environment variable ‘QUIMB_NUMBA_PARALLEL’.

expec(a, b)

Alias for expectation().

qu(data[, qtype, normalized, chopped, …])

Alias of quimbify().

ket(data, *[, qtype, normalized, chopped, …])

Convert an object into a ket.

bra(data, *[, qtype, normalized, chopped, …])

Convert an object into a bra.

dop(data, *[, qtype, normalized, chopped, …])

Convert an object into a density operator.

sparse(data[, qtype, normalized, chopped, …])

Convert an object into sparse form.

eye(d[, sparse, stype, dtype])

Alias for identity().

speye(d, *[, sparse, stype, dtype])

Sparse identity.

nmlz(qob[, inplace])

Alias for normalize().

tr(mat)

Alias for trace().

ptr(p, dims, keep)

Alias for partial_trace().

Functions

allclose_sparse(A, B, **allclose_opts)

chop(qob[, tol, inplace])

Set small values of a dense or sparse array to zero.

common_type(*arrays)

Quick compute the minimal dtype sufficient for arrays.

complex_array(real, imag)

Accelerated creation of complex array.

csr_mulvec_wrap(fn)

Dispatch sparse csr-vector multiplication to parallel method.

dag(qob)

Conjugate transpose.

dim_compress(dims, inds)

Compress neighbouring subsytem dimensions.

dim_map(dims, coos[, cyclic, trim])

Flatten 2d+ dimensions and coordinates.

divide_update_(X, c, out)

Accelerated computation of X / c into out.

dot(a, b)

Matrix multiplication, dispatched to dense or sparse functions.

dot_sparse(a, b)

Dot product for sparse matrix, dispatching to parallel for v large nnz.

dynal(x, bases)

Generate ‘dynamic decimal’ for x given dims.

ensure_qarray(fn)

Decorator that wraps output as a qarray.

expec(a, b)

Alias for expectation().

expectation(a, b)

‘Expectation’ between a vector/operator and another vector/operator.

eye(d[, sparse, stype, dtype])

Alias for identity().

gen_matching_dynal(ri, rf, dims)

Return the matching dynal part of ri and rf, plus the first pair that don’t match.

gen_ops_maybe_sliced(ops, ix)

Take ops and slice the first few, according to the length of ix and with ix, and leave the rest.

identity(d[, sparse, stype, dtype])

Return identity of size d in complex format, optionally sparse.

ikron(ops, dims, inds[, sparse, stype, …])

Tensor an operator into a larger space by padding with identities.

ind_complement(inds, n)

Return the indices below n not contained in inds.

infer_size(p[, base])

Infer the size, i.e. number of ‘sites’ in a state.

isbra(qob)

Checks if qob is in bra form – an array row.

isdense(qob)

Checks if qob is explicitly dense.

isherm(qob, **allclose_opts)

Checks if qob is hermitian.

isket(qob)

Checks if qob is in ket form – an array column.

isop(qob)

Checks if qob is an operator.

ispos(qob[, tol])

Checks if the dense hermitian qob is approximately positive semi-definite, using the cholesky decomposition.

isreal(qob, **allclose_opts)

Checks if qob is approximately real.

issparse(qob)

Checks if qob is explicitly sparse.

isvec(qob)

Checks if qob is row-vector, column-vector or one-dimensional.

itrace(a[, axes])

General tensor trace, i.e. multiple contractions, for a dense array.

kron(*ops[, stype, coo_build, parallel, …])

Tensor (kronecker) product of variable number of arguments.

kron_dense(a, b[, par_thresh])

kron_dispatch(a, b[, stype])

Kronecker product of two arrays, dispatched based on dense/sparse and also size of product.

kron_sparse(a, b[, stype])

Sparse tensor (kronecker) product,

kronpow(a, p, **kron_opts)

Returns a tensored with itself p times

l_diag_dot_dense(diag, mat)

Dot product of diagonal matrix (with only diagonal supplied) and dense matrix.

l_diag_dot_sparse(diag, mat)

Dot product of digonal matrix (with only diagonal supplied) and sparse matrix.

ldmul(diag, mat)

Accelerated left diagonal multiplication.

make_immutable(mat)

Make array read only, in-place.

mul(x, y)

Element-wise multiplication, dispatched to correct dense or sparse function.

mul_dense(a, b)

nmlz(qob[, inplace])

Alias for normalize().

normalize(qob[, inplace])

Normalize a quantum object.

outer(a, b)

Outer product between two vectors (no conjugation).

par_dot_csr_matvec(A, x)

Parallel sparse csr-matrix vector dot product.

par_reduce(fn, seq[, nthreads])

Parallel reduce.

partial_trace(p, dims, keep)

Partial trace of a dense or sparse state.

permute(p, dims, perm)

Permute the subsytems of state or opeator.

pkron(op, dims, inds, **ikron_opts)

Advanced, padded tensor product.

prod(xs)

Product (as in multiplication) of an iterable.

ptr(p, dims, keep)

Alias for partial_trace().

qu(data[, qtype, normalized, chopped, …])

Alias of quimbify().

quimbify(data[, qtype, normalized, chopped, …])

Converts data to ‘quantum’ i.e. complex matrices, kets being columns.

r_diag_dot_dense(mat, diag)

Dot product of dense matrix and digonal matrix (with only diagonal supplied).

r_diag_dot_sparse(mat, diag)

Dot product of sparse matrix and digonal matrix (with only diagonal supplied).

rdmul(mat, diag)

Accelerated left diagonal multiplication.

rdot(a, b)

realify(fn[, imag_tol])

Decorator that drops fn’s output imaginary part if very small.

realify_scalar(x[, imag_tol])

sp_mulvec_wrap(fn)

Scipy sparse doesn’t call __array_finalize__ so need to explicitly make sure qarray input -> qarray output.

sparse_matrix(data[, stype, dtype])

Construct a sparse matrix of a particular format.

subtract_update_(X, c, Y)

Accelerated inplace computation of X -= c * Y.

tr(mat)

Alias for trace().

trace(mat)

Trace of a dense or sparse operator.

upcast(fn)

Decorator to make sure the types of two numpy arguments match.

vdot(a, b)

Accelerated ‘Hermitian’ inner product of two arrays.

zeroify(fn[, tol])

Decorator that rounds fn’s output to zero if very small.

Classes

CacheThreadPool(func)

qarray(data[, dtype, order])

Thin subclass of numpy.ndarray with some convenient quantum linear algebra related methods and attributes (.H, &, etc.), and matrix-like preservation of at least 2-dimensions so as to distiguish kets and bras.

class quimb.core.CacheThreadPool(func)[source]
quimb.core.bra(data, *, qtype='bra', normalized=False, chopped=False, sparse=None, stype=None, dtype=<class 'complex'>)

Convert an object into a bra.

quimb.core.chop(qob, tol=1e-15, inplace=True)[source]

Set small values of a dense or sparse array to zero.

Parameters
  • qob (dense or sparse vector or operator) – Quantum object to chop.

  • tol (float, optional) – Fraction of max(abs(qob)) to chop below.

  • inplace (bool, optional) – Whether to act on input array or return copy.

Returns

Chopped quantum object.

Return type

dense or sparse vector or operator

quimb.core.chop_(qob, tol=1e-15, *, inplace=True)

Set small values of a dense or sparse array to zero.

Parameters
  • qob (dense or sparse vector or operator) – Quantum object to chop.

  • tol (float, optional) – Fraction of max(abs(qob)) to chop below.

  • inplace (bool, optional) – Whether to act on input array or return copy.

Returns

Chopped quantum object.

Return type

dense or sparse vector or operator

quimb.core.common_type(*arrays)[source]

Quick compute the minimal dtype sufficient for arrays.

quimb.core.complex_array(real, imag)[source]

Accelerated creation of complex array.

quimb.core.csr_mulvec_wrap(fn)[source]

Dispatch sparse csr-vector multiplication to parallel method.

quimb.core.dag(qob)[source]

Conjugate transpose.

quimb.core.dim_compress(dims, inds)[source]

Compress neighbouring subsytem dimensions.

Take some dimensions and target indices and compress both, i.e. merge adjacent dimensions that are both either in dims or not. For example, if tensoring an operator onto a single site, with many sites the identity, treat these as single large identities.

Parameters
  • dims (tuple of int) – List of system’s dimensions - 1d or flattened (e.g. with dim_map).

  • inds (tuple of int) – List of target indices, i.e. dimensions not to merge.

Returns

  • dims (tuple of int) – New compressed dimensions.

  • inds (tuple of int) – New indexes corresponding to the compressed dimensions. These are guaranteed to now be alternating i.e. either (0, 2, …) or (1, 3, …).

Examples

>>> dims = [2] * 10
>>> inds = [3, 4]
>>> compressed_dims, compressed_inds = dim_compress(dims, inds)
>>> compressed_dims
(8, 4, 32)
>>> compressed_inds
(1,)
quimb.core.dim_map(dims, coos, cyclic=False, trim=False)[source]

Flatten 2d+ dimensions and coordinates.

Maps multi-dimensional coordinates and indices to flat arrays in a regular way. Wraps or deletes coordinates beyond the system size depending on parameters cyclic and trim.

Parameters
  • dims (nested tuple of int) – Multi-dim array of systems’ internal dimensions.

  • coos (list of tuples of int) – Array of coordinate tuples to convert

  • cyclic (bool, optional) – Whether to automatically wrap coordinates beyond system size or delete them.

  • trim (bool, optional) – If True, any coordinates beyond dimensions will be deleted, overidden by cyclic.

Returns

  • flat_dims (tuple) – Flattened version of dims.

  • inds (tuple) – Indices corresponding to the original coordinates.

Examples

>>> dims = [[2, 3], [4, 5]]
>>> coords = [(0, 0), (1, 1)]
>>> flat_dims, inds = dim_map(dims, coords)
>>> flat_dims
(2, 3, 4, 5)
>>> inds
(0, 3)
>>> dim_map(dims, [(2, 0), (-1, 1)], cyclic=True)
((2, 3, 4, 5), (0, 3))
quimb.core.divide_update_(X, c, out)[source]

Accelerated computation of X / c into out.

quimb.core.dop(data, *, qtype='dop', normalized=False, chopped=False, sparse=None, stype=None, dtype=<class 'complex'>)

Convert an object into a density operator.

quimb.core.dot(a, b)[source]

Matrix multiplication, dispatched to dense or sparse functions.

Parameters
  • a (dense or sparse operator) – First array.

  • b (dense or sparse operator) – Second array.

Returns

Dot product of a and b.

Return type

dense or sparse operator

quimb.core.dot_sparse(a, b)[source]

Dot product for sparse matrix, dispatching to parallel for v large nnz.

quimb.core.dynal(x, bases)[source]

Generate ‘dynamic decimal’ for x given dims.

Examples

>>> dims = [13, 2, 7, 3, 10]
>>> prod(dims)  # total hilbert space size
5460
>>> x = 3279
>>> drep = list(dyn_bin(x, dims))  # dyn bases repr
>>> drep
[7, 1, 4, 0, 9]
>>> bs_szs = [prod(dims[i + 1:]) for i in range(len(dims))]
>>> bs_szs
[420, 210, 30, 10, 1]
>>> # reconstruct x
>>> sum(d * b for d, b in zip(drep, bs_szs))
3279
quimb.core.ensure_qarray(fn)[source]

Decorator that wraps output as a qarray.

quimb.core.expec(a, b)

Alias for expectation().

quimb.core.expectation(a, b)[source]

‘Expectation’ between a vector/operator and another vector/operator.

The ‘operator’ inner product between a and b, but also for vectors. This means that for consistency:

  • for two vectors it will be the absolute expec squared |<a|b><b|a>|, not <a|b>.

  • for a vector and an operator its will be <a|b|a>

  • for two operators it will be the Hilbert-schmidt inner product tr(A @ B)

In this way expectation(a, b) == expectation(dop(a), b) == expectation(dop(a), dop(b)).

Parameters
  • a (vector or operator) – First state or operator - assumed to be ket if vector.

  • b (vector or operator) – Second state or operator - assumed to be ket if vector.

Returns

x – ‘Expectation’ of a with b.

Return type

float

quimb.core.eye(d, sparse=False, stype='csr', dtype=<class 'complex'>)

Alias for identity().

quimb.core.gen_matching_dynal(ri, rf, dims)[source]

Return the matching dynal part of ri and rf, plus the first pair that don’t match.

quimb.core.gen_ops_maybe_sliced(ops, ix)[source]

Take ops and slice the first few, according to the length of ix and with ix, and leave the rest.

quimb.core.identity(d, sparse=False, stype='csr', dtype=<class 'complex'>)[source]

Return identity of size d in complex format, optionally sparse.

Parameters
  • d (int) – Dimension of identity.

  • sparse (bool, optional) – Whether to output in sparse form.

  • stype (str, optional) – If sparse, what format to use.

Returns

id – Identity operator.

Return type

qarray or sparse matrix

quimb.core.ikron(ops, dims, inds, sparse=None, stype=None, coo_build=False, parallel=False, ownership=None)[source]

Tensor an operator into a larger space by padding with identities.

Automatically placing a large operator over several dimensions is allowed and a list of operators can be given which are then placed cyclically.

Parameters
  • op (operator or sequence of operators) – Operator(s) to place into the tensor space. If more than one, these are cyclically placed at each of the dims specified by inds.

  • dims (sequence of int or nested sequences of int) – The subsystem dimensions. If treated as an array, should have the same number of dimensions as the system.

  • inds (tuple of int, or sequence of tuple of int) – Indices, or coordinates, of the dimensions to place operator(s) on. Each dimension specified can be smaller than the size of op (as long as it factorizes it).

  • sparse (bool, optional) – Whether to construct the new operator in sparse form.

  • stype (str, optional) – If sparse, which format to use for the output.

  • coo_build (bool, optional) – Whether to build the intermediary matrices using the 'coo' format - can be faster to build sparse in this way, then convert to chosen format, including dense.

  • parallel (bool, optional) – Whether to build the operator in parallel using threads (only good for big (d > 2**16) operators).

  • ownership ((int, int), optional) – If given, only construct the rows in range(*ownership). Such that the final operator is actually X[slice(*ownership), :]. Useful for constructing operators in parallel, e.g. for MPI.

Returns

Operator such that ops act on dims[inds].

Return type

qarray or sparse matrix

See also

kron(), pkron()

Examples

Place an operator between two identities:

>>> IZI = ikron(pauli('z'), [2, 2, 2], 1)
>>> np.allclose(IZI, eye(2) & pauli('z') & eye(2))
True

Overlay a large operator on several sites:

>>> rho_ab = rand_rho(4)
>>> rho_abc = ikron(rho_ab, [5, 2, 2, 7], [1, 2])  # overlay both 2s
>>> rho_abc.shape
(140, 140)

Place an operator at specified sites, regardless of size:

>>> A = rand_herm(5)
>>> ikron(A, [2, -1, 2, -1, 2, -1], [1, 3, 5]).shape
(1000, 1000)
quimb.core.ind_complement(inds, n)[source]

Return the indices below n not contained in inds.

quimb.core.infer_size(p, base=2)[source]

Infer the size, i.e. number of ‘sites’ in a state.

Parameters
  • p (vector or operator) – An array representing a state with a shape attribute.

  • base (int, optional) – Size of the individual states that p is composed of, e.g. this defauts 2 for qubits.

Returns

Number of composite systems.

Return type

int

Examples

>>> infer_size(singlet() & singlet())
4
>>> infersize(rand_rho(5**3), base=5)
3
quimb.core.isbra(qob)[source]

Checks if qob is in bra form – an array row.

quimb.core.isdense(qob)[source]

Checks if qob is explicitly dense.

quimb.core.isherm(qob, **allclose_opts)[source]

Checks if qob is hermitian.

Parameters

qob (dense or sparse operator) – Matrix to check.

Returns

Return type

bool

quimb.core.isket(qob)[source]

Checks if qob is in ket form – an array column.

quimb.core.isop(qob)[source]

Checks if qob is an operator.

quimb.core.ispos(qob, tol=1e-15)[source]

Checks if the dense hermitian qob is approximately positive semi-definite, using the cholesky decomposition.

Parameters

qob (dense operator) – Matrix to check.

Returns

Return type

bool

quimb.core.isreal(qob, **allclose_opts)[source]

Checks if qob is approximately real.

quimb.core.issparse(qob)[source]

Checks if qob is explicitly sparse.

quimb.core.isvec(qob)[source]

Checks if qob is row-vector, column-vector or one-dimensional.

quimb.core.itrace(a, axes=0, 1)[source]

General tensor trace, i.e. multiple contractions, for a dense array.

Parameters
  • a (numpy.ndarray) – Tensor to trace.

  • axes ((2,) int or (2,) array of int) –

    • (2,) int: Perform trace on the two indices listed.

    • (2,) array of int: Trace out first sequence of indices with second sequence indices.

Returns

The tensor remaining after tracing out the specified axes.

Return type

numpy.ndarray

Examples

Trace out a single pair of dimensions:

>>> a = randn(2, 3, 4, 2, 3, 4)
>>> itrace(a, axes=(0, 3)).shape
(3, 4, 3, 4)

Trace out multiple dimensions:

>>> itrace(a, axes=([1, 2], [4, 5])).shape
(2, 2)
quimb.core.ket(data, *, qtype='ket', normalized=False, chopped=False, sparse=None, stype=None, dtype=<class 'complex'>)

Convert an object into a ket.

quimb.core.kron(*ops, stype=None, coo_build=False, parallel=False, ownership=None)[source]

Tensor (kronecker) product of variable number of arguments.

Parameters
  • ops (sequence of vectors or matrices) – Objects to be tensored together.

  • stype (str, optional) – Desired output format if resultant object is sparse. Should be one of {'csr', 'bsr', 'coo', 'csc'}. If None, infer from input matrices.

  • coo_build (bool, optional) – Whether to force sparse construction to use the 'coo' format (only for sparse matrices in the first place.).

  • parallel (bool, optional) – Perform a parallel reduce on the operators, can be quicker.

  • ownership ((int, int), optional) – If given, only construct the rows in range(*ownership). Such that the final operator is actually X[slice(*ownership), :]. Useful for constructing operators in parallel, e.g. for MPI.

Returns

X – Tensor product of ops.

Return type

dense or sparse vector or operator

Notes

  1. The product is performed as (a & (b & (c & ...)))

Examples

Simple example:

>>> a = np.array([[1, 2], [3, 4]])
>>> b = np.array([[1., 1.1], [1.11, 1.111]])
>>> kron(a, b)
qarray([[1.   , 1.1  , 2.   , 2.2  ],
        [1.11 , 1.111, 2.22 , 2.222],
        [3.   , 3.3  , 4.   , 4.4  ],
        [3.33 , 3.333, 4.44 , 4.444]])

Partial construction of rows:

>>> ops = [rand_matrix(2, sparse=True) for _ in range(10)]
>>> kron(*ops, ownership=(256, 512))
<256x1024 sparse matrix of type '<class 'numpy.complex128'>'
        with 13122 stored elements in Compressed Sparse Row format>
quimb.core.kron_dispatch(a, b, stype=None)[source]

Kronecker product of two arrays, dispatched based on dense/sparse and also size of product.

quimb.core.kron_sparse(a, b, stype=None)[source]

Sparse tensor (kronecker) product,

Output format can be specified or will be automatically determined.

quimb.core.kronpow(a, p, **kron_opts)[source]

Returns a tensored with itself p times

Equivalent to reduce(lambda x, y: x & y, [a] * p).

Parameters
  • a (dense or sparse vector or operator) – Object to tensor power.

  • p (int) – Tensor power.

  • kron_opts – Supplied to kron().

Returns

Return type

dense or sparse vector or operator

quimb.core.l_diag_dot_dense(diag, mat)[source]

Dot product of diagonal matrix (with only diagonal supplied) and dense matrix.

quimb.core.l_diag_dot_sparse(diag, mat)[source]

Dot product of digonal matrix (with only diagonal supplied) and sparse matrix.

quimb.core.ldmul(diag, mat)[source]

Accelerated left diagonal multiplication. Equivalent to numpy.diag(diag) @ mat, but faster than numpy.

Parameters
  • diag (vector or 1d-array) – Vector representing the diagonal of a matrix.

  • mat (dense or sparse matrix) – A normal (non-diagonal) matrix.

Returns

Dot product of the matrix whose diagonal is diag and mat.

Return type

dense or sparse matrix

quimb.core.make_immutable(mat)[source]

Make array read only, in-place.

Parameters

mat (sparse or dense array) – Matrix to make immutable.

quimb.core.mul(x, y)[source]

Element-wise multiplication, dispatched to correct dense or sparse function.

Parameters
  • x (dense or sparse operator) – First array.

  • y (dense or sparse operator) – Second array.

Returns

Element wise product of x and y.

Return type

dense or sparse operator

quimb.core.njit(*args, **kws)

Numba no-python jit, but obeying cache setting.

quimb.core.nmlz(qob, inplace=True)

Alias for normalize().

quimb.core.normalize(qob, inplace=True)[source]

Normalize a quantum object.

Parameters
  • qob (dense or sparse vector or operator) – Quantum object to normalize.

  • inplace (bool, optional) – Whether to act inplace on the given operator.

Returns

Normalized quantum object.

Return type

dense or sparse vector or operator

quimb.core.normalize_(qob, *, inplace=True)

Normalize a quantum object.

Parameters
  • qob (dense or sparse vector or operator) – Quantum object to normalize.

  • inplace (bool, optional) – Whether to act inplace on the given operator.

Returns

Normalized quantum object.

Return type

dense or sparse vector or operator

quimb.core.outer(a, b)[source]

Outer product between two vectors (no conjugation).

quimb.core.par_dot_csr_matvec(A, x)[source]

Parallel sparse csr-matrix vector dot product.

Parameters
Returns

Result of A @ x.

Return type

dense vector

Notes

The main bottleneck for sparse matrix vector product is memory access, as such this function is only beneficial for pretty large matrices.

quimb.core.par_reduce(fn, seq, nthreads=1)[source]

Parallel reduce.

Parameters
  • fn (callable) – Two argument function to reduce with.

  • seq (sequence) – Sequence to reduce.

  • nthreads (int, optional) – The number of threads to reduce with in parallel.

Returns

Return type

depends on fn and seq.

Notes

This has a several hundred microsecond overhead.

quimb.core.partial_trace(p, dims, keep)[source]

Partial trace of a dense or sparse state.

Parameters
  • p (ket or density operator) – State to perform partial trace on - can be sparse.

  • dims (sequence of int or nested sequences of int) – The subsystem dimensions. If treated as an array, should have the same number of dimensions as the system.

  • keep (int, sequence of int or sequence of tuple[int]) – Index or indices of subsytem(s) to keep. If a sequence of integer tuples, each should be a coordinate such that the length matches the number of dimensions of the system.

Returns

rho – Density operator of subsytem dimensions dims[keep].

Return type

qarray

See also

itrace()

Examples

Trace out single subsystem of a ket:

>>> psi = bell_state('psi-')
>>> ptr(psi, [2, 2], keep=0)  # expect identity
qarray([[ 0.5+0.j,  0.0+0.j],
        [ 0.0+0.j,  0.5+0.j]])

Trace out multiple subsystems of a density operator:

>>> rho_abc = rand_rho(3 * 4 * 5)
>>> rho_ab = partial_trace(rho_abc, [3, 4, 5], keep=[0, 1])
>>> rho_ab.shape
(12, 12)

Trace out qutrits from a 2D system:

>>> psi_abcd = rand_ket(3 ** 4)
>>> dims = [[3, 3],
...         [3, 3]]
>>> keep = [(0, 0), (1, 1)]
>>> rho_ac = partial_trace(psi_abcd, dims, keep)
>>> rho_ac.shape
(9, 9)
quimb.core.permute(p, dims, perm)[source]

Permute the subsytems of state or opeator.

Parameters
  • p (vector or operator) – State or operator to permute.

  • dims (tuple of int) – Internal dimensions of the system.

  • perm (tuple of int) – New order of indexes range(len(dims)).

Returns

pp – Permuted state or operator.

Return type

vector or operator

See also

pkron()

Examples

>>> IX = speye(2) & pauli('X', sparse=True)
>>> XI = permute(IX, dims=[2, 2], perm=[1, 0])
>>> np.allclose(XI.A, pauli('X') & eye(2))
True
quimb.core.pkron(op, dims, inds, **ikron_opts)[source]

Advanced, padded tensor product.

Construct an operator such that op acts on dims[inds], and allow it to be arbitrarily split and reversed etc., in other words, permute and then tensor it into a larger space.

Parameters
  • ops (matrix-like or tuple of matrix-like) – Operator to place into the tensor space.

  • dims (tuple of int) – Dimensions of tensor space.

  • inds (tuple of int) – Indices of the dimensions to place operators on. If multiple operators are specified, inds[1] corresponds to ops[1] and so on.

  • sparse (bool, optional) – Whether to construct the new operator in sparse form.

  • stype (str, optional) – If sparse, which format to use for the output.

  • coo_build (bool, optional) – Whether to build the intermediary matrices using the 'coo' format - can be faster to build sparse in this way, then convert to chosen format, including dense.

Returns

Operator such that ops act on dims[inds].

Return type

operator

See also

ikron(), permute()

Examples

Here we take an operator that acts on spins 0 and 1 with X and Z, and transform it to act on spins 2 and 0 – i.e. reverse it and sandwich an identity between the two sites it acts on.

>>> XZ = pauli('X') & pauli('Z')
>>> ZIX = pkron(XZ, dims=[2, 3, 2], inds=[2, 0])
>>> np.allclose(ZIX, pauli('Z') & eye(3) & pauli('X'))
True
quimb.core.pnjit(*args, **kws)

Numba no-python jit, but obeying cache setting, with optional parallel target, depending on environment variable ‘QUIMB_NUMBA_PARALLEL’.

quimb.core.prod(xs)[source]

Product (as in multiplication) of an iterable.

quimb.core.ptr(p, dims, keep)

Alias for partial_trace().

quimb.core.pvectorize(ftylist_or_function=(), **kws)

Numba vectorize, but obeying cache setting, with optional parallel target, depending on environment variable ‘QUIMB_NUMBA_PARALLEL’.

class quimb.core.qarray(data, dtype=None, order=None)[source]

Thin subclass of numpy.ndarray with some convenient quantum linear algebra related methods and attributes (.H, &, etc.), and matrix-like preservation of at least 2-dimensions so as to distiguish kets and bras.

quimb.core.qu(data, qtype=None, normalized=False, chopped=False, sparse=None, stype=None, dtype=<class 'complex'>)

Alias of quimbify().

quimb.core.quimbify(data, qtype=None, normalized=False, chopped=False, sparse=None, stype=None, dtype=<class 'complex'>)[source]

Converts data to ‘quantum’ i.e. complex matrices, kets being columns.

Parameters
  • data (dense or sparse array_like) – Array describing vector or operator.

  • qtype ({'ket', 'bra' or 'dop'}, optional) – Quantum object type output type. Note that if an operator is given as data and 'ket' or 'bra' as qtype, the operator will be unravelled into a column or row vector.

  • sparse (bool, optional) – Whether to convert output to sparse a format.

  • normalized (bool, optional) – Whether to normalise the output.

  • chopped (bool, optional) – Whether to trim almost zero entries of the output.

  • stype ({'csr', 'csc', 'bsr', 'coo'}, optional) – Format of output matrix if sparse, defaults to 'csr'.

Returns

Return type

dense or sparse vector or operator

Notes

  1. Will unravel an array if 'ket' or 'bra' given.

  2. Will conjugate if 'bra' given.

  3. Will leave operators as is if 'dop' given, but construct one if vector given with the assumption that it was a ket.

Examples

Create a ket (column vector):

>>> qu([1, 2j, 3])
qarray([[1.+0.j],
        [0.+2.j],
        [3.+0.j]])

Create a single precision bra (row vector):

>>> qu([1, 2j, 3], qtype='bra', dtype='complex64')
qarray([[1.-0.j, 0.-2.j, 3.-0.j]], dtype=complex64)

Create a density operator from a vector:

>>> qu([1, 2j, 3], qtype='dop')
qarray([[1.+0.j, 0.-2.j, 3.+0.j],
        [0.+2.j, 4.+0.j, 0.+6.j],
        [3.+0.j, 0.-6.j, 9.+0.j]])

Create a sparse density operator:

>>> qu([1, 0, 0], sparse=True, qtype='dop')
<3x3 sparse matrix of type '<class 'numpy.complex128'>'
    with 1 stored elements in Compressed Sparse Row format>
quimb.core.r_diag_dot_dense(mat, diag)[source]

Dot product of dense matrix and digonal matrix (with only diagonal supplied).

quimb.core.r_diag_dot_sparse(mat, diag)[source]

Dot product of sparse matrix and digonal matrix (with only diagonal supplied).

quimb.core.rdmul(mat, diag)[source]

Accelerated left diagonal multiplication.

Equivalent to mat @ numpy.diag(diag), but faster.

Parameters
  • mat (dense or sparse matrix) – A normal (non-diagonal) matrix.

  • diag (vector or 1d-array) – Vector representing the diagonal of a matrix.

Returns

Dot product of mat and the matrix whose diagonal is diag.

Return type

dense or sparse matrix

quimb.core.realify(fn, imag_tol=1e-12)[source]

Decorator that drops fn’s output imaginary part if very small.

quimb.core.sp_mulvec_wrap(fn)[source]

Scipy sparse doesn’t call __array_finalize__ so need to explicitly make sure qarray input -> qarray output.

quimb.core.sparse(data, qtype=None, normalized=False, chopped=False, *, sparse=True, stype=None, dtype=<class 'complex'>)

Convert an object into sparse form.

quimb.core.sparse_matrix(data, stype='csr', dtype=<class 'complex'>)[source]

Construct a sparse matrix of a particular format.

Parameters
  • data (array_like) – Fed to scipy.sparse constructor.

  • stype ({'csr', 'csc', 'coo', 'bsr'}, optional) – Sparse format.

Returns

Of format stype.

Return type

scipy sparse matrix

quimb.core.speye(d, *, sparse=True, stype='csr', dtype=<class 'complex'>)

Sparse identity.

quimb.core.subtract_update_(X, c, Y)[source]

Accelerated inplace computation of X -= c * Y. This is mainly for Lanczos iteration.

quimb.core.tr(mat)

Alias for trace().

quimb.core.trace(mat)[source]

Trace of a dense or sparse operator.

Parameters

mat (operator) – Operator, dense or sparse.

Returns

x – Trace of mat

Return type

float

quimb.core.upcast(fn)[source]

Decorator to make sure the types of two numpy arguments match.

quimb.core.vdot(a, b)[source]

Accelerated ‘Hermitian’ inner product of two arrays. In other words, b here will be conjugated by the function.

quimb.core.vectorize(ftylist_or_function=(), target='cpu', identity=None, **kws)

Numba vectorize, but obeying cache setting.

quimb.core.zeroify(fn, tol=1e-14)[source]

Decorator that rounds fn’s output to zero if very small.