git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@13826 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp
2015-08-04 18:54:11 +00:00
parent 57ed470003
commit b0dec93ffb
618 changed files with 136681 additions and 0 deletions

View File

@ -0,0 +1,23 @@
Description:
This is a simulation of pyramid-shaped objects resting on an immobile surface
(resembling graphene). Each pyramid is built from spherical particles stacked
like cannon-balls, or stacked fruit. Ordinarily, the stack does not move
because the particles at the ground layer are immobilized. However,
given an initial (small) perturbation the pyramids collapse in an avalanche.
(In this example, the perturbation is due to shock because we (intentionally)
did not minimize the system before starting the simulation. This shock
causes an avalanche to occur approximately 5000 timesteps later.)
The particles roll down the pyramid and bounce off the "ground". The bouncing
is due to a repulsive external force which is added artificially.
(See the "run.in" file.) The simulation looks weird without something
to bounce off of. So I added a graphene surface at the bottom as scenery.
(It does not exert any force on the atoms.)
(Random comment: This could be a fun example to illustrate the Boltzmann
distribution. Because there is no damping, in a small region, I'm guessing
the particle heights should eventually approach the Boltzmann distribution
for some temperature consistent with the initial potential energy of the
system.)

View File

@ -0,0 +1,28 @@
------- A note on building the graphene sheet in VMD: ------
Probably you can ignore these instructions.
These instructions are not necessary for this example to run.
This example contains several pyramid shaped objects resting on a surface
made of graphene. The instructions in this file explain how to build the
graphene (representing the "ground") using VMD instead of with moltemplate.
Why do this?
VMD can create graphene sheets with bonds connecting neighboring carbon atoms,
(which looks more pretty). However, as of 2013-4-29, moltemplate currently
can not generate these bonds. It does not matter physically in this case,
because the graphene sheet used here does not move. It is only used as
scenery, to graphically represent the ground surface.
Select "Extensions"->"Modeling"->"Carbon Nanotube Builder"
Build a graphene sheet of size 39.8 x 39.8 (units: nm)
400.3358398 399.876008
(try to use a size compatible with the periodic boundaries)
Select "Extensions"->"Tk Console", and type
display backgroundgradient on
Note: If you want to do this, before you run moltemplate, you may want to delete
the sections of the "system.lt" file (located in "moltemplate_files")
which define the graphene wall. Instead create the graphene data file
in VMD. You will have to manually merge the data file for graphene
with the data file for the pyramids created by moltemplate,
(taking care to avoid overlapping atom-id numbers).

View File

@ -0,0 +1,23 @@
# Use these commands to generate the LAMMPS input script and data file
# (and other auxilliary files):
# Create LAMMPS input files this way:
cd moltemplate_files
# run moltemplate
moltemplate.sh system.lt
# This will generate various files with names ending in *.in* and *.data.
# These files are the input files directly read by LAMMPS. Move them to
# the parent directory (or wherever you plan to run the simulation).
mv -f system.in* system.data ../
# Optional:
# The "./output_ttree/" directory is full of temporary files generated by
# moltemplate. They can be useful for debugging, but are usually thrown away.
rm -rf output_ttree/
cd ../

View File

@ -0,0 +1,76 @@
------- To view a lammps trajectory in VMD --------
1) Build a PSF file for use in viewing with VMD.
This step works with VMD 1.9 and topotools 1.2.
(Older versions, like VMD 1.8.6, don't support this.)
a) Start VMD
b) Menu Extensions->Tk Console
c) Enter:
(I assume that the the DATA file is called "system.data")
topo readlammpsdata system.data full
animate write psf system.psf
2)
Later, to Load a trajectory in VMD:
Start VMD
Select menu: File->New Molecule
-Browse to select the PSF file you created above, and load it.
(Don't close the window yet.)
-Browse to select the trajectory file.
If necessary, for "file type" select: "LAMMPS Trajectory"
Load it.
---- A note on trajectory format: -----
If the trajectory is a DUMP file, then make sure the it contains the
information you need for pbctools (see below. I've been using this
command in my LAMMPS scripts to create the trajectories:
dump 1 all custom 5000 DUMP_FILE.lammpstrj id mol type x y z ix iy iz
It's a good idea to use an atom_style which supports molecule-ID numbers
so that you can assign a molecule-ID number to each atom. (I think this
is needed to wrap atom coordinates without breaking molecules in half.)
Of course, you don't have to save your trajectories in DUMP format,
(other formats like DCD work fine) I just mention dump files
because these are the files I'm familiar with.
3) ----- Wrap the coordinates to the unit cell
(without cutting the molecules in half)
a) Start VMD
b) Load the trajectory in VMD (see above)
c) Menu Extensions->Tk Console
d) Try entering these commands:
pbc wrap -compound res -all
pbc box
----- Optional ----
To shift the box by a fraction in the x direction (for example)
do this:
pbc wrap -compound res -all -shiftcenterrel {-0.50 -0.52 0.0 }
pbc box -shiftcenterrel {-0.50 -0.52 0.0 }
# Alternately if you have a solute whose atoms are all of type 1,
# then you can also try this to center the box around it:
pbc wrap -sel type=1 -all -centersel type=2 -center com
4) Optional: If you like, change the atom types in the PSF file so
that VMD recognizes the atom types, use something like:
sed -e 's/ 1 1 / C C /g' < system.psf > temp1.psf
sed -e 's/ 2 2 / H H /g' < temp1.psf > temp2.psf
sed -e 's/ 3 3 / P P /g' < temp2.psf > system.psf
(If you do this, it might effect step 2 above.)

Binary file not shown.

After

Width:  |  Height:  |  Size: 86 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 89 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 96 KiB

View File

@ -0,0 +1,15 @@
# This directory contains moltemplate files for the "Pyramids of Giza" example.
# (Note: the ground lattice work that appears in some images was not generated
# by moltemplate. Moltemplate can not currently create bonded periodic
# structures as of 2013-4-04. Those were generated by topotools.)
#
# To run moltemplate, use:
moltemplate.sh system.lt
# This will generate:system.data, system.in, system.in.init, system.in.settings
#
# The output_ttree/ directory will contain files like "Data Atoms", "Data Bonds"
# which contain the corresponding structures of the system.data file.
# (This might make it slightly easier to combine them with atom data and
# bond data generated by other programs, such as topotools, for example.)

View File

@ -0,0 +1,61 @@
# This file contains a unit cell for building graphene and nanotubes
#
#
# The 4AtomRectCellXY "molecule" defined below is a reactangular unit cell
# for hexagonal tesselations in 2-dimensions. (See "graphene_unit_cell.jpg")
# Surfaces constructed with this unit cell can be flat or curved into tubes.
# The distance between nearest-neighbor carbon atoms (ie the length of a
# carbon-carbon bond) is equal to "d" which I set to 1.42 Angstroms.
#
# d = length of each hexagon's side = 1.42 Angstroms
# L = length of each hexagon = 2*d = 2.84 Angstroms
# W = width of each hexagon = 2*d*sqrt(3)/2 = 2.4595121467478056 Angstroms
# 2w = width of hexagon rows = 3*l = 4.26 Angstroms
#
# Consequently, the Lattice-cell vectors for singe-layer graphene are:
# (2.4595121467478, 0, 0) (aligned with X axis)
# (0, 4.26, 0) (aligned with Y axis)
# So, to build a sheet of graphite, you could use:
# sheet = new Graphene/4AtomRectCellXY [10].move(2.4595121467478, 0, 0)
# [10].move(0, 4.26, 0)
Graphene {
4AtomRectCellXY
{
# atomID molID atomType charge x y z
write("Data Atoms") {
$atom:C11 $mol:... @atom:../C 0.0 0.61487803668695 0.71 0.0
$atom:C21 $mol:... @atom:../C 0.0 1.84463411006085 1.42 0.0
$atom:C12 $mol:... @atom:../C 0.0 0.61487803668695 3.55 0.0
$atom:C22 $mol:... @atom:../C 0.0 1.84463411006085 2.84 0.0
}
}
# Now define properties of the Carbon graphene atom
write_once("In Init") {
pair_style hybrid lj/cut 9.0
}
write_once("Data Masses") {
@atom:C 12.0
}
write_once("In Settings") {
# i j epsilon sigma
pair_coeff @atom:C @atom:C lj/cut 0.068443 3.407
# These Lennard-Jones parameters come from
# R. Saito, R. Matsuo, T. Kimura, G. Dresselhaus, M.S. Dresselhaus,
# Chem Phys Lett, 348:187 (2001)
# Define a group consisting of only carbon atoms in graphene molecules
group gGraphene type @atom:C
}
} # Graphene

View File

@ -0,0 +1,21 @@
import "graphene.lt"
# -------------- graphene sheet -----------------
# Notes:
# Hexagonal lattice with:
# l = length of each hexagonal side = 1.42 Angstroms
# L = length of each hexagon = 2*l = 2.84 Angstroms
# W = width of each hexagon = 2*l*sqrt(3)/2 ~= 2.4595121467478 Angstroms
# 2w = width of hexagon rows = 3.0*l = 4.26 Angstroms
GrapheneWall {
unitcells = new Graphene/4AtomRectCellXY [163].move(2.456, 0, 0)
[94].move(0, 4.254, 0)
# (Note: I fudged the spacing slightly to make it line up better with the
# lattice spacing for graphene generated by VMD's graphene builder.)
}

View File

@ -0,0 +1,283 @@
# Brick is a very simple molecule containing one "atom".
# "ImmobileBrick" and "GoldBrick" are identical to "Brick" but are
# given different atom types. (This makes it easier to put them in
# different groups and apply different LAMMPS "fixes" to them.)
Brick {
# atomID molID atomType charge x y z
write("Data Atoms") {
$atom $mol @atom 0.0 0.0 0.0 0.0
}
write_once("Data Masses") {
@atom 1.0
}
write_once("In Settings") {
# U(r) = 4*epsilon_ij*((sigma_ij/r)^12 - (sigma_ij/r)^6)
#
# i j eps sig
pair_coeff @atom @atom lj/cut 1.0 0.8908987181403393
}
write_once("In Settings") {
group gMobile type @atom
# (Atoms of this type belong to the "gMobile" group)
}
write_once("In Init") {
atom_style full
units lj
pair_style lj/cut 1.0
}
}
#We want to stack "Brick"s the same way a green-grocer sometimes stack apples:
#Place the apples at the base an square lattice of apples at the base.
#The apples in the next layer up are placed in between the 4 apples beneath them.
#Each new layer is smaller and placed above the previous layer at a height
#of sigma / sqrt(2), where "sigma" is the diameter of each spherical "Brick".
#We will artificially hold the apples at the base in place
#(to keep the entire stack from collapsing).
#
# The lines below were generated from the following python loop:
#
#from math import *
#N=50
#s=1.0
#for i in range(0,N):
# print(" layer"+str(i)+" = new Brick.move("+str(-(N-(i+1))*s*0.5)+","+
# str(-(N-(i+1))*s*0.5)+","+str(i*s/sqrt(2))+") ["+str(N-i)+"].move("+
# str(s)+",0,0) ["+str(N-i)+"].move(0,"+str(s)+",0)")
PyramidKhufu {
layer0 = new ImmobileBrick.move(-24.5,-24.5,0.0) [50].move(1.0,0,0) [50].move(0,1.0,0)
layer1 = new Brick.move(-24.0,-24.0,0.707106781187) [49].move(1.0,0,0) [49].move(0,1.0,0)
layer2 = new Brick.move(-23.5,-23.5,1.41421356237) [48].move(1.0,0,0) [48].move(0,1.0,0)
layer3 = new Brick.move(-23.0,-23.0,2.12132034356) [47].move(1.0,0,0) [47].move(0,1.0,0)
layer4 = new Brick.move(-22.5,-22.5,2.82842712475) [46].move(1.0,0,0) [46].move(0,1.0,0)
layer5 = new Brick.move(-22.0,-22.0,3.53553390593) [45].move(1.0,0,0) [45].move(0,1.0,0)
layer6 = new Brick.move(-21.5,-21.5,4.24264068712) [44].move(1.0,0,0) [44].move(0,1.0,0)
layer7 = new Brick.move(-21.0,-21.0,4.94974746831) [43].move(1.0,0,0) [43].move(0,1.0,0)
layer8 = new Brick.move(-20.5,-20.5,5.65685424949) [42].move(1.0,0,0) [42].move(0,1.0,0)
layer9 = new Brick.move(-20.0,-20.0,6.36396103068) [41].move(1.0,0,0) [41].move(0,1.0,0)
layer10 = new Brick.move(-19.5,-19.5,7.07106781187) [40].move(1.0,0,0) [40].move(0,1.0,0)
layer11 = new Brick.move(-19.0,-19.0,7.77817459305) [39].move(1.0,0,0) [39].move(0,1.0,0)
layer12 = new Brick.move(-18.5,-18.5,8.48528137424) [38].move(1.0,0,0) [38].move(0,1.0,0)
layer13 = new Brick.move(-18.0,-18.0,9.19238815543) [37].move(1.0,0,0) [37].move(0,1.0,0)
layer14 = new Brick.move(-17.5,-17.5,9.89949493661) [36].move(1.0,0,0) [36].move(0,1.0,0)
layer15 = new Brick.move(-17.0,-17.0,10.6066017178) [35].move(1.0,0,0) [35].move(0,1.0,0)
layer16 = new Brick.move(-16.5,-16.5,11.313708499) [34].move(1.0,0,0) [34].move(0,1.0,0)
layer17 = new Brick.move(-16.0,-16.0,12.0208152802) [33].move(1.0,0,0) [33].move(0,1.0,0)
layer18 = new Brick.move(-15.5,-15.5,12.7279220614) [32].move(1.0,0,0) [32].move(0,1.0,0)
layer19 = new Brick.move(-15.0,-15.0,13.4350288425) [31].move(1.0,0,0) [31].move(0,1.0,0)
layer20 = new Brick.move(-14.5,-14.5,14.1421356237) [30].move(1.0,0,0) [30].move(0,1.0,0)
layer21 = new Brick.move(-14.0,-14.0,14.8492424049) [29].move(1.0,0,0) [29].move(0,1.0,0)
layer22 = new Brick.move(-13.5,-13.5,15.5563491861) [28].move(1.0,0,0) [28].move(0,1.0,0)
layer23 = new Brick.move(-13.0,-13.0,16.2634559673) [27].move(1.0,0,0) [27].move(0,1.0,0)
layer24 = new Brick.move(-12.5,-12.5,16.9705627485) [26].move(1.0,0,0) [26].move(0,1.0,0)
layer25 = new Brick.move(-12.0,-12.0,17.6776695297) [25].move(1.0,0,0) [25].move(0,1.0,0)
layer26 = new Brick.move(-11.5,-11.5,18.3847763109) [24].move(1.0,0,0) [24].move(0,1.0,0)
layer27 = new Brick.move(-11.0,-11.0,19.091883092) [23].move(1.0,0,0) [23].move(0,1.0,0)
layer28 = new Brick.move(-10.5,-10.5,19.7989898732) [22].move(1.0,0,0) [22].move(0,1.0,0)
layer29 = new Brick.move(-10.0,-10.0,20.5060966544) [21].move(1.0,0,0) [21].move(0,1.0,0)
layer30 = new Brick.move(-9.5,-9.5,21.2132034356) [20].move(1.0,0,0) [20].move(0,1.0,0)
layer31 = new Brick.move(-9.0,-9.0,21.9203102168) [19].move(1.0,0,0) [19].move(0,1.0,0)
layer32 = new Brick.move(-8.5,-8.5,22.627416998) [18].move(1.0,0,0) [18].move(0,1.0,0)
layer33 = new Brick.move(-8.0,-8.0,23.3345237792) [17].move(1.0,0,0) [17].move(0,1.0,0)
layer34 = new Brick.move(-7.5,-7.5,24.0416305603) [16].move(1.0,0,0) [16].move(0,1.0,0)
layer35 = new Brick.move(-7.0,-7.0,24.7487373415) [15].move(1.0,0,0) [15].move(0,1.0,0)
layer36 = new Brick.move(-6.5,-6.5,25.4558441227) [14].move(1.0,0,0) [14].move(0,1.0,0)
layer37 = new Brick.move(-6.0,-6.0,26.1629509039) [13].move(1.0,0,0) [13].move(0,1.0,0)
layer38 = new Brick.move(-5.5,-5.5,26.8700576851) [12].move(1.0,0,0) [12].move(0,1.0,0)
layer39 = new Brick.move(-5.0,-5.0,27.5771644663) [11].move(1.0,0,0) [11].move(0,1.0,0)
layer40 = new GoldBrick.move(-4.5,-4.5,28.2842712475) [10].move(1.0,0,0) [10].move(0,1.0,0)
layer41 = new GoldBrick.move(-4.0,-4.0,28.9913780286) [9].move(1.0,0,0) [9].move(0,1.0,0)
layer42 = new GoldBrick.move(-3.5,-3.5,29.6984848098) [8].move(1.0,0,0) [8].move(0,1.0,0)
layer43 = new GoldBrick.move(-3.0,-3.0,30.405591591) [7].move(1.0,0,0) [7].move(0,1.0,0)
layer44 = new GoldBrick.move(-2.5,-2.5,31.1126983722) [6].move(1.0,0,0) [6].move(0,1.0,0)
layer45 = new GoldBrick.move(-2.0,-2.0,31.8198051534) [5].move(1.0,0,0) [5].move(0,1.0,0)
layer46 = new GoldBrick.move(-1.5,-1.5,32.5269119346) [4].move(1.0,0,0) [4].move(0,1.0,0)
layer47 = new GoldBrick.move(-1.0,-1.0,33.2340187158) [3].move(1.0,0,0) [3].move(0,1.0,0)
layer48 = new GoldBrick.move(-0.5,-0.5,33.941125497) [2].move(1.0,0,0) [2].move(0,1.0,0)
layer49 = new GoldBrick.move(0.0,0.0,34.6482322781) [1].move(1.0,0,0) [1].move(0,1.0,0)
}
PyramidKhafre {
layer0 = new ImmobileBrick.move(-23.5,-23.5,0.0) [48].move(1.0,0,0) [48].move(0,1.0,0)
layer1 = new Brick.move(-23.0,-23.0,0.707106781187) [47].move(1.0,0,0) [47].move(0,1.0,0)
layer2 = new Brick.move(-22.5,-22.5,1.41421356237) [46].move(1.0,0,0) [46].move(0,1.0,0)
layer3 = new Brick.move(-22.0,-22.0,2.12132034356) [45].move(1.0,0,0) [45].move(0,1.0,0)
layer4 = new Brick.move(-21.5,-21.5,2.82842712475) [44].move(1.0,0,0) [44].move(0,1.0,0)
layer5 = new Brick.move(-21.0,-21.0,3.53553390593) [43].move(1.0,0,0) [43].move(0,1.0,0)
layer6 = new Brick.move(-20.5,-20.5,4.24264068712) [42].move(1.0,0,0) [42].move(0,1.0,0)
layer7 = new Brick.move(-20.0,-20.0,4.94974746831) [41].move(1.0,0,0) [41].move(0,1.0,0)
layer8 = new Brick.move(-19.5,-19.5,5.65685424949) [40].move(1.0,0,0) [40].move(0,1.0,0)
layer9 = new Brick.move(-19.0,-19.0,6.36396103068) [39].move(1.0,0,0) [39].move(0,1.0,0)
layer10 = new Brick.move(-18.5,-18.5,7.07106781187) [38].move(1.0,0,0) [38].move(0,1.0,0)
layer11 = new Brick.move(-18.0,-18.0,7.77817459305) [37].move(1.0,0,0) [37].move(0,1.0,0)
layer12 = new Brick.move(-17.5,-17.5,8.48528137424) [36].move(1.0,0,0) [36].move(0,1.0,0)
layer13 = new Brick.move(-17.0,-17.0,9.19238815543) [35].move(1.0,0,0) [35].move(0,1.0,0)
layer14 = new Brick.move(-16.5,-16.5,9.89949493661) [34].move(1.0,0,0) [34].move(0,1.0,0)
layer15 = new Brick.move(-16.0,-16.0,10.6066017178) [33].move(1.0,0,0) [33].move(0,1.0,0)
layer16 = new Brick.move(-15.5,-15.5,11.313708499) [32].move(1.0,0,0) [32].move(0,1.0,0)
layer17 = new Brick.move(-15.0,-15.0,12.0208152802) [31].move(1.0,0,0) [31].move(0,1.0,0)
layer18 = new Brick.move(-14.5,-14.5,12.7279220614) [30].move(1.0,0,0) [30].move(0,1.0,0)
layer19 = new Brick.move(-14.0,-14.0,13.4350288425) [29].move(1.0,0,0) [29].move(0,1.0,0)
layer20 = new Brick.move(-13.5,-13.5,14.1421356237) [28].move(1.0,0,0) [28].move(0,1.0,0)
layer21 = new Brick.move(-13.0,-13.0,14.8492424049) [27].move(1.0,0,0) [27].move(0,1.0,0)
layer22 = new Brick.move(-12.5,-12.5,15.5563491861) [26].move(1.0,0,0) [26].move(0,1.0,0)
layer23 = new Brick.move(-12.0,-12.0,16.2634559673) [25].move(1.0,0,0) [25].move(0,1.0,0)
layer24 = new Brick.move(-11.5,-11.5,16.9705627485) [24].move(1.0,0,0) [24].move(0,1.0,0)
layer25 = new Brick.move(-11.0,-11.0,17.6776695297) [23].move(1.0,0,0) [23].move(0,1.0,0)
layer26 = new Brick.move(-10.5,-10.5,18.3847763109) [22].move(1.0,0,0) [22].move(0,1.0,0)
layer27 = new Brick.move(-10.0,-10.0,19.091883092) [21].move(1.0,0,0) [21].move(0,1.0,0)
layer28 = new Brick.move(-9.5,-9.5,19.7989898732) [20].move(1.0,0,0) [20].move(0,1.0,0)
layer29 = new Brick.move(-9.0,-9.0,20.5060966544) [19].move(1.0,0,0) [19].move(0,1.0,0)
layer30 = new Brick.move(-8.5,-8.5,21.2132034356) [18].move(1.0,0,0) [18].move(0,1.0,0)
layer31 = new Brick.move(-8.0,-8.0,21.9203102168) [17].move(1.0,0,0) [17].move(0,1.0,0)
layer32 = new Brick.move(-7.5,-7.5,22.627416998) [16].move(1.0,0,0) [16].move(0,1.0,0)
layer33 = new Brick.move(-7.0,-7.0,23.3345237792) [15].move(1.0,0,0) [15].move(0,1.0,0)
layer34 = new Brick.move(-6.5,-6.5,24.0416305603) [14].move(1.0,0,0) [14].move(0,1.0,0)
layer35 = new GoldBrick.move(-6.0,-6.0,24.7487373415) [13].move(1.0,0,0) [13].move(0,1.0,0)
layer36 = new GoldBrick.move(-5.5,-5.5,25.4558441227) [12].move(1.0,0,0) [12].move(0,1.0,0)
layer37 = new GoldBrick.move(-5.0,-5.0,26.1629509039) [11].move(1.0,0,0) [11].move(0,1.0,0)
layer38 = new GoldBrick.move(-4.5,-4.5,26.8700576851) [10].move(1.0,0,0) [10].move(0,1.0,0)
layer39 = new GoldBrick.move(-4.0,-4.0,27.5771644663) [9].move(1.0,0,0) [9].move(0,1.0,0)
layer40 = new GoldBrick.move(-3.5,-3.5,28.2842712475) [8].move(1.0,0,0) [8].move(0,1.0,0)
layer41 = new GoldBrick.move(-3.0,-3.0,28.9913780286) [7].move(1.0,0,0) [7].move(0,1.0,0)
layer42 = new GoldBrick.move(-2.5,-2.5,29.6984848098) [6].move(1.0,0,0) [6].move(0,1.0,0)
layer43 = new GoldBrick.move(-2.0,-2.0,30.405591591) [5].move(1.0,0,0) [5].move(0,1.0,0)
layer44 = new GoldBrick.move(-1.5,-1.5,31.1126983722) [4].move(1.0,0,0) [4].move(0,1.0,0)
layer45 = new GoldBrick.move(-1.0,-1.0,31.8198051534) [3].move(1.0,0,0) [3].move(0,1.0,0)
layer46 = new GoldBrick.move(-0.5,-0.5,32.5269119346) [2].move(1.0,0,0) [2].move(0,1.0,0)
layer47 = new GoldBrick.move(0.0,0.0,33.2340187158) [1].move(1.0,0,0) [1].move(0,1.0,0)
}
PyramidMenkaure {
layer0 = new ImmobileBrick.move(-9.0,-9.0,0.0) [19].move(1.0,0,0) [19].move(0,1.0,0)
layer1 = new Brick.move(-8.5,-8.5,0.707106781187) [18].move(1.0,0,0) [18].move(0,1.0,0)
layer2 = new Brick.move(-8.0,-8.0,1.41421356237) [17].move(1.0,0,0) [17].move(0,1.0,0)
layer3 = new Brick.move(-7.5,-7.5,2.12132034356) [16].move(1.0,0,0) [16].move(0,1.0,0)
layer4 = new Brick.move(-7.0,-7.0,2.82842712475) [15].move(1.0,0,0) [15].move(0,1.0,0)
layer5 = new Brick.move(-6.5,-6.5,3.53553390593) [14].move(1.0,0,0) [14].move(0,1.0,0)
layer6 = new Brick.move(-6.0,-6.0,4.24264068712) [13].move(1.0,0,0) [13].move(0,1.0,0)
layer7 = new Brick.move(-5.5,-5.5,4.94974746831) [12].move(1.0,0,0) [12].move(0,1.0,0)
layer8 = new Brick.move(-5.0,-5.0,5.65685424949) [11].move(1.0,0,0) [11].move(0,1.0,0)
layer9 = new Brick.move(-4.5,-4.5,6.36396103068) [10].move(1.0,0,0) [10].move(0,1.0,0)
layer10 = new Brick.move(-4.0,-4.0,7.07106781187) [9].move(1.0,0,0) [9].move(0,1.0,0)
layer11 = new Brick.move(-3.5,-3.5,7.77817459305) [8].move(1.0,0,0) [8].move(0,1.0,0)
layer12 = new Brick.move(-3.0,-3.0,8.48528137424) [7].move(1.0,0,0) [7].move(0,1.0,0)
layer13 = new Brick.move(-2.5,-2.5,9.19238815543) [6].move(1.0,0,0) [6].move(0,1.0,0)
layer14 = new Brick.move(-2.0,-2.0,9.89949493661) [5].move(1.0,0,0) [5].move(0,1.0,0)
layer15 = new Brick.move(-1.5,-1.5,10.6066017178) [4].move(1.0,0,0) [4].move(0,1.0,0)
layer16 = new Brick.move(-1.0,-1.0,11.313708499) [3].move(1.0,0,0) [3].move(0,1.0,0)
layer17 = new Brick.move(-0.5,-0.5,12.0208152802) [2].move(1.0,0,0) [2].move(0,1.0,0)
layer18 = new Brick.move(0.0,0.0,12.7279220614) [1].move(1.0,0,0) [1].move(0,1.0,0)
}
PyramidQueens1 {
layer0 = new ImmobileBrick.move(-3.5,-3.5,0.0) [8].move(1.0,0,0) [8].move(0,1.0,0)
layer1 = new ImmobileBrick.move(-3.0,-3.0,0.707106781187) [7].move(1.0,0,0) [7].move(0,1.0,0)
layer2 = new ImmobileBrick.move(-2.0,-2.0,1.707106781187) [5].move(1.0,0,0) [5].move(0,1.0,0)
layer3 = new Brick.move(-1.5,-1.5,2.41421356237) [4].move(1.0,0,0) [4].move(0,1.0,0)
layer4 = new Brick.move(-0.5,-0.5,3.41421356237) [2].move(1.0,0,0) [2].move(0,1.0,0)
layer5 = new Brick.move(0.0,0.0,4.12132034356) [1].move(1.0,0,0) [1].move(0,1.0,0)
}
PyramidQueens2 {
layer0 = new ImmobileBrick.move(-3.5,-3.5,0.0) [8].move(1.0,0,0) [8].move(0,1.0,0)
layer1 = new ImmobileBrick.move(-3.0,-3.0,0.707106781187) [7].move(1.0,0,0) [7].move(0,1.0,0)
layer2 = new ImmobileBrick.move(-2.0,-2.0,1.707106781187) [5].move(1.0,0,0) [5].move(0,1.0,0)
layer3 = new Brick.move(-1.5,-1.5,2.41421356237) [4].move(1.0,0,0) [4].move(0,1.0,0)
layer4 = new Brick.move(-0.5,-0.5,3.41421356237) [2].move(1.0,0,0) [2].move(0,1.0,0)
layer5 = new Brick.move(0.0,0.0,4.12132034356) [1].move(1.0,0,0) [1].move(0,1.0,0)
}
PyramidQueens3 {
layer0 = new ImmobileBrick.move(-3.5,-3.5,0.0) [8].move(1.0,0,0) [8].move(0,1.0,0)
layer1 = new Brick.move(-3.0,-3.0,0.707106781187) [7].move(1.0,0,0) [7].move(0,1.0,0)
layer2 = new Brick.move(-2.5,-2.5,1.41421356237) [6].move(1.0,0,0) [6].move(0,1.0,0)
layer3 = new Brick.move(-2.0,-2.0,2.12132034356) [5].move(1.0,0,0) [5].move(0,1.0,0)
layer4 = new Brick.move(-1.5,-1.5,2.82842712475) [4].move(1.0,0,0) [4].move(0,1.0,0)
layer5 = new Brick.move(-1.0,-1.0,3.53553390593) [3].move(1.0,0,0) [3].move(0,1.0,0)
layer6 = new Brick.move(-0.5,-0.5,4.24264068712) [2].move(1.0,0,0) [2].move(0,1.0,0)
layer7 = new Brick.move(0.0,0.0,4.94974746831) [1].move(1.0,0,0) [1].move(0,1.0,0)
}
# "ImmobileBrick"s are identical to "Brick"s,
# except that they have a different atom type.
# We can define groups based on atom type
# and apply fixes to them.
ImmobileBrick {
# atomID molID atomType charge x y z
write("Data Atoms") {
$atom $mol @atom 0.0 0.0 0.0 0.0
}
write_once("Data Masses") {
@atom 1.0
}
write_once("In Settings") {
# U(r) = 4*epsilon_ij*((sigma_ij/r)^12 - (sigma_ij/r)^6)
#
# i j eps sig
pair_coeff @atom @atom lj/cut 1.0 0.8908987181403393
}
write_once("In Settings") {
group gImmobile type @atom
# (Atoms of this type belong to the "gImmobile" group)
}
write_once("In Init") {
atom_style full
units lj
pair_style hybrid lj/cut 1.0
}
}
GoldBrick {
# atomID molID atomType charge x y z
write("Data Atoms") {
$atom $mol @atom 0.0 0.0 0.0 0.0
}
write_once("Data Masses") {
@atom 1.0
}
write_once("In Settings") {
# U(r) = 4*epsilon_ij*((sigma_ij/r)^12 - (sigma_ij/r)^6)
#
# i j eps sig
pair_coeff @atom @atom lj/cut 1.0 0.8908987181403393
}
write_once("In Settings") {
group gMobile type @atom
# (Atoms of this type belong to the "gMobile" group)
}
write_once("In Init") {
atom_style full
units lj
pair_style lj/cut 1.0
}
}

View File

@ -0,0 +1,80 @@
# Description.
# This is a simulation of pyramid-like objects made of particles stacked
# and arranged like cannon-balls, or fruit-stands. Ordinarilly, the stack
# does not collapse because the particles at the ground layer are immobilized.
# However given an initial perterbation the pyramids collapse in an avalanche.
# (This can happen, for example when you do not minimize the system beforehand.)
# The particles roll down the pyramid and bounce off the "ground". The bouncing
# is due to a repulsive external force which is added artificially.
# (See the "run.in" file.) The simulation looks weird without something
# to bounce off of. So I added a graphene surface at the bottom as scenery.
# The ground does not serve any purpose except to look pretty.
#
# (Because there is no damping, I suspect that the distribution of heights of
# the particles in a small area should approach the Boltzmann distribution,
# if you run the simulation long enough.)
# ----------------- Pyramids: -----------------
import "pyramids.lt"
# Move the pyramids into their locations in Giza (approximate)
pyramidKhufu = new PyramidKhufu.move(210, 215, 1)
pyramidKhafre = new PyramidKhafre.move(150, 150, 1)
pyramidMenkaure = new PyramidMenkaure.move(105, 082, 1)
PyramidQueens1 = new PyramidQueens1.move(089, 059, 1)
PyramidQueens2 = new PyramidQueens2.move(100, 059, 1)
PyramidQueens3 = new PyramidQueens3.move(111, 059, 1)
# --------------- Scenery: --------------------
import "graphene_wall.lt"
graphene_wall = new GrapheneWall
write_once("In Settings") {
# Turn off all interactions with the graphene atoms by setting epsilon to 0.
# (We will use a different repulsive barrier to represent the ground instead.)
# These atoms are just "for show". epsilon sigma
pair_coeff @atom:Graphene/C @atom:Graphene/C lj/cut 0.00000 3.407
# Optional: Add the graphene atoms to the "gImmobile" group. Later freeze them
group gImmobile type @atom:Graphene/C
}
# Unfortunately, the ground still looks kind of ugly because moltemplate does
# not yet know how to automatically connect nearby carbon atoms with C-C bonds
# (based on distance). (As of 2013-4-29, moltemplate is not good at
# generating crystalline objects containing explicit bonds.)
# If you want bonds between atoms, use VMD's "carbon-nanotube-builder plugin"
# (which creates data files with bonds) and then merge the two data files
# manually later. (This is not done here.)
# -------- override earlier settings ----------
write_once("In Init") {
# Override any earlier style settings
atom_style full
units lj
pair_style hybrid lj/cut 1.0
bond_style none
angle_style none
dihedral_style none
improper_style none
pair_modify mix arithmetic
special_bonds lj 0.0 0.0 0.0
}
# ------------ boundary conditions ------------
write_once("Data Boundary") {
-1.842033 398.493813 xlo xhi
-0.708994 399.167013 ylo yhi
0.0 400.0 zlo zhi
}
# ---------------------------------------------

View File

@ -0,0 +1,64 @@
# -- Init Section --
include system.in.init
boundary p p f
# -- Atom Definition Section --
read_data system.data
# -- Settings Section --
include system.in.settings
# -- Run Section --
timestep 0.0025
dump 1 all custom 200 traj_nvt.lammpstrj id mol type x y z ix iy iz
thermo_style custom step temp pe etotal
thermo 100 # time interval for printing out "thermo" data
# ---- Set up the physical environment ----
# Add gravity:
fix fxGrav gMobile gravity 0.05 vector 0 0 -1
# Create a "ground" surface.
# This is a repulsive "wall" which particles can bounce off of:
fix fxWall gMobile wall/lj126 zlo EDGE 1.0 0.8908987181403393 1.0
# ---- Evolve the system: ----
# Evolve the (mobile) atoms using ordinary Newton's laws (NVE)
fix fxNVE gMobile nve
# IF YOU WANT TO ADD DAMPING, THEN UNCOMMENT THE NEXT LINE:
#fix fxLan gMobile langevin 0.001 0.001 10000.0 48279
# To use Langevin dynamics in LAMMPS you need both "fix langevin" and "fix nve"
# (See http://lammps.sandia.gov/doc/fix_langevin.html for details.)
# This was not tested.
# For efficient simulation in parallel, try using "fix balance":
# (This will adjust the spatial decomposition as the distribution of
# particles changes over time.)
fix fxBalance gMobile balance 1000 xy 20 1.3
# Optional: Improve efficiency by omitting the calcuation of interactions
# between immobile atoms:
neigh_modify exclude group gImmobile gImmobile
restart 50000 restart_nvt
run 200000
write_data system_after_nvt.data