Ab initio calculation of polarons: algorithms and benchmarks*

Calculations of polarons using density-functional theory and supercell approaches face two key challenges: (i) large computational cells are required to describe intermediate and large polarons, (ii) the formation energy and localization of the polaron wavefunction are sensitive to the exchange-correlation functional. We developed a new approach where the polaron is expressed as a superposition of Bloch states, and the calculation of wavefunctions and energies is cast into the solution of a nonlinear system involving Kohn-Sham energies, phonon frequencies, and el-ph matrix elements from density-functional perturbation theory [PRL 122, 246403 (2019)]. Here we report on further optimization, specifically improvements of the iterative eigensolver, parallelism, memory management, and increased modularity using EPW and Quantum ESPRESSO. We analyze the performance of the method in terms of Brillouin-zone sampling and wavefunction initialization and report benchmark on prototypical polaronic systems, from ionic to covalent insulators.