Match

• MATCH
• I. Introduction
  • I.1. The physics
  • I.2. Software package match
  • I.3. Program match
  • I.4. Program incom
  • I.5. Program outgo
  • I.6. Program domatch
• II. Input file
  • II.1. Structure of the file
  • II.2. Program calcleis
• III. Configuration file
  • III.1. Structure of the file
  • III.2. Program makelattice
• IV. Angular tables file
  • IV.1. Structure of the file
  • IV.2. Program maketable
  • IV.3. Program showtable
  • IV.4. Program table2asc
  • IV.5. Program table2bin
• V. Output file
  • V.1. Structure of the file
  • V.2. Program showscan
  • V.3. Program sumscan
• VI. Specular reflection
  • VI.1. Program spec
  • VI.2. Program showspec
  • VI.3. Program makesurface

MATCH

I. Introduction

I.1. The physics

1. The difference between the energies of projectile and scattered or recoiled particle has to correspond to the kinematics of the collision. The tolerance is determined by the chosen width of the energy bins, the kinematic factor of the collision, and the energy loss between this collision and the detector.
   
   
2. The incoming and outgoing trajectories cross each other within 0.05 angstrom. As a matter of fact, the trajectories should intersect. This means that matching can be accomplished by shifting both incoming and outgoing trajectory (i.e. all atoms at wich incoming and outgoing particles have been scattered) by half the distance between the two trajectories. This approach is assumed to be valid as long as the shift is appreciably smaller than the thermal vibration amplitude of the atoms.

1. The probability of having the atom at the proper position that would have caused the scattering/recoiling event. This probability is calculated from the thermal vibration distributions. The position of the atom considered for the scattering/recoiling event is derived from the intersection point of the incoming and the outgoing trajectories, and the impact parameter correponding to the scattering/recoiling angle.
   
   
2. The differential scattering or recoiling cross-section of the collision.

    See: M.H. Langelaar, M. Breeman, A.V. Mijiritskii and D.O. Boerma
         Nucl. Instr. and Meth. in Phys. Res. B132 (1997) 578-584

I.2. Software package match

I.3. Program match

• Input file
  
  This ASCII file contains global information that is partly passed to the programs incom, outgo and domatch. Some input parameters, like scattering angle and energy of the outgoing particle, are conveniently calculated with program calcleis. The file must be specified on the command line.
  
  
• Configuration file
  
  This ASCII file contains information about the slab size, lattice structure, constituing elements, atom positions, and thermal lattice vibrations. The file may be created with program makelattice. The file must be specified on the command line.
  
  
• Angular tables file
  
  These tables contain values for different scattering and recoiling angles, scattering cross-sections, and electronic energy loss. The tables are created with program maketable. The tables can be viewed with program showtable. Watch for irregularities in the colour patterns as they may indicate computational errors ! Note that the angular tables contain unformatted information. Use programs table2asc and table2bin to port the tables to other computer systems.

• Output file
  
  This ASCII file contains the calculated detector yield. There is a separate file for each energy bin. They may be viewed with program showscan. The results for the lowest energy bin can be displayed while match is running.
  
  
• <stout>
  
  Error messages are printed to standard output.

• incom
  
  Create input for program incom, start this program and check the status every 5 seconds. The output files have the extention .<sym1>.<sym2>. In the scattering case <sym1> is the incoming projectile, whereas in the recoiling case <sym1> is either the incoming projectile or a recoiled target atom; in all cases <sym2> is one of the target atoms. Next, for each element to be detected and for each energy bin:
  
  
• outgo
  
  Create input for program outgo, start this program and check the status every 5 seconds. The output files have the extention .<sym1>.<sym2>. In the scattering case <sym1> is the scattered projectile and <sym2> is a target atom, whereas in the recoiling case both <sym1> and <sym2> are recoiled target atoms. For each combination of incom and outgo output files:
  
  
• domatch
  
  Create input for program domatch, start this program and check the status every 5 seconds. Update the detector yield. Update the display window if outputmode is TRUE. Note that only the results for the last detected element and the lowest energy bin are displayed while the program is running.

I.4. Program incom

    LEIS-ELEMENT   Ar              # primary projectile
    MASS           39.98
    #SURF_ELEMENTS 2               # number of elements in the slab
    HOEKTABEL_1    tabs/ArFe95.ang # name of the angular tables file and
    MAXIMPACT_1    0.6             # the maximum impact parameter; there
    HOEKTABEL_2    tabs/ArN80.ang  # is a table for each combination of
    MAXIMPACT_2    0.4             # projectile and slab element, and for
    HOEKTABEL_3    tabs/FeFe90.ang # each combination of two slab elements
    MAXIMPACT_3    0.7
    HOEKTABEL_4    tabs/FeN80.ang
    MAXIMPACT_4    0.4
    HOEKTABEL_5    tabs/NFe80.ang
    MAXIMPACT_5    0.6
    HOEKTABEL_6    tabs/NN98.ang
    MAXIMPACT_6    0.4
    MIN_ENERGY_0   2.30            # minimum energy till which element i
    MIN_ENERGY_1   2.30            # is followed in the calculation
    MIN_ENERGY_2   2.30
    START_ENERGY   5.00            # energy of the projectile
    INCOMING_ANGLE 168.0           # incoming angle
    OUTGOING_ANGLE 136.0           # outgoing angle
    TILT           0.0             # tilt angle
    AZIMUTH        -74             # azimuthal slab orientation
    #FROZ.LAT      100             # number of frozen lattices
    DOSE_PER_FROZ. 100             # dose per frozen lattice
    OUTPUT_MODE    1               # 0 = no messages, no graphics

• layertab[Natoms]
  
  index of each layer in the slab; 0 = surface layer
  
  
• surftab[Natoms]
  
  index of each surface atom
  
  
• lattice[Natoms]
  
  actual position (in angstrom) of each atom in the slab
  
  
• equilat[Natoms]
  
  equilibrium position (in angstrom) of each atom in the slab
  
  
• Nneighbour[Natoms]
  
  number of neighbours of each atom in the slab
  
  
• neighbour[Natoms][MAXNEIGHBOUR]
  
  indices of the neighbours of each atom in the slab
  
  
• maximpact[MAXANGTAB]
  
  maximum tabulated impact parameter in each angular table
  
  
• ELSTOP[MAXANGTAB][maxs][maxe]
  
  electronic energy loss as a function of the impact parameter and the energy of the incoming atom, for each angular table
  
  
• KTAB[MAXANGTAB][maxs][maxe]
  
  
• ANGTAB[MAXANGTAB][maxs][maxe]
  
  the kinematic factor used to calculate the outgoing energy from the incoming energy and the impact parameter, and the scattering angle in laboratory system, for each angular table; the following formula are used:

      ksi = (pi-ph)/2
      Er/E = [ 4mu/(1+mu)**2 ]*cos(ksi)**2
      Es/E = 1 - Er/E
      tg(th) = sin(ph) / [ 1/mu + cos(ph) ]
  
      see: D.O. Boerma
           Lecture Notes on Atomic Collisions (1995), formula (2.11),
           (2.14) and (2.22)
  

1. Create a frozen lattice by adding random displacements to the equilibrium lattice positions.
   
   
2. Construct a sorted neighbour list for each of the atoms in the slab.
   
   
3. Emit a primary projectile into the incoming direction, at a random position on the slab and at a height of 2.2 angstrom above the surface.
   
   
4. Store all relevant information in local variables, which will be put on the stack on a recursive call of the scattering function.
   
   
5. Select the next scattering partner from either surftab (first time) or the neighbour list of the current scattering partner. Calculate energy, position and velocity vector of the primary projectile after the scattering event, and store the information in a trajectory structure.
   
   
6. Calculate energy, positon and velocity vector for the recoiling atom, and restart the simulation at point 4.
   
   
7. Restart the simulation at point 3 while the number of emitted projectiles is less than the dose specified in the input file.
   
   
8. Restart the simulation at point 1 while the number of frozen lattices is smaller than the value specified in the input file.

1. The energy of the projectile has dropped below the minimum energy specified in the input file.
   
   
2. The primary atom has moved to a depth below the bottom of the slab.
   
   
3. The primary atom is scattered away from the slab, to a height of 2.0 angstrom above the surface.

I.5. Program outgo

    DETEC-ELEMENT   N              # outgoing particle
    DETEC-ELEMENTNO 2
    MASS            14.00
    #SURF_ELEMENTS  2              # number of elements in the slab
    HOEKTABEL_4     tabs/NFe80.ang # name of the angular tables file and
    MAXIMPACT_4     0.6            # the maximum impact parameter; there
    HOEKTABEL_5     tabs/NN98.ang  # is a table for each combination of
    MAXIMPACT_5     0.4            # outgoing particle and slab element
    MAX_ENERGY      2.30           # maximum energy of outgoing particle
    DET_ENERGY      1.40           # current energy bin
    INCOMING_ANGLE  168.0          # incoming angle
    OUTGOING_ANGLE  136.0          # outgoing angle
    TILT            0.0            # tilt angle
    AZIMUTH         -74            # azimuthal slab orientation
    #FROZ.LAT       100            # number of frozen lattices
    DOSE_PER_FROZ.  100            # dose per frozen lattice
    OUTPUT_MODE     1              # 0 = no messages, no graphics

• layertab[Natoms]
• surftab[Natoms]
• lattice[Natoms]
• equilat[Natoms]
• Nneighbour[Natoms]
• neighbour[Natoms][MAXNEIGHBOUR]
• maximpact[MAXANGTAB]
• ELSTOP[MAXANGTAB][maxs][maxe]
• KTAB[MAXANGTAB][maxs][maxe]
• ANGTAB[MAXANGTAB][maxs][maxe]

The outline of program incom is as follows:

1. Create a frozen lattice by adding random displacements to the equilibrium lattice positions.
   
   
2. Construct a sorted neighbour list for each of the atoms in the slab.
   
   
3. Emit a primary projectile in a rondom direction within a cone with opening (half) angle of 0.25 degrees about the incoming direction, at a random position on the slab, and at a height of 2.2 angstrom above the surface.
   
   
4. Select the next scattering partner from either surftab (first time) or the neighbour list of the current scattering partner. Calculate energy, position and velocity vector of the primary projectile after the scattering event, and store the information in a trajectory structure.
   
   
5. Restart the simulation at point 3 while the number of emitted projectiles is less than the dose specified in the input file.
   
   
6. Restart the simulation at point 1 while the number of frozen lattices is smaller than the value specified in the input file.

1. The energy of the projectile is more than the maximum energy specified in the input file.
   
   
2. The primary atom has moved to a depth below the bottom of the slab.
   
   
3. The primary atom is scattered away from the slab, to a height of 2.0 angstrom above the surface.

I.6. Program domatch

    INCOMFILE      dum/incom.output.N.N # location of output from INCOM
    OUTGOFILE      dum/outgo.output.N.N # location of output from OUTGO
    FSTMATCHLAYER  1                    # first and last slab layer used
    LASTMATCHLAYER 4                    # in the calculation
    MU             0.350175             # mass ratio of target and incoming atom
    MAIN_SCAT_NO   2                    # index of detected element (scattering)
                                        # or scattering element(recoil)
    SCAT_RECOIL    1                    # 0 = scattering, 1 = recoiling
    DIST_TOL       0.050000             # distance tolerance in angstrom
    E_TOL          0.035714             # relative energy tolerance
    ANGTAB         tabs/NN98.ang        # location of angtab file
    CONFIGFILE     FeN100.cfg           # location of config file
    OUTPUT_MODE    1                    # 0 = no messages, no graphics
    SIZEX          3.795000             # slab size in angstrom
    SIZEY          3.795000
    5 1.897500 1.897500 1.897500 2.000000 0.120000 0.120000 0.120000
    3 1.897500 1.897500 5.692500 2.000000 0.120000 0.120000 0.120000
    1 1.897500 1.897500 9.487500 2.000000 0.120000 0.120000 0.150000
    ^ <------------------------> <------> <------------------------>
    |       atom position      MAIN_SCAT_NO   vibration amplitude
    |       (in angstrom)                     (in lattice units)
    number of layer

• atomtab[MAXATOMS]
  
  maps each atom in the slab on the equivalent atom within the unit cell; the index is set to -1 if the element differs from the one for which the yield is currently calculated
  
  
• impartab[maxangle*2][maxe]
• crosstab[maxangle*2][maxe]
  
  impact parameter and cross-section in the laboratory system, as a function of the scattering angle in the laboratory system and the energy of the incoming atom
  
  
• recimpartab[maxangle][maxe]
• reccrosstab[maxangle][maxe]
  
  impact parameter and cross-section in the laboratory system, as a function of the recoiling angle in the laboratory system and the energy of the incoming atom
  
  
• Ktab[MAXK]
• Ktab2[MAXK]
  
  kinematic factors used to calculate the outgoing energy from the incoming energy and the scattering/recoiling angle in the laboratory system; the following formula are used:

      Es/Ep = [ cos(th) +- sqrt( mu**2 - sin(th)**2 ) ]**2 / [ 1+mu ]**2
        for mu >= 1 only the + sign holds
        for mu < 1 and th < asin(mu) both + sign and - sign
        should be considered
      Er/Ep = [ 4mu/(1+mu)**2 ]*cos(ksi)**2
  
      see: D.O. Boerma
           Lecture Notes on Atomic Collisions (1995), formula (2.16),
           (2.18) and (2.22)
  

  data for the + sign are stored in Ktab, those for the - sign in Ktab2

1. The outgoing particle is of the element for which the yield is currently calculated.
   
   
2. The incoming particle is of the element for which the yield is currently calculated.
   
   
3. The scattering centres for incoming and outgoing particles both map on the same atom in the unit cell.
   
   
4. The scattering angle is smaller than pi if mu > 1, or smaller than asin(mu) if mu < 1.
   
   
5. The recoiling angle is smaller than pi/2.
   
   
6. The outgoing energy calculated from the energy of the incoming trajectory, the scattering angle and the kinematic factor doesn't deviate from the energy of the outgoing trajectory by more than the energy tolerance.

1. The distance between incoming and outgoing trajectories is calculated; the matching procedure is terminated if the distance exceeds the distance tolerance, else the two trajectories are shifted over half the distance to make them intersect.
   
   
2. The intersection point is calculated; the matching procedure is terminated if the intersection point is more than 0.3 angstrom outside the layer corresponding to the incoming trajectory, or if the intersection point is more than 1.0 angstrom away from the equilibrium position of the scattering atom corresponding to the outgoing trajectory.
   
   
3. The impact parameter b for the scattering/recoiling event is obtained from table impartab or recimpartab.
   
   
4. The actual location of the scattering centre is calculated using the property that in a central-force field the motion of a particle is symmetrical about the point of closest approach to the force centre:

                                               scattering
                               .                 .
                                .              .
                                 R...........Q
                                  .        .  .
                                   . ksi .     .
                                ksi .  . phi    .
       incoming .....................A...........P....
                                |     . ksi
                                |      .
                               b|       .
                                |        .           phi = pi-2ksi
                            ..............B          b = Rm sin(ksi)
                                           .           = Rm cos(phi/2)
                                            .
                                             .
                                              .
                                               S
                                                .
                                              recoiling
   
       A       position of incoming atom at closest approach
       B       position of target atom at closest approach
       b       impact parameter
       |A->B|  distance Rm of closest approach
       A       intersection of incoming and outgoing trajectories
       A->P    (unity) velocity vector of incoming atom
       A->Q    (unity) velocity vector of scattering atom
       B->S    (unity) velocity vector of recoiling atom
   

   Since |A->P| = |A->Q| = 1, we have angle(AQP) = angle(APQ) = ksi, and thus Q->P // R->A and R->Q // A->P. Therefore, both Vin-Vout(scattering) = R->A and Vout(recoil) = B->S point into the direction of the vector A->B. The distance between A and B is equal to Rm = b/cos(alpha/2) in the scattering case, and Rm/sin(alpha) in the recoiling case, where alpha is the angle between incoming and outgoing trajectories.

       see: Jerry B. Marion
            Classical Dynamics, figure 11-9
            (Academic Press, 1965)
            D.O. Boerma
            Lecture Notes on Atomic Collisions, figure 3.3
   

5. The probability P of finding the scattering atom at the actual location of the scattering centre is calculated from the vibration amplitudes, assuming a gaussian probability distribution
   
   
6. The scattering/recoiling cross-section sigma of the scattering/recoiling event is obtained from table crosstab or reccrosstab.
   
   
7. Finally, the detector yield is incremented by an amount equal to the product of probability P and cross-section sigma.

II. Input file

II.1. Structure of the file

    ** PARTICLES
    CONFIG_FILE    FeN100.cgf      # name of configuration file
    PROJECTIEL     Ar              # name of element 0 (projectile)
    MASS           39.98           # mass of element 0
    #ELEMENTEN     2               # number of elements in the slab
    ELEMENT_1      Fe              # name of element 1
    MASS_1         56.94           # mass of element 1
    ELEMENT_2      N               # name of element 2
    MASS_2         14.00           # mass of element 2
    #DETECELEM     1               # number of detected elements
    DETECTIENO_0   0               # index of the detected element
    HOEKTABEL_0_1  ArFe95.ang      # name of the angular tables file and
    MAXIMPACT_0_1  0.6             # the maximum impact parameter; there
    HOEKTABEL_0_2  ArN80.ang       # is a table for each combination of
    MAXIMPACT_0_2  0.4             # projectile and slab element, and for
    HOEKTABEL_1_1  FeFe90.ang      # each combination of two slab elements
    MAXIMPACT_1_1  0.7
    HOEKTABEL_1_2  FeN80.ang
    MAXIMPACT_1_2  0.4
    HOEKTABEL_2_1  NFe80.ang
    MAXIMPACT_2_1  0.6
    HOEKTABEL_2_2  NN98.ang
    MAXIMPACT_2_2  0.4

    ** ENERGIES
    PRIM_ENERGY    5.00            # energy of the projectile
    MINE_PROJ_0    2.30            # mimimum energy till which element i
    MINE_PROJ_1    2.30            # is followed in the calculation
    MINE_PROJ_2    2.30
    MINE_DETECT_0  2.40            # lower and upper energy windows,
    MAXE_DETECT_0  3.60            # and maximum outgoing energy, for
    MAXE_BACKW_0   3.70            # detected element i

    ** GEOMETRY
    INCOMING_ANGLE 168.0           # incoming angle
    OUTGOING_ANGLE 136.0           # outgoing angle
    TILT           0.0             # tilt angle
    AZIMUTH        0.0             # azimuthal slab orientation
    SCANTYPE       1               # 0 = polar, 1 = azimuthal
    FIRST_ANGLE    -90.0           # begin of scan
    LAST_ANGLE     90.0            # end of scan
    INCR_ANGLE     2.0             # angular step
    #FROZEN_LATT   100             # number of frozen lattices
    DOSE_PER_FROZ  100             # dose per frozen lattice

    ** MATCHING
    FSTMATCHLAYER  1               # first and last slab layer used
    LASTMATCHLAYER 4               # in the calculation
    DIST_TOL       0.05            # tolerance in distance
    E_TOL          0.10            # tolerance in energy 

    ** FILES
    INCOM_DIR      /usr/local/bin/ # location of program incom
    OUTGO_DIR      /usr/local/bin/ # location of program outgo
    DOMATCH_DIR    /usr/local/bin/ # location of program domatch
    ANGTAB_DIR     tabs/           # location of angtab files
    INPUT_INCOM    dum/incom.in    # names of input and output files
    OUTPUT_INCOM   dum/incom.out   # of programs incom, outgo and
    INPUT_OUTGO    dum/outgo.in    # domatch
    OUTPUT_OUTGO   dum/outgo.out
    INPUT_DOMATCH  dum/domatch.in
    OUTPUT_DOMATCH dum/domatch.out

    ** VERBOSE
    OUTPUT_MODE    1               # 0 = no messages, no graphics

II.2. Program calcleis

    nele                            # number of listed elements
    sym[1]     m[1]       z[1]      # element symbol, mass (in amu)
    ...                             # and atomic number
    sym[nele]  m[nele]    z[nele]
    ndet                            # number of listed detectors
    name[1]    path[1]    ang[1]    # name, distance from target (in m)
    ...                             # and angle (in degree)
    name[ndet] path[ndet] ang[ndet]
    path0                           # distance (in m) between beam chopper
                                    # and target

• Input parameters
  
  The following input parameters can be set in the left part of the graphics menu:

      Projectile                     # pop-up menu
      Target                         # pop-up menu
      Energy (in keV)                # edit line
      Detection angle (in degree)    # edit line
      Impact parameter (in angstrom) # edit line
      Detector                       # pop-up menu
      Screening correction           # edit line
  

• Actions
  
  The following actions may be selected by clicking on the proper button in the right part of the graphics menu:
  • Angle
    
    The scattering and recoiling angles in the laboratory system are calculated for the specified combination of elements, energy and impact parameter. The scattering angle phi in the centre-of-mass system is found by calculating the scattering integral:

        phi = pi-2*b*Int[ f(r)dr ]
        1/f(r) = r*r*sqrt( 1 - (b/r)**2 - V(r)/E )
    

    where b is the impact parameter, V(r) is the interatomic potential, and E is the energy of the incoming particle in the centre-of-mass-system. The scattering angle theta and the recoiling angle ksi in the laboratory system follow from:

        tg(theta) = mu*sin(phi) / [ 1 + mu*cos(phi) ]
        cos(ksi) = sin(phi/2)
    
        mu = m2/m1
        m1   mass of the incoming atom (in amu)
        m2   mass of the target atom (in amu)
    

    Note that scattering is not possible if mu = 1 and phi = pi, because the incoming atom comes to a full stop and never will reach the detector.

        See: D.O. Boerma
             Lecture Notes on Atomic Collisions (1995), formula (2.11),
             (2.14) and (3.13)
    

  • Impact parameter

    The impact parameters for scattering and recoiling processes are calculated for the specified combination of elements, energy and detection angle. The following formula are used to calculate the scattering angle phi in the centre-of-mass system from the specified detection angles theta and ksi, respectively:

        mu*cos(phi) = -sin(theta)**2
                      +- cos(theta)*sqrt[ mu**2 - sin(theta)**2 ]
          for mu >= 1 only the + sign holds
          for mu < 1 both + sign and - sign should be considered
        cos(phi/2) = sin(ksi)
    

    Once phi is known, the impact parameter b is obtained by solving the equation

        CMSangle(b) = phi
    

    by successive approximation. Scattering is assumed only to be possible for an impact parameter in te range from 0.001 to 2.0 angstrom.

        See: D.O. Boerma
             Lecture Notes on Atomic Collisions (1995), formula (2.14)
             and (2.17)
    

  • Cross-section
    
    The scattering and recoiling cross-sections in the laboratory system are calculated for the specified combination of elements, energy and detection angle. First the scattering angle phi in the centre-of-mass system is calculated from the specified detection angles theta and ksi, respectively (see Impact parameter). Then the cross-sections are obtained from the formula:

        sigma(phi) = [ b/sin(phi) ]*| db/dphi |
        sigma(theta) = F(phi,mu)*sigma(phi)
        sigma(ksi) = G(phi)*sigma(phi)
        F(phi,mu) = | [ 1 + 2*cos(phi)/mu + 1/mu**2 ]*
                      mu/[ mu + cos(phi) ] |
        G(phi) = 4*sin(phi/2)
    
        mu = m2/m1
        m1      mass of the incoming atom (in amu)
        m2      mass of the target atom (in amu)
    

    Note that sigma(theta) has a singularity for cos(phi) = -mu.

        See: D.O. Boerma
             Lecture Notes on Atomic Collisions (1995), formula (3.23),
             (3.34) and (3.35)
    

  • Energy

    The energies of scattering and recoiling atoms are calculated for the specified combination of elements, energy and detection angle. The energy ene is calculated from the specified detection angles theta and ksi, respectively, using the formula:

        sqrt( ene/E ) = [ cos(theta) +- sqrt( mu**2 - sin(theta)**2 ] /
                        [ 1 + mu ]
            for mu >= 1 only the + sign holds
            for mu < 1 both + sign and - sign should be considered
        sqrt( ene/E ) =  2*sqrt(mu)*cos(ksi) / [ 1 + mu ]
    

    where E is the specified energy of the incoming atom. Note that scattering or recoiling is not possible if the energy is too small, because the atom never will reach the detector.

        See: D.O. Boerma
             Lecture Notes on Atomic Collisions (1995), formula (2.18)
             and (2.22)
    

  • Time_of_flight
    
    The flight time of incoming, scattering and recoiling atoms are calculated for the specified combination of elements, energy and detection angle. First the energy ene is calculated from the specified detection angles theta and ksi, respectively (see Energy). Then the propagation time t is obtained from:

        t = path / velocity
          = path*sqrt( 0.5*mass/ene )
    

    Multiply by sqrt(1.6605655e-27/1.6021892e-16) to change from kg/J to amu/keV, and multiply by 1.0e+9 to change from s to ns.

III. Configuration file

III.1. Structure of the file

    ORIENTATION     FeN(100)
    LATT_CONSTANT   3.795    # lattice constant in angstrom
    #ELEMENTS       2        # number of elements in slab
    Fe 56.94                 # first element and mass
    N  14.00                 # second element and mass
    4 4 4                    # atom rows on the slab
    2 2 4                    # unit cells on the slab
    2.000000 2.000000        # slab size / lattice constant
    40                       # number of atoms in the slab
    0.000000 0.000000 0.000000 1 0.120000 0.120000 0.120000 4
    1.000000 0.000000 0.000000 1 0.120000 0.120000 0.120000 4
    0.500000 0.500000 0.000000 1 0.120000 0.120000 0.120000 4
    1.500000 0.500000 0.000000 1 0.120000 0.120000 0.120000 4
    0.000000 1.000000 0.000000 1 0.120000 0.120000 0.120000 4
    1.000000 1.000000 0.000000 1 0.120000 0.120000 0.120000 4
    0.500000 1.500000 0.000000 1 0.120000 0.120000 0.120000 4
    1.500000 1.500000 0.000000 1 0.120000 0.120000 0.120000 4
    1.500000 0.500000 0.500000 2 0.120000 0.120000 0.120000 3
    ............................. ................. .........
    ............................. ................. .........
    0.000000 1.000000 1.000000 1 0.120000 0.120000 0.120000 2
    0.500000 0.500000 1.500000 2 0.120000 0.120000 0.150000 1
    1.500000 0.000000 1.500000 1 0.120000 0.120000 0.150000 1
    1.000000 0.500000 1.500000 1 0.120000 0.120000 0.150000 1
    1.500000 1.000000 1.500000 1 0.120000 0.120000 0.150000 1
    1.000000 1.500000 1.500000 1 0.120000 0.120000 0.150000 1
    0.500000 0.000000 1.500000 1 0.120000 0.120000 0.150000 1
    0.000000 1.500000 1.500000 1 0.120000 0.120000 0.150000 1
    1.500000 0.500000 1.500000 2 0.120000 0.120000 0.150000 1
    0.500000 1.500000 1.500000 2 0.120000 0.120000 0.150000 1
    0.500000 1.000000 1.500000 1 0.120000 0.120000 0.150000 1
    1.500000 1.500000 1.500000 2 0.120000 0.120000 0.150000 1
    0.000000 0.500000 1.500000 1 0.120000 0.120000 0.150000 1
    <------------------------> ^ <------------------------> ^
          atom position        |     vibration amplitude    |
                            element                        layer
                           0=vacancy

III.2. Program makelattice

• Input parameters
  
  The program reads the following data (this example creates the configuration file shown in the previous section):

      STRUCTURE      sc                   # lattice structure
      ORIENTATION    (100)                # slab orientation
      LATT_CONSTANT  3.795                # lattice constant in angstrom
      #SUBLATTICES   5                    # number of sublattices
      ROWS_ON_SLAB   2 2 2                # unit cells on the slab
      ROWS_PER_CELL  2 2 2                # atom rows per unit cell
      SUBLATTICE_0   Fe 56.94 0.0 0.0 0.0 # element symbol, element mass,
                                          # and shift of the sublattice
      SUBLATTICE_1   Fe 56.94 0.5 0.5 0.0
      SUBLATTICE_2   Fe 56.94 0.5 0.0 0.5
      SUBLATTICE_3   Fe 56.94 0.0 0.5 0.5
      SUBLATTICE_4   N  14.00 0.5 0.5 0.5
      THVIB_LAYER_0  0.12 0.12 0.15       # vibration amplitudes
      THVIB_LAYER_1  0.12 0.12 0.12
      THVIB_LAYER_2  0.12 0.12 0.12
      THVIB_LAYER_3  0.12 0.12 0.12
      OUTFILE        FeN100.cfg           # name of the output file
  

  The structure of the input file is checked. The input parameters must fulfil the following constraints:

      structure                     fcc, bcc, hcp, sc (= simple cubic)
      orientation                   (100), (110), (111), (0001)
                                    for hcp only 0001 is allowed, for
                                    others structures (0001) is not
                                    allowed                              
      maximum number of sublattices 10
      maximum number of rows        10
      rows on the slab              integer multiple of the number of
                                    rows per unit cell
      mass                          any value from 1.0 till 300.0
      shift                         any value from 0.0 till 1.0
      vibration amplitude           any value from 0.0 till 1.0
  

  
  
  
• Unit cell
  
  The unit cells used in the program are described by the following vectors (the latter are defined in the crystallographical unit cell):

                        x              y               z
      =========================================================
      fcc(100)    (1/2,1/2,0  )  (-1/2,1/2,0  )  (0,   0,  1/2)
      fcc(110)    (1/2,0  ,1/2)  (0   ,1  ,0  )  (-1/4,0  ,1/4)
      fcc(111)    (1/2,1/2,0  )  (-1/4,1/4,1/2)  (1/3,-1/3,1/3)
  
      bcc(100)    (1  ,0  ,0  )  (0   ,1  ,0  )  (0   ,0  ,1/2)
      bcc(110)    (1  ,0  ,0  )  (0   ,1/2,1/2)  (0  ,-1/2,1/2)
      bcc(111)    (1  ,1  ,0  )  (-1/2,1/2,1  )  (1/6,-1/6,1/6)
  
      hcp(0001)
  
      sc(100)     (1  ,0  ,0  )  (0   ,1  ,0  )  (0   ,0  ,1  )
      sc(110)     (1  ,0  ,0  )  (0   ,1  ,1  )  (0  ,-1/2,1/2)
      sc(111)     (1  ,1  ,0  )  (-1/2,1/2,1  )  (1/3,-1/3,1/3)
  

• Atom positions
  
  The atom positions with respect to the previously described unit cell are calculated for each individual sublattice, properly shifted and tagged with the element index. Then the rows are sorted according to the x, y and z coordinate, respectively, and the slab is truncated if the number of rows exceeds the slab size specified in the input file. The result is written to the named output file.

IV. Angular tables file

IV.1. Structure of the file

• float angtab[maxs][maxe]
  
  scattering angle phi in the centre-of-mass system as a function of the impact parameter b and the energy Ep of the incoming atom
  
  
• float elstop[maxs][maxe]
  
  electronic energy loss QsE as a function of the impact parameter b and the energy Ep of the incoming atom
  
  
• float revangtab[maxs][maxe]
  
  scattering angle phi in the centre-of-mass system as a function of the impact parameter b and the energy Es of the scattered atom
  
  
• float revelstop[maxs][maxe]
  
  electronic energy loss QsE as a function of the impact parameter b and the energy Es of the scattered atom
  
  
• float bth[maxangle*2][maxe]
  
  impact parameter b as a function of the scattering angle theta in the laboratory system and the energy Ep of the incoming atom
  
  
• float sig[maxangle*2][maxe]
  
  scattering cross-section sigma in the laboratory system as a function of the scattering angle theta in the laboratory system and the energy Ep of the incoming atom
  
  
• float bthrec[maxangle][maxe]
  
  impact parameter b as a function of the recoiling angle ksi in the laboratory system and the energy Ep of the incoming atom
  
  
• float sigrec[maxangle][maxe]
  
  recoiling cross-section sigma in the laboratory system as a function of the recoiling angle ksi in the laboratory system and the energy Ep of the incoming atom

IV.2. Program maketable

• Input parameters
  
  The program reads the following data:

      m1        mass of the incoming atom (in amu)
      Z1        atomic number of the incoming atom
      m2        mass of the target atom (in amu)
      Z2        atomic number of the target atom
      fcorr     screening length reduction factor
      outfile   name of the output file
  

  The structure of the input file is checked.
  
  
• Universal potential
  
  The interatomic potential V(r) as a function of the interatomic distance r is:

      V(r) = c[0]*SUM[ c[i]*exp( a[i]*r ) ] / r, i = 1..4
  
      a[0] = [ 0.88534*a0*fcorr ] / [ Z1**0.23 + Z2**0.23 ]
      a[1] = -3.2/a[0]
      a[2] = -0.9423/a[0]
      a[3] = -0.4029/a[0]
      a[4] = -0.2016/a[0]
      c[0] = [ Z1*Z2*Q*1.0e-3 ] / [ 4*pi*eps0*1.0e-10 ]
      c[1] = 0.1818
      c[2] = 0.5099
      c[3] = 0.2802
      c[4] = 0.02817
  

  where a0 = 0.52917706 angstrom is the Bohr radius, Q = 1.6021892e-19 C is the electronic charge, and eps0 = 8.85418782e-12 F/m is the permittivity of the vacuum. The factor 1.0e-3 accounts for the conversion from eV to keV, and the factor 1.0e-10 accounts for the conversion from m to angstrom.

      See: D.J. O'Connor and J.P. Biersack
           Nucl. Instr. and Meth. in Phys. Res. B15 91986) 14-19
  

• Angtab
  
  The scattering angle phi in the centre-of-mass system is found by calculating the scattering integral:

      phi = pi - 2*b*Int[ f(r)dr ]
      1/f(r) = r*r*sqrt( 1 - (b/r)**2 - V(r)/E )
  

  Here b is the impact parameter, V(r) is the interatomic potential, and E is the energy of the incoming atom in the centre-of-mass system. The integral runs from rm till infinity, where rm is the distance of closest approach.

      See: D.O. Boerma
           Lecture Notes on Atomic Collisions (1995), formula (3.12) and
           (3.13)
  

  The value of rm is obtained by using Newton's method to solve the equation

      1 - (b/rm)**2 - V(rm)/E = 0
  

  The variable x = sqrt(1-(rm/r)) is substituted for the integration variable r, and the integral is evaluated using Gaussian quadrature with 10 nodes. The integration interval is split in two sub-intervals x = [0,0.5] and x = [0.5,1].

      See: Numerical Recipes in C
           (Cambridge University Press, 1992), sect. 4.5 
  

• Elstop
  
  The electronic energy loss is calculated using the impact parameter dependency introduced by Oen and Robinson:

      Q(s,E) = C*LSS(E)*exp( -0.3*rm(s,E)/a12 )
  

  where s is the impact parameter, E is the energy of the incoming atom in the laboratory system, LSS(E) is the electronic stopping power according to the LSS-theory, rm is the distance of closest approach, and a12 is the screening length according to the Moliere theory (in this program replaced by the reduced screening length a[0]; see Universal potential).

      See: O.S. Oen and M.T. Robinson
           Nucl. Instr. and Meth 132 (1976) 647-653
  

  The normalization constant C follows from the condition that the electronic stopping power should be equal to the LSS-value, and the assumption that rm = b (only valid if V(r) << E):

      LSS(E) = 2*pi*Int[ Q(s,E)sds ] -> C = 0.045 / [ pi*a12*a12 ]
  

  The LSS-stopping power (in eV*m**2) is:

      LSS(E) = [ 8.8e-25*ksi*Z1*Z2*v0 ] / [ Z1**(2/3) + Z2**(2/3) ]**3/2
  
      ksi = Z**(1/6)
      Z1   atomic number of the incoming atom
      Z2   atomic number of the target atom
      v0   velocity of the incoming atom in the laboratory system (in m/s)
  

  Substitute sqrt(2*E/m1) for v0, multiply by sqrt(1000Q/amu), where Q = 1.6021892e-19 C is the electronic charge and amu = 1.6605655e-27 kg, to account for the conversion from J/kg to keV/amu, multiply by (1.0e+10)**2 to account for the change from m to angstrom, and multiply by 1.0-3 to account for the conversion from eV to keV.

      See: J. Lindhard, M. Scharff and H.E. Schiott
           Mat. Fys. Medd. Dan. Vid. Selsk. 33, no.14 (1963) 1-42,
           formula (2.10)
           D.O. Boerma
           Lecture Notes on Atomic Collisions (1995), formula (6.34)   
  

• RevAngtab
• RevElstop
  
  The scattering angle phi and the recoiling angle ksi in the centre-of-mass system are related by the expression:

      phi = pi - 2*ksi
  

  The energy Er of the recoiling atom is:

      Er = Ep*[ 4*mu/(1+mu)**2 ]*cos(ksi)**2
  
      mu = m2/m1
      m1   mass of the incoming atom (in amu)
      m2   mass of the target atom (in amu)
      Ep   energy of the incoming atom in the centre-of-mass system
      Er   energy of the recoiling atom in the centre-of-mass system
  

  For each value of b and Es, the scattering angle phi is obtained in the following iterative way. First the energy of the incoming particle is assumed to be x, and phi and QsE are obtained by linear interpolation from Angtab and Elstop, respectively. Once phi is known, ksi and Er can be calculated. If Ep-Er differs from Es, a better approximation of x is chosen and a new iteration is started.

      See: D.O. Boerma
           Lecture Notes on Atomic Collisions (1995), formula (2.14) and
           (2.22)
  

• Bth
• Sig
  
  The following formula are used:

      mu*cos(phi) = -sin(theta)**2
                    +- cos(theta)*sqrt( mu**2 - sin(theta)**2 )
      dphi/dtheta = [ mu + 2*cos(phi) + (1/mu) ] / [ mu + cos(phi) ]
      sigma(phi) = [ b/sin(phi) ]*| db/dphi |
      sigma(theta) = sigma(phi)*[ sin(phi)/sin(theta) ]*| dphi/dtheta |
                   = [ b/sin(theta) ]*| db/dphi |*| dphi/dtheta |
  
      mu = m2/m1
      m1      mass of the incoming atom (in amu)
      m2      mass of the target atom (in amu)
      phi     scattering angle in the centre-of-mass system
      theta   scattering angle in the laboratory system
  
      See: D.O. Boerma
           Lecture Notes on Atomic Collisions (1995), formula (2.15),
           (2.17), (3.23) and (3.33)
  

  The scattering angle phi is obtained by interpolation from Angtab. For this purpose, a polynomial poly(b) is least-squares fitted to log(phi) as a function of the impact parameter b, for a fixed value of Ep:

      log(phi) = SUM[ coeff[i]*b**i ], i = 0..7
  

  Thus:

      db/dphi = [ db/dlog(phi) ]*[ dlog(phi)/dphi ]
              = [ db/dlog(phi) ] / phi
              = 1 / [ phi*dlog(phi)/db ]
              = 1 / [ phi*dpoly()/db ]
  

  For each value of theta and Ep, the scatering angle phi, the derivatives dphi/dtheta and db/dphi, and the cross-section sigma(theta) are calculated. Finally, the impact parameter b is obtained by using Newton's method to solve the equation

      poly(b) = log(phi)
  

• BthRec
• SigRec
  
  The following formula are used:

      phi = pi - 2*ksi
      sigma(ksi) = 4*sin(phi/2)*sigma(phi)
      sigma(phi) = [ b/sin(phi) ]*| db/dphi |
      db/dphi = 1 / [ phi*dpoly()/db ]
  
      phi     scattering angle in the centre-of-mass system
      ksi     recoiling angle in the laboratory system
  

  For each value of ksi and Ep, the scattering angle phi is calculated. Then the impact parameter b is obtained by using Newton's method to solve the equation

      poly(b) = log(phi)
  

  Finally, the derivative db/dphi, and the cross-sections sigma(phi) and sigma(ksi) are calculated.

      See: D.O. Boerma
           Lecture Notes on Atomic Collisions (1995), formula (2.14), (3.23)
           and (3.25)
  

IV.3. Program showtable

IV.4. Program table2asc

IV.5. Program table2bin

V. Output file

V.1. Structure of the file

    <ncol>
    <ang[1]> <yield[1]> <yield[2]> ............. <yield[ncol-1]> <totyield>
    <ang[1]> <yield[1]> <yield[2]> ............. <yield[ncol-1]> <totyield>
    .......................................................................
    .......................................................................
             <---------------><--------------->
                 element 1        element 2    ................

V.2. Program showscan

V.3. Program sumscan

VI. Specular reflection

VI.1. Program spec

1. If a lattice site in the surface layer is empty, neighbours of this site are added to the list of surface atoms, provided the neighbour is not already a member of this list and the neighbour is located in the second layer of the slab.
   
   
2. Recoil atoms are not followed.
   
   
3. At each scattering event the total number of projectiles, the number of scattering events in the current scattering sequence and the current trajectory structure are written to file <outfile>.<sym>, where <outfile> is the name specified on the command line and <sym> is the element symbol of the projectile.

1. The energy has dropped below the minimum value specified in the input file.
   
   
2. The primary atom has moved to a depth below the bottom of the slab.
   
   
3. The primary atom is scattered away from the slab, to a height of 2.0 angstrom above the surface.

VI.2. Program showspec

VI.3. Program makesurface

    PERCENTAGE   <number>              # percentage vacancies to be added
    CONFIG_FILE  FeFCC.cfg             # original configuration file
    OUTPUT_FILE  dum/FeFCC_<N>%vac.cfg # modified configuration file