10_Quantum_Hubbard

=Quantum statistics and Hamiltonians= Let us first note that the Hamiltonian in the Schrödinger equation

math H \Psi = E \Psi math

is independent of the symmetry of the wavefunction. For a many-particle system, the distinguishable case makes no symmetry assumptions on Ψ, and it follows that the spectrum of the distinguishable-particle system contains the spectrum of the symmetric (bosonic) and the anti-symmetric (fermionic) system. It contains also states of mixed symmetry, neither bosonic nor fermionic. In an earlier lecture, we saw how the symmetry requirement was introduced into the density matrix via a sum over permutations. Here, we see how the symmetry requirement is realized in first and second quantization for lattice systems of distinguishable particles, of bosons, and of fermions.

NB: In this lecture, particles (even fermions) are spinless.

An exact diagonalization algorithm for distinguishable quantum lattice particles
Here is an exact diagonalization algorithm for the one-dimensional Hubbard model of //N// spin-less distinguishable particles on //L// lattice sites with periodic boundary conditions, with a potential energy math \sum_k U n_k(n_k-1) math and a hopping term -//t//, if there is one particle hopping from one site to one of its neighbors. code format="python" import numpy def general_count(i_min,i_max,N): I_vec = [i_min] * N   states=[] condition = False while condition == False: states.append(I_vec[:]) for k in range(N-1,-2,-1): if k == -1: condition = True if I_vec[k] == i_max: I_vec[k] = i_min else: I_vec[k] += 1 break return states
 * 1) enumerate states: each particle can sit on any of the N sites of the lattice
 * 1) enumerate states: each particle can sit on any of the N sites of the lattice

U = 4 t = 0.75 L = 4 N = 3 n_states = L**N a_mat = numpy.zeros((n_states,n_states)) nbr = {}

nbr[0] = [L-1,1] nbr[L-1] = [L-2,0] for i in range(1, L - 1): nbr[i] = (i - 1, i + 1) print nbr states = general_count(0, L - 1, N) print states for i in states: ii = states.index(i) a_mat[ii][ii] = U * sum(i.count(k) * (i.count(k) - 1) for k in range(L)) for k in range(len(i)): for l in nbr[i[k]]: j= i[:] j[k] = l           print i, j            jj = states.index(j) a_mat[ii][jj] = -t print a_mat w, v = numpy.linalg.eigh(a_mat) print " Spectrum of H" for x in w:   print x code

Here is the (partial) list of the 64 states generated by the above program: code format="python" [0, 0, 0] [0, 0, 1] [0, 0, 2] [0, 0, 3] [0, 1, 0] [0, 1, 1] [0, 1, 2] ... code

Here is the (partial) 64 × 64 Hamiltonian matrix: code format="python" 24.   -0.75  0.   ..., 0. 0. 0. ] [  -0.75   8.   -0.75 ..., 0. 0. 0. ] [   0.    -0.75  8.   ..., 0. 0. 0. ] ..., [ 0. 0. 0. ..., 8. -0.75 0. ] [ 0. 0. 0. ..., -0.75 8. -0.75] [ 0. 0. 0. ..., 0. -0.75 24.  code

Here, the element //H// 00 = 24, because three particles are on one site (0), and math H_{01} = -3/4 math because particle no 3 hops from site 0 to site 1.

An exact diagonalization algorithm for quantum lattice bosons
Here is an equivalent exact diagonalization algorithm for the one-dimensional Hubbard model of spin-less bosonic particles. Second quantized notations are used, so that the hamiltonian is now

math H = - t \sum_{\langle i,j\rangle}c_i^+ c_j + \text{h.c.} + U \sum_i n_i (n_i-1) math We note the combinatorial factors math c_i^+ |n_i\rangle = \sqrt{n_i+1} |(n_i+1)\rangle math and math c_i |n_i\rangle = \sqrt{n_i} |(n_i-1)\rangle math (see R. P. Feynman Statistical mechanics: a set of lectures (Benjamin Cummings, 1972) to understand the square roots). Besides the diagonal terms, the Hamiltonian connects only states with correspond to a single hop of a particle from site //i// to a neighboring site //j//.

code format="python" import numpy from math import sqrt def general_count(i_min, i_max, L,N): I_vec = [i_min] * L   states = [] condition = False while condition == False: if sum(I_vec) == N: states.append(I_vec[ : ]) for k in range(L-1, -2, -1): if k == -1: condition = True if I_vec[k] == i_max: I_vec[k] = i_min else: I_vec[k] += 1 break return states U = 4.0 t = 0.75 L = 4 N = 3 nbr = {} nbr[0] = [L - 1, 1] nbr[L - 1] = [L - 2, 0] for i in range(1, L - 1): nbr[i] = (i-1, i + 1) print nbr states = general_count(0, N, L, N) a_mat = numpy.zeros((len(states), len(states))) for i in states: ii = states.index(i) a_mat[ii][ii] = U * sum(k * (k-1) for k in i)   for k in range(L): for l in nbr[k]: j = i[ : ] if (j[l] > 0): hopping = sqrt((i[k] + 1) * i[l]) j[k] = i[k] + 1 j[l] = i[l] - 1 print j               jj = states.index(j) a_mat[ii][jj] = -hopping * t print a_mat w, v = numpy.linalg.eigh(a_mat) for x in w:   print x code

Here, the states have length 4, and there are 20 of them ( see W. Krauth Statistical mechanics: Algorithms and Computations (Oxford University Press, 2006) See also the book's website eq. 4.12) code format="python" [0, 0, 0, 3] [0, 0, 1, 2] [0, 0, 2, 1] [0, 0, 3, 0] [0, 1, 0, 2] [0, 1, 1, 1] [0, 1, 2, 0] ...

code

and the matrix is code format="python" [ 24.       -1.29903811 0.0         0.0        0.0        0.0 ... [ -1.29903811 8.0       -1.5         0.0       -0.75       0.0 ... [  0.0       -1.5        8.0        -1.29903811 0.0       -1.06066017 ... [  0.0        0.0       -1.29903811 24.0        0.0        0.0 ... [  0.0       -0.75       0.0         0.0        8.0       -1.06066017 ... [  0.0        0.0       -1.06066017  0.0       -1.06066017 0.0 ... ... code

here, we note that math 1.299038.. = \frac{3}{4} \sqrt{3} math

Exact diagonalization algorithm for (spinless) lattice fermions
The Hubbard hamiltonian for lattice fermions is usually written for two spin species, and For fermions, the program is much simpler, because there are only 4 states (for 3 particles on 4 sites). Here is the program computing the Hamiltonian and the spectrum:

code format="python" import numpy import numpy.linalg U = 4.0 t = 0.75 L = 4 n_states = 4 a_mat = 0., -t, 0., -t], [-t, 0., -t, 0.], [0., -t, 0., -t], [-t, 0., -t, 0. print a_mat w, v = numpy.linalg.eigh(a_mat) for dummy in w:   print dummy code

Fermionic states are code format="python" (1, 1, 1, 0) (1, 1, 0, 1) (1, 0, 1, 1) (0, 1, 1, 1) code

The Hamiltonian matrix is code format="python" 0.0, -0.75, 0.0, -0.75], [ -0.75, 0.0, -1.0, 0.0], [  0.0, -0.75, 0.0, -0.75], [ -0.75, 0.0, -0.75, 0.0 code

The eigenvalues of this matrix are -1.5, 0, 0, 1.5

Comparing spectra


600px|left|visualizing two spectra with vimdiff Here we visualize the spectra of hubbard_bose.py (to the left) and of hubbard_dist.py (to the right), using an editor (vimdiff) capable of comparing two files. We see that the groundstate is of bosonic symmetry, but that not all of the states of the distinguishable hamiltonian appear either for the bosons or for the fermions.

Checking states
We discussed before that the groundstate of the distinguishable-particle problem had bosonic statistics, so let us check that this is indeed true: code format="python" -2.28331340813 Energy... (this is the distinguishable-particle                         and the bosonic groundstate) (0.0067261131197429158, [1, 1, 1]) (0.0067261131197430841, [3, 3, 3]) (0.0067261131197431639, [0, 0, 0]) (0.006726113119743274, [2, 2, 2]) cc (0.039285453143280247, [0, 1, 0]) (0.039285453143280691, [0, 0, 3]) (0.03928545314328085, [3, 0, 0]) (0.039285453143280878, [1, 0, 0]) (0.039285453143280892, [0, 3, 0]) (0.03928545314328092, [1, 1, 0]) (0.03928545314328092, [3, 0, 3]) (0.039285453143280927, [0, 3, 3]) (0.039285453143280927, [3, 2, 3]) (0.039285453143280968, [1, 1, 2]) (0.039285453143280975, [2, 1, 1]) (0.039285453143280982, [2, 3, 3]) (0.039285453143280989, [1, 0, 1]) (0.039285453143280989, [3, 3, 0]) (0.039285453143280989, [3, 3, 2]) (0.039285453143280996, [1, 2, 1]) (0.039285453143281024, [0, 1, 1]) (0.039285453143281059, [2, 2, 3]) (0.039285453143281086, [2, 1, 2]) (0.039285453143281086, [2, 2, 1]) (0.0392854531432811, [1, 2, 2]) (0.039285453143281107, [3, 2, 2]) (0.039285453143281128, [2, 3, 2]) (0.039285453143281177, [0, 0, 1]) ...etc code

and let's check, for fun, the fermionic statistics of the state with energy -1.5: code format="python" == -1.5 (this is the fermionic groundstate, and an excited state of the distinguishable-particle Hamiltonian) == (-0.20412414523193678, [3, 0, 2]) (-0.20412414523193539, [1, 2, 3]) (-0.20412414523193423, [3, 1, 2]) (-0.20412414523193395, [1, 2, 0]) (-0.20412414523193284, [1, 3, 0]) (-0.20412414523193179, [3, 0, 1]) (-0.20412414523193165, [0, 1, 2]) (-0.20412414523192998, [0, 2, 3]) (-0.20412414523192959, [0, 1, 3]) (-0.20412414523192821, [2, 0, 1]) (-0.20412414523192746, [2, 3, 0]) (-0.20412414523192479, [2, 3, 1]) (0.20412414523191447, [1, 0, 2]) (0.20412414523191599, [1, 0, 3]) (0.20412414523192099, [1, 3, 2]) (0.20412414523192221, [2, 0, 3]) (0.20412414523193054, [2, 1, 3]) (0.2041241452319319, [0, 3, 2]) (0.20412414523193434, [2, 1, 0]) (0.20412414523193689, [3, 1, 0]) (0.20412414523193881, [3, 2, 0]) (0.20412414523194106, [0, 3, 1]) (0.2041241452319435, [3, 2, 1]) (0.20412414523194503, [0, 2, 1]) code

Conclusion
Specifically, we have studied the Hubbard model (for 3 particles on 4 sites with periodic boundary conditions)
 * 1) for distinguishable particles (64 × 64 Hamiltonian matrix, first program)
 * 2) for bosons (20 × 20 Hamiltonian using creation and annihilation operators, second program)
 * 3) for fermions (4 × 4 Hamiltonian matrix, final program).

The most important observation was that the spectrum of the distinguishable-particle Hamiltonian contained all possible states, the symmetric ones, the antisymmetric ones, but also the states with mixed symmetry. We also notice that the ground state of the Fermionic problem, in fact, is an excited state of the distinguishable-particle Hamiltonian. This, in a sense, is the origin of the fermion sign problem. It is interesting to see how three different formulations yield different aspects of the same system: the quantum many-body problem. We restricted ourselves to a very small number of sites and of particles, because of the exponential increase in the size of the Hamiltonian matrix.

More specifically, we have performed a first socalled exact diagonalization of a quantum hamiltonian. Exact diagonalization is a mature field of research which allows to compute the spectrum, or parts of the spectrum (for example the ground state only) of lattice hamiltonians. For thermodynamics, it is mostly used for fermions. Bosons can be handled very successfully using Monte Carlo algorithms, so that exact diagonalization methods are less common. The main problem in exact diagonalization is the exponential growth of the Hilbert space, which seriously limits the system sizes which can be taken into account. The main advantage is, as the name indicates, that results are numerically exact, without error bars, so that extrapolation can be very powerful.

=Lattice Monte Carlo algorithms=

In earlier lectures, we treated the case of Monte Carlo algorithms in the continuum. Key elements were the convolution property, the explicit solution for the free density matrix, and the Trotter formula. These elements remain valid for the Monte Carlo

Lattice bosons
For lattice bosons, the procedure is the same as for the continuum system,

Lattice Fermions
For lattice fermions, we will treat a successful algorithm next time...

Quantum spins
... unfortunately not treated in this course...Quantum Monte Carlo, similar to what we learned earlier, is possible...

Computation of Fermion determinant
The formula of this section can be derived using Grassmann algebras. We prefer a pedestrian derivation, provided by J. E. Hirsch in the appendix of his famous paper J. E. Hirsch Two-dimensional Hubbard model: Numerical simulation study Physical Review B 31, 4403

math Z = \sum_{\text{states} \alpha } \langle \alpha |\exp(- \beta H) | \alpha \rangle math We first show how to treat bilinear terms, such as

math \text{Tr} \exp( - c_i^+ B_{ij} c_j) = \text{Tr} \prod_\mu \exp( - c_\mu^+ b_\mu c_\mu) math For this simple case, we can expand the exponential into

math \exp( - c_\mu^+ b_\mu c_\mu) = 1 - c_\mu^+ b_\mu c_\mu + \frac{1}{2}c_\mu^+ b_\mu c_\mu c_\mu^+ b_\mu c_\mu - \frac{1}{3!}..... math Using c_\mu^+ c_\mu + c_\mu c_\mu^+ = 1, we reach

math \exp( - c_\mu^+ b_\mu c_\mu) = 1 [+ c_\mu^+ b_\mu c_\mu ] + \frac{1}{2}c_\mu^+ b_\mu^2 c_\mu - \frac{1}{3!}c_\mu^+ b_\mu^3 c_\mu .... [-c_\mu^+ b_\mu c_\mu ] = \prod_\mu [ 1 + ( \exp (- b_\mu) - 1) c_\mu^+ c_\mu math This expression is easily evaluated. It yields the expression math \det[1 + \exp(-B)] math More generally, we find the expression:

math \text{Tr} \exp( - c_i^+ A_{ij} c_j)\exp( - c_i^+ B_{ij} c_j) = \det[1 + \exp(-A) \exp(-B) ] math

which can be derived from the action on single-particle states, or via Grassmann algebras. This non-trivial formula is the key to determinantal fermion methods.

=References=