B. Civalleri 
Overview and goals 
The simplest type of calculation in CRYSTAL is a
single point energy calculation. It consists in the calculation of the
wave function and energy
for a given crystalline system with a wellspecified geometric structure, i.e. at a single,
fixed point on the potential energy surface. The computed energy is the total energy,
sum of the electronic energy and nuclear repulsion energy.
To perform a single point
energy calculation a welldefined level of calculation must be specified. A level of
calculation is uniquely defined by the combination of a theoretical method with a basis set.
The analysis of the computed wave function, an
important step in predicting properties, is presented in the tutorial "One
electron properties".
Many aspects of the crystalline behavior can be understood when the total energy of the system is known, or more often, when the energy difference between two or more electronic or nuclear configurations is known as a function of some parameters. Other important properties are related to the derivatives of the total energy. For instance, first derivatives with respect to the nuclear atomic coordinates give the forces acting on the atoms and allow to obtain the equilibrium geometry when a geometry optimization is performed. Other important properties are related to second derivatives of the total energy of the system with respect to, for instance: volume, parameters of the unit cell, and the fractional coordinates of the atoms. The related observables are the bulk modulus, the elastic constants and vibrational properties, respectively.
Some of the energy differences of interest are listed below:
Computed energy  System 1  System2  Example 
Cohesive  Bulk  Isolated atoms  Ionic, covalent 
Interaction  Bulk  Isolated molecules  Molecular crystals 
Stability  Bulk  Bulk  Polymorphism 
Surface  Bulk  Slab  MgO(100) 
Surface stability  Slab  Slab  MgO(100) vs MgO(110) 
Adsorption  Slab + molecule  Slab, molecule  CO on MgO(100) 
Substitution  Defective solid  Perfect system, defect  C in Si 
Super exchange  AFM bulk  FM bulk  NiO 
Solid state reaction  Bulk1,Blk2  Bulk12  MgO + Al_{2}O_{3} MgAl_{2}O_{4} 
We will focus mainly on the total energy calculation and a few related energy differences for test systems: MgO, Si, Be, and Urea. They are simple examples of the most important categories of solid compounds: ionic, covalent, metallic and molecular, respectively.
For a list of references on CRYSTAL applications to solidstate systems and processes see the CRYSTAL web site: http://www.crystal.unito.it
This tutorial consists of a worked example and three exercises. The example concerns the calculation of the energy of formation and a few other related properties of an ionic crystal, MgO. The same strategy as for the example holds for the proposed exercises, in which a covalent system, Silicon, a metallic solid, Be, and a molecular crystal, urea, are considered. Since example and exercises cover different kinds of chemical bond a direct comparison among the computed properties can be performed showing how CRYSTAL can be a powerful tool in computational solidstate chemistry.
Example 1. A simple inorganic system: MgO 
MgO has been selected as a case system for ionic compounds.
Technical notes
Input decks and some output files for the proposed example and exercises
can be found in mgo.zip.
HartreeFock (HF) calculation
The total energy of the system, computed at HF level, is given by these lines
of the output:
=== SCF ENDED CONVERGENCE ON ENERGY E(AU) 2.7466415145499E+02 CYCLES 10 TOTAL ENERGY(HF)(AU)( 10) 2.7466415145499E+02 DE2.6E06 DP 3.0E04 PX 3.2E03 
Total energy is in hartree. The number of SCF cycles as well as SCF convergence criteria are printed
Density functional theory (DFT) calculation
In CRYSTAL09, the total energy of the system, at DFT level, is computed
at each cycle of the SCF iteration.
At the end, the DFT total energy is given by these lines of the output:
=== SCF ENDED CONVERGENCE ON ENERGY E(AU) 2.7402721224954E+02 CYCLES 7 ENERGY EXPRESSION=HARTREE+FOCK EXCH*0.00000+(LDA EXCH)*1.00000+PZ CORR TOTAL ENERGY(DFT)(AU)( 7) 2.7402721224954E+02 DE4.5E07 DP 1.5E04 
The value is in hartree. The adopted exchangecorrelation functional is also reported.
Run the MgO calculation at both HF and DFT level of theory and locate the lines
above in the output.
Note: For DFT calculation use
the LEVSHIFT and FMIXING options.
For instance, for a LDA (LDAPZ)
calculation the DFT input is (see the
CRYSTAL09 User's Manual
for details):
DFT  DFT input block  
EXCHANGE
LDA CORRELAT PZ ENDDFT 
Keyword to define the exchange functional
Selected exchange functional Keyword to define the correlation functional Selected correlation functional End of the DFT input block 

END  End of method input section  
LEVSHIFT
6 1 FMIXING 30 
Keyword to apply an energy shift
Energy shift (6*0.1 Hartrees), locking active Keyword to mix the Fock/KS matrix percent of Fock/KS matrix mixing 

END  End of SCF input section 
Draw the total energy versus the lattice parameter (xy chart) and locate the equilibrium lattice parameter. By starting from the experimental geometry, repeat it at different levels of theory: HF, LDA (LDAPZ), GGA (PWGGA) and with a hybrid (B3LYP) method.
Hints: Span the lattice parameter in the region: 4.12  4.26 (15 pts) at HF, GGA and B3LYP level, and 4.07  4.21 (15 pts) at LDA level. Use the xmgr (or gnuplot) program to plot the curve.
Here are the equilibrium lattice parameters for the adopted level of theory:
MgO  HF  LDA  GGA  B3LYP  Exp. 
a  4.191  4.139  4.208  4.205  4.195 
HF gives a lattice parameter in good agreement with the experimental value.
DFT results show that LDA gives a lattice parameter which is underestimated
with respect to the experimental data, whereas gradientcorrected and hybrid
methods lattice parameters are in better agreement although slightly overestimated.
where E_{a} is the atomic energy for each atom (or ion) belonging to
the crystal unit cell and N is the number of atoms in the unit cell.
The choice of an atomic reference energy is a fundamental and controversial point
in the evaluation of the cohesive energy.
As the basis set is not complete, two different basis sets are used for isolated
atoms
and for bulk. The same basis set is used for inner electrons, whereas extra
functions must be added to the bulk basis set for an accurate description of the valence
electron of the isolated atoms.
Three calculations are necessary to compute
cohesive energy of MgO: MgO bulk, oxygen atom, and magnesium atom.
At HF level, calculate the atomic energies of isolated atoms by crystal
(ATOMHF keyword: atomic wave function only is computed, no periodic
calculation). Test one
more sp function for better describing the atoms. Repeat this calculation
with a 8511G basis set for Mg and a 8411G basis set for O. The required basis sets
can be found at the CRYSTAL web site:
http://www.crystal.unito.it/Basis_Sets/ptable.html
CRYSTAL does not allow calculation of atomic wave function with DFT Hamiltonian
a) Modify the MgO input deck to create super cells containing 8, 16 and 32 atoms. Run each calculation at HF level and then compute the total energy per formula unit. Compare the results with the MgO bulk calculation. Hint: To speed up the calculation for the 32 atoms super cell reduce the value of the shrinking factor from 8 to 4.
Here is the HF total energy for the MgO super cell as well as
the total energy per formula unit (hartree):
2  8  16  32  
Total energy (TE)  274.66415  1098.65661  2197.31322  4394.62644 
TE per formula unit  274.66415  274.66415  274.66415  274.66415 
Exercise: Repeat such calculations at LDA level. Is the energy stable also for DFT methods?
b) Have a look at the CPU time and at the elapsed time. The total CPU time should be
separated into two contributions: on one hand, the CPU time to compute the integrals,
on the other hand, the CPU time to perform the SCF iteration.
To locate the two contributions, digit the command:
grep 'TTT MONMAD' *.out
and:
grep 'TTT END' *.out
The former is the CPU time for the integrals calculation, the latter is the total CPU time.
To obtain the SCF CPU time, then, subtract the CPU time for the integrals calculation
to the total CPU time.
How does the cost of the calculation increase with the size of the system?
Extensions
Many extensions of the proposed example are possible in which the knowledge of the total energy is a fundamental step.
For instance, in the field of surface science, different surfaces (slab models) can be cut out from the bulk (e.g. MgO(100) and MgO(110)). By comparing the total energy of bulk and slab models the surface energy can be easily computed. When the surface energy is available the different stability of surfaces can be obtained as well as the dependence of the surface energy on the thickness of the slab model. See "Surfaces and adsorption".
Another important application in materials science is the study of
defective systems, in which a defect is present in the bulk structure. See
"Defects".
Exercise 1. A simple covalent system: Si
Silicon is an important material which finds many application as semiconductor in electronics.
Repeat the calculation previously proposed for MgO and compare the results with the MgO ones. Results can be interpreted in the same way.
Technical notes
Input decks and some output files for the proposed exercise can be found
in si.zip.
Exercise 2. A metallic solid: Be
Be is a simple metallic solid. The Fermi energy must be computed by integration: the value of the Gilat shrinking factor, ISP (see "SCF") must be set to 2*IS (shrinking factor PackMonkhorst net).
Repeat the calculation proposed for MgO and Si and compare results with the previous ones. Note that Be bulk structure has two lattice parameters to be optimized. In this case, the geometry optimization of the cell may be performed by fixing the c/a ratio at its experimental value (1.5677) and by plotting the total energy vs the volume of the cell.
Technical notes
Input decks and some output files for the proposed exercise can be found
in be.zip
Exercise 3. A molecular crystal: Urea
The formation energy of a molecular crystal is obtained as difference between the energy of the bulk and the sum of the energies of isolated molecules.
Two calculations are required: urea crystal and urea molecule (there are two equivalent molecules in the crystallographic cell).
By starting at HF/321G level and adopting the experimental structure
determined at 12 K (see "Geometry input"
for the structural data):
1) Compute the total energy for the urea crystal;
2) For the molecular calculation, use the keyword MOLECULE to extract
a molecule from the molecular crystal and run CRYSTAL
Then, calculate the formation energy of the urea crystal.
Check the computed formation energy for the basis set dependence. Repeat the HF calculation with a 621G(d,p), a 631G(d,p), and a 6311G(d,p) basis set.
Compare the HF results with DFT calculations. Repeat the calculations above by using the following exchangecorrelation functionals: SVWN, PWGGA, B3LYP.
To
check the importance of the Basis Set Superposition Error (BSSE) use the
MOLEBSSE keyword: a molecule is extracted from the bulk, surrounded by
the basis functions of the neighboring molecule, up to a given distance from the
outermost atoms.
Technical notes
Input decks and some output files for the proposed exercise can be found
in urea.zip