Temperature and metallic systems

In this example we consider the modeling of a magnesium lattice as a simple example for a metallic system. For our treatment we will use the PBE exchange-correlation functional. First we import required packages and setup the lattice. Again notice that DFTK uses the convention that lattice vectors are specified column by column.

In [1]:
using DFTK
using Plots

a = 3.01794  # bohr
b = 5.22722  # bohr
c = 9.77362  # bohr
lattice = [[-a -a  0]; [-b  b  0]; [0   0 -c]]
Mg = ElementPsp(:Mg, psp=load_psp("hgh/pbe/Mg-q2"))
atoms = [Mg => [[2/3, 1/3, 1/4], [1/3, 2/3, 3/4]]];
┌ Warning: Package DFTK does not have Plots in its dependencies:
│ - If you have DFTK checked out for development and have
│   added Plots as a dependency but haven't updated your primary
│   environment's manifest file, try `Pkg.resolve()`.
│ - Otherwise you may need to report an issue with DFTK
│ Loading Plots into DFTK from project dependency, future warnings for DFTK are suppressed.
└ @ nothing nothing:909

Next we build the PBE model and discretize it. Since magnesium is a metal we apply a small smearing temperature to ease convergence using the Fermi-Dirac smearing scheme. Note that both the Ecut is too small as well as the minimal k-point spacing kspacing far too large to give a converged result. These have been selected to obtain a fast execution time. By default PlaneWaveBasis chooses a kspacing of 2π * 0.022 inverse Bohrs, which is much more reasonable.

In [2]:
kspacing = 0.5      # Minimal spacing of k-points,
#                    in units of wavevectors (inverse Bohrs)
Ecut = 5            # kinetic energy cutoff in Hartree
temperature = 0.01  # Smearing temperature in Hartree

model = model_DFT(lattice, atoms, [:gga_x_pbe, :gga_c_pbe];
                  temperature=temperature,
                  smearing=DFTK.Smearing.FermiDirac())
kgrid = kgrid_size_from_minimal_spacing(lattice, kspacing)
basis = PlaneWaveBasis(model, Ecut, kgrid=kgrid);

Finally we run the SCF. Two magnesium atoms in our pseudopotential model result in four valence electrons being explicitly treated. Nevertheless this SCF will solve for eight bands by default in order to capture partial occupations beyond the Fermi level due to the employed smearing scheme.

In [3]:
scfres = self_consistent_field(basis);
n     Free energy       Eₙ-Eₙ₋₁     ρout-ρin   Diag
---   ---------------   ---------   --------   ----
  1   -1.761876675688         NaN   5.05e-02    4.3 
  2   -1.762189266484   -3.13e-04   9.85e-03    2.7 
  3   -1.762240544347   -5.13e-05   6.95e-04    2.0 
  4   -1.762241608044   -1.06e-06   8.64e-05    3.7 
  5   -1.762241620417   -1.24e-08   5.50e-06    2.7 
In [4]:
scfres.occupation[1]
Out[4]:
9-element Array{Float64,1}:
 1.9999999999077946
 1.9999975862873673
 0.004016894311510368
 3.000471271318729e-15
 1.1116102630424774e-18
 1.1115354537046484e-18
 7.977986029335707e-19
 7.977289382637639e-19
 3.2639285292055033e-22
In [5]:
scfres.energies
Out[5]:
Energy breakdown:
    Kinetic             0.7180638 
    AtomicLocal         0.3145369 
    AtomicNonlocal      0.3265786 
    Ewald               -2.1544222
    PspCorrection       -0.1026056
    Hartree             0.0055002 
    Xc                  -0.8610488
    Entropy             -0.0088446

    total               -1.762241620417

The fact that magnesium is a metal is confirmed by plotting the density of states around the Fermi level.

In [6]:
plot_dos(scfres)
Out[6]: