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.916345639962                   -0.69    4.2    681ms
  2   -7.921179253813       -2.32       -1.52    1.0    279ms
  3   -7.921407946927       -3.64       -2.48    1.4    329ms
  4   -7.921440887483       -4.48       -2.84    2.3    433ms
  5   -7.921441758865       -6.06       -3.16    1.0    318ms
  6   -7.921442007420       -6.60       -4.50    1.1    290ms
  7   -7.921442022075       -7.83       -5.01    2.6    474ms
  8   -7.921442022136      -10.22       -5.31    1.7    341ms
  9   -7.921442022142      -11.22       -5.80    1.1    298ms
 10   -7.921442022144      -11.82       -6.18    1.2    307ms
 11   -7.921442022144      -12.66       -7.28    1.3    319ms
 12   -7.921442022144      -13.88       -7.33    2.2    453ms
 13   -7.921442022144      -14.75       -8.09    1.0    291ms

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.92s  Non-interacting
   1        0       1        1.03     15.2   13.1s  
   2        0       2        0.53     14.5   2.11s  
   3        0       3       -0.01     14.1   1.12s  
   4        0       4       -0.87     13.7   1.17s  
   5        0       5       -1.74     11.8   1.03s  
   6        0       6       -2.43     10.0   890ms  
   7        0       7       -3.33      8.9   845ms  
   8        0       8       -4.02      7.8   764ms  
   9        0       9       -5.04      6.4   666ms  
  10        0      10       -6.06      4.0   497ms  
  11        0      11       -7.08      2.3   390ms  
  12        0      12       -8.27      1.0   320ms  
  13        1       1       -7.52     17.8   1.86s  Restart
  14        1       2       -7.89      3.1   406ms  
  15        1       3       -8.72      2.1   324ms  
                                      19.2   2.16s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.43s  Non-interacting
   1        0       1        1.03     15.2   1.36s  
   2        0       2        0.53     14.5   1.22s  
   3        0       3       -0.01     14.1   1.20s  
   4        0       4       -0.87     13.7   1.17s  
   5        0       5       -1.74     11.8   1.06s  
   6        0       6       -2.43     10.0   914ms  
   7        0       7       -3.33      8.9   856ms  
   8        0       8       -4.02      7.8   778ms  
   9        0       9       -5.04      6.4   671ms  
  10        0      10       -6.06      3.9   496ms  
  11        0      11       -7.08      2.2   393ms  
  12        0      12       -8.27      1.0   310ms  
  13        1       1       -7.54     17.7   1.67s  Restart
  14        1       2       -7.90      3.1   396ms  
  15        1       3       -8.70      2.0   345ms  
                                      19.2   1.54s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   2.31s  Non-interacting
   1        0       1        1.03     15.2   1.27s  
   2        0       2        0.53     14.6   1.20s  
   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   923ms  
   7        0       7       -3.33      8.9   837ms  
   8        0       8       -4.02      7.8   759ms  
   9        0       9       -5.04      6.4   663ms  
  10        0      10       -6.06      3.9   507ms  
  11        0      11       -7.08      2.3   396ms  
  12        0      12       -8.27      1.0   312ms  
  13        1       1       -7.52     17.8   1.68s  Restart
  14        1       2       -7.89      3.1   398ms  
  15        1       3       -8.69      2.1   341ms  
                                      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.22s  
   3        0       3       -0.01     14.1   1.18s  
   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   912ms  
   7        0       7       -3.33      8.9   828ms  
   8        0       8       -4.02      7.8   773ms  
   9        0       9       -5.04      6.4   662ms  
  10        0      10       -6.06      4.0   492ms  
  11        0      11       -7.08      2.3   390ms  
  12        0      12       -8.27      1.0   309ms  
  13        1       1       -7.52     17.8   1.66s  Restart
  14        1       2       -7.89      3.1   393ms  
  15        1       3       -8.72      2.1   347ms  
                                      19.2   1.63s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   2.14s  Non-interacting
   1        0       1        1.03     15.2   1.29s  
   2        0       2        0.53     14.5   1.19s  
   3        0       3       -0.01     14.1   1.16s  
   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   911ms  
   7        0       7       -3.33      8.9   831ms  
   8        0       8       -4.02      7.8   774ms  
   9        0       9       -5.04      6.4   665ms  
  10        0      10       -6.06      3.9   486ms  
  11        0      11       -7.08      2.2   386ms  
  12        0      12       -8.27      1.0   316ms  
  13        1       1       -7.54     17.7   1.66s  Restart
  14        1       2       -7.90      3.1   389ms  
  15        1       3       -8.70      2.1   339ms  
                                      19.2   1.55s  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.33s  
   2        0       2        0.53     14.6   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.05s  
   6        0       6       -2.43     10.0   913ms  
   7        0       7       -3.33      8.9   828ms  
   8        0       8       -4.02      7.8   767ms  
   9        0       9       -5.04      6.4   656ms  
  10        0      10       -6.06      3.9   491ms  
  11        0      11       -7.08      2.3   383ms  
  12        0      12       -8.27      1.0   309ms  
  13        1       1       -7.53     17.8   2.50s  Restart
  14        1       2       -7.89      3.1   364ms  
  15        1       3       -8.69      2.1   309ms  
                                      19.2   1.47s  Final orbitals
138.339944 seconds (102.85 M allocations: 30.881 GiB, 4.33% gc time, 26.12% 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