Preconditioning #98
Preconditioning #98
Conversation
|
Can you add a test? Instructions here: https://github.com/wannier-developers/wannier90/blob/develop/test-suite/README You can just copy an existing test (e.g. one from the examples), give it a reasonable name, add the precond=true flag, and adapt the reference output. Better if you can have two tests for both types of optimisation - try however to choose a short test (<10secs ideally, doesn't need to be converged) and also that does not consume too much memory, as the testing machines might not have many GBs... |
|
Done. I just used gaas with a bad input guess, to make it a bit challenging to converge. Preconditioning does not do much in this case (2x2x2 grid), but it doesn't change the convergence very much. I put in two tests: one for the GEMM based implementation, one for the loop: results are identical to machine precision. |
This introduces a new keyword, precond, that turns on preconditioning of the SD/CG algorithm to optimize the spread. It is based on the same idea as the Tetter-Payne-Allan preconditioner of DFT: the functional is nearly diagonal in R-space, and the gradients put too much emphasis on the large-R contributions (plotting the gradients as a function of k shows a characteristic checkerboard pattern). Using this preconditioner speeds up the iteration for fine k-point grids. In my tests I got about x2 on 8x8x8 Silicon and x3 on 14x14 Arsenene, (exact speedups are sensitive to other input parameters such as tolerance, trial step, CG parameters...). I have not observed any slowdown for coarser grids.
In terms of implementation, I use a slow Fourier transform: I really should be using a FFT, but this requires a lot more code than this barebone implementation, that reuses the machinery in hamiltonian.F90 to construct a R-space grid (thanks to G. Pizzi for the idea!). On the systems I tried, the additional cost of preconditioning is nearly zero. It does require more RAM (sizeof(complex double) times the number of kpoints squared, independently of the number of Wannier functions), although modestly so: at 10k kpoints (22 x 22 x 22), it's 1.6 GB. When the "optimisation" flag is < 3, I turn off the GEMM-based implementation to save RAM: this is slow, but still yields a speedup it because of the gain in iterations.
There are a number of tunable parameters in the algorithm, which I do not expose: the default values seem good enough not to bother.
Any comments welcome!