Numerical Calculation of Berry Curvature I am trying to calculate some berry curvature (BC) in a 2D lattice and I have some things I am getting lost with.
In the 2D lattice, we set up the eigenvalue problem $H|u_1\rangle = \epsilon_i|u_i\rangle$. Numerically, I can find $\epsilon,u$ by diagonalizing some matrix.  This gives my energy bands E(qx,qy).
When I want to calculate the BC, the formula is 
$$\Omega_i(q)=i\left(\langle \partial_{qx}u_i|\partial_{qy}u_i \rangle-\langle \partial_{qy}u_i|\partial_{qx}u_i \rangle\right)$$
I am kind of confused as to what $u_i$ I should be using.  I had just been using $u_i=E_i$, so taking the derivatives of the energy.  Should I be using the eigenvectors.  If I assume bloch waves, then I know there is the $e^{iqr}$ factor, but I want the derivatives in k space, so I am just not even sure what some pseudo code would look like.
As an aside, I am following this thesis which does the Harper model for a flux of 1/5.
So the prodedure is:
for each qx,qy:


*

*Find E(q), $\vec{u(q)}$
from here I just dont know what to use for u.

 A: In some cases the problem can be found analytically, but otherwise you don't really need to think about Bloch waves - especially if you're looking for spectra and properties within the first Brillouin zone. For numerical Berry curvature calculations, it usually makes most sense to first introduce a grid or lattice to discretize the relevant momentum space (i.e. the first Brillouin zone). Then, at each point in this lattice you can think of the Hamiltonian as a matrix $H(q_x,q_y)$, allowing you to obtain the energy eigenvalue for the $n$:th band $\epsilon_n(q_x,q_y)$ and eigenvector $u_n(q_x,q_y)$. (Which band is the $n$:th one can be decided e.g. by sorting the bands by energy.) Next, you can just introduce numerical derivatives using finite differences.
Unfortunately, the above method often gives bad answers. The reason is you need a smooth gauge for the wave function to evaluate this difference, and the numerical eigensolver often fails to produce this. Even when it works, it tends to require a very fine grid to approach the limit of $\delta q\rightarrow 0$, especially for complicated Hamiltonians. That's why you'll likely want to use a manifestly gauge invariant expression for the Berry curvature, and an efficient sampling method. Both of these problems are avoided if you use the so-called Fukui method, introduced in 

Fukui, Hatsugai, and Suzuki: Chern Numbers in Discretized Brillouin Zone: Efficient Method of Computing (Spin) Hall Conductances J. Phys. Soc. Jpn. 74, pp 1674-1677 (2005).

The published version is available here and a non-paywalled arXiv version can be found here.
On further investigation, it seems that in the dissertation you linked, M. Aidelsburger uses this very method.
