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.916332322535                   -0.69    4.1    640ms
  2   -7.921175427663       -2.31       -1.52    1.0    297ms
  3   -7.921406933331       -3.64       -2.48    1.4    332ms
  4   -7.921440898817       -4.47       -2.84    2.3    441ms
  5   -7.921441757977       -6.07       -3.15    1.0    296ms
  6   -7.921442006333       -6.60       -4.44    1.0    299ms
  7   -7.921442022044       -7.80       -5.03    2.5    460ms
  8   -7.921442022137      -10.03       -5.38    1.8    343ms
  9   -7.921442022143      -11.18       -6.07    1.1    274ms
 10   -7.921442022144      -12.23       -6.05    1.9    354ms
 11   -7.921442022144      -12.60       -7.03    1.0    364ms
 12   -7.921442022144      -13.80       -7.55    2.0    1.23s
 13   -7.921442022144   +    -Inf       -8.34    1.2    304ms

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.99s  Non-interacting
   1        0       1        1.03     15.2   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.7   1.15s  
   5        0       5       -1.74     11.8   1.05s  
   6        0       6       -2.43     10.0   929ms  
   7        0       7       -3.33      8.8   844ms  
   8        0       8       -4.02      7.8   783ms  
   9        0       9       -5.04      6.4   657ms  
  10        0      10       -6.06      3.9   490ms  
  11        0      11       -7.08      2.3   390ms  
  12        0      12       -8.28      1.0   319ms  
  13        1       1       -7.52     17.9   1.87s  Restart
  14        1       2       -7.87      3.2   425ms  
  15        1       3       -8.68      2.1   366ms  
                                      19.2   2.18s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.41s  Non-interacting
   1        0       1        1.03     15.2   1.35s  
   2        0       2        0.53     14.5   1.23s  
   3        0       3       -0.01     14.1   1.21s  
   4        0       4       -0.87     13.7   1.17s  
   5        0       5       -1.74     11.8   1.99s  
   6        0       6       -2.43     10.0   893ms  
   7        0       7       -3.33      8.9   810ms  
   8        0       8       -4.02      7.8   718ms  
   9        0       9       -5.04      6.4   657ms  
  10        0      10       -6.06      3.9   518ms  
  11        0      11       -7.08      2.2   383ms  
  12        0      12       -8.27      1.0   318ms  
  13        1       1       -7.54     17.9   1.63s  Restart
  14        1       2       -7.88      3.2   394ms  
  15        1       3       -8.66      2.1   346ms  
                                      19.2   1.51s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.40s  Non-interacting
   1        0       1        1.03     15.2   1.31s  
   2        0       2        0.53     14.5   1.21s  
   3        0       3       -0.01     14.1   1.17s  
   4        0       4       -0.87     13.8   1.16s  
   5        0       5       -1.74     11.8   1.04s  
   6        0       6       -2.43     10.0   917ms  
   7        0       7       -3.33      8.9   830ms  
   8        0       8       -4.02      7.8   760ms  
   9        0       9       -5.04      6.4   648ms  
  10        0      10       -6.06      3.9   494ms  
  11        0      11       -7.08      2.3   395ms  
  12        0      12       -8.27      1.0   310ms  
  13        1       1       -7.52     17.9   1.65s  Restart
  14        1       2       -7.87      3.2   412ms  
  15        1       3       -8.68      2.1   343ms  
                                      19.2   1.52s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.37s  Non-interacting
   1        0       1        1.03     15.2   1.34s  
   2        0       2        0.53     14.5   1.21s  
   3        0       3       -0.01     14.1   1.20s  
   4        0       4       -0.87     13.7   1.18s  
   5        0       5       -1.74     11.8   1.05s  
   6        0       6       -2.43     10.0   914ms  
   7        0       7       -3.33      8.9   837ms  
   8        0       8       -4.02      7.8   776ms  
   9        0       9       -5.04      6.4   1.63s  
  10        0      10       -6.06      3.9   464ms  
  11        0      11       -7.08      2.3   373ms  
  12        0      12       -8.27      1.0   283ms  
  13        1       1       -7.53     17.8   1.63s  Restart
  14        1       2       -7.88      3.1   374ms  
  15        1       3       -8.69      2.1   348ms  
                                      19.2   1.51s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.41s  Non-interacting
   1        0       1        1.03     15.2   1.34s  
   2        0       2        0.53     14.5   1.20s  
   3        0       3       -0.01     14.1   1.15s  
   4        0       4       -0.87     13.7   1.16s  
   5        0       5       -1.74     11.8   1.05s  
   6        0       6       -2.43     10.0   928ms  
   7        0       7       -3.33      8.9   825ms  
   8        0       8       -4.02      7.8   749ms  
   9        0       9       -5.04      6.4   656ms  
  10        0      10       -6.06      3.9   499ms  
  11        0      11       -7.08      2.2   398ms  
  12        0      12       -8.27      1.0   314ms  
  13        1       1       -7.54     17.9   1.68s  Restart
  14        1       2       -7.88      3.2   403ms  
  15        1       3       -8.66      2.1   344ms  
                                      19.2   1.54s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.42s  Non-interacting
   1        0       1        1.03     15.2   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.16s  
   5        0       5       -1.74     11.8   1.04s  
   6        0       6       -2.43     10.0   911ms  
   7        0       7       -3.33      8.9   827ms  
   8        0       8       -4.02      7.8   761ms  
   9        0       9       -5.04      6.5   660ms  
  10        0      10       -6.06      3.9   486ms  
  11        0      11       -7.08      2.3   395ms  
  12        0      12       -8.27      1.0   317ms  
  13        1       1       -7.52     17.9   2.57s  Restart
  14        1       2       -7.87      3.2   372ms  
  15        1       3       -8.68      2.1   313ms  
                                      19.2   1.49s  Final orbitals
138.950646 seconds (102.38 M allocations: 30.840 GiB, 3.86% gc time, 26.71% 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