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.916320546596                   -0.69    4.2    659ms
  2   -7.921166431814       -2.31       -1.52    1.0    274ms
  3   -7.921408085869       -3.62       -2.48    1.5    303ms
  4   -7.921440930301       -4.48       -2.84    2.3    413ms
  5   -7.921441778633       -6.07       -3.18    1.0    359ms
  6   -7.921442008318       -6.64       -4.42    1.1    1.03s
  7   -7.921442022053       -7.86       -5.17    2.4    428ms
  8   -7.921442022139      -10.06       -5.52    2.1    367ms
  9   -7.921442022143      -11.35       -6.04    1.2    294ms
 10   -7.921442022144      -12.27       -6.18    1.3    308ms
 11   -7.921442022144      -12.74       -7.14    1.0    270ms
 12   -7.921442022144      -13.77       -7.94    2.1    416ms
 13   -7.921442022144   +  -15.05       -8.58    2.0    382ms

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.83s  Non-interacting
   1        0       1        1.03     15.2   13.0s  
   2        0       2        0.53     14.6   1.22s  
   3        0       3       -0.01     14.1   1.16s  
   4        0       4       -0.87     13.7   1.12s  
   5        0       5       -1.74     11.8   1.00s  
   6        0       6       -2.43     10.0   877ms  
   7        0       7       -3.33      8.9   799ms  
   8        0       8       -4.02      7.8   739ms  
   9        0       9       -5.04      6.4   636ms  
  10        0      10       -6.06      3.9   476ms  
  11        0      11       -7.08      2.3   374ms  
  12        0      12       -8.28      1.0   294ms  
  13        1       1       -7.53     17.8   1.78s  Restart
  14        1       2       -7.89      3.1   393ms  
  15        1       3       -8.71      2.1   340ms  
                                      19.2   2.10s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   2.22s  Non-interacting
   1        0       1        1.03     15.2   1.25s  
   2        0       2        0.53     14.6   1.14s  
   3        0       3       -0.01     14.1   1.13s  
   4        0       4       -0.87     13.8   1.12s  
   5        0       5       -1.74     11.8   1.00s  
   6        0       6       -2.43     10.0   885ms  
   7        0       7       -3.33      8.9   797ms  
   8        0       8       -4.02      7.8   746ms  
   9        0       9       -5.04      6.4   638ms  
  10        0      10       -6.06      3.9   468ms  
  11        0      11       -7.08      2.3   376ms  
  12        0      12       -8.27      1.0   296ms  
  13        1       1       -7.53     17.8   1.60s  Restart
  14        1       2       -7.89      3.1   383ms  
  15        1       3       -8.71      2.1   327ms  
                                      19.2   1.49s  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.6   1.16s  
   3        0       3       -0.01     14.1   1.14s  
   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   882ms  
   7        0       7       -3.33      8.9   801ms  
   8        0       8       -4.02      7.8   737ms  
   9        0       9       -5.04      6.5   634ms  
  10        0      10       -6.06      3.9   470ms  
  11        0      11       -7.08      2.3   372ms  
  12        0      12       -8.27      1.0   294ms  
  13        1       1       -7.51     17.9   1.60s  Restart
  14        1       2       -7.86      3.2   377ms  
  15        1       3       -8.66      2.1   328ms  
                                      19.2   2.32s  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.6   1.18s  
   3        0       3       -0.01     14.1   1.13s  
   4        0       4       -0.87     13.7   1.12s  
   5        0       5       -1.74     11.8   1.01s  
   6        0       6       -2.43     10.0   877ms  
   7        0       7       -3.33      8.9   806ms  
   8        0       8       -4.02      7.8   737ms  
   9        0       9       -5.04      6.4   636ms  
  10        0      10       -6.06      3.9   473ms  
  11        0      11       -7.08      2.3   375ms  
  12        0      12       -8.27      1.0   294ms  
  13        1       1       -7.52     17.8   1.61s  Restart
  14        1       2       -7.88      3.1   374ms  
  15        1       3       -8.71      2.1   326ms  
                                      19.2   1.48s  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.6   1.16s  
   3        0       3       -0.01     14.1   1.13s  
   4        0       4       -0.87     13.8   1.13s  
   5        0       5       -1.74     11.8   995ms  
   6        0       6       -2.43     10.0   881ms  
   7        0       7       -3.33      8.9   795ms  
   8        0       8       -4.02      7.8   736ms  
   9        0       9       -5.04      6.4   636ms  
  10        0      10       -6.06      3.9   467ms  
  11        0      11       -7.08      2.3   371ms  
  12        0      12       -8.27      1.0   293ms  
  13        1       1       -7.53     17.8   1.58s  Restart
  14        1       2       -7.90      3.1   380ms  
  15        1       3       -8.71      2.1   325ms  
                                      19.2   2.33s  Final orbitals
Iter  Restart  Krydim  log10(res)  avg(CG)  Δtime   Comment
----  -------  ------  ----------  -------  ------  ---------------
                                      20.0   1.31s  Non-interacting
   1        0       1        1.03     15.2   1.25s  
   2        0       2        0.53     14.6   1.15s  
   3        0       3       -0.01     14.1   1.14s  
   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   877ms  
   7        0       7       -3.33      8.9   807ms  
   8        0       8       -4.02      7.8   734ms  
   9        0       9       -5.04      6.5   636ms  
  10        0      10       -6.06      3.9   463ms  
  11        0      11       -7.08      2.3   377ms  
  12        0      12       -8.27      1.0   295ms  
  13        1       1       -7.52     17.9   1.61s  Restart
  14        1       2       -7.87      3.1   375ms  
  15        1       3       -8.67      2.1   327ms  
                                      19.2   1.48s  Final orbitals
133.530743 seconds (102.04 M allocations: 30.822 GiB, 3.82% gc time, 26.70% 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