quimb.gen.states

Functions for generating quantum states.

Attributes

qu

Alias of quimbify().

eye

Alias for identity().

eigh

zplus

zminus

xplus

xminus

Functions

dag(qob)

Conjugate transpose.

ldmul(diag, mat)

Accelerated left diagonal multiplication. Equivalent to

dot(a, b)

Matrix multiplication, dispatched to dense or sparse functions.

make_immutable(mat)

Make array read only, in-place.

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

Tensor (kronecker) product of variable number of arguments.

ikron(ops, dims, inds[, sparse, stype, coo_build, ...])

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

kronpow(a, p, **kron_opts)

Returns a tensored with itself p times

pauli(xyz[, dim])

Generates the pauli operators for dimension 2 or 3.

controlled(s[, dtype, sparse])

Construct a controlled pauli gate for two qubits.

basis_vec(i, dim[, ownership])

Constructs a unit vector ket:

up(**kwargs)

Returns up-state, aka. |0>, +Z eigenstate:

down(**kwargs)

Returns down-state, aka. |1>, -Z eigenstate:

plus(**kwargs)

Returns plus-state, aka. |+>, +X eigenstate:

minus(**kwargs)

Returns minus-state, aka. |->, -X eigenstate:

yplus(**kwargs)

Returns yplus-state, aka. |y+>, +Y eigenstate:

yminus(**kwargs)

Returns yplus-state, aka. |y->, -Y eigenstate:

bloch_state(ax, ay, az[, purified])

Construct qubit density operator from bloch vector:

bell_state(s, **kwargs)

One of the four bell-states.

singlet(**kwargs)

Alias for the 'psi-' bell-state.

thermal_state(ham, beta[, precomp_func])

Generate a thermal state of a Hamiltonian.

computational_state(binary, **kwargs)

Generate the qubit computational state with binary.

neel_state(n[, down_first])

Construct Neel state for n spins, i.e. alternating up/down.

singlet_pairs(n, **kwargs)

Construct fully dimerised spin chain.

werner_state(p, **kwargs)

Construct Werner State, i.e. fractional mix of identity with singlet.

ghz_state(n, **kwargs)

Construct GHZ state of n spins, i.e. equal superposition of all up

w_state(n, **kwargs)

Construct W-state: equal superposition of all single spin up states.

levi_civita(perm)

Compute the generalised levi-civita coefficient for a permutation.

perm_state(ps)

Construct the anti-symmetric state which is the +- sum of all

graph_state_1d(n[, cyclic, sparse])

Graph State on a line.

Module Contents

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

Conjugate transpose.

quimb.gen.states.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.gen.states.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.states.make_immutable(mat)[source]

Make array read only, in-place.

Parameters:

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

quimb.gen.states.qu[source]

Alias of quimbify().

quimb.gen.states.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.states.eye[source]

Alias for identity().

quimb.gen.states.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)

Create a two site interaction (note the coefficient jx we only need to multiply into a single input operator):

>>> Sx = spin_operator('X')
>>> jx = 0.123
>>> jSxSx = ikron([jx * Sx, Sx], [2, 2, 2, 2], [0, 3])
>>> np.allclose(jSxSx, jx * (Sx & eye(2) & eye(2) & Sx))
True
quimb.gen.states.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().

Return type:

dense or sparse vector or operator

quimb.gen.states.eigh
quimb.gen.states.pauli(xyz, dim=2, **kwargs)[source]

Generates the pauli operators for dimension 2 or 3.

Parameters:
  • xyz (str) – Which spatial direction, upper or lower case from {'I', 'X', 'Y', 'Z'}.

  • dim (int, optional) – Dimension of spin operator (e.g. 3 for spin-1), defaults to 2 for spin half.

  • kwargs – Passed to quimbify.

Returns:

P – The pauli operator.

Return type:

immutable operator

See also

spin_operator

quimb.gen.states.controlled(s, dtype=complex, sparse=False)[source]

Construct a controlled pauli gate for two qubits.

Parameters:
  • s (str) – Which pauli to use, including ‘not’ aliased to ‘x’.

  • sparse (bool, optional) – Whether to construct a sparse operator.

Returns:

C – The controlled two-qubit gate operator.

Return type:

qarray

quimb.gen.states.basis_vec(i, dim, ownership=None, **kwargs)[source]

Constructs a unit vector ket:

\[\begin{split}|i\rangle = \begin{pmatrix} 0 \\ \vdots \\ 1 \\ \vdots \\ 0 \end{pmatrix}\end{split}\]
Parameters:
  • i (int) – Which index should the single non-zero, unit entry.

  • dim (int) – Total size of hilbert space.

  • sparse (bool, optional) – Return vector as sparse matrix.

  • kwargs – Supplied to qu.

Returns:

The basis vector.

Return type:

vector

quimb.gen.states.up(**kwargs)[source]

Returns up-state, aka. |0>, +Z eigenstate:

\[\begin{split}|0\rangle = \begin{pmatrix} 1 \\ 0 \end{pmatrix}\end{split}\]
quimb.gen.states.zplus[source]
quimb.gen.states.down(**kwargs)[source]

Returns down-state, aka. |1>, -Z eigenstate:

\[\begin{split}|1\rangle = \begin{pmatrix} 0 \\ 1 \end{pmatrix}\end{split}\]
quimb.gen.states.zminus[source]
quimb.gen.states.plus(**kwargs)[source]

Returns plus-state, aka. |+>, +X eigenstate:

\[\begin{split}|+\rangle = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 \\ 1 \end{pmatrix}\end{split}\]
quimb.gen.states.xplus[source]
quimb.gen.states.minus(**kwargs)[source]

Returns minus-state, aka. |->, -X eigenstate:

\[\begin{split}|-\rangle = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 \\ -1 \end{pmatrix}\end{split}\]
quimb.gen.states.xminus[source]
quimb.gen.states.yplus(**kwargs)[source]

Returns yplus-state, aka. |y+>, +Y eigenstate:

\[\begin{split}|y+\rangle = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 \\ i \end{pmatrix}\end{split}\]
quimb.gen.states.yminus(**kwargs)[source]

Returns yplus-state, aka. |y->, -Y eigenstate:

\[\begin{split}|y-\rangle = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 \\ -i \end{pmatrix}\end{split}\]
quimb.gen.states.bloch_state(ax, ay, az, purified=False, **kwargs)[source]

Construct qubit density operator from bloch vector:

\[\rho = \frac{1}{2} \left( I + a_x X + a_y Y + a_z Z \right)\]
Parameters:
  • ax (float) – X component of bloch vector.

  • ay (float) – Y component of bloch vector.

  • az (float) – Z component of bloch vector.

  • purified – Whether to map vector to surface of bloch sphere.

Returns:

Density operator of qubit ‘pointing’ in (ax, ay, az) direction.

Return type:

Matrix

quimb.gen.states.bell_state(s, **kwargs)[source]

One of the four bell-states.

If n = 2**-0.5, they are:

  1. 'psi-' : n * ( |01> - |10> )

  2. 'psi+' : n * ( |01> + |10> )

  3. 'phi-' : n * ( |00> - |11> )

  4. 'phi+' : n * ( |00> + |11> )

They can be enumerated in this order.

Parameters:
  • s (str or int) – String of number of state corresponding to above.

  • kwargs – Supplied to qu called on state.

Returns:

p – The bell-state s.

Return type:

immutable vector

quimb.gen.states.singlet(**kwargs)[source]

Alias for the ‘psi-’ bell-state.

quimb.gen.states.thermal_state(ham, beta, precomp_func=False)[source]

Generate a thermal state of a Hamiltonian.

Parameters:
  • ham (operator or (1d-array, 2d-array)) – Hamiltonian, either full or tuple of (evals, evecs).

  • beta (float) – Inverse temperature of state.

  • precomp_func (bool, optional) – If True, return a function that takes beta only and is closed over the solved hamiltonian.

Returns:

Density operator of thermal state, or function to generate such given a temperature.

Return type:

operator or callable

quimb.gen.states.computational_state(binary, **kwargs)[source]

Generate the qubit computational state with binary.

Parameters:

binary (sequence of 0s and 1s) – The binary of the computation state.

Examples

>>> computational_state('101'):
qarray([[0.+0.j],
        [0.+0.j],
        [0.+0.j],
        [0.+0.j],
        [0.+0.j],
        [1.+0.j],
        [0.+0.j],
        [0.+0.j]])
>>> qu.computational_state([0, 1], qtype='dop')
qarray([[0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j],
        [0.+0.j, 1.+0.j, 0.+0.j, 0.+0.j],
        [0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j],
        [0.+0.j, 0.+0.j, 0.+0.j, 0.+0.j]])

See also

MPS_computational_state, basic_vec

quimb.gen.states.neel_state(n, down_first=False, **kwargs)[source]

Construct Neel state for n spins, i.e. alternating up/down.

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

  • down_first (bool, optional) – Whether to start with ‘1’ or ‘0’ first.

  • kwargs – Supplied to qu called on state.

See also

computational_state, MPS_neel_state

quimb.gen.states.singlet_pairs(n, **kwargs)[source]

Construct fully dimerised spin chain.

I.e. bell_state('psi-') & bell_state('psi-') & ...

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

  • kwargs – Supplied to qu called on state.

Return type:

vector

quimb.gen.states.werner_state(p, **kwargs)[source]

Construct Werner State, i.e. fractional mix of identity with singlet.

Parameters:
  • p (float) – Singlet Fraction.

  • kwargs – Supplied to qu() called on state.

Return type:

qarray

quimb.gen.states.ghz_state(n, **kwargs)[source]

Construct GHZ state of n spins, i.e. equal superposition of all up and down.

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

  • kwargs – Supplied to qu called on state.

Return type:

vector

quimb.gen.states.w_state(n, **kwargs)[source]

Construct W-state: equal superposition of all single spin up states.

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

  • kwargs – Supplied to qu called on state.

Return type:

vector

quimb.gen.states.levi_civita(perm)[source]

Compute the generalised levi-civita coefficient for a permutation.

Parameters:

perm (sequence of int) – The permutation, a re-arrangement of range(n).

Returns:

Either -1, 0 or 1.

Return type:

int

quimb.gen.states.perm_state(ps)[source]

Construct the anti-symmetric state which is the +- sum of all tensored permutations of states ps.

Parameters:

ps (sequence of states) – The states to combine.

Returns:

The permutation state, dimension same as kron(*ps).

Return type:

vector or operator

Examples

A singlet is the perm_state of up and down.

>>> states = [up(), down()]
>>> pstate = perm_state(states)
>>> expec(pstate, singlet())
1.0
quimb.gen.states.graph_state_1d(n, cyclic=False, sparse=False)[source]

Graph State on a line.

Parameters:
  • n (int) – The number of spins.

  • cyclic (bool, optional) – Whether to use cyclic boundary conditions for the graph.

  • sparse (bool, optional) – Whether to return a sparse state.

Returns:

The 1d-graph state.

Return type:

vector