Phonon computations

This is a quick sketch how to run a simple phonon calculation using DFTK.

Preliminary implementation

Practical phonon computations have only seen rudimentary testing as of now. As of now we do not yet recommend relying on this feature for production calculations. Some of the limitations are:

  • symmetries must be disabled (pass symmetries=false to the model)
  • only LDA functionals are supported
  • non-linear core corrections from the pseudopotentials are not supported
  • MPI parallelization over k-points is not supported (due to $k$ and $k+q$ interactions)

We appreciate any issues, bug reports or PRs.

First we run an SCF calculation.

using AtomsBuilder
using DFTK
using Printf
using PseudoPotentialData

pseudopotentials = PseudoFamily("cp2k.nc.sr.lda.v0_1.semicore.gth")
# Make sure to disable symmetries:
model  = model_DFT(bulk(:Si); pseudopotentials, functionals=LDA(), symmetries=false)
basis  = PlaneWaveBasis(model; Ecut=10, kgrid=[4, 4, 4])
scfres = self_consistent_field(basis, tol=1e-8);
nothing  # hide
n     Energy            log10(ΔE)   log10(Δρ)   Diag   Δtime 
---   ---------------   ---------   ---------   ----   ------
  1   -7.916327480915                   -0.69    4.2    693ms
  2   -7.921174353403       -2.31       -1.52    1.0    292ms
  3   -7.921408184453       -3.63       -2.48    1.5    319ms
  4   -7.921440926448       -4.48       -2.85    2.2    425ms
  5   -7.921441774803       -6.07       -3.17    1.0    284ms
  6   -7.921442007903       -6.63       -4.49    1.1    295ms
  7   -7.921442022051       -7.85       -4.98    2.5    471ms
  8   -7.921442022134      -10.08       -5.28    1.6    340ms
  9   -7.921442022143      -11.03       -6.07    1.1    292ms
 10   -7.921442022144      -12.00       -6.29    2.0    476ms
 11   -7.921442022144      -12.90       -7.31    1.0    1.19s
 12   -7.921442022144      -14.15       -7.79    2.2    433ms
 13   -7.921442022144      -14.57       -8.61    1.2    287ms

Next we compute the phonon modes at the q-point [1/4, 1/4, 1/4].

phret_q0 = @time DFTK.phonon_modes(scfres; q=[0.25, 0.25, 0.25]);
nothing  # hide
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.88s  Non-interacting
   1        0       1        1.03     15.1   13.4s  
   2        0       2        0.53     14.5   1.26s  
   3        0       3       -0.01     14.1   1.19s  
   4        0       4       -0.87     13.8   1.17s  
   5        0       5       -1.74     11.8   1.05s  
   6        0       6       -2.43     10.0   916ms  
   7        0       7       -3.33      8.9   834ms  
   8        0       8       -4.02      7.8   768ms  
   9        0       9       -5.04      6.4   660ms  
  10        0      10       -6.06      3.9   491ms  
  11        0      11       -7.08      2.3   387ms  
  12        0      12       -8.28      1.0   303ms  
  13        1       1       -7.55     17.8   1.84s  Restart
  14        1       2       -7.89      3.1   408ms  
  15        1       3       -8.71      2.1   351ms  
                                      19.2   2.18s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.44s  Non-interacting
   1        0       1        1.03     15.1   1.32s  
   2        0       2        0.53     14.5   1.21s  
   3        0       3       -0.01     14.0   1.19s  
   4        0       4       -0.87     13.7   1.17s  
   5        0       5       -1.74     11.8   1.04s  
   6        0       6       -2.43     10.0   907ms  
   7        0       7       -3.33      8.9   927ms  
   8        0       8       -4.02      7.8   1.58s  
   9        0       9       -5.04      6.4   635ms  
  10        0      10       -6.06      3.9   457ms  
  11        0      11       -7.08      2.3   365ms  
  12        0      12       -8.27      1.0   305ms  
  13        1       1       -7.53     17.8   1.63s  Restart
  14        1       2       -7.89      3.1   390ms  
  15        1       3       -8.69      2.1   325ms  
                                      19.2   1.57s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.44s  Non-interacting
   1        0       1        1.03     15.1   1.33s  
   2        0       2        0.53     14.6   1.22s  
   3        0       3       -0.01     14.1   1.19s  
   4        0       4       -0.87     13.8   1.18s  
   5        0       5       -1.74     11.8   1.05s  
   6        0       6       -2.43     10.0   911ms  
   7        0       7       -3.33      8.9   833ms  
   8        0       8       -4.02      7.8   769ms  
   9        0       9       -5.04      6.4   660ms  
  10        0      10       -6.06      3.9   483ms  
  11        0      11       -7.08      2.3   387ms  
  12        0      12       -8.27      1.0   304ms  
  13        1       1       -7.53     17.8   1.66s  Restart
  14        1       2       -7.89      3.1   385ms  
  15        1       3       -8.69      2.1   344ms  
                                      19.2   1.56s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.42s  Non-interacting
   1        0       1        1.03     15.1   1.32s  
   2        0       2        0.53     14.5   1.21s  
   3        0       3       -0.01     14.1   1.18s  
   4        0       4       -0.87     13.8   1.17s  
   5        0       5       -1.74     11.8   1.04s  
   6        0       6       -2.43     10.0   916ms  
   7        0       7       -3.33      8.9   825ms  
   8        0       8       -4.02      7.8   861ms  
   9        0       9       -5.04      6.4   1.47s  
  10        0      10       -6.06      3.9   482ms  
  11        0      11       -7.08      2.3   369ms  
  12        0      12       -8.28      1.0   286ms  
  13        1       1       -7.54     17.9   1.60s  Restart
  14        1       2       -7.87      3.1   380ms  
  15        1       3       -8.70      2.1   327ms  
                                      19.2   1.53s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.43s  Non-interacting
   1        0       1        1.03     15.1   1.33s  
   2        0       2        0.53     14.5   1.21s  
   3        0       3       -0.01     14.0   1.19s  
   4        0       4       -0.87     13.7   1.18s  
   5        0       5       -1.74     11.8   1.04s  
   6        0       6       -2.43     10.0   907ms  
   7        0       7       -3.33      8.9   826ms  
   8        0       8       -4.02      7.8   764ms  
   9        0       9       -5.04      6.4   658ms  
  10        0      10       -6.06      3.9   478ms  
  11        0      11       -7.08      2.3   393ms  
  12        0      12       -8.27      1.0   304ms  
  13        1       1       -7.53     17.8   1.66s  Restart
  14        1       2       -7.89      3.1   392ms  
  15        1       3       -8.70      2.1   337ms  
                                      19.2   1.57s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.44s  Non-interacting
   1        0       1        1.03     15.1   1.32s  
   2        0       2        0.53     14.5   1.22s  
   3        0       3       -0.01     14.1   1.19s  
   4        0       4       -0.87     13.8   1.17s  
   5        0       5       -1.74     11.8   1.04s  
   6        0       6       -2.43     10.0   917ms  
   7        0       7       -3.33      8.9   831ms  
   8        0       8       -4.02      7.8   765ms  
   9        0       9       -5.04      6.4   658ms  
  10        0      10       -6.06      3.9   482ms  
  11        0      11       -7.08      2.3   393ms  
  12        0      12       -8.27      1.0   299ms  
  13        1       1       -7.53     17.8   2.51s  Restart
  14        1       2       -7.88      3.1   369ms  
  15        1       3       -8.70      2.1   319ms  
                                      19.2   1.56s  Final orbitals
138.108581 seconds (102.04 M allocations: 30.809 GiB, 3.87% gc time, 26.32% compilation time)

These are the final phonon frequencies:

for (i, ω) in enumerate(phret_q0.frequencies)
    @printf("Mode %2d: %8.3f cm-1\n", i, ω .* DFTK.hartree_to_cm⁻¹)
end
Mode  1:   93.616 cm-1
Mode  2:   93.616 cm-1
Mode  3:  223.147 cm-1
Mode  4:  478.527 cm-1
Mode  5:  484.617 cm-1
Mode  6:  484.617 cm-1