Calculate the Band Structure of CopperΒΆ
This section shows how we can use KSSOLV2.0 to compute the band structure of a FCC copper crystal using a small unit cell and several k-points.
We first perform a ground state SCF calculation to obtain the electron density and the corresponding Hamiltonian associated with the ground state of the copper crystal.
The atomic species and the unit cell are set up as follows
a1 = Atom('Cu');
atomlist = [a1];
C = [-3.411 0.000 3.411
0.000 3.411 3.411
-3.411 3.411 0.000];
The unit cell contains one copper atom at the corner of the cubic unit cell.
xyzlist = [0.0 0.0 0.0];
The cu_scf_k.m script in the Examples directory shows how one can either using 4x4x4 Monkhorst k-points or user specified k-points and weights to perform the SCF calculation. In the latter case, we set
kpts = [
0.0000 0.0000 0.0000
0.0000 0.2500 0.0000
0.0000 0.5000 0.0000
0.0000 0.7500 0.0000
0.2500 0.0000 0.0000
0.2500 0.2500 0.0000
0.2500 0.5000 0.0000
0.2500 0.7500 0.0000
0.5000 0.0000 0.0000
0.5000 0.2500 0.0000
0.5000 0.5000 0.0000
0.5000 0.7500 0.0000
0.7500 0.0000 0.0000
0.7500 0.2500 0.0000
0.7500 0.5000 0.0000
0.7500 0.7500 0.0000
0.0000 0.0000 0.2500
0.0000 0.2500 0.2500
0.0000 0.5000 0.2500
0.0000 0.7500 0.2500
0.2500 0.0000 0.2500
0.2500 0.2500 0.2500
0.2500 0.5000 0.2500
0.2500 0.7500 0.2500
0.5000 0.0000 0.2500
0.5000 0.2500 0.2500
0.5000 0.5000 0.2500
0.5000 0.7500 0.2500
0.7500 0.0000 0.2500
0.7500 0.2500 0.2500
0.7500 0.5000 0.2500
0.7500 0.7500 0.2500
0.0000 0.0000 0.5000
0.0000 0.2500 0.5000
0.0000 0.5000 0.5000
0.0000 0.7500 0.5000
0.2500 0.0000 0.5000
0.2500 0.2500 0.5000
0.2500 0.5000 0.5000
0.2500 0.7500 0.5000
0.5000 0.0000 0.5000
0.5000 0.2500 0.5000
0.5000 0.5000 0.5000
0.5000 0.7500 0.5000
0.7500 0.0000 0.5000
0.7500 0.2500 0.5000
0.7500 0.5000 0.5000
0.7500 0.7500 0.5000
0.0000 0.0000 0.7500
0.0000 0.2500 0.7500
0.0000 0.5000 0.7500
0.0000 0.7500 0.7500
0.2500 0.0000 0.7500
0.2500 0.2500 0.7500
0.2500 0.5000 0.7500
0.2500 0.7500 0.7500
0.5000 0.0000 0.7500
0.5000 0.2500 0.7500
0.5000 0.5000 0.7500
0.5000 0.7500 0.7500
0.7500 0.0000 0.7500
0.7500 0.2500 0.7500
0.7500 0.5000 0.7500
0.7500 0.7500 0.7500
];
wks = 0.0156250 * ones(size(kpts,1),1);
cry = Crystal('supercell',C, 'atomlist',atomlist, 'xyzlist' ,xyzlist, ...
'ecutwfc',30, 'name','Cu', 'temperature', 300, 'kpts', kpts, ...
'wks', wks);
After setting some optional parameters
opt = setksopt;
opt.maxscfiter = 30;
We perform the calculation by passing cry and opt to the scf function
[cry,H,X,info] = scf(cry,opt);
After the SCF iteration converges (up to \(10^{-8}\) default error tolerance in potential) in 15 iterations.
SCF iter 15:
Rel Vtot Err = 2.147e-09
Total Energy = -1.8198831034619e+02
Convergence is reached!
Etot = -1.8198831034619e+02
Entropy = -4.6058235803353e-02
Ekin = -3.3554142343570e+01
Eewald = -1.1961236762130e+02
Ealphat = 0.0000000000000e+00
Ecor = -4.9871469579648e+01
Ecoul = 3.2102334213766e+01
Exc = -1.1006606779632e+01
Efermi = -7.3851302930231e-01
--------------------------------------
Total time used = 3.840e+02
||HX-XD||_F = 9.365e-09
We can plot the eigenvalues of the converged Hamiltonian and show the band structure by using the following commands (contained in cu_nscf_k.m)
endkpts = { 0.00 0.00 0.00 '$\Gamma$'
0.00 0.50 0.50 'X'
0.25 0.75 0.50 'W'
0.50 0.50 0.50 'L'
0.375 0.75 0.375 'K'};
plotband(cry,H.rho,info.efermi,15,endkpts,[0.4777 1]);