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.916337246028                   -0.69    4.2    632ms
  2   -7.921175051279       -2.32       -1.52    1.0    287ms
  3   -7.921406983061       -3.63       -2.48    1.4    296ms
  4   -7.921440848661       -4.47       -2.83    2.2    411ms
  5   -7.921441751755       -6.04       -3.15    1.0    270ms
  6   -7.921442006970       -6.59       -4.45    1.1    272ms
  7   -7.921442022030       -7.82       -4.94    2.6    443ms
  8   -7.921442022133       -9.98       -5.25    1.7    316ms
  9   -7.921442022143      -11.03       -5.93    1.0    266ms
 10   -7.921442022144      -11.87       -6.10    1.8    347ms
 11   -7.921442022144      -12.58       -7.17    1.0    275ms
 12   -7.921442022144      -13.73       -7.37    2.1    410ms
 13   -7.921442022144      -15.05       -8.18    1.0    271ms

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.84s  Non-interacting
   1        0       1        1.03     15.1   14.5s  
   2        0       2        0.53     14.5   1.15s  
   3        0       3       -0.01     14.1   1.17s  
   4        0       4       -0.87     13.7   1.18s  
   5        0       5       -1.74     11.8   1.03s  
   6        0       6       -2.43     10.0   909ms  
   7        0       7       -3.33      8.8   823ms  
   8        0       8       -4.02      7.8   767ms  
   9        0       9       -5.04      6.4   651ms  
  10        0      10       -6.06      3.9   476ms  
  11        0      11       -7.08      2.3   386ms  
  12        0      12       -8.27      1.0   305ms  
  13        1       1       -7.54     17.8   1.84s  Restart
  14        1       2       -7.88      3.2   411ms  
  15        1       3       -8.71      2.1   352ms  
                                      19.2   2.15s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.41s  Non-interacting
   1        0       1        1.03     15.1   1.31s  
   2        0       2        0.53     14.5   1.22s  
   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.03s  
   6        0       6       -2.43     10.0   902ms  
   7        0       7       -3.33      8.9   820ms  
   8        0       8       -4.02      7.8   776ms  
   9        0       9       -5.04      6.5   653ms  
  10        0      10       -6.06      3.9   471ms  
  11        0      11       -7.08      2.3   383ms  
  12        0      12       -8.27      1.0   301ms  
  13        1       1       -7.54     17.9   1.64s  Restart
  14        1       2       -7.88      3.2   380ms  
  15        1       3       -8.68      2.1   329ms  
                                      19.2   1.52s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.40s  Non-interacting
   1        0       1        1.03     15.1   2.25s  
   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.15s  
   5        0       5       -1.74     11.8   1.01s  
   6        0       6       -2.43     10.0   898ms  
   7        0       7       -3.33      8.9   835ms  
   8        0       8       -4.02      7.8   761ms  
   9        0       9       -5.04      6.4   664ms  
  10        0      10       -6.06      3.9   495ms  
  11        0      11       -7.08      2.3   384ms  
  12        0      12       -8.27      1.0   308ms  
  13        1       1       -7.54     17.9   1.66s  Restart
  14        1       2       -7.88      3.2   385ms  
  15        1       3       -8.67      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.1   1.34s  
   2        0       2        0.53     14.5   1.20s  
   3        0       3       -0.01     14.1   1.17s  
   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   907ms  
   7        0       7       -3.33      8.9   822ms  
   8        0       8       -4.02      7.8   761ms  
   9        0       9       -5.04      6.4   668ms  
  10        0      10       -6.06      3.9   488ms  
  11        0      11       -7.08      2.3   384ms  
  12        0      12       -8.27      1.0   302ms  
  13        1       1       -7.54     17.8   1.64s  Restart
  14        1       2       -7.89      3.2   383ms  
  15        1       3       -8.70      2.1   332ms  
                                      19.2   1.50s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.39s  Non-interacting
   1        0       1        1.03     15.1   1.30s  
   2        0       2        0.53     14.5   1.19s  
   3        0       3       -0.01     14.1   1.27s  
   4        0       4       -0.87     13.8   2.00s  
   5        0       5       -1.74     11.8   999ms  
   6        0       6       -2.43     10.0   881ms  
   7        0       7       -3.33      8.9   786ms  
   8        0       8       -4.02      7.8   764ms  
   9        0       9       -5.04      6.4   643ms  
  10        0      10       -6.06      3.9   471ms  
  11        0      11       -7.08      2.3   380ms  
  12        0      12       -8.27      1.0   301ms  
  13        1       1       -7.53     17.9   1.64s  Restart
  14        1       2       -7.87      3.2   391ms  
  15        1       3       -8.66      2.1   339ms  
                                      19.2   1.53s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.41s  Non-interacting
   1        0       1        1.03     15.1   1.32s  
   2        0       2        0.53     14.5   1.20s  
   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   901ms  
   7        0       7       -3.33      8.9   834ms  
   8        0       8       -4.02      7.8   760ms  
   9        0       9       -5.04      6.4   657ms  
  10        0      10       -6.06      3.9   480ms  
  11        0      11       -7.08      2.3   385ms  
  12        0      12       -8.27      1.0   303ms  
  13        1       1       -7.53     17.9   1.66s  Restart
  14        1       2       -7.88      3.2   380ms  
  15        1       3       -8.67      2.1   336ms  
                                      19.2   1.53s  Final orbitals
138.978715 seconds (96.50 M allocations: 30.453 GiB, 4.78% gc time, 26.85% 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