quimb.gen.operators

Functions for generating quantum operators.

Attributes

qu

Alias of quimbify().

eye

Alias for identity().

Rx

Ry

Rz

cswap

fredkin

toffoli

Classes

qarray

Thin subclass of numpy.ndarray with some convenient quantum

Functions

make_immutable(mat)

Make array read only, in-place.

get_thread_pool([num_workers])

par_reduce(fn, seq[, nthreads])

Parallel reduce.

isreal(qob, **allclose_opts)

Checks if qob is approximately real.

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.

spin_operator(label[, S])

Generate a general spin-operator.

pauli(xyz[, dim])

Generates the pauli operators for dimension 2 or 3.

hadamard([dtype, sparse])

The Hadamard gate.

phase_gate(phi[, dtype, sparse])

The generalized qubit phase-gate, which adds phase phi to the

T_gate([dtype, sparse])

The T-gate (pi/8 gate).

S_gate([dtype, sparse])

The S-gate (phase gate).

rotation(phi[, xyz, dtype, sparse])

The single qubit rotation gate.

U_gate(theta, phi, lamda[, dtype, sparse])

Arbitrary unitary single qubit gate.

Xsqrt(**qu_opts)

Rx(pi / 2).

Ysqrt(**qu_opts)

Ry(pi / 2).

Zsqrt(**qu_opts)

Rz(pi / 2).

Wsqrt(**qu_opts)

R[X + Y](pi / 2).

swap([dim, dtype])

The SWAP operator acting on subsystems of dimension dim.

fsim(theta, phi[, dtype])

The 'fermionic simulation' gate:

fsimg(theta, zeta, chi, gamma, phi[, dtype])

The 'fermionic simulation' gate, with:

iswap([dtype])

ncontrolled_gate(ncontrol, gate[, dtype, sparse])

Build an n-qubit controlled gate. The control qubits are the

controlled(s[, dtype, sparse])

Construct a controlled pauli gate for two qubits.

CNOT([dtype, sparse])

The controlled-not gate.

cX([dtype, sparse])

The controlled-X gate.

cY([dtype, sparse])

The controlled-Y gate.

cZ([dtype, sparse])

The controlled-Z gate.

ccX([dtype, sparse])

The double controlled X gate, or Toffoli gate.

ccY([dtype, sparse])

The double controlled Y gate.

ccZ([dtype, sparse])

The double controlled Z gate.

controlled_swap([dtype, sparse])

The controlled swap or Fredkin gate. The control qubit is the first

hamiltonian_builder(fn)

Wrap a function to perform some generic postprocessing and take the

ham_heis(n[, j, b, cyclic, parallel, nthreads, ownership])

Constructs the nearest neighbour 1d heisenberg spin-1/2 hamiltonian.

ham_ising(n[, jz, bx])

Generate the quantum transverse field ising model hamiltonian. This is a

ham_XY(n, jxy, bz, **ham_opts)

Generate the quantum transverse field XY model hamiltonian. This is a

ham_XXZ(n, delta[, jxy])

Generate the XXZ-model hamiltonian. This is a

ham_j1j2(n[, j1, j2, bz, cyclic, ownership])

Generate the j1-j2 hamiltonian, i.e. next nearest neighbour

_gen_mbl_random_factors(n, dh, dh_dim, dh_dist[, ...])

ham_mbl(n, dh[, j, bz, cyclic, seed, dh_dist, dh_dim, ...])

Constructs a heisenberg hamiltonian with isotropic coupling and

ham_heis_2D(n, m[, j, bz, cyclic, parallel, ownership])

Construct the 2D spin-1/2 heisenberg model hamiltonian:

uniq_perms(xs)

Generate all the unique permutations of sequence xs.

zspin_projector(n[, sz, stype, dtype])

Construct the projector onto spin-z subpspaces.

create([n])

The creation operator acting on an n-level system.

destroy([n])

The annihilation operator acting on an n-level system.

num(n, **qu_opts)

The number operator acting on an n-level system.

ham_hubbard_hardcore(n[, t, V, mu, cyclic, parallel, ...])

Generate the spinless fermion hopping hamiltonian.

Module Contents

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

Make array read only, in-place.

Parameters:

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

quimb.gen.operators.get_thread_pool(num_workers=None)
quimb.gen.operators.par_reduce(fn, seq, nthreads=_NUM_THREAD_WORKERS)[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.

Return type:

depends on fn and seq.

Notes

This has a several hundred microsecond overhead.

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

Checks if qob is approximately real.

quimb.gen.operators.qu[source]

Alias of quimbify().

quimb.gen.operators.eye[source]

Alias for identity().

quimb.gen.operators.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.operators.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.operators.spin_operator(label, S=1 / 2, **kwargs)[source]

Generate a general spin-operator.

Parameters:
  • label (str) –

    The type of operator, can be one of six options:

    • {'x', 'X'}, x-spin operator.

    • {'y', 'Y'}, y-spin operator.

    • {'z', 'Z'}, z-spin operator.

    • {'+', 'p'}, Raising operator.

    • {'-', 'm'}, Lowering operator.

    • {'i', 'I'}, identity operator.

  • S (float, optional) – The spin of particle to act on, default to spin-1/2.

  • kwargs – Passed to quimbify().

Returns:

S – The spin operator.

Return type:

immutable operator

See also

pauli

Examples

>>> spin_operator('x')
qarray([[0. +0.j, 0.5+0.j],
        [0.5+0.j, 0. +0.j]])
>>> qu.spin_operator('+', S=1)
qarray([[0.        +0.j, 1.41421356+0.j, 0.        +0.j],
        [0.        +0.j, 0.        +0.j, 1.41421356+0.j],
        [0.        +0.j, 0.        +0.j, 0.        +0.j]])
>>> qu.spin_operator('Y', sparse=True)
<2x2 sparse matrix of type '<class 'numpy.complex128'>'
    with 2 stored elements in Compressed Sparse Row format>
quimb.gen.operators.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.operators.hadamard(dtype=complex, sparse=False)[source]

The Hadamard gate.

quimb.gen.operators.phase_gate(phi, dtype=complex, sparse=False)[source]

The generalized qubit phase-gate, which adds phase phi to the |1> state.

quimb.gen.operators.T_gate(dtype=complex, sparse=False)[source]

The T-gate (pi/8 gate).

quimb.gen.operators.S_gate(dtype=complex, sparse=False)[source]

The S-gate (phase gate).

quimb.gen.operators.rotation(phi, xyz='Z', dtype=complex, sparse=False)[source]

The single qubit rotation gate.

quimb.gen.operators.Rx
quimb.gen.operators.Ry
quimb.gen.operators.Rz
quimb.gen.operators.U_gate(theta, phi, lamda, dtype=complex, sparse=False)[source]

Arbitrary unitary single qubit gate.

\[\begin{split}U_3(\theta, \phi, \lambda) = \begin{bmatrix} \cos(\theta / 2) & - e^{i \lambda} \sin(\theta / 2) \\ e^{i \phi} \sin(\theta / 2) & e^{i(\lambda + \phi)}\cos(\theta / 2) \end{bmatrix}\end{split}\]
Parameters:
  • theta (float) – Angle between 0 and pi.

  • phi (float) – Angle between 0 and 2 pi.

  • lamba (float) – Angle between 0 and 2 pi.

Returns:

U – The unitary matrix, cached.

Return type:

(2, 2) array

quimb.gen.operators.Xsqrt(**qu_opts)[source]

Rx(pi / 2).

\[\begin{split}X^{\frac{1}{2}} = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 & - i \\ - i & 1 \end{bmatrix}\end{split}\]
quimb.gen.operators.Ysqrt(**qu_opts)[source]

Ry(pi / 2).

\[\begin{split}Y^{\frac{1}{2}} = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 & - 1 \\ 1 & 1 \end{bmatrix}\end{split}\]
quimb.gen.operators.Zsqrt(**qu_opts)[source]

Rz(pi / 2).

\[\begin{split}Z^{\frac{1}{2}} = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 - i & 0 \\ 0 & 1 + i \end{bmatrix}\end{split}\]
quimb.gen.operators.Wsqrt(**qu_opts)[source]

R[X + Y](pi / 2).

\[\begin{split}W^{\frac{1}{2}} = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 & -\sqrt{i} \\ \sqrt{-i} & 1 \end{bmatrix}\end{split}\]
quimb.gen.operators.swap(dim=2, dtype=complex, **kwargs)[source]

The SWAP operator acting on subsystems of dimension dim.

quimb.gen.operators.fsim(theta, phi, dtype=complex, **kwargs)[source]

The ‘fermionic simulation’ gate:

\[\begin{split}\mathrm{fsim}(\theta, \phi) = \begin{bmatrix} 1 & 0 & 0 & 0\\ 0 & \cos(\theta) & -i sin(\theta) & 0\\ 0 & -i sin(\theta) & \cos(\theta) & 0\\ 0 & 0 & 0 & \exp(-i \phi) \end{bmatrix}\end{split}\]

Note that theta and phi should be specified in radians and the sign convention with this gate varies. Here for example, fsim(- pi / 2, 0) == iswap().

quimb.gen.operators.fsimg(theta, zeta, chi, gamma, phi, dtype=complex, **kwargs)[source]

The ‘fermionic simulation’ gate, with:

  • \(\theta\) is the iSWAP angle

  • \(\phi\) is the controlled-phase angle

  • \(\zeta, \chi, \gamma\) are single-qubit phase angles.

\[\begin{split}\mathrm{fsimg}(\theta, \zeta, \chi, \gamma, \phi) = \begin{bmatrix} 1 & 0 & 0 & 0\\ 0 & \exp(-i(\gamma +\zeta )) \cos(\theta) & -i \exp(-i(\gamma - \chi )) sin(\theta) & 0\\ 0 & -i \exp(-i(\gamma + \chi )) sin(\theta) & \exp(-i(\gamma - \zeta )) \cos(\theta) & 0\\ 0 & 0 & 0 & \exp(-i (\phi +2 \gamma)) \end{bmatrix}\end{split}\]

See Equation 18 of https://arxiv.org/abs/2010.07965. Note that theta, phi, zeta, chi, gamma should be specified in radians and the sign convention with this gate varies. Here for example, fsimg(- pi / 2, 0, 0, 0,0) == iswap().

quimb.gen.operators.iswap(dtype=complex, **kwargs)[source]
quimb.gen.operators.ncontrolled_gate(ncontrol, gate, dtype=complex, sparse=False)[source]

Build an n-qubit controlled gate. The control qubits are the first ncontrol qubits.

Parameters:
  • ncontrol (int) – The number of control qubits.

  • gate (array_like) – The gate to apply to the controlled qubit(s).

  • dtype (str) – The data type of the returned matrix.

  • sparse (bool) – Whether to return a sparse matrix.

Returns:

C – The n-qubit controlled gate.

Return type:

qarray

quimb.gen.operators.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.operators.CNOT(dtype=complex, sparse=False)[source]

The controlled-not gate.

quimb.gen.operators.cX(dtype=complex, sparse=False)[source]

The controlled-X gate.

quimb.gen.operators.cY(dtype=complex, sparse=False)[source]

The controlled-Y gate.

quimb.gen.operators.cZ(dtype=complex, sparse=False)[source]

The controlled-Z gate.

quimb.gen.operators.ccX(dtype=complex, sparse=False)[source]

The double controlled X gate, or Toffoli gate.

quimb.gen.operators.ccY(dtype=complex, sparse=False)[source]

The double controlled Y gate.

quimb.gen.operators.ccZ(dtype=complex, sparse=False)[source]

The double controlled Z gate.

quimb.gen.operators.controlled_swap(dtype=complex, sparse=False)[source]

The controlled swap or Fredkin gate. The control qubit is the first qubit, if in state |1> a swap is performed on the last two qubits.

quimb.gen.operators.cswap[source]
quimb.gen.operators.fredkin[source]
quimb.gen.operators.toffoli[source]
quimb.gen.operators.hamiltonian_builder(fn)[source]

Wrap a function to perform some generic postprocessing and take the kwargs stype and sparse. This assumes the core function always builds the hamiltonian in sparse form. The wrapper then:

  1. Checks if the operator is real and, if so, discards imaginary part if no explicity dtype was given

  2. Converts the operator to dense or the correct sparse form

  3. Makes the operator immutable so it can be safely cached

quimb.gen.operators.ham_heis(n, j=1.0, b=0.0, cyclic=False, parallel=False, nthreads=None, ownership=None)[source]

Constructs the nearest neighbour 1d heisenberg spin-1/2 hamiltonian.

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

  • j (float or tuple(float, float, float), optional) – Coupling constant(s), with convention that positive = antiferromagnetic. Can supply scalar for isotropic coupling or vector (jx, jy, jz).

  • b (float or tuple(float, float, float), optional) – Magnetic field, defaults to z-direction only if tuple not given.

  • cyclic (bool, optional) – Whether to couple the first and last spins.

  • sparse (bool, optional) – Whether to return the hamiltonian in sparse form.

  • stype (str, optional) – What format of sparse operator to return if sparse.

  • parallel (bool, optional) – Whether to build the operator in parallel. By default will do this for n > 16.

  • nthreads (int optional) – How mny threads to use in parallel to build the operator.

  • ownership ((int, int), optional) – If given, which range of rows to generate.

  • kwargs – Supplied to quimbify().

Returns:

H – The Hamiltonian.

Return type:

immutable operator

quimb.gen.operators.ham_ising(n, jz=1.0, bx=1.0, **ham_opts)[source]

Generate the quantum transverse field ising model hamiltonian. This is a simple alias for ham_heis() with Z-interactions and an X-field.

quimb.gen.operators.ham_XY(n, jxy, bz, **ham_opts)[source]

Generate the quantum transverse field XY model hamiltonian. This is a simple alias for ham_heis() with X- and Y-interactions and a Z-field.

quimb.gen.operators.ham_XXZ(n, delta, jxy=1.0, **ham_opts)[source]

Generate the XXZ-model hamiltonian. This is a simple alias for ham_heis() with matched X- and Y-interactions and delta Z coupling.

quimb.gen.operators.ham_j1j2(n, j1=1.0, j2=0.5, bz=0.0, cyclic=False, ownership=None)[source]

Generate the j1-j2 hamiltonian, i.e. next nearest neighbour interactions.

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

  • j1 (float, optional) – Nearest neighbour coupling strength.

  • j2 (float, optional) – Next nearest neighbour coupling strength.

  • bz (float, optional) – B-field strength in z-direction.

  • cyclic (bool, optional) – Cyclic boundary conditions.

  • sparse (bool, optional) – Return hamiltonian as sparse-csr operator.

  • ownership ((int, int), optional) – If given, which range of rows to generate.

  • kwargs – Supplied to quimbify().

Returns:

H – The Hamiltonian.

Return type:

immutable operator

quimb.gen.operators._gen_mbl_random_factors(n, dh, dh_dim, dh_dist, seed=None, beta=None)[source]
quimb.gen.operators.ham_mbl(n, dh, j=1.0, bz=0.0, cyclic=False, seed=None, dh_dist='s', dh_dim=1, beta=None, ownership=None)[source]

Constructs a heisenberg hamiltonian with isotropic coupling and random fields acting on each spin - the many-body localized (MBL) spin hamiltonian.

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

  • dh (float or (float, float, float)) – Strength of random fields (stdev of gaussian distribution), can be scalar (isotropic noise) or 3-vector for (x, y, z) directions.

  • j (float or (float, float, float), optional) – Coupling strength, can be scalar (isotropic) or 3-vector.

  • bz (float, optional) – Global magnetic field (in z-direction).

  • cyclic (bool, optional) – Whether to use periodic boundary conditions.

  • seed (int, optional) – Number to seed random number generator with.

  • dh_dist ({'g', 's', 'qr'}, optional) –

    Type of random distribution for the noise:

    • ”s”: square, with bounds (-dh, dh)

    • ”g”: gaussian, with standard deviation dh

    • ”qp”: quasi periodic, with amplitude dh and ‘wavenumber’ beta so that the field at site i is dh * cos(2 * pi * beta * i + delta) with delta a random offset between (0, 2 * pi), possibly seeded by seed.

  • dh_dim ({1, 2, 3} or str, optional) – The number of dimensions the noise acts in, or string specifier like 'yz'.

  • beta (float, optional) – The wave number if dh_dist='qr', defaults to the golden ratio``(5**0.5 - 1) / 2``.

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

  • stype ({'csr', 'csc', 'coo'}, optional) – The sparse format.

  • ownership ((int, int), optional) – If given, which range of rows to generate.

  • kwargs – Supplied to quimbify().

Returns:

H – The MBL hamiltonian for spin-1/2.

Return type:

operator

See also

MPO_ham_mbl

quimb.gen.operators.ham_heis_2D(n, m, j=1.0, bz=0.0, cyclic=False, parallel=False, ownership=None)[source]

Construct the 2D spin-1/2 heisenberg model hamiltonian:

\[\hat{H} = \sum_{<i, j>} J_X S^X_i S^X_j + J_Y S^Y_i S^Y_j + J_Z S^Z_i S^Z_j\]

where the sum runs over pairs \(<i,j>\) on a 2D square lattice.

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

  • m (int) – The number of columns.

  • j (float or (float, float, float), optional) – The coupling strength(s). Isotropic if scalar else if vector (Jx, Jy, Jz) = j.

  • bz (float, optional) – The z direction magnetic field.

  • cyclic (bool, optional) – Whether to use periodic boundary conditions.

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

  • stype ({'csr', 'csc', 'coo'}, optional) – The sparse format.

  • parallel (bool, optional) – Construct the hamiltonian in parallel. Faster but might use more memory.

  • ownership ((int, int), optional) – If given, which range of rows to generate.

  • kwargs – Supplied to quimbify().

Returns:

H – The hamiltonian.

Return type:

operator

quimb.gen.operators.uniq_perms(xs)[source]

Generate all the unique permutations of sequence xs.

Examples

>>> list(uniq_perms('0011'))
[('0', '0', '1', '1'),
 ('0', '1', '0', '1'),
 ('0', '1', '1', '0'),
 ('1', '0', '0', '1'),
 ('1', '0', '1', '0'),
 ('1', '1', '0', '0')]
quimb.gen.operators.zspin_projector(n, sz=0, stype='csr', dtype=float)[source]

Construct the projector onto spin-z subpspaces.

Parameters:
  • n (int) – Total size of spin system.

  • sz (float or sequence of floats) – Spin-z value(s) subspace(s) to find projector for.

  • stype (str) – Sparse format of the output operator.

  • dtype ({float, complex}, optional) – The data type of the operator to generate.

Returns:

prj – The (non-square) projector onto the specified subspace(s). The subspace size D is given by n choose (n / 2 + s) for each s specified in sz.

Return type:

immutable sparse operator, shape (2**n, D)

Examples

>>> zspin_projector(n=2, sz=0).A
array([[0., 0.],
       [1., 0.],
       [0., 1.],
       [0., 0.]]

Project a 9-spin Heisenberg-Hamiltonian into its spin-1/2 subspace:

>>> H = ham_heis(9, sparse=True)
>>> H.shape
(512, 512)
>>> P = zspin_projector(n=9, sz=1 / 2)
>>> H0 = P.T @ H @ P
>>> H0.shape
(126, 126)
quimb.gen.operators.create(n=2, **qu_opts)[source]

The creation operator acting on an n-level system.

quimb.gen.operators.destroy(n=2, **qu_opts)[source]

The annihilation operator acting on an n-level system.

quimb.gen.operators.num(n, **qu_opts)[source]

The number operator acting on an n-level system.

quimb.gen.operators.ham_hubbard_hardcore(n, t=0.5, V=1.0, mu=1.0, cyclic=False, parallel=False, ownership=None)[source]

Generate the spinless fermion hopping hamiltonian.

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

  • t (float, optional) – The hopping energy.

  • V (float, optional) – The interaction energy.

  • mu (float, optional) – The chemical potential - defaults to half-filling.

  • cyclic (bool, optional) – Whether to use periodic boundary conditions.

  • parallel (bool, optional) – Construct the hamiltonian in parallel. Faster but might use more memory.

  • ownership ((int, int), optional) – If given, which range of rows to generate.

  • kwargs – Supplied to quimbify().

Returns:

H – The hamiltonian.

Return type:

operator