# Crystals

A crystal is a bulk structure consisting of repeated units periodically extended in one to three dimensions. An infinitely long wire is an 1d crystal; an infinitely large film is a 2d crystal; and an infinitely large 3d structure is a 3d crystal. Many DFT packages are for calculating these periodic systems. The Bloch theorem allows one to focus on a single repeating unit of an infinitely large periodic crystal by the k-sampling technique.

The smallest unit in a periodic crystal is a primitive cell. For
instance the primitive cell of Si has two Si atoms; the primitive cell
of Fe has one atom, *etc*. Once the atomic positions of the atoms inside
the primitive cell is specified, lattice vectors are used to
periodically extend the primitive cell to infinity forming a crystal.
Therefore one has to define these lattice vectors. Since the primitive
cell is the smallest unit, it is desirable to use primitive cell in DFT
calculations to save some computation time. It is true that the smaller
the cell is, the larger the number of k-sampling is needed. But even
using large enough k-sampling will still save some time. On the other
hand, since we are dealing with a periodic crystal, we can actually use
any unit-cell (called supercell) as our basic unit to extend to the full
infinitely large crystal. These supercells typically have more atoms
than that contained in the primitive cell, thus the computation will
cost more time. For cubic crystals, it is some times simpler to set up
atomic positions of a supercell having a cubic shape.

Similar to that discussed in section Molecules, **NanoDCAL**
solves the Kohn-Sham equation for the supercell inside a computation box
of linear dimensions \((L_x,L_y,L_z)\). Because we are dealing with
a periodic structure, the values of \((L_x,L_y,L_z)\) must guarantee
periodicity when extended to form the crystal.
Fig. 5 is a 2d schematic of a bulk
crystal.

## Basic input file

The basic input file for bulk calculation is very similar to that of
molecule calculation discussed section Molecules. There are
some special points to be noted and we will focus on them in the
following discussions. There are several examples of bulk calculations
and let’s go to `/examples/Group2_Fe_Bulk/`

which is for computing
electronic structure of Fe crystal.
Fe scf input shows the `scf.input`

for the
bulk Fe calculation.

Since Fe is a magnetic metal, line-3 of
Fe scf input indicates that this is a
collinear spin calculation. Line-4 uses Matlab function eye(3) to set up
a \(3\times 3\) identity matrix, and then multiply it by 5.22 au.
This is to set up the supercell vectors to be 5.22 au in all three
dimensions. Here, the atomBlock has two Fe atoms forming the supercell.
For BCC structure, one Fe atom sits at the origin and the other at the
body centred location. A new parameter appears in line-6, *i.e.* the
keyWord SpinPolarization. Spin polarization is defined as,

where \(Q_{\uparrow,\downarrow}\) are charges in the spin-up and -down channels. For spin dependent calculations, it is usually helpful to give some initial value to \(\eta\) and here \(\eta = 0.2\). During the DFT iterations, \(\eta\) will be determined self-consistently. Lines 10-11 indicate the units to be used, thus the numbers in the atomBlock are in au. Please refer to Section Default units for the default units. Here we override the default by lines 10-11.

Finally, lines 12-14 tells nanodcal to do self-consistent DFT iterations using the multisecant mixing method [ML08]. Refer to Input reference section for more details.

The input file scf.input for bulk Fe calculation:

```
% The input file for the self-consistent field (hamiltonian) calculation.
% The file Fe_DZP.nad is needed for the calculation.
% type "nanodcal -parameter ? \.SCF\." for a complete list of input parameters
% which are used for the SCF calculation.
% no default value for the calculation name; must be given.
calculation.name = scf
% type "nanodcal -parameter ? calculation.name" for a complete list of
% calculation names.
% some parameters to define the system
system.name = ’FeBulk’
system.spinType = ’CollinearSpin’
system.centralCellVectors = eye(3)*5.2202
system.atomBlock = 2
AtomType X
Y
Z orbitaltype SpinPolarization
Fe 0.000000 0.000000 0.000000 DZP
0.200000
Fe 2.610100 2.610100 2.610100 DZP
0.200000
end
% type "nanodcal -parameter ? system\." for a complete list of input
% parameters which are used to define the system.
% the length and evergy units used for the input and output.
calculation.control.lengthUnit = ’au’
calculation.control.energyUnit = ’au’
% some parameters for SCF calculation
calculation.SCF.mixer.class = cMixerMultiSecant
calculation.SCF.mixer.parameter.beta = 2e-1
calculation.SCF.mixer.parameter.maxhistory = 3
% type "nanodcal -parameter ? \.SCF\." for a complete list of input parameters
% which are used for the SCF calculation.
```

Standard output of the Fe bulk run:

```
System Summary:
System Name: FeBulk
System Type: periodic-3D, collinear spin
System Symmetry: Oh
# of Atoms in Central Cell: 2
# of Electrons of Central Cell Atoms: 16
# of Basis of Central Cell Atoms: 30
Parameters Summary:
XC functional type: LDA_PZ81
central cell base vectors:
v1 = (5.220200e+000, 0.000000e+000, 0.000000e+000)
v2 = (0.000000e+000, 5.220200e+000, 0.000000e+000)
v3 = (0.000000e+000, 0.000000e+000, 5.220200e+000)
realspace grid numbers: [17 17 17]
realspace grid volume: 0.028954 Bohr3̂
k-space grid numbers: [15 15 15]
electronic temperature: 100 K
....
```

## Basic standard output

As already demonstrated in Section Basic standard output, we run the Fe bulk calculation by typing:

```
nanodcal scf.input
```

Many lines of output appear on the computer screen similar to those
shown in Benzene output.
Fe scf input (b) shows the system and parameter
summary blocks of the standard output. In particular, since we have
given a name to this run in `scf.input`

, line-2 of
Fe scf output (b) presents it. Line-3
indicates that it was a periodic 3d bulk run with collinear spin
configuration. Lines 4-6 are similar to that in
Benzene output.

In the parameter block, only line-14 is new in comparison to
Listing 2. Since
Listing 2 is for a molecular calculation,
there is no k-sampling. But for the bulk Fe run, it is a periodic system
and appropriate k-sampling must be done. Thus, line-14 of
Fe scf output indicates that we did k-sampling using
a k-grid of \([15,15,15]\) points. How did **NanoDCAL** decide the
k-grids? Again, as discussed in Section
Basic standard output, **NanoDCAL** sets up k-grid by a
precision control parameter:

```
calculation.control.precision = low, or normal, or high
```

These precision parameters tell **NanoDCAL** the numerical accuracy to
be enforced. In Section Basic standard output we discussed
how this affects the energy cutoff that determines the real space grid
(see (1) and associated discussions), here we
see how the k-grid is determined. In particular, if `scf.input`

does not
specify what precision to use, nanodcal applies the default which is
‘normal’. For k-grids, low, normal or high precisions mean 20 Bohr, 40
Bohr and 80 Bohr upper cutoff for 1/k values. As a result, it was
decides that k-grid should be \(15^3\) points for normal precision
as shown in Fe scf output. On the
other hand, the user can enforce his own k-grid by adding the following
keyWord in `scf.input`

:

```
calculation.k_spacegrids.L_cutoff = 55
```

this lines enforces a higher cutoff of 55 Bohr than the default normal
precision of 40 Bohr. As a result, the k-grid will become \(21^3\)
instead of \(15^3\) for the Fe example. One can also directly tell
**NanoDCAL** what k-grid to use by adding the following line to
`scf.input`

:

```
calculation.k_spacegrids.number = [8 12 14]
```

This would set up a k-grid of 8 points in \(k_x\), 12 points in \(k_y\) and 14 in \(k_z\).

Finally, the actual self-consistent iteration proceeds just like that in Listing 2 discussed in Section Basic standard output, and will not be repeated here. In addition, for bulk systems the essential output files of nanodcal are the same as that shown in Table 5.

## Calculating band structure

Once the SCF computation is converged, one may wish to obtain the band
structure of the crystal. In the `/examples/Group2_Fe_Bulk/`

directory,
you can find a four line input file `band.input`

for this purpose, as
shown below:

```
% The input file for the electronic structure (band structure) calculation.
% The crystal structure is automatically determined by the code.
% The file NanodcalObject.mat is needed for the calculation.
% type "nanodcal -parameter ? calculation.bandstructure" for a complete list
% of input parameters which are used for the band structure calculation.
% the hamiltonian of the cNA_LCAO object saved in the file
% NanodcalObject.mat will be used for the calculation.
system.object = NanodcalObject.mat
% Another way to define the hamiltonian is through the input parameter
% system.hamiltonian.
% type "nanodcal -parameter ? system\." for a complete list of input
% parameters which are used to define the system.
% no default value for the calculation name; must be given.
calculation.name = bandstructure
% type "nanodcal -parameter ? calculation.name" for a complete list of
% calculation names.
% some input parameters for the calculation.
calculation.control.energyUnit = eV
calculation.bandStructure.plot = 1
% type "nanodcal -parameter ? calculation.bandstructure" for a complete list
% of input parameters which are used for the band structure calculation.
```

Here, line-1 tells **NanoDCAL** to start from an object file obtained
from a previous self-consistent DFT run. Line-2 requires a calculation
of band structure. Line-3 requires the energy unit to be eV in the
present calculation. Note that even though we used atomic unit for the
SCF run (see Fe scf input), line-3 here will
convert it to eV during the present band structure calculation and,
hence, line-4 shall plot a figure using eV as energy unit. Type:

```
nanodcal band.input
```

a band structure of the Fe crystal is produced and plotted. The band
structure calculation uses the same control parameters stored in the
system object output file `NanodcalObject.mat`

. Note that NanoDCAL
does band unfolding automatically so that the plotted bands are already
unfolded.

Now, we return to entry 5 of
Table 5, namely
the output file `CalculatedResults.mat`

. Because we have calculated band
structure of Fe, this file now contains the information about the band
structure. By loading this file into Matlab, you can examine the stored
data by typing the following in the Matlab command window:

```
load CalculatedResults.mat
data
data.bandStructure
data.bandStructure.data
...
```

Here, each line provides increasingly more information.

## Calculating complex band structure

In typical electronic calculations, one is after the band structure of the crystal. For quantum transport, one may be interested in the complex band structure. A complex band structure is \(E=E(k_1,k_2,k_3)\) where any of the \(k_i\) is a complex value. Why complex band structure is interesting? For example, suppose an electron tunnels into the band gap of an insulator: its wave function becomes evanescent and it’s \(k\) becomes an imaginary hence complex number. Complex bands can connect real bands across the band gap.

In the `/examples/Group2_Fe_Bulk/`

directory, there is an input file
`cband.input`

which essentially has the following lines:

```
% The input file for the complex-band structure calculation.
% The file NanodcalObject.mat is needed for the calculation.
% type "nanodcal -parameter ? complexBandStructure" for a complete list
% of input parameters which are used for the complex-band structure
% calculation.
% the hamiltonian of the cNA_LCAO object saved in the file
% NanodcalObject.mat will be used for the calculation.
system.object = NanodcalObject.mat
% Another way to define the hamiltonian is through the input parameter
% system.hamiltonian.
% type "nanodcal -parameter ? system\." for a complete list of input
% parameters which are used to define the system.
% no default value for the calculation name; must be given.
calculation.name = cbs
% type "nanodcal -parameter ? calculation.name" for a complete list of
% calculation names.
% some input parameters for the calculation.
calculation.control.energyUnit = eV
calculation.complexBandStructure.plot = true
calculation.complexBandStructure.energyPoints = -20:0.1:20
calculation.complexBandStructure.numberOfKpoints = [1 1 1]
% type "nanodcal -parameter ? complexBandStructure" for a complete list
% of input parameters which are used for the complex-band structure
% calculation.
```

Running **NanoDCAL** with this input produces a complex band structure
having energy range \(E=-20\)eV to \(+20\)eV with steps of
\(0.1\)eV, versus real and complex \(k\) values along the
default direction. The default direction is the transport direction, see
Input reference section for more details.

## Calculating total energy

To calculate total energy of the system, let’s open the input file
`energy.input`

in the `/examples/Group2_Fe_Bulk/`

directory which
essentially has the following lines:

```
% The input file for the (decomposed) energy calculation.
% The file NanodcalObject.mat is needed for the calculation.
% type "nanodcal -parameter ? calculation.totalEnergy" for a complete list
% of input parameters which are used for the energy calculation.
% the hamiltonian and realspace potential of the cNA_LCAO object saved in
% the file NanodcalObject.mat will be used for the calculation.
system.object = NanodcalObject.mat
% no default value for the calculation name; must be given.
calculation.name = totalenergy
% type "nanodcal -parameter ? calculation.name" for a complete list of
% calculation names.
% some input parameters for the calculation.
calculation.control.energyUnit = eV
calculation.totalEnergy.decomposed = true
calculation.totalEnergy.plot = true
% type "nanodcal -parameter ? calculation.totalEnergy" for a complete list
% of input parameters which are used for the energy calculation.
```

Here, line-4 tells **NanoDCAL** to also calculate and display individual
terms contributing to the total energy. If this line is absent, only
total energy will be given.

## Calculating transmission channel

Next, we can determine the number of Bloch states per unit cell in the Fe crystal. This information is quite useful for transport analysis, for instance when Fe is used as device electrodes. Suppose we wish to consider transport in the z direction. Because our Fe electrode is a perfect crystal, every Bloch state - called transmission channel of the electrode, goes through the crystal from \(z=-\infty\) to \(z=+\infty\) without scattering. Hence for any given energy \(E\) and transverse momentum \(k_1,k_2\), the number of transmission channel \(N_{trc}(E,k_1,k_2)\) is an integer.

In the `/examples/Group2_Fe_Bulk/`

directory, there is an input file
transChannel.input which essentially has the following six lines:

```
% The input file for the transmission channel calculation.
% The file NanodcalObject.mat is needed for the calculation.
% type "nanodcal -parameter ? transChannel" for a complete list
% of input parameters which are used for the transmission channel
% calculation.
% the hamiltonian of the cNA_LCAO object saved in the file
% NanodcalObject.mat will be used for the calculation.
system.object = NanodcalObject.mat
% Another way to define the hamiltonian is through the input parameter
% system.hamiltonian.
% type "nanodcal -parameter ? system\." for a complete list of input
% parameters which are used to define the system.
% no default value for the calculation name; must be given.
calculation.name = trc
% type "nanodcal -parameter ? calculation.name" for a complete list of
% calculation names.
% some input parameters for the calculation.
calculation.transChannel.plot = true
calculation.transChannel.energyPoints = [0]
calculation.transChannel.kSpaceGridNumber = [30 30 1]
calculation.control.energyUnit = eV
% type "nanodcal -parameter ? transChannel" for a complete list
% of input parameters which are used for the transmission channel
% calculation.
```

By typing `nanodcal transChannel.input`

, **NanoDCAL** is told to
calculate transmission channels at \(E=0\) (Fermi level of the Fe
crystal), with \(30\times 30 =900\) k-points in the two transverse
directions. A plot of \(N_{trc}(E,k_1,k_2)\) at \(E=0\) versus
\(k_1,k_2\) is produced for both spin-up and spin-down channels,
shown in Fig. 6 is for the spin-up channel
produced by **NanoDCAL**.

## Calculating density of states

As detailed in the Input reference section, density of states (DOS) can
be expressed in different ways. In the `/examples/Group2_Fe_Bulk/`

directory, there is an input file dos.input which essentially has the
following lines:

```
% The input file for the density of states calculation.
% The file NanodcalObject.mat is needed for the calculation.
% type "nanodcal -parameter ? calculation.densityOfStates" for a complete list
% of input parameters which are used for the density of states calculation.
% the hamiltonian of the cNA_LCAO object saved in the file
% NanodcalObject.mat will be used for the calculation.
system.object = NanodcalObject.mat
% Another way to define the hamiltonian is through the input parameter
% system.hamiltonian.
% type "nanodcal -parameter ? system\." for a complete list of input
% parameters which are used to define the system.
% no default value for the calculation name; must be given.
calculation.name = dos
% type "nanodcal -parameter ? calculation.name" for a complete list of
% calculation names.
% some input parameters for the calculation.
calculation.densityOfStates.isIntegrated = true
calculation.densityOfStates.energyRange = [-0.1, 0.1]
calculation.densityOfStates.kSpaceGridNumber = [1 1 1]’*8
calculation.densityOfStates.whatProjected = ’Ell’
calculation.densityOfStates.indexProjected = {0 1 2}
calculation.densityOfStates.whatLocalized = ’Point’
calculation.densityOfStates.regionPosition = [0,0,0]
calculation.densityOfStates.regionVectors = eye(3)*5.2202
calculation.densityOfStates.regionGridNumber = [50 50 1]
calculation.densityOfStates.plot = true
calculation.control.lengthUnit = ang
calculation.control.energyUnit = eV
% type "nanodcal -parameter ? calculation.densityOfStates" for a complete list
% of input parameters which are used for the density of states calculation.
```

Lines 2-4 tell **NanoDCAL** to construct DOS from \(E=-5eV\) to
\(+5eV\) (Fermi level is shifted at \(E=0\)) using 200 points,
and line-4 plots the DOS versus energy. This is the typical DOS plot.

Lines 8-16 were commented out. Suppose we instead comment out lines 3-4
and uncomment lines 8-16, then we would obtain DOS on a real space grid
(the local DOS) and integrated over energies. So, line-8 tells
**NanoDCAL** to integrate over \(E\) in a range shown in line-9. How
about k? As there is nothing said not to integrate it, k is integrated
by default using a k-grid in line-10. But if you wish to resolve the k,
you must explicitly specify it by adding a line to the above:

```
calculation.densityOfStates.isResolved = true
```

Line-11 tells **NanoDCAL** to project DOS on orbitals with the given
angular momentum specified in line-12. These angular momentum states are
called projectors which will appear in the DOS plots. Line-13 tells
nanodcal to calculate DOS on real space grid points and lines 14-16
specify the real space region and points. These parameters can take
other values, refer to the Input reference section. There are many possible
projections of DOS, and some practice is needed to get familiar with
them.

## Calculating potential

One can inspect the potential terms of the calculated Hamiltonian. In
the `/examples/Group2_Fe_Bulk/`

directory, there is an input file
potential.input which essentially has the following lines:

```
% The input file for the realspace coulomb (delta-hartree) potential.
% The file NanodcalObject.mat is needed for the calculation.
% type "nanodcal -parameter ? calculation.potential" for a complete list
% of input parameters which are used for the potential calculation.
% the realspace potential of the cNA_LCAO object saved in the file
% NanodcalObject.mat will be used for the calculation.
system.object = NanodcalObject.mat
% no default value for the calculation name; must be given.
calculation.name = potential
% type "nanodcal -parameter ? calculation.name" for a complete list of
% calculation names.
% some input parameters for the calculation.
calculation.control.energyUnit = eV
calculation.potential.whatPotential = ’Vdh’
calculation.potential.plot = 1
% type "nanodcal -parameter ? calculation.potential" for a complete list
% of input parameters which are used for the potential calculation.
```

Line-3 tells **NanoDCAL** to calculate and plot the Hartree potential
\(V_{\delta H}\) which satisfies the following Poisson equation,

where

is the neutral atom density obtained by summing over individual neutral atom densities \(\rho^{{NA,I}}(\textbf{r}-\textbf{R}_{I})\). We refer to Theory section for more details.

By replacing ‘Vdh’ in line-3 with other values, one can calculate and display other terms of the Hamiltonian. Please consulate Input reference for more values.

## Calculating charge

One can inspect the charge in the system using the input file
`charge.input`

in the `/examples/Group2_Fe_Bulk/`

directory. This input
essentially has the following lines:

```
% The input file for the analysis on the electronic charge.
% The file NanodcalObject.mat is needed for the calculation.
% type "nanodcal -parameter ? calculation.charge" for a complete list
% of input parameters which are used for the charge calculation.
% the hamiltonian of the cNA_LCAO object saved in the file
% NanodcalObject.mat will be used for the calculation.
system.object = NanodcalObject.mat
% Another way to define the hamiltonian is through the input parameter
% system.hamiltonian.
% type "nanodcal -parameter ? system\." for a complete list of input
% parameters which are used to define the system.
% no default value for the calculation name; must be given.
calculation.name = charge
% type "nanodcal -parameter ? calculation.name" for a complete list of
% calculation names.
% some parameters for the charge analysis
calculation.charge.whatProjected = Ell
calculation.charge.indexProjected = {0 1 2}
% type "nanodcal -parameter ? calculation.charge" for a complete list
% of input parameters which are used for the charge calculation.
```

Here, line-3 tells **NanoDCAL** to calculate the number of electrons
with the angular momentum specified in line-4, namely project the charge
onto the s, p and d states (angular momentum 0, 1 and 2, respectively).

The results are saved in a text file `Charge.txt`

or MAT file `Charge.mat`

.
Opening these files, one notices that the spin polarization of the
charge is also saved. There are several ways of computing charge
information, for instance we may want to know the charge on a particular
atom or the total charge of the system. Please go to Input reference section
for these and other options.

## Summary

While calculation of crystal is important on its own right, for
**NanoDCAL** crystals also serve as electrodes of devices. For quantum
transport calculations, a first step is actually determining the
electronic structure and properties of the electrodes where crystal
calculations must be done. Many other control and system parameters can
be set by **NanoDCAL** so that very flexible calculations of crystal can
be done. The Input reference section lists all these possibilities.

Marks, L. D., & Luke, D. R. (2008). Robust mixing for ab initio quantum mechanical calculations. Physical Review B - Condensed Matter and Materials Physics, 78(7), 075114.