quimb.gen.operators =================== .. py:module:: quimb.gen.operators .. autoapi-nested-parse:: Functions for generating quantum operators. Attributes ---------- .. autoapisummary:: quimb.gen.operators.Rx quimb.gen.operators.Ry quimb.gen.operators.Rz quimb.gen.operators.cswap quimb.gen.operators.fredkin quimb.gen.operators.toffoli Functions --------- .. autoapisummary:: quimb.gen.operators.spin_operator quimb.gen.operators.pauli quimb.gen.operators.hadamard quimb.gen.operators.phase_gate quimb.gen.operators.T_gate quimb.gen.operators.S_gate quimb.gen.operators.rotation quimb.gen.operators.U_gate quimb.gen.operators.Xsqrt quimb.gen.operators.Ysqrt quimb.gen.operators.Zsqrt quimb.gen.operators.Wsqrt quimb.gen.operators.swap quimb.gen.operators.fsim quimb.gen.operators.fsimg quimb.gen.operators.iswap quimb.gen.operators.ncontrolled_gate quimb.gen.operators.controlled quimb.gen.operators.CNOT quimb.gen.operators.cX quimb.gen.operators.cY quimb.gen.operators.cZ quimb.gen.operators.ccX quimb.gen.operators.ccY quimb.gen.operators.ccZ quimb.gen.operators.controlled_swap quimb.gen.operators.hamiltonian_builder quimb.gen.operators.ham_heis quimb.gen.operators.ham_ising quimb.gen.operators.ham_XY quimb.gen.operators.ham_XXZ quimb.gen.operators.ham_j1j2 quimb.gen.operators._gen_mbl_random_factors quimb.gen.operators.ham_mbl quimb.gen.operators.ham_heis_2D quimb.gen.operators.uniq_perms quimb.gen.operators.zspin_projector quimb.gen.operators.create quimb.gen.operators.destroy quimb.gen.operators.num quimb.gen.operators.ham_hubbard_hardcore Module Contents --------------- .. py:function:: spin_operator(label, S=1 / 2, **kwargs) Generate a general spin-operator. :param label: 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. :type label: str :param S: The spin of particle to act on, default to spin-1/2. :type S: float, optional :param kwargs: Passed to :func:`quimbify`. :returns: **S** -- The spin operator. :rtype: immutable operator .. seealso:: :py:obj:`pauli` .. rubric:: 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 '' with 2 stored elements in Compressed Sparse Row format> .. py:function:: pauli(xyz, dim=2, **kwargs) Generates the pauli operators for dimension 2 or 3. :param xyz: Which spatial direction, upper or lower case from ``{'I', 'X', 'Y', 'Z'}``. :type xyz: str :param dim: Dimension of spin operator (e.g. 3 for spin-1), defaults to 2 for spin half. :type dim: int, optional :param kwargs: Passed to ``quimbify``. :returns: **P** -- The pauli operator. :rtype: immutable operator .. seealso:: :py:obj:`spin_operator` .. py:function:: hadamard(dtype=complex, sparse=False) The Hadamard gate. .. py:function:: phase_gate(phi, dtype=complex, sparse=False) The generalized qubit phase-gate, which adds phase ``phi`` to the ``|1>`` state. .. py:function:: T_gate(dtype=complex, sparse=False) The T-gate (pi/8 gate). .. py:function:: S_gate(dtype=complex, sparse=False) The S-gate (phase gate). .. py:function:: rotation(phi, xyz='Z', dtype=complex, sparse=False) The single qubit rotation gate. .. py:data:: Rx .. py:data:: Ry .. py:data:: Rz .. py:function:: U_gate(theta, phi, lamda, dtype=complex, sparse=False) Arbitrary unitary single qubit gate. .. math:: 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} :param theta: Angle between 0 and pi. :type theta: float :param phi: Angle between 0 and 2 pi. :type phi: float :param lamba: Angle between 0 and 2 pi. :type lamba: float :returns: **U** -- The unitary matrix, cached. :rtype: (2, 2) array .. py:function:: Xsqrt(**qu_opts) Rx(pi / 2). .. math:: X^{\frac{1}{2}} = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 & - i \\ - i & 1 \end{bmatrix} .. py:function:: Ysqrt(**qu_opts) Ry(pi / 2). .. math:: Y^{\frac{1}{2}} = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 & - 1 \\ 1 & 1 \end{bmatrix} .. py:function:: Zsqrt(**qu_opts) Rz(pi / 2). .. math:: Z^{\frac{1}{2}} = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 - i & 0 \\ 0 & 1 + i \end{bmatrix} .. py:function:: Wsqrt(**qu_opts) R[X + Y](pi / 2). .. math:: W^{\frac{1}{2}} = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 & -\sqrt{i} \\ \sqrt{-i} & 1 \end{bmatrix} .. py:function:: swap(dim=2, dtype=complex, **kwargs) The SWAP operator acting on subsystems of dimension `dim`. .. py:function:: fsim(theta, phi, dtype=complex, **kwargs) The 'fermionic simulation' gate: .. math:: \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} 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()``. .. py:function:: fsimg(theta, zeta, chi, gamma, phi, dtype=complex, **kwargs) The 'fermionic simulation' gate, with: * :math:`\theta` is the iSWAP angle * :math:`\phi` is the controlled-phase angle * :math:`\zeta, \chi, \gamma` are single-qubit phase angles. .. math:: \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} 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()``. .. py:function:: iswap(dtype=complex, **kwargs) .. py:function:: ncontrolled_gate(ncontrol, gate, dtype=complex, sparse=False) Build an n-qubit controlled gate. The control qubits are the first ``ncontrol`` qubits. :param ncontrol: The number of control qubits. :type ncontrol: int :param gate: The gate to apply to the controlled qubit(s). :type gate: array_like :param dtype: The data type of the returned matrix. :type dtype: str :param sparse: Whether to return a sparse matrix. :type sparse: bool :returns: **C** -- The n-qubit controlled gate. :rtype: qarray .. py:function:: controlled(s, dtype=complex, sparse=False) Construct a controlled pauli gate for two qubits. :param s: Which pauli to use, including 'not' aliased to 'x'. :type s: str :param sparse: Whether to construct a sparse operator. :type sparse: bool, optional :returns: **C** -- The controlled two-qubit gate operator. :rtype: qarray .. py:function:: CNOT(dtype=complex, sparse=False) The controlled-not gate. .. py:function:: cX(dtype=complex, sparse=False) The controlled-X gate. .. py:function:: cY(dtype=complex, sparse=False) The controlled-Y gate. .. py:function:: cZ(dtype=complex, sparse=False) The controlled-Z gate. .. py:function:: ccX(dtype=complex, sparse=False) The double controlled X gate, or Toffoli gate. .. py:function:: ccY(dtype=complex, sparse=False) The double controlled Y gate. .. py:function:: ccZ(dtype=complex, sparse=False) The double controlled Z gate. .. py:function:: controlled_swap(dtype=complex, sparse=False) 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. .. py:data:: cswap .. py:data:: fredkin .. py:data:: toffoli .. py:function:: hamiltonian_builder(fn) 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 .. py:function:: ham_heis(n, j=1.0, b=0.0, cyclic=False, parallel=False, nthreads=None, ownership=None) Constructs the nearest neighbour 1d heisenberg spin-1/2 hamiltonian. :param n: Number of spins. :type n: int :param j: Coupling constant(s), with convention that positive = antiferromagnetic. Can supply scalar for isotropic coupling or vector ``(jx, jy, jz)``. :type j: float or tuple(float, float, float), optional :param b: Magnetic field, defaults to z-direction only if tuple not given. :type b: float or tuple(float, float, float), optional :param cyclic: Whether to couple the first and last spins. :type cyclic: bool, optional :param sparse: Whether to return the hamiltonian in sparse form. :type sparse: bool, optional :param stype: What format of sparse operator to return if ``sparse``. :type stype: str, optional :param parallel: Whether to build the operator in parallel. By default will do this for n > 16. :type parallel: bool, optional :param nthreads: How mny threads to use in parallel to build the operator. :type nthreads: int optional :param ownership: If given, which range of rows to generate. :type ownership: (int, int), optional :param kwargs: Supplied to :func:`~quimb.core.quimbify`. :returns: **H** -- The Hamiltonian. :rtype: immutable operator .. py:function:: ham_ising(n, jz=1.0, bx=1.0, **ham_opts) Generate the quantum transverse field ising model hamiltonian. This is a simple alias for :func:`~quimb.gen.operators.ham_heis` with Z-interactions and an X-field. .. py:function:: ham_XY(n, jxy, bz, **ham_opts) Generate the quantum transverse field XY model hamiltonian. This is a simple alias for :func:`~quimb.gen.operators.ham_heis` with X- and Y-interactions and a Z-field. .. py:function:: ham_XXZ(n, delta, jxy=1.0, **ham_opts) Generate the XXZ-model hamiltonian. This is a simple alias for :func:`~quimb.gen.operators.ham_heis` with matched X- and Y-interactions and ``delta`` Z coupling. .. py:function:: ham_j1j2(n, j1=1.0, j2=0.5, bz=0.0, cyclic=False, ownership=None) Generate the j1-j2 hamiltonian, i.e. next nearest neighbour interactions. :param n: Number of spins. :type n: int :param j1: Nearest neighbour coupling strength. :type j1: float, optional :param j2: Next nearest neighbour coupling strength. :type j2: float, optional :param bz: B-field strength in z-direction. :type bz: float, optional :param cyclic: Cyclic boundary conditions. :type cyclic: bool, optional :param sparse: Return hamiltonian as sparse-csr operator. :type sparse: bool, optional :param ownership: If given, which range of rows to generate. :type ownership: (int, int), optional :param kwargs: Supplied to :func:`~quimb.core.quimbify`. :returns: **H** -- The Hamiltonian. :rtype: immutable operator .. py:function:: _gen_mbl_random_factors(n, dh, dh_dim, dh_dist, seed=None, beta=None) .. py:function:: ham_mbl(n, dh, j=1.0, bz=0.0, cyclic=False, seed=None, dh_dist='s', dh_dim=1, beta=None, ownership=None) Constructs a heisenberg hamiltonian with isotropic coupling and random fields acting on each spin - the many-body localized (MBL) spin hamiltonian. :param n: Number of spins. :type n: int :param dh: Strength of random fields (stdev of gaussian distribution), can be scalar (isotropic noise) or 3-vector for (x, y, z) directions. :type dh: float or (float, float, float) :param j: Coupling strength, can be scalar (isotropic) or 3-vector. :type j: float or (float, float, float), optional :param bz: Global magnetic field (in z-direction). :type bz: float, optional :param cyclic: Whether to use periodic boundary conditions. :type cyclic: bool, optional :param seed: Number to seed random number generator with. :type seed: int, optional :param dh_dist: 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``. :type dh_dist: {'g', 's', 'qr'}, optional :param dh_dim: The number of dimensions the noise acts in, or string specifier like ``'yz'``. :type dh_dim: {1, 2, 3} or str, optional :param beta: The wave number if ``dh_dist='qr'``, defaults to the golden ratio``(5**0.5 - 1) / 2``. :type beta: float, optional :param sparse: Whether to construct the hamiltonian in sparse form. :type sparse: bool, optional :param stype: The sparse format. :type stype: {'csr', 'csc', 'coo'}, optional :param ownership: If given, which range of rows to generate. :type ownership: (int, int), optional :param kwargs: Supplied to :func:`~quimb.core.quimbify`. :returns: **H** -- The MBL hamiltonian for spin-1/2. :rtype: operator .. seealso:: :py:obj:`MPO_ham_mbl` .. py:function:: ham_heis_2D(n, m, j=1.0, bz=0.0, cyclic=False, parallel=False, ownership=None) Construct the 2D spin-1/2 heisenberg model hamiltonian: .. math:: \hat{H} = \sum_{} 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 :math:`` on a 2D square lattice. :param n: The number of rows. :type n: int :param m: The number of columns. :type m: int :param j: The coupling strength(s). Isotropic if scalar else if vector ``(Jx, Jy, Jz) = j``. :type j: float or (float, float, float), optional :param bz: The z direction magnetic field. :type bz: float, optional :param cyclic: Whether to use periodic boundary conditions. :type cyclic: bool, optional :param sparse: Whether to construct the hamiltonian in sparse form. :type sparse: bool, optional :param stype: The sparse format. :type stype: {'csr', 'csc', 'coo'}, optional :param parallel: Construct the hamiltonian in parallel. Faster but might use more memory. :type parallel: bool, optional :param ownership: If given, which range of rows to generate. :type ownership: (int, int), optional :param kwargs: Supplied to :func:`~quimb.core.quimbify`. :returns: **H** -- The hamiltonian. :rtype: operator .. py:function:: uniq_perms(xs) Generate all the unique permutations of sequence ``xs``. .. rubric:: 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')] .. py:function:: zspin_projector(n, sz=0, stype='csr', dtype=float) Construct the projector onto spin-z subpspaces. :param n: Total size of spin system. :type n: int :param sz: Spin-z value(s) subspace(s) to find projector for. :type sz: float or sequence of floats :param stype: Sparse format of the output operator. :type stype: str :param dtype: The data type of the operator to generate. :type dtype: {float, complex}, optional :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``. :rtype: immutable sparse operator, shape (2**n, D) .. rubric:: Examples >>> zspin_projector(n=2, sz=0).toarray() 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) .. py:function:: create(n=2, **qu_opts) The creation operator acting on an n-level system. .. py:function:: destroy(n=2, **qu_opts) The annihilation operator acting on an n-level system. .. py:function:: num(n, **qu_opts) The number operator acting on an n-level system. .. py:function:: ham_hubbard_hardcore(n, t=0.5, V=1.0, mu=1.0, cyclic=False, parallel=False, ownership=None) Generate the spinless fermion hopping hamiltonian. :param n: The number of sites. :type n: int :param t: The hopping energy. :type t: float, optional :param V: The interaction energy. :type V: float, optional :param mu: The chemical potential - defaults to half-filling. :type mu: float, optional :param cyclic: Whether to use periodic boundary conditions. :type cyclic: bool, optional :param parallel: Construct the hamiltonian in parallel. Faster but might use more memory. :type parallel: bool, optional :param ownership: If given, which range of rows to generate. :type ownership: (int, int), optional :param kwargs: Supplied to :func:`~quimb.core.quimbify`. :returns: **H** -- The hamiltonian. :rtype: operator