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.916320812030                   -0.69    4.1    655ms
  2   -7.921171547687       -2.31       -1.52    1.0    270ms
  3   -7.921408023880       -3.63       -2.48    1.4    292ms
  4   -7.921440882261       -4.48       -2.85    2.3    404ms
  5   -7.921441769329       -6.05       -3.17    1.1    371ms
  6   -7.921442007274       -6.62       -4.55    1.0    1.17s
  7   -7.921442022091       -7.83       -5.03    2.6    451ms
  8   -7.921442022138      -10.33       -5.32    1.7    334ms
  9   -7.921442022142      -11.45       -5.64    1.1    279ms
 10   -7.921442022144      -11.71       -6.83    1.0    267ms
 11   -7.921442022144      -12.90       -7.30    2.4    421ms
 12   -7.921442022144      -14.45       -7.50    1.3    291ms
 13   -7.921442022144   +    -Inf       -8.42    1.0    280ms

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.86s  Non-interacting
   1        0       1        1.03     15.2   13.3s  
   2        0       2        0.53     14.5   1.23s  
   3        0       3       -0.01     14.1   1.15s  
   4        0       4       -0.87     13.7   1.12s  
   5        0       5       -1.74     11.8   977ms  
   6        0       6       -2.43     10.0   853ms  
   7        0       7       -3.33      8.9   776ms  
   8        0       8       -4.02      7.8   737ms  
   9        0       9       -5.04      6.4   618ms  
  10        0      10       -6.06      3.9   462ms  
  11        0      11       -7.08      2.2   364ms  
  12        0      12       -8.27      1.0   290ms  
  13        1       1       -7.52     17.8   1.76s  Restart
  14        1       2       -7.87      3.2   394ms  
  15        1       3       -8.71      2.1   338ms  
                                      19.2   2.05s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   2.32s  Non-interacting
   1        0       1        1.03     15.2   1.22s  
   2        0       2        0.53     14.6   1.13s  
   3        0       3       -0.01     14.1   1.11s  
   4        0       4       -0.87     13.8   1.09s  
   5        0       5       -1.74     11.8   978ms  
   6        0       6       -2.43     10.0   863ms  
   7        0       7       -3.33      8.9   792ms  
   8        0       8       -4.02      7.8   725ms  
   9        0       9       -5.04      6.5   634ms  
  10        0      10       -6.06      3.9   470ms  
  11        0      11       -7.08      2.3   366ms  
  12        0      12       -8.28      1.0   293ms  
  13        1       1       -7.54     17.8   1.59s  Restart
  14        1       2       -7.90      3.1   368ms  
  15        1       3       -8.69      2.0   329ms  
                                      19.2   1.46s  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.27s  
   2        0       2        0.53     14.5   1.16s  
   3        0       3       -0.01     14.1   1.12s  
   4        0       4       -0.87     13.8   1.12s  
   5        0       5       -1.74     11.8   982ms  
   6        0       6       -2.43     10.0   856ms  
   7        0       7       -3.33      8.9   818ms  
   8        0       8       -4.02      7.8   726ms  
   9        0       9       -5.04      6.4   640ms  
  10        0      10       -6.06      3.9   456ms  
  11        0      11       -7.08      2.2   383ms  
  12        0      12       -8.27      1.0   290ms  
  13        1       1       -7.52     17.9   1.58s  Restart
  14        1       2       -7.87      3.2   380ms  
  15        1       3       -8.66      2.1   326ms  
                                      19.2   2.48s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.32s  Non-interacting
   1        0       1        1.03     15.2   1.25s  
   2        0       2        0.53     14.5   1.14s  
   3        0       3       -0.01     14.1   1.11s  
   4        0       4       -0.87     13.7   1.10s  
   5        0       5       -1.74     11.8   996ms  
   6        0       6       -2.43     10.0   863ms  
   7        0       7       -3.33      8.9   787ms  
   8        0       8       -4.02      7.8   758ms  
   9        0       9       -5.04      6.4   627ms  
  10        0      10       -6.06      3.9   458ms  
  11        0      11       -7.08      2.2   365ms  
  12        0      12       -8.27      1.0   300ms  
  13        1       1       -7.51     17.9   1.61s  Restart
  14        1       2       -7.87      3.2   376ms  
  15        1       3       -8.69      2.2   325ms  
                                      19.2   1.48s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.36s  Non-interacting
   1        0       1        1.03     15.2   1.27s  
   2        0       2        0.53     14.6   1.16s  
   3        0       3       -0.01     14.1   1.14s  
   4        0       4       -0.87     13.8   1.11s  
   5        0       5       -1.74     11.8   992ms  
   6        0       6       -2.43     10.0   870ms  
   7        0       7       -3.33      8.9   791ms  
   8        0       8       -4.02      7.8   727ms  
   9        0       9       -5.04      6.5   640ms  
  10        0      10       -6.06      3.9   463ms  
  11        0      11       -7.08      2.3   370ms  
  12        0      12       -8.28      1.0   291ms  
  13        1       1       -7.54     17.8   1.58s  Restart
  14        1       2       -7.90      3.1   385ms  
  15        1       3       -8.68      2.0   323ms  
                                      19.2   2.53s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.32s  Non-interacting
   1        0       1        1.03     15.2   1.27s  
   2        0       2        0.53     14.5   1.13s  
   3        0       3       -0.01     14.1   1.12s  
   4        0       4       -0.87     13.8   1.13s  
   5        0       5       -1.74     11.8   1.00s  
   6        0       6       -2.43     10.0   865ms  
   7        0       7       -3.33      8.9   803ms  
   8        0       8       -4.02      7.8   726ms  
   9        0       9       -5.04      6.4   637ms  
  10        0      10       -6.06      3.9   472ms  
  11        0      11       -7.08      2.2   374ms  
  12        0      12       -8.27      1.0   295ms  
  13        1       1       -7.51     17.9   1.60s  Restart
  14        1       2       -7.86      3.2   382ms  
  15        1       3       -8.65      2.1   322ms  
                                      19.2   1.47s  Final orbitals
135.588719 seconds (102.15 M allocations: 30.820 GiB, 4.64% gc time, 27.76% 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