Unverified Commit 795cdf45 authored by Axel Kohlmeyer's avatar Axel Kohlmeyer
Browse files

simplify example and skip the system generation step, so it gives consistent...

simplify example and skip the system generation step, so it gives consistent results in serial and parallel
parent e481c5f3
Loading
Loading
Loading
Loading
+2 −6
Original line number Diff line number Diff line
The input script in.lammps simulates bulk water using the 2015 E3B potential.
It can be modified to use the 2011 E3B potential by
   1) using a tip4p molecule file instead of tip4p2005.mol
   2) changing the qdist parameter for lj/cut/tip4p/long in the pair_style
   hybrid/overlay command
   3) using the 2011 "preset" command in the e3b pair_coeff command

This script also demonstrates the use of compute pe/e3b to calculate the
potential energy contribution of the e3b pair style. These potential energy
contributions can be found in the output file e3b.txt. See the LAMMPS
+1973 −0

File added.

Preview size limit exceeded, changes collapsed.

+10 −30
Original line number Diff line number Diff line
@@ -2,7 +2,6 @@
#to simulate bulk E3B3 water model

#####################################################################
clear

variable	samp_rate 	equal 10
variable	thermo_rate 	equal 10
@@ -27,20 +26,11 @@ boundary p p p

#############################################################################
#setup box

lattice	  	sc ${Wlat}
region		simbox block -$L $L -$L $L -$L $L units lattice
read_data       e3b_box.data

#############################################################################
#set up potential

create_box	2 simbox bond/types 1 angle/types 1 extra/bond/per/atom 2 &
		extra/special/per/atom 2 extra/angle/per/atom 1

molecule 	h2o tip4p2005.mol
mass            1 15.9994 #oxygen
mass            2 1.008   #hydrogen

pair_style 	hybrid/overlay e3b 1 lj/cut/tip4p/long 1 2 1 1 0.1546 8.5
pair_modify 	table 0 table/disp 0 shift yes

@@ -53,9 +43,9 @@ pair_coeff * * lj/cut/tip4p/long 0.0 0.0
pair_coeff	1 1 lj/cut/tip4p/long 0.1852 3.1589
pair_coeff	* * e3b preset 2015

#intramolecular bond/angle coeffs very stiff for minimization
bond_coeff	1 100000 0.9572
angle_coeff	1 100000 104.52
#potential coeffs aren't too important since will be rigid anyways
bond_coeff	1 554.13 0.9572
angle_coeff	1 45.769 104.52

#############################################################################
#setup for run
@@ -69,25 +59,15 @@ neighbor 2.0 bin
neigh_modify	every 1 delay 3 check yes

#############################################################################
#make atoms and rough minimze
create_atoms	0 box mol h2o 15856

#dump positions only in first batch run
dump 		7 all custom ${samp_rate} dump.lammpstrj id x y z
dump_modify	7 sort id

min_style 	cg
minimize	1.0e-4 1.0e-4 10000 100000

#potential coeffs aren't too important since will be rigid anyways
bond_coeff	1 554.13 0.9572
angle_coeff	1 45.769 104.52

#dump 		7 all custom ${samp_rate} dump.lammpstrj id x y z
#dump_modify	7 sort id

#############################################################################
#initialize velocity and rigid constraint

fix 		rigid all shake 1.0e-8 100 0 b 1 a 1 t 1 2 mol h2o
fix 		rigid all shake 1.0e-8 100 0 b 1 a 1 t 1 2
velocity 	all create ${myT} 15856 dist gaussian rot yes mom yes

#scale velocity
@@ -107,9 +87,9 @@ run ${equil}
#############################################################################
#run at NVT

dump 		1 all custom ${samp_rate} dump.lammpstrj id x y z type
dump_modify	1 sort id
#dump 		1 all custom ${samp_rate} dump.lammpstrj id x y z type
#dump_modify	1 sort id

run 		${run}

write_restart 	lammps.restart
# write_restart 	lammps.restart
+332 −0
Original line number Diff line number Diff line
LAMMPS (29 Mar 2019)
#LAMMPS input file
#to simulate bulk E3B3 water model

#####################################################################

variable	samp_rate 	equal 10
variable	thermo_rate 	equal 10
variable	Wlat		equal 3.10744  #for water density 0.997g/mL
variable        L		equal 3  #(L*2)^3 = Nmolec, L=3 -> N=216

variable	equil		equal 100
variable	run		equal 100

variable	ts		equal 2.0
variable	Tdamp		equal 100*${ts}
variable	Tdamp		equal 100*2
variable 	Pdamp		equal 1000*${ts}
variable 	Pdamp		equal 1000*2
variable	myT		equal 298.0
variable	myP		equal 1.0

units 		real
atom_style	full

dimension 	3

boundary 	p p p

#############################################################################
#setup box
read_data       e3b_box.data
  orthogonal box = (-9.32232 -9.32232 -9.32232) to (9.32232 9.32232 9.32232)
  1 by 1 by 1 MPI processor grid
  reading atoms ...
  648 atoms
  reading velocities ...
  648 velocities
  scanning bonds ...
  2 = max bonds/atom
  scanning angles ...
  1 = max angles/atom
  reading bonds ...
  432 bonds
  reading angles ...
  216 angles
  2 = max # of 1-2 neighbors
  1 = max # of 1-3 neighbors
  1 = max # of 1-4 neighbors
  2 = max # of special neighbors
  special bonds CPU = 0.000257254 secs
  read_data CPU = 0.00286555 secs

#############################################################################
#set up potential

pair_style 	hybrid/overlay e3b 1 lj/cut/tip4p/long 1 2 1 1 0.1546 8.5
pair_modify 	table 0 table/disp 0 shift yes

bond_style  	harmonic
angle_style 	harmonic

kspace_style 	pppm/tip4p 1.0e-6

pair_coeff 	* * lj/cut/tip4p/long 0.0 0.0
pair_coeff	1 1 lj/cut/tip4p/long 0.1852 3.1589
pair_coeff	* * e3b preset 2015

#potential coeffs aren't too important since will be rigid anyways
bond_coeff	1 554.13 0.9572
angle_coeff	1 45.769 104.52

#############################################################################
#setup for run
thermo 		${thermo_rate}
thermo 		10
thermo_style 	custom step vol temp epair pe etotal press density

timestep 	${ts}
timestep 	2
run_style 	verlet

neighbor	2.0 bin
neigh_modify	every 1 delay 3 check yes

#############################################################################

#dump positions only in first batch run
#dump 		7 all custom ${samp_rate} dump.lammpstrj id x y z
#dump_modify	7 sort id

#############################################################################
#initialize velocity and rigid constraint

fix 		rigid all shake 1.0e-8 100 0 b 1 a 1 t 1 2
  0 = # of size 2 clusters
  0 = # of size 3 clusters
  0 = # of size 4 clusters
  216 = # of frozen angles
  find clusters CPU = 0.000185728 secs
velocity 	all create ${myT} 15856 dist gaussian rot yes mom yes
velocity 	all create 298 15856 dist gaussian rot yes mom yes

#scale velocity
run		0
PPPM initialization ...
  extracting TIP4P info from pair style
  using polynomial approximation for long-range coulomb (../kspace.cpp:319)
  G vector (1/distance) = 0.409658
  grid = 36 36 36
  stencil order = 5
  estimated absolute RMS force accuracy = 0.000341883
  estimated relative force accuracy = 1.02957e-06
  using double precision FFTs
  3d grid and FFT values/proc = 91125 46656
Neighbor list info ...
  update every 1 steps, delay 3 steps, check yes
  max neighbors/atom: 2000, page size: 100000
  master list distance cutoff = 10.8092
  ghost atom cutoff = 10.8092
  binsize = 5.4046, bins = 4 4 4
  2 neighbor lists, perpetual/occasional/extra = 2 0 0
  (1) pair e3b, perpetual, skip from (2)
      attributes: half, newton on
      pair build: skip
      stencil: none
      bin: none
  (2) pair lj/cut/tip4p/long, perpetual
      attributes: half, newton on
      pair build: half/bin/newton
      stencil: half/bin/3d/newton
      bin: standard
Per MPI rank memory allocation (min/avg/max) = 16.4 | 16.4 | 16.4 Mbytes
Step Volume Temp E_pair PotEng TotEng Press Density 
       0    6481.2982          298   -512.62432   -512.62432   -129.77504    14088.322   0.99697602 
Loop time of 1.19209e-06 on 1 procs for 0 steps with 648 atoms

251.7% CPU use with 1 MPI tasks x no OpenMP threads

MPI task timing breakdown:
Section |  min time  |  avg time  |  max time  |%varavg| %total
---------------------------------------------------------------
Pair    | 0          | 0          | 0          |   0.0 |  0.00
Bond    | 0          | 0          | 0          |   0.0 |  0.00
Kspace  | 0          | 0          | 0          |   0.0 |  0.00
Neigh   | 0          | 0          | 0          |   0.0 |  0.00
Comm    | 0          | 0          | 0          |   0.0 |  0.00
Output  | 0          | 0          | 0          |   0.0 |  0.00
Modify  | 0          | 0          | 0          |   0.0 |  0.00
Other   |            | 1.192e-06  |            |       |100.00

Nlocal:    648 ave 648 max 648 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost:    5943 ave 5943 max 5943 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs:    18984 ave 18984 max 18984 min
Histogram: 1 0 0 0 0 0 0 0 0 0

Total # of neighbors = 18984
Ave neighs/atom = 29.2963
Ave special neighs/atom = 2
Neighbor list builds = 0
Dangerous builds = 0
velocity	all scale ${myT}
velocity	all scale 298

compute    	e3b all pe/e3b
fix		e3b all ave/time 1 1 ${thermo_rate} c_e3b c_e3b[*] 		file e3b.txt title2 "step pe_e3b pe_ea pe_eb pe_ec pe_e2"
fix		e3b all ave/time 1 1 10 c_e3b c_e3b[*] 		file e3b.txt title2 "step pe_e3b pe_ea pe_eb pe_ec pe_e2"

#############################################################################
#equilibrate bulk water at NVT

fix 		1 all nvt temp ${myT} ${myT} ${Tdamp}
fix 		1 all nvt temp 298 ${myT} ${Tdamp}
fix 		1 all nvt temp 298 298 ${Tdamp}
fix 		1 all nvt temp 298 298 200
run 		${equil}
run 		100
PPPM initialization ...
  extracting TIP4P info from pair style
  using polynomial approximation for long-range coulomb (../kspace.cpp:319)
  G vector (1/distance) = 0.409658
  grid = 36 36 36
  stencil order = 5
  estimated absolute RMS force accuracy = 0.000341883
  estimated relative force accuracy = 1.02957e-06
  using double precision FFTs
  3d grid and FFT values/proc = 91125 46656
Neighbor list info ...
  update every 1 steps, delay 3 steps, check yes
  max neighbors/atom: 2000, page size: 100000
  master list distance cutoff = 10.8092
  ghost atom cutoff = 10.8092
  binsize = 5.4046, bins = 4 4 4
  2 neighbor lists, perpetual/occasional/extra = 2 0 0
  (1) pair e3b, perpetual, skip from (2)
      attributes: half, newton on
      pair build: skip
      stencil: none
      bin: none
  (2) pair lj/cut/tip4p/long, perpetual
      attributes: half, newton on
      pair build: half/bin/newton
      stencil: half/bin/3d/newton
      bin: standard
Per MPI rank memory allocation (min/avg/max) = 16.88 | 16.88 | 16.88 Mbytes
Step Volume Temp E_pair PotEng TotEng Press Density 
       0    6481.2982          298   -512.33264   -512.33264   -129.48336    14091.383   0.99697602 
      10    6481.2982    905.58028   -1442.0378   -1442.0378   -278.61245   -798.38952   0.99697602 
      20    6481.2982    816.39844   -1363.5999   -1363.5999   -314.74903    3023.9064   0.99697602 
      30    6481.2982     783.3897   -1370.6594   -1370.6594   -364.21587    7095.4765   0.99697602 
      40    6481.2982    793.12519   -1425.8404   -1425.8404   -406.88933     7030.242   0.99697602 
      50    6481.2982    810.90264   -1495.4822   -1495.4822   -453.69195    6944.2325   0.99697602 
      60    6481.2982    766.64937   -1491.2317   -1491.2317   -506.29493    9062.0151   0.99697602 
      70    6481.2982    761.77292   -1538.7368   -1538.7368   -560.06492    7693.2197   0.99697602 
      80    6481.2982    730.44938   -1554.1818   -1554.1818   -615.75215    7345.9601   0.99697602 
      90    6481.2982    695.46244   -1563.9869   -1563.9869   -670.50605    7809.0685   0.99697602 
     100    6481.2982     691.5674   -1613.2754   -1613.2754   -724.79866     7143.062   0.99697602 
Loop time of 1.94577 on 1 procs for 100 steps with 648 atoms

Performance: 8.881 ns/day, 2.702 hours/ns, 51.393 timesteps/s
99.7% CPU use with 1 MPI tasks x no OpenMP threads

MPI task timing breakdown:
Section |  min time  |  avg time  |  max time  |%varavg| %total
---------------------------------------------------------------
Pair    | 0.91572    | 0.91572    | 0.91572    |   0.0 | 47.06
Bond    | 6.7234e-05 | 6.7234e-05 | 6.7234e-05 |   0.0 |  0.00
Kspace  | 0.92654    | 0.92654    | 0.92654    |   0.0 | 47.62
Neigh   | 0.087331   | 0.087331   | 0.087331   |   0.0 |  4.49
Comm    | 0.0054724  | 0.0054724  | 0.0054724  |   0.0 |  0.28
Output  | 0.0003047  | 0.0003047  | 0.0003047  |   0.0 |  0.02
Modify  | 0.0093319  | 0.0093319  | 0.0093319  |   0.0 |  0.48
Other   |            | 0.001007   |            |       |  0.05

Nlocal:    648 ave 648 max 648 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost:    5911 ave 5911 max 5911 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs:    19136 ave 19136 max 19136 min
Histogram: 1 0 0 0 0 0 0 0 0 0

Total # of neighbors = 19136
Ave neighs/atom = 29.5309
Ave special neighs/atom = 2
Neighbor list builds = 15
Dangerous builds = 0

#############################################################################
#run at NVT

#dump 		1 all custom ${samp_rate} dump.lammpstrj id x y z type
#dump_modify	1 sort id

run 		${run}
run 		100
PPPM initialization ...
  extracting TIP4P info from pair style
  using polynomial approximation for long-range coulomb (../kspace.cpp:319)
  G vector (1/distance) = 0.409658
  grid = 36 36 36
  stencil order = 5
  estimated absolute RMS force accuracy = 0.000341883
  estimated relative force accuracy = 1.02957e-06
  using double precision FFTs
  3d grid and FFT values/proc = 91125 46656
Neighbor list info ...
  update every 1 steps, delay 3 steps, check yes
  max neighbors/atom: 2000, page size: 100000
  master list distance cutoff = 10.8092
  ghost atom cutoff = 10.8092
  binsize = 5.4046, bins = 4 4 4
  2 neighbor lists, perpetual/occasional/extra = 2 0 0
  (1) pair e3b, perpetual, skip from (2)
      attributes: half, newton on
      pair build: skip
      stencil: none
      bin: none
  (2) pair lj/cut/tip4p/long, perpetual
      attributes: half, newton on
      pair build: half/bin/newton
      stencil: half/bin/3d/newton
      bin: standard
Per MPI rank memory allocation (min/avg/max) = 16.88 | 16.88 | 16.88 Mbytes
Step Volume Temp E_pair PotEng TotEng Press Density 
     100    6481.2982     691.5674   -1613.2754   -1613.2754   -724.79866    7131.3961   0.99697602 
     110    6481.2982    668.27004   -1635.5867   -1635.5867   -777.04068    6965.8705   0.99697602 
     120    6481.2982    646.18686   -1659.5025   -1659.5025   -829.32738    6196.7432   0.99697602 
     130    6481.2982     650.7802   -1716.9557   -1716.9557   -880.87943    5626.5466   0.99697602 
     140    6481.2982      586.262   -1681.8052   -1681.8052   -928.61728     6103.747   0.99697602 
     150    6481.2982    615.88299   -1767.8614   -1767.8614   -976.61859     3897.648   0.99697602 
     160    6481.2982    585.23516   -1773.2038   -1773.2038   -1021.3352    4821.7742   0.99697602 
     170    6481.2982    558.77885   -1782.2817   -1782.2817   -1064.4022    5092.7248   0.99697602 
     180    6481.2982    564.01576   -1830.0301   -1830.0301   -1105.4226    4316.1636   0.99697602 
     190    6481.2982    526.53776   -1821.3122   -1821.3122   -1144.8538    4529.5062   0.99697602 
     200    6481.2982    537.81273   -1873.6662   -1873.6662   -1182.7226     4244.313   0.99697602 
Loop time of 1.9157 on 1 procs for 100 steps with 648 atoms

Performance: 9.020 ns/day, 2.661 hours/ns, 52.200 timesteps/s
99.6% CPU use with 1 MPI tasks x no OpenMP threads

MPI task timing breakdown:
Section |  min time  |  avg time  |  max time  |%varavg| %total
---------------------------------------------------------------
Pair    | 0.90464    | 0.90464    | 0.90464    |   0.0 | 47.22
Bond    | 6.9141e-05 | 6.9141e-05 | 6.9141e-05 |   0.0 |  0.00
Kspace  | 0.92769    | 0.92769    | 0.92769    |   0.0 | 48.43
Neigh   | 0.067551   | 0.067551   | 0.067551   |   0.0 |  3.53
Comm    | 0.0051386  | 0.0051386  | 0.0051386  |   0.0 |  0.27
Output  | 0.00029993 | 0.00029993 | 0.00029993 |   0.0 |  0.02
Modify  | 0.0092504  | 0.0092504  | 0.0092504  |   0.0 |  0.48
Other   |            | 0.001062   |            |       |  0.06

Nlocal:    648 ave 648 max 648 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost:    5901 ave 5901 max 5901 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs:    18971 ave 18971 max 18971 min
Histogram: 1 0 0 0 0 0 0 0 0 0

Total # of neighbors = 18971
Ave neighs/atom = 29.2762
Ave special neighs/atom = 2
Neighbor list builds = 12
Dangerous builds = 0

# write_restart 	lammps.restart

Please see the log.cite file for references relevant to this simulation

Total wall time: 0:00:03
+332 −0

File added.

Preview size limit exceeded, changes collapsed.

Loading