A number of utility functions are provided to set up the most common surfaces, to add vacuum layers, and to add adsorbates to a surface. In general, all surfaces can be set up with the modules described in the section General crystal structures and surfaces, but these utility functions make common tasks easier.
To setup an Al(111) surface with a hydrogen atom adsorbed in an on-top position:
from ase.lattice.surface import fcc111
slab = fcc111('Al', size=(2,2,3), vacuum=10.0)
This will produce a slab 2x2x3 times the minimal possible size, with a (111) surface in the z direction. A 10 Å vacuum layer is added on each side.
To set up the same surface with with a hydrogen atom adsorbed in an on-top position 1.5 Å above the top layer:
from ase.lattice.surface import fcc111, add_adsorbate
slab = fcc111('Al', size=(2,2,3))
add_adsorbate(slab, 'H', 1.5, 'ontop')
slab.center(vacuum=10.0, axis=2)
Note that in this case is is probably not meaningful to use the vacuum keyword to fcc111, as we want to leave 10 Å of vacuum after the adsorbate has been added. Instead, the center() method of the Atoms is used to add the vacuum and center the system.
The atoms in the slab will have tags set to the layer number: First layer atoms will have tag=1, second layer atoms will have tag=2, and so on. Adsorbates get tag=0:
>>> print atoms.get_tags()
[3 3 3 3 2 2 2 2 1 1 1 1 0]
This can be useful for setting up ase.constraints (see Diffusion of gold atom on Al(100) surface (NEB)).
All the functions setting up surfaces take the same arguments.
Each function defines a number of standard adsorption sites that can later be used when adding an adsorbate with ase.lattice.surface.add_adsorbate().
These allways give orthorhombic cells:
fcc100 | ![]() |
fcc110 | ![]() |
bcc100 | ![]() |
hcp10m10 | ![]() |
diamond100 | ![]() |
These can give both non-orthorhombic and orthorhombic cells:
fcc111 | ![]() |
![]() |
fcc111_root | varies on root | not_implemented |
fcc211 | not implemented | ![]() |
bcc110 | ![]() |
![]() |
bcc111 | ![]() |
![]() |
hcp0001 | ![]() |
![]() |
diamond111 | ![]() |
not implemented |
The adsorption sites are marked with:
ontop | hollow | fcc | hcp | bridge | shortbridge | longbridge |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
This can be used for MX2 2D structures such as MoS2:
Create three-layer 2D materials with hexagonal structure.
For metal dichalcogenites, ect.
The kind argument accepts ‘2H’, which gives a mirror plane symmetry and ‘1T’, which gives an inversion symmetry.
If you need other surfaces than the ones above, the easiest is to look in the source file surface.py, and adapt one of the existing functions. Please contribute any such function that you make either by checking it into SVN or by mailing it to the developers.
After a slab has been created, a vacuum layer can be added. It is also possible to add one or more adsorbates.
Add an adsorbate to a surface.
This function adds an adsorbate to a slab. If the slab is produced by one of the utility functions in ase.lattice.surface, it is possible to specify the position of the adsorbate by a keyword (the supported keywords depend on which function was used to create the slab).
If the adsorbate is a molecule, the atom indexed by the mol_index optional argument is positioned on top of the adsorption position on the surface, and it is the responsibility of the user to orient the adsorbate in a sensible way.
This function can be called multiple times to add more than one adsorbate.
Parameters:
slab: The surface onto which the adsorbate should be added.
height: Height above the surface.
Note position is given in absolute xy coordinates (or as a keyword), whereas offset is specified in unit cells. This can be used to give the positions in units of the unit cell by using offset instead.
New in version 3.5.2.
In addition to the most normal surfaces, a function has been constructed to create more uncommon surfaces that one could be interested in. It is constructed upon the Miller Indices defining the surface and can be used for both fcc, bcc and hcp structures. The theory behind the implementation can be found here: general_surface.pdf.
Create surface from a given lattice and Miller indices.
To setup a Au(211) surface with 9 layers and 10 Å of vacuum:
from ase.lattice.surface import surface
s1 = surface('Au', (2, 1, 1), 9)
s1.center(vacuum=10, axis=2)
This is the easy way, where you use the experimental lattice constant for gold bulk structure. You can write:
from ase.visualize import view
view(s1)
or simply s1.edit() if you want to see and rotate the structure.
Next example is a molybdenum bcc(321) surface where we decide what lattice constant to use:
from ase.lattice import bulk
Mobulk = bulk('Mo', 'bcc', a=3.16, cubic=True)
s2 = surface(Mobulk, (3, 2, 1), 9)
s2.center(vacuum=10, axis=2)
As the last example, creation of alloy surfaces is also very easily carried out with this module. In this example, two Pt3Rh fcc(211) surfaces will be created:
a = 4.0
from ase import Atoms
Pt3Rh = Atoms('Pt3Rh',
scaled_positions=[(0, 0, 0),
(0.5, 0.5, 0),
(0.5, 0, 0.5),
(0, 0.5, 0.5)],
cell=[a, a, a],
pbc=True)
s3 = surface(Pt3Rh, (2, 1, 1), 9)
s3.center(vacuum=10, axis=2)
Pt3Rh.set_chemical_symbols('PtRhPt2')
s4 = surface(Pt3Rh , (2, 1, 1), 9)
s4.center(vacuum=10, axis=2)