This script calculates the bare Coulomb matrix elements using real space Wannier functions (xsf files) $W_i(\mathbf{r})$ from wannier90 via Monte Carlo sampling:

$$U = \int d\mathbf{r} d\mathbf{r'} \frac{W_1^2(\mathbf{r}) W_1^2(\mathbf{r'})}{|\mathbf{r} - \mathbf{r'}|}$$ $$V = \int d\mathbf{r} d\mathbf{r'} \frac{W_1^2(\mathbf{r}) W_2^2(\mathbf{r'})}{|\mathbf{r} - \mathbf{r'}|}$$ $$J = \int d\mathbf{r} d\mathbf{r'} \frac{W_1(\mathbf{r}) W_2(\mathbf{r}) W_1(\mathbf{r'}) W_2(\mathbf{r'}) }{|\mathbf{r} - \mathbf{r'}|}$$

Dependencies

python version requires numba
c++ version can be complied with OpenMP via CMake:
mkdir build
cd build/
cmake ..
make

Usage

Run python Wannier_Coulomb.py or Wannier_Coulomb.x in the same directory with W1.xsf and W2.xsf files from wannier90. Play with mc_steps parameter to reach the required accuracy. Optionally set r_center and r_cut to reduce the size for Monte Carlo sampling, but keeping the norms of Wannier functions close to 1.