Groundstate energy of the Heisenberg S=1/2 chain, with momentum conservation

Computing the ground state energy of the S=1/2 chain

We compute the ground state energy of the \(S=1/2\) Heisenberg chain with periodic boundary conditions in each momentum sector \(k\). The Hamiltonian is given by

\[H = J\sum_{\langle i,j \rangle} \mathbf{S}_i \cdot \mathbf{S}_j\]

where \(\mathbf{S}_i = (S_i^x, S_i^y, S_i^z)\) are the spin \(S=1/2\) operators and \(\langle i,j \rangle\) denotes summation over nearest-meighbor sites \(i\) and \(j\). We set \(J=1\).

We start by defining the Hamiltonian.

  BondList bonds;
  for (int i = 0; i < n_sites; ++i) {
    bonds << Bond("HB", "J", {i, (i + 1) % n_sites});
  }
  bonds["J"] = 1.0;

Next, we need to create the symmetry group, which consists of all translations of the lattice. To create this group, we first define a Permutation object, which is the generating translation \(T: i \rightarrow i+1\). This is done by:

  std::vector<int> translation;
  for (int s = 0; s < n_sites; ++s) {
    translation.push_back((s + 1) % n_sites);
  }
  Permutation perm(translation);

To create the full group of all translations out of this generator, we use the function generated_group

  auto group = generated_group(perm);

Next we can loop over all momenta. To define the irreducible representation (irrep) of the momentum \(k\) we use the function generated_irrep. It returns an object of type Representation.

    complex phase = exp(2i * pi * k / n_sites);
    auto irrep = generated_irrep(perm, phase);

Finally, we simply compute the groundstate energy in this block using

    auto block = Spinhalf(n_sites, nup, group, irrep);
    double e0 = eig0(bonds, block);

The full working example is given by this code:

#include <hydra/all.h>

int main() {
  using namespace hydra;

  int n_sites = 16;
  int nup = n_sites / 2;

  // Define the nearest-neighbor Heisenberg Hamiltonian
  BondList bonds;
  for (int i = 0; i < n_sites; ++i) {
    bonds << Bond("HB", "J", {i, (i + 1) % n_sites});
  }
  bonds["J"] = 1.0;

  // Create the permutation group
  std::vector<int> translation;
  for (int s = 0; s < n_sites; ++s) {
    translation.push_back((s + 1) % n_sites);
  }
  Permutation perm(translation);
  auto group = generated_group(perm);

  // Loop over momenta k
  for (int k = 0; k < n_sites; ++k) {

    // Create irrep with momentum k
    complex phase = exp(2i * pi * k / n_sites);
    auto irrep = generated_irrep(perm, phase);

    // Compute the groundstate energy
    auto block = Spinhalf(n_sites, nup, group, irrep);
    double e0 = eig0(bonds, block);
    Log("k: {}, e0: {}", k, e0);
  }

  return EXIT_SUCCESS;
}