.. _examples: Usage examples ============== Holstein model -------------- This example shows how to construct Hamiltonian of the `Holstein model `_ on a square lattice and to check that the total number of electrons is a conserved quantity of the model. The Hamiltonian considered here is a sum of three terms. * An electronic tight-binding model on a square lattice with only nearest-neighbour hopping allowed, .. math:: \hat H_\text{e} = -t \sum_\sigma \sum_{\langle i,j\rangle} c^\dagger_{i,\sigma} c_{j,\sigma}. * A harmonic oscillator at each lattice site (a localized phonon), .. math:: \hat H_\text{ph} = \omega_0 \sum_i a^\dagger_i a_i. * A coupling between the electrons and the phonons. .. math:: \hat H_\text{e-ph} = g \sum_\sigma \sum_i n_{i,\sigma}(a^\dagger_i + a_i). Instead of writing the sums over lattice sites explicitly, we call library functions :func:`tight_binding() `, :func:`dispersion() ` and :func:`holstein_int() `. We also make use of the `NetworkX `_ package to easily generate the adjacency matrix of the periodic square lattice. .. literalinclude:: examples/holstein.py :language: python :lines: 11- :linenos: Spin-1/2 Heisenberg chain ------------------------- The spin-1/2 Heisenberg chain is a textbook example of an integrable quantum system. Its Hamiltonian .. math:: \hat H = g \sum_i \mathbf{S}_i \cdot \mathbf{S}_{i+1} conserves three projections of the total spin .. math:: \mathbf{S} = \sum_i \mathbf{S}_i as well as a series of higher order charges :math:`Q_n`. Existence of these charges can be derived from the transfer matrix theory. Explicit expressions for :math:`Q_n` were obtained in [GM94]_. The following script constructs Hamiltonian of the Heisenberg chain with periodic boundary conditions (see :func:`heisenberg() `) and checks that :math:`[\hat H, \mathbf{S}] = 0`, :math:`[\hat H, Q_n] = 0` and :math:`[Q_n, Q_m] = 0` for :math:`m,n = 3,4,5`. .. literalinclude:: examples/heisenberg_chain.py :language: python :lines: 11- :linenos: Spectrum of Tavis-Cummings model -------------------------------- In this example we demonstrate how to use tools from the :mod:`loperator ` module along with NumPy's eigensolvers to diagonalize Hamiltonians of simple quantum models. We will consider a two-atom Jaynes-Cummings (Tavis-Cummings) Hamiltonian [TC68]_, which models a system of two qubits coherently coupled via a bosonic degree of freedom (quantum oscillator), .. math:: \hat H = \epsilon (\hat S_{z,0} + \hat S_{z,1}) + \omega \hat a^\dagger \hat a + g_0 (\hat a^\dagger \hat S_{-,0} + \hat a \hat S_{+,0}) + g_1 (\hat a^\dagger \hat S_{-,1} + \hat a \hat S_{+,1}). An expression of this form can be conveniently generated by a single call to :func:`jaynes_cummings() `. Since the occupation number of the bosonic mode can formally grow to :math:`+\infty`, we have to artificially restrict the dimension of the bosonic Hilbert space in order to be able to construct a finite matrix representation of our problem. .. note:: Performance of repeated calls to :class:`LOperatorR ` or to :class:`LOperatorC ` in order to construct a matrix representation of a linear operator is limited by Python method call overhead. For large-scale problems it is advised to call the utility function :func:`make_matrix() `. .. literalinclude:: examples/tavis_cummings.py :language: python :lines: 11- :linenos: Sectors of a multi-orbital interaction Hamiltonian -------------------------------------------------- When implementing an exact diagonalization algorithm, it is usually advantageous to take an extra preparatory step and split the Hilbert space of the problem into subspaces (sectors) characterized by a set of quantum numbers. In many situations, it is possible to find the sectors even without knowing the quantities conserved by the Hamiltonian (see, for instance, [SKFP16]_). In the script shown below we construct an electron-electron interaction Hamiltonian of an atomic :math:`d`-shell using :func:`Slater ` parametrization. Then, we employ the :class:`SpacePartition ` and :class:`BasisMapper ` utility classes (and optionally the :func:`make_matrix() ` function) to independently diagonalize the Hamiltonian within each sector. .. literalinclude:: examples/space_partition.py :language: python :lines: 11- :linenos: N-fermion sector views ---------------------- This example shows how to partially diagonalize the single band Fermi-Hubbard model on a large 2D square lattice with periodic boundary conditions. .. math:: \hat H = -t\sum_{\langle i,j \rangle, \sigma} c^\dagger_{i,\sigma} c_{j,\sigma} -\mu \sum_{i, \sigma} n_{i,\sigma} + U \sum_i n_{i,\uparrow} n_{i,\downarrow}. The number of atoms in the considered lattice (cluster) is set to 16, which makes for an intractably large Hilbert space of dimension :math:`2^{32}`. Since the Hamiltonian of the Hubbard model preserves total numbers of spin-up and spin-down electrons, one can diagonalize the Hamiltonian within a single :math:`N`-fermion sector or within a :math:`(N_\uparrow, N_\downarrow)`-multisector. An :math:`N`-fermion :py:class:`sector ` is a subspace of a full Hilbert space, which is spanned by all basis states with a fixed total occupation of fermionic degrees of freedom (FDOF). Similarly, a :py:class:`multisector ` is spanned by all basis states with a fixed occupation :math:`N_1` of a subset of the FDOF :math:`\{S_1\}`, occupation :math:`N_2` of another subset :math:`\{S_2\}` and so on. There can be any number of pairs :math:`(\{S_i\}, N_i)` (sectors contributing to the multisector) as long as all the subsets :math:`\{S_i\}` are disjoint. In our example we consider a moderately sized sector with :math:`N = 2` (:math:`\dim = {32 \choose 2} = 496`) and a multisector with :math:`N_\uparrow = 1, N_\downarrow = 1` (:math:`\dim = {16 \choose 1}{16 \choose 1} = 256`). .. note:: In general, the N-fermion (multi)sector views and functions do not require a purely fermionic system. The definitions of a sector and a multisector stand in presence of bosons or spins. .. literalinclude:: examples/n_fermion_sectors.py :language: python :lines: 16- :linenos: References ---------- .. [GM94] "Quantum Integrals of Motion for the Heisenberg Spin Chain", M. P. Grabowski and P. Mathieu, Mod. Phys. Lett. A, Vol. 09, No. 24, pp. 2197-2206 (1994), https://doi.org/10.1142/S0217732394002057 .. [TC68] "Exact Solution for an N-Molecule-Radiation-Field Hamiltonian", M. Tavis and F. W. Cummings Phys. Rev. 170, 379 (1968) https://doi.org/10.1103/PhysRev.170.379 .. [SKFP16] "TRIQS/CTHYB: A continuous-time quantum Monte Carlo hybridisation expansion solver for quantum impurity problems", P. Seth, I. Krivenko, M. Ferrero and O. Parcollet, Comp. Phys. Comm. 200, March 2016, 274-284, http://dx.doi.org/10.1016/j.cpc.2015.10.023 (section 4.2)