git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@10359 f3b2605a-c512-4ea7-a41b-209d697bcdaa
This commit is contained in:
sjplimp
33
tools/moltemplate/examples/all_atom_examples/ice_crystal/README_run.sh
Executable file
33
tools/moltemplate/examples/all_atom_examples/ice_crystal/README_run.sh
Executable file
@ -0,0 +1,33 @@
|
||||
# --- Running LAMMPS ---
|
||||
# -- Prerequisites: --
|
||||
# The 2 files "run.in.npt", and "run.in.nvt" are LAMMPS
|
||||
# input scripts which link to the input scripts and data files
|
||||
# you hopefully have created earlier with moltemplate.sh:
|
||||
# system.in.init, system.in.settings, system.data
|
||||
# If not, carry out the instructions in "README_setup.sh".
|
||||
#
|
||||
# -- Instructions: --
|
||||
# If "lmp_linux" is the name of the command you use to invoke lammps,
|
||||
# then you would run lammps on these files this way:
|
||||
|
||||
|
||||
lmp_linux -i run.in.npt # minimization and simulation at constant pressure
|
||||
|
||||
# or
|
||||
|
||||
lmp_linux -i run.in.nvt # minimization and simulation at constant volume
|
||||
|
||||
#(Note: The constant volume simulation lacks pressure equilibration. These are
|
||||
# completely separate simulations. The results of the constant pressure
|
||||
# simulation are ignored when beginning the simulation at constant volume.
|
||||
# This can be fixed. Read "run.in.nvt" for equilibration instructions.)
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
# If you have compiled the MPI version of lammps, you can run lammps in parallel
|
||||
#mpirun -np 4 lmp_linux -i run.in.npt
|
||||
# or
|
||||
#mpirun -np 4 lmp_linux -i run.in.nvt
|
||||
# (assuming you have 4 processors available)
|
||||
23
tools/moltemplate/examples/all_atom_examples/ice_crystal/README_setup.sh
Executable file
23
tools/moltemplate/examples/all_atom_examples/ice_crystal/README_setup.sh
Executable 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 -atomstyle full 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 ../
|
||||
@ -0,0 +1,87 @@
|
||||
|
||||
------- 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 500 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 ----
|
||||
Sometimes the solvent or membrane obscures the view of the solute.
|
||||
It can help to shift the location of the periodic boundary box
|
||||
To shift the box in the y direction (for example) do this:
|
||||
|
||||
pbc wrap -compound res -all -shiftcenterrel {0.0 0.15 0.0}
|
||||
pbc box -shiftcenterrel {0.0 0.15 0.0}
|
||||
|
||||
Distances are measured in units of box-length fractions, not Angstroms.
|
||||
|
||||
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)
|
||||
You should check if your periodic boundary conditions are too small.
|
||||
To do that:
|
||||
select Graphics->Representations menu option
|
||||
click on the "Periodic" tab, and
|
||||
click on the "+x", "-x", "+y", "-y", "+z", "-z" checkboxes.
|
||||
|
||||
5) 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: 25 KiB |
@ -0,0 +1,52 @@
|
||||
# file "spce.lt"
|
||||
#
|
||||
# H1 H2
|
||||
# \ /
|
||||
# O
|
||||
|
||||
SPCE {
|
||||
|
||||
write_once("In Init") {
|
||||
# -- Default styles (for solo "SPCE" water) --
|
||||
units real
|
||||
atom_style full
|
||||
# (Hybrid force fields were not necessary but are used for portability.)
|
||||
pair_style hybrid lj/charmm/coul/long 9.0 10.0 10.0
|
||||
bond_style hybrid harmonic
|
||||
angle_style hybrid harmonic
|
||||
kspace_style pppm 0.0001
|
||||
pair_modify mix arithmetic
|
||||
}
|
||||
|
||||
write("Data Atoms") {
|
||||
$atom:O $mol:. @atom:O -0.8476 0.0000000 0.00000 0.000000
|
||||
$atom:H1 $mol:. @atom:H 0.4238 0.8164904 0.00000 0.5773590
|
||||
$atom:H2 $mol:. @atom:H 0.4238 -0.8164904 0.00000 0.5773590
|
||||
}
|
||||
|
||||
write_once("Data Masses") {
|
||||
@atom:O 15.9994
|
||||
@atom:H 1.008
|
||||
}
|
||||
|
||||
write("Data Bonds") {
|
||||
$bond:OH1 @bond:OH $atom:O $atom:H1
|
||||
$bond:OH2 @bond:OH $atom:O $atom:H2
|
||||
}
|
||||
|
||||
write("Data Angles") {
|
||||
$angle:HOH @angle:HOH $atom:H1 $atom:O $atom:H2
|
||||
}
|
||||
|
||||
write_once("In Settings") {
|
||||
bond_coeff @bond:OH harmonic 1000.0 1.0
|
||||
angle_coeff @angle:HOH harmonic 1000.0 109.47
|
||||
pair_coeff @atom:O @atom:O lj/charmm/coul/long 0.1553 3.166
|
||||
pair_coeff @atom:H @atom:H lj/charmm/coul/long 0.0 2.058
|
||||
group spce type @atom:O @atom:H
|
||||
fix fSHAKE spce shake 0.0001 10 100 b @bond:OH a @angle:HOH
|
||||
# (Remember to "unfix" fSHAKE during minimization.)
|
||||
}
|
||||
|
||||
} # end of definition of "SPCE" water molecule type
|
||||
|
||||
@ -0,0 +1,57 @@
|
||||
# This ice (1h) unit cell is rectangular and contains 8 water molecules.
|
||||
# (Coordinates and cell dimensions converted were from a PDB file.)
|
||||
# The dimensions of the unit cell (in Angstroms) are: 4.521 7.832 7.362
|
||||
|
||||
|
||||
import "spce.lt" # <-- define the "SPCE" molecule
|
||||
|
||||
SpceIceRect8 {
|
||||
|
||||
# Create a 3-dimensional array of 8 water molecules
|
||||
|
||||
wat = new SPCE[2][2][2]
|
||||
|
||||
# Array indices will be correlated with position [xindex][yindex][zindex]
|
||||
|
||||
# You can overwrite coordinates of atoms after they were created this way:
|
||||
# (Order is not important)
|
||||
# atom-ID molecule-ID atomType charge newX newY newZ
|
||||
|
||||
write("Data Atoms") {
|
||||
$atom:wat[1][0][0]/O $mol:wat[1][0][0] @atom:SPCE/O -0.8476 3.391 1.305 1.381
|
||||
$atom:wat[1][0][0]/H1 $mol:wat[1][0][0] @atom:SPCE/H 0.4238 3.391 0.370 1.710
|
||||
$atom:wat[1][0][0]/H2 $mol:wat[1][0][0] @atom:SPCE/H 0.4238 2.582 1.772 1.710
|
||||
$atom:wat[1][0][1]/O $mol:wat[1][0][1] @atom:SPCE/O -0.8476 3.391 1.305 5.981
|
||||
$atom:wat[1][0][1]/H1 $mol:wat[1][0][1] @atom:SPCE/H 0.4238 3.391 1.305 6.970
|
||||
$atom:wat[1][0][1]/H2 $mol:wat[1][0][1] @atom:SPCE/H 0.4238 4.200 1.772 5.652
|
||||
$atom:wat[0][0][0]/O $mol:wat[0][0][0] @atom:SPCE/O -0.8476 1.131 2.611 2.300
|
||||
$atom:wat[0][0][0]/H1 $mol:wat[0][0][0] @atom:SPCE/H 0.4238 1.131 2.611 3.289
|
||||
$atom:wat[0][0][0]/H2 $mol:wat[0][0][0] @atom:SPCE/H 0.4238 0.320 2.143 1.971
|
||||
$atom:wat[0][0][1]/O $mol:wat[0][0][1] @atom:SPCE/O -0.8476 1.131 2.611 5.061
|
||||
$atom:wat[0][0][1]/H1 $mol:wat[0][0][1] @atom:SPCE/H 0.4238 1.940 2.143 5.391
|
||||
$atom:wat[0][0][1]/H2 $mol:wat[0][0][1] @atom:SPCE/H 0.4238 1.131 3.546 5.391
|
||||
$atom:wat[0][1][0]/O $mol:wat[0][1][0] @atom:SPCE/O -0.8476 1.131 5.221 1.381
|
||||
$atom:wat[0][1][0]/H1 $mol:wat[0][1][0] @atom:SPCE/H 0.4238 1.131 4.286 1.710
|
||||
$atom:wat[0][1][0]/H2 $mol:wat[0][1][0] @atom:SPCE/H 0.4238 0.320 5.688 1.710
|
||||
$atom:wat[0][1][1]/O $mol:wat[0][1][1] @atom:SPCE/O -0.8476 1.131 5.221 5.981
|
||||
$atom:wat[0][1][1]/H1 $mol:wat[0][1][1] @atom:SPCE/H 0.4238 1.131 5.221 6.970
|
||||
$atom:wat[0][1][1]/H2 $mol:wat[0][1][1] @atom:SPCE/H 0.4238 1.940 5.688 5.652
|
||||
$atom:wat[1][1][0]/O $mol:wat[1][1][0] @atom:SPCE/O -0.8476 3.391 6.526 2.300
|
||||
$atom:wat[1][1][0]/H1 $mol:wat[1][1][0] @atom:SPCE/H 0.4238 3.391 6.526 3.289
|
||||
$atom:wat[1][1][0]/H2 $mol:wat[1][1][0] @atom:SPCE/H 0.4238 2.582 6.058 1.971
|
||||
$atom:wat[1][1][1]/O $mol:wat[1][1][1] @atom:SPCE/O -0.8476 3.391 6.526 5.061
|
||||
$atom:wat[1][1][1]/H1 $mol:wat[1][1][1] @atom:SPCE/H 0.4238 4.200 6.058 5.391
|
||||
$atom:wat[1][1][1]/H2 $mol:wat[1][1][1] @atom:SPCE/H 0.4238 3.391 7.462 5.391
|
||||
}
|
||||
} # IceRect8
|
||||
|
||||
# Credit goes to Martin Chaplin.
|
||||
# These coordinates were orignally downloaded from Martin Chaplin's
|
||||
# website: http://www.btinternet.com/~martin.chaplin/ice1h.html
|
||||
# ... and then they were stretched independently in the xy and z
|
||||
# directions in order to match the lattice parameters measured by
|
||||
# Rottger et al.,
|
||||
# "Lattice constants and thermal expansion of H2O and D2O ice Ih"
|
||||
# between 10 and 265K", Acta Crystallogr. B, 50 (1994) 644-648
|
||||
# I am using the lattice constants measured at temperature 265K
|
||||
# (and pressure=100Torr).
|
||||
@ -0,0 +1,11 @@
|
||||
import "spce_ice_rect8.lt"
|
||||
|
||||
cells = new SpceIceRect8 [3].move(4.521, 0.0, 0.0)
|
||||
[2].move( 0.0, 7.832, 0.0)
|
||||
[2].move( 0.0, 0.0, 7.362)
|
||||
|
||||
write_once("Data Boundary") {
|
||||
0 13.563 xlo xhi
|
||||
0 15.664 ylo yhi
|
||||
0 14.724 zlo zhi
|
||||
}
|
||||
@ -0,0 +1,38 @@
|
||||
# ------------------------------- Initialization Section --------------------
|
||||
|
||||
include system.in.init
|
||||
|
||||
# ------------------------------- Atom Definition Section -------------------
|
||||
|
||||
read_data system.data
|
||||
|
||||
# ------------------------------- Settings Section --------------------------
|
||||
|
||||
include system.in.settings
|
||||
|
||||
# ------------------------------- Run Section -------------------------------
|
||||
|
||||
|
||||
# -- minimization protocol --
|
||||
|
||||
# Note: The minimization step is not necessary in this example. However
|
||||
# in general, it's always a good idea to minimize the system beforehand.
|
||||
# fSHAKE was defined in system.in.settings. It is incompatible with "minimize".
|
||||
unfix fSHAKE
|
||||
minimize 1.0e-5 1.0e-7 100000 400000
|
||||
# Now read "system.in.settings" in order to redefine fSHAKE again:
|
||||
include system.in.settings
|
||||
|
||||
# -- simulation protocol --
|
||||
|
||||
|
||||
timestep 2.0
|
||||
dump 1 all custom 200 traj_npt.lammpstrj id mol type x y z ix iy iz
|
||||
fix fxnpt all npt temp 400.0 400.0 100.0 iso 1.0 1.0 1000.0 drag 1.0
|
||||
|
||||
thermo 100
|
||||
#thermo_modify flush yes
|
||||
|
||||
run 20000
|
||||
|
||||
write_restart system_after_npt.rst
|
||||
@ -0,0 +1,61 @@
|
||||
# ------------------------------- Initialization Section --------------------
|
||||
|
||||
include system.in.init
|
||||
|
||||
# ------------------------------- Atom Definition Section -------------------
|
||||
|
||||
|
||||
# Normally, I would minimize the system and equilibrate the system at constant
|
||||
# pressure and temperature beforehand. If you run lammps with "run.in.npt",
|
||||
# it will generate a restart file "system_after_npt.rst" with reasonable
|
||||
# coordinates at that temperature and pressure. Then we could load it now:
|
||||
#
|
||||
#read_restart system_after_npt.rst
|
||||
#
|
||||
# Unfortunately the LAMMPS "read_restart" command has been undependable over
|
||||
# the past year (2012), and I feel it is safer to remove it from the examples.
|
||||
# Instead, for this example, I just read the raw coordinates generated by
|
||||
# moltemplate (and the default volume). (I get fewer questions this way.)
|
||||
# However you should never run any liquid simulations at constant volume without
|
||||
# pressure equilibration first. Hopefully in the future "read_restart" will
|
||||
# work. Until then, try "read_dump", "dump2data.py", or "restart2data".
|
||||
|
||||
read_data system.data
|
||||
|
||||
# ------------------------------- Settings Section --------------------------
|
||||
|
||||
include system.in.settings
|
||||
|
||||
# ------------------------------- Run Section -------------------------------
|
||||
|
||||
|
||||
# -- minimization protocol --
|
||||
|
||||
# Note: The minimization step is not necessary in this example. However
|
||||
# in general, it's always a good idea to minimize the system beforehand.
|
||||
# fSHAKE was defined in system.in.settings. It is incompatible with "minimize".
|
||||
unfix fSHAKE
|
||||
minimize 1.0e-5 1.0e-7 100000 400000
|
||||
# Now read "system.in.settings" in order to redefine fSHAKE again:
|
||||
include system.in.settings
|
||||
|
||||
# -- simulation protocol --
|
||||
|
||||
|
||||
|
||||
timestep 2.0
|
||||
dump 1 all custom 500 traj_nvt.lammpstrj id mol type x y z ix iy iz
|
||||
fix fxnvt all nvt temp 300.0 300.0 500.0 tchain 1
|
||||
|
||||
|
||||
thermo 500
|
||||
#thermo_modify flush yes
|
||||
|
||||
|
||||
#restart 100000 restart_nvt
|
||||
|
||||
run 50000
|
||||
|
||||
write_restart system_after_nvt.rst
|
||||
|
||||
|
||||
Reference in New Issue
Block a user