12. Mulliken charge calculation

This section shows how to obtain the Mulliken charges for an AlP crystal.

12.1. Example: Mulliken populations of AlP

We first obtain the ground state density which determines the Hamiltonian to be used in the Mulliken charge calculation.

LCAO.status = true
info.savepath = 'results/alp_lcao_scf'
info.outfile = 'falp_lcao_scf.out'
info.calculationType = 'self-consistent'
atom.element = [1 2]
atom.fracxyz = 0.25*[0 0 0; 1 1 1]
domain.latvec = 5.4510*0.5*[0 1 1;1 0 1;1 1 0]
domain.lowres = 0.3
units.latvec = 'A'
element(1).species = 'Al'
element(1).path = './Al_TM_LDA.mat'
element(2).species = 'P'
element(2).path = './P_TM_LDA.mat'
functional.list = {'XC_LDA_X','XC_LDA_C_PW'}
kpoint.gridn = [4,4,4]

Let’s call the previous input file alp_lcao_scf. We pass this standard self-consistent calculation input to RESCU

rescu -i alp_lcao_scf

and get the ground state density.

The next step is calculating the Mulliken charges. This is done using the following input file

info.savepath = 'results/alp_lcao_mul'
info.outfile = 'falp_lcao_mul.out'
info.calculationType = 'mulliken'
kpoint.gridn = [8 8 8]
rho.in = 'results/alp_lcao_scf'

named alp_lcao_mul. RESCU performs a Mulliken calculation if the keyword info.calculationType is set to mulliken. The path to the converged density is passed using the keyword rho.in. System data like atomic and structural information are also read from rho.in. Remember to modify the output path, otherwise the self-consistent results will be overwritten. Finally, the keyword kpoint.gridn can take the same value as in the ground state calculation, but we choose a grid twice as fine here. We can now run RESCU

rescu -i alp_lcao_mul

After completing the calculation, the results are written to results/alp_lcao_mul.h5 in /LCAO/mullikenPop1. When doing a spin-dependent calculation, the spin-up populations are saved in /LCAO/mullikenPop1 and the spin-down populations are saved in /LCAO/mullikenPop2.

12.2. Analysis of the Mulliken populations

Table 12.2.1 LCAO.orbInfo field and the corresponding indexing.

Field	Index
Aorb	atom
Lorb	\(L\)
Morb	\(L_z\)
Oorb	basis
Rorb	cutoff radius
Sorb	species
Zorb	zeta

The Mulliken population matrix is imported as follows. In the MATLAB command window, type

h5path = './results/alp_lcao_mul.h5';
mul = loadDistArray(h5path,'/LCAO/mullikenPop1');
mul = mul.data;

and obtain the \(n_{nao}\times n_{nao}\) mul matrix.

Suppose we want to know the total charge of the first atom, then we need to know which orbitals are associated with it. The orbital information for each row or column is stored in LCAO.orbInfo in results/alp_lcao_mul.mat. The various fields of LCAO.orbInfo are summarized in Table 2.7.1. The field LCAO.orbInfo.Aorb contains the atomic indices, and hence we calculate the first atom’s charge as follows

Aorb = LCAO.orbInfo.Aorb
sum(sum(mul(:,Aorb == 1)))
ans =

   2.6670 - 0.0020i

Note that a partial sum such as the above is generally complex, but the sum of all entries in the Mulliken population matrix is real and equal to the number of valence electrons. The last point is not strictly satisfied if the calculation is performed on real space grids since the atomic orbitals do not form a complete basis.