CLAssical Many-Particle Simulator
Full description of Keywords
This is the keyword which defines a type of identical particles labeled with name.
Number of particles, default 1.
Mass of particle default 1. Only used for Molecular dynamics.
Charge of this particle in units of E1/2L1/2 . Ewald sums will be used if this is nonzero but note comment about "couple" below. The charge-charge interaction is in addition to the other potentials given by the POTENTIAL keyword. The charge parameter is also used to compute:
Sk = <|rk|2 where rk = Si qi exp(ikri)
The inverse of the friction coefficient used in Langevin dynamics and Brownian dynamics. Default is the mass.
Spatial dimensionality. Range 1 to 3. Default is 3.
Size of the simulation box in the x, y and z directions. Actual coordinates are in the range (- L2 to L/2). Default for ly is lx. Default for lz is ly.
Total number of particles per unit volume (or area in 2D). Only used if BOX_SIZE is absent.
The interaction between a pair of particles of types namei and namej separated by a distance r= |ri-rj| will be:
epsilon*(F(type,r/sigma)-F(type,cutoff/sigma))
type = potential function number. Default is 1
sigma = length scale of potential. Default is 1.
epsilon = energy scale of potential. Default is 1.
cutoff is the distance outside of which the potential is set to zero (smoothly). Default is half the box size.
If cutoff is positive, a tail correction will be added to the potential energy to account for the particles outside of the potential. If you do not want this tail correction to be added set cutoff to a negative value=- ABS(cutoff). If cutoff=0 (default) it will be set to L/2 for that pair of types which first use a particular function number. The succeeding pairs of types using that function and having cutoff=0 will use the same potential tables and hence have the same dimensionless cutoff. These should be checked on the output to see if they are less than L/2. If not the defaults should be overridden.
The various potential functions are:
1. 4*(1/x12-1/x6) (default)
2. 1/x
3. ERFC(x)/x
4. e-x/x
5. (e-2x- e-x)/x
6-8. refer to the Rahman, Stillinger and Lemberg central force model for water J. Chem. Phys. 63,5225(1975). They are respectively the HH, OH and OO interactions.
Value of temperature (units either C or K) as specified on units. No default. You must specify the temperature.
The program uses units for the following purposes
By default these conversions will not be done. This means that the user must calculate the correct relationship amongst the input quantities so that energy=mass*(length/time)**2 and temperature is in units of energy.
Creates bonding potentials between the "named" particles. The interaction between two bonded particles on a chain will be:
v(r) = -SPRING r02LOG(l.-(r/r0)2)/2
spring constant is in units of energy/length**2 with a default of 1.0.
r0 is maximum separation of two bonded particles. Default is 105
CLAMPS will stop if the total number of monomers of some chemical type exceeds that given in the type declaration. As currently implemented, both r0 and spring_constant are independent of the chemical type and only particles of the same type can be connected by a polymer potential.
This tells CLAMPS how many bins in which to divide the simulation box to speed the computation for systems with more than several hundred particles. The bins are used to construct a neighbor list and thus avoid computing pair distances for particles outside the potential cutoff. The program can select the optimal number of bins itself so this keyword can be absent. However, a user may wish to compute all pair distances, or for a very large system a user may wish to restrict the amount of memory the bins use. The program will examine dividing into between min and max divisions in each dimension. The defaults are 1 to 15. When it sees the computer time per step increasing it stops searching and uses the previous division. MD uses only odd bin numbers.
Sizeof the force and potential tables with a default of 1000. One should make some tests runs to determine how it effects the potential energy, particularly with Molecular Dynamics (e.g., a value of up to 100,000 may be needed). The memory used by these tables is 2#entries ntypes2.
dxpot is the increment in r of potential and force tables printed out at the beginning of the program. If dxpot is absent there will be no printing of these tables.
The EWALD keyword controls the accuracy of the image potential for charged systems. #k_shells is the distance in k space at which to truncate the fourier space coefficients for the Ewald sums and for which S(k) is computed.The largest value will be of order (k_shells,0,0) The Ewald summation technique used in CLAMPS is described in D. Ceperley, Phys. Rev. B18, 3126 (1978) appendix. Note that CPU time will have a factor roughly k_shellsndim, so it should not be made too large.
couple is a parameter which multiplies the contribution of the fourier coefficients to the energy. Default is 1. If set to zero then the k_shells are computed for S(k) but do not change the energy.
LATTICE will initialize the particle coordinates to a lattice.
The crystal_type is an integer rang (1-4) corresponding to:
1 Simple cubic lattice 3D/ square 2D
2 body centered cubic lattice 3D / square rotated by 45o 2D
3 hexagonal close packed lattice 3D / triangular lattice 2D Note that the sides of the simulation box will not be equal.
4 face centered cubic lattice 3D / not defined 2D
Default for crystal type is to choose the lattice type leading to the fewest number of vacancies. The vacancies, if any, are randomly distributed among the crystal sites without changing the random number sequence.
ncellx is the number of units cells in the x direction. CLAMPS computes ncellx*ncelly*ncelly*NXTAL. If this is less than the total number of particles, ncell is increased for NXTAL =1,2,4 and roughly cubic for NXTAL=3.
This keyword can be used to select a noncubic box. For example, the input line
LATTICE 2 5 5 10
will give a BCC lattice of 2*5*5*10 = 500 points as long as the number of particles is less than or equal to 500 The box side in the z direction will be twice the sides in the x and y direction. Note that sides are always selected so that the number of particles per unit volume equals the density (which must have been specified).
Sets the positions of the particles to those in file_name (see format.) extension of file_name must be .crd Default file_name is id.crd. Velocities will be initialized to a Maxwell-Boltzmann distribution. No checking of particle type or number is done.
Sets the internal state (positions, velocities, random number seed, ..) to those in filename. Default file_name is id.chk. This will supercede all initializations.
Begin with last configuration generated by previous run and continue as if job had not concluded. Restarting is appropriate if run has been interrupted before completion or if you want to extend the length of a run. The only changes allowed to id.in are to increase the number of steps or to add more commands "RUN xx" after the last executed step. CHECKPOINT files must have been written for RESTART to be successful. The RESTART command can be placed anywhere in the .in file and it will supercede all other state initializations. The difference between RESTART and READ_STATE is the RESTART will skip steps already performed.
Random number seed. (integer) See the SPRNG documentation. The default is a random seed constructed from the date and time. If you wish to repeat a run, and you chose the default for the initial run use the number printed on the output listing. To change the generator type, the Makefile must be changed.
Writes the particle coordinates onto id.crd with this period. The total number of writes generated by this command will be nsteps/period (each containing natoms lines) period=0 will turn off writes (default).
Period with which to write scalars for analysis onto the file id.sca. The default is not to write out scalars. This file contains the energy, acceptance ratio and computer time/step.
Period with which to output the structure factor to file id.sk. See TYPE, EWALD and READ_K for the definition of S(k) and how to input the k-values. Default is not to write out S(k).
Reads filename for additional k points to compute S(k) for. This can be used for monitoring the melting of a crystal. The first line of the file are the number of distinct shells and number of k-vectors. Then for the first shell one has the number of k-vectors in that shell. Following that are the xyz coordinates of the k-vector. Examples shown are for a 43 fcc lattice and a triangle lattice in 2D.
Period to spill out coordinates and averages onto file id.chk The checkpoint file is always written at the end of a "RUN" command, but if this command is absent, it will be impossible to restart from a partially complete "RUN".
Will cause more printouts on the id.out file.
These keywords cause the simulation to begin. They execute in the order in the input file. If RESTART is present, the run will begin in the last CHECKPOINTed state.
nsteps refers to the total number of steps performed.
Molecular Dynamics; Newton's equations of motion.
tau is the time step
With tbath nonzero a random force and viscous drag, characteristic of a heat bath will be added to the force (Langevin equation). The friction coefficient for a particle will be 1/(tbath xinv). (not currently implemented)
The velocities for molecular dynamics are rescaled so that the kinetic energy equals (1/2)KBT.
If rate is present, then the temperature is also scaled by Tnew=rate T. rate<1 cools the system, rate1 heats up the system. rate=1 stabilizes T. Default is rate=1.
Applies a Nose-Hoover thermostat to the molecular dynamics. Inverse_mass is related to how quickly the kinetic energy is forced to the desired temperature. If inverse_mass=0 then the microcanonical (constant energy) MD is recovered. If inverse_mass is much larger than 1 then it will quickly quench the system, but there may be additional timestep errors. For inverse_mass not equal to zero, the usual energy fluctuates. The extended energy (sum of system and heat bath) recorded in MD_ext should be checked for conservation.
delta is side of the cube in which a particle is moved. Adjust delta so the the acceptance ratio is reasonable.
Brownian dynamics. More precisely a Monte Carlo algorithm for simulating the diffusion equation (Smoluchowski equation) applicable when viscous forces are much larger than inertial forces.
tau is the timestep.
Polymer reptation the slithering snake algorithm in continuous space.
* For these keywords order is important: they have an effect after their occurrence.
Variable types: string, integer, real, KEYWORD, required keyword, required choices.
D. Ceperley 8/29/98