1%
2% $Id$
3%
4\label{sec:nwmd}
5%\newcommand{\mc}[3]{\multicolumn{#1}{#2}{#3}}
6\newcommand{\mc}[3]{\mbox{\bf #3}}
7\newcommand{\vb}[1]{\mbox{\verb.#1.}}
8\newcommand{\none}{\multicolumn{2}{|c|}{ }}
9%%%%%%%\renewcommand{\thetable}{\Roman{table}}
10\newcommand{\mcc}[1]{\multicolumn{2}{c}{#1}}
11\def\bmu{\mbox{\boldmath $\mu$}}
12\def\bE{\mbox{\bf E}}
13\def\br{\mbox{\bf r}}
14\def\tT{\tilde{T}}
15\def\t{\tilde{1}}
16\def\ip{i\prime}
17\def\jp{j\prime}
18\def\ipp{i\prime\prime}
19\def\jpp{j\prime\prime}
20\def\etal{{\sl et al.}}
21\def\nwchem{{\bf NWChem}}
22\def\prepare{{\bf prepare}}
23\def\nwtop{{\bf nwtop}}
24\def\nwrst{{\bf nwrst}}
25\def\nwsgm{{\bf nwsgm}}
26
27\section{Introduction}
28
29\subsection{Spacial decomposition}
30The molecular dynamics module of \nwchem\ uses a distribution of data
31based on a spacial decomposition of the molecular system, offering
32an efficient parallel implementation in terms of both memory
33requirements and communication costs, especially for simulations of
34large molecular systems.
35
36Inter-processor communication using the global array tools and the
37design of a data structure allowing distribution based on spacial
38decomposition are the key elements in taking advantage of
39the distribution of memory requirements and computational work with
40minimal communication.
41
42In the spacial decomposition approach, the physical simulation
43volume is divided into rectangular cells, each of which is
44assigned to a processor. Depending on the conditions of the
45calculation and the number of available processors, each processor
46contains one or more of these spacially grouped cells.
47The most important aspects of this decomposition are the dependence
48of the cell sizes and communication cost on the number of processors
49and the shape of the cells, the frequent reassignment of atoms to
50cells leading to a fluctuating number of atoms per cell, and the
51locality of communication which is the main reason for the efficiency
52of this approach for very large molecular systems.
53
54To improve efficiency, molecular systems are broken up into separately
55treated solvent and solute parts. Solvent molecules are assigned to
56the domains according to their center of geometry and are always owned
57by a one node. This avoids solvent--solvent bonded interactions
58crossing node boundaries.  Solute molecules are broken up into
59segments, with each segment assigned to a processor based on its
60center of geometry.  This limits the number of solute bonded
61interactions that cross node boundaries.  The processor to which a
62particular cell is assigned is responsible for the calculation of all
63interactions between atoms within that cell.  For the calculation of
64forces and energies in which atoms in cells assigned to different
65processors are involved, data are exchanged between processors. The
66number of neighboring cells is determined by the size and shape of the
67cells and the range of interaction. The data exchange that takes place
68every simulation time step represents the main communication
69requirements.  Consequently, one of the main efforts is to design
70algorithms and data structures to minimize the cost of this
71communication. However, for very large molecular systems, memory
72requirements also need to be taken into account.
73
74To compromise between these requirements exchange of data is performed
75in successive point to point communications rather than using the
76shift algorithm which reduces the number of communication calls
77for the same amount of communicated data.
78
79For inhomogeneous systems, the computational load of evaluating
80atomic interactions will generally differ between cell pairs.
81This will lead to load imbalance between processors.
82Two algorithms have been implemented that allow for dynamically
83balancing the workload of each processor.
84One method is the dynamic resizing of cells such that cells gradually
85become smaller on the busiest node, thereby reducing the computational
86load of that node. Disadvantages of this method are that the
87efficiency depends on the solute distribution in the simulation volume
88and the redistribution of work depends on the number of nodes which
89could lead to results that depend on the number of nodes used.
90The second method is based on the dynamic redistribution of intra-node
91cell-cell interactions. This method represents a more coarse load
92balancing scheme, but does not have the disadvantages of the cell
93resizing algorithm. For most molecular systems the cell pair
94redistribution is the more efficient and preferred  method.
95
96The description of a molecular system consists of static and dynamic
97information. The static information does not change during a
98simulation and includes items such as connectivity, excluded and third
99neighbor lists, equilibrium values and force constants for all
100bonded and non-bonded interactions. The static information is called
101the topology of the molecular system, and is kept on a separate
102topology file. The dynamic information includes coordinates and
103velocities for all atoms in the molecular system, and is kept in a
104so-called restart file.
105
106\subsection{Topology}
107\label{sec:nwatopology}
108The static information about a molecular system that is needed for
109a molecular simulation is provided to the simulation module in a
110topology file.
111Items in this file include, among many other things,
112a list of atoms, their non-bonded parameters for van der Waals and
113electrostatic interactions, and the complete connectivity in terms
114of bonds, angles and dihedrals.
115
116In molecular systems, a distinction is made between
117{\it solvent} and {\it solute}, which are treated separately.
118A solvent molecule is defined only once in the topology file,
119even though many solvent molecules usually are included in the
120actual molecular system. In the current implementation only one
121solvent can be defined. Everything that is not solvent in the
122molecular system is solute. Each solute atom in the system must
123be explicitly defined in the topology.
124
125Molecules are defined in terms of one or more {\it segment}s.
126Typically, repetitive parts of a molecule are each defined as a single
127segment, such as the amino acid residues in a protein.
128Segments  can be quite complicated to define and are, therefore,
129collected in a set of database files.
130The definition of a molecular system in terms of segments is a
131{\it sequence}.
132
133Topology files are created using the \prepare\ module.
134
135\subsection{Files}
136\label{sec:nwafilenames}
137
138File names used have the form \verb+$system$_$calc$.$ext$+, with
139exception of the topology file (Section \ref{sec:nwatopology}), which is named
140\verb+$system$.top+.
141Anything that refers to the definition of the chemical system can be used
142for \verb+$system$+, as long as no periods or underlines are used.
143The identifier \verb+$calc$+ can be anything that refers to the type of
144calculation to be performed for the system with the topology defined.
145This file naming convention allows for the creation of a single
146topology file \verb+$system$.top+ that can be used for a number of
147different calculations, each identified with a different \verb+$calc$+.
148For example, if {\tt crown.top} is the name of the topology file for
149a crown ether, {\tt crown\_em}, {\tt crown\_md}, {\tt crown\_ti} could
150be used with appropriate extensions for the filenames for energy
151minimization, molecular dynamics simulation and multi-configuration
152thermodynamic integration, respectively. All of these calculations
153would use the same topology file {\tt crown.top}.
154
155\label{sec:nwaextensions}
156
157The extensions \verb+<ext>+ identify the kind of information on a file,
158and are pre-determined.
159\begin{table}[htbp]
160\begin{center}
161\begin{tabular}{ll}
162%{\bf acf} & free energy correlation data file\\
163%{\bf cnv} & free energy convergence file\\
164{\bf dbg} & debug file\\
165%{\bf dre} & distance restraints file\\
166%{\bf fet} & free energy step contribution file\\
167{\bf frg} & fragment file\\
168{\bf gib} & free energy data file\\
169{\bf mri} & free energy multiple run input file\\
170{\bf mro} & free energy multiple run output file\\
171{\bf nw}  & \nwchem\ input file\\
172{\bf nwout}  & \nwchem\ output file\\
173{\bf out} & molecular dynamics output file\\
174{\bf pdb} & PDB formatted coordinate file\\
175{\bf prp} & property file\\
176{\bf qrs} & quenched restart file, resulting from an energy minimization\\
177%{\bf rdf} & radial distribution function output file\\
178%{\bf rdi} & radial distribution function input file\\
179{\bf rst} & restart file, used to start or restart a simulation \\
180{\bf seq} & sequence file, describing the system in segments\\
181{\bf sgm} & segment file, describing segments\\
182{\bf syn} & synchronization time file\\
183{\bf tst} & test file\\
184{\bf tim} & timing analysis file\\
185{\bf top} & topology file, contains the static description of a system\\
186{\bf trj} & trajectory file\\
187\end{tabular}
188\end{center}
189\caption{List of file extensions for nwchem chemical system files.}
190\end{table}
191
192
193\subsection{Databases}
194Database file supplied with NWChem and used by the \prepare\ module are found in
195directories with name \verb+$ffield$_$level$+, \\
196where \verb+$ffield$+ is any of the
197supported force fields (Section \ref{sec:nwaforcefields}).
198The source of the data is identified by \verb+$level$+, and can be
199\begin{center}
200\begin{tabular}{ll}
201\hline
202level   & Description                \\
203{\bf s} & original published data    \\
204{\bf x} & additional published data  \\
205{\bf q} & contributed data           \\
206{\bf u} & user preferred data        \\
207{\bf t} & user defined temporary data\\
208\hline
209\end{tabular}
210\end{center}
211
212The user is can replace these directories or add additional database files
213by specifying them in the $.nwchemrc$ file. or in the \verb+prepare+ input file.
214
215The extension \verb+1-9+ defines the priority of database file.
216
217\begin{table}[htbp]
218\begin{center}
219\begin{tabular}{ll}
220{\bf frg} & fragments\\
221{\bf par} & parameters\\
222{\bf seq} & sequences\\
223{\bf sgm} & segments\\
224\end{tabular}
225\end{center}
226\caption{List of database file extensions.}
227\end{table}
228
229
230The paths of the different database directories should be defined in a file
231{\tt .nwchemrc} in a user's home directory, and provides the user the
232option to select which database files are scanned.
233
234\subsection{Force fields}
235\label{sec:nwaforcefields}
236Force fields recognized are
237\begin{center}
238\begin{tabular}{lll}
239\hline
240Keyword       & Force field   & Status \\
241{\tt amber}   & AMBER99       & AMBER95,GLYCAM also available\\
242{\tt charmm}  & CHARMM        & \\
243%{\tt cvff}   & CVFF          & \\
244%{\tt gromos} & GROMOS87      & \\
245%{\tt oplsa}  & OPLS/AMBER3.0 & \\
246%{\tt oplsg}  & OPLS/GROMOS87 & \\
247\hline
248\end{tabular}
249\end{center}
250The units for lengths, angles, and energies are correspondingly
251nanometers, radians, and kJ/mol.
252
253\section{Format of fragment files}
254Fragment files contain the basic information needed to specify all
255interactions that need to be considered in a molecular simulation.
256The format of the fragment files is described in Table \ref{tbl:nwmdfrg}.
257Normally these files are created by the \prepare\ module. Manual
258editing is needed when, for example, the \prepare\ module could not
259complete atom typing, or when modified charges are required.
260
261\section{Creating segment files}
262\label{sec:nwanwsgm}
263The \prepare\ module is used to generate segment files
264from corresponding fragment files. A segment file contains all
265information for the calculation of bonded and non-bonded interactions
266for a given chemical system using a specific force field.
267
268Which atoms form a fragment is specified in the coordinate file,
269currently only in PDB format.
270the restriction is that bonded interactions may only involve atoms on at
271most two segments does no longer exist as of \nwchem\ release 3.2.1.
272The segment entries define three sets of parameters
273for each interaction.
274
275Free energy perturbations can be performed using set 1 for the
276generation of the ensemble while using sets 2 and/or 3
277as perturbations. Free energy multiconfiguration thermodynamic
278integration and multistep thermodynamic perturbation calculations are
279performed by gradually changing the interactions in the system from
280parameter set 2 to parameter set 3. These modifications can be
281edited into the segment files manually, or introduced directly into
282the topology file using the \verb+modify+ commands in the input for
283the \prepare\ module.
284
285The format of a segment is
286described in Tables \ref{tbl:nwmdseg1}--\ref{tbl:nwmdseg6}.
287
288\section{Creating sequence files}
289A sequence file describes a molecular system in terms of segments. This
290file is generated by the \prepare\ module for the molecular system
291provided on a PDB-formatted coordinate file.
292The file format is given in Table \ref{tbl:nwmdseq}
293
294\section{Creating topology files}
295\label{sec:nwanwtop}
296
297The topology (Section \ref{sec:nwatopology}) describes all static information
298that describes a molecular system. This includes the connectivity in
299terms of bond-stretching, angle-bending and torsional interactions, as well as
300the non-bonded van der Waals and Coulombic interactions.
301
302The topology of a molecular system is generated by the \prepare\ module
303from the sequence in terms of segments as specified on the PDB file.
304For each unique segment specified in this file the
305segment database directories are searched for the segment definition.
306For segments not found in one of the database directories a segment definition
307is generated in the temporary directory if a fragment file was found.
308If a fragment file could not be found, it is generated by the \prepare\ module
309base on what is found on the PDB file.
310
311When all segments are found or created, the parameter substitutions are
312performed, using force field parameters taken from the parameter
313databases. After all lists have been generated the
314topology is written to a local topology file \verb+$system$.top+.
315
316\section{Creating restart files}
317\label{sec:nwanwrst}
318
319Restart files contain all dynamical information about a molecular
320system and are created by the \prepare\ module if a topology file
321is available. The \prepare\ module will automatically generate
322coordinates for hydrogen atoms and monatomic counter ions
323not found on the PDB formatted coordinate file, if no fragment or
324segment files were generated using that PDB file.
325
326The \prepare\ module has a number of other optional input command,
327including solvation.
328
329\section{Molecular simulations}
330The type of molecular dynamics simulation is specified by the
331\nwchem\ task directive.
332\begin{verbatim}
333task md [ energy | optimize | dynamics | thermodynamics ]
334\end{verbatim}
335where the theory keyword {\tt md} specifies use of the molecular
336dynamics module, and the operation keyword is one of
337\begin{description}
338\item
339{\tt energy} for single configuration energy evaluation
340\item
341{\tt optimize} for energy minimization
342\item
343{\tt dynamics} for molecular dynamics simulations and single step
344thermodynamic perturbation free energy molecular dynamics simulations
345\item
346{\tt thermodynamics} for combined multi-configuration thermodynamic
347integration and multiple step thermodynamic perturbation free
348energy molecular dynamics simulations.
349\end{description}
350
351\section{System specification}
352
353The chemical system for a calculation is specified in the topology
354and restart files. These files should be created using the utilities
355\nwtop\ and \nwrst\ before a simulation can be performed.
356The names of these files are determined from the required \verb+system+
357directive.
358
359\begin{verbatim}
360system <string systemid>_<string calcid>
361\end{verbatim}
362
363where the strings \verb+systemid+ and \verb+calcid+ are user defined names
364for the chemical system and the type of calculation to ber performed,
365respectively. These names are used to derive the filenames used for the
366calculation. The topoly file used will be \verb+systemid.top+, while all
367other files are named \verb+systemid_calcid.ext+.
368
369\section{Restarting and continuing simulations}
370\begin{description}
371\item
372\begin{verbatim}
373finish
374\end{verbatim}
375specifies that the current job will finish a previous, incomplete
376simulation, using the input data that have been recorded by that
377previous run in the restart file. Most of the input in the current
378md input block will be ignored.
379
380\item
381\begin{verbatim}
382resume
383\end{verbatim}
384specifies that the current job will be an extension of a previous
385simulation, using most of the  input data that have been recorded by that
386previous run in the restart file. Typically the input in the current
387md input block defines a larger number of steps than the previous job.
388
389\section{Parameter set}
390
391\item
392\begin{verbatim}
393set <integer iset>
394\end{verbatim}
395specifies the use of parameter set \verb+<iset>+ for the
396molecular dynamics simulation.
397The topology file contains three separate parameters sets that can
398be used. The default for \verb+<iset>+ is 1.
399
400\item
401\begin{verbatim}
402lambda <integer ilambda> <integer ilambda>
403\end{verbatim}
404specifies the use of parameter set for the \verb+ilambda+-th
405of \verb+mlambda+ steps.
406
407\item
408\begin{verbatim}
409pset <integer isetp1> [<integer isetp2>]
410\end{verbatim}
411specifies the parameter sets to be used as perturbation potentials
412in single step thermodynamic perturbation free energy evaluations,
413where \verb+<isetp1>+ specifies the first perturbation parameter set and
414\verb+<isetp2>+ specifies the second perturbation parameter set. Legal
415values for \verb+<isetp1>+ are 2 and 3. Legal value for \verb+<isetp2>+ is
4163, in which case \verb+<isetp1>+ can only be 2. If specified, \verb+<iset>+
417is automatically set to 1.
418
419\item
420\begin{verbatim}
421pmf [ equilharm <integer npmfc> | scale <real facpmf>]
422\end{verbatim}
423specifies that any potential of mean force functions defined in the
424topology files are to be used. If \verb+equilharm+ is specified, the
425first \verb+npmfc+ dynamics steps will use a harmonic potential
426in stead of any pmf constraint. If \verb+scale+ is specified, all
427pmf force constants are scaled by a factor \verb+facpmf+.
428
429\item
430\begin{verbatim}
431distar [draver [<integer ndaver default 1>]]
432   [scale <real drsscl>]
433   [after <integer nfdrss>]
434\end{verbatim}
435specifies that any distance restraint functions defined in the
436topology files are to be used.
437
438%Optional keywords are \verb+draver+
439%to specify that average distances over \verb+ndaver+ steps should be
440%used, \verb+scale+ to scale the specified force constants, and
441%\verb+after+
442
443\item
444\begin{verbatim}
445qhop [<integer nfhop default 10>]
446   [<real rhop default 0.35>]
447   [<real thop default 0.02>]
448\end{verbatim}
449specifies that a Q-HOP simulation is to be carried out with attempted
450proton hops every $nfhop$ steps, a cutoff for the donor-acceptor pair
451distance of $rhop$ nm, and a minimum time before back hopping can occur
452of $thop$ ps.
453
454
455\end{description}
456
457\section{Energy minimization algorithms}
458The energy minimization of the system as found in the restart file
459is performed with the following directives. If both are specified,
460steepest descent energy minimization precedes conjugate gradient
461minimization.
462
463%\begin{description}
464\begin{description}
465\item
466\begin{verbatim}
467sd <integer msdit> [init <real dx0sd>] [min <real dxsdmx>] \
468                   [max <real dxmsd>]
469\end{verbatim}
470specifies the variables for steepest descent energy minimizations,
471where \verb+<msdit>+ is the maximum number of steepest descent steps taken,
472for which the default is 100, \verb+<dx0sd>+ is the initial step size in nm
473for which the default is 0.001, \verb+<dxsdmx>+ is the threshold for the
474step size in nm for which the default is 0.0001, and \verb+<dxmsd>+ is the
475maximum allowed step size in nm for which the default is 0.05.
476\item
477\begin{verbatim}
478cg <integer mcgit> [init <real dx0cg>] [min <real dxcgmx>] \
479                   [cy <integer ncgcy>]
480\end{verbatim}
481specifies the variables for conjugate gradient energy minimizations,
482where \verb+<mcgit>+ is the maximum number of conjugate gradient steps
483taken, for which the default is 100, \verb+<dx0cg>+ is the initial search
484interval size in nm for which the default is 0.001, \verb+<dxcgmx>+ is the
485threshold for the step size in nm for which the default is 0.0001, and
486\verb+<ncgcy>+ is the number of conjugate gradient steps after which the
487gradient history is discarded for which the default is 10. If conjugate
488gradient energy minimization is preceded by steepest descent energy
489minimization, the search interval is set to twice the final step of the
490steepest descent energy minimization.
491\end{description}
492%\end{description}
493
494\section{Multi-configuration thermodynamic integration}
495The following keywords control free energy difference simulations.
496Multi-configuration thermodynamic integrations are always combined
497with multiple step thermodynamic perturbations.
498\begin{description}
499
500\item
501\begin{verbatim}
502(forward | reverse) [[<integer mrun> of] <integer maxlam>]
503\end{verbatim}
504specifies the direction and number of integration steps in free
505energy evaluations, with {\tt forward} being the default direction.
506\verb+<mrun>+ is the number of ensembles that will be generated in
507this calculation, and \verb+<maxlam>+ is the total number of ensembles
508to complete the thermodynamic integration. The default value for
509\verb+<maxlam>+ is 21. The default value of \verb+<mrun>+ is the
510value of \verb+<maxlam>+.
511
512\item
513\begin{verbatim}
514error <real edacq>
515\end{verbatim}
516specifies the maximum allowed statistical error in each generated
517ensemble, where \verb+<edacq>+ is the maximum error allowed in the
518ensemble average derivative of the Hamiltonian with respect to
519$\lambda$ with a default of 5.0 kJ~mol$^{-1}$.
520
521\item
522\begin{verbatim}
523drift <real ddacq>
524\end{verbatim}
525specifies the maximum allowed drift in the free energy result,
526where \verb+<ddacq>+ is the maximum drift allowed in the
527ensemble average derivative of the Hamiltonian with respect to
528$\lambda$with a default of 5.0 kJ~mol$^{-1}$ps$^{-1}$.
529
530\item
531\begin{verbatim}
532factor <real fdacq>
533\end{verbatim}
534specifies the maximum allowed change in ensemble size
535where \verb+<fdacq>+ is the minimum size of an ensemble relative to the
536previous ensemble in the calculation with a default value of 0.75.
537
538\item
539\begin{verbatim}
540decomp
541\end{verbatim}
542specifies that a free energy decomposition is to be carried out.
543Since free energy contributions are path dependent, results from a
544decomposition analysis can no be unambiguously interpreted, and
545the default is not to perform this decomposition.
546
547\item
548\begin{verbatim}
549sss [delta <real delta>]
550\end{verbatim}
551specifies that atomic non-bonded interactions describe a dummy atom
552in either the initial or final state of the thermodynamic calculation
553will be calculated using separation-shifted scaling, where \verb+<delta>+
554is the separation-shifted scaling factor with a default of 0.075 nm$^2$.
555This scaling method prevents problems associated with singularities in
556the interaction potentials.
557
558\item
559\begin{verbatim}
560new | renew | extend
561\end{verbatim}
562specifies the initial conditions for thermodynamic calculations.
563{\tt new} indicates that this is an initial mcti calculation, which
564is the default. {\tt renew} instructs to obtain the initial
565conditions for each $\lambda$ from the {\bf mro}-file from a previous
566mcti calculation, which has to be renamed to an {\bf mri}-file. The
567keyword {\tt extend} will extend a previous mcti calculation from the
568data read from an {\bf mri}-file.
569\end{description}
570
571\section{Time and integration algorithm directives}
572Following directives control the integration of the equations of motion.
573\begin{description}
574
575\item
576\begin{verbatim}
577leapfrog | leapfrog_bc
578\end{verbatim}
579specifies the integration algorithm,
580where {\tt leapfrog} specifies the default leap frog integration, and
581{\tt leapfrog\_bc} specifies the Brown-Clarke leap frog integrator.
582
583\item
584\begin{verbatim}
585guided [<real fguide default 0.2> [<real tguide default 0.2>]]
586\end{verbatim}
587specifies the use of the guided molecular dynamics simulation
588technique. Variable $fguide$ defines the fraction of the averaged
589forces $g$ to be added to the forces $f^{f}$ evaluated using the force
590field functions to obtain the forces $f$ used to advance the coordinates.
591\begin{equation}
592f_i=f^{f}_i+fguide * g_{i-1}
593\end{equation}
594Variable $tguide$ defines the length of the averaging relative to the
595timestep $\Delta t$.
596\begin{equation}
597g_i = {\Delta t\over tguide} f_i + \left(1- {\Delta t\over tguide}\right)
598g_{i-1}
599\end{equation}
600The current implementation is still under
601development.
602
603\item
604\begin{verbatim}
605equil <integer mequi>
606\end{verbatim}
607specifies the number of equilibration steps \verb+<mequi>+, with a default
608of 100.
609
610\item
611\begin{verbatim}
612data <integer mdacq> [over <integer ldacq>]]
613\end{verbatim}
614specifies the number of data gathering steps \verb+<mdacq>+ with a
615default of 500. In multi-configuration thermodynamic integrations
616\verb+<mequi>+ and \verb+<mdacq>+ are for each of the ensembles, and
617variable \verb+<ldacq>+ specifies the minimum number of data gathering steps
618in each ensemble. In regular molecular dynamics simulations \verb+<ldacq>+
619is not used. The default value for \verb+<ldacq>+ is the value of \verb+<mdacq>+.
620
621\item
622\begin{verbatim}
623time <real stime>
624\end{verbatim}
625specifies the initial time \verb+<stime>+ of a molecular simulation in ps,
626with a default of 0.0.
627
628\item
629\begin{verbatim}
630step <real tstep>
631\end{verbatim}
632specifies the time step \verb+<tstep>+ in ps, with 0.001 as the default value.
633\end{description}
634
635\section{Ensemble selection}
636Following directives control the ensemble type.
637
638\begin{description}
639
640\item
641\begin{verbatim}
642isotherm [<real tmpext> [<real tmpext2>]] [trelax <real tmprlx> [<real tmsrlx>]] \
643         [anneal [<real tann1>] <real tann2>]
644\end{verbatim}
645specifies a constant temperature ensemble using Berendsen's thermostat,
646where \verb+<tmpext>+ is the external temperature with a default of 298.15~K,
647and \verb+<tmprlx>+ and \verb+<tmsrlx>+ are temperature relaxation times in ps
648with a default of 0.1. If only \verb+<tmprlx>+ is given the complete system
649is coupled to the heat bath with relaxation time \verb+<tmprlx>+. If both
650relaxation times are supplied, solvent and solute are independently coupled
651to the heat bath with relaxation times \verb+<tmprlx>+ and \verb+<tmsrlx>+,
652respectively. If keyword \verb+anneal+ is specified, the external temperature
653will change from \verb+tmpext+ to \verb+tempext2+ between simulation time
654\verb+tann1+ and \verb+tann2+
655
656\item
657\begin{verbatim}
658isobar [<real prsext>] [trelax <real prsrlx> ] \
659       [compress <real compr>] [anisotropic] [xy | z | xy-z]
660\end{verbatim}
661specifies a constant pressure ensemble using Berendsen's piston,
662where \verb+<prsext>+ is the external pressure with a default of 1.025~10$^5$ Pa,
663\verb+<prsrlx>+ is the pressure relaxation time in ps with a default of 0.5, and
664\verb+<compr>+ is the system compressibility in m$^2$N$^{-1}$ with a
665default of 4.53E-10. Optional keywords \verb+xy+, \verb+z+ and \verb+xy-z+
666may be used to specify that pressure scaling is to be applied in
667the \verb+x+ and \verb+y+ dimension only, the \verb+z+ dimension only, or,
668in all three dimensions with identical scaling in the \verb+x+ and \verb+y+
669dimension. The last option requires that \verb+anisotropic+ is also specified.
670
671\end{description}
672
673\section{Velocity reassignments}
674Velocities can be periodically reassigned to reflect a certain temperature.
675\begin{description}
676\item
677\begin{verbatim}
678vreass <integer nfgaus> <real tgauss>
679       [fraction [<real frgaus default 0.5]]
680       [once]
681       [(first | initial)] [(last | final)]
682\end{verbatim}
683specifies that velocities will be reassigned every \verb+<nfgaus>+ molecular
684dynamics steps, reflecting a temperature of \verb+<tgauss>+~K. The default
685is not to reassign velocities, i.e.\ \verb+<nfgaus>+ is 0. Keyword
686\verb+fraction+ allows the specification of the fraction of the new
687velocities are random. Keyword \verb+once+ specifies that velocity
688reassignment only should be done in the first step. Keywords \verb+first+
689or \verb+initial+ and \verb+last+ or \verb+final+ specify that
690velocity reassigment should only be applied in the first and last
691window of multiple run simulations.
692
693\end{description}
694
695\section{Cutoff radii}
696Cutoff radii can be specified for short range and long range interactions.
697\begin{description}
698\item
699\begin{verbatim}
700cutoff [short] <real rshort> [long <real rlong>] \
701       [qmmm <real rqmmm>]
702\end{verbatim}
703specifies the short range cutoff radius \verb+<rshort>+, and the long range
704cutoff radius \verb+<rlong>+ in nm. If the long range cutoff radius
705is larger than the short range cutoff radius the twin range method will
706be used, in which short range forces and energies are evaluated every
707molecular dynamics step, and long range forces and energies with a
708frequency of \verb+<nflong>+ molecular dynamics steps. Keyword
709\verb+qmmm+ specifies the radius of the zone around quantum atoms
710defining the QM/MM bare charges.
711The default value for \verb+<rshort>+, \verb+<rlong>+ and \verb+<rqmmm>+
712is 0.9~nm.
713\end{description}
714
715\section{Polarization}
716First order and self consistent electronic polarization models have
717been implemented.
718\begin{description}
719\item
720\begin{verbatim}
721polar (first | scf [[<integer mpolit>] <real ptol>])
722\end{verbatim}
723specifies the use of polarization potentials,
724where the keyword {\tt first} specifies the first order polarization
725model, and {\tt scf} specifies the self consistent polarization field
726model, iteratively determined with a maximum of \verb+<mpolit>+
727iterations to within a tolerance of \verb+<ptol>+ D in the generated
728induced dipoles. The default is not to use polarization models.
729\end{description}
730
731\section{External electrostatic field}
732\begin{description}
733
734\item
735\begin{verbatim}
736field <real xfield> [freq <real xffreq>] [vector <real xfvect(1:3)>]
737\end{verbatim}
738specifies an external electrostatic field,
739where \verb+<xfield>+ is the field strength, \verb+<xffreq>+ is the
740frequency in MHz and \verb+<xfvect>+ is the external field vector.
741\end{description}
742
743\section{Constraints}
744Constraints are satisfied using the SHAKE
745coordinate resetting procedure.
746\begin{description}
747
748\item
749\begin{verbatim}
750shake [<integer mshitw> [<integer mshits>]]  \
751      [<real tlwsha> [<real tlssha>]]
752\end{verbatim}
753specifies the use of SHAKE constraints,
754where \verb+<mshitw>+ is the maximum number of solvent SHAKE iterations,
755and \verb+<mshits>+ is the maximum number of solute SHAKE iterations. If
756only \verb+<mshitw>+ is specified, the value will also be used for \verb+<mshits>+.
757The default maximum number of iterations is 100 for both.
758\verb+<tlwsha>+ is the solvent SHAKE tolerance in nm, and \verb+<tlssha>+ is
759the solute SHAKE tolerance in nm. If only \verb+<tlwsha>+ is specified, the
760value given will also be used for \verb+<tlssha>+. The default tolerance
761is 0.001~nm for both.
762
763\item
764\begin{verbatim}
765noshake (solvent | solute)
766\end{verbatim}
767disables SHAKE and treats the bonded interaction according to the force
768field.
769
770\end{description}
771
772\section{Long range interaction corrections}
773Long range electrostatic interactions are implemented using the
774smooth particle mesh Ewald technique, for neutral periodic cubic systems in
775the constant volume ensemble, using pair interaction potentials. Particle-mesh
776Ewald long range interactions can only be used in molecular dynamics simulations
777using effective pair potentials, and not in free energy simulations, QMD or
778QM/MM simulations.
779
780\begin{description}
781\item
782\begin{verbatim}
783pme [grid <integer ng>] [alpha <real ealpha>] \
784    [order <integer morder>] [fft <integer imfft>]\
785    [procs <integer nprocs>] [solvent]
786\end{verbatim}
787specifies the use of smooth particle-mesh Ewald long range
788interaction treatment,
789where \verb+ng+ is the number of grid points per dimension,
790\verb+ealpha+ is the Ewald coefficient in nm$^{-1}$, with a default
791that leads to a tolerance of $10^{-4}$ at the short range cutoff radius,
792and \verb+morder+ is order of the Cardinal B-spline
793interpolation which must be an even number and at least 4 (default
794value). A platform specific 3D fast Fourier transform is used, if
795available, when \verb+imfft+ is set to 2. \verb+nprocs+ can be used to
796define a subset of processors to be used to do the FFT calculations.
797If \verb+solvent+ is specified, the charge grid will be calculated from
798the solvent charges only.
799
800\item
801\begin{verbatim}
802react [<real dielec default 80.0>]
803\end{verbatim}
804specifies that a simple reaction field correction is used with a
805dielectric constant \verb+dielec+. This is an
806experimental option that has not been well tested.
807
808\end{description}
809
810\section{Fixing coordinates}
811Solvent or solute may be fixed using the following keywords.
812\begin{description}
813\item
814\begin{verbatim}
815( fix | free )
816 solvent ( [<integer idfirst> [<integer idlast>]] |
817            ( within | beyond) <real rfix> <string atomname> ) | \
818 solute  ( [<integer idfirst> [<integer idlast>]] [ heavy | {<string atomname>}] |
819            ( within | beyond) <real rfix> <string atomname> )
820 [permanent]
821\end{verbatim}
822For solvent the molecule numbers \verb+idfirst+ and \verb+idlast+may be
823specified to be the first and last molecule to which the directive
824applies. If omitted, the directive applies to all molecules. For solute,
825the segment numbers \verb+idfirst+ and \verb+idlast+may be
826specified to be the first and last segment to which the directive
827applies. If omitted, the directive applies to all segments. In addition,
828the keyword \verb+heavy+ may be specified to apply to all non hydrogen
829atoms in the solute, or a set of atom names may be specified in which
830a wildcard character \verb+?+ may be used. Keyword \verb+permanent+
831is used to keep the specification on the restart file for subsequent
832simulations.
833\end{description}
834
835\section{Special options}
836\begin{description}
837
838\item
839\begin{verbatim}
840import [<integer impfr default 1> [<integer impto default impfr> \
841       [<integer nftri default 1>]]]
842\end{verbatim}
843specifies the import of frames \verb+impfr+ to \verb+impto+ with
844frequency \verb+nftri+ from a trajectory file with extension
845\verb+tri+ for which energies and forces are to be recalculated.
846This option only applied to \verb+task md energy+.
847
848\item
849\begin{verbatim}
850detail
851\end{verbatim}
852specifies that moments of inertia and radii of gyration will be part of the
853recorded properties.
854
855\item
856\begin{verbatim}
857profile
858\end{verbatim}
859specifies that execution time profiling data will be part of the
860recorded properties.
861
862\item
863\begin{verbatim}
864scale <real scaleq>
865\end{verbatim}
866specifies that all charges will be scaled by the factro \verb+scaleq+.
867
868\item
869\begin{verbatim}
870collapse [<real fcoll default 10.0> [ segment | z | xy ]
871\end{verbatim}
872specifies that additional forces directed to the origin of the
873simulation cell with strength \verb+fcoll+ will be
874applied to all solute molecules. If \verb+z+ or \verb+xy+ is
875specified, these forces will only apply in the specified dimension(s).
876
877
878\item
879\begin{verbatim}
880include fixed
881\end{verbatim}
882specifies that energies will be evaluated between fixed atoms.
883Normally these interactions are excluded from the pairlists.
884
885
886\item
887\begin{verbatim}
888eqm <real eqm>
889\end{verbatim}
890specifies the zero point of energy in QMD simulations.
891
892\item
893\begin{verbatim}
894atomlist
895\end{verbatim}
896specifies that pairlists will be atom based. Normally pairlist
897are charge group based.
898
899\end{description}
900
901\section{Autocorrelation function}
902For the evaluation of the statistical error of multi-configuration
903thermodynamic integration free energy results a correlated data
904analysis is carried out, involving the calculation of the
905autocorrelation function of the derivative of the Hamiltonian with
906respect to the control variable $\lambda$.
907\begin{description}
908
909\item
910\begin{verbatim}
911auto <integer lacf> [fit <integer nfit>] [weight <real weight>]
912\end{verbatim}
913controls the calculation of the autocorrelation,
914where \verb+<lacf>+ is the length of the autocorrelation function, with
915a default of 1000, \verb+<nfit>+ is the number of functions used in the
916fit of the autocorrelation function, with a default of 15, and
917\verb+<weight>+ is the weight factor for the autocorrelation function,
918with a default value of 0.0.
919\end{description}
920
921\section{Print options}
922Keywords that control print to the output file, with extension {\bf out}.
923Print directives may be combined to a single directive.
924
925\begin{description}
926
927\item
928\begin{verbatim}
929print [topol [nonbond] [solvent] [solute]]     \
930      [step <integer nfoutp> [extra] [energy]] \
931      [stat <integer nfstat>]                  \
932      [energies [<integer nfener>]]            \
933      [forces [<integer nfforce>]]             \
934      [matrix]                                 \
935      [expect <integer npxpct>]                \
936      [timing]                                 \
937      [pmf [<integer iprpmf>]]                 \
938      [out6]                                   \
939      [dayout]
940\end{verbatim}
941
942Keyword \verb+topol+ specifies printing the topology information,
943where {\tt nonbond} refers to the non-bonded interaction parameters,
944{\tt solvent} to the solvent bonded parameters, and {\tt solute} to the
945solute bonded parameters. If only {\tt topol} is specified, all
946topology information will be printed to the output file.
947
948Keyword \verb+step+
949specifies the frequency \verb+nfoutp+ of printing molecular dynamics step
950information to the output file. If the keyword {\tt extra} is specified
951additional energetic data are printed for solvent and solute separately.
952If the keyword {\tt energy} is specified, information is printed for
953all bonded solute interactions.
954The default for \verb+nfoutp+ is 0. For molecular dynamics simulations
955this frequency is in time steps, and for multi-configuration thermodynamic
956integration in $\lambda$-steps.
957
958Keyword \verb+stat+
959specifies the frequency \verb+<nfstat>+ of printing statistical information
960of properties that are calculated during the simulation.
961For molecular dynamics simulation
962this frequency is in time steps, for multi-configuration thermodynamic
963integration in $\lambda$-steps.
964
965Keyword \verb+energies+
966specifies the frequency \verb+nfener+ of printing solute bonded energies
967the output file for energy/import calculations.
968The default for \verb+nfener+ is 0.
969
970Keyword \verb+forces+
971specifies the frequency \verb+nfforc+ of printing solute forces
972the output file for energy/import calculations.
973The default for \verb+nfforc+ is 0.
974
975Keyword \verb+matrix+ specifies that a solute distance matrix is to
976be printed.
977
978Keyword \verb+expect+ is obsolete.
979
980Keyword \verb+timing+ specifies that timing data is printed.
981
982Keyword \verb+pmf+ specifies that pmf data is printed every \verb+iprpmf+
983steps.
984Keyword \verb+out6+ specifies that output is written to standard out in stead
985of the output file with extension \verb+out+.
986
987Keyword \verb+dayout+ is obsolete.
988
989
990\end{description}
991
992\section{Periodic updates}
993Following keywords control periodic events during a molecular
994dynamics or thermodynamic integration simulation.
995Update directives may be combined to a single directive.
996
997\begin{description}
998\item
999\begin{verbatim}
1000update [pairs <integer nfpair default 1>]                  \
1001       [long <integer nflong default 1>]                   \
1002       [center <integer nfcntr default 0> [zonly | xyonly] \
1003               [fraction <integer idscb(1:5)>]             \
1004       [motion <integer nfslow default 0>]                 \
1005       [analysis <integer nfanal default 0>]               \
1006       [rdf <integer nfrdf default 0>                      \
1007               [range <real rrdf>] [bins <integer ngl>]    \
1008\end{verbatim}
1009
1010Keyword \verb+pairs+
1011specifies the frequency \verb+<nfpair>+ in molecular dynamics steps of
1012updating the pair lists. The default for the frequency is 1.
1013In addition, pair lists are also updated after each step in which
1014recording of the restart or trajectory files is performed. Updating
1015the pair lists includes the redistribution of atoms that changed
1016domain and load balancing, if specified.
1017
1018Keyword \verb+long+
1019specifies the frequency \verb+<nflong>+ in molecular dynamics steps
1020of updating the long range forces. The default frequency is 1.
1021The distinction of short range and long range forces is only
1022made if the long range cutoff radius was specified to be larger
1023than the short range cutoff radius. Updating the long range forces
1024is also done in every molecular dynamics step in which the
1025pair lists are regenerated.
1026
1027Keywrod \verb+center+
1028specifies the frequency \verb+<nfcntr>+ in molecular dynamics steps in
1029which the center of geometry of the solute(s) is translated to the
1030center of the simulation volume. Optional keyword \verb+zonly+ or
1031\verb+xyonly+ can be used to specify that centering will take place in
1032the z-direction or in the xy-plane only.
1033The solute fractions determining the
1034solutes that will be centered are specified by the keyword
1035{\tt fraction} and the vector \verb+<idscb>+, with a maximum of 5 entries.
1036This translation is implemented such that it has no effect on any
1037aspect of the simulation. The default is not to center, i.e. nfcntr is
10380. The default fraction used to center solute is 1.
1039
1040Keyword \verb+motion+
1041specifies the frequency \verb+<nfslow>+ in molecular dynamics steps of
1042removing the overall rotational and center of mass translational motion.
1043
1044Keyword \verb+analysis+
1045specifies the frequency \verb+<nfanal>+ in molecular dynamics steps of
1046invoking the analysis module. This option is obsolete.
1047
1048Keyword \verb+rdf+
1049specifies the frequency \verb+<nfrdf>+ in molecular dynamics steps of
1050calculating contributions to the radial distribution functions.
1051The default is 0. The range of the radial distribution
1052functions is given by \verb+<rrdf>+ in nm, with a default of the short
1053range cutoff radius. Note that radial distribution functions are not
1054evaluated beyond the short range cutoff radius. The number of
1055bins in each radial distribution function is given by \verb+<ngl>+, with
1056a default of 1000. This option is no longer supported.
1057If radial distribution function are to be
1058calculated, a {\bf rdi} files needs to be available in which the
1059contributions are specified as follows.
1060\begin{center}
1061\begin{tabular}{lll}
1062\hline\hline
1063Card & Format & Description \\ \hline
1064I-1  & i & Type, 1=solvent-solvent, 2=solvent-solute,
10653-solute-solute\\
1066I-2  & i & Number of the rdf for this contribution\\
1067I-3  & i & First atom number \\
1068I-4  & i & Second atom number \\
1069\hline
1070\end{tabular}
1071\end{center}
1072\end{description}
1073
1074\section{Recording}
1075The following keywords control recording data to file.
1076Record directives may be combined to a single directive.
1077
1078\begin{description}
1079
1080\item
1081\begin{verbatim}
1082record [rest <integer nfrest> [keep]]     \
1083       [coord <integer nfcoor default 0>] \
1084       [wcoor <integer nfwcoo default 0>] \
1085       [scoor <integer nfscoo default 0>] \
1086       [veloc <integer nfvelo default 0>] \
1087       [wvelo <integer nfwvel default 0>] \
1088       [svelo <integer nfsvel default 0>] \
1089       [force <integer nfvelo default 0>] \
1090       [wforc <integer nfwvel default 0>] \
1091       [sforc <integer nfsvel default 0>] \
1092       [(prop | prop_average) <integer nfprop default 0>]  \
1093       [free <integer nffree default 1>]  \
1094       [sync <integer nfsync default 0>]  \
1095       [times <integer nftime default 0>] \
1096       [acf] [cnv] [fet]
1097       [binary] [ascii] [ecce] [argos]
1098\end{verbatim}
1099
1100Keyword \verb+rest+
1101specifies the frequency \verb+<nfrest>+ in molecular dynamics steps
1102of rewriting the restart file, with extension \verb+rst+.
1103For multi-configuration
1104thermodynamic integration simulations the frequency is in
1105steps in $\lambda$. The default is not to record. The restart
1106file is used to start or restart simulations. The keyword {\tt keep}
1107causes all restart files written to be kept on disk, rather than
1108to be overwritten.
1109
1110Keyword \verb+coord+
1111specifies the frequency \verb+<nfcoor>+ in molecular dynamics steps
1112of writing coordinates to the trajectory file. This directive redefines
1113previous \verb+coord+, \verb+wcoor+ and \verb+scoor+
1114directives. The default is not to record.
1115
1116Keyword \verb+wcoor+
1117specifies the frequency \verb+<nfcoor>+ in molecular dynamics steps
1118of writing solvent coordinates to the trajectory file. This keyword
1119takes precedent over \verb+coord+. This directive redefines
1120previous \verb+coord+, \verb+wcoor+ and \verb+scoor+
1121directives.
1122The default is not to record.
1123
1124Keyword \verb+scoor+
1125specifies the frequency \verb+<nfscoo>+ in molecular dynamics steps
1126of writing solute coordinates to the trajectory file. This keyword
1127takes precedent over \verb+coord+. This directive redefines
1128previous \verb+coord+, \verb+wcoor+ and \verb+scoor+
1129directives.
1130The default is not to record.
1131
1132Keyword \verb+veloc+
1133specifies the frequency \verb+<nfvelo>+ in molecular dynamics steps
1134of writing velocities to the trajectory file. This directive redefines
1135previous \verb+veloc+, \verb+wvelo+ and \verb+svelo+
1136directives.
1137The default is not to record.
1138
1139Keyword \verb+wvelo+
1140specifies the frequency \verb+<nfvelo>+ in molecular dynamics steps
1141of writing solvent velocitiesto the trajectory file. This keyword
1142takes precedent over \verb+veloc+. This directive redefines
1143previous \verb+veloc+, \verb+wvelo+ and \verb+svelo+
1144directives.
1145The default is not to record.
1146
1147Keyword \verb+svelo+
1148specifies the frequency \verb+<nfsvel>+ in molecular dynamics steps
1149of writing solute velocities to the trajectory file. This keyword
1150takes precedent over \verb+veloc+. This directive redefines
1151previous \verb+veloc+, \verb+wvelo+ and \verb+svelo+
1152directives.
1153The default is not to record.
1154
1155Keyword \verb+force+
1156specifies the frequency \verb+<nfvelo>+ in molecular dynamics steps
1157of writing forces to the trajectory file. This directive redefines
1158previous \verb+vforce+, \verb+wforc+ and \verb+sforc+
1159directives.
1160The default is not to record.
1161
1162Keyword \verb+wforc+
1163specifies the frequency \verb+<nfvelo>+ in molecular dynamics steps
1164of writing solvent forcesto the trajectory file. This keyword
1165takes precedent over \verb+force+. This directive redefines
1166previous \verb+vforce+, \verb+wforc+ and \verb+sforc+
1167directives.
1168The default is not to record.
1169
1170Keyword \verb+sforc+
1171specifies the frequency \verb+<nfsvel>+ in molecular dynamics steps
1172of writing solute forces to the trajectory file. This keyword
1173takes precedent over \verb+force+. This directive redefines
1174previous \verb+vforce+, \verb+wforc+ and \verb+sforc+
1175directives.
1176The default is not to record.
1177
1178Keyword \verb+prop+
1179specifies the frequency \verb+<nfprop>+ in molecular dynamics steps
1180of writing information to the property file, with extension
1181\verb+prp+. The default is not to record.
1182
1183Keyword \verb+prop_average+
1184specifies the frequency \verb+<nfprop>+ in molecular dynamics steps
1185of writing average information to the property file, with extension
1186\verb+prp+.
1187The default is not to record.
1188
1189Keyword \verb+free+
1190specifies the frequency \verb+<nffree>+ in multi-configuration
1191thermodynamic integration steps to record data to the
1192free energy data file, with extension \verb+gib+.
1193The default is 1, i.e.\ to record at every $\lambda$.
1194This option is obsolete. All data are required to do the
1195final analysis.
1196
1197%\item
1198%\begin{verbatim}
1199%record cnv
1200%\end{verbatim}
1201%specifies that free energy convergence data will be written to the
1202%free energy convergence file, with extension \verb+cnv+.
1203%
1204%\item
1205%\begin{verbatim}
1206%record acf
1207%\end{verbatim}
1208%specifies that free energy derivative autocorrelation data will be
1209%written to the free energy autocorrelation file, with extension
1210%\verb+acf+.
1211%
1212%\item
1213%\begin{verbatim}
1214%record fet
1215%\end{verbatim}
1216%that free energy vs.\ time data will be recorded to the free energy
1217%data file, with extension \verb+fet+.
1218%
1219
1220Keyword \verb+sync+
1221specifies the frequency  \verb+<nfsync>+ in molecular dynamics steps
1222of writing information to the synchronization file, with extension
1223\verb+syn+.
1224The default is not to record.
1225The information written is the simulation time, the wall clock time
1226of the previous MD step, the wall clock time of the previous force
1227evaluation, the total synchronization time, the largest
1228synchronization time and the node on which the largest synchronization
1229time was found. The recording of synchronization times is part of the
1230load balancing algorithm. Since load balancing is only performed when
1231pair-lists are updated, the frequency \verb+<nfsync>+ is correlated
1232with the frequency of pair-list updates \verb+<nfpair>+. This directive
1233is only needed for analysis of the load balancing performance. For
1234normal use this directive is not used.
1235
1236Keyword \verb+times+
1237specifies the frequency  \verb+<nfsync>+ in molecular dynamics steps
1238of writing information to the timings file, with extension
1239\verb+tim+.
1240The default is not to record.
1241The information written is wall clock time used by each of the
1242processors for the different components in the force evaluation.
1243This directive is only needed for analysis of the wall clock time
1244distribution. For normal use this directive is not used.
1245
1246Keywords \verb+acf+, \verb+cnv+ and \verb+fet+ are obsolete.
1247
1248Keywords \verb+binary+, \verb+ascii+, \verb+ecce+ and \verb+argos+ are obsolete.
1249
1250\end{description}
1251
1252\section{Program control options}
1253\begin{description}
1254\item
1255\begin{verbatim}
1256load [reset]
1257     ( none |
1258       size [<real factld>] |
1259       sizez [<real factld>] | pairs |
1260      (pairs [<integer ldpair>] size [<real factld>]) )
1261     [last]
1262     [minimum]
1263     [average]
1264     [combination]
1265     [iotime]
1266     [experimental]
1267
1268\end{verbatim}
1269determines the type of dynamic load balancing performed,
1270where the default is {\tt none}. Load balancing option {\tt size}
1271is resizing cells on a node, and {\tt pairs} redistributes the
1272cell-cell interactions over nodes. Keyword \verb+reset+ will reset the
1273load balancing read from the restart file. The level of cell resizing
1274can be influenced with $factld$. The cells on the busiest node are
1275resized with a factor
1276\begin{equation}
1277\left( 1 - factld * { {T_{sync} \over n_p} - t^{min}_{sync} \over t_{wall}}
1278\right)^{1\over 3}
1279\end{equation}
1280where $T_{sync}$ is the accumulated synchronization time of all nodes,
1281$n_p$ is the total number of nodes, $t^{min}_{sync}$ is the synchronization
1282time of the busiest node, and $t_{wall}$ is the wall clock time of the
1283molecular dynamics step.
1284For the combined load balancing, \verb+ldpair+ is the number of successive pair
1285redistribution load balancing steps in which the accumulated synchronization
1286time increases, before a resizing load balancing step will be attempted.
1287Load balancing is only performed in molecular dynamics steps in which the
1288pair-list is updated. The default load balancing is equivalent to specifying\\
1289\begin{verbatim}
1290load pairs 10 size 0.75
1291\end{verbatim}
1292Keyword \verb+last+ specifies that the load balancing is based on the
1293synchronization times of the last step. This is the default.
1294Keyword \verb+average+ specifies that the load balancing is based on the
1295average synchronization times since the last load balancing step.
1296Keyword \verb+minimum+ specifies that the load balancing is based on the
1297minimum synchronization times since the last load balancing step.
1298Keywords \verb+combination+, \verb+iotime+ and \verb+experimental+ are
1299experimental load balancing options that should not be used in
1300production runs.
1301
1302\item
1303\begin{verbatim}
1304(pack | nopack)
1305\end{verbatim}
1306specifies if data are communicated in packed or unpacked form. The
1307default is \verb+pack+.
1308
1309\item
1310\begin{verbatim}
1311procs <integer npx> <integer npy> <integer npz>
1312\end{verbatim}
1313specifies the distribution of the available processors over the three
1314Cartesian dimensions. The default distribution is chosen such that,
1315\verb+<npx>+$*$\verb+<npy>+$*$\verb+<npz>+=\verb+<np>+
1316and \verb+<npx>+ $<=$ \verb+<npy>+ $<=$ \verb+<npz>+,
1317where \verb+<npx>+, \verb+<npy>+ and \verb+<npz>+ are the processors in the
1318\verb+x+, \verb+y+ and \verb+z+ dimension respectively, and \verb+<np>+ is the number of processors
1319allocated for the calculation. Where more than one combination
1320of \verb+<npx>+, \verb+<npy>+ and \verb+<npz>+ are possible, the
1321combination is chosen with the minimum value of
1322\verb+<npx>+$+$\verb+<npy>+$+$\verb+<npz>+. To change the default setting
1323the following optional input option is provided.
1324
1325\item
1326\begin{verbatim}
1327cells <integer nbx> <integer nby> <integer nbz>
1328\end{verbatim}
1329specifies the distribution of cells,
1330where \verb+<nbx>+, \verb+<nby>+ and \verb+<nbz>+ are the number of
1331cells in \verb+x+, \verb+y+ and \verb+z+ direction, respectively.
1332The molecular system is decomposed into cells that form the smallest
1333unit for communication of atomic data between nodes. The size of the
1334cells is per default set to the short-range cutoff radius. If
1335long-range cutoff radii  are used the cell size is set to half the
1336long-range cutoff radius if it is larger than the short-range cutoff.
1337If the number of cells in a dimension is less than the number of
1338processors in that dimension, the number of cells is set to the number
1339of processors.
1340
1341\item
1342\begin{verbatim}
1343extra <integer madbox>
1344\end{verbatim}
1345sets the number of additional cells for which memory is allocated.
1346In rare events the amount of memory set aside per node is insufficient
1347to hold all atomic coordinates assigned to that node. This leads to
1348execution which aborts with the message that {\tt mwm} or {\tt msa} is too
1349small. Jobs may be restarted with additional space allocated by
1350where \verb+<madbox>+ is the number of additional cells that are allocated
1351on each node. The default for \verb+<madbox>+ is 6.
1352In some cases \verb+<madbox>+ can be reduced to 4 if memory usage is a
1353concern. Values of 2 or less will almost certainly result in memory
1354shortage.
1355
1356\item
1357\begin{verbatim}
1358mwm <integer mwmreq>
1359\end{verbatim}
1360sets the maximum number of solvent molecules \verb+<mwmreq>+ per node,
1361allowing increased memory to be allocated for solvent molecules. This
1362option can be used if execution aborted because \verb+mwm+ was too
1363small.
1364
1365\item
1366\begin{verbatim}
1367msa <integer msareq>
1368\end{verbatim}
1369sets the maximum number of solute atoms \verb+<msareq>+ per node,
1370allowing increased memory to be allocated for solute atoms. This
1371option can be used if execution aborted because \verb+msa+ was too
1372small.
1373
1374\item
1375\begin{verbatim}
1376mcells <integer mbbreq>
1377\end{verbatim}
1378sets the maximum number of cell pairs \verb+<mbbreq>+ per node,
1379allowing increased memory to be allocated for the cell pair lists.
1380This option can be used if execution aborted because \verb+mbbl+ was too
1381small.
1382
1383\item
1384\begin{verbatim}
1385boxmin <real rbox>
1386\end{verbatim}
1387sets the minimum size of a cell. This directive is obsolete. The
1388use of mcells is preferred.
1389
1390\item
1391\begin{verbatim}
1392segmentsize <real rsgm>
1393\end{verbatim}
1394sets the maximum size of a segment. This value is used to determine
1395which segments at the boundary of the cutoff radius should be considered
1396in the generation of the pairlists. This value is also determined by the
1397prepare module and written to the restart file. Use of this directive
1398is not needed for simulations that use the current prepare module to
1399generate the restart file.
1400
1401\item
1402\begin{verbatim}
1403memory <integer memlim>
1404\end{verbatim}
1405sets a limit \verb+<memlim>+ in kB on the allocated amount of memory used by
1406the molecular dynamics module.
1407Per default all available memory is allocated. Use of this command
1408is required for QM/MM simulations only.
1409
1410\item
1411\begin{verbatim}
1412expert
1413\end{verbatim}
1414enables the use of certain combinations of features that are considered
1415unsafe. This directive should not be used for production runs.
1416
1417\item
1418\begin{verbatim}
1419develop <integer idevel>
1420\end{verbatim}
1421enables the use of certain development options specified by the
1422integer \verb+idevel+. This option is for development purposes only,
1423and should not be used for production runs.
1424
1425\item
1426\begin{verbatim}
1427control <integer icntrl>
1428\end{verbatim}
1429enables the use of certain development options specified by the
1430integer \verb+icntrl+. This option is for development purposes only,
1431and should not be used for production runs.
1432
1433\item
1434\begin{verbatim}
1435numerical
1436\end{verbatim}
1437writes out analytical and finite difference forces for test purposes.
1438
1439\item
1440\begin{verbatim}
1441server <string servername> <integer serverport>
1442\end{verbatim}
1443allows monitoring over a socket connection to the specified port on the
1444named server of basic data as a simulation is running.
1445
1446\item
1447For development purposes debug information can be written to the debug
1448file with extension {\bf dbg} with
1449
1450\begin{verbatim}
1451debug <integer idebug>
1452\end{verbatim}
1453where $idebug$ specifies the type of debug information being written.
1454
1455\item
1456For testing purposes test information can be written to the test
1457file with extension {\bf tst} with
1458
1459\begin{verbatim}
1460test <integer itest>
1461\end{verbatim}
1462where $itest$ specifies the number of steps test information is written.
1463
1464\item
1465On some platforms prefetching of data can improve the efficiency. This
1466feature can be turned on using
1467
1468\begin{verbatim}
1469prefetch [<integer nbget>]
1470\end{verbatim}
1471where $nbget$ is the number of outstanding communication operations.
1472
1473\item
1474Application of periodic boundary conditions for the evaluation of
1475forces can be controlled with
1476
1477\begin{verbatim}
1478pbc ( atom | residue | molecule )
1479\end{verbatim}
1480
1481This option rarely needs to be used.
1482
1483\item
1484Autocorrelation functions for error analysis are controlled using
1485
1486\begin{verbatim}
1487auto [ fit <integer iapprx> | weight <real weight> ]
1488\end{verbatim}
1489
1490This option is disabled in the current release.
1491
1492\item
1493Membrane system equilibration can be made more efficient using
1494
1495\begin{verbatim}
1496membrane [ rotations ]
1497\end{verbatim}
1498
1499\item
1500Constraining the center of mass of solute molecules in the xy plane is
1501accomplished using
1502
1503\begin{verbatim}
1504scmxy [<integer icmopt default 1>]
1505\end{verbatim}
1506
1507where $icmopt$ determines if the constraint is mass weighted (2).
1508
1509\item
1510Radius of gyration calculations are enabled using
1511
1512\begin{verbatim}
1513radius_gyration
1514\end{verbatim}
1515
1516\item
1517Calculations of diffusion coefficients is enabled using
1518
1519\begin{verbatim}
1520diffusion
1521\end{verbatim}
1522This option is disabled in the current release.
1523
1524\item
1525
1526\begin{verbatim}
1527comlim ( on | off )
1528\end{verbatim}
1529is disabled
1530
1531\item
1532
1533To limit the size of recoding files, new files are opened every $nfnewf$ md steps using
1534
1535\begin{verbatim}
1536batch <integer nfnewf>
1537\end{verbatim}
1538
1539
1540\end{description}
1541