CRYSTAL home page CRYSTAL home page

CO on MgO(001) (loading=1:4)

Surfaces and adsorption:
Modeling the surface - The slab model


Modeling the surface - The slab model


The bulk crystal (MgO)

We take MgO (rock salt structure) as a simple example.

The bulk structure is specified as:

CRYSTAL     Select a 3D periodic structure
0 0 0     Default options for space group
225     The space group
4.21     The lattice parameter in Angstrom
2     2 symmetry inequivalent atoms
12   0.0   0.0   0.0     The Mg atom at the origin
8     0.5   0.5   0.5     The O atom at ½ ½ ½ (fractional coordinates)

In the CRYSTAL03 output a description of the primitive cell contents and details on the direct space cell in the adopted Cartesian frame are printed, as follows:

CRYSTAL CALCULATION
 (INPUT ACCORDING TO THE INTERNATIONAL TABLES FOR X-RAY CRYSTALLOGRAPHY)
 CRYSTAL FAMILY                       :  CUBIC
 CRYSTAL CLASS  (GROTH - 1921)        :  CUBIC HEXAKISOCTAHEDRAL

 SPACE GROUP (CENTROSYMMETRIC)        :  F M 3 M

 LATTICE PARAMETERS  (ANGSTROMS AND DEGREES) - CONVENTIONAL CELL
        A           B           C        ALPHA        BETA       GAMMA
     4.21000     4.21000     4.21000    90.00000    90.00000    90.00000

 NUMBER OF IRREDUCIBLE ATOMS IN THE CONVENTIONAL CELL:    2

 INPUT COORDINATES

 ATOM AT. N.              COORDINATES
   1  12     0.000000000000E+00  0.000000000000E+00  0.000000000000E+00
   2   8     5.000000000000E-01  5.000000000000E-01  5.000000000000E-01
 *******************************************************************************
 << INFORMATION >>: FROM NOW ON, ALL COORDINATES REFER TO THE PRIMITIVE CELL
 *******************************************************************************

 LATTICE PARAMETERS  (ANGSTROMS AND DEGREES) - PRIMITIVE CELL
       A          B          C         ALPHA      BETA     GAMMA        VOLUME
    2.97692    2.97692    2.97692     60.0000   60.0000   60.0000      18.65462

 COORDINATES OF THE EQUIVALENT ATOMS (FRACTIONARY UNITS)

 N. ATOM EQUIV AT. N.          X                  Y                  Z

   1   1   1   12 MG    0.00000000000E+00  0.00000000000E+00  0.00000000000E+00
   2   2   1    8 O    -5.00000000000E-01 -5.00000000000E-01 -5.00000000000E-01

 NUMBER OF SYMMETRY OPERATORS         :   48
 ******************************************************************************   

This part of output reports the transformation of the crystallographic cell given in input (Figure 1) into a more convenient primitive unit cell (Figure 2).

Figure 1 The Crystallographic cell of MgO Figure 2 The Primitive Cell of MgO


Specifying the surface plane – Miller indices

Any crystallographic plane can be specified by three integers (h,k,l) – the Miller indices. The indices of any given plane can be worked out by considering the intersection of the plane with the crystallographic unit cell  vectors a, b, c. If x, y, z are the fractional coordinates of the intersections of the plane with a, b, c,  Miller indices are the smallest integers in the same ratio as (1/x, 1/y, 1/z). Here are some examples for cubic lattices:
 
 



Miller indices required in CRYSTAL03 input refer to the crystallographic cell.

Not all crystalline surfaces are physically stable or worthy of investigation. This is specially true for ionic or semi-ionic crystals. For example, in a polar surface consisting of charged layers alternating there is a net dipole moment per unit area, normal to the surface. Such a surface is unstable and can only be prepared with substantial reconstruction or with the adsorption of charged species. 

Exercise: See, for instance, the MgO (111) surface.

Choosing the slab thickness and the surface termination

In general having selected the Miller indices of a surface there will be a number of possible choices for how to terminate the solid at the surface. In the case of MgO we have selected the particularly simple 100 surface and there is only one symmetry distinct layer at Z=0. To select this layer, (LAYER 1), as the terminating layer for the surface and to create a 2D-periodic slab 3 layers thick enter the following at the end of the geometry section:

SLABCUT     Transform cell to have two lattice vectors in the chosen surface plane and create 2D translational symmetry
1 0 0     Miller indices of the surface plane
1 3     Select Layer 1 as the terminating layer and request a slab of 3 layers


SLABCUT cuts a 2D-periodic slab of the required thickness. Wave function and energy of 2D system are computed, if requested. The geometry of the surface layers and the surface unit cell for a slab with a single layer are displayed in the figure below.

Figure 3 Four primitive cells of the MgO(100) surface

CRYSTAL chooses to place the origin of the surface unit cell at the centre of the slab and has rotated the lattice vectors by 45° relative to the bulk crystallographic cell. This produces the minimum number of atoms and maximum symmetry description of the 2D slab - both factors are exploited to reduce the computer resources required.

The calculation of the wave function (if complete input is supplied) will take a few minutes to run. Using the slab directive above and the computational parameters in the file mgo100_3 (extension d12) the energy of the slab will be -823.9326404 Hartree. Note that the tolerances on the integrals (set by the TOLINTEG keyword) have been raised from the default values - a rather accurate calculation is required to converge the surface energy.

Exercise: Construct a MgO (110) 5 layers slab. Compare your result with:

 

 ATOMS IN THE ASYMMETRIC UNIT    6 - ATOMS IN THE UNIT CELL:   10
     ATOM              X/A                 Y/B             Z(ANGSTROM)
   1 T  12 MG   -1.110223024625E-16  5.000000000000E-01  2.981869296264E+00
   2 T   8 O    -1.110223024625E-16  0.000000000000E+00  2.981869296264E+00
   3 T  12 MG    5.000000000000E-01  0.000000000000E+00  1.490934648132E+00
   4 T   8 O     5.000000000000E-01  5.000000000000E-01  1.490934648132E+00
   5 T  12 MG    0.000000000000E+00  5.000000000000E-01  0.000000000000E+00
   6 T   8 O     0.000000000000E+00  6.884683795283E-17  0.000000000000E+00
   7 F  12 MG    5.000000000000E-01 -2.220446049250E-16 -1.490934648132E+00
   8 F   8 O     5.000000000000E-01  5.000000000000E-01 -1.490934648132E+00
   9 F  12 MG    1.110223024625E-16  5.000000000000E-01 -2.981869296264E+00
  10 F   8 O     1.110223024625E-16  1.110223024625E-16 -2.981869296264E+00

Note: by default, the coordinates of atoms in 2-D systems are fractionary in the directions with translational symmetry,  and Å in the non periodic direction.

Choosing the surface termination

MgO (100) surface is very simple case, there is 1 distinct layer only. When the structure is more complex, there are different possible terminations of the slab. The keyword REDEFINE, followed (next record) by the crystallographic (Miller) indices  of the surface, defines a new 3D cell, with  two lattice vectors perpendicular to the [hkl] direction. The indices refer to the Bravais lattice of the crystal; the hexagonal lattice is used for the rhombohedral systems, the cubic lattice for cubic systems (non primitive).

The atoms in the new reference cell are re-ordered according to their z coordinate, in order to recognize the layered structure, parallel to the (hkl) plane. The layers of atoms are numbered. 

Example: surface (001) of corundum (test 24 of CRYSTAL test cases):

TEST24 - CORUNDUM 001 2 LAYERS (3D-->2D) - SEE TEST04
CRYSTAL
0 0 0
167
4.7602  12.9933
2
13  0.        0.  0.35216
8   0.30624   0.  0.25
REDEFINE
0 0 1
TESTGEOM
END

The input defines a bulk corundum crystal, and the analysis of the layered structure parallel to the surface of crystallographic indices is requested. The program stops after the geometry step (TESTGEOM), no further input data required. The following output is obtained: 

CRYSTAL CALCULATION
 (INPUT ACCORDING TO THE INTERNATIONAL TABLES FOR X-RAY CRYSTALLOGRAPHY)
 CRYSTAL FAMILY                       :  HEXAGONAL
 CRYSTAL CLASS  (GROTH - 1921)        :  DITRIGONAL SCALENOHEDRAL

 SPACE GROUP (CENTROSYMMETRIC)        :  R -3 C

 LATTICE PARAMETERS  (ANGSTROMS AND DEGREES) - CONVENTIONAL CELL
        A           B           C        ALPHA        BETA       GAMMA
     4.76020     4.76020    12.99330    90.00000    90.00000   120.00000


 NUMBER OF IRREDUCIBLE ATOMS IN THE CONVENTIONAL CELL:    2

 INPUT COORDINATES

 ATOM AT. N.              COORDINATES
   1  13     0.000000000000E+00  0.000000000000E+00  3.521600000000E-01
   2   8     3.062400000000E-01  0.000000000000E+00  2.500000000000E-01
   
 LATTICE PARAMETERS  (ANGSTROMS AND DEGREES) - PRIMITIVE CELL
       A          B          C         ALPHA      BETA     GAMMA        VOLUME
    5.12948    5.12948    5.12948     55.2915   55.2915   55.2915      84.99223

 COORDINATES OF THE EQUIVALENT ATOMS (FRACTIONARY UNITS)

 N. ATOM EQUIV AT. N.          X                  Y                  Z

   1   1   1   13 AL    3.52160000000E-01  3.52160000000E-01  3.52160000000E-01
   2   1   2   13 AL    1.47840000000E-01  1.47840000000E-01  1.47840000000E-01
   3   1   3   13 AL   -3.52160000000E-01 -3.52160000000E-01 -3.52160000000E-01
   4   1   4   13 AL   -1.47840000000E-01 -1.47840000000E-01 -1.47840000000E-01

   5   2   1    8 O    -4.43760000000E-01 -5.62400000000E-02  2.50000000000E-01
   6   2   2    8 O     2.50000000000E-01 -4.43760000000E-01 -5.62400000000E-02
   7   2   3    8 O    -5.62400000000E-02  2.50000000000E-01 -4.43760000000E-01
   8   2   4    8 O     4.43760000000E-01  5.62400000000E-02 -2.50000000000E-01
   9   2   5    8 O    -2.50000000000E-01  4.43760000000E-01  5.62400000000E-02
  10   2   6    8 O     5.62400000000E-02 -2.50000000000E-01  4.43760000000E-01

 NUMBER OF SYMMETRY OPERATORS         :   12
 ******************************************************************************
 DIRECT-LATTICE FUNDAMENTAL VECTORS
             X         Y         Z
 A1       2.7483    0.0000    4.3311
 A2      -1.3742    2.3801    4.3311
 A3      -1.3742   -2.3801    4.3311
 

There are 10 atoms in the primitive cell of bulk corundum. The shape of the cell is modified in order to have two lattice vectors perpendicular to the [001] direction. The indices (001) given in input refer to the crystallographic cell, they are transformed into the equivalent indices (111) referred to the primitive cell:

**********************************************************************
 *            DEFINITION OF THE NEW LATTICE VECTORS                   *
 *            TWO OF WHICH BELONG TO THE SELECTED PLANE               *
 **********************************************************************

                       PLANE INDICES =   0   0   1
                   PRIMITIVE INDICES =   1   1   1

 NEW FUNDAMENTAL DIRECT-LATTICE VECTORS B1,B2,B3
 B1=      0 A1   -1 A2    1 A3
 B2=     -1 A1    1 A2    0 A3
 B3=      0 A1    0 A2   -1 A3

The coordinate of the atoms are computed in the new reference system:

 UNIT CELL ATOM COORDINATES :
 LAB AT.NO.       CARTESIAN (ANG)               CRYSTALLOGRAPHIC
   1  13       0.000    -2.748     4.086     0.295680  -0.352160   0.943520
   2  13      -2.380    -1.374     2.410    -0.295680  -0.147840   0.556480
   3  13      -2.380     1.374     0.245    -0.295680   0.352160   0.056480
   4  13       0.000     0.000     1.921     0.295680   0.147840   0.443520
   5   8       0.729     1.486     1.083     0.500000   0.443760   0.250000
   6   8       0.922    -1.374     1.083     0.193760  -0.250000   0.250000
   7   8      -1.651    -0.112     1.083    -0.193760   0.056240   0.250000
   8   8      -3.109    -2.860     3.248    -0.500000  -0.443760   0.750000
   9   8      -3.302     0.000     3.248    -0.193760   0.250000   0.750000
  10   8      -0.729    -1.262     3.248     0.193760  -0.056240   0.750000

the sequence of the atoms is reordered according to the z coordinate. Each value of z corresponds to a "layer". The top layer contains an Aluminum:

 ATOMS CLASSIFIED ACCORDING TO THE Z COORDINATE :

  LAYER   1 Z=   4.0865; LABEL AT.NO.              X,Y,Z (ANG.)          
                           1    13        0.00000    -2.74830     4.08648
  LAYER   2 Z=   3.2483; LABEL AT.NO.              X,Y,Z (ANG.)          
                           8     8       -3.10898    -2.85999     3.24832
                           9     8       -3.30244     0.00000     3.24832
                          10     8       -0.72888    -1.26246     3.24832
  LAYER   3 Z=   2.4102; LABEL AT.NO.              X,Y,Z (ANG.)          
                           2    13       -2.38010    -1.37415     2.41017
  LAYER   4 Z=   1.9209; LABEL AT.NO.              X,Y,Z (ANG.)          
                           4    13        0.00000     0.00000     1.92093
  LAYER   5 Z=   1.0828; LABEL AT.NO.              X,Y,Z (ANG.)          
                           5     8        0.72888     1.48584     1.08278
                           6     8        0.92234    -1.37415     1.08278
                           7     8       -1.65122    -0.11169     1.08278
  LAYER   6 Z=   0.2446; LABEL AT.NO.              X,Y,Z (ANG.)          
                           3    13       -2.38010     1.37415     0.24462
 

The repetitive unit contains 5 atoms, Al - 3 Oxygens - Al, in 3 layers.

A model of the surface (001) of corundum could be a slab of 6 layers (2 repetitive units), with Aluminum in the outer layer is::

TEST24 - CORUNDUM 001 2 LAYERS (3D-->2D) - SEE TEST04
CRYSTAL
0 0 0
167
4.7602  12.9933
2
13  0.        0.  0.35216
8   0.30624   0.  0.25
SLABCUT
0 0 1
1 6             the top layer is the layer # 1; the slab is 6 layers thick.
END
........ basis set, .... scf input follows

 

The figure shows top view and lateral view of 12 layers (001) slab of  a-alumina parallel to the (001) plane.

Alumina slab (001) - lateral view Alumina slab (001) - top view

A slab can be defined as 2D structure (see CRYSTAL03 User's Manual). 
A 2D input for the 6 layer slab defined above (starting from the 3D structure) is as follows:

 
TEST04 - CORUNDUM 111 (0001) 2 LAYERS- 2D
SLAB
66
4.7602
3
13 0. 0. 1.9209
8 0.333333333 -0.027093 1.0828
13 -0.333333333 0.333333333 0.2446
TESTGEOM
END


The output is as follows:

TEST04 - CORUNDUM 111 (0001) 2 LAYERS- 2D

 SLAB CALCULATION
 SYSTEM AND LATTICE                   :  HEXAGONAL
 TWO-SIDED PLANE GROUP N. 66          :  P -3
 CORRESPONDING SPACE GROUP N. 147

 LATTICE PARAMETERS  (ANGSTROMS AND DEGREES) - CONVENTIONAL CELL
           A              B            GAMMA
        4.76020        4.76020      120.00000

 NUMBER OF IRREDUCIBLE ATOMS IN THE CONVENTIONAL CELL:    3

 INPUT COORDINATES (X AND Y IN FRACTIONARY UNITS, Z IN ANGSTROMS)

 ATOM AT. N.              COORDINATES
   1  13     0.000000000000E+00  0.000000000000E+00  1.920900000000E+00
   2   8     3.333333330000E-01 -2.709300000000E-02  1.082800000000E+00
   3  13    -3.333333330000E-01  3.333333330000E-01  2.446000000000E-01


 LATTICE PARAMETERS  (ANGSTROMS AND DEGREES) - PRIMITIVE CELL
       A          B          C         ALPHA      BETA     GAMMA        VOLUME
    4.76020    4.76020  500.00000     90.0000   90.0000  120.0000      19.62371

 COORDINATES OF THE EQUIVALENT ATOMS (X AND Y IN FRACT. UNITS, Z IN ANGSTROMS)

 N. ATOM EQUIV AT. N.          X                  Y                  Z

   1   1   1   13 AL    0.00000000000E+00  0.00000000000E+00  1.92090000000E+00
   2   1   2   13 AL    0.00000000000E+00  0.00000000000E+00 -1.92090000000E+00

   3   2   1    8 O     3.33333333000E-01 -2.70930000000E-02  1.08280000000E+00
   4   2   2    8 O     2.70930000000E-02  3.60426333000E-01  1.08280000000E+00
   5   2   3    8 O    -3.60426333000E-01 -3.33333333000E-01  1.08280000000E+00
   6   2   4    8 O    -3.33333333000E-01  2.70930000000E-02 -1.08280000000E+00
   7   2   5    8 O    -2.70930000000E-02 -3.60426333000E-01 -1.08280000000E+00
   8   2   6    8 O     3.60426333000E-01  3.33333333000E-01 -1.08280000000E+00

   9   3   1   13 AL   -3.33333333000E-01  3.33333333000E-01  2.44600000000E-01
  10   3   2   13 AL    3.33333333000E-01 -3.33333333000E-01 -2.44600000000E-01

To model the surface (110), analysis of the layered structure orthogonal to the [110] direction:

TEST27 - CORUNDUM - SLAB (3D-->2D) - SEE TEST07
CRYSTAL
0 0 0
167
4.7602  12.9933
2
8     0.30624   0.   0.25
13    0.        0.   0.35216
REDEFINE
1 1 0
TESTGEOM
END


gives:


ATOMS CLASSIFIED ACCORDING TO THE Z COORDINATE :
  LAYER   1 Z=   1.6512; LABEL AT.NO.              X,Y,Z (ANG.)
                           1     8       -0.85440     0.67444     1.65122
                           2     8       -3.91909    -4.15525     1.65122
  LAYER   2 Z=   1.4578; LABEL AT.NO.              X,Y,Z (ANG.)
                           3     8       -2.38675    -1.74041     1.45776
  LAYER   3 Z=   0.9223; LABEL AT.NO.              X,Y,Z (ANG.)
                           6     8        0.17800    -1.74041     0.92234
  LAYER   4 Z=   0.7289; LABEL AT.NO.              X,Y,Z (ANG.)
                           4     8       -1.35435    -4.15525     0.72888
                           5     8        1.71034     0.67444     0.72888
  LAYER   5 Z=   0.0000; LABEL AT.NO.              X,Y,Z (ANG.)
                           7    13       -1.26595    -2.45160     0.00000
                           8    13        1.62194    -1.02921     0.00000
                           9    13        1.26595     2.45160     0.00000
                          10    13       -1.62194     1.02921     0.00000

Exercise: Model the (110) surface of corundum (solution1, solution2).

 

Computing the surface energy

If E(n) is the energy of an n-layer slab then our estimate of the energy per unit area required to form the surface from the bulk crystal based on an n-layer slab (Ens) is given by:

              Ens = [E(n) - nEbulk] / 2A

where A is the area of the primitive surface unit cell (8.862 Å2 as reported in the output) and Ebulk is the energy of a single layers worth of bulk material – in this case the energy of the bulk primitive unit cell which can be obtained by running mgobulk (extension d12): -274.6641058 Hartree. 

Thus, E3s = 3.367 mHa/ Å2 = 1.4679 J/m-2 
Note: 1 Ha/ Å2 = 435.9748 J/m-2. 

As more layers are used in the calculation Ens will converge to the surface energy if numerical errors in the calculation of E(n) and Ebulk cancel exactly. As CRYSTAL constructs the Hamiltonian matrix in direct space consistently in the SLAB and CRYSTAL calculations the main source of numerical differences between E(n) and Ebulk is due to k-space sampling.

Exercise: Compute E(n) and thus Ens for n=3, 4, 5, 6, 7 … layers and confirm that for MgO(100) the surface energy converges to  1.4679 Jm-2.

MgO is a wide band gap insulator and the convergence of Ebulk and E(n) with respect to k-space sampling is excellent .
Ens converges rapidly and the surface energy is well defined. In metals or small band gap semiconductors this may not be the case. A more numerically stable definition of the surface energy is then:


              Ens = {E(n) - n [E(n) - E(n-1)]} / 2A

In this expression Ebulk has been replaced by E(n)-E(n-1). If you consider each additional layer in the slab to be added as the central layer it is clear that E(n)-E(n-1) should converge to the energy of a single layer in the bulk crystal. However, in this case the bulk energy has effectively been computed under exactly the same numerical conditions as the surface layers leading to a cancellation of any systematic error.

Exercise: Plot (E(n)-E(n-1))  as a function of n –confirm that it converges to Ebulk 

Note: E(n) is tabulated E_n.dat

# MgO(100), Tols 7 7 7 7 14, Basis 861/851 Eb=--2.746641058E+02, confac = 435.9748

# nlayer        E(n)           Es(n)          Es(n)      E(n)-E(n-1)     Es(n)_method2
             hartree           Ha/A^2        J/m^2        hartree            J/m^2
    1     -274.60487715     0.00334172     1.45690559   -274.60487715     0.00000000
    2     -549.26850201     0.00336885     1.46873559   -274.66362486     1.44507559
    3     -823.93264035     0.00336702     1.46793476   -274.66413835     1.47033723
    4    -1098.59674765     0.00336693     1.46789765   -274.66410730     1.46804610
    5    -1373.26085590     0.00336679     1.46783731   -274.66410824     1.46813904
    6    -1647.92496516     0.00336660     1.46775187   -274.66410926     1.46826449
    7    -1922.58907607     0.00336631     1.46762597   -274.66411091     1.46850727
    8    -2197.25318892     0.00336591     1.46745235   -274.66411285     1.46884131
    9    -2471.91730400     0.00336539     1.46722393   -274.66411508     1.46927975
   10    -2746.58142159     0.00336472     1.46693352   -274.66411760     1.46983763

 

Properties of the surface

Many of the interesting properties of surfaces (optical, chemical, electronic, magnetic…) are related to the change in electronic structure which occurs in the near surface region. CRYSTAL provides a variety of tools for analyzing the electronic structure. These will be briefly reviewed here and discussed more extensively in the second tutorial on surfaces.

Charge distributions (Mulliken, charge density plots)

The Mulliken atomic charges are  printed at each SCF cycle, and on top of  properties module output.  To obtain full Mulliken population analysis insert the keyword PPAN in SCF input block or in properties input.

Mulliken charges for the 6-layer slab (keyword PPAN):

 

 MULLIKEN POPULATION ANALYSIS - NO. OF ELECTRONS  120.000000
  
  ATOM    Z CHARGE  SHELL POPULATION

   1 MG  12 10.038  2.000  7.986  0.052  
   2 O    8  9.957  2.007  4.597  3.353
   3 MG  12 10.025  2.000  7.985  0.039  
   4 O    8  9.981  2.007  4.617  3.357
   5 MG  12 10.021  2.000  7.986  0.035  
   6 O    8  9.979  2.007  4.617  3.355
   7 MG  12 10.021  2.000  7.986  0.035  
   8 O    8  9.979  2.007  4.617  3.355
   9 MG  12 10.025  2.000  7.985  0.039  
  10 O    8  9.981  2.007  4.617  3.357
  11 MG  12 10.038  2.000  7.986  0.052  
  12 O    8  9.957  2.007  4.597  3.353


 OVERLAP POPULATION CONDENSED TO ATOMS FOR FIRST   6 NEIGHBORS

 ATOM A   1 MG   ATOM B     CELL      R(AB)/AU R(AB)/ANG    OVPOP(AB)
                  2 O    (  0  0  0)     3.978     2.105    -0.006
                  1 MG   ( -1  0  0)     5.626     2.977     0.000
                  4 O    (  1  0  0)     6.890     3.646     0.000
                  1 MG   ( -1 -1  0)     7.956     4.210     0.000
                  2 O    ( -1  0  0)     8.895     4.707     0.000
                  3 MG   ( -1  0  0)     9.744     5.156     0.000

Bulk
MgO electronic charges: Mg 10.021 electrons ; O 9.979 electrons.

Measured in this way the bulk crystal is highly ionic, very nearly Mg2+O2-. At the surface the ionicity is very slightly reduced – by about 0.02|e|. The perturbation in the charges is rapidly screened and they  converge rapidly reaching bulk values at the third layer. 

Note: the slab is symmetric – layer 6 is identical to layer 1. It is often more convenient to study systems with symmetric slabs – even when surface adsorption is being studied.

ECHG computes the charge density in a plane orthogonal to the surface. Run properties with the following input (3 layer slab):


ECHG        keyword to compute charge density                  A----------D
0           order of the charge density derivative             |          |
65          number of points along the B-A segment             |          |
ATOMS       point A, B, C at atomic sites                      |          |
1 0 0 0     point A: atom 1 (Mg top layer) in cell (0 0 0)     B----------C
5 0 0 0     point B: atom 5 (Mg bottom layer) in cell (0 0 0)
6 0 0 0     point C: atom 6 (O bottom layer) in cell (0 0 0)   Mg  
MARGINS                                                        |
4. 4. 6. 6. margins along AB, CD, AD, BC                       |
END         end of definition of the map                       | 
END         end of properties input                            Mg---------O


Reference to the atom coordinates in a three layer slab

 

 *******************************************************************************
 *      ATOM          X(ANGSTROM)         Y(ANGSTROM)         Z(ANGSTROM)
 *******************************************************************************
   1    12 MG   -3.310539949135E-16  0.000000000000E+00  2.108500000000E+00
   2     8 O     1.490934648132E+00  1.490934648132E+00  2.108500000000E+00
   3    12 MG    1.490934648132E+00  1.490934648132E+00  0.000000000000E+00
   4     8 O     1.441895329408E-16  1.026461361182E-16  0.000000000000E+00
   5    12 MG    3.310069232598E-16  0.000000000000E+00 -2.108500000000E+00
   6     8 O     1.490934648132E+00  1.490934648132E+00 -2.108500000000E+00

The script runprop09 executes the program to draw contour lines, and visualize the postscript file generated.

Electrostatic potential

POTM computes the electrostatic potential in a plane orthogonal to the surface with the following input:.

POLI              multipole calculation
4 0 0             maximum order of the poles
POTM              electrostatic potential 
0 5               see CRYSTAL03 User's Manual
65                number of points along the B-A segment              
ATOMS             point A, B, C at atomic sites                       
1 0 0 0           point A: atom 1 (Mg top layer) in cell (0 0 0)     
5 0 0 0           point B: atom 5 (Mg bottom layer) in cell (0 0 0)
6 0 0 0           point C: atom 6 (O bottom layer) in cell (0 0 0)    
MARGINS                                                               
4. 4. 6. 6.       margins along AB, CD, AD, BC                       
END                                                                  
END               end of properties input

Electrostatic potential is positive on Mg ion and negative on Oxygen.

Charge density of 3-layer slab of MgO

Electrostatic potential of 3-layer slab of MgO

 

Surface reconstructions and supercell terminology

In general the atoms near the surface will move slightly away from their bulk positions; this is referred to as surface relaxation. In some systems, surfaces have a tendency to reconstruct – that is the periodicity of the surface layer changes from that implied by the bulk termination. 

If the primitive cell of the surface is defined by the lattice vectors a, b then a reconstruction introducing  a new periodicity involving two cells in a and three cells in b is called a (2x3) reconstruction. 
CRYSTAL provides a facility for creating supercells which can be used to describe any reconstruction which is commensurate with the underlying lattice. The new surface lattice vectors (a'b') are generated from the primitive ones by introducing a SUPERCEL matrix,



Our (2x3) reconstruction can then be entered at the end of the geometry section as:

SUPERCEL   

Transform cell to have two lattice vectors in the chosen surface plane

2    0   

S11 and S21

0    3   

S12 and S22


Exercise: Figure out the supercell which corresponds to the bulk crystallographic unit cell terminating at the surface – as shown in the figure below. How does the energy of a three layer slab for this cell compare to that of the primitive surface unit cell ?


 

A supercell of the MgO(100) surface.

 

Vicinal Surfaces – Modeling Steps and Kinks


Real surfaces are rarely atomically flat. Typically a surface will be covered in plateaux, with edges and corners. Explicit modeling of the low co-ordination sites associated with these structures involves very large unit cells and expensive calculations. 

A useful trick is to use vicinal surfaces (those cut at an angle slightly different to that of the low energy surface) to model steps and kinks. For example the figure shows a few unit cells of the (501) surface of MgO which contains a sequence of steps in the [100] direction.



The input to generate this slab is:

SLABCUT     Transform cell to have two lattice vectors in the chosen surface plane
5   0   1   Miller indices of the surface plane
1      10  Select Layer 1 as the terminating layer and request a slab of 10 layers

The primitive unit cell of this surface is 3.714 Angstrom thick and contains 20 atoms in 10 layers at distance of 0.413 Angstrom.  The system has only 4 symmetry operators.

 NEW FUNDAMENTAL DIRECT-LATTICE VECTORS B1,B2
             X         Y         Z
 B1       4.2100    0.0000    0.0000
 B2       0.0000   10.7334    0.0000

 AREA OF THE 2D CELL        45.188

 COORDINATES OF THE ATOMS BELONGING TO THE SLAB
 LAB AT.N.     CARTESIAN (ANG)          FRAC (2D)     Z(ANG)
   1  12  1.05250 -1.44489  1.85771   0.250 -0.135   1.85771    layer # 1
   2   8 -1.05250 -1.44489  1.85771  -0.250 -0.135   1.85771
   3  12 -1.05250 -3.50901  1.44489  -0.250 -0.327   1.44489    layer # 2
   4   8  1.05250 -3.50901  1.44489   0.250 -0.327   1.44489
   5  12  1.05250  5.16031  1.03206   0.250  0.481   1.03206    layer # 3
   6   8 -1.05250  5.16031  1.03206  -0.250  0.481   1.03206
   7  12 -1.05250  3.09618  0.61924  -0.250  0.288   0.61924    layer # 4
   8   8  1.05250  3.09618  0.61924   0.250  0.288   0.61924
   9  12  1.05250  1.03206  0.20641   0.250  0.096   0.20641    layer # 5
  10   8 -1.05250  1.03206  0.20641  -0.250  0.096   0.20641
  11  12 -1.05250 -1.03206 -0.20641  -0.250 -0.096  -0.20641    layer # 6
  12   8  1.05250 -1.03206 -0.20641   0.250 -0.096  -0.20641
  13  12  1.05250 -3.09618 -0.61924   0.250 -0.288  -0.61924    layer # 7
  14   8 -1.05250 -3.09618 -0.61924  -0.250 -0.288  -0.61924
  15  12 -1.05250 -5.16031 -1.03206  -0.250 -0.481  -1.03206    layer # 8
  16   8  1.05250 -5.16031 -1.03206   0.250 -0.481  -1.03206
  17  12  1.05250  3.50901 -1.44489   0.250  0.327  -1.44489    layer # 9
  18   8 -1.05250  3.50901 -1.44489  -0.250  0.327  -1.44489
  19  12 -1.05250  1.44489 -1.85771  -0.250  0.135  -1.85771    layer #10
  20   8  1.05250  1.44489 -1.85771   0.250  0.135  -1.85771

CHARGE NORMALIZATION FACTOR   1.00000000; TOTAL ATOMIC CHARGES: (last SCF cycle)

  10.0620154   9.9290463  10.0472099   9.9609569  10.0412329   9.9588219
  10.0411186   9.9598752  10.0454597   9.9542632  10.0454597   9.9542632
  10.0411186   9.9598752  10.0412329   9.9588219  10.0472099   9.9609569
  10.0620154   9.9290463

It takes a little longer to run in CRYSTAL (see 501_10.d12 and 501_10.out).

Exercise: Run the (501) 10-layer slab in CRYSTAL. How do the Mulliken charges on the low coordinate sites compare to those on the flat surface ? It should be clear from the charges that this slab is not thick enough to converge the electronic structure or the surface energy.

Exercise: If you have the patience try running a 20 layer version of this slab. Check the disk space necessary to write mono- and bi- electron integrals by inserting TESTRUN in the third input block. 
If the space occupied by the integrals exceed 2Gb, check if extended files are allowed by your installation. Otherwise split the integrals in several files using the keyword BIESPLIT (third input block). If the space is not available, run in SCF direct mode by inserting the keyword SCFDIR in the third input block (remove TESTRUN). This is a reasonable model of the [100] step and a starting point for considering the adsorption of small molecules on the stepped surface.

Exercise:
Modify the number of layers of the slab MgO (001) and see how the number of symmetry operators change.

Exercise:
Define a 3-layers slab parallel to the MgO (110) crystallographic plane.

Exercise:
Define a MgO surface with higher indices, e.g. (103). Such planes can be used to simulate a step on the surface. Note: play with the number of atomic layers in order to define a suitable surface model.

 

Acknowledgements

We thank Professor N.M. Harrison, who prepared the first version of this tutorial for MSSC2000.


References.

    J.W. He, C.A. Estrada, J.S. Corneille, M.C. Wu and D.W. Goodman, Surf. Sci., 261 (1992) 164.
    R. Wichtendhal, M. Rodriguez-Rodrigo, U. Härtel, H. Kuhlenbeck and H.J. Freund, Phys. Stat. Sol. (a), 173 (1999) 93.

    A. Zecchina, D. Scarano, S. Bordiga, G. Ricchiardi, G. Spoto, F. Geobaldo, Catalysis Today, 27 (1996) 403.

    G. Pacchioni, T. Minerva and P.S. Bagus, Surf. Sci., 275 (1992) 450.

    M.A. Nygren, L.G.M. Petterson, Z. Barandiaran, L. Seijo, J. Chem. Phys., 100 (1994) 2010.

    J.A. Meijas, A.M. Marquez, J. Fernandez-Sanz, M. Fernandez-Garcia, J.M. Ricart, C. Sousa and F. Illas, Surf. Sci., 327 (1995) 59.

    A.M. Ferrari and G. Pacchioni, Int. J. Quantum Chem., 58 (1996) 241.

    G. Pacchioni, A.M. Ferrari, A.M. Marquez and F. Illas, Int. J. Quantum Chem., 18 (1997) 617.

    J. L. Pascual, L.G.M. Petterson, Chem. Phys. Lett., 270 (1997) 351.

    J.A. Snyder, D. R. Alfonso, J.E. Jaffe, Zijing Lin, A.C. Hess and M. Gutowski, J. Phys. Chem. B, 104 (2000) 4717.

    S.F. Boys and F. Bernardi, Mol. Phys., 19 (1970) 553.

    A. Damin, R. Dovesi, A. Zecchina, P. Ugliengo, Surf. Sci., 479 (2001) 255.



CRYSTAL home page CRYSTAL home page