quimb.gen.rand

Functions for generating random quantum objects and states.

Attributes

_NUM_THREAD_WORKERS

nmlz

Alias for normalize().

ptr

Alias for partial_trace().

pvectorize

Numba vectorize, but obeying cache setting, with optional parallel

qu

Alias of quimbify().

vectorize

Numba vectorize, but obeying cache setting.

_RG_HANDLER

choice

_phase_sigs

_phase_to_complex_seq

Turn array of phases into unit circle complex numbers - sequential.

_phase_to_complex_par

Turn array of phases into unit circle complex numbers - parallel.

rand_mps

Classes

qarray

Thin subclass of numpy.ndarray with some convenient quantum

_RGenHandler

Private object that handles pool of random number generators for

Functions

complex_array(real, imag)

Accelerated creation of complex array.

dag(qob)

Conjugate transpose.

dot(a, b)

Matrix multiplication, dispatched to dense or sparse functions.

get_thread_pool([num_workers])

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

Tensor (kronecker) product of variable number of arguments.

prod(iterable)

rdmul(mat, diag)

Accelerated left diagonal multiplication.

seed_rand(seed)

See the random number generators, by instantiating a new set of bit

set_rand_bitgen(bitgen)

Set the core bit generator type to use, from either numpy or

_get_rgens(num_threads)

randn([shape, dtype, scale, loc, num_threads, seed, dist])

Fast multithreaded generation of random normally distributed data.

rand(*args, **kwargs)

random_seed_fn(fn)

Modify fn to take a seed argument (so as to seed the random

_randint(*args, **kwargs)

_choice(*args, **kwargs)

rand_rademacher(shape[, scale, loc, dtype])

_phase_to_complex_base(x)

phase_to_complex(x)

rand_phase(shape[, scale, dtype])

Generate random complex numbers distributed on the unit sphere.

rand_matrix(d[, scaled, sparse, stype, density, ...])

Generate a random matrix of order d with normally distributed

rand_herm(d[, sparse, density, dtype])

Generate a random hermitian operator of order d with normally

rand_pos(d[, sparse, density, dtype])

Generate a random positive operator of size d, with normally

rand_rho(d[, sparse, density, dtype])

Generate a random positive operator of size d with normally

rand_uni(d[, dtype])

Generate a random unitary operator of size d, distributed according to

rand_ket(d[, sparse, stype, density, dtype])

Generates a ket of length d with normally distributed entries.

rand_haar_state(d[, dtype])

Generate a random state of dimension d according to the Haar

gen_rand_haar_states(d, reps[, dtype])

Generate many random Haar states, recycling a random unitary operator

rand_mix(d[, tr_d_min, tr_d_max, mode, dtype])

Constructs a random mixed state by tracing out a random ket

rand_product_state(n[, qtype, dtype])

Generates a ket of n many random pure qubits.

rand_matrix_product_state(n, bond_dim[, phys_dim, ...])

Generate a random matrix product state (in dense form, see

rand_seperable(dims[, num_mix, dtype])

Generate a random, mixed, seperable state. E.g rand_seperable([2, 2])

rand_iso(n, m[, dtype])

Generate a random isometry of shape (n, m).

rand_mera(n[, invariant, dtype])

Generate a random mera state of n qubits, which must be a power

Module Contents

quimb.gen.rand._NUM_THREAD_WORKERS
quimb.gen.rand.complex_array(real, imag)[source]

Accelerated creation of complex array.

quimb.gen.rand.dag(qob)[source]

Conjugate transpose.

quimb.gen.rand.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.gen.rand.get_thread_pool(num_workers=None)
quimb.gen.rand.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.gen.rand.nmlz[source]

Alias for normalize().

quimb.gen.rand.prod(iterable)
quimb.gen.rand.ptr[source]

Alias for partial_trace().

quimb.gen.rand.pvectorize

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

class quimb.gen.rand.qarray(shape, dtype=float, buffer=None, offset=0, strides=None, order=None)[source]

Bases: numpy.ndarray

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.

property H
property A
__array__()[source]
__and__(other)[source]
normalize(inplace=True)[source]
nmlz(inplace=True)[source]
chop(inplace=True)[source]
tr()[source]
partial_trace(dims, keep)[source]
ptr(dims, keep)[source]
__str__()[source]

Return str(self).

__repr__()[source]

Return repr(self).

quimb.gen.rand.qu[source]

Alias of quimbify().

quimb.gen.rand.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.gen.rand.vectorize

Numba vectorize, but obeying cache setting.

class quimb.gen.rand._RGenHandler(initial_seed=None, initial_bitgen=None)[source]

Private object that handles pool of random number generators for parallel number generation - seeding them & changing the underlying bit generators.

set_bitgen(bitgen)[source]

Set the core underlying bit-generator.

Parameters:

bitgen ({None, str}) – Which bit generator to use, either from numpy or randomgen - https://bashtage.github.io/randomgen/bit_generators/index.html.

set_seed(seed=None)[source]

Set the seed for the bit generators.

Parameters:

seed ({None, int}, optional) – Seed supplied to numpy.random.SeedSequence. None will randomly the generators (default).

get_rgens(num_threads)[source]

Get a list of the numpy.random.Generator instances, having made sure there are at least num_threads.

Parameters:

num_threads (int) – The number of generators to return.

Return type:

list[numpy.random.Generator]

quimb.gen.rand._RG_HANDLER
quimb.gen.rand.seed_rand(seed)[source]

See the random number generators, by instantiating a new set of bit generators with a ‘seed sequence’.

quimb.gen.rand.set_rand_bitgen(bitgen)[source]

Set the core bit generator type to use, from either numpy or randomgen.

Parameters:

bitgen ({'PCG64', 'SFC64', 'MT19937', 'Philox', str}) – Which bit generator to use.

quimb.gen.rand._get_rgens(num_threads)[source]
quimb.gen.rand.randn(shape=(), dtype=float, scale=1.0, loc=0.0, num_threads=None, seed=None, dist='normal')[source]

Fast multithreaded generation of random normally distributed data.

Parameters:
  • shape (tuple[int]) – The shape of the output random array.

  • dtype ({'complex128', 'float64', 'complex64' 'float32'}, optional) – The data-type of the output array.

  • scale (float, optional) – A multiplicative scale for the random numbers.

  • loc (float, optional) – An additive location for the random numbers.

  • num_threads (int, optional) – How many threads to use. If None, decide automatically.

  • dist ({'normal', 'uniform', 'rademacher', 'exp'}, optional) – Type of random number to generate.

quimb.gen.rand.rand(*args, **kwargs)[source]
quimb.gen.rand.random_seed_fn(fn)[source]

Modify fn to take a seed argument (so as to seed the random generators once-only at beginning of function not every randn call).

quimb.gen.rand._randint(*args, **kwargs)[source]
quimb.gen.rand._choice(*args, **kwargs)[source]
quimb.gen.rand.choice[source]
quimb.gen.rand.rand_rademacher(shape, scale=1, loc=0.0, dtype=float)[source]
quimb.gen.rand._phase_to_complex_base(x)[source]
quimb.gen.rand._phase_sigs = ['complex64(float32)', 'complex128(float64)']
quimb.gen.rand._phase_to_complex_seq[source]

Turn array of phases into unit circle complex numbers - sequential.

quimb.gen.rand._phase_to_complex_par

Turn array of phases into unit circle complex numbers - parallel.

quimb.gen.rand.phase_to_complex(x)[source]
quimb.gen.rand.rand_phase(shape, scale=1, dtype=complex)[source]

Generate random complex numbers distributed on the unit sphere.

quimb.gen.rand.rand_matrix(d, scaled=True, sparse=False, stype='csr', density=None, dtype=complex, seed=None)[source]

Generate a random matrix of order d with normally distributed entries. If scaled is True, then in the limit of large d the eigenvalues will be distributed on the unit complex disk.

Parameters:
  • d (int) – Matrix dimension.

  • scaled (bool, optional) – Whether to scale the matrices values such that its spectrum approximately lies on the unit disk (for dense matrices).

  • sparse (bool, optional) – Whether to produce a sparse matrix.

  • stype ({'csr', 'csc', 'coo', ...}, optional) – The type of sparse matrix if sparse=True.

  • density (float, optional) – Target density of non-zero elements for the sparse matrix. By default aims for about 10 entries per row.

  • dtype ({complex, float}, optional) – The data type of the matrix elements.

Returns:

mat – Random matrix.

Return type:

qarray or sparse matrix

quimb.gen.rand.rand_herm(d, sparse=False, density=None, dtype=complex)[source]

Generate a random hermitian operator of order d with normally distributed entries. In the limit of large d the spectrum will be a semi-circular distribution between [-1, 1].

quimb.gen.rand.rand_pos(d, sparse=False, density=None, dtype=complex)[source]

Generate a random positive operator of size d, with normally distributed entries. In the limit of large d the spectrum will lie between [0, 1].

quimb.gen.rand.rand_rho(d, sparse=False, density=None, dtype=complex)[source]

Generate a random positive operator of size d with normally distributed entries and unit trace.

quimb.gen.rand.rand_uni(d, dtype=complex)[source]

Generate a random unitary operator of size d, distributed according to the Haar measure.

quimb.gen.rand.rand_ket(d, sparse=False, stype='csr', density=0.01, dtype=complex)[source]

Generates a ket of length d with normally distributed entries.

quimb.gen.rand.rand_haar_state(d, dtype=complex)[source]

Generate a random state of dimension d according to the Haar distribution.

quimb.gen.rand.gen_rand_haar_states(d, reps, dtype=complex)[source]

Generate many random Haar states, recycling a random unitary operator by using all of its columns (not a good idea?).

quimb.gen.rand.rand_mix(d, tr_d_min=None, tr_d_max=None, mode='rand', dtype=complex)[source]

Constructs a random mixed state by tracing out a random ket where the composite system varies in size between 2 and d. This produces a spread of states including more purity but has no real meaning.

quimb.gen.rand.rand_product_state(n, qtype=None, dtype=complex)[source]

Generates a ket of n many random pure qubits.

quimb.gen.rand.rand_matrix_product_state(n, bond_dim, phys_dim=2, dtype=complex, cyclic=False, trans_invar=False)[source]

Generate a random matrix product state (in dense form, see MPS_rand_state() for tensor network form).

Parameters:
  • n (int) – Number of sites.

  • bond_dim (int) – Dimension of the bond (virtual) indices.

  • phys_dim (int, optional) – Physical dimension of each local site, defaults to 2 (qubits).

  • cyclic (bool (optional)) – Whether to impose cyclic boundary conditions on the entanglement structure.

  • trans_invar (bool (optional)) – Whether to generate a translationally invariant state, requires cyclic=True.

Returns:

ket – The random state, with shape (phys_dim**n, 1)

Return type:

qarray

quimb.gen.rand.rand_mps[source]
quimb.gen.rand.rand_seperable(dims, num_mix=10, dtype=complex)[source]

Generate a random, mixed, seperable state. E.g rand_seperable([2, 2]) for a mixed two qubit state with no entanglement.

Parameters:
  • dims (tuple of int) – The local dimensions across which to be seperable.

  • num_mix (int, optional) – How many individual product states to sum together, each with random weight.

Returns:

Mixed seperable state.

Return type:

qarray

quimb.gen.rand.rand_iso(n, m, dtype=complex)[source]

Generate a random isometry of shape (n, m).

quimb.gen.rand.rand_mera(n, invariant=False, dtype=complex)[source]

Generate a random mera state of n qubits, which must be a power of 2. This uses quimb.tensor.