2. Core Hole Localisation: Azulene
In this part, you will learn how to perform corehole calculations for highly symmetrical systems. In such cases, the corestate wave functions are not localised on a single atom but rather form linear combinations of the single coreorbitals. However, The experimental evidence shows that the corehole is localised on one atom. Therefore, we must also force the corehole of our excited wave function to localise on one single atom.
As an example, to show the principal approach to this problem, we chose the azulene molecule. This molecule only contains carbon and hydrogen atoms, but due to its nonalternant topology, the chemical shift between the carbon atoms is sufficiently large to see differences. Azulene possesses C_{2v} point symmetry, and several pairs of atoms are symmetry equivalent. There is experimental data for the multilayer of azulene, to which we can compare our results later.^{1}
Scheme 3. Structure of the azulene molecule and numbering scheme used in this tutorial.
The following geometry.in
file can be used for setting up the azulene calculations:
atom 0.00000000 0.00000000 2.15610643 C
atom 1.15055599 0.00000000 1.34809461 C
atom 1.15055599 0.00000000 1.34809461 C
atom 0.74965762 0.00000000 0.00057279 C
atom 0.74965762 0.00000000 0.00057279 C
atom 1.59583767 0.00000000 1.10633748 C
atom 1.59583767 0.00000000 1.10633748 C
atom 1.26699680 0.00000000 2.46635052 C
atom 1.26699680 0.00000000 2.46635052 C
atom 0.00000000 0.00000000 3.05928669 C
atom 0.00000000 0.00000000 3.24440762 H
atom 2.18076818 0.00000000 1.69493655 H
atom 2.18076818 0.00000000 1.69493655 H
atom 2.66541089 0.00000000 0.87579776 H
atom 2.66541089 0.00000000 0.87579776 H
atom 2.11174769 0.00000000 3.15756315 H
atom 2.11174769 0.00000000 3.15756315 H
atom 0.00000000 0.00000000 4.15206631 H
In general, we follow the same steps as in Part 1.
Ground State Calculation
First, we perform the single point ground state calculation:
 Use the functional
PBE
and add atighttier2
basis set for carbon and hydrogen.  Add the basis set augmentation to the basis set definition of carbon
 Take care to use the keyword
spin collinear
 Once again, add
output cube eigenstate n
with n = 1 to 10 to the input file to check for the localisation of the wave functions  Run the calculation
By looking at the exported cube files, we can see that the wave functions of the core states are now mixed and not localised properly at single atoms anymore.
Using deltascf_projector
and Breaking the Symmetry
We can again use the deltascf_projector
command, but due to the symmetry equivalent atoms in azulene, we must now break the symmetry of the wave function first to force a proper localisation.
To break the symmetry, we need to use a trick: We differentiate the atom, where the wavefunction is to be localised, by increasing its nuclear charge in an intermittent calculation and carrying the localised wave function over to a later calculation with the physically correct core charge. This later calculation then keeps to the localised wave function but yields the soughtafter physically meaningful energy.
It turns out that the best results in the localisation are obtained if we perform two socalled initialisation calculations before the final hole calculation, yielding a threestep procedure:

The
init_1
calculation with the nuclear charge increased by +0.1 and the system charge set to +0.1. A full SCF cycle is performed and arestart_file
is written) 
The
init_2
calculation with the nuclear charge increased by +0.1 and the system charge set to +1.1.deltascf_projector
is used to initialise the core hole at the correct atom. It is enough to calculate 1 SCF step, therestart_file
from theinita
calculation is read and overwritten) 
The
hole
calculation with the nuclear charge back to the original value and the system charge set to 1.0.deltascf_projector
is used to deoccupy the corehole, a full SCF cycle is performed and therestart_file
from theinit_2
calculation is read. This calculation produces the corehole energy used to determine the binding energy according to equation (3)
To determine which atom the corehole is to be localised at, it has to be differentiated in the geometry.in
file, e.g. by calling it C1. As described above, the same basis set is used as for all the other atoms, but the control.in
file must now contain a second definition for the carbon basis set, called C1. We will use the C1 basis set definition to change the nuclear charge only for this specific atom.
Note on Boys Localisation: We are currently investigating using Boys' Localisation here to converge the core hole rather than using an additional partial nuclear charge. Until this is thoroughly tested, the partial core hole should suffice as a method in most circumstances.
Initialisation and Hole Calculations for deltascf_projector
Following the procedure laid out in the last section, now we perform the initialisation and hole calculations:
 Create 3 folders for each atom:
init_1
,init_2
, andhole
 Create a special
geometry.in
file for each atom where that atom is labelled as C1
Running the init_1
calculation requires the following settings:
 In the C1 basis set information, the nuclear charge is increased by 0.1 by setting
nucleus
from 6 to 6.1  To compensate this increased nuclear charge, the valence electron count of the C1 basis set is also increased by 0.1 by setting
valence 2 p
from 2 to 2.1.  In the general part of the control.in file, the keyword
charge = 0.1
is added to ensure charge neutrality is maintained  add
restart_write_only restart_file
andforce_single_restartfile .true.
to write out a restart file  Run all ten
init_1
calculations
Now we move on to the init_2
calculations
 Copy the
restart file
from eachinit_1
calculation into the correspondinginit_2
folder  the changes in the C1 basis set definition remain the same
 we change the system charge to
charge = 1.1
to remove one complete electron in addition to the compensation charge  we add
deltascf_projector 1 1 0.0 1 10
to introduce the core hole. Note that we do not need to change theKS_state
value, because for this system, the added nuclear charge will always force the corresponding 1s state to the lowest energy.  set
sc_iter_limit
to 1, because we only need to initialise the correct wave function  add
restart restart_file
andforce_single_restartfile .true.
to both load the restart file from the previousinit_1
calculation and overwrite it for use in the next calculation  Run all ten
init_2
calculations
Finally, we can perform the hole
calculations
 Copy the
restart_file
from theinit_2
calculations into the correspondinghole
folder  Remove the changes in the C1 basis set definition.
nucleus
goes back to 6 andvalence 2 p
goes back to 2  Change the system charge to
charge = 1.0
to have a regularly charged system deltascf_projector 1 1 0.0 1 10
to introduce the core hole remains the same. Set
sc_iter_limit
back to a normal value, e.g. 500, to perform a normal SCF cycle  Add
restart_read_only restart_file
andforce_single_restartfile .true.
to load the restart file from the previousinit_2
calculation  To check the localisation of the core hole after calculation, we can again use the
output cube spin_density
keyword  Run all ten
hole
calculations
With the total energies of these hole calculations and the total energy of the reference ground state single point calculation, we again calculate the core electron binding energies using equation (3).
deltascf_basis
Command
As an alternative to the deltascf_projector
procedure described above, there is an additional method implemented in FHIaims to force localisation of the corehole: the deltascf_basis
command. This command also uses the maximum overlap method (MOM) methodology to track the state with the occupation constraint. However, deltascf_basis
tries to force the occupation to a state similar to a specified atomic basis function in the system rather than a previously determined KS state.
The deltascf_basis
keyword follows the following syntax:
deltascf_basis atom spin basis_type basis_n basis_l basis_m occ max_KS_state
atom
indicates the number of the specified atom (in thegeometry.in
file), on which the basis function will be centredspin
is again arbitrarily set to 1basis_type
is chosen to beatomic
, because we want to constrain the hole to an atomic orbital like holebasis_n
,basis_l
, andbasis_m
now define the basis function via quantum numbers, to specify our 1s state we choose accordinglybasis_n = 1
,basis_l = 0
andbasis_m = 0
occ
is again set to 0.0 to introduce a full core holemax_KS_state
is again set to the highest number of 1s states present, which is 10 for the azulene molecule
The deltascf_basis
command is often more convenient than the deltascf_projector
keyword because it needs less calculations and attention to detail when setting up the input files. However, deltascf_basis
is only available for the aperiodic case and can be problematic when used with hybrid functionals.
Hole Calculations using deltascf_basis
Performing the calculations is much simpler with the deltascf_basis
keyword:
 Generate one
hole
folder for each carbon atom  It is not required to specify a restart file
 The control.in file should contain:
charge = 1.0
to remove one electrondeltascf_basis N 1 atomic 1 0 0 0.0 10
where N is chosen from 1 to 10 for each hole calculationoutput cube spin_density
to check for the correct localisation of the core hole Run all corehole calculations
The total energies of these hole calculations and the total energy of the reference ground state single point calculation are again used to calculate the core electron binding energies via equation (3).
Comparison of deltascf_basis
and deltascf_projector
The following table contains the results for the absolute and relative binding energies for the ten atoms of azulene calculated both with the deltascf_basis
and deltascf_projector
commands.
Table 3. Absolute and relative C 1s coreelectron binding energies for azulene were calculated with the keywords deltascf_basis
and deltascf_projector
(PBE/tighttier2). All values are given in eV.
Catom  E_{B} (proj.)  ΔE_{B} (proj.)  E_{B} (basis)  ΔE_{B} (basis) 

C1  289.03  0.00  289.03  0.00 
C2  288.70  0.33  288.70  0.33 
C3  288.70  0.33  288.70  0.33 
C4  289.44  0.41  289.44  0.41 
C5  289.44  0.41  289.44  0.41 
C6  289.64  0.62  289.64  0.62 
C7  289.64  0.62  289.64  0.62 
C8  289.29  0.26  289.29  0.26 
C9  289.29  0.26  289.29  0.26 
C10  289.54  0.51  289.54  0.51 
We calculated the binding energies of all ten carbon atoms in azulene. As is visible in the resulting energies, there are two unique atoms and four symmetryequivalent pairs of atoms. Overall, only six hole calculations were necessary after all. However, in this case, it was good practice to calculate the binding energies for all atoms in order to use the symmetryenforced equivalence of binding energies as an internal sanity check.
We can also clearly see that deltascf_basis
and deltascf_projector
provide equivalent results when performed correctly.
Comparison with Experimental Data
There exists no experimental XPS data for azulene in the gas phase. However, experimental XPS data for a thin film multilayer of azulene was published by Klein et al. ^{1} Assuming that the intermolecular interactions in the multilayer only include weak dispersion interactions, which do not influence the corehole binding energies much, we can still compare our "gasphase" DFT data to these experiments.
Figure 2 shows this comparison; again, we used the pseudoVoigt scheme above to broaden the simulated ΔSCF data and simulate the experimental resolution. We used an FWHM of 0.7 eV and a GLratio of 0.4 for the azulene data. Then, we summed up the ten subpeaks to generate the whole spectrum.
Figure 2. Comparison of the calculated and experimental C 1s XP spectrum of azulene. The experimental data was measured for a thin film multilayer. The calculated data are for a gasphase molecule. The calculated data was shifted by 4.8 eV to match the experimental energy scale.
In the figure above, the experimental data is taken from Klein et al.^{1}
In Figure 2, we used the absolute energy scale of the experimental data. We needed to rigidly shift the calculated data by 4.8 eV to match this energy scale. Comparing the absolute energy scales of experiment and simulation would necessitate a cautious treatment of the work function for the experimental and simulated data. At the same time, systematic approximations in the DFT calculation would remain. The comparison in Figure 2, however, shows that the peak shape and, therefore, the relative binding energies of carbon atoms in different chemical environments are already in much better agreement than the size of the empirical shift parameter (4.8 eV) would suggest.
The experimental resolution in these condensed phase measurements allows only a comparison of the overall shape. However, a closer inspection is possible in highresolution XPS experiments in the gas phase. For example, in the gas phase spectra for naphthalene by Minkov et al., ^{2} two subpeaks with a difference in binding energy of 0.36 eV are resolved. This relative binding energy shift is in excellent agreement with the value of 0.32 eV obtained by ΔSCF (PBEtighttier2).
References

B.P. Klein, N. J. van der Heijden, S. R. Kachel, M. Franke, C. K. Krug, K. K. Greulich, L. Ruppenthal, P. Müller, P. Rosenow, S. Parhizkar, F. C. Bocquet, M. Schmid, W. Hieringer, R. J. Maurer, R. Tonner, C. Kumpf, I. Swart, and J. M. Gottfried, Phys. Rev. X 2019, 9, 011030. https://doi.org/10.1103/PhysRevX.9.011030 ↩↩↩

I. Minkov, F. Gel’mukhanov, R. Friedlein, W. Osikowicz, C. Suess, G. Öhrwall, S. L. Sorensen, S. Braun, R. Murdey, W. R. Salaneck, and H. Ågren, J. Chem. Phys., 2004, 121, 5733. https://doi.org/10.1063/1.1784450 ↩