17 Kasım 2009 Salı

LAMMPS Howto's : How to simulate non-orthogonal systems in LAMMPS?

Primitive vectors of diamond crystal is:
A1=0, 1/2 , 1/2
A2=1/2 , 0 , 1/2
A3=1/2 , 1/2 , 0
To simulate triclinic systems in lammps lets check LAMMPS manual first:
If the xy xz yz line does not appear, LAMMPS will set up an axis-aligned (orthogonal) simulation box. If the line does appear, LAMMPS creates a non-orthogonal simulation domain shaped as a parallelepiped with triclinic symmetry. See the region prism command for a description of how the extent of the parallelepiped is defined. The parallelepiped has its "origin" at (xlo,ylo,zlo) and 3 edge vectors starting from the origin given by a = (xhi-xlo,0,0); b = (xy,yhi-ylo,0); c = (xz,yz,zhi-zlo). Note that if your simulation will tilt the box, e.g. via the fix deform command, the simulation box must be triclinic, even if the tilt factors are initially 0.0.

So what we need is a set of primitive vectors with first one is aligned with -x axis , second one is at xy plane.
If we apply a rotation for the primitive vectors of diamond to align A1 to X axis and A2 to xy plane our rotated primitive vectors will be:
A1=1/2, 0 , 0
A2=1/2 , sqrt(3)/2 , 0
A3=1/2 , sqrt(3)/6 , sqrt(2/3)
and basis points are:
B1=0 , 0 , 0
B2= 1/4 A1 + 1/4 A2 + 1/4 A3
If we set lattice constant to 5.43 A for Silicon
we will have the following data file for lammps:
Note that
For Lsi stands for lattice constant of Silicon
xhi-xlo=Lsi/sqrt(2)
yhi-ylo=Lsi*sqrt(3/2)*1/2
zhi-zlo=Lsi/sqrt(3)
and tilting factors are:
xy=xz=Lsi/(2*sqrt(2))
yz=Lsi*1/6*sqrt(3/2)


2 atoms
6 atom types

0 3.83959 xlo xhi
0 3.32518 ylo yhi
0 3.1350 zlo zhi
1.91979 1.91979 1.1084 xy xz yz
Masses

1 12.0000
2 1.0080
3 15.9990
4 14.0000
5 14.0000
6 28.0000

Atoms

1 6 0 0.0 0.0 0.0
2 6 0 1.91979 1.10839 0.78375

This data file is prepared to use ReaxFFsio with LAMMPS. However LAMMPS distribution does not contain Si parameters. Older distribution does. You may use those too.
Dundar Yilmaz.