Lammps -- My Documentation

This is my documentation about LAMMPS (Large-scale Atomic/Molecular Massively Parallel Simulator). Here I'm putting some very especific info for a very especific application. The official documentation cab be found in LAMMPS-OFF-DOC.

Input Script

  1. Initialization:
  2. Commands: units, dimension, processors, boundary, atom_style, atom_modify

  3. Atom definition:
  4. Create atoms on a lattice using lattice, region, create_box, create_atoms. Atoms cab be defined also froam a data or restart file.

  5. Settings:
  6. Once atoms are defined , we need to set:

  7. Run simulation:
  8. Run simulation using run. Energy minimization with minimize command.

Simple Script

# 2d Born potential gas
# U = A exp(\frac{\sigma-r}{\rho} - \frac{C}{r^6} + \frac{D}{r^8}

units		lj
dimension	2
boundary	p p p
atom_style	atomic

# create geometry
lattice		sq 0.8
region		box block 0 35 0 35 -0.25 0.25
create_box	1 box
create_atoms	1 box

mass		1 1.0
# velocity	all create temperature random_seed
velocity	all create 1.0 200
velocity	all zero linear

# Potential
pair_style	born 15.0
pair_coeff	* * 1.0 1.0 0.0 0.0 0.0

# fix
fix		1 all nve

# Run
timestep	0.001
thermo		1000
thermo_style	custom step etotal ke pe temp press
dump		1 all custom 1000 dump.vel id vx vy vz

compute         myRDF all rdf 50
# fix   ID group-ID ave/time Nevery Nrepeat Nfreq 
fix             2 all ave/time 1 100 1000 c_myRDF[1] c_myRDF[2] file tmp.rdf mode vector

run		10000
    
  1. Initialization:
  2. Units are Lennard-Jones (defalult value): all quantities are unitless. Mass, sigma, epsilon, and Boltzmann contant are set to 1. Other options are real, metal, si, cgs, electron, micro, and nano.

    We are defined a two dimensional simulation with the command: dimension 2.

    For our simulation we are using periodic boundary condition on all coordinates (defalult values). Other options are f: non-periodic and fixed, s: non-periodic and shrink-wrapped, and m: non-periodic and shrink-wrapped witn minimum value.

    Atom style is atomic (defalult value), which is used coarse-grain liquids, solids and metals.

  3. Atom definition:
  4. lattice sq 0.8: This is a 2D lattice with reduced density equals to 0.8 (arbitarily choosen), where $$ \rho^{*}=\frac{N}{V^{*}}=\frac{N}{{r^{*}}^d} $$ with \(N\) the number of basis molucules, \(V^{*}\) is the reduce volume, and \(d\) is the dimension.

    region box block 0 35 0 35 -0.25 0.25: this simulation will create 1225 atoms, separated by 1 lj-unit, so in one of the directions we will have \(\sqrt{1225}=35\) atoms separated one cell unit. Then this command will create a box with coordinates in x from 0 to 35, in y from 0 to 35 and from -0.25 to 0.25 in z. The input box is just an ID asssigned by the user, and block is a style of region (others options are cone, cylinder, etc).

    create_box 1 box creates simulation box for N types of atoms (1 in this case) for the region with ID called "box".

    create_atoms 1 box creates N types of atoms (1 in this case) in the region with ID called "box". The simulation box must already exists.

    mass 1 1.0 Set the mass for the atom type N (in this case mass=1.0 for atom type 1).

    velocity all create 1.0 200 creates ensamble of velocities for 'all' atoms using random number generator (with seed=200) at a specific temperature (T*=1.0).

    velocity all zero linear set the total linear momentum of the atomic system to zero.

  5. Settings:
  6. pair_style born 15.0

    pair_coeff * * 1.0 1.0 0.0 0.0 0.0

    Here we will use the Born-Mayer-Huggins potential, with cutoff radius of \(r_c=15.0\) lj-units. The potential is given by $$ E = A \exp\left( \frac{\sigma - r}{\rho} \right) - \frac{C}{r^6} + \frac{D}{r^8}, ~~r < r_c $$ For now, we just need the exponential part of the potential so we set its parameters with pair_coeff: first two entries specified the atoms with this interaction (in this case the interaction is for all atoms, so we used * *), then we specified the numerical values of the constants A=1.0, \(\rho=1.0\), \( \sigma=0.0 \), C=0.0, and D=0.0

    fix 1 all nve The fix 'nve' is applied to a group of atoms ('all' in this case). nve performs a constant NVE integration, this creates a system trajectory consistent with the microcanonical ensamble.

  7. Run:
  8. timestep 0.001 Set timestep

    thermo 1000

    thermo_style custom step etotal ke pe temp press Output themodynamics info every N steps (1000 in this case-- If 0 just gives info in the initial and final timestep), In thermo_style we specified the info we want to output (custom). In this case we are getting the number of timestep, total energy, kinetic and potential enegies, temperature and pressure.

    dump 1 all custom 1000 dump.vel id vx vy vz Dump a snapshot of atoms N (1 in this example) quantities (velocities components in this case). The snapshots are stored in the file 'dump.vel'

    compute myRDF all rdf 50

    Compute the Radial Distribution Function for all particles with N=50 bins

    fix 2 all ave/time 10 100 1000 c_myRDF[1] c_myRDF[2] file tmp.rdf mode vector

    The averaged RDF is generated on timesteps multiple of Nfreq. The average is calculated over Nrepeat quantities, computed using samples taken every Nevery.

    Nfreq must be a multiple of Nevery. Nrepeat*Nevery can not exceed Nfreq. If Nrepeat=1 and Nfreq=100, there are not time averaging, values are generated on timesteps 100, 200, etc

    run 10000 Run the simulation for 10,000 timesteps.