author Giacomo Fiorin Mon, 15 Oct 2018 18:08:11 +0000 (14:08 -0400) committer Giacomo Fiorin Mon, 15 Oct 2018 18:30:24 +0000 (14:30 -0400)
Addresses the latex2html issue for the NAMD user's guide, and also improves
the learning curve by defining the variable functions earlier in the chapter.

The resulting diff is fairly large, mostly due to having merged one chapter
with another chapter that was not consecutive.

https://github.com/Colvars/colvars/pull/183

Change-Id: I1e6a2005c62efc863bf47c98f1a81ae951535883

 ug/Makefile patch | blob | history ug/ug_colvars-cv.tex [deleted file] patch | blob | history ug/ug_colvars.tex patch | blob | history

index 0a169e7..a78c8ea 100644 (file)
@@ -29,7 +29,7 @@ namd_title.tex namd_version.tex ug.tex ug_accel.tex ug_avail.tex \
ug_dynamics.tex ug_forcefield.tex ug_performance.tex ug_userdef.tex \
ug_intro.tex ug_io.tex ug_macros.tex ug_runit.tex ug_sample.tex ug_start.tex \
ug_xplor.tex ug_fenergy.tex ug_psfgen.tex psfgen_macros.tex ug_alchemy.tex \
-ug_colvars-cv.tex colvars-cv.tex ug_colvars.bib ug_colvars_macros.tex \
+ug_colvars.bib ug_colvars_macros.tex \
ug_colvars.tex ug_analysis.tex ug_gbis.tex ug_constantph.tex

PDFFIGS = figures/fmaOn.pdf figures/pairlistdist.pdf \
diff --git a/ug/ug_colvars-cv.tex b/ug/ug_colvars-cv.tex
deleted file mode 100644 (file)
index 136d2bb..0000000
+++ /dev/null
@@ -1,108 +0,0 @@
-% This file is part of the Collective Variables module (Colvars).
-% The original version of Colvars and its updates are located at:
-% https://github.com/colvars/colvars
-% Please update all Colvars source files before making any changes.
-% If you wish to distribute your changes, please submit them to the
-% Colvars repository at GitHub.
-
-Collective variables and biases can be added, queried and deleted through the scripting command \texttt{cv}, with the following syntax: \texttt{cv~$<$subcommand$>$ [args...]}.
-For example, to query the value of a collective variable named \texttt{myVar},
-use the following syntax: \texttt{set value [cv colvar myVar value]}.
-All subcommands of \texttt{cv} are documented below.
-
-% The NAMD and VMD manuals have this section at different levels
-\newcommand{\cvparagraph}[2]{\cvsubsec{#1}{#2}}
-\cvvmdonly{
-\renewcommand{\cvparagraph}[2]{\cvsubsubsec{#1}{#2}}
-}
-
-\cvvmdonly{
-\cvparagraph{Example use of the \texttt{cv} command: analyze a trajectory}{sec:cv_command_example}
-
-By far the most typical use of Colvars in VMD is computing the values of one or more variables along an existing trajectory:
-\bigskip
-{\noindent\ttfamily\\
-\-{\bfseries\# Activate the module on the current VMD molecule}\\
-\-cv molid top\\
-\-{\bfseries\# Load a Colvars config file}\\
-\-cv configfile test.in\\
-\-set out [open "test.colvars.traj" "w"]\\
-\-{\bfseries\# Write the labels to the file}\\
-\-puts -nonewline \$\{out\} [cv printframelabels]\\ -\-for \{ set fr 0 \} \{ \$\{fr\} < [molinfo top get numframes] \} \{ incr fr \} \{\\
-\-~~{\bfseries\# Point Colvars to this trajectory frame}\\
-\-~~cv frame \$\{fr\}\\ -\-~~{\bfseries\# Recompute variables and biases}\\ -\-~~cv update\\ -\-~~{\bfseries\# Print variables and biases to the file}\\ -\-~~puts -nonewline \$\{out\} [cv printframe]\\
-\-\}\\
-\-close \$\{out\}\\ -} - -The following sections explain the syntax of the various sub-commands of \texttt{cv}. -} - -\cvparagraph{Managing the Colvars module}{sec:cv_command_general} - -\begin{itemize} -\cvvmdonly{ -\item \texttt{molid} \emph{$<$molid$>$}: setup the Colvars module for the given molecule; \emph{$<$molid$>$} is the identifier of a valid VMD molecule, either an integer or \texttt{top} for the top molecule; -\item \texttt{delete}: delete the Colvars module; after this, the module can be setup again using \texttt{cv molid}; -} -\item \texttt{configfile} \emph{$<$file name$>$}: read configuration from a file; -\item \texttt{config} \emph{$<$string$>$}: read configuration from the given string; both \texttt{config} and \texttt{configfile} subcommands may be invoked multiple times; -\cvvmdonly{see \ref{sec:colvars_config} for the configuration syntax;} -\item \texttt{reset}: delete all internal configuration of the Colvars module; -\item \texttt{version}: return the version of the colvars code. -\end{itemize} - -\cvparagraph{Input and output}{sec:cv_command_io} - -\begin{itemize} -\item \texttt{list}: return a list of all currently defined variables; -\cvvmdonly{See \ref{sec:colvar}, \ref{sec:cvc} and \ref{sec:colvar_atom_groups} for how to configure a collective variable;} -\item \texttt{list biases}: return a list of all currently defined biases (i.e.~sampling and analysis algorithms); -\item \texttt{load} \emph{$<$file name$>$}: load a collective variables state file, typically produced during a previous simulation; \cvvmdonly{this parameter requires that the corresponding configuration has already been loaded by \texttt{configfile} or \texttt{config}; see \ref{sec:colvars_output} for a brief description of this file; the contents of this file are not required for as long as the VMD molecule has valid coordinates and \texttt{cv update} is used;} -\item \texttt{load} \emph{$<$prefix$>$}: same as above, but without the \texttt{.colvars.state} suffix; -\item \texttt{save} \emph{$<$prefix$>$}: save the current state in a file whose name begins with the given argument; if any of the biases have additional output files defined, those are saved as well using the same prefix; -\item \texttt{update}: recalculate all colvars and biases based on the current atomic coordinates; -\item \texttt{addenergy} \emph{$<$E$>$}: add value \emph{E} to the total bias energy; this must be used within \texttt{calc\_colvar\_forces}; -\item \texttt{printframe}: return a summary of the current frame, in a format equivalent to a line of the collective variables trajectory file; -\item \texttt{printframelabels}: return text labels for the columns of \texttt{printframe}'s output; -\cvvmdonly{ -\item \texttt{frame}: return the current frame number, from which colvar values are calculated; -\item \texttt{frame} \emph{$<$new frame$>$}: set the frame number; returns 0 if successful, nonzero if the requested frame does not exist. -} -\end{itemize} - -\cvparagraph{Accessing collective variables}{sec:cv_command_colvar} - -\begin{itemize} -\item \texttt{colvar} \emph{$<$name$>$} \texttt{value}: return the current value of colvar \emph{$<$name$>$}; -\item \texttt{colvar} \emph{$<$name$>$} \texttt{update}: recalculate colvar \emph{$<$name$>$}; -\item \texttt{colvar} \emph{$<$name$>$} \texttt{type}: return the type of colvar \emph{$<$name$>$}; -\cvvmdonly{ -\item \texttt{colvar} \emph{$<$name$>$} \texttt{delete}: delete colvar \emph{$<$name$>$}; -} -\cvnamdonly{ -\item \texttt{colvar} \emph{$<$name$>$} \texttt{addforce} \emph{$<$F$>$}: apply given force on colvar \emph{$<$name$>$}; -\item \texttt{colvar} \emph{$<$name$>$} \texttt{getappliedforce}: get the force currently being applied on colvar \emph{$<$name$>$}; -\item \texttt{colvar} \emph{$<$name$>$} \texttt{gettotalforce}: get the total force acting on colvar \emph{$<$name$>$} at the previous step (see \ref{sec:cvc_sys_forces}); -} -\item \texttt{colvar} \emph{$<$name$>$} \texttt{getconfig}: return config string of colvar \emph{$<$name$>$}. -\item \texttt{colvar} \emph{$<$name$>$} \texttt{cvcflags} \emph{$<$flags$>$}: for a colvar with several CVCs (numbered according to their name string order), set which CVCs are enabled or disabled in subsequent evaluations according to a list of 0/1 flags (one per CVC). -\item \texttt{colvar} \emph{$<$name$>$} \texttt{modifycvcs} \emph{$<$strings$>$}: for a colvar with one or more CVCs (numbered according to their name string order), pass a list of new configuration strings to modify each CVC without needing to delete the colvar. -This option is currently limited to changing the values of \refkey{componentCoeff}{sec:cvc_superp} and \refkey{componentExp}{sec:cvc_superp} (e.g.{} to update the polynomial superposition parameters on the fly), of \refkey{period}{sec:cvc_periodic} and \refkey{wrapAround}{sec:cvc_periodic} (e.g.{} to update the period of variables such as \texttt{distanceZ}), and of the \texttt{forceNoPBC} option for those CVCs that support it. -Changes in configuration done by this command are not saved to state files, and are lost when restarting a simulation, deleting the parent colvar, or resetting the module with \cvvmdonly{\texttt{cv delete} or} \texttt{cv reset}. - -\end{itemize} - -\cvparagraph{Accessing biases}{sec:cv_command_bias} - -\begin{itemize} -\item \texttt{bias} \emph{$<$name$>$} \texttt{energy}: return the current energy of the bias \emph{$<$name$>$}; -\item \texttt{bias} \emph{$<$name$>$} \texttt{update}: recalculate the bias \emph{$<$name$>$}; -\item \texttt{bias} \emph{$<$name$>$} \texttt{delete}: delete the bias \emph{$<$name$>$}; -\item \texttt{bias} \emph{$<$name$>$} \texttt{getconfig}: return config string of bias \emph{$<$name$>$}. -\end{itemize} index 9e0f0fc..b1c0097 100644 (file) @@ -43,28 +43,25 @@ Please ask any usage questions through the VMD mailing list, and development que In molecular dynamics simulations, it is often useful to reduce the large number of degrees of freedom of a physical system into few parameters whose statistical distributions can be analyzed individually, or used to define biasing potentials to alter the dynamics of the system in a controlled manner. These have been called order parameters', collective variables', (surrogate) reaction coordinates', and many other terms. -Here we use primarily the term collective variable' (shortened to \textit{colvar}), which indicates any differentiable function of atomic Cartesian coordinates,$\bm{x}_{i}$, with$i$between$1$and$N$, the total +Here we use primarily the term collective variable', often shortened to \textit{colvar}, to indicate any differentiable function of atomic Cartesian coordinates,$\bm{x}_{i}$, with$i$between$1$and$N$, the total number of atoms: \label{eq:colvar_basic} - \xi(t) \; = \xi\left(\bm{x}_{i}(t), \bm{x}_{j}(t), \bm{x}_{k}(t), + \xi(t) \; = \xi(\bm{X}(t)) \; = \xi\left(\bm{x}_{i}(t), \bm{x}_{j}(t), \bm{x}_{k}(t), \ldots \right)\;, \;\; 1 \leq i,j,k\ldots \leq N -\cvnamdugonly{ -The Colvars module in NAMD may be used in both MD simulations and energy minimization runs. -} \cvvmdugonly{% The Colvars module in VMD may be used to calculate these functions over a molecular structure, and to analyze the results of previous simulations.} \cvrefmanonly{% -This manual documents the collective variables module (\textbf{Colvars}), a portable software that interfaces multiple MD simulation programs, with a focus on flexibility, robustness and high performance.} +This manual documents the collective variables module (\textbf{Colvars}), a software that provides an implementation for the functions$\xi(\bm{X})$with a focus on flexibility, robustness and high performance.} The module is designed to perform multiple tasks concurrently during or after a simulation, the most common of which are: \begin{itemize} -\item apply restraints or biasing potentials to multiple colvars, tailored on the system by choosing from a wide set of basis functions, without limitations on their number or on the number of atoms involved; \cvnamdonly{while this can in principle be done through a TclForces script, using the Colvars module is both easier and computationally more efficient;} +\item apply restraints or biasing potentials to multiple variables, tailored on the system by choosing from a wide set of basis functions, without limitations on their number or on the number of atoms involved; \cvnamdonly{while this can in principle be done through a TclForces script, using the Colvars module is both easier and computationally more efficient;} -\item calculate potentials of mean force (PMFs) along any set of colvars, using different enhanced sampling methods, such as Adaptive Biasing Force (ABF), metadynamics, steered MD and umbrella sampling; variants of these methods that make use of an ensemble of replicas are supported as well; +\item calculate potentials of mean force (PMFs) along any set of variables, using different enhanced sampling methods, such as Adaptive Biasing Force (ABF), metadynamics, steered MD and umbrella sampling; variants of these methods that make use of an ensemble of replicas are supported as well; -\item calculate statistical properties of the colvars, such as running averages and standard deviations, correlation functions of pairs of colvars, and multidimensional histograms: this can be done either at run-time without the need to save very large trajectory files, or after a simulation has been completed using VMD and the \texttt{cv} command\cvnamdugonly{ or NAMD and the \texttt{coorfile read} command as illustrated in \ref{section:sample}}. +\item calculate statistical properties of the variables, such as running averages and standard deviations, correlation functions of pairs of variables, and multidimensional histograms: this can be done either at run-time without the need to save very large trajectory files, or after a simulation has been completed using VMD and the \texttt{cv} command\cvnamdugonly{ or NAMD and the \texttt{coorfile read} command as illustrated in \ref{section:sample}}. \end{itemize} @@ -72,13 +69,14 @@ The module is designed to perform multiple tasks concurrently during or after a Detailed explanations of the design of the Colvars module are provided in reference~\cite{Fiorin2013}. Please cite this reference whenever publishing work that makes use of this module. + \cvsec{A crash course}{sec:colvars_crash_course} % FIXME: the example uses NAMD syntax to select Calphas, which may confuse % LAMMPS users -Suppose that we want to run a steered MD experiment where a small molecule is pulled away from a binding site. +Suppose that we want to run a steered MD experiment where a small molecule is pulled away from a protein binding site. In Colvars terms, this is done by applying a moving restraint to the distance between the two objects. -The configuration will contain two blocks, one defining the distance variable (see sections \ref{sec:colvar} and \ref{sec:cvc}), and the other the moving harmonic restraint (\ref{sec:colvarbias_harmonic}). +The configuration will contain two blocks, one defining the distance variable (see section \ref{sec:colvar} and \ref{sec:cvc_distance}), and the other the moving harmonic restraint (\ref{sec:colvarbias_harmonic}). \bigskip % verbatim can't appear within commands @@ -109,12 +107,11 @@ The atom selection syntax is detailed in section \ref{sec:colvar_atom_groups}. To the \emph{dist}'' variable, we apply a \texttt{harmonic} potential of force constant 20~\cvnamdonly{kcal/mol}\cvvmdonly{kcal/mol}\cvlammpsonly{energy units}/\AA$^2$, initially centered around a value of 4~\AA, which will increase to 15~\AA{} over 500,000 simulation steps. -\cvsec{General parameters and input/output files}{sec:colvarmodule} +\cvsec{General parameters}{sec:colvarmodule} Here, we document the syntax of the commands and parameters used to set up and use the Colvars module in \MDENGINE{}. -One of these parameters is the configuration file or the configuration text for the module itself, whose syntax is described in \ref{sec:colvars_config} and in the following sections. +One of these parameters is the configuration file or the configuration text for the module itself, whose syntax is described in \ref{sec:colvars_config_syntax} and in the following sections. -% Removed VMD section by hand to make latex2html work correctly \cvnamdonly{ \cvsubsec{NAMD parameters}{sec:colvars_mdengine_params} @@ -145,24 +142,6 @@ An optional third parameter, \texttt{colvarsInput}, can be used to continue a pr This file contains the definition of all collective variables and their biasing or analysis methods. This file can also be provided by the Tcl command \texttt{cv configfile}; alternatively, the contents of the file itself can be given as an argument to the command \texttt{cv config}. - % Parameters within the configuration file can be controlled from - % a NAMD config file using Tcl variables in the following way: - - % {\ttfamily - % colvars on\\ - % colvarsConfig colvars\_subst.tmp\\ - % set myParameter someValue\\ - % \# Parse template and create specific config file on the fly\\ - % set infile [open colvars\_template.in r] \\ - % set outfile [open colvars\_subst.tmp w+] \\ - % puts \$outfile [subst [read \$infile]] \\ - % close \$infile \\
-    % close \$outfile} - - % In this example, the string \texttt{\$myParameter} will be replaced
-    % with the value \texttt{someValue} wherever it appears in the file
-    % \texttt{colvars\_template.in}. This value will then be read in by
-    % the Colvars module when it parses its input.
}

\item %
@@ -173,12 +152,14 @@ An optional third parameter, \texttt{colvarsInput}, can be used to continue a pr
UNIX filename}{%
When continuing a previous simulation run, this file contains the current state of all collective variables and of their associated algorithms.
It is written automatically at the end of any simulation with collective variables.
-    This file can also be provided by the Tcl command \texttt{cv load}.
+    The step number contained by this file will be used internally by Colvars to control time-dependent biases, unless \texttt{firstTimestep} is given, in which case that value will be used.
+    This file can also be loaded with the Tcl command \texttt{cv load}.
}

\end{itemize}
}

+
\cvlammpsonly{
\subsection{LAMMPS keywords}
\label{sec:colvars_mdengine_parameters}
@@ -198,7 +179,9 @@ where \emph{ID} is a string that uniquely identifies this fix command inside a L
Name or prefix of the input state file}{%
string}{%
If a value is provided, it is interpreted as either the name of the input state file, or as the prefix of the file named \emph{input}\texttt{.colvars.state}.
-This allows to continue a previous collective variables-based calculation when a regular binary LAMMPS restart file is not available (see \ref{sec:colvars_input}).}
+    This file contains information needed to continue a previous collective variables-based calculation, including the number of the last computed step (useful for time-dependent biases).
+    The same information is also stored in the binary restart files written by LAMMPS, so this option is not needed when continuing a calculation from a LAMMPS restart.
+  }

\item %
\keydef
@@ -239,10 +222,126 @@ Setting this to \emph{no} will use the current local coordinates that are wrappe
This keyword provides the \emph{ID} of an applicable thermostating fix command. This will be used to provide the Colvars module with the current thermostat target temperature when using a method that needs this information.}

\end{itemize}
+
+}
+
+
+\cvscriptonly{
+\cvsubsec{Using the \texttt{cv} command}{sec:cv_scripting}
+\cvvmdugonly{\index{Colvars!\texttt{cv} command}}
+
+At any moment during the execution of \MDENGINE{}, several options in the Colvars module can be read or modified by the command \texttt{cv} with the following syntax:\\
+\noindent\texttt{cv~$<$subcommand$>$ [args...]}\\
+\noindent{}For example, to record the value of a collective variable named \texttt{myVar} into the Tcl variable \texttt{value}, use the following syntax:\\
+\noindent\texttt{set value [cv colvar myVar value]}\\
+\noindent{}All subcommands of \texttt{cv} are documented below.
+
+\cvvmdonly{
+\textbf{Note:} in VMD, Colvars must be attached to one molecule (system).
+Therefore, the \texttt{cv} command must be used for the first time as \texttt{cv molid }\emph{$<$molid$>$} to set up the Colvars module for the molecule identified by \emph{$<$molid$>$}.
+In all following invocations, the \texttt{cv} command will continue operating on the same molecule, regardless of its top'' status.
+To use the \texttt{cv} command on a different molecule, use \texttt{cv delete} first and then \texttt{cv molid }\emph{$<$molid$>$}.
+Invoking the \texttt{cv} command with no arguments prints a help screen.
+}
+
+\cvvmdonly{
+\cvsubsubsec{Example use of the \texttt{cv} command: analyze a trajectory}{sec:cv_command_example}
+
+By far the most typical use of Colvars in VMD is computing the values of one or more variables along an existing trajectory:
+\bigskip
+{\noindent\ttfamily\\
+\-{\bfseries\# Activate the module on the current VMD molecule}\\
+\-cv molid top\\
+\-{\bfseries\# Load a Colvars config file}\\
+\-cv configfile test.in\\
+\-set out [open "test.colvars.traj" "w"]\\
+\-{\bfseries\# Write the labels to the file}\\
+\-puts -nonewline \$\{out\} [cv printframelabels]\\ +\-for \{ set fr 0 \} \{ \$\{fr\} < [molinfo top get numframes] \} \{ incr fr \} \{\\
+\-~~{\bfseries\# Point Colvars to this trajectory frame}\\
+\-~~cv frame \$\{fr\}\\ +\-~~{\bfseries\# Recompute variables and biases}\\ +\-~~cv update\\ +\-~~{\bfseries\# Print variables and biases to the file}\\ +\-~~puts -nonewline \$\{out\} [cv printframe]\\
+\-\}\\
+\-close \$\{out\}\\ +} + +The following sections explain the syntax of the various sub-commands of \texttt{cv}. +} + + +\cvsubsubsec{Managing the Colvars module}{sec:cv_command_general} + +\begin{itemize} +\cvvmdonly{ +\item \texttt{molid} \emph{$<$molid$>$}: setup the Colvars module for the given molecule; \emph{$<$molid$>$} is the identifier of a valid VMD molecule, either an integer or \texttt{top} for the top molecule; +\item \texttt{delete}: delete the Colvars module; after this, the module can be setup again using \texttt{cv molid}; +} +\item \texttt{configfile} \emph{$<$file name$>$}: read configuration from a file; +\item \texttt{config} \emph{$<$string$>$}: read configuration from the given string; both \texttt{config} and \texttt{configfile} subcommands may be invoked multiple times; +\cvvmdonly{see \ref{sec:colvars_config_syntax} for the configuration syntax;} +\item \texttt{reset}: delete all internal configuration and atom selections of the Colvars module; +\item \texttt{version}: return the version of the Colvars module (see \ref{sec:colvars_config_changes} for any changes to older keywords). +\end{itemize} + +\cvsubsubsec{Input and output commands}{sec:cv_command_io} + +\begin{itemize} +\item \texttt{list}: return a list of all currently defined variables; +\cvvmdonly{See \ref{sec:colvar} and \ref{sec:colvar_atom_groups} for how to configure a collective variable;} +\item \texttt{list biases}: return a list of all currently defined biases (i.e.~sampling and analysis algorithms); +\item \texttt{load} \emph{$<$file name$>$}: load a collective variables state file, typically produced during a previous simulation; + \cvvmdonly{this parameter requires that the corresponding configuration has already been loaded by \texttt{configfile} or \texttt{config}; see \ref{sec:colvars_output} for a brief description of this file; the contents of this file are not required for as long as the VMD molecule has valid coordinates and \texttt{cv update} is used;} + \cvnamdonly{the step number contained by this file will be used internally by Colvars to control time-dependent biases, unless \texttt{firstTimestep} is given, in which case that value will be used;} +\item \texttt{load} \emph{$<$prefix$>$}: same as above, but without the \texttt{.colvars.state} suffix; +\item \texttt{save} \emph{$<$prefix$>$}: save the current state in a file whose name begins with the given argument; if any of the biases have additional output files defined, those are saved as well using the same prefix; +\item \texttt{update}: recalculate all variables and biases based on the current atomic coordinates\cvnamdonly{ (this is typically not needed in NAMD, as these functions are called already during a timestep)}; +\item \texttt{addenergy} \emph{$<$E$>$}: add value \emph{E} to the total bias energy; this can be used within \texttt{calc\_colvar\_forces}; +\item \texttt{printframe}: return a summary of the current frame, in a format equivalent to a line of the collective variables trajectory file; +\item \texttt{printframelabels}: return text labels for the columns of \texttt{printframe}'s output; +\cvvmdonly{ +\item \texttt{frame}: return the current frame number, from which colvar values are calculated; +\item \texttt{frame} \emph{$<$new frame$>$}: set the frame number; returns 0 if successful, nonzero if the requested frame does not exist. +} +\end{itemize} + +\cvsubsubsec{Managing collective variables}{sec:cv_command_colvar} + +\begin{itemize} +\item \texttt{colvar} \emph{$<$name$>$} \texttt{value}: return the current value of colvar \emph{$<$name$>$}; +\item \texttt{colvar} \emph{$<$name$>$} \texttt{update}: recalculate colvar \emph{$<$name$>$}; +\item \texttt{colvar} \emph{$<$name$>$} \texttt{type}: return the type of colvar \emph{$<$name$>$}; +\cvvmdonly{ +\item \texttt{colvar} \emph{$<$name$>$} \texttt{delete}: delete colvar \emph{$<$name$>$}; +} +\cvnamdonly{ +\item \texttt{colvar} \emph{$<$name$>$} \texttt{addforce} \emph{$<$F$>$}: apply given force on colvar \emph{$<$name$>$}; this can be used within \texttt{calc\_colvar\_forces} +\item \texttt{colvar} \emph{$<$name$>$} \texttt{getappliedforce}: get the force currently being applied on colvar \emph{$<$name$>$}; +\item \texttt{colvar} \emph{$<$name$>$} \texttt{gettotalforce}: get the total force acting on colvar \emph{$<$name$>$} at the previous step (see \ref{sec:cvc_sys_forces}); +} +\item \texttt{colvar} \emph{$<$name$>$} \texttt{getconfig}: return config string of colvar \emph{$<$name$>$}. +\item \texttt{colvar} \emph{$<$name$>$} \texttt{cvcflags} \emph{$<$flags$>$}: for a colvar with several CVCs (numbered according to their name string order), set which CVCs are enabled or disabled in subsequent evaluations according to a list of 0/1 flags (one per CVC). +\item \texttt{colvar} \emph{$<$name$>$} \texttt{modifycvcs} \emph{$<$strings$>$}: for a colvar with one or more CVCs (numbered according to their name string order), pass a list of new configuration strings to modify each CVC without needing to delete the colvar. +This option is currently limited to changing the values of \refkey{componentCoeff}{sec:cvc_superp} and \refkey{componentExp}{sec:cvc_superp} (e.g.{} to update the polynomial superposition parameters on the fly), of \refkey{period}{sec:cvc_periodic} and \refkey{wrapAround}{sec:cvc_periodic} (e.g.{} to update the period of variables such as \texttt{distanceZ}), and of the \texttt{forceNoPBC} option for those CVCs that support it. +Changes in configuration done by this command are not saved to state files, and are lost when restarting a simulation, deleting the parent colvar, or resetting the module with \cvvmdonly{\texttt{cv delete} or} \texttt{cv reset}. + +\end{itemize} + +\cvsubsubsec{Managing biases}{sec:cv_command_bias} + +\begin{itemize} +\item \texttt{bias} \emph{$<$name$>$} \texttt{energy}: return the current energy of the bias \emph{$<$name$>$}; +\item \texttt{bias} \emph{$<$name$>$} \texttt{update}: recalculate the bias \emph{$<$name$>$}; +\item \texttt{bias} \emph{$<$name$>$} \texttt{delete}: delete the bias \emph{$<$name$>$}; +\item \texttt{bias} \emph{$<$name$>$} \texttt{getconfig}: return config string of bias \emph{$<$name$>$}. +\end{itemize} + } -\cvsubsec{Configuration syntax for the Colvars module}{sec:colvars_config} +\cvsubsec{Configuration syntax}{sec:colvars_config_syntax} \cvnamdonly{All the parameters defining variables and their biasing or analysis algorithms are read from the file specified by the configuration option \texttt{colvarsConfig}, or by the Tcl commands \texttt{cv config} and \texttt{cv configfile}. \emph{None of the keywords described in the remainder of this manual are recognized directly in the NAMD configuration file, unless as arguments of \texttt{cv config}.}} @@ -273,6 +372,9 @@ The following rules apply: \end{itemize} + +\cvsubsec{Global keywords}{sec:colvars_global} + The following keywords are available in the global context of the Colvars configuration, i.e.~they are not nested inside other keywords: \begin{itemize} @@ -345,7 +447,7 @@ The following keywords are available in the global context of the Colvars config To illustrate the flexibility of the Colvars module, a non-trivial setup is represented in Figure~\ref{fig:colvars_diagram}. -The corresponding configuration is given below. The options within the \texttt{colvar} blocks are described in \ref{sec:colvar} and \ref{sec:cvc}, those within the \texttt{harmonic} and \texttt{histogram} blocks in \ref{sec:colvarbias}. +The corresponding configuration is given below. The options within the \texttt{colvar} blocks are described in \ref{sec:colvar}, those within the \texttt{harmonic} and \texttt{histogram} blocks in \ref{sec:colvarbias}. \textbf{Note:} \emph{except }\texttt{colvar}\emph{, none of the keywords shown is mandatory}.\\ \begin{figure}[!ht] @@ -415,75 +517,65 @@ harmonic \{\\ %\cvlammpsonly{\textbf{Note:} \emph{currently, the$\alpha$-helical content colvar is unavailable in LAMMPS, as it requires a name-based topology; future releases will overcome this limitation.}} Section \ref{sec:colvar} explains how to define a colvar and its behavior, regardless of its specific functional form. -To define colvars that are appropriate to a specific physical system, Section \ref{sec:colvar_atom_groups} documents how to select atoms, and section \ref{sec:cvc} lists all of the available functional forms, which we call colvar components''. +To define colvars that are appropriate to a specific physical system, Section \ref{sec:colvar_atom_groups} documents how to select atoms, and section \ref{sec:colvar} lists all of the available functional forms, which we call colvar components''. Finally, section \ref{sec:colvarbias} lists the available methods and algorithms to perform biased simulations and multidimensional analysis of colvars. -\cvsubsec{Input state file (optional)}{sec:colvars_input} - -Aside from the Colvars configuration, an optional input state file may be provided to load the relevant data from a previous simulation. -\cvnamdonly{The name of this file is provided as a value to the keyword \texttt{colvarsInput}.} -\cvvmdonly{The name of this file is provided through \texttt{cv load }\emph{$<$file name$>$}.} -\cvlammpsonly{The name of this file is provided as the argument to the \texttt{input} keyword of the \texttt{fix ID all colvars} command. The same information is stored in the binary restart files of LAMMPS, so it not needed when continuing a calculation from such a restart.} -This file contains information about the current run, including the number of the last computed step. \cvnamdonly{Within Colvars, this number will be carried over to the new simulation, unless \texttt{firstTimestep} is given: in this case, the value of \texttt{firstTimestep} will always be used as first step by Colvars as well}. - - \cvsubsec{Output files}{sec:colvars_output} During a simulation with collective variables defined, the following three output files are written: \begin{itemize} -\item a \emph{state file}, named \outputName\texttt{.colvars.state}; this file is in ASCII (plain text) format\cvnamdonly{, regardless of the value of \texttt{binaryOutput} in the NAMD configuration}, and is meant to be used as input when continuing a simulation (see \ref{sec:colvars_input}). +\item A \emph{state file}, named \outputName\texttt{.colvars.state}; this file is in ASCII (plain text) format\cvnamdonly{, regardless of the value of \texttt{binaryOutput} in the NAMD configuration}. This file is written at the end of the specified run\cvscriptonly{, but can also be written at any time with the command \texttt{cv save} (\ref{sec:cv_command_io})}.\\ + \emph{This is the only Colvars output file needed to continue a simulation.} -\item if \cvnamdonly{the NAMD parameter \texttt{restartFreq} or } the parameter \texttt{colvarsRestartFrequency} is larger than zero, a \emph{restart file} named \restartName\texttt{.colvars.state} is written every that many steps: this file is equivalent to the final state file; +\item If the parameter \refkey{colvarsRestartFrequency}{Colvars-global|colvarsRestartFrequency} is larger than zero, a \emph{restart file} is written every that many steps: this file is fully equivalent to the final state file. + The name of this file is \restartName\texttt{.colvars.state}. -\item if the parameter \texttt{colvarsTrajFrequency} is greater than 0 (default: 100), a \emph{trajectory file} is written during the simulation: its name is \outputName\texttt{.colvars.traj}; unlike the state file, it is not needed to restart a simulation, but can be used later for post-processing and analysis. +\item If the parameter \refkey{colvarsTrajFrequency}{Colvars-global|colvarsTrajFrequency} is greater than 0 (default: 100), a \emph{trajectory file} is written during the simulation: its name is \outputName\texttt{.colvars.traj}; unlike the state file, it is not needed to restart a simulation, but can be used later for post-processing and analysis. \end{itemize} -Other output files may be written by specific methods applied to the colvars (e.g.~by the ABF method, see \ref{sec:colvarbias_abf}, or the metadynamics method, see \ref{sec:colvarbias_meta}). -Like the colvar trajectory file, they are needed only for analyzing, not continuing a simulation. -All such files' names also begin with the prefix \outputName\texttt{}. +Other output files may also be written by specific methods, e.g.{} the ABF or metadynamics methods (\ref{sec:colvarbias_abf}, \ref{sec:colvarbias_meta}). +Like the trajectory file, they are needed only for analyzing, not continuing a simulation. +All such files' names also begin with the prefix \outputName. -\cvnamdonly{Finally, the total energy of all biases or restraints applied to the colvars appears under the NAMD standard output, under the MISC column.} +\cvnamdonly{Lastly, the total energy of all biases or restraints applied to the colvars appears under the NAMD standard output, under the MISC column.} -\cvsec{Defining and setting up collective variables}{sec:colvar} +\cvsec{Defining collective variables}{sec:colvar} A collective variable is defined by the keyword \texttt{colvar} followed by its configuration options contained within curly braces: {\bigskip\noindent\ttfamily colvar \{\\ \-~~name xi\\ -\-~~<other options>\\ -\-~~\ldots\\ +\-~~$<$other options$>$\\ +\-~~function\_name \{\\ +\-~~~~$<$parameters$>$\\ +\-~~~~$<$atom selection$>$\\ +\-~~\}\\ \}\\ } -There are multiple ways of defining a variable: +\noindent{}There are multiple ways of defining a variable: \begin{itemize} -\item The \emph{simplest and most common way} way is using one of the precompiled functions (here called components''), which are listed in section~\ref{sec:cvc}. For example, using the keyword \texttt{rmsd} (section \ref{sec:cvc_rmsd}) defines the variable as the root mean squared deviation (RMSD) of the selected atoms. -\item A new variable may also be constructed as a linear or polynomial combination of the components listed in section~\ref{sec:cvc} (see \ref{sec:cvc_superp} for details). +\item The \emph{simplest and most common way} way is using one of the precompiled functions (here called components''), which are listed in section~\ref{sec:cvc_list}. For example, using the keyword \texttt{rmsd} (section \ref{sec:cvc_rmsd}) defines the variable as the root mean squared deviation (RMSD) of the selected atoms. +\item A new variable may also be constructed as a linear or polynomial combination of the components listed in section~\ref{sec:cvc_list} (see \ref{sec:cvc_superp} for details). \cvleptononly{ -\item A user-defined mathematical function of the existing components (see list in section~\ref{sec:cvc}), or of the atomic coordinates directly (see the \texttt{cartesian} keyword in \ref{sec:cvc_cartesian}). +\item A user-defined mathematical function of the existing components (see list in section~\ref{sec:cvc_list}), or of the atomic coordinates directly (see the \texttt{cartesian} keyword in \ref{sec:cvc_cartesian}). The function is defined through the keyword \refkey{customFunction}{colvar|customFunction} (see \ref{sec:colvar_custom_function} for details). } \cvscriptonly{ -\item A user-defined Tcl function of the existing components (see list in section~\ref{sec:cvc}), or of the atomic coordinates directly (see the \texttt{cartesian} keyword in \ref{sec:cvc_cartesian}). +\item A user-defined Tcl function of the existing components (see list in section~\ref{sec:cvc_list}), or of the atomic coordinates directly (see the \texttt{cartesian} keyword in \ref{sec:cvc_cartesian}). The function is provided by a separate Tcl script, and referenced through the keyword \refkey{scriptedFunction}{colvar|scriptedFunction} (see \ref{sec:colvar_scripted} for details). } \end{itemize} - -All components listed in section~\ref{sec:cvc}, including the \texttt{cartesian} component make use of atom selection keyword described in section~\ref{sec:colvar_atom_groups}. - -The remainder of this section lists options to control the behavior of a variable, regardless of how it was defined. -Only the options \texttt{name}, \texttt{width}, \texttt{lowerBoundary} and \texttt{upperBoundary} are the most commonly used (\ref{sec:colvar_general}). - -\cvsubsec{General options for a collective variable}{sec:colvar_general} - +Choosing a component (function) is the only parameter strictly required to define a collective variable. +It is also highly recommended to specify a name for the variable: \begin{itemize} - +\label{sec:colvar_general} \item % \labelkey{colvar|name} \keydef Only the options \texttt{name}, \texttt{width}, \texttt{lowerBoundary} and \text used in the trajectory file to label to the columns corresponding to this colvar.} -\item % - \labelkey{colvar|width} - \keydef - {width}{% - \texttt{colvar}}{% - Colvar fluctuation scale, or resolution for grid-based methods}{% - positive decimal}{% - 1.0}{% - This number has the same physical unit as the colvar value and defines an effective colvar unit. -% Biasing algorithms use it for different purposes. - Harmonic restraints (\ref{sec:colvarbias_harmonic}) use it to set the physical unit of the force constant, which is useful for multidimensional restraints involving variables with very different units (for examples,$\AA$or degrees$^\circ$) with a single, scaled force constant. - Histograms (\ref{sec:colvarbias_histogram}), ABF (\ref{sec:colvarbias_abf}) and metadynamics (\ref{sec:colvarbias_meta}) all use this number as the initial choice for the grid spacing along this variable: for this reason, \texttt{width} should generally be no larger than the standard deviation of the colvar in an unbiased simulation. -% When a non-unity width is required by the application, the optimal value is application-dependent, but can often be thought of as a user-provided estimate of the fluctuation amplitude for the colvar. - Unless it is required to control the spacing, it is usually simplest to keep the default value of 1, so that restraint force constants are provided in more intuitive units. - } +\end{itemize} -\item % - \labelkey{colvar|lowerBoundary} - \key - {lowerBoundary}{% - \texttt{colvar}}{% - Lower boundary of the colvar}{% - decimal}{% - Defines the lowest end of the interval of relevant'' values for the colvar. - This number can be either a true physical boundary, or a user-defined number. - Together with \texttt{upperBoundary} and \texttt{width}, it is used to define a grid of values along the colvar (not available for colvars based on \texttt{distanceDir}, \texttt{distanceVec}, and \texttt{orientation}). - This option does not affect dynamics: to confine a colvar within a certain interval, use a \texttt{harmonicWalls} bias. + +\cvsubsec{Choosing a function}{sec:cvc_list} + +In this context, the function that computes a colvar is called a \emph{component}. +A component's choice and definition consists of including in the variable's configuration a keyword indicating the type of function (e.g.{} \texttt{rmsd}), followed by a definition block specifying the atoms involved (see \ref{sec:colvar_atom_groups}) and any additional parameters (cutoffs, reference'' values, \ldots). +\emph{At least one component must be chosen to define a variable:} if none of the keywords listed below is found, an error is raised. + +The following components implement functions with a scalar value (i.e.{} a real number): +\begin{itemize} +\item \refkey{distance}{colvar|distance}: distance between two groups; +\item \refkey{distanceZ}{colvar|distanceZ}: projection of a distance vector on an axis; +\item \refkey{distanceXY}{colvar|distanceXY}: projection of a distance vector on a plane; +\item \refkey{distanceInv}{colvar|distanceInv}: mean distance between two groups of atoms (e.g.~NOE-based distance); +\item \refkey{angle}{colvar|angle}: angle between three groups; +\item \refkey{dihedral}{colvar|dihedral}: torsional (dihedral) angle between four groups; +\item \refkey{dipoleAngle}{colvar|dipoleAngle}: angle between two groups and dipole of a third group; +\item \refkey{polarTheta}{colvar|polarTheta}: polar angle of a group in spherical coordinates; +\item \refkey{polarPhi}{colvar|polarPhi}: azimuthal angle of a group in spherical coordinates; +\item \refkey{coordNum}{colvar|coordNum}: coordination number between two groups; +\item \refkey{selfCoordNum}{colvar|selfCoordNum}: coordination number of atoms within a + group; +\item \refkey{hBond}{colvar|hBond}: hydrogen bond between two atoms; +\item \refkey{rmsd}{colvar|rmsd}: root mean square deviation (RMSD) from a set of + reference coordinates; +\item \refkey{eigenvector}{colvar|eigenvector}: projection of the atomic coordinates on a + vector; +\item \refkey{orientationAngle}{colvar|orientationAngle}: angle of the best-fit rotation from + a set of reference coordinates; +\item \refkey{orientationProj}{colvar|orientationProj}: cosine of \refkey{orientationProj}{colvar|orientationProj}; +\item \refkey{spinAngle}{colvar|spinAngle}: projection orthogonal to an axis of the best-fit rotation + from a set of reference coordinates; +\item \refkey{tilt}{colvar|tilt}: projection on an axis of the best-fit rotation + from a set of reference coordinates; +\item \refkey{gyration}{colvar|gyration}: radius of gyration of a group of atoms; +\item \refkey{inertia}{colvar|inertia}: moment of inertia of a group of atoms; +\item \refkey{inertiaZ}{colvar|inertiaZ}: moment of inertia of a group of atoms around a chosen axis; +\cvnamebasedonly{ +\item \refkey{alpha}{colvar|alpha}:$\alpha$-helix content of a protein segment. +\item \refkey{dihedralPC}{colvar|dihedralPC}: projection of protein backbone dihedrals onto a dihedral principal component. +} +\end{itemize} + +Some components do not return scalar, but vector values: +\begin{itemize} +\item \refkey{distanceVec}{colvar|distanceVec}: distance vector between two groups (length: 3); +\item \refkey{distanceDir}{colvar|distanceDir}: unit vector parallel to distanceVec (length: 3); +\item \refkey{cartesian}{colvar|cartesian}: vector of atomic Cartesian coordinates (length:$N$times the number of Cartesian components requested, X, Y or Z); +\item \refkey{distancePairs}{colvar|distancePairs}: vector of mutual distances (length:$N_{\mathrm{1}}\times{}N_{\mathrm{2}}$); +\item \refkey{orientation}{colvar|orientation}: best-fit rotation, expressed as a unit quaternion (length: 4). +\end{itemize} + +The types of components used in a colvar (scalar or not) determine the +properties of that colvar, and particularly which biasing or analysis methods +can be applied. + +\textbf{What if X'' is not listed?} If a function type is not available on this list, it may be possible to define it as a polynomial superposition of existing ones (see \ref{sec:cvc_superp})\cvleptononly{, a custom function (see \ref{sec:colvar_custom_function})}\cvscriptonly{, or a scripted function (see \ref{sec:colvar_scripted})}. + +In the rest of this section, all available component types are listed, along with their physical units and the ranges of values, if limited. +Such limiting values can be used to define \refkey{lowerBoundary}{colvar|lowerBoundary} and \refkey{upperBoundary}{colvar|upperBoundary} in the parent colvar. + +For each type of component, the available configurations keywords are listed: +when two components share certain keywords, the second component references to +the documentation of the first one that uses that keyword. +The very few keywords that are available for all types of components are listed in a separate section \ref{sec:cvc_common}. + +\newenvironment{cvcoptions}% + {\noindent\textbf{List of keywords} (see also \ref{sec:cvc_superp} for additional options): + \begin{itemize}} + { + \end{itemize} } +\cvsubsubsec{\texttt{distance}: center-of-mass distance between two groups.}{sec:cvc_distance} +\labelkey{colvar|distance} + +The \texttt{distance \{...\}} block defines a distance component between the two atom groups, \texttt{group1} and \texttt{group2}. + +\begin{cvcoptions} \item % - \labelkey{colvar|upperBoundary} % + \labelkey{colvar|distance|group1} \key - {upperBoundary}{% - \texttt{colvar}}{% - Upper boundary of the colvar}{% - decimal}{% - Similarly to \texttt{lowerBoundary}, defines the highest possible or allowed value.} + {group1}{% + \texttt{distance}}{% + First group of atoms}{% + Block \texttt{group1 \{...\}}}{% + First group of atoms.} \item % - \keydef - {hardLowerBoundary}{% - \texttt{colvar}}{% - Whether the lower boundary is the physical lower limit}{% - boolean}{% - \texttt{off}}{% - This option does not affect simulation results, but enables some internal optimizations. - Depending on its mathematical definition, a colvar may have natural'' boundaries: for example, a \texttt{distance} colvar has a natural'' lower boundary at 0. Setting this option instructs the Colvars module that the user-defined lower boundary is natural''. -See Section~\ref{sec:cvc} for the physical ranges of values of each component.} + \labelkey{colvar|distance|group2} + \simkey{group2}{\texttt{distance}}{group1} \item % + \labelkey{colvar|distance|forceNoPBC} \keydef - {hardUpperBoundary}{% - \texttt{colvar}}{% - Whether the upper boundary is the physical upper limit of the colvar's values}{% + {forceNoPBC}{% + \texttt{distance}}{% + Calculate absolute rather than minimum-image distance?}{% boolean}{% - \texttt{off}}{% - Analogous to \texttt{hardLowerBoundary}.} + \texttt{no}}{% + By default, in calculations with periodic boundary conditions, the + \texttt{distance} component returns the distance according to the + minimum-image convention. If this parameter is set to \texttt{yes}, + PBC will be ignored and the distance between the coordinates as maintained + internally will be used. This is only useful in a limited number of + special cases, e.g. to describe the distance between remote points + of a single macromolecule, which cannot be split across periodic cell + boundaries, and for which the minimum-image distance might give the + wrong result because of a relatively small periodic cell.} \item % - \labelkey{colvar|expandBoundaries} % + \labelkey{colvar|distance|oneSiteTotalForce} \keydef - {expandBoundaries}{% - \texttt{colvar}}{% - Allow to expand the two boundaries if needed}{% + {oneSiteTotalForce}{% + \texttt{angle}, \texttt{dipoleAngle}, \texttt{dihedral}}{% + Measure total force on group 1 only?}{% boolean}{% - \texttt{off}}{% - If defined, biasing and analysis methods may keep their own copies - of \texttt{lowerBoundary} and \texttt{upperBoundary}, and expand - them to accommodate values that do not fit in the initial range. - Currently, this option is used by the metadynamics bias - (\ref{sec:colvarbias_meta}) to keep all of its hills fully within - the grid. This option cannot be used when - the initial boundaries already span the full period of a periodic - colvar.} + \texttt{no}}{% + If this is set to \texttt{yes}, the total force is measured along + a vector field (see equation~(\ref{eq:gradient_vector}) in + section~\ref{sec:colvarbias_abf}) that only involves atoms of + \texttt{group1}. This option is only useful for ABF, or custom + biases that compute total forces. See + section~\ref{sec:colvarbias_abf} for details.} -\item % - \keydef - {subtractAppliedForce}{% - \texttt{colvar}}{% - Do not include biasing forces in the total force for this colvar}{% - boolean}{% - \texttt{off}}{% - If the colvar supports total force calculation (see \ref{sec:cvc_sys_forces}), all forces applied to this colvar by biases will be removed from the total force. - This keyword allows to recover some of the system force'' calculation available in the Colvars module before version 2016-08-10. - Please note that removal of \emph{all} other external forces (including biasing forces applied to a different colvar) is \emph{no longer supported}, due to changes in the underlying simulation engines (primarily NAMD). - This option may be useful when continuing a previous simulation where the removal of external/applied forces is essential. - \emph{For all new simulations, the use of this option is not recommended.} -} +\end{cvcoptions} -\end{itemize} +The value returned is a positive number (in \AA), ranging from$0$+to the largest possible interatomic distance within the chosen +boundary conditions (with PBCs, the minimum image convention is used +unless the \texttt{forceNoPBC} option is set). -\cvsubsec{Trajectory output}{sec:colvar_traj_output} +\cvsubsubsec{\texttt{distanceZ}: projection of a distance vector on an axis.}{sec:cvc_distanceZ} +\labelkey{colvar|distanceZ} -\begin{itemize} -\item % - \keydef - {outputValue}{% - \texttt{colvar}}{% - Output a trajectory for this colvar}{% - boolean}{% - \texttt{on}}{% - If \texttt{colvarsTrajFrequency} is non-zero, the value of this - colvar is written to the trajectory file every - \texttt{colvarsTrajFrequency} steps in the column labeled - $<$\texttt{name}$>$''.} +The \texttt{distanceZ~\{...\}} block defines a distance projection +component, which can be seen as measuring the distance between two +groups projected onto an axis, or the position of a group along such +an axis. The axis can be defined using either one reference group and +a constant vector, or dynamically based on two reference groups. +One of the groups can be set to a dummy atom to allow the use of an absolute Cartesian coordinate. +\begin{cvcoptions} \item % - \keydef - {outputVelocity}{% - \texttt{colvar}}{% - Output a velocity trajectory for this colvar}{% - boolean}{% - \texttt{off}}{% - If \texttt{colvarsTrajFrequency} is defined, the - finite-difference calculated velocity of this colvar are written - to the trajectory file under the label - \texttt{v\_}$<$\texttt{name}$>$''.} + \labelkey{colvar|distanceZ|main} + \key + {main}{% + \texttt{distanceZ}}{% + Main group of atoms}{% + Block \texttt{main \{...\}}}{% + Group of atoms whose position$\bm{r}$is measured.} \item % - \keydef - {outputEnergy}{% - \texttt{colvar}}{% - Output an energy trajectory for this colvar}{% - boolean}{% - \texttt{off}}{% - This option applies only to extended Lagrangian colvars. If - \texttt{colvarsTrajFrequency} is defined, the kinetic energy of - the extended degree and freedom and the potential energy of the - restraining spring are are written to the trajectory file under - the labels \texttt{Ek\_}$<$\texttt{name}$>$'' and - \texttt{Ep\_}$<$\texttt{name}$>$''.} + \labelkey{colvar|distanceZ|ref} + \key + {ref}{% + \texttt{distanceZ}}{% + Reference group of + atoms}{% + Block \texttt{ref \{...\}}}{% + Reference group of atoms. The position of its center of mass is + noted$\bm{r}_1$below.} \item % - \labelkey{colvar|outputTotalForce} + \labelkey{colvar|distanceZ|ref2} \keydef - {outputTotalForce}{% - \texttt{colvar}}{% - Output a total force trajectory for this - colvar}{% - boolean}{% - \texttt{off}}{% - If \texttt{colvarsTrajFrequency} is defined, the total force on this - colvar (i.e.~the projection of all atomic total forces - onto this colvar --- see - equation~(\ref{eq:gradient_vector}) in - section~\ref{sec:colvarbias_abf}) are written to the trajectory - file under the label \texttt{fs\_}$<$\texttt{name}$>$''. - For extended Lagrangian colvars, the total force'' felt by the extended degree of freedom - is simply the force from the harmonic spring. - \textbf{Note:} not all components support this option. The - physical unit for this force is \cvnamdonly{kcal/mol}\cvvmdonly{kcal/mol}\cvlammpsonly{the unit of energy specified by \texttt{units}}, divided by the colvar unit U.} + {ref2}{% + \texttt{distanceZ}}{% + Secondary reference + group}{% + Block \texttt{ref2 \{...\}}}{% + none}{% + Optional group of reference atoms, whose position$\bm{r}_2$can + be used to define a dynamic projection axis:$\bm{e}=(\| \bm{r}_2
+    - \bm{r}_1\|)^{-1} \times (\bm{r}_2 - \bm{r}_1)$. In this case, + the origin is$\bm{r}_m = 1/2 (\bm{r}_1+\bm{r}_2)$, and the value + of the component is$\bm{e} \cdot (\bm{r}-\bm{r}_m)$.} \item % + \labelkey{colvar|distanceZ|axis} \keydef - {outputAppliedForce}{% - \texttt{colvar}}{% - Output an applied force trajectory for this - colvar}{% - boolean}{% - \texttt{off}}{% - If \texttt{colvarsTrajFrequency} is defined, the total force - applied on this colvar by Colvars biases are - written to the trajectory under the label - \texttt{fa\_}$<$\texttt{name}$>$''. - For extended Lagrangian colvars, this force is actually applied to the - extended degree of freedom rather than the geometric colvar itself. - The physical unit for this - force is \cvnamdonly{kcal/mol}\cvvmdonly{kcal/mol}\cvlammpsonly{the unit of energy specified by \texttt{units}} divided by the colvar unit.} + {axis}{% + \texttt{distanceZ}}{% + Projection axis (\AA{})}{% + \texttt{(x, y, z)} triplet}{% + \texttt{(0.0, 0.0, 1.0)}}{% + The three components of this vector define a + projection axis$\bm{e}$for the distance vector$\bm{r} -
+    \bm{r}_1$joining the centers of groups \texttt{ref} and + \texttt{main}. The value of the component is then$\bm{e} \cdot
+    (\bm{r}-\bm{r}_1)$. The vector should be written as three + components separated by commas and enclosed in parentheses.} -\end{itemize} +\item % + \dupkey{forceNoPBC}{\texttt{distanceZ}}{colvar|distance|forceNoPBC}{\texttt{distance} component} +\item % + \dupkey{oneSiteTotalForce}{\texttt{distanceZ}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} +\end{cvcoptions} +This component returns a number (in \AA{}) whose range is determined +by the chosen boundary conditions. For instance, if the$z$axis is +used in a simulation with periodic boundaries, the returned value ranges +between$-b_{z}/2$and$b_{z}/2$, where$b_{z}$is the box length +along$z$(this behavior is disabled if \texttt{forceNoPBC} is set). -\cvsubsec{Extended Lagrangian.}{sec:colvar_extended} -The following options enable extended-system -dynamics, where a colvar is coupled to an additional degree of freedom -(fictitious particle) by a harmonic spring. -All biasing and confining forces are then applied to the extended degree -of freedom. The actual'' geometric colvar (function of Cartesian -coordinates) only feels the force from the harmonic spring. -This is particularly useful when combined with an ABF bias (\ref{sec:colvarbias_abf}) -to perform eABF simulations (\ref{sec:eABF}). +\cvsubsubsec{\texttt{distanceXY}: modulus of the projection of a distance vector on a plane.}{sec:cvc_distanceXY} +\labelkey{colvar|distanceXY} -\begin{itemize} -\item % - \keydef - {extendedLagrangian}{% - \texttt{colvar}}{% - Add extended degree of freedom}{% - boolean}{% - \texttt{off}}{% - Adds a fictitious particle to be coupled to the colvar by a harmonic - spring. The fictitious mass and the force constant of the coupling - potential are derived from the parameters \texttt{extendedTimeConstant} - and \texttt{extendedFluctuation}, described below. Biasing forces on the - colvar are applied to this fictitious particle, rather than to the - atoms directly. This implements the extended Lagrangian formalism - used in some metadynamics simulations~\cite{Iannuzzi2003}. - \cvnamdonly{The energy associated with the extended degree of freedom is reported - under the MISC title in NAMD's energy output.} - } +The \texttt{distanceXY~\{...\}} block defines a distance projected on +a plane, and accepts the same keywords as the component \texttt{distanceZ}, i.e. +\texttt{main}, \texttt{ref}, either \texttt{ref2} or \texttt{axis}, +and \texttt{oneSiteTotalForce}. It returns the norm of the +projection of the distance vector between \texttt{main} and +\texttt{ref} onto the plane orthogonal to the axis. The axis is +defined using the \texttt{axis} parameter or as the vector joining +\texttt{ref} and \texttt{ref2} (see \texttt{distanceZ} above). +\begin{cvcoptions} \item % - \key - {extendedFluctuation}{% - \texttt{colvar}}{% - Standard deviation between the colvar and the fictitious - particle (colvar unit)}{% - positive decimal}{% - Defines the spring stiffness for the \texttt{extendedLagrangian} - mode, by setting the typical deviation between the colvar and the extended - degree of freedom due to thermal fluctuation. - The spring force constant is calculated internally as$k_B T / \sigma^2$, - where$\sigma$is the value of \texttt{extendedFluctuation}.} - + \dupkey{main}{\texttt{distanceXY}}{colvar|distanceZ|main}{\texttt{distanceZ} component} \item % - \keydef - {extendedTimeConstant}{% - \texttt{colvar}}{% - Oscillation period of the fictitious particle (fs)}{% - positive decimal}{% - \texttt{200}}{% - Defines the inertial mass of the fictitious particle, by setting the - oscillation period of the harmonic oscillator formed by the fictitious - particle and the spring. The period - should be much larger than the MD time step to ensure accurate integration - of the extended particle's equation of motion. - The fictitious mass is calculated internally as$k_B T (\tau/2 \pi \sigma)^2$, - where$\tau$is the period and$\sigma$is the typical fluctuation (see above).} - + \dupkey{ref}{\texttt{distanceXY}}{colvar|distanceZ|ref}{\texttt{distanceZ} component} \item % - \keydef - {extendedTemp}{% - \texttt{colvar}}{% - Temperature for the extended degree of freedom (K)}{% - positive decimal}{% - thermostat temperature}{% - Temperature used for calculating the coupling force constant of the - extended variable (see \texttt{extendedFluctuation}) and, if needed, as a - target temperature for extended Langevin dynamics (see - \texttt{extendedLangevinDamping}). This should normally be left at its - default value.} - + \dupkey{ref2}{\texttt{distanceXY}}{colvar|distanceZ|ref2}{\texttt{distanceZ} component} \item % - \keydef - {extendedLangevinDamping}{% - \texttt{colvar}}{% - Damping factor for extended Langevin dynamics - (ps$^{-1}$)}{% - positive decimal}{% - \texttt{1.0}}{% - If this is non-zero, the extended degree of freedom undergoes Langevin dynamics - at temperature \texttt{extendedTemp}. The friction force is minus - \texttt{extendedLangevinDamping} times the velocity. This is useful because - the extended dynamics coordinate may heat up in the transient - non-equilibrium regime of ABF. Use moderate damping values, to limit - viscous friction (potentially slowing down diffusive sampling) and stochastic - noise (increasing the variance of statistical measurements). In - doubt, use the default value.} -\end{itemize} + \dupkey{axis}{\texttt{distanceXY}}{colvar|distanceZ|axis}{\texttt{distanceZ} component} +\item % + \dupkey{forceNoPBC}{\texttt{distanceXY}}{colvar|distance|forceNoPBC}{\texttt{distance} component} +\item % + \dupkey{oneSiteTotalForce}{\texttt{distanceZ}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} +\end{cvcoptions} -\cvsubsec{Statistical analysis of collective variables}{sec:colvar_acf} -When the global keyword \texttt{analysis} is defined in the -configuration file, run-time calculations of statistical properties for -individual colvars can be performed. At the moment, several types of -time correlation functions, running averages and running standard -deviations are available. +\cvsubsubsec{\texttt{distanceVec}: distance vector between two groups.}{sec:cvc_distanceVec} +\labelkey{colvar|distanceVec} -\begin{itemize} +The \texttt{distanceVec~\{...\}} block defines +a distance vector component, which accepts the same keywords as +the component \texttt{distance}: \texttt{group1}, \texttt{group2}, and +\texttt{forceNoPBC}. Its value is the 3-vector joining the centers +of mass of \texttt{group1} and \texttt{group2}. +\begin{cvcoptions} \item % - \keydef - {corrFunc}{% - \texttt{colvar}}{% - Calculate a time correlation function?}{% - boolean}{% - \texttt{off}}{% - Whether or not a time correlaction function should be calculated - for this colvar.} - + \dupkey{group1}{\texttt{distanceVec}}{colvar|distance|group1}{\texttt{distance} component} \item % - \key - {corrFuncWithColvar}{% - \texttt{colvar}}{% - Colvar name for the correlation function}{% - string}{% - By default, the auto-correlation function (ACF) of this colvar, -$\xi_{i}$, is calculated. When this option is specified, the - correlation function is calculated instead with another colvar, -$\xi_{j}$, which must be of the same type (scalar, vector, or - quaternion) as$\xi_{i}$.} + \simkey{group2}{\texttt{distanceVec}}{group1} +\item % + \dupkey{forceNoPBC}{\texttt{distanceVec}}{colvar|distance|forceNoPBC}{\texttt{distance} component} +\item % + \dupkey{oneSiteTotalForce}{\texttt{distanceVec}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} +\end{cvcoptions} -\item% - \keydef - {corrFuncType}{% - \texttt{colvar}}{% - Type of the correlation function}{% - \texttt{velocity}, \texttt{coordinate} or - \texttt{coordinate\_p2}}{% - \texttt{velocity}}{% - With \texttt{coordinate} or \texttt{velocity}, the correlation - function$C_{i,j}(t)$~=$\left\langle \Pi\left(\xi_{i}(t_{0}),
-        \xi_{j}(t_{0}+t)\right) \right\rangle$is calculated between - the variables$\xi_{i}$and$\xi_{j}$, or their velocities. -$\Pi(\xi_{i}, \xi_{j})$is the scalar product when calculated - between scalar or vector values, whereas for quaternions it is the - cosine between the two corresponding rotation axes. With - \texttt{coordinate\_p2}, the second order Legendre polynomial, -$(3\cos(\theta)^{2}-1)/2$, is used instead of the cosine.} -\item % - \keydef - {corrFuncNormalize}{% - \texttt{colvar}}{% - Normalize the time correlation function?}{% - boolean}{% - \texttt{on}}{% - If enabled, the value of the correlation function at$t$~= 0 - is normalized to 1; otherwise, it equals to$\left\langle
-      O\left(\xi_{i}, \xi_{j}\right) \right\rangle$.} -\item % - \keydef - {corrFuncLength}{% - \texttt{colvar}}{% - Length of the time correlation function}{% - positive integer}{% - \texttt{1000}}{% - Length (in number of points) of the time correlation function.} +\cvsubsubsec{\texttt{distanceDir}: distance unit vector between two groups.}{sec:cvc_distanceDir} +\labelkey{colvar|distanceDir} -\item % - \keydef - {corrFuncStride}{% - \texttt{colvar}}{% - Stride of the time correlation function}{% - positive integer}{% - \texttt{1}}{% - Number of steps between two values of the time correlation function.} +The \texttt{distanceDir~\{...\}} block defines +a distance unit vector component, which accepts the same keywords as +the component \texttt{distance}: \texttt{group1}, \texttt{group2}, and +\texttt{forceNoPBC}. It returns a +3-dimensional unit vector$\mathbf{d} = (d_{x}, d_{y}, d_{z})$, with +$|\mathbf{d}| = 1$. +\begin{cvcoptions} \item % - \keydef - {corrFuncOffset}{% - \texttt{colvar}}{% - Offset of the time correlation function}{% - positive integer}{% - \texttt{0}}{% - The starting time (in number of steps) of the time correlation - function (default:$t$~= 0). \textbf{Note:} \emph{the value at$t$~= 0 is always - used for the normalization}.} - + \dupkey{group1}{\texttt{distanceDir}}{colvar|distance|group1}{\texttt{distance} component} \item % - \keydef - {corrFuncOutputFile}{% - \texttt{colvar}}{% - Output file for the time correlation function}{% - UNIX filename}{% - \texttt{$<$name$>$.corrfunc.dat}}{% - The time correlation function is saved in this file.} - + \simkey{group2}{\texttt{distanceDir}}{group1} \item % - \keydef - {runAve}{% - \texttt{colvar}}{% - Calculate the running average and standard deviation}{% - boolean}{% - \texttt{off}}{% - Whether or not the running average and standard deviation should - be calculated for this colvar.} - + \dupkey{forceNoPBC}{\texttt{distanceDir}}{colvar|distance|forceNoPBC}{\texttt{distance} component} \item % - \keydef - {runAveLength}{% - \texttt{colvar}}{% - Length of the running average window}{% - positive integer}{% - \texttt{1000}}{% - Length (in number of points) of the running average window.} + \dupkey{oneSiteTotalForce}{\texttt{distanceDir}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} +\end{cvcoptions} -\item % - \keydef - {runAveStride}{% - \texttt{colvar}}{% - Stride of the running average window values}{% - positive integer}{% - \texttt{1}}{% - Number of steps between two values within the running average window.} -\item % - \keydef - {runAveOutputFile}{% - \texttt{colvar}}{% - Output file for the running average and standard deviation}{% - UNIX filename}{% - \texttt{$<$name$>$.runave.dat}}{% - The running average and standard deviation are saved in this file.} +\cvsubsubsec{\texttt{distanceInv}: mean distance between two groups of atoms.}{sec:cvc_distanceInv} +\labelkey{colvar|distanceInv} -\end{itemize} +The \texttt{distanceInv~\{...\}} block defines a generalized mean distance between two groups of atoms 1 and 2, weighted with exponent$1/n$: + + \label{eq:distanceInv} + d_{\mathrm{1,2}}^{[n]} \; = \; \left(\frac{1}{N_{\mathrm{1}}N_{\mathrm{2}}}\sum_{i,j} \left(\frac{1}{\Vert\mathbf{d}^{ij}\Vert}\right)^{n} \right)^{-1/n} + +where$\Vert\mathbf{d}^{ij}\Vert$is the distance between atoms$i$and$j$in groups 1 and 2 respectively, and$n$is an even integer. +\begin{cvcoptions} +\item % + \dupkey{group1}{\texttt{distanceInv}}{colvar|distance|group1}{\texttt{distance} component} +\item % + \simkey{group2}{\texttt{distanceInv}}{group1} +\item % + \dupkey{oneSiteTotalForce}{\texttt{distanceInv}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} +\item % + \keydef + {exponent}{% + \texttt{distanceInv}}{% + Exponent$n$in equation~\ref{eq:distanceInv}}{% + positive even integer}{% + 6}{Defines the exponent to which the individual distances are elevated before averaging. The default value of 6 is useful for example to applying restraints based on NOE-measured distances.} +\end{cvcoptions} +This component returns a number in \AA{}, ranging from$0$to the largest possible distance within the chosen boundary conditions. -\cvsec{Selecting atoms}{sec:colvar_atom_groups} -To define collective variables, atoms are usually selected as groups. Each group is defined using an identifier that is unique in the context of the specific colvar component (e.g.~for a distance component, the two groups are \texttt{group1} and \texttt{group2}). -The identifier is followed by a brace-delimited block containing selection keywords and other parameters, including an optional \texttt{name}: +\cvsubsubsec{\texttt{distancePairs}: set of pairwise distances between two groups.}{sec:cvc_distancePairs} +\labelkey{colvar|distancePairs} -\begin{itemize} -\item \key - {name}{% - atom group}{% - Unique name for the atom group}{% - string}{% - This parameter defines a unique name for this atom group, which can be referred to - in the definition of other atom groups (including in other colvars) by invoking - \texttt{atomsOfGroup} as a selection keyword.} -\end{itemize} +The \texttt{distancePairs~\{...\}} block defines a$N_{\mathrm{1}}\times{}N_{\mathrm{2}}$-dimensional variable that includes all mutual distances between the atoms of two groups. +This can be useful, for example, to develop a new variable defined over two groups, by using the \texttt{scriptedFunction} feature. +\begin{cvcoptions} +\item % + \dupkey{group1}{\texttt{distancePairs}}{colvar|distance|group1}{\texttt{distance} component} +\item % + \simkey{group2}{\texttt{distancePairs}}{group1} +\item % + \dupkey{forceNoPBC}{\texttt{distancePairs}}{colvar|distance|forceNoPBC}{\texttt{distance} component} +\end{cvcoptions} +This component returns a$N_{\mathrm{1}}\times{}N_{\mathrm{2}}$-dimensional vector of numbers, each ranging from$0$to the largest possible distance within the chosen boundary conditions. -\cvsubsec{Atom selection keywords}{sec:colvar_atom_groups_sel} -Selection keywords may be used individually or in combination with each other, and each can be repeated any number of times. -Selection is incremental: each keyword adds the corresponding atoms to the selection, so that different sets of atoms can be combined. -However, atoms included by multiple keywords are only counted once. -Below is an example configuration for an atom group called \texttt{atoms}''. -\textbf{Note: }\emph{this is an unusually varied combination of selection keywords, demonstrating how they can be combined together: most simulations only use one of them.}\\ +\cvsubsubsec{\texttt{cartesian}: vector of atomic Cartesian coordinates.}{sec:cvc_cartesian} +\labelkey{colvar|cartesian} -{% -% verbatim can't appear within commands -\noindent\ttfamily atoms \{\\ -\\ -\-~~\# add atoms 1 and 3 to this group (note: the first atom in the system is 1)\\ -\-~~atomNumbers \{ \\ -\-~~~~1 3\\ -\-~~\}\\ -\\ -\-~~\# add atoms starting from 20 up to and including 50\\ -\-~~atomNumbersRange 20-50\\ -} -\cvnamebasedonly{{% -\noindent\ttfamily\\ -\-~~\# add all the atoms with occupancy 2 in the file atoms.pdb\\ -\-~~atomsFile atoms.pdb\\ -\-~~atomsCol O\\ -\-~~atomsColValue 2.0\\ -\\ -\-~~\# add all the C-alphas within residues 11 to 20 of segments "PR1" and "PR2"\\ -\-~~psfSegID PR1 PR2\\ -\-~~atomNameResidueRange CA 11-20\\ -\-~~atomNameResidueRange CA 11-20\\ -}} -{\noindent\ttfamily\\ -\-~~\# add index group (requires a .ndx file to be provided globally)\\ -\-~~indexGroup Water\\ -\}\\} +The \texttt{cartesian~\{...\}} block defines a component returning a flat vector containing +the Cartesian coordinates of all participating atoms, in the order +$(x_1, y_1, z_1, \cdots, x_n, y_n, z_n)$. +\begin{cvcoptions} +\item % + \key + {atoms}{% + \texttt{cartesian}}{% + Group of atoms}{% + Block \texttt{atoms \{...\}}}{% + Defines the atoms whose coordinates make up the value of the component. + If \texttt{rotateReference} or \texttt{centerReference} are defined, coordinates + are evaluated within the moving frame of reference.} +\end{cvcoptions} -The resulting selection includes atoms 1 and 3, those between 20 and 50,\cvnamebasedonly{ the$\mathrm{C}_{\alpha}$atoms between residues 11 and 20 of the two segments \texttt{PR1} and \texttt{PR2},} and those in the index group called Water''. -The indices of this group are read from the file provided by the global keyword \refkey{indexFile}{Colvars-global|indexFile}. -\cvvmdonly{In the current version, the Colvars module does not manipulate VMD atom selections directly: however, these can be converted to atom groups within the Colvars configuration string, using selection keywords such as \texttt{atomNumbers}.} -The complete list of selection keywords available in \MDENGINE{} is: +\cvsubsubsec{\texttt{angle}: angle between three groups.}{sec:cvc_angle} +\labelkey{colvar|angle} -\begin{itemize} +The \texttt{angle~\{...\}} block defines an angle, and contains the +three blocks \texttt{group1}, \texttt{group2} and \texttt{group3}, defining +the three groups. It returns an angle (in degrees) within the +interval$[0:180]$. +\begin{cvcoptions} \item % - \key - {atomNumbers}{% - atom group}{% - List of atom numbers}{% - space-separated list of positive integers}{% - This option adds to the group all the atoms whose numbers are in - the list. \emph{The number of the first atom in the system is 1: to convert from a VMD selection, use atomselect get serial''.} - } - + \dupkey{group1}{\texttt{angle}}{colvar|distance|group1}{\texttt{distance} component} \item % - \key - {indexGroup}{% - atom group}{% - Name of index group to be used (GROMACS format)}{% - string}{% - If the name of an index file has been provided by \texttt{indexFile}, this option allows to select one index group from that file: the atoms from that index group will be used to define the current group.} - + \simkey{group2}{\texttt{angle}}{group1} \item % - \key - {atomsOfGroup}{% - atom group}{% - Name of group defined previously}{% - string}{% - Refers to a group defined previously using its user-defined \texttt{name}. - This adds all atoms of that named group to the current group.} - + \simkey{group3}{\texttt{angle}}{group1} \item % - \key - {atomNumbersRange}{% - atom group}{% - Atoms within a number range}{% -$<$Starting number$>$-$<$Ending number$>$}{% - This option includes in the group all atoms whose numbers are within the range specified. \emph{The number of the first atom in the system is 1.} - } - -\cvnamebasedonly{ + \dupkey{forceNoPBC}{\texttt{angle}}{colvar|distance|forceNoPBC}{\texttt{distance} component} \item % - \key - {atomNameResidueRange}{% - atom group}{% - Named atoms within a range of residue numbers}{% -$<$Atom name$><$Starting residue$>$-$<$Ending residue$>$}{% - This option adds to the group all the atoms with the provided - name, within residues in the given range.} + \dupkey{oneSiteTotalForce}{\texttt{angle}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} +\end{cvcoptions} -\item % - \key - {psfSegID}{% - atom group}{% - PSF segment identifier}{% - space-separated list of strings (max 4 characters)}{% - This option sets the PSF segment identifier for - \texttt{atomNameResidueRange}. Multiple values may be provided, - which correspond to multiple instances of - \texttt{atomNameResidueRange}, in order of their occurrence. - This option is only necessary if a PSF topology file is used.} -\item % - \key - {atomsFile}{% - atom group}{% - PDB file name for atom selection}{% - UNIX filename}{% - This option selects atoms from the PDB file provided and adds them - to the group according to numerical flags in the column - \texttt{atomsCol}. \textbf{Note:} \emph{the sequence of atoms in the PDB file - provided must match that in the system's topology}.} -\item % - \key - {atomsCol}{% - atom group}{% - PDB column to use for atom selection flags}{% - \texttt{O}, \texttt{B}, \texttt{X}, \texttt{Y}, or \texttt{Z}}{% - This option specifies which PDB column in \texttt{atomsFile} is used to determine which atoms are to be included in the group. - } +\cvsubsubsec{\texttt{dipoleAngle}: angle between two groups and dipole of a third group.}{sec:cvc_dipoleAngle} +\labelkey{colvar|dipoleAngle} -\item % - \key - {atomsColValue}{% - atom group}{% - Atom selection flag in the PDB column}{% - positive decimal}{% - If defined, this value in \texttt{atomsCol} identifies atoms in \texttt{atomsFile} that are included in the group. - If undefined, all atoms with a non-zero value in \texttt{atomsCol} are included.} -} +The \texttt{dipoleAngle~\{...\}} block defines an angle, and contains the +three blocks \texttt{group1}, \texttt{group2} and \texttt{group3}, defining +the three groups, being \texttt{group1} the group where dipole is calculated. +It returns an angle (in degrees) within the interval$[0:180]$. +\begin{cvcoptions} \item % - \key - {dummyAtom}{% - atom group}{% - Dummy atom position (\AA{})}{% - \texttt{(x, y, z)} triplet}{% - Instead of selecting any atom, this option makes the group a virtual particle at a fixed position in space. This is useful e.g.~to replace a group's center of geometry with a user-defined position.} + \dupkey{group1}{\texttt{dipoleAngle}}{colvar|distance|group1}{\texttt{distance} component} +\item % + \simkey{group2}{\texttt{dipoleAngle}}{group1} +\item % + \simkey{group3}{\texttt{dipoleAngle}}{group1} +\item % + \dupkey{forceNoPBC}{\texttt{dipoleAngle}}{colvar|distance|forceNoPBC}{\texttt{distance} component} +\item % + \dupkey{oneSiteTotalForce}{\texttt{dipoleAngle}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} +\end{cvcoptions} -\end{itemize} -\cvsubsec{Moving frame of reference.}{sec:colvar_atom_groups_ref_frame} +\cvsubsubsec{\texttt{dihedral}: torsional angle between four groups.}{sec:cvc_dihedral} +\labelkey{colvar|dihedral} -The following options define an automatic calculation of an optimal translation (\texttt{centerReference}) or optimal rotation (\texttt{rotateReference}), that superimposes the positions of this group to a provided set of reference coordinates. -This can allow, for example, to effectively remove from certain colvars the effects of molecular tumbling and of diffusion. -Given the set of atomic positions$\mathbf{x}_{i}$, the colvar$\xi$can be defined on a set of roto-translated positions$\mathbf{x}_{i}' = R(\mathbf{x}_{i} - \mathbf{x}^{\mathrm{C}}) + \mathbf{x}^{\mathrm{ref}}$. -$\mathbf{x}^{\mathrm{C}}$is the geometric center of the$\mathbf{x}_{i}$,$R$is the optimal rotation matrix to the reference positions and$\mathbf{x}^{\mathrm{ref}}$is the geometric center of the reference positions. +The \texttt{dihedral~\{...\}} block defines a torsional angle, and +contains the blocks \texttt{group1}, \texttt{group2}, \texttt{group3} +and \texttt{group4}, defining the four groups. It returns an angle +(in degrees) within the interval$[-180:180]$. The Colvars module +calculates all the distances between two angles taking into account +periodicity. For instance, reference values for restraints or range +boundaries can be defined by using any real number of choice. -Components that are defined based on pairwise distances are naturally invariant under global roto-translations. -Other components are instead affected by global rotations or translations: however, they can be made invariant if they are expressed in the frame of reference of a chosen group of atoms, using the \texttt{centerReference} and \texttt{rotateReference} options. -Finally, a few components are defined by convention using a roto-translated frame (e.g. the minimal RMSD): for these components, \texttt{centerReference} and \texttt{rotateReference} are enabled by default. -In typical applications, the default settings result in the expected behavior. +\begin{cvcoptions} +\item % + \dupkey{group1}{\texttt{dihedral}}{colvar|distance|group1}{\texttt{distance} component} +\item % + \simkey{group2}{\texttt{dihedral}}{group1} +\item % + \simkey{group3}{\texttt{dihedral}}{group1} +\item % + \simkey{group4}{\texttt{dihedral}}{group1} +\item % + \dupkey{forceNoPBC}{\texttt{dihedral}}{colvar|distance|forceNoPBC}{\texttt{distance} component} +\item % + \dupkey{oneSiteTotalForce}{\texttt{dihedral}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} +\end{cvcoptions} -\paragraph*{Warning on rotating frames of reference and periodic boundary conditions.} -\texttt{rotateReference} affects coordinates that depend on minimum-image distances in periodic boundary conditions (PBC). -After rotation of the coordinates, the periodic cell vectors become irrelevant: the rotated system is effectively non-periodic. -A safe way to handle this is to ensure that the relevant inter-group distance vectors remain smaller than the half-size of the periodic cell. -If this is not desirable, one should avoid the rotating frame of reference, and apply orientational restraints to the reference group instead, in order to keep the orientation of the reference group consistent with the orientation of the periodic cell. -\paragraph*{Warning on rotating frames of reference and ABF.} -Note that \texttt{centerReference} and \texttt{rotateReference} may affect the Jacobian derivative of colvar components in a way that is not taken into account by default. -Be careful when using these options in ABF simulations or when using total force values. +\cvsubsubsec{\texttt{polarTheta}: polar angle in spherical coordinates.}{sec:cvc_polarTheta} +\labelkey{colvar|polarTheta} -\begin{itemize} +The \texttt{polarTheta~\{...\}} block defines the polar angle in +spherical coordinates, for the center of mass of a group of atoms +described by the block \texttt{atoms}. It returns an angle +(in degrees) within the interval$[0:180]. +To obtain spherical coordinates in a frame of reference tied to +another group of atoms, use the \texttt{fittingGroup} (\ref{sec:colvar_atom_groups_ref_frame}) option +within the \texttt{atoms} block. +An example is provided in file \texttt{examples/11\_polar\_angles.in} of the Colvars public repository. +\begin{cvcoptions} \item % - \keydef - {centerReference}{% - atom group}{% - Implicitly remove translations for this group}{% - boolean}{% - \texttt{off}}{% - If this option is \texttt{on}, the center of geometry of the group will be aligned with that of the reference positions provided by \cvnamebasedonly{either} \texttt{refPositions} or \texttt{refPositionsFile}. - Colvar components will only have access to the aligned positions. -\textbf{Note}: unless otherwise specified, \texttt{rmsd} and \texttt{eigenvector} set this option to \texttt{on} \emph{by default}. -} + \labelkey{colvar|polarTheta|atoms} + \key + {atoms}{% + \texttt{polarPhi}}{% + Atom group}{% + \texttt{atoms~\{...\}} block}{% + Defines the group of atoms for the COM of which the angle should be calculated. + } +\end{cvcoptions} -\item % - \keydef - {rotateReference}{% - atom group}{% - Implicitly remove rotations for this group}{% - boolean}{% - \texttt{off}}{% - If this option is \texttt{on}, the coordinates of this group will be optimally superimposed to the reference positions provided by \cvnamebasedonly{either} \texttt{refPositions} or \texttt{refPositionsFile}. - The rotation will be performed around the center of geometry if \texttt{centerReference} is \texttt{on}, or around the origin otherwise. - The algorithm used is the same employed by the \texttt{orientation} colvar component~\cite{Coutsias2004}. - Forces applied to the atoms of this group will also be implicitly rotated back to the original frame. - \textbf{Note}: unless otherwise specified, \texttt{rmsd} and \texttt{eigenvector} set this option to \texttt{on} \emph{by default}. -} -\item % - \labelkey{atom-group|refPositions} - \key - {refPositions}{% - atom group}{% - Reference positions for fitting (\AA)}{% - space-separated list of \texttt{(x, y, z)} triplets}{% - \label{key:colvars:atom_group:refPositions} - This option provides a list of reference coordinates for \texttt{centerReference} and/or \texttt{rotateReference}, and is mutually exclusive with \texttt{refPositionsFile}. - If only \texttt{centerReference} is \texttt{on}, the list may contain a single (x, y, z) triplet; if also \texttt{rotateReference} is \texttt{on}, the list should be as long as the atom group, and \emph{its order must match the order in which atoms were defined}. -} +\cvsubsubsec{\texttt{polarPhi}: azimuthal angle in spherical coordinates.}{sec:cvc_polarPhi} +\labelkey{colvar|polarPhi} -\item % - \labelkey{atom-group|refPositionsFile} - \key - {refPositionsFile}{% - atom group}{% - File containing the reference positions for fitting}{% - UNIX filename}{% - \label{key:colvars:atom_group:refPositionsFile} - This option provides a list of reference coordinates for \texttt{centerReference} and/or \texttt{rotateReference}, and is mutually exclusive with \texttt{refPositions}. - The acceptable file format is XYZ, which is read in double precision\cvnamebasedonly{, or PDB; \emph{the latter is discouraged if the precision of the reference coordinates is a concern}}. - Atomic positions are read differently depending on the following scenarios: - \textbf{(i)} the file contains exactly as many records as the atoms in the group: all positions are read in sequence; - \textbf{(ii)} (most common case) the file contains coordinates for the entire system: only the positions corresponding to the numeric indices of the atom group are read\cvnamebasedonly{; - \textbf{(iii)} if the file is a PDB file and \texttt{refPositionsCol} is specified, positions are read according to the value of the column \texttt{refPositionsCol} (which may be the same as \texttt{atomsCol})}. - In each case, atoms are read from the file \emph{in order of increasing number}. -} +The \texttt{polarPhi~\{...\}} block defines the azimuthal angle in +spherical coordinates, for the center of mass of a group of atoms +described by the block \texttt{atoms}. It returns an angle +(in degrees) within the interval[-180:180]$. The Colvars module +calculates all the distances between two angles taking into account +periodicity. For instance, reference values for restraints or range +boundaries can be defined by using any real number of choice. +To obtain spherical coordinates in a frame of reference tied to +another group of atoms, use the \texttt{fittingGroup} (\ref{sec:colvar_atom_groups_ref_frame}) option +within the \texttt{atoms} block. +An example is provided in file \texttt{examples/11\_polar\_angles.in} of the Colvars public repository. -\cvnamebasedonly{ +\begin{cvcoptions} \item % + \labelkey{colvar|polarPhi|atoms} \key - {refPositionsCol}{% - atom group}{% - PDB column containing atom flags}{% - \texttt{O}, \texttt{B}, \texttt{X}, \texttt{Y}, or \texttt{Z}}{% - Like \texttt{atomsCol} for \texttt{atomsFile}, indicates which column to use to identify the atoms in \texttt{refPositionsFile} (if this is a PDB file).} + {atoms}{% + \texttt{polarPhi}}{% + Atom group}{% + \texttt{atoms~\{...\}} block}{% + Defines the group of atoms for the COM of which the angle should be calculated. + } +\end{cvcoptions} -\item % - \key - {refPositionsColValue}{% - atom group}{% - Atom selection flag in the PDB column}{% - positive decimal}{% - Analogous to \texttt{atomsColValue}, but applied to \texttt{refPositionsCol}.} -} -\item % - \labelkey{atom-group|fittingGroup} % - \keydef - {fittingGroup}{% - atom group}{% - Use an alternate set of atoms to define the roto-translation}{% - Block \texttt{fittingGroup \{ ... \}}}{% - This group itself}{% - If either \texttt{centerReference} or \texttt{rotateReference} is defined, this keyword defines an alternate atom group to calculate the optimal roto-translation. - Use this option to define a continuous rotation if the structure of the group involved changes significantly (a typical symptom would be the message Warning: discontinuous rotation!''). +\cvsubsubsec{\texttt{coordNum}: coordination number between two groups.}{sec:cvc_coordNum} +\labelkey{colvar|coordNum} -\cvnamebasedonly{ - The following example illustrates the syntax of \texttt{fittingGroup}: a group called \texttt{atoms}'' is defined, including 8 C$_{\alpha}$atoms of a protein of 100 residues. - An optimal roto-translation is calculated automatically by fitting the C$_{\alpha}$trace of the rest of the protein onto the coordinates provided by a PDB file.} +The \texttt{coordNum \{...\}} block defines +a coordination number (or number of contacts), which calculates the +function$(1-(d/d_0)^{n})/(1-(d/d_0)^{m})$, where$d_0$is the +cutoff'' distance, and$n$and$m$are exponents that can control +its long range behavior and stiffness \cite{Iannuzzi2003}. This +function is summed over all pairs of atoms in \texttt{group1} and +\texttt{group2}: + + \label{eq:cvc_coordNum} + C (\mathtt{group1}, \mathtt{group2}) \; = \; + \sum_{i\in\mathtt{group1}}\sum_{j\in\mathtt{group2}} { + \frac{1 - (|\mathbf{x}_{i}-\mathbf{x}_{j}|/d_{0})^{n}}{ + 1 - (|\mathbf{x}_{i}-\mathbf{x}_{j}|/d_{0})^{m} } + } + -{% -\noindent\ttfamily -\# Example: defining a group "atoms", with its coordinates expressed \\ -\# on a roto-translated frame of reference defined by a second group\\ -atoms \{\\ -\\ -\-~~psfSegID PROT\\ -\-~~atomNameResidueRange CA 41-48\\ -\\ -\-~~centerReference yes\\ -\-~~rotateReference yes\\ -\-~~fittingGroup \{\\ -\-~~~~\# define the frame by fitting the rest of the protein\\ -\-~~~~psfSegID PROT PROT\\ -\-~~~~atomNameResidueRange CA 1-40\\ -\-~~~~atomNameResidueRange CA 49-100\\ -\-~~\} \\ -\-~~refPositionsFile all.pdb \# can be the entire system\\ -\}\\} -} -\end{itemize} +\begin{cvcoptions} -The following two options have default values appropriate for the vast majority of applications, and are only provided to support rare, special cases. -\begin{itemize} +\item % + \labelkey{colvar|coordNum|group1} + \dupkey{group1}{\texttt{coordNum}}{colvar|distance|group1}{\texttt{distance} component} \item % - \keydef - {enableFitGradients}{% - atom group}{% - Include the roto-translational contribution to colvar gradients}{% - boolean}{% - \texttt{on}}{% - When either \texttt{centerReference} or \texttt{rotateReference} is on, - the gradients of some colvars include terms proportional to -$\partial{}R/\partial\mathbf{x}_{i}$(rotational gradients) and -$\partial\mathbf{x}^{\mathrm{C}}/\partial\mathbf{x}_{i}$(translational gradients). - By default, these terms are calculated and included in the total gradients; - if this option is set to \texttt{off}, they are neglected. - In the case of a minimum RMSD component, this flag is automatically disabled - because the contributions of those derivatives to the gradients cancel out. -} + \labelkey{colvar|coordNum|group2} + \simkey{group2}{\texttt{coordNum}}{group1} \item % + \labelkey{colvar|coordNum|cutoff} \keydef - {enableForces}{% - atom group}{% - Apply forces from this colvar to this group}{% - boolean}{% - \texttt{on}}{% - If this option is \texttt{off}, no forces are applied the atoms in the group. - Other forces are not affected (i.e. those - from the MD engine, from other colvars, and other external forces). - For dummy atoms, this option is \texttt{off} by default. - } - -\end{itemize} - + {cutoff}{% + \texttt{coordNum}}{% + Interaction'' distance (\AA)}{% + positive decimal}{% + 4.0}{% + This number defines the switching distance to define an + interatomic contact: for$d \ll d_0$, the switching function +$(1-(d/d_0)^{n})/(1-(d/d_0)^{m})$is close to 1, at$d = d_0$it + has a value of$n/m$($1/2$with the default$n$and$m$), and at +$d \gg d_0$it goes to zero approximately like$d^{m-n}$. Hence, + for a proper behavior,$m$must be larger than$n$.} -\cvsubsec{Treatment of periodic boundary conditions.}{sec:colvar_atom_groups_wrapping} +\item % + \labelkey{colvar|coordNum|cutoff3} + \keydef + {cutoff3}{% + \texttt{coordNum}}{% + Reference distance vector (\AA)}{% + \texttt{(x, y, z)}'' triplet of positive decimals}{% + \texttt{(4.0, 4.0, 4.0)}}{% + The three components of this vector define three different cutoffs +$d_{0}$for each direction. This option is mutually exclusive with + \texttt{cutoff}.} -\cvnamdonly{ - In simulations with periodic boundary conditions, NAMD maintains - the coordinates of all the atoms within a molecule contiguous to - each other (i.e.~there are no spurious jumps'' in the molecular - bonds). The Colvars module relies on this when calculating a group's - center of geometry, but this condition may fail if the group spans - different molecules. In that case, writing the NAMD output and restart files - using \texttt{wrapAll} or \texttt{wrapWater} could produce wrong results - when a simulation run is continued from a previous one. - The user should then determine, according to which - type of colvars are being calculated, whether \texttt{wrapAll} or - \texttt{wrapWater} can be enabled. +\item % + \labelkey{colvar|coordNum|expNumer} + \keydef + {expNumer}{% + \texttt{coordNum}}{% + Numerator exponent}{% + positive even integer}{% + 6}{% + This number defines the$n$exponent for the switching function.} - In general, internal coordinate wrapping by NAMD does not affect the calculation of colvars if each atom group satisfies one or more of the following: -} -\cvlammpsonly{ -In simulations with periodic boundary conditions, many of the implemented colvar components rely on the fact that each position within a group of atoms is at the nearest periodic image from the center of geometry of the group itself. -However, due to the internal wrapping of individual atomic positions done by LAMMPS, this assumption is broken if the group straddles one of the unit cell's boundaries. -For this reason, within the Colvars module all coordinates are unwrapped by default to avoid discontinuities (see \texttt{unwrap} keyword in \ref{sec:colvars_mdengine_parameters}). +\item % + \labelkey{colvar|coordNum|expDenom} + \keydef + {expDenom}{% + \texttt{coordNum}}{% + Denominator exponent}{% + positive even integer}{% + 12}{% + This number defines the$m$exponent for the switching function.} -The user should determine whether maintaining the default value of \texttt{unwrap}, depending on the specifics of each system. -In general, internal coordinate wrapping by LAMMPS does not affect the calculation of colvars if each atom group satisfies one or more of the following: -} -\cvvmdonly{ - When periodic boundary conditions are defined, the Colvars module requires that the coordinates of each molecular fragment are contiguous, without jumps'' when a fragment is partially wrapped near a periodic boundary. - The Colvars module relies on this assumption when calculating a group's center of geometry, but the condition may fail if the group spans different molecules. - In general, coordinate wrapping does not affect the calculation of colvars if each atom group satisfies one or more of the following: +\item % + \labelkey{colvar|coordNum|group2CenterOnly} + \keydef + {group2CenterOnly}{% + \texttt{coordNum}}{% + Use only \texttt{group2}'s center of + mass}{% + boolean}{% + \texttt{off}}{% + If this option is \texttt{on}, only contacts between each atoms in \texttt{group1} and the center of mass of \texttt{group2} are calculated (by default, the sum extends over all pairs of atoms in \texttt{group1} and \texttt{group2}). +If \texttt{group2} is a \texttt{dummyAtom}, this option is set to \texttt{yes} by default. } -\begin{enumerate} - \item[\emph{i)}] it is composed by only one atom; - \item[\emph{ii)}] it is used by a colvar component which does not make use of its center of geometry, but only of pairwise distances (\texttt{distanceInv}, \texttt{coordNum}, \texttt{hBond}, \texttt{alpha}, \texttt{dihedralPC}); - \item[\emph{iii)}] it is used by a colvar component that ignores the ill-defined Cartesian components of its center of mass (such as the$x$and$y$components of a membrane's center of mass modeled with \texttt{distanceZ})\cvnamdonly{; - \item[\emph{iv)}] it has all of its atoms within the same molecular fragment% -}. -\end{enumerate} -\cvvmdonly{If none of these conditions are met, wrapping may affect the calculation of collective variables: a possible solution is to use \texttt{pbc wrap} or \texttt{pbc unwrap} prior to processing a trajectory with the Colvars module.} - -\cvsubsec{Performance of a Colvars calculation based on group size.}{sec:colvar_atom_groups_scaling} - -In simulations performed with message-passing programs (such as NAMD or LAMMPS), the calculation of energy and forces is distributed (i.e., parallelized) across multiple nodes, as well as over the processor cores of each node. -When Colvars is enabled, certain atomic coordinates are collected on a single node, where the calculation of collective variables and of their biases is executed. -This means that for simulations over large numbers of nodes, a Colvars calculation may produce a significant overhead, coming from the costs of transmitting atomic coordinates to one node and of processing them. -\cvnamdonly{The latency-tolerant design and dynamic load balancing of NAMD may alleviate both factors, but a noticeable performance impact may be observed.} - -Performance can be improved in multiple ways: -\begin{itemize} -\item The calculation of variables, components and biases can be distributed over the processor cores of the node where the Colvars module is executed. - Currently, an equal weight is assigned to each colvar, or to each component of those colvars that include more than one component. - The performance of simulations that use many colvars or components is improved automatically. - For simulations that use a single large colvar, it may be advisable to partition it in multiple components, which will be then distributed across the available cores. - \cvnamdonly{In NAMD, this feature is enabled in all binaries compiled using SMP builds of Charm++ with the CkLoop extension.} - \cvlammpsonly{In LAMMPS, this feature is supported automatically when LAMMPS is compiled with OpenMP support.} - If printed, the message SMP parallelism is available.'' indicates the availability of the option\cvvmdonly{ (will be supported in a future relase of VMD)}. - If available, the option is turned on by default, but may be disabled using the keyword \refkey{smp}{Colvars-global|smp} if required for debugging. - -\cvnamdonly{ - % Use the following command to identify them: - % grep -B10 'provide(f_cvc_com_based' * |grep '\:\:'|grep '(std::string const &conf)' -\item NAMD also offers a parallelized calculation of the centers of mass of groups of atoms. - This option is on by default for all components that are simple functions of centers of mass, and is controlled by the keyword \refkey{scalable}{sec:cvc_common}. - When supported, the message Will enable scalable calculation for group \ldots'' is printed for each group. +\item % + \labelkey{colvar|coordNum|tolerance} + \keydef + {tolerance}{% + \texttt{coordNum}}{% + Pairlist control}{% + decimal}{% + 0.0}{This controls the pairlist feature, dictating the minimum value for each summation element in Eq.~\ref{eq:cvc_coordNum} such that the pair that contributed the summation element is included in subsequent simulation timesteps until the next pairlist recalculation. For most applications, this value should be small (eg. 0.001) to avoid missing important contributions to the overall sum. Higher values will improve performance, although values above 1 will exclude all possible pair interactions. Similarly, values below 0 will never exclude a pair from consideration. } -\item As a general rule, the size of atom groups should be kept relatively small (up to a few thousands of atoms, depending on the size of the entire system in comparison). -To gain an estimate of the computational cost of a large colvar, one can use a test calculation of the same colvar in VMD (hint: use the \texttt{time} Tcl command to measure the cost of running \texttt{cv update}). -\end{itemize} - - -\cvsec{Collective variable types (available functions)}{sec:cvc} - -In this context, the function that computes a colvar is here called a \emph{component}. -A component's choice and definition consists of a keyword indicating the -type of function (e.g.{} \texttt{rmsd}), followed by a definition block -specifying the atoms involved (see \ref{sec:colvar_atom_groups}) and any additional parameters (cutoffs, reference'' values, \ldots). -At least one component must be chosen: if none of the keywords listed below is found, an error is raised. - -Most components return a scalar value (i.e.{} a real number): -\begin{itemize} -\item \texttt{distance}: distance between two groups; -\item \texttt{distanceZ}: projection of a distance vector on an axis; -\item \texttt{distanceXY}: projection of a distance vector on a plane; -\item \texttt{distanceInv}: mean distance between two groups of atoms (e.g.~NOE-based distance); -\item \texttt{angle}: angle between three groups; -\item \texttt{dihedral}: torsional (dihedral) angle between four groups; -\item \texttt{dipoleAngle}: angle between two groups and dipole of a third group; -\item \texttt{polarTheta}: polar angle of a group in spherical coordinates; -\item \texttt{polarPhi}: azimuthal angle of a group in spherical coordinates; -\item \texttt{coordNum}: coordination number between two groups; -\item \texttt{selfCoordNum}: coordination number of atoms within a - group; -\item \texttt{hBond}: hydrogen bond between two atoms; -\item \texttt{rmsd}: root mean square deviation (RMSD) from a set of - reference coordinates; -\item \texttt{eigenvector}: projection of the atomic coordinates on a - vector; -\item \texttt{orientationAngle}: angle of the best-fit rotation from - a set of reference coordinates; -\item \texttt{orientationProj}: cosine of \texttt{orientationProj}; -\item \texttt{spinAngle}: projection orthogonal to an axis of the best-fit rotation - from a set of reference coordinates; -\item \texttt{tilt}: projection on an axis of the best-fit rotation - from a set of reference coordinates; -\item \texttt{gyration}: radius of gyration of a group of atoms; -\item \texttt{inertia}: moment of inertia of a group of atoms; -\item \texttt{inertiaZ}: moment of inertia of a group of atoms around a chosen axis; -\cvnamebasedonly{ -\item \texttt{alpha}:$\alpha$-helix content of a protein segment. -\item \texttt{dihedralPC}: projection of protein backbone dihedrals onto a dihedral principal component. +\item % + \labelkey{colvar|coordNum|pairListFrequency} + \keydef + {pairListFrequency}{% + \texttt{coordNum}}{% + Pairlist regeneration frequency}{% + positive integer}{% + 100}{This controls the pairlist feature, dictating how many steps are taken between regenerating pairlists if the tolerance is greater than 0. At this interval, the colvar defined will be exact, as though it was the all-to-all pair summation defined in Eq.~\ref{eq:cvc_coordNum}. All other steps will report only values and gradients from pairs in the pairlist. } -\end{itemize} +\end{cvcoptions} -Some components do not return scalar, but vector values: -\begin{itemize} -\item \texttt{distanceVec}: distance vector between two groups (length: 3); -\item \texttt{distanceDir}: unit vector parallel to distanceVec (length: 3); -\item \texttt{cartesian}: vector of atomic Cartesian coordinates (length:$N$); -\item \texttt{distancePairs}: vector of mutual distances (length:$N_{\mathrm{1}}\times{}N_{\mathrm{2}}$); -\item \texttt{orientation}: best-fit rotation, expressed as a unit quaternion (length: 4). -\end{itemize} +This component returns a dimensionless number, which ranges from +approximately 0 (all interatomic distances are much larger than the +cutoff) to$N_{\mathtt{group1}} \times N_{\mathtt{group2}}$(all distances +are less than the cutoff), or$N_{\mathtt{group1}}$if +\texttt{group2CenterOnly} is used. For performance reasons, at least +one of \texttt{group1} and \texttt{group2} should be of limited size or \texttt{group2CenterOnly} should be used: the cost of the loop over all pairs grows as$N_{\mathtt{group1}} \times N_{\mathtt{group2}}$. -The types of components used in a colvar (scalar or not) determine the -properties of that colvar, and particularly which biasing or analysis methods -can be applied. -\textbf{What if X'' is not listed?} If a function type is not available on this list, it may be possible to define it as a polynomial superposition of existing ones (see \ref{sec:cvc_superp})\cvleptononly{, a custom function (see \ref{sec:colvar_custom_function})}\cvscriptonly{, or a scripted function (see \ref{sec:colvar_scripted})}. +\cvsubsubsec{\texttt{selfCoordNum}: coordination number between atoms within a group.}{sec:cvc_selfCoordNum} +\labelkey{colvar|selfCoordNum} +The \texttt{selfCoordNum \{...\}} block defines +a coordination number similarly to the component \texttt{coordNum}, +but the function is summed over atom pairs within \texttt{group1}: + + \label{eq:cvc_selfCoordNum} + C (\mathtt{group1}) \; = \; + \sum_{i\in\mathtt{group1}}\sum_{j > i} { + \frac{1 - (|\mathbf{x}_{i}-\mathbf{x}_{j}|/d_{0})^{n}}{ + 1 - (|\mathbf{x}_{i}-\mathbf{x}_{j}|/d_{0})^{m} } + } + +The keywords accepted by \texttt{selfCoordNum} are a subset of +those accepted by \texttt{coordNum}, namely \texttt{group1} +(here defining \emph{all} of the atoms to be considered), +\texttt{cutoff}, \texttt{expNumer}, and \texttt{expDenom}. -\cvsubsec{List of available colvar components}{sec:cvc_list} +\begin{cvcoptions} +\item % + \dupkey{group1}{\texttt{selfCoordNum}}{colvar|coordNum|group1}{\texttt{coordNum} component} +\item % + \dupkey{cutoff}{\texttt{selfCoordNum}}{colvar|coordNum|cutoff}{\texttt{coordNum} component} +\item % + \dupkey{cutoff3}{\texttt{selfCoordNum}}{colvar|coordNum|cutoff3}{\texttt{coordNum} component} +\item % + \dupkey{expNumer}{\texttt{selfCoordNum}}{colvar|coordNum|expNumer}{\texttt{coordNum} component} +\item % + \dupkey{expDenom}{\texttt{selfCoordNum}}{colvar|coordNum|expDenom}{\texttt{coordNum} component} +\item % + \dupkey{tolerance}{\texttt{selfCoordNum}}{colvar|coordNum|tolerance}{\texttt{coordNum} component} +\item % + \dupkey{pairListFrequency}{\texttt{selfCoordNum}}{colvar|coordNum|pairListFrequency}{\texttt{coordNum} component} +\end{cvcoptions} -In this section, all available component types are listed, along -with their physical units and the ranges of values, if limited. Such -limiting values can be used to define \texttt{lowerBoundary} and -\texttt{upperBoundary} in the parent colvar. +This component returns a dimensionless number, which ranges from +approximately 0 (all interatomic distances much larger than the +cutoff) to$N_{\mathtt{group1}} \times (N_{\mathtt{group1}} - 1) / 2$(all +distances within the cutoff). For performance reasons, +\texttt{group1} should be of limited size, because the cost of the +loop over all pairs grows as$N_{\mathtt{group1}}^2$. -For each type of component, the available configurations keywords are listed: -when two components share certain keywords, the second component references to -the documentation of the first one that uses that keyword. -The very few keywords that are available for all types of components are listed in a separate section \ref{sec:cvc_common}. -\newenvironment{cvcoptions}% - {\noindent\textbf{List of keywords} (see also \ref{sec:cvc_superp} for additional options): - \begin{itemize}} - { - \end{itemize} -} -\cvsubsubsec{\texttt{distance}: center-of-mass distance between two groups.}{sec:cvc_distance} +\cvsubsubsec{\texttt{hBond}: hydrogen bond between two atoms.}{sec:cvc_hBond} +\labelkey{colvar|hBond} -The \texttt{distance \{...\}} block defines a distance component between the two atom groups, \texttt{group1} and \texttt{group2}. +The \texttt{hBond \{...\}} block defines a hydrogen +bond, implemented as a coordination number (eq.~\ref{eq:cvc_coordNum}) +between the donor and the acceptor atoms. Therefore, it accepts the +same options \texttt{cutoff} (with a different default value of +3.3~\AA{}), \texttt{expNumer} (with a default value of 6) and +\texttt{expDenom} (with a default value of 8). Unlike +\texttt{coordNum}, it requires two atom numbers, \texttt{acceptor} and +\texttt{donor}, to be defined. It returns an adimensional number, +with values between 0 (acceptor and donor far outside the cutoff +distance) and 1 (acceptor and donor much closer than the cutoff). \begin{cvcoptions} \item % - \labelkey{colvar|distance|group1} \key - {group1}{% - \texttt{distance}}{% - First group of atoms}{% - Block \texttt{group1 \{...\}}}{% - First group of atoms.} - + {acceptor}{% + \texttt{hBond}}{% + Number of the acceptor atom}{% + positive integer}{% + Number that uses the same convention as \texttt{atomNumbers}.} \item % - \labelkey{colvar|distance|group2} - \simkey{group2}{\texttt{distance}}{group1} - + \simkey{donor}{\texttt{hBond}}{acceptor} \item % - \labelkey{colvar|distance|forceNoPBC} - \keydef - {forceNoPBC}{% - \texttt{distance}}{% - Calculate absolute rather than minimum-image distance?}{% - boolean}{% - \texttt{no}}{% - By default, in calculations with periodic boundary conditions, the - \texttt{distance} component returns the distance according to the - minimum-image convention. If this parameter is set to \texttt{yes}, - PBC will be ignored and the distance between the coordinates as maintained - internally will be used. This is only useful in a limited number of - special cases, e.g. to describe the distance between remote points - of a single macromolecule, which cannot be split across periodic cell - boundaries, and for which the minimum-image distance might give the - wrong result because of a relatively small periodic cell.} - + \dupkey{cutoff}{\texttt{hBond}}{colvar|coordNum|cutoff}{\texttt{coordNum} component}\\ + \textbf{Note:} default value is 3.3~\AA. \item % - \labelkey{colvar|distance|oneSiteTotalForce} - \keydef - {oneSiteTotalForce}{% - \texttt{angle}, \texttt{dipoleAngle}, \texttt{dihedral}}{% - Measure total force on group 1 only?}{% - boolean}{% - \texttt{no}}{% - If this is set to \texttt{yes}, the total force is measured along - a vector field (see equation~(\ref{eq:gradient_vector}) in - section~\ref{sec:colvarbias_abf}) that only involves atoms of - \texttt{group1}. This option is only useful for ABF, or custom - biases that compute total forces. See - section~\ref{sec:colvarbias_abf} for details.} - + \dupkey{expNumer}{\texttt{hBond}}{colvar|coordNum|expNumer}{\texttt{coordNum} component}\\ + \textbf{Note:} default value is 6. +\item % + \dupkey{expDenom}{\texttt{hBond}}{colvar|coordNum|expDenom}{\texttt{coordNum} component}\\ + \textbf{Note:} default value is 8. \end{cvcoptions} -The value returned is a positive number (in \AA), ranging from$0$-to the largest possible interatomic distance within the chosen -boundary conditions (with PBCs, the minimum image convention is used -unless the \texttt{forceNoPBC} option is set). - -\cvsubsubsec{\texttt{distanceZ}: projection of a distance vector on an axis.}{sec:cvc_distanceZ} +\cvsubsubsec{\texttt{rmsd}: root mean square displacement (RMSD) from reference positions.}{sec:cvc_rmsd} +\labelkey{colvar|rmsd} -The \texttt{distanceZ~\{...\}} block defines a distance projection -component, which can be seen as measuring the distance between two -groups projected onto an axis, or the position of a group along such -an axis. The axis can be defined using either one reference group and -a constant vector, or dynamically based on two reference groups. -One of the groups can be set to a dummy atom to allow the use of an absolute Cartesian coordinate. +The block \texttt{rmsd~\{...\}} defines the root mean square replacement +(RMSD) of a group of atoms with respect to a reference structure. For +each set of coordinates$\{ \mathbf{x}_1(t), \mathbf{x}_2(t), \ldots
+\mathbf{x}_N(t) \}$, the colvar component \texttt{rmsd} calculates the +optimal rotation +$U^{\{\mathbf{x}_{i}(t)\}\rightarrow\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}}$+that best superimposes the coordinates$\{\mathbf{x}_{i}(t)\}$onto a +set of reference coordinates$\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}$. +Both the current and the reference coordinates are centered on their +centers of geometry,$\mathbf{x}_{\mathrm{cog}}(t)$and +$\mathbf{x}_{\mathrm{cog}}^{\mathrm{(ref)}}$. The root mean square +displacement is then defined as: + + \label{eq:cvc_rmsd} + { \mathrm{RMSD}(\{\mathbf{x}_{i}(t)\}, + \{\mathbf{x}_{i}^{\mathrm{(ref)}}\}) } \; = \; \sqrt{ + \frac{1}{N} \sum_{i=1}^{N} \left| + U + \left(\mathbf{x}_{i}(t) - \mathbf{x}_{\mathrm{cog}}(t)\right) - + \left(\mathbf{x}_{i}^{\mathrm{(ref)}} - + \mathbf{x}_{\mathrm{cog}}^{\mathrm{(ref)}} \right) \right|^{2} } + +The optimal rotation +$U^{\{\mathbf{x}_{i}(t)\}\rightarrow\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}}$+is calculated within the formalism developed in +reference~\cite{Coutsias2004}, which guarantees a continuous +dependence of +$U^{\{\mathbf{x}_{i}(t)\}\rightarrow\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}}$+with respect to$\{\mathbf{x}_{i}(t)\}$. \begin{cvcoptions} -\item % - \labelkey{colvar|distanceZ|main} - \key - {main}{% - \texttt{distanceZ}}{% - Main group of atoms}{% - Block \texttt{main \{...\}}}{% - Group of atoms whose position$\bm{r}$is measured.} \item % - \labelkey{colvar|distanceZ|ref} + \labelkey{colvar|rmsd|atoms} \key - {ref}{% - \texttt{distanceZ}}{% - Reference group of - atoms}{% - Block \texttt{ref \{...\}}}{% - Reference group of atoms. The position of its center of mass is - noted$\bm{r}_1$below.} + {atoms}{% + \texttt{rmsd}}{% + Atom group}{% + \texttt{atoms~\{...\}} block}{% + Defines the group of atoms of which the RMSD should be calculated. + Optimal fit options (such as \texttt{refPositions} and + \texttt{rotateReference}) should typically NOT be set within this + block. Exceptions to this rule are the special cases discussed in + the \emph{Advanced usage} paragraph below. + } \item % - \labelkey{colvar|distanceZ|ref2} - \keydef - {ref2}{% - \texttt{distanceZ}}{% - Secondary reference - group}{% - Block \texttt{ref2 \{...\}}}{% - none}{% - Optional group of reference atoms, whose position$\bm{r}_2$can - be used to define a dynamic projection axis:$\bm{e}=(\| \bm{r}_2
-    - \bm{r}_1\|)^{-1} \times (\bm{r}_2 - \bm{r}_1)$. In this case, - the origin is$\bm{r}_m = 1/2 (\bm{r}_1+\bm{r}_2)$, and the value - of the component is$\bm{e} \cdot (\bm{r}-\bm{r}_m)$.} + \labelkey{colvar|rmsd|refPositions} + \key + {refPositions}{% + \texttt{rmsd}}{% + Reference coordinates}{% + space-separated list of \texttt{(x, y, z)} triplets}{% + This option (mutually exclusive with \texttt{refPositionsFile}) sets the reference coordinates for RMSD calculation, and uses these to compute the roto-translational fit. + It is functionally equivalent to the option \refkey{refPositions}{atom-group|refPositions} in the atom group definition, which also supports more advanced fitting options. + } \item % - \labelkey{colvar|distanceZ|axis} - \keydef - {axis}{% - \texttt{distanceZ}}{% - Projection axis (\AA{})}{% - \texttt{(x, y, z)} triplet}{% - \texttt{(0.0, 0.0, 1.0)}}{% - The three components of this vector define a - projection axis$\bm{e}$for the distance vector$\bm{r} -
-    \bm{r}_1$joining the centers of groups \texttt{ref} and - \texttt{main}. The value of the component is then$\bm{e} \cdot
-    (\bm{r}-\bm{r}_1)$. The vector should be written as three - components separated by commas and enclosed in parentheses.} + \labelkey{colvar|rmsd|refPositionsFile} + \key + {refPositionsFile}{% + \texttt{rmsd}}{% + Reference coordinates file}{% + UNIX filename}{% + This option (mutually exclusive with \texttt{refPositions}) sets the reference coordinates for RMSD calculation, and uses these to compute the roto-translational fit. + It is functionally equivalent to the option \refkey{refPositionsFile}{atom-group|refPositionsFile} in the atom group definition, which also supports more advanced fitting options. + } +\cvnamebasedonly{ \item % - \dupkey{forceNoPBC}{\texttt{distanceZ}}{colvar|distance|forceNoPBC}{\texttt{distance} component} + \labelkey{colvar|rmsd|refPositionsCol} + \key + {refPositionsCol}{% + \texttt{rmsd}}{% + PDB column containing atom flags}{% + \texttt{O}, \texttt{B}, \texttt{X}, \texttt{Y}, or \texttt{Z}}{% + If \texttt{refPositionsFile} is a PDB file that contains all the atoms in the topology, this option may be provided to set which PDB field is used to flag the reference coordinates for \texttt{atoms}. + } \item % - \dupkey{oneSiteTotalForce}{\texttt{distanceZ}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} + \labelkey{colvar|rmsd|refPositionsColValue} + \key + {refPositionsColValue}{% + \texttt{rmsd}}{% + Atom selection flag in the PDB column}{% + positive decimal}{% + If defined, this value identifies in the PDB column + \texttt{refPositionsCol} of the file \texttt{refPositionsFile} + which atom positions are to be read. Otherwise, all positions + with a non-zero value are read. + } +} \end{cvcoptions} -This component returns a number (in \AA{}) whose range is determined -by the chosen boundary conditions. For instance, if the$z$axis is -used in a simulation with periodic boundaries, the returned value ranges -between$-b_{z}/2$and$b_{z}/2$, where$b_{z}$is the box length -along$z$(this behavior is disabled if \texttt{forceNoPBC} is set). +This component returns a positive real number (in \AA). -\cvsubsubsec{\texttt{distanceXY}: modulus of the projection of a distance vector on a plane.}{sec:cvc_distanceXY} -The \texttt{distanceXY~\{...\}} block defines a distance projected on -a plane, and accepts the same keywords as the component \texttt{distanceZ}, i.e. -\texttt{main}, \texttt{ref}, either \texttt{ref2} or \texttt{axis}, -and \texttt{oneSiteTotalForce}. It returns the norm of the -projection of the distance vector between \texttt{main} and -\texttt{ref} onto the plane orthogonal to the axis. The axis is -defined using the \texttt{axis} parameter or as the vector joining -\texttt{ref} and \texttt{ref2} (see \texttt{distanceZ} above). +\cvsubsubsec{Advanced usage of the \texttt{rmsd} component.}{sec:cvc_rmsd_advanced} +In the standard usage as described above, the \texttt{rmsd} component +calculates a minimum RMSD, that is, current coordinates are optimally +fitted onto the same reference coordinates that are used to +compute the RMSD value. The fit itself is handled by the atom group +object, whose parameters are automatically set by the \texttt{rmsd} +component. +For very specific applications, however, it may be +useful to control the fitting process separately from the definition +of the reference coordinates, to evaluate various types of +non-minimal RMSD values. This can be achieved by setting the +related options (\texttt{refPositions}, etc.) explicitly in the +atom group block. This allows for the following non-standard cases: -\begin{cvcoptions} -\item % - \dupkey{main}{\texttt{distanceXY}}{colvar|distanceZ|main}{\texttt{distanceZ} component} -\item % - \dupkey{ref}{\texttt{distanceXY}}{colvar|distanceZ|ref}{\texttt{distanceZ} component} -\item % - \dupkey{ref2}{\texttt{distanceXY}}{colvar|distanceZ|ref2}{\texttt{distanceZ} component} -\item % - \dupkey{axis}{\texttt{distanceXY}}{colvar|distanceZ|axis}{\texttt{distanceZ} component} -\item % - \dupkey{forceNoPBC}{\texttt{distanceXY}}{colvar|distance|forceNoPBC}{\texttt{distance} component} -\item % - \dupkey{oneSiteTotalForce}{\texttt{distanceZ}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} -\end{cvcoptions} +\begin{enumerate} +\item applying the optimal translation, but no rotation +(\texttt{rotateReference off}), to bias or restrain the shape and +orientation, but not the position of the atom group; +\item applying the optimal rotation, but no translation +(\texttt{translateReference off}), to bias or restrain the shape and +position, but not the orientation of the atom group; +\item disabling the application of optimal roto-translations, which +lets the RMSD component decribe the deviation of atoms +from fixed positions in the laboratory frame: this allows for custom +positional restraints within the Colvars module; +\item fitting the atomic positions to different reference coordinates +than those used in the RMSD calculation itself; +\item applying the optimal rotation and/or translation from a separate +atom group, defined through \texttt{fittingGroup}: the RMSD then +reflects the deviation from reference coordinates in a separate, moving +reference frame. +\end{enumerate} +\cvscriptonly{ +\cvsubsubsec{Path collective variables}{sec:pathcv} -\cvsubsubsec{\texttt{distanceVec}: distance vector between two groups.}{sec:cvc_distanceVec} -The \texttt{distanceVec~\{...\}} block defines -a distance vector component, which accepts the same keywords as -the component \texttt{distance}: \texttt{group1}, \texttt{group2}, and -\texttt{forceNoPBC}. Its value is the 3-vector joining the centers -of mass of \texttt{group1} and \texttt{group2}. +An application of the \texttt{rmsd} component is "path collective variables",\cite{Branduardi2007} +which are implemented as Tcl-scripted combinations or RMSDs. +The implementation is available as file \texttt{colvartools/pathCV.tcl}, and +an example is provided in file \texttt{examples/10\_pathCV.namd} of the Colvars public repository.} + + +\cvsubsubsec{\texttt{eigenvector}: projection of the atomic coordinates on a vector.}{sec:cvc_eigenvector} +\labelkey{colvar|eigenvector} + +The block \texttt{eigenvector~\{...\}} defines the projection of the coordinates +of a group of atoms (or more precisely, their deviations from the +reference coordinates) onto a vector in$\mathbb{R}^{3n}$, where$n$is the +number of atoms in the group. The computed quantity is the +total projection: + + \label{eq:cvc_eigenvector} + { p(\{\mathbf{x}_{i}(t)\}, + \{\mathbf{x}_{i}^{\mathrm{(ref)}}\}) } \; = \; { + \sum_{i=1}^{n} \mathbf{v}_{i} \cdot + \left(U(\mathbf{x}_{i}(t) - \mathbf{x}_{\mathrm{cog}}(t)) - + (\mathbf{x}_{i}^{\mathrm{(ref)}} - + \mathbf{x}_{\mathrm{cog}}^{\mathrm{(ref)}}) \right)\mathrm{,} } + +where, as in the \texttt{rmsd} component,$U$is the optimal rotation +matrix,$\mathbf{x}_{\mathrm{cog}}(t)$and +$\mathbf{x}_{\mathrm{cog}}^{\mathrm{(ref)}}$are the centers of +geometry of the current and reference positions respectively, and +$\mathbf{v}_{i}$are the components of the vector for each atom. +Example choices for$(\mathbf{v}_{i})$are an eigenvector +of the covariance matrix (essential mode), or a normal +mode of the system. It is assumed that$\sum_{i}\mathbf{v}_{i} = 0$: +otherwise, the Colvars module centers the$\mathbf{v}_{i}$+automatically when reading them from the configuration. \begin{cvcoptions} \item % - \dupkey{group1}{\texttt{distanceVec}}{colvar|distance|group1}{\texttt{distance} component} + \dupkey{atoms}{\texttt{eigenvector}}{colvar|rmsd|atoms}{\texttt{rmsd} component} \item % - \simkey{group2}{\texttt{distanceVec}}{group1} + \dupkey{refPositions}{\texttt{eigenvector}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} \item % - \dupkey{forceNoPBC}{\texttt{distanceVec}}{colvar|distance|forceNoPBC}{\texttt{distance} component} + \dupkey{refPositionsFile}{\texttt{eigenvector}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} +\cvnamebasedonly{ \item % - \dupkey{oneSiteTotalForce}{\texttt{distanceVec}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} -\end{cvcoptions} - - + \dupkey{refPositionsCol}{\texttt{eigenvector}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} +\item % + \dupkey{refPositionsColValue}{\texttt{eigenvector}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} +} -\cvsubsubsec{\texttt{distanceDir}: distance unit vector between two groups.}{sec:cvc_distanceDir} -The \texttt{distanceDir~\{...\}} block defines -a distance unit vector component, which accepts the same keywords as -the component \texttt{distance}: \texttt{group1}, \texttt{group2}, and -\texttt{forceNoPBC}. It returns a -3-dimensional unit vector$\mathbf{d} = (d_{x}, d_{y}, d_{z})$, with -$|\mathbf{d}| = 1$. +\item % + \key + {vector}{% + \texttt{eigenvector}}{% + Vector components}{% + space-separated list of \texttt{(x, y, z)} triplets}{% + This option (mutually exclusive with \texttt{vectorFile}) sets the values of the vector components.} -\begin{cvcoptions} \item % - \dupkey{group1}{\texttt{distanceDir}}{colvar|distance|group1}{\texttt{distance} component} + \key + {vectorFile}{% + \texttt{eigenvector}}{% + file containing vector components}{% + UNIX filename}{% + This option (mutually exclusive with \texttt{vector}) sets the name of a coordinate file containing the vector components; the file is read according to the same format used for \texttt{refPositionsFile}. + \cvnamebasedonly{For a PDB file specifically, the components are read from the X, Y and Z fields. + \textbf{Note:} \emph{The PDB file has limited precision and fixed-point numbers: in some cases, the vector components may not be accurately represented; a XYZ file should be used instead, containing floating-point numbers.}} + } + +\cvnamebasedonly{ \item % - \simkey{group2}{\texttt{distanceDir}}{group1} + \key + {vectorCol}{% + \texttt{eigenvector}}{% + PDB column used to flag participating atoms}{% + \texttt{O} or \texttt{B}}{% + Analogous to \texttt{atomsCol}.} + \item % - \dupkey{forceNoPBC}{\texttt{distanceDir}}{colvar|distance|forceNoPBC}{\texttt{distance} component} + \key + {vectorColValue}{% + \texttt{eigenvector}}{% + Value used to flag participating atoms in the PDB file}{% + positive decimal}{% + Analogous to \texttt{atomsColValue}.} +} + \item % - \dupkey{oneSiteTotalForce}{\texttt{distanceDir}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} + \keydef + {differenceVector}{% + \texttt{eigenvector}}{% + The$3n$-dimensional vector is the difference between \texttt{vector} and \texttt{refPositions}}{% + boolean}{% + \texttt{off}}{% + If this option is \texttt{on}, the numbers provided by \texttt{vector}\cvnamebasedonly{ or \texttt{vectorFile}} are interpreted as another set of positions,$\mathbf{x}_{i}'$: the vector$\mathbf{v}_{i}$is then defined as$\mathbf{v}_{i} = \left(\mathbf{x}_{i}' - \mathbf{x}_{i}^{\mathrm{(ref)}}\right)$. +This allows to conveniently define a colvar$\xi$as a projection on the linear transformation between two sets of positions, A'' and B''. +For convenience, the vector is also normalized so that$\xi = 0$when the atoms are at the set of positions A'' and$\xi = 1$at the set of positions B''. +} \end{cvcoptions} +This component returns a number (in \AA), whose value ranges between +the smallest and largest absolute positions in the unit cell during +the simulations (see also \texttt{distanceZ}). Due to the +normalization in eq.~\ref{eq:cvc_eigenvector}, this range does not +depend on the number of atoms involved. -\cvsubsubsec{\texttt{distanceInv}: mean distance between two groups of atoms.}{sec:cvc_distanceInv} -The \texttt{distanceInv~\{...\}} block defines a generalized mean distance between two groups of atoms 1 and 2, weighted with exponent$1/n$: +\cvsubsubsec{\texttt{gyration}: radius of gyration of a group of atoms.}{sec:cvc_gyration} +\labelkey{colvar|gyration} + +The block \texttt{gyration~\{...\}} defines the +parameters for calculating the radius of gyration of a group of atomic +positions$\{ \mathbf{x}_1(t), \mathbf{x}_2(t), \ldots \mathbf{x}_N(t)
+\}$with respect to their center of geometry, +$\mathbf{x}_{\mathrm{cog}}(t)$: - \label{eq:distanceInv} - d_{\mathrm{1,2}}^{[n]} \; = \; \left(\frac{1}{N_{\mathrm{1}}N_{\mathrm{2}}}\sum_{i,j} \left(\frac{1}{\Vert\mathbf{d}^{ij}\Vert}\right)^{n} \right)^{-1/n} + \label{eq:colvar_gyration} + R_{\mathrm{gyr}} \; = \; \sqrt{ \frac{1}{N} + \sum_{i=1}^{N} \left|\mathbf{x}_{i}(t) - + \mathbf{x}_{\mathrm{cog}}(t)\right|^{2} } -where$\Vert\mathbf{d}^{ij}\Vert$is the distance between atoms$i$and$j$in groups 1 and 2 respectively, and$n$is an even integer. +This component must contain one \texttt{atoms~\{...\}} block to +define the atom group, and returns a positive number, expressed in +\AA{}. \begin{cvcoptions} \item % - \dupkey{group1}{\texttt{distanceInv}}{colvar|distance|group1}{\texttt{distance} component} -\item % - \simkey{group2}{\texttt{distanceInv}}{group1} -\item % - \dupkey{oneSiteTotalForce}{\texttt{distanceInv}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} -\item % - \keydef - {exponent}{% - \texttt{distanceInv}}{% - Exponent$n$in equation~\ref{eq:distanceInv}}{% - positive even integer}{% - 6}{Defines the exponent to which the individual distances are elevated before averaging. The default value of 6 is useful for example to applying restraints based on NOE-measured distances.} + \dupkey{atoms}{\texttt{gyration}}{colvar|rmsd|atoms}{\texttt{rmsd} component} \end{cvcoptions} -This component returns a number in \AA{}, ranging from$0$to the largest possible distance within the chosen boundary conditions. -\cvsubsubsec{\texttt{distancePairs}: set of pairwise distances between two groups.}{sec:cvc_distancePairs} -The \texttt{distancePairs~\{...\}} block defines a$N_{\mathrm{1}}\times{}N_{\mathrm{2}}$-dimensional variable that includes all mutual distances between the atoms of two groups. -This can be useful, for example, to develop a new variable defined over two groups, by using the \texttt{scriptedFunction} feature. +\cvsubsubsec{\texttt{inertia}: total moment of inertia of a group of atoms.}{sec:cvc_inertia} +\labelkey{colvar|inertia} + +The block \texttt{inertia~\{...\}} defines the +parameters for calculating the total moment of inertia of a group of atomic +positions$\{ \mathbf{x}_1(t), \mathbf{x}_2(t), \ldots \mathbf{x}_N(t)
+\}$with respect to their center of geometry, +$\mathbf{x}_{\mathrm{cog}}(t)$: + + \label{eq:colvar_inertia} + I \; = \; \sum_{i=1}^{N} \left|\mathbf{x}_{i}(t) - + \mathbf{x}_{\mathrm{cog}}(t)\right|^{2} + +\emph{Note that all atomic masses are set to 1 for simplicity.} +This component must contain one \texttt{atoms~\{...\}} block to +define the atom group, and returns a positive number, expressed in +\AA{}$^{2}$. \begin{cvcoptions} \item % - \dupkey{group1}{\texttt{distancePairs}}{colvar|distance|group1}{\texttt{distance} component} -\item % - \simkey{group2}{\texttt{distancePairs}}{group1} -\item % - \dupkey{forceNoPBC}{\texttt{distancePairs}}{colvar|distance|forceNoPBC}{\texttt{distance} component} + \dupkey{atoms}{\texttt{inertia}}{colvar|rmsd|atoms}{\texttt{rmsd} component} \end{cvcoptions} -This component returns a$N_{\mathrm{1}}\times{}N_{\mathrm{2}}$-dimensional vector of numbers, each ranging from$0$to the largest possible distance within the chosen boundary conditions. -\cvsubsubsec{\texttt{cartesian}: vector of atomic Cartesian coordinates.}{sec:cvc_cartesian} -The \texttt{cartesian~\{...\}} block defines a component returning a flat vector containing -the Cartesian coordinates of all participating atoms, in the order -$(x_1, y_1, z_1, \cdots, x_n, y_n, z_n)$. +\cvsubsubsec{\texttt{inertiaZ}: total moment of inertia of a group of atoms around a chosen axis.}{sec:cvc_inertiaZ} +\labelkey{colvar|inertiaZ} + +The block \texttt{inertiaZ~\{...\}} defines the +parameters for calculating the component along the axis$\mathbf{e}$of the moment of inertia of a group of atomic +positions$\{ \mathbf{x}_1(t), \mathbf{x}_2(t), \ldots \mathbf{x}_N(t)
+\}$with respect to their center of geometry, +$\mathbf{x}_{\mathrm{cog}}(t)$: + + \label{eq:colvar_inertia_z} + I_{\mathbf{e}} \; = \; \sum_{i=1}^{N} \left(\left(\mathbf{x}_{i}(t) - + \mathbf{x}_{\mathrm{cog}}(t)\right)\cdot\mathbf{e}\right)^{2} + +\emph{Note that all atomic masses are set to 1 for simplicity.} +This component must contain one \texttt{atoms~\{...\}} block to +define the atom group, and returns a positive number, expressed in +\AA{}$^{2}$. \begin{cvcoptions} \item % - \key - {atoms}{% - \texttt{cartesian}}{% - Group of atoms}{% - Block \texttt{atoms \{...\}}}{% - Defines the atoms whose coordinates make up the value of the component. - If \texttt{rotateReference} or \texttt{centerReference} are defined, coordinates - are evaluated within the moving frame of reference.} + \dupkey{atoms}{\texttt{inertiaZ}}{colvar|rmsd|atoms}{\texttt{rmsd} component} +\item % + \keydef + {axis}{% + \texttt{inertiaZ}}{% + Projection axis (\AA{})}{% + \texttt{(x, y, z)} triplet}{% + \texttt{(0.0, 0.0, 1.0)}}{% + The three components of this vector define (when normalized) the + projection axis$\mathbf{e}$.} \end{cvcoptions} -\cvsubsubsec{\texttt{angle}: angle between three groups.}{sec:cvc_angle} -The \texttt{angle~\{...\}} block defines an angle, and contains the -three blocks \texttt{group1}, \texttt{group2} and \texttt{group3}, defining -the three groups. It returns an angle (in degrees) within the -interval$[0:180]$. +\cvsubsubsec{\texttt{orientation}: orientation from reference coordinates.}{sec:cvc_orientation} +\labelkey{colvar|orientation} + +The block \texttt{orientation~\{...\}} returns the +same optimal rotation used in the \texttt{rmsd} component to +superimpose the coordinates$\{\mathbf{x}_{i}(t)\}$onto a set of +reference coordinates$\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}$. Such +component returns a four dimensional vector$\mathsf{q} = (q_0, q_1,
+q_2, q_3)$, with$\sum_{i} q_{i}^{2} = 1$; this \emph{quaternion} +expresses the optimal rotation$\{\mathbf{x}_{i}(t)\} \rightarrow
+\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}$according to the formalism in +reference~\cite{Coutsias2004}. The quaternion$(q_0, q_1, q_2, q_3)$+can also be written as$\left(\cos(\theta/2), \,
+  \sin(\theta/2)\mathbf{u}\right)$, where$\theta$is the angle and +$\mathbf{u}$the normalized axis of rotation; for example, a rotation +of 90$^{\circ}$around the$z$axis is expressed as +\texttt{(0.707, 0.0, 0.0, 0.707)}''. The script +\texttt{quaternion2rmatrix.tcl} provides Tcl functions for converting +to and from a$4\times{}4$rotation matrix in a format suitable for +usage in VMD. + +As for the component \texttt{rmsd}, the available options are \texttt{atoms}, \texttt{refPositionsFile}\cvnamebasedonly{, \texttt{refPositionsCol} and \texttt{refPositionsColValue}, } and \texttt{refPositions}. + +\textbf{Note:} \texttt{refPositions}and \texttt{refPositionsFile} define the set of positions \emph{from which} the optimal rotation is calculated, but this rotation is not applied to the coordinates of the atoms involved: it is used instead to define the variable itself. \begin{cvcoptions} \item % - \dupkey{group1}{\texttt{angle}}{colvar|distance|group1}{\texttt{distance} component} + \dupkey{atoms}{\texttt{orientation}}{colvar|rmsd|atoms}{\texttt{rmsd} component} \item % - \simkey{group2}{\texttt{angle}}{group1} + \dupkey{refPositions}{\texttt{orientation}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} \item % - \simkey{group3}{\texttt{angle}}{group1} + \dupkey{refPositionsFile}{\texttt{orientation}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} + +\cvnamebasedonly{ \item % - \dupkey{forceNoPBC}{\texttt{angle}}{colvar|distance|forceNoPBC}{\texttt{distance} component} + \dupkey{refPositionsCol}{\texttt{orientation}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} \item % - \dupkey{oneSiteTotalForce}{\texttt{angle}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} + \dupkey{refPositionsColValue}{\texttt{orientation}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} +} + +\item % + \keydef + {closestToQuaternion}{% + \texttt{orientation}}{% + Reference rotation}{% + \texttt{(q0, q1, q2, q3)}'' quadruplet}{% + \texttt{(1.0, 0.0, 0.0, 0.0)} (null'' rotation)}{% + Between the two equivalent quaternions$(q_0, q_1, q_2, q_3)$and +$(-q_0, -q_1, -q_2, -q_3)$, the closer to \texttt{(1.0, 0.0, 0.0, + 0.0)} is chosen. This simplifies the visualization of the + colvar trajectory when sampled values are a smaller subset of all + possible rotations. \textbf{Note:} \emph {this only affects the + output, never the dynamics}.} + \end{cvcoptions} +\textbf{Tip: stopping the rotation of a protein.} To stop the +rotation of an elongated macromolecule in solution (and use an +anisotropic box to save water molecules), it is possible to define a +colvar with an \texttt{orientation} component, and restrain it throuh +the \texttt{harmonic} bias around the identity rotation, \texttt{(1.0, + 0.0, 0.0, 0.0)}. Only the overall orientation of the macromolecule +is affected, and \emph{not} its internal degrees of freedom. The user +should also take care that the macromolecule is composed by a single +chain, or disable \texttt{wrapAll} otherwise. + -\cvsubsubsec{\texttt{dipoleAngle}: angle between two groups and dipole of a third group.}{sec:cvc_dipoleAngle} -The \texttt{dipoleAngle~\{...\}} block defines an angle, and contains the -three blocks \texttt{group1}, \texttt{group2} and \texttt{group3}, defining -the three groups, being \texttt{group1} the group where dipole is calculated. -It returns an angle (in degrees) within the interval$[0:180]$. +\cvsubsubsec{\texttt{orientationAngle}: angle of rotation from reference coordinates.}{sec:cvc_orientationAngle} +\labelkey{colvar|orientationAngle} + +The block \texttt{orientationAngle~\{...\}} accepts the same base options as +the component \texttt{orientation}: \texttt{atoms}, \texttt{refPositions}, \texttt{refPositionsFile}\cvnamebasedonly{, \texttt{refPositionsCol} and \texttt{refPositionsColValue}}. +The returned value is the angle of rotation$\theta$between the current and the reference positions. +This angle is expressed in degrees within the range [0$^{\circ}$:180$^{\circ}$]. \begin{cvcoptions} \item % - \dupkey{group1}{\texttt{dipoleAngle}}{colvar|distance|group1}{\texttt{distance} component} + \dupkey{atoms}{\texttt{orientationAngle}}{colvar|rmsd|atoms}{\texttt{rmsd} component} \item % - \simkey{group2}{\texttt{dipoleAngle}}{group1} + \dupkey{refPositions}{\texttt{orientationAngle}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} \item % - \simkey{group3}{\texttt{dipoleAngle}}{group1} + \dupkey{refPositionsFile}{\texttt{orientationAngle}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} + +\cvnamebasedonly{ \item % - \dupkey{forceNoPBC}{\texttt{dipoleAngle}}{colvar|distance|forceNoPBC}{\texttt{distance} component} + \dupkey{refPositionsCol}{\texttt{orientationAngle}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} \item % - \dupkey{oneSiteTotalForce}{\texttt{dipoleAngle}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} + \dupkey{refPositionsColValue}{\texttt{orientationAngle}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} +} + \end{cvcoptions} -\cvsubsubsec{\texttt{dihedral}: torsional angle between four groups.}{sec:cvc_dihedral} -The \texttt{dihedral~\{...\}} block defines a torsional angle, and -contains the blocks \texttt{group1}, \texttt{group2}, \texttt{group3} -and \texttt{group4}, defining the four groups. It returns an angle -(in degrees) within the interval$[-180:180]$. The Colvars module -calculates all the distances between two angles taking into account -periodicity. For instance, reference values for restraints or range -boundaries can be defined by using any real number of choice. +\cvsubsubsec{\texttt{orientationProj}: cosine of the angle of rotation from reference coordinates.} {sec:cvc_orientationProj} +\labelkey{colvar|orientationProj} + +The block \texttt{orientationProj~\{...\}} accepts the same base options as +the component \texttt{orientation}: \texttt{atoms}, \texttt{refPositions}, \texttt{refPositionsFile}\cvnamebasedonly{, \texttt{refPositionsCol} and \texttt{refPositionsColValue}}. +The returned value is the cosine of the angle of rotation$\theta$between the current and the reference positions. +The range of values is [-1:1]. \begin{cvcoptions} \item % - \dupkey{group1}{\texttt{dihedral}}{colvar|distance|group1}{\texttt{distance} component} -\item % - \simkey{group2}{\texttt{dihedral}}{group1} + \dupkey{atoms}{\texttt{orientationProj}}{colvar|rmsd|atoms}{\texttt{rmsd} component} \item % - \simkey{group3}{\texttt{dihedral}}{group1} + \dupkey{refPositions}{\texttt{orientationProj}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} \item % - \simkey{group4}{\texttt{dihedral}}{group1} + \dupkey{refPositionsFile}{\texttt{orientationProj}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} +\cvnamebasedonly{ \item % - \dupkey{forceNoPBC}{\texttt{dihedral}}{colvar|distance|forceNoPBC}{\texttt{distance} component} + \dupkey{refPositionsCol}{\texttt{orientationProj}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} \item % - \dupkey{oneSiteTotalForce}{\texttt{dihedral}}{colvar|distance|oneSiteTotalForce}{\texttt{distance} component} + \dupkey{refPositionsColValue}{\texttt{orientationProj}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} +} \end{cvcoptions} -\cvsubsubsec{\texttt{polarTheta}: polar angle in spherical coordinates.}{sec:cvc_polarTheta} -The \texttt{polarTheta~\{...\}} block defines the polar angle in -spherical coordinates, for the center of mass of a group of atoms -described by the block \texttt{atoms}. It returns an angle -(in degrees) within the interval$[0:180]$. -To obtain spherical coordinates in a frame of reference tied to -another group of atoms, use the \texttt{fittingGroup} (\ref{sec:colvar_atom_groups_ref_frame}) option -within the \texttt{atoms} block. -An example is provided in file \texttt{examples/11\_polar\_angles.in} of the Colvars public repository. +\cvsubsubsec{\texttt{spinAngle}: angle of rotation around a given axis.}{sec:cvc_spinAngle} +\labelkey{colvar|spinAngle} + +The complete rotation described by \texttt{orientation} can optionally be decomposed into two sub-rotations: one is a \emph{spin}'' rotation around \textbf{e}, and the other a \emph{tilt}'' rotation around an axis orthogonal to \textbf{e}. +The component \texttt{spinAngle} measures the angle of the spin'' sub-rotation around \textbf{e}. \begin{cvcoptions} \item % - \labelkey{colvar|polarTheta|atoms} - \key - {atoms}{% - \texttt{polarPhi}}{% - Atom group}{% - \texttt{atoms~\{...\}} block}{% - Defines the group of atoms for the COM of which the angle should be calculated. - } + \dupkey{atoms}{\texttt{spinAngle}}{colvar|rmsd|atoms}{\texttt{rmsd} component} +\item % + \dupkey{refPositions}{\texttt{spinAngle}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} +\item % + \dupkey{refPositionsFile}{\texttt{spinAngle}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} +\cvnamebasedonly{ +\item % + \dupkey{refPositionsCol}{\texttt{spinAngle}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} +\item % + \dupkey{refPositionsColValue}{\texttt{spinAngle}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} +} +\item % + \labelkey{colvar|spinAngle|axis} + \keydef + {axis}{% + \texttt{tilt}}{% + Special rotation axis (\AA{})}{% + \texttt{(x, y, z)} triplet}{% + \texttt{(0.0, 0.0, 1.0)}}{% + The three components of this vector define (when normalized) the special rotation axis used to calculate the \texttt{tilt} and \texttt{spinAngle} components.} \end{cvcoptions} +The component \texttt{spinAngle} returns an angle (in degrees) within the periodic interval$[-180:180]$. +\textbf{Note:} the value of \texttt{spinAngle} is a continuous function almost everywhere, with the exception of configurations with the corresponding tilt'' angle equal to 180$^\circ$(i.e.~the \texttt{tilt} component is equal to$-1$): in those cases, \texttt{spinAngle} is undefined. If such configurations are expected, consider defining a \texttt{tilt} colvar using the same axis \textbf{e}, and restraining it with a lower wall away from$-1$. -\cvsubsubsec{\texttt{polarPhi}: azimuthal angle in spherical coordinates.}{sec:cvc_polarPhi} -The \texttt{polarPhi~\{...\}} block defines the azimuthal angle in -spherical coordinates, for the center of mass of a group of atoms -described by the block \texttt{atoms}. It returns an angle -(in degrees) within the interval$[-180:180]$. The Colvars module -calculates all the distances between two angles taking into account -periodicity. For instance, reference values for restraints or range -boundaries can be defined by using any real number of choice. -To obtain spherical coordinates in a frame of reference tied to -another group of atoms, use the \texttt{fittingGroup} (\ref{sec:colvar_atom_groups_ref_frame}) option -within the \texttt{atoms} block. -An example is provided in file \texttt{examples/11\_polar\_angles.in} of the Colvars public repository. +\cvsubsubsec{\texttt{tilt}: cosine of the rotation orthogonal to a given axis.}{sec:cvc_tilt} +\labelkey{colvar|tilt} + +The component \texttt{tilt} measures the cosine of the angle of the tilt'' sub-rotation, which combined with the spin'' sub-rotation provides the complete rotation of a group of atoms. +The cosine of the tilt angle rather than the tilt angle itself is implemented, because the latter is unevenly distributed even for an isotropic system: consider as an analogy the angle$\theta$in the spherical coordinate system. +The component \texttt{tilt} relies on the same options as \texttt{spinAngle}, including the definition of the axis \textbf{e}. +The values of \texttt{tilt} are real numbers in the interval$[-1:1]$: the value$1$represents an orientation fully parallel to \textbf{e} (tilt angle = 0$^\circ$), and the value$-1$represents an anti-parallel orientation. \begin{cvcoptions} \item % - \labelkey{colvar|polarPhi|atoms} - \key - {atoms}{% - \texttt{polarPhi}}{% - Atom group}{% - \texttt{atoms~\{...\}} block}{% - Defines the group of atoms for the COM of which the angle should be calculated. - } + \dupkey{atoms}{\texttt{tilt}}{colvar|rmsd|atoms}{\texttt{rmsd} component} +\item % + \dupkey{refPositions}{\texttt{tilt}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} +\item % + \dupkey{refPositionsFile}{\texttt{tilt}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} +\cvnamebasedonly{ +\item % + \dupkey{refPositionsCol}{\texttt{tilt}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} +\item % + \dupkey{refPositionsColValue}{\texttt{tilt}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} +} +\item % + \dupkey{axis}{\texttt{tilt}}{colvar|spinAngle|axis}{\texttt{spinAngle} component} \end{cvcoptions} +\cvnamebasedonly{ -\cvsubsubsec{\texttt{coordNum}: coordination number between two groups.}{sec:cvc_coordNum} -The \texttt{coordNum \{...\}} block defines -a coordination number (or number of contacts), which calculates the -function$(1-(d/d_0)^{n})/(1-(d/d_0)^{m})$, where$d_0$is the -cutoff'' distance, and$n$and$m$are exponents that can control -its long range behavior and stiffness \cite{Iannuzzi2003}. This -function is summed over all pairs of atoms in \texttt{group1} and -\texttt{group2}: +\cvsubsubsec{\texttt{alpha}:$\alpha$-helix content of a protein segment.}{sec:cvc_alpha} +\labelkey{colvar|alpha} + +The block \texttt{alpha~\{...\}} defines the +parameters to calculate the helical content of a segment of protein +residues. The$\alpha$-helical content across the$N+1$residues +$N_{0}$to$N_{0}+N$is calculated by the formula: +\begin{eqnarray} + \label{eq:colvars_alpha} + { + \alpha\left( + \mathrm{C}_{\alpha}^{(N_{0})}, + \mathrm{O}^{(N_{0})}, + \mathrm{C}_{\alpha}^{(N_{0}+1)}, + \mathrm{O}^{(N_{0}+1)}, + \ldots + \mathrm{N}^{(N_{0}+5)}, + \mathrm{C}_{\alpha}^{(N_{0}+5)}, + \mathrm{O}^{(N_{0}+5)}, + \ldots + \mathrm{N}^{(N_{0}+N)}, + \mathrm{C}_{\alpha}^{(N_{0}+N)} + \right) + } \; = \; \; \; \; \\ \; \; \; \; { + \nonumber + \frac{1}{2(N-2)} + \sum_{n=N_{0}}^{N_{0}+N-2} + \mathrm{angf}\left( + \mathrm{C}_{\alpha}^{(n)}, + \mathrm{C}_{\alpha}^{(n+1)}, + \mathrm{C}_{\alpha}^{(n+2)}\right) + } \; + \; { + \frac{1}{2(N-4)} + \sum_{n=N_{0}}^{N_{0}+N-4} + \mathrm{hbf}\left( + \mathrm{O}^{(n)}, + \mathrm{N}^{(n+4)}\right) \mathrm{,} + } \\ +\end{eqnarray} +where the score function for the$\mathrm{C}_{\alpha} -
+\mathrm{C}_{\alpha} - \mathrm{C}_{\alpha}$angle is defined as: - \label{eq:cvc_coordNum} - C (\mathtt{group1}, \mathtt{group2}) \; = \; - \sum_{i\in\mathtt{group1}}\sum_{j\in\mathtt{group2}} { - \frac{1 - (|\mathbf{x}_{i}-\mathbf{x}_{j}|/d_{0})^{n}}{ - 1 - (|\mathbf{x}_{i}-\mathbf{x}_{j}|/d_{0})^{m} } + \label{eq:colvars_alpha_Calpha} + { + \mathrm{angf}\left( + \mathrm{C}_{\alpha}^{(n)}, + \mathrm{C}_{\alpha}^{(n+1)}, + \mathrm{C}_{\alpha}^{(n+2)}\right) + } \; = \; { + \frac{1 - \left(\theta( + \mathrm{C}_{\alpha}^{(n)}, + \mathrm{C}_{\alpha}^{(n+1)}, + \mathrm{C}_{\alpha}^{(n+2)}) - + \theta_{0}\right)^{2} / + \left(\Delta\theta_{\mathrm{tol}}\right)^{2}}{ + 1 - \left(\theta( + \mathrm{C}_{\alpha}^{(n)}, + \mathrm{C}_{\alpha}^{(n+1)}, + \mathrm{C}_{\alpha}^{(n+2)}) - + \theta_{0}\right)^{4} / + \left(\Delta\theta_{\mathrm{tol}}\right)^{4}} \mathrm{,} } +and the score function for the$\mathrm{O}^{(n)} \leftrightarrow
+\mathrm{N}^{(n+4)}$hydrogen bond is defined through a \texttt{hBond} +colvar component on the same atoms. \begin{cvcoptions} \item % - \labelkey{colvar|coordNum|group1} - \dupkey{group1}{\texttt{coordNum}}{colvar|distance|group1}{\texttt{distance} component} + \labelkey{colvar|alpha|residueRange} + \key + {residueRange}{% + \texttt{alpha}}{% + Potential$\alpha$-helical residues}{% + $<$Initial residue number$>$-$<$Final residue number$>$''}{% + This option specifies the range of residues on which this + component should be defined. The Colvars module looks for the + atoms within these residues named \texttt{CA}'', \texttt{N}'' + and \texttt{O}'', and raises an error if any of those atoms is + not found.} \item % - \labelkey{colvar|coordNum|group2} - \simkey{group2}{\texttt{coordNum}}{group1} + \labelkey{colvar|alpha|psfSegID} + \key + {psfSegID}{% + \texttt{alpha}}{% + PSF segment identifier}{% + string (max 4 characters)}{% + This option sets the PSF segment identifier for the residues + specified in \texttt{residueRange}. This option is only required + when PSF topologies are used.} + \item % - \labelkey{colvar|coordNum|cutoff} + \labelkey{colvar|alpha|hBondCoeff} \keydef - {cutoff}{% - \texttt{coordNum}}{% - Interaction'' distance (\AA)}{% - positive decimal}{% - 4.0}{% - This number defines the switching distance to define an - interatomic contact: for$d \ll d_0$, the switching function -$(1-(d/d_0)^{n})/(1-(d/d_0)^{m})$is close to 1, at$d = d_0$it - has a value of$n/m$($1/2$with the default$n$and$m$), and at -$d \gg d_0$it goes to zero approximately like$d^{m-n}$. Hence, - for a proper behavior,$m$must be larger than$n$.} + {hBondCoeff}{% + \texttt{alpha}}{% + Coefficient for the hydrogen bond term}{% + positive between 0 and 1}{% + 0.5}{% + This number specifies the contribution to the total value from the + hydrogen bond terms. 0 disables the hydrogen bond terms, 1 + disables the angle terms.} \item % - \labelkey{colvar|coordNum|cutoff3} + \labelkey{colvar|alpha|angleRef} \keydef - {cutoff3}{% - \texttt{coordNum}}{% - Reference distance vector (\AA)}{% - \texttt{(x, y, z)}'' triplet of positive decimals}{% - \texttt{(4.0, 4.0, 4.0)}}{% - The three components of this vector define three different cutoffs -$d_{0}$for each direction. This option is mutually exclusive with - \texttt{cutoff}.} + {angleRef}{% + \texttt{alpha}}{% + Reference$\mathrm{C}_{\alpha} -
+    \mathrm{C}_{\alpha} - \mathrm{C}_{\alpha}$angle}{% + positive decimal}{% + 88$^{\circ}$}{% + This option sets the reference angle used in the score function + (\ref{eq:colvars_alpha_Calpha}).} \item % - \labelkey{colvar|coordNum|expNumer} + \labelkey{colvar|alpha|angleTol} \keydef - {expNumer}{% - \texttt{coordNum}}{% - Numerator exponent}{% - positive even integer}{% - 6}{% - This number defines the$n$exponent for the switching function.} + {angleTol}{% + \texttt{alpha}}{% + Tolerance in the$\mathrm{C}_{\alpha} -
+    \mathrm{C}_{\alpha} - \mathrm{C}_{\alpha}$angle}{% + positive decimal}{% + 15$^{\circ}$}{% + This option sets the angle tolerance used in the score function + (\ref{eq:colvars_alpha_Calpha}).} \item % - \labelkey{colvar|coordNum|expDenom} + \labelkey{colvar|alpha|hBondCutoff} \keydef - {expDenom}{% - \texttt{coordNum}}{% - Denominator exponent}{% - positive even integer}{% - 12}{% - This number defines the$m$exponent for the switching function.} + {hBondCutoff}{% + \texttt{alpha}}{% + Hydrogen bond cutoff}{% + positive decimal}{% + 3.3~\AA{}}{% + Equivalent to the \texttt{cutoff} option in the \texttt{hBond} + component.} \item % - \labelkey{colvar|coordNum|group2CenterOnly} + \labelkey{colvar|alpha|hBondExpNumer} \keydef - {group2CenterOnly}{% - \texttt{coordNum}}{% - Use only \texttt{group2}'s center of - mass}{% - boolean}{% - \texttt{off}}{% - If this option is \texttt{on}, only contacts between each atoms in \texttt{group1} and the center of mass of \texttt{group2} are calculated (by default, the sum extends over all pairs of atoms in \texttt{group1} and \texttt{group2}). -If \texttt{group2} is a \texttt{dummyAtom}, this option is set to \texttt{yes} by default. -} - -\item % - \labelkey{colvar|coordNum|tolerance} - \keydef - {tolerance}{% - \texttt{coordNum}}{% - Pairlist control}{% - decimal}{% - 0.0}{This controls the pairlist feature, dictating the minimum value for each summation element in Eq.~\ref{eq:cvc_coordNum} such that the pair that contributed the summation element is included in subsequent simulation timesteps until the next pairlist recalculation. For most applications, this value should be small (eg. 0.001) to avoid missing important contributions to the overall sum. Higher values will improve performance, although values above 1 will exclude all possible pair interactions. Similarly, values below 0 will never exclude a pair from consideration. -} + {hBondExpNumer}{% + \texttt{alpha}}{% + Hydrogen bond numerator exponent}{% + positive integer}{% + 6}{% + Equivalent to the \texttt{expNumer} option in the \texttt{hBond} + component.} \item % - \labelkey{colvar|coordNum|pairListFrequency} - \keydef - {pairListFrequency}{% - \texttt{coordNum}}{% - Pairlist regeneration frequency}{% + \labelkey{colvar|alpha|hBondExpDenom} + \keydef + {hBondExpDenom}{% + \texttt{alpha}}{% + Hydrogen bond denominator exponent}{% positive integer}{% - 100}{This controls the pairlist feature, dictating how many steps are taken between regenerating pairlists if the tolerance is greater than 0. At this interval, the colvar defined will be exact, as though it was the all-to-all pair summation defined in Eq.~\ref{eq:cvc_coordNum}. All other steps will report only values and gradients from pairs in the pairlist. -} + 8}{% + Equivalent to the \texttt{expDenom} option in the \texttt{hBond} + component.} + \end{cvcoptions} -This component returns a dimensionless number, which ranges from -approximately 0 (all interatomic distances are much larger than the -cutoff) to$N_{\mathtt{group1}} \times N_{\mathtt{group2}}$(all distances -are less than the cutoff), or$N_{\mathtt{group1}}$if -\texttt{group2CenterOnly} is used. For performance reasons, at least -one of \texttt{group1} and \texttt{group2} should be of limited size or \texttt{group2CenterOnly} should be used: the cost of the loop over all pairs grows as$N_{\mathtt{group1}} \times N_{\mathtt{group2}}$. +This component returns positive values, always comprised between 0 +(lowest$\alpha$-helical score) and 1 (highest$\alpha$-helical +score). +\cvsubsubsec{\texttt{dihedralPC}: protein dihedral pricipal component}{sec:cvc_dihedralPC} +\labelkey{colvar|dihedralPC} -\cvsubsubsec{\texttt{selfCoordNum}: coordination number between atoms within a group.}{sec:cvc_selfCoordNum} -The \texttt{selfCoordNum \{...\}} block defines -a coordination number similarly to the component \texttt{coordNum}, -but the function is summed over atom pairs within \texttt{group1}: +The block \texttt{dihedralPC~\{...\}} defines the +parameters to calculate the projection of backbone dihedral angles within +a protein segment onto a \emph{dihedral principal component}, following +the formalism of dihedral principal component analysis (dPCA) proposed by +Mu et al.\cite{Mu2005} and documented in detail by Altis et +al.\cite{Altis2007}. +Given a peptide or protein segment of$N$residues, each with Ramachandran +angles$\phi_i$and$\psi_i$, dPCA rests on a variance/covariance analysis +of the$4(N-1)$variables$\cos(\psi_1), \sin(\psi_1), \cos(\phi_2), \sin(\phi_2)
+\cdots \cos(\phi_N), \sin(\phi_N)$. Note that angles$\phi_1$and$\psi_N$+have little impact on chain conformation, and are therefore discarded, +following the implementation of dPCA in the analysis software Carma.\cite{Glykos2006} + +For a given principal component (eigenvector) of coefficients +$(k_i)_{1 \leq i \leq 4(N-1)}$, +the projection of the current backbone conformation is: - \label{eq:cvc_selfCoordNum} - C (\mathtt{group1}) \; = \; - \sum_{i\in\mathtt{group1}}\sum_{j > i} { - \frac{1 - (|\mathbf{x}_{i}-\mathbf{x}_{j}|/d_{0})^{n}}{ - 1 - (|\mathbf{x}_{i}-\mathbf{x}_{j}|/d_{0})^{m} } - } +\xi = \sum_{n=1}^{N-1} k_{4n-3} \cos(\psi_n) + k_{4n-2} \sin (\psi_n) ++ k_{4n-1} \cos (\phi_{n+1}) + k_{4n} \sin(\phi_{n+1}) -The keywords accepted by \texttt{selfCoordNum} are a subset of -those accepted by \texttt{coordNum}, namely \texttt{group1} -(here defining \emph{all} of the atoms to be considered), -\texttt{cutoff}, \texttt{expNumer}, and \texttt{expDenom}. -\begin{cvcoptions} -\item % - \dupkey{group1}{\texttt{selfCoordNum}}{colvar|coordNum|group1}{\texttt{coordNum} component} -\item % - \dupkey{cutoff}{\texttt{selfCoordNum}}{colvar|coordNum|cutoff}{\texttt{coordNum} component} -\item % - \dupkey{cutoff3}{\texttt{selfCoordNum}}{colvar|coordNum|cutoff3}{\texttt{coordNum} component} -\item % - \dupkey{expNumer}{\texttt{selfCoordNum}}{colvar|coordNum|expNumer}{\texttt{coordNum} component} -\item % - \dupkey{expDenom}{\texttt{selfCoordNum}}{colvar|coordNum|expDenom}{\texttt{coordNum} component} -\item % - \dupkey{tolerance}{\texttt{selfCoordNum}}{colvar|coordNum|tolerance}{\texttt{coordNum} component} -\item % - \dupkey{pairListFrequency}{\texttt{selfCoordNum}}{colvar|coordNum|pairListFrequency}{\texttt{coordNum} component} -\end{cvcoptions} +\texttt{dihedralPC} expects the same parameters as the \texttt{alpha} +component for defining the relevant residues (\texttt{residueRange} +and \texttt{psfSegID}) in addition to the following: -This component returns a dimensionless number, which ranges from -approximately 0 (all interatomic distances much larger than the -cutoff) to$N_{\mathtt{group1}} \times (N_{\mathtt{group1}} - 1) / 2$(all -distances within the cutoff). For performance reasons, -\texttt{group1} should be of limited size, because the cost of the -loop over all pairs grows as$N_{\mathtt{group1}}^2$. +\begin{cvcoptions} +\item % + \dupkey{residueRange}{\texttt{dihedralPC}}{colvar|alpha|residueRange}{\texttt{alpha} component} +\item % + \dupkey{psfSegID}{\texttt{dihedralPC}}{colvar|alpha|psfSegID}{\texttt{alpha} component} -\cvsubsubsec{\texttt{hBond}: hydrogen bond between two atoms.}{sec:cvc_hBond} -The \texttt{hBond \{...\}} block defines a hydrogen -bond, implemented as a coordination number (eq.~\ref{eq:cvc_coordNum}) -between the donor and the acceptor atoms. Therefore, it accepts the -same options \texttt{cutoff} (with a different default value of -3.3~\AA{}), \texttt{expNumer} (with a default value of 6) and -\texttt{expDenom} (with a default value of 8). Unlike -\texttt{coordNum}, it requires two atom numbers, \texttt{acceptor} and -\texttt{donor}, to be defined. It returns an adimensional number, -with values between 0 (acceptor and donor far outside the cutoff -distance) and 1 (acceptor and donor much closer than the cutoff). +\item % + \key + {vectorFile}{% + \texttt{dihedralPC}}{% + File containing dihedral PCA eigenvector(s)}{% + file name}{% + A text file containing the coefficients of dihedral PCA eigenvectors on the + cosine and sine coordinates. The vectors should be arranged in columns, + as in the files output by Carma.\cite{Glykos2006}} -\begin{cvcoptions} \item % \key - {acceptor}{% - \texttt{hBond}}{% - Number of the acceptor atom}{% + {vectorNumber}{% + \texttt{dihedralPC}}{% + File containing dihedralPCA eigenvector(s)}{% positive integer}{% - Number that uses the same convention as \texttt{atomNumbers}.} -\item % - \simkey{donor}{\texttt{hBond}}{acceptor} -\item % - \dupkey{cutoff}{\texttt{hBond}}{colvar|coordNum|cutoff}{\texttt{coordNum} component}\\ - \textbf{Note:} default value is 3.3~\AA. -\item % - \dupkey{expNumer}{\texttt{hBond}}{colvar|coordNum|expNumer}{\texttt{coordNum} component}\\ - \textbf{Note:} default value is 6. -\item % - \dupkey{expDenom}{\texttt{hBond}}{colvar|coordNum|expDenom}{\texttt{coordNum} component}\\ - \textbf{Note:} default value is 8. + Number of the eigenvector to be used for this component.} \end{cvcoptions} +} % end of \cvnamebasedonly -\cvsubsubsec{\texttt{rmsd}: root mean square displacement (RMSD) from reference positions.}{sec:cvc_rmsd} -The block \texttt{rmsd~\{...\}} defines the root mean square replacement -(RMSD) of a group of atoms with respect to a reference structure. For -each set of coordinates$\{ \mathbf{x}_1(t), \mathbf{x}_2(t), \ldots
-\mathbf{x}_N(t) \}$, the colvar component \texttt{rmsd} calculates the -optimal rotation -$U^{\{\mathbf{x}_{i}(t)\}\rightarrow\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}}$-that best superimposes the coordinates$\{\mathbf{x}_{i}(t)\}$onto a -set of reference coordinates$\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}$. -Both the current and the reference coordinates are centered on their -centers of geometry,$\mathbf{x}_{\mathrm{cog}}(t)$and -$\mathbf{x}_{\mathrm{cog}}^{\mathrm{(ref)}}$. The root mean square -displacement is then defined as: - - \label{eq:cvc_rmsd} - { \mathrm{RMSD}(\{\mathbf{x}_{i}(t)\}, - \{\mathbf{x}_{i}^{\mathrm{(ref)}}\}) } \; = \; \sqrt{ - \frac{1}{N} \sum_{i=1}^{N} \left| - U - \left(\mathbf{x}_{i}(t) - \mathbf{x}_{\mathrm{cog}}(t)\right) - - \left(\mathbf{x}_{i}^{\mathrm{(ref)}} - - \mathbf{x}_{\mathrm{cog}}^{\mathrm{(ref)}} \right) \right|^{2} } - -The optimal rotation -$U^{\{\mathbf{x}_{i}(t)\}\rightarrow\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}}$-is calculated within the formalism developed in -reference~\cite{Coutsias2004}, which guarantees a continuous -dependence of -$U^{\{\mathbf{x}_{i}(t)\}\rightarrow\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}}$-with respect to$\{\mathbf{x}_{i}(t)\}$. - -\begin{cvcoptions} -\item % - \labelkey{colvar|rmsd|atoms} - \key - {atoms}{% - \texttt{rmsd}}{% - Atom group}{% - \texttt{atoms~\{...\}} block}{% - Defines the group of atoms of which the RMSD should be calculated. - Optimal fit options (such as \texttt{refPositions} and - \texttt{rotateReference}) should typically NOT be set within this - block. Exceptions to this rule are the special cases discussed in - the \emph{Advanced usage} paragraph below. - } +\cvsubsec{Shared keywords for all components}{sec:cvc_common} +The following options can be used for any of the above colvar components in order to obtain a polynomial combination\cvscriptonly{ or any user-supplied function provided by \refkey{scriptedFunction}{sec:cvc_superp}}. +\begin{itemize} \item % - \labelkey{colvar|rmsd|refPositions} - \key - {refPositions}{% - \texttt{rmsd}}{% - Reference coordinates}{% - space-separated list of \texttt{(x, y, z)} triplets}{% - This option (mutually exclusive with \texttt{refPositionsFile}) sets the reference coordinates for RMSD calculation, and uses these to compute the roto-translational fit. - It is functionally equivalent to the option \refkey{refPositions}{atom-group|refPositions} in the atom group definition, which also supports more advanced fitting options. - } + \keydef + {name}{% + any component}{% + Name of this component}{% + string}{% + type of component + numeric id}{% + The name is an unique case-sensitive string which allows the + Colvars module to identify this component. This is useful, for example, + when combining multiple components via a \texttt{scriptedFunction}. + \cvnamdonly{It also defines the variable name representing the component's value in a \refkey{customFunction}{colvar|customFunction} expression.}} \item % - \labelkey{colvar|rmsd|refPositionsFile} - \key - {refPositionsFile}{% - \texttt{rmsd}}{% - Reference coordinates file}{% - UNIX filename}{% - This option (mutually exclusive with \texttt{refPositions}) sets the reference coordinates for RMSD calculation, and uses these to compute the roto-translational fit. - It is functionally equivalent to the option \refkey{refPositionsFile}{atom-group|refPositionsFile} in the atom group definition, which also supports more advanced fitting options. - } - -\cvnamebasedonly{ -\item % - \labelkey{colvar|rmsd|refPositionsCol} - \key - {refPositionsCol}{% - \texttt{rmsd}}{% - PDB column containing atom flags}{% - \texttt{O}, \texttt{B}, \texttt{X}, \texttt{Y}, or \texttt{Z}}{% - If \texttt{refPositionsFile} is a PDB file that contains all the atoms in the topology, this option may be provided to set which PDB field is used to flag the reference coordinates for \texttt{atoms}. + \keydef + {scalable}{% + any component}{% + Attempt to calculate this component in parallel?}{% + boolean}{% + \texttt{on}, if available}{% + If set to \texttt{on} (default), the Colvars module will attempt to calculate this component in parallel to reduce overhead. + Whether this option is available depends on the type of component: currently supported are \texttt{distance}, \texttt{distanceZ}, \texttt{distanceXY}, \texttt{distanceVec}, \texttt{distanceDir}, \texttt{angle} and \texttt{dihedral}. + This flag influences computational cost, but does not affect numerical results: therefore, it should only be turned off for debugging or testing purposes. } +\end{itemize} + +\cvsubsec{Periodic components}{sec:cvc_periodic} +The following components returns +real numbers that lie in a periodic interval: +\begin{itemize} +\item \texttt{dihedral}: torsional angle between four groups; +\item \texttt{spinAngle}: angle of rotation around a predefined axis + in the best-fit from a set of reference coordinates. +\end{itemize} +In certain conditions, \texttt{distanceZ} can also be periodic, namely +when periodic boundary conditions (PBCs) are defined in the simulation +and \texttt{distanceZ}'s axis is parallel to a unit cell vector. + +In addition, a custom \cvscriptonly{or scripted} scalar colvar may be periodic +depending on its user-defined expression. It will only be treated as such by +the Colvars module if the period is specified using the \texttt{period} keyword, +while \texttt{wrapAround} is optional. + +The following keywords can be used within periodic components\cvleptononly{, or within custom variables (\ref{sec:colvar_custom_function})}\cvscriptonly{, or wthin scripted variables \ref{sec:colvar_scripted}}). + +\begin{itemize} \item % - \labelkey{colvar|rmsd|refPositionsColValue} - \key - {refPositionsColValue}{% - \texttt{rmsd}}{% - Atom selection flag in the PDB column}{% + \keydef + {period}{% + \texttt{distanceZ}, custom colvars}{% + Period of the component}{% positive decimal}{% - If defined, this value identifies in the PDB column - \texttt{refPositionsCol} of the file \texttt{refPositionsFile} - which atom positions are to be read. Otherwise, all positions - with a non-zero value are read. - } -} -\end{cvcoptions} -This component returns a positive real number (in \AA). + 0.0}{% + Setting this number enables the treatment of \texttt{distanceZ} as + a periodic component: by default, \texttt{distanceZ} is not + considered periodic. The keyword is supported, but irrelevant + within \texttt{dihedral} or \texttt{spinAngle}, because their + period is always 360~degrees.} -\cvsubsubsec{Advanced usage of the \texttt{rmsd} component.}{sec:cvc_rmsd_advanced} -In the standard usage as described above, the \texttt{rmsd} component -calculates a minimum RMSD, that is, current coordinates are optimally -fitted onto the same reference coordinates that are used to -compute the RMSD value. The fit itself is handled by the atom group -object, whose parameters are automatically set by the \texttt{rmsd} -component. -For very specific applications, however, it may be -useful to control the fitting process separately from the definition -of the reference coordinates, to evaluate various types of -non-minimal RMSD values. This can be achieved by setting the -related options (\texttt{refPositions}, etc.) explicitly in the -atom group block. This allows for the following non-standard cases: +\item % + \keydef + {wrapAround}{% + \texttt{distanceZ}, \texttt{dihedral}, \texttt{spinAngle}, custom colvars}{% + Center of the wrapping interval for periodic variables}{% + decimal}{% + 0.0}{% + By default, values of the periodic components are centered around zero, ranging from$-P/2$to$P/2$, where$P$is the period. + Setting this number centers the interval around this value. + This can be useful for convenience of output, or to set the walls for a \texttt{harmonicWalls} in an order that would not otherwise be allowed.} +\end{itemize} -\begin{enumerate} -\item applying the optimal translation, but no rotation -(\texttt{rotateReference off}), to bias or restrain the shape and -orientation, but not the position of the atom group; -\item applying the optimal rotation, but no translation -(\texttt{translateReference off}), to bias or restrain the shape and -position, but not the orientation of the atom group; -\item disabling the application of optimal roto-translations, which -lets the RMSD component decribe the deviation of atoms -from fixed positions in the laboratory frame: this allows for custom -positional restraints within the Colvars module; -\item fitting the atomic positions to different reference coordinates -than those used in the RMSD calculation itself; -\item applying the optimal rotation and/or translation from a separate -atom group, defined through \texttt{fittingGroup}: the RMSD then -reflects the deviation from reference coordinates in a separate, moving -reference frame. -\end{enumerate} +Internally, all differences between two values of a periodic colvar +follow the minimum image convention: they are calculated based on +the two periodic images that are closest to each other. -\cvscriptonly{ -\cvsubsubsec{Path collective variables}{sec:pathcv} +\emph{Note: linear or polynomial combinations of periodic components (see \ref{sec:cvc_superp}) may become meaningless when components cross the periodic boundary. Use such combinations carefully: estimate the range of possible values of each component in a given simulation, and make use of \texttt{wrapAround} to limit this problem whenever possible.} -An application of the \texttt{rmsd} component is "path collective variables",\cite{Branduardi2007} -which are implemented as Tcl-scripted combinations or RMSDs. -The implementation is available as file \texttt{colvartools/pathCV.tcl}, and -an example is provided in file \texttt{examples/10\_pathCV.namd} of the Colvars public repository.} -\cvsubsubsec{\texttt{eigenvector}: projection of the atomic coordinates on a vector.}{sec:cvc_eigenvector} -The block \texttt{eigenvector~\{...\}} defines the projection of the coordinates -of a group of atoms (or more precisely, their deviations from the -reference coordinates) onto a vector in$\mathbb{R}^{3n}$, where$n$is the -number of atoms in the group. The computed quantity is the -total projection: +\cvsubsec{Non-scalar components}{sec:cvc_non_scalar} + +When one of the following components are used, the defined colvar returns a value that is not a scalar number: +\begin{itemize} +\item \texttt{distanceVec}: 3-dimensional vector of the distance + between two groups; +\item \texttt{distanceDir}: 3-dimensional unit vector of the distance + between two groups; +\item \texttt{orientation}: 4-dimensional unit quaternion representing + the best-fit rotation from a set of reference coordinates. +\end{itemize} +The distance between two 3-dimensional unit vectors is computed as the +angle between them. The distance between two quaternions is computed +as the angle between the two 4-dimensional unit vectors: because the +orientation represented by$\mathsf{q}$is the same as the one +represented by$-\mathsf{q}$, distances between two quaternions are +computed considering the closest of the two symmetric images. + +Non-scalar components carry the following restrictions: +\begin{itemize} +\item Calculation of total forces (\texttt{outputTotalForce} option) + is currently not implemented. +\item Each colvar can only contain one non-scalar component. +\item Binning on a grid (\texttt{abf}, \texttt{histogram} and + \texttt{metadynamics} with \texttt{useGrids} enabled) is currently + not implemented for colvars based on such components. +\end{itemize} + +\emph{Note: while these restrictions apply to individual colvars based + on non-scalar components, no limit is set to the number of scalar + colvars. To compute multi-dimensional histograms and PMFs, use sets + of scalar colvars of arbitrary size.} + + +\cvsubsubsec{Calculating total forces}{sec:cvc_sys_forces} +In addition to the restrictions due to the type of value computed (scalar or non-scalar), +a final restriction can arise when calculating total force +(\texttt{outputTotalForce} option or application of a \texttt{abf} +bias). total forces are available currently only for the following +components: \texttt{distance}, \texttt{distanceZ}, +\texttt{distanceXY}, \texttt{angle}, \texttt{dihedral}, \texttt{rmsd}, +\texttt{eigenvector} and \texttt{gyration}. + + + +\cvsubsec{Linear and polynomial combinations of components}{sec:cvc_superp} + +To extend the set of possible definitions of colvars$\xi(\mathbf{r})$, multiple components +$q_i(\mathbf{r})$can be summed with the formula: - \label{eq:cvc_eigenvector} - { p(\{\mathbf{x}_{i}(t)\}, - \{\mathbf{x}_{i}^{\mathrm{(ref)}}\}) } \; = \; { - \sum_{i=1}^{n} \mathbf{v}_{i} \cdot - \left(U(\mathbf{x}_{i}(t) - \mathbf{x}_{\mathrm{cog}}(t)) - - (\mathbf{x}_{i}^{\mathrm{(ref)}} - - \mathbf{x}_{\mathrm{cog}}^{\mathrm{(ref)}}) \right)\mathrm{,} } + \label{eq:colvar_combination} + \xi(\mathbf{r}) = \sum_i c_i [q_i(\mathbf{r})]^{n_i} -where, as in the \texttt{rmsd} component,$U$is the optimal rotation -matrix,$\mathbf{x}_{\mathrm{cog}}(t)$and -$\mathbf{x}_{\mathrm{cog}}^{\mathrm{(ref)}}$are the centers of -geometry of the current and reference positions respectively, and -$\mathbf{v}_{i}$are the components of the vector for each atom. -Example choices for$(\mathbf{v}_{i})$are an eigenvector -of the covariance matrix (essential mode), or a normal -mode of the system. It is assumed that$\sum_{i}\mathbf{v}_{i} = 0$: -otherwise, the Colvars module centers the$\mathbf{v}_{i}$-automatically when reading them from the configuration. +where each component appears with a unique coefficient$c_i$(1.0 by +default) the positive integer exponent$n_i$(1 by default). -\begin{cvcoptions} -\item % - \dupkey{atoms}{\texttt{eigenvector}}{colvar|rmsd|atoms}{\texttt{rmsd} component} +Any set of components can be combined within a colvar, provided that +they return the same type of values (scalar, unit vector, vector, or +quaternion). By default, the colvar is the sum of its components. +Linear or polynomial combinations (following +equation~(\ref{eq:colvar_combination})) can be obtained by setting the +following parameters, which are common to all components: +\begin{itemize} \item % - \dupkey{refPositions}{\texttt{eigenvector}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} + \keydef + {componentCoeff}{% + any component}{% + Coefficient of this component in the colvar}{% + decimal}{% + \texttt{1.0}}{% + Defines the coefficient by which this component is multiplied + (after being raised to \texttt{componentExp}) before being added + to the sum.} + \item % - \dupkey{refPositionsFile}{\texttt{eigenvector}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} -\cvnamebasedonly{ + \keydef + {componentExp}{% + any component}{% + Exponent of this component in the colvar}{% + integer}{% + \texttt{1}}{% + Defines the power at which the value of this component is raised + before being added to the sum. When this exponent is + different than 1 (non-linear sum), total forces and the Jacobian + force are not available, making the colvar unsuitable for ABF calculations.} +\end{itemize} + +\textbf{Example:} To define the \emph{average} of a colvar across +different parts of the system, simply define within the same colvar +block a series of components of the same type (applied to different +atom groups), and assign to each component a \texttt{componentCoeff} +of$1/N$. + + +\cvleptononly{ +\cvsubsec{Custom functions}{sec:colvar_custom_function} + +Collective variables may be defined by specifying a custom function as an analytical +expression such as \texttt{cos(x) + y\^{}2}. +The expression is parsed by the Lepton expression parser (written by Peter Eastman), +which produces efficient evaluation routines for the function itself as well as its derivatives. +The expression may use the collective variable components as variables, refered to as their \texttt{name} string. +Scalar elements of vector components may be accessed by appending a 1-based index to their \texttt{name}. +When implementing generic functions of Cartesian coordinates rather +than functions of existing components, the \texttt{cartesian} component +may be particularly useful. +A scalar-valued custom variable may be manually defined as periodic by providing +the keyword \texttt{period}, and the optional keyword \texttt{wrapAround}, with the +same meaning as in periodic components (see \ref{sec:cvc_periodic} for details). +A vector variable may be defined by specifying the \texttt{customFunction} parameter several times: each expression defines one scalar element of the vector colvar. +This is illustrated in the example below. + +\bigskip +{% +% verbatim can't appear within commands +\noindent\ttfamily colvar \{\\ +\-~~name custom\\ +\\ +\-~~\# A 2-dimensional vector function of a scalar x and a 3-vector r\\ +\-~~customFunction cos(x) * (r1 + r2 + r3)\\ +\-~~customFunction sqrt(r1 * r2)\\ +\\ +\-~~distance \{\\ +\-~~~~name x\\ +\-~~~~group1 \{ atomNumbers 1 \}\\ +\-~~~~group2 \{ atomNumbers 50 \}\\ +\-~~\}\\ +\-~~distanceVec \{\\ +\-~~~~name r\\ +\-~~~~group1 \{ atomNumbers 10 11 12 \}\\ +\-~~~~group2 \{ atomNumbers 20 21 22 \}\\ +\-~~\}\\ +\}} + +\begin{itemize} \item % - \dupkey{refPositionsCol}{\texttt{eigenvector}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} + \labelkey{colvar|customFunction} + \key + {customFunction}{% + \texttt{colvar}}{% + Compute colvar as a custom function of its components}{% + string}{% + Defines the colvar as a scalar expression of its colvar components. Multiple mentions can + be used to define a vector variable (as in the example above).} + \item % - \dupkey{refPositionsColValue}{\texttt{eigenvector}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} + \keydef + {customFunctionType}{% + \texttt{colvar}}{% + Type of value returned by the scripted colvar}{% + string}{% + \texttt{scalar}}{% + With this flag, the user may specify whether the + colvar is a scalar or one of the following vector types: \texttt{vector3} + (a 3D vector), \texttt{unit\_vector3} (a normalized 3D vector), or + \texttt{unit\_quaternion} (a normalized quaternion), or \texttt{vector}. + Note that the scalar and vector cases are not necessary, as they are detected automatically.} +\end{itemize} } -\item % - \key - {vector}{% - \texttt{eigenvector}}{% - Vector components}{% - space-separated list of \texttt{(x, y, z)} triplets}{% - This option (mutually exclusive with \texttt{vectorFile}) sets the values of the vector components.} -\item % - \key - {vectorFile}{% - \texttt{eigenvector}}{% - file containing vector components}{% - UNIX filename}{% - This option (mutually exclusive with \texttt{vector}) sets the name of a coordinate file containing the vector components; the file is read according to the same format used for \texttt{refPositionsFile}. - \cvnamebasedonly{For a PDB file specifically, the components are read from the X, Y and Z fields. - \textbf{Note:} \emph{The PDB file has limited precision and fixed-point numbers: in some cases, the vector components may not be accurately represented; a XYZ file should be used instead, containing floating-point numbers.}} +\cvscriptonly{ +\cvsubsec{Scripted functions}{sec:colvar_scripted} +When scripting is supported\cvnamdonly{ (default in NAMD)}\cvvmdonly{ (default in VMD)}, +a colvar may be defined as a scripted function of its components, +rather than a linear or polynomial combination. +When implementing generic functions of Cartesian coordinates rather +than functions of existing components, the \texttt{cartesian} component +may be particularly useful. +A scalar-valued scripted variable may be manually defined as periodic by providing +the keyword \texttt{period}, and the optional keyword \texttt{wrapAround}, with the +same meaning as in periodic components (see \ref{sec:cvc_periodic} for details). + +An example of elaborate scripted colvar is given in example 10, in the +form of path-based collective variables as defined by Branduardi et al\cite{Branduardi2007} +(\ref{sec:pathcv}). + +\begin{itemize} + \item % + \labelkey{colvar|scriptedFunction} + \key + {scriptedFunction}{% + \texttt{colvar}}{% + Compute colvar as a scripted function of its components}{% + string}{% + If this option is specified, the colvar will be computed as a + scripted function of the values of its components. + To that effect, the user should define two Tcl procedures: + \texttt{calc\_$<$scriptedFunction$>$} and \texttt{calc\_$<$scriptedFunction$>$\_gradient}, + both accepting as many parameters as the colvar has components. + Values of the components will be passed to those procedures in the + order defined by their sorted \texttt{name} strings. Note that if all + components are of the same type, their default names are sorted in the + order in which they are defined, so that names need only be specified + for combinations of components of different types. + \texttt{calc\_$<$scriptedFunction$>$} should return one value of + type$<$scriptedFunctionType$>$, corresponding to the colvar value. + \texttt{calc\_$<$scriptedFunction$>$\_gradient} should return a Tcl list + containing the derivatives of the function with respect to each + component. + If both the function and some of the components are vectors, the gradient + is really a Jacobian matrix that should be passed as a linear vector in + row-major order, i.e. for a function$f_i(x_j)$:$\nabla_x f_1 \nabla_x f_2 \cdots$. } -\cvnamebasedonly{ -\item % - \key - {vectorCol}{% - \texttt{eigenvector}}{% - PDB column used to flag participating atoms}{% - \texttt{O} or \texttt{B}}{% - Analogous to \texttt{atomsCol}.} +\item% + \keydef + {scriptedFunctionType}{% + \texttt{colvar}}{% + Type of value returned by the scripted colvar}{% + string}{% + \texttt{scalar}}{% + If a colvar is defined as a scripted function, its type is not constrained by + the types of its components. With this flag, the user may specify whether the + colvar is a scalar or one of the following vector types: \texttt{vector3} + (a 3D vector), \texttt{unit\_vector3} (a normalized 3D vector), or + \texttt{unit\_quaternion} (a normalized quaternion), or \texttt{vector} + (a vector whose size is specified by \texttt{scriptedFunctionVectorSize}). + Non-scalar values should be passed as space-separated lists.} -\item % + \item % \key - {vectorColValue}{% - \texttt{eigenvector}}{% - Value used to flag participating atoms in the PDB file}{% - positive decimal}{% - Analogous to \texttt{atomsColValue}.} -} + {scriptedFunctionVectorSize}{% + \texttt{colvar}}{% + Dimension of the vector value of a scripted colvar}{% + positive integer}{% + This parameter is only valid when \texttt{scriptedFunctionType} is + set to \texttt{vector}. It defines the vector length of the colvar value + returned by the function.} +\end{itemize} +} % \cvscriptonly -\item % - \keydef - {differenceVector}{% - \texttt{eigenvector}}{% - The$3n$-dimensional vector is the difference between \texttt{vector} and \texttt{refPositions}}{% - boolean}{% - \texttt{off}}{% - If this option is \texttt{on}, the numbers provided by \texttt{vector}\cvnamebasedonly{ or \texttt{vectorFile}} are interpreted as another set of positions,$\mathbf{x}_{i}'$: the vector$\mathbf{v}_{i}$is then defined as$\mathbf{v}_{i} = \left(\mathbf{x}_{i}' - \mathbf{x}_{i}^{\mathrm{(ref)}}\right)$. -This allows to conveniently define a colvar$\xi$as a projection on the linear transformation between two sets of positions, A'' and B''. -For convenience, the vector is also normalized so that$\xi = 0$when the atoms are at the set of positions A'' and$\xi = 1$at the set of positions B''. -} -\end{cvcoptions} -This component returns a number (in \AA), whose value ranges between -the smallest and largest absolute positions in the unit cell during -the simulations (see also \texttt{distanceZ}). Due to the -normalization in eq.~\ref{eq:cvc_eigenvector}, this range does not -depend on the number of atoms involved. -\cvsubsubsec{\texttt{gyration}: radius of gyration of a group of atoms.}{sec:cvc_gyration} -The block \texttt{gyration~\{...\}} defines the -parameters for calculating the radius of gyration of a group of atomic -positions$\{ \mathbf{x}_1(t), \mathbf{x}_2(t), \ldots \mathbf{x}_N(t)
-\}$with respect to their center of geometry, -$\mathbf{x}_{\mathrm{cog}}(t)$: - - \label{eq:colvar_gyration} - R_{\mathrm{gyr}} \; = \; \sqrt{ \frac{1}{N} - \sum_{i=1}^{N} \left|\mathbf{x}_{i}(t) - - \mathbf{x}_{\mathrm{cog}}(t)\right|^{2} } - -This component must contain one \texttt{atoms~\{...\}} block to -define the atom group, and returns a positive number, expressed in -\AA{}. +\cvsubsec{Defining grid parameters}{sec:colvar_grid_params} -\begin{cvcoptions} -\item % - \dupkey{atoms}{\texttt{gyration}}{colvar|rmsd|atoms}{\texttt{rmsd} component} -\end{cvcoptions} +Many algorithms require the definition of boundaries and/or characteristic spacings that can be used to define discrete states'' in the collective variable, or to combine variables with very different units. +The parameters described below offer a way to specify these parameters only once for each variable, while using them multiple times in restraints, time-dependent biases or analysis methods. +\begin{itemize} -\cvsubsubsec{\texttt{inertia}: total moment of inertia of a group of atoms.}{sec:cvc_inertia} -The block \texttt{inertia~\{...\}} defines the -parameters for calculating the total moment of inertia of a group of atomic -positions$\{ \mathbf{x}_1(t), \mathbf{x}_2(t), \ldots \mathbf{x}_N(t)
-\}$with respect to their center of geometry, -$\mathbf{x}_{\mathrm{cog}}(t)$: - - \label{eq:colvar_inertia} - I \; = \; \sum_{i=1}^{N} \left|\mathbf{x}_{i}(t) - - \mathbf{x}_{\mathrm{cog}}(t)\right|^{2} - -\emph{Note that all atomic masses are set to 1 for simplicity.} -This component must contain one \texttt{atoms~\{...\}} block to -define the atom group, and returns a positive number, expressed in -\AA{}$^{2}$. +\item % + \labelkey{colvar|width} + \keydef + {width}{% + \texttt{colvar}}{% + Colvar fluctuation scale, or resolution for grid-based methods}{% + positive decimal}{% + 1.0}{% + This number has the same physical unit as the colvar value and defines an effective colvar unit. +% Biasing algorithms use it for different purposes. + Harmonic restraints (\ref{sec:colvarbias_harmonic}) use it to set the physical unit of the force constant, which is useful for multidimensional restraints involving variables with very different units (for examples,$\AA$or degrees$^\circ$) with a single, scaled force constant. + Histograms (\ref{sec:colvarbias_histogram}), ABF (\ref{sec:colvarbias_abf}) and metadynamics (\ref{sec:colvarbias_meta}) all use this number as the initial choice for the grid spacing along this variable: for this reason, \texttt{width} should generally be no larger than the standard deviation of the colvar in an unbiased simulation. +% When a non-unity width is required by the application, the optimal value is application-dependent, but can often be thought of as a user-provided estimate of the fluctuation amplitude for the colvar. + Unless it is required to control the spacing, it is usually simplest to keep the default value of 1, so that restraint force constants are provided in more intuitive units. + } -\begin{cvcoptions} \item % - \dupkey{atoms}{\texttt{inertia}}{colvar|rmsd|atoms}{\texttt{rmsd} component} -\end{cvcoptions} + \labelkey{colvar|lowerBoundary} + \key + {lowerBoundary}{% + \texttt{colvar}}{% + Lower boundary of the colvar}{% + decimal}{% + Defines the lowest end of the interval of relevant'' values for the colvar. + This number can be either a true physical boundary, or a user-defined number. + Together with \texttt{upperBoundary} and \texttt{width}, it is used to define a grid of values along the colvar (not available for colvars based on \texttt{distanceDir}, \texttt{distanceVec}, and \texttt{orientation}). + This option does not affect dynamics: to confine a colvar within a certain interval, use a \texttt{harmonicWalls} bias. +} +\item % + \labelkey{colvar|upperBoundary} % + \key + {upperBoundary}{% + \texttt{colvar}}{% + Upper boundary of the colvar}{% + decimal}{% + Similarly to \texttt{lowerBoundary}, defines the highest possible or allowed value.} -\cvsubsubsec{\texttt{inertiaZ}: total moment of inertia of a group of atoms around a chosen axis.}{sec:cvc_inertiaZ} -The block \texttt{inertiaZ~\{...\}} defines the -parameters for calculating the component along the axis$\mathbf{e}$of the moment of inertia of a group of atomic -positions$\{ \mathbf{x}_1(t), \mathbf{x}_2(t), \ldots \mathbf{x}_N(t)
-\}$with respect to their center of geometry, -$\mathbf{x}_{\mathrm{cog}}(t)$: - - \label{eq:colvar_inertia_z} - I_{\mathbf{e}} \; = \; \sum_{i=1}^{N} \left(\left(\mathbf{x}_{i}(t) - - \mathbf{x}_{\mathrm{cog}}(t)\right)\cdot\mathbf{e}\right)^{2} - -\emph{Note that all atomic masses are set to 1 for simplicity.} -This component must contain one \texttt{atoms~\{...\}} block to -define the atom group, and returns a positive number, expressed in -\AA{}$^{2}$. +\item % + \keydef + {hardLowerBoundary}{% + \texttt{colvar}}{% + Whether the lower boundary is the physical lower limit}{% + boolean}{% + \texttt{off}}{% + This option does not affect simulation results, but enables some internal optimizations. + Depending on its mathematical definition, a colvar may have natural'' boundaries: for example, a \texttt{distance} colvar has a natural'' lower boundary at 0. Setting this option instructs the Colvars module that the user-defined lower boundary is natural''. +See Section~\ref{sec:cvc_list} for the physical ranges of values of each component.} -\begin{cvcoptions} \item % - \dupkey{atoms}{\texttt{inertiaZ}}{colvar|rmsd|atoms}{\texttt{rmsd} component} + \keydef + {hardUpperBoundary}{% + \texttt{colvar}}{% + Whether the upper boundary is the physical upper limit of the colvar's values}{% + boolean}{% + \texttt{off}}{% + Analogous to \texttt{hardLowerBoundary}.} + \item % + \labelkey{colvar|expandBoundaries} % \keydef - {axis}{% - \texttt{inertiaZ}}{% - Projection axis (\AA{})}{% - \texttt{(x, y, z)} triplet}{% - \texttt{(0.0, 0.0, 1.0)}}{% - The three components of this vector define (when normalized) the - projection axis$\mathbf{e}$.} -\end{cvcoptions} + {expandBoundaries}{% + \texttt{colvar}}{% + Allow to expand the two boundaries if needed}{% + boolean}{% + \texttt{off}}{% + If defined, biasing and analysis methods may keep their own copies + of \texttt{lowerBoundary} and \texttt{upperBoundary}, and expand + them to accommodate values that do not fit in the initial range. + Currently, this option is used by the metadynamics bias + (\ref{sec:colvarbias_meta}) to keep all of its hills fully within + the grid. This option cannot be used when + the initial boundaries already span the full period of a periodic + colvar.} +\end{itemize} -\cvsubsubsec{\texttt{orientation}: orientation from reference coordinates.}{sec:cvc_orientation} -The block \texttt{orientation~\{...\}} returns the -same optimal rotation used in the \texttt{rmsd} component to -superimpose the coordinates$\{\mathbf{x}_{i}(t)\}$onto a set of -reference coordinates$\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}$. Such -component returns a four dimensional vector$\mathsf{q} = (q_0, q_1,
-q_2, q_3)$, with$\sum_{i} q_{i}^{2} = 1$; this \emph{quaternion} -expresses the optimal rotation$\{\mathbf{x}_{i}(t)\} \rightarrow
-\{\mathbf{x}_{i}^{\mathrm{(ref)}}\}$according to the formalism in -reference~\cite{Coutsias2004}. The quaternion$(q_0, q_1, q_2, q_3)$-can also be written as$\left(\cos(\theta/2), \,
-  \sin(\theta/2)\mathbf{u}\right)$, where$\theta$is the angle and -$\mathbf{u}$the normalized axis of rotation; for example, a rotation -of 90$^{\circ}$around the$z$axis is expressed as -\texttt{(0.707, 0.0, 0.0, 0.707)}''. The script -\texttt{quaternion2rmatrix.tcl} provides Tcl functions for converting -to and from a$4\times{}4$rotation matrix in a format suitable for -usage in VMD. -As for the component \texttt{rmsd}, the available options are \texttt{atoms}, \texttt{refPositionsFile}\cvnamebasedonly{, \texttt{refPositionsCol} and \texttt{refPositionsColValue}, } and \texttt{refPositions}. -\textbf{Note:} \texttt{refPositions}and \texttt{refPositionsFile} define the set of positions \emph{from which} the optimal rotation is calculated, but this rotation is not applied to the coordinates of the atoms involved: it is used instead to define the variable itself. +\cvsubsec{Trajectory output}{sec:colvar_traj_output} -\begin{cvcoptions} -\item % - \dupkey{atoms}{\texttt{orientation}}{colvar|rmsd|atoms}{\texttt{rmsd} component} +\begin{itemize} \item % - \dupkey{refPositions}{\texttt{orientation}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} + \keydef + {outputValue}{% + \texttt{colvar}}{% + Output a trajectory for this colvar}{% + boolean}{% + \texttt{on}}{% + If \texttt{colvarsTrajFrequency} is non-zero, the value of this + colvar is written to the trajectory file every + \texttt{colvarsTrajFrequency} steps in the column labeled + $<$\texttt{name}$>$''.} + \item % - \dupkey{refPositionsFile}{\texttt{orientation}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} + \keydef + {outputVelocity}{% + \texttt{colvar}}{% + Output a velocity trajectory for this colvar}{% + boolean}{% + \texttt{off}}{% + If \texttt{colvarsTrajFrequency} is defined, the + finite-difference calculated velocity of this colvar are written + to the trajectory file under the label + \texttt{v\_}$<$\texttt{name}$>$''.} -\cvnamebasedonly{ \item % - \dupkey{refPositionsCol}{\texttt{orientation}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} + \keydef + {outputEnergy}{% + \texttt{colvar}}{% + Output an energy trajectory for this colvar}{% + boolean}{% + \texttt{off}}{% + This option applies only to extended Lagrangian colvars. If + \texttt{colvarsTrajFrequency} is defined, the kinetic energy of + the extended degree and freedom and the potential energy of the + restraining spring are are written to the trajectory file under + the labels \texttt{Ek\_}$<$\texttt{name}$>$'' and + \texttt{Ep\_}$<$\texttt{name}$>$''.} + \item % - \dupkey{refPositionsColValue}{\texttt{orientation}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} -} + \labelkey{colvar|outputTotalForce} + \keydef + {outputTotalForce}{% + \texttt{colvar}}{% + Output a total force trajectory for this + colvar}{% + boolean}{% + \texttt{off}}{% + If \texttt{colvarsTrajFrequency} is defined, the total force on this + colvar (i.e.~the projection of all atomic total forces + onto this colvar --- see + equation~(\ref{eq:gradient_vector}) in + section~\ref{sec:colvarbias_abf}) are written to the trajectory + file under the label \texttt{fs\_}$<$\texttt{name}$>$''. + For extended Lagrangian colvars, the total force'' felt by the extended degree of freedom + is simply the force from the harmonic spring. + \textbf{Note:} not all components support this option. The + physical unit for this force is \cvnamdonly{kcal/mol}\cvvmdonly{kcal/mol}\cvlammpsonly{the unit of energy specified by \texttt{units}}, divided by the colvar unit U.} \item % \keydef - {closestToQuaternion}{% - \texttt{orientation}}{% - Reference rotation}{% - \texttt{(q0, q1, q2, q3)}'' quadruplet}{% - \texttt{(1.0, 0.0, 0.0, 0.0)} (null'' rotation)}{% - Between the two equivalent quaternions$(q_0, q_1, q_2, q_3)$and -$(-q_0, -q_1, -q_2, -q_3)$, the closer to \texttt{(1.0, 0.0, 0.0, - 0.0)} is chosen. This simplifies the visualization of the - colvar trajectory when samples values are a smaller subset of all - possible rotations. \textbf{Note:} \emph {this only affects the - output, never the dynamics}.} + {outputAppliedForce}{% + \texttt{colvar}}{% + Output an applied force trajectory for this + colvar}{% + boolean}{% + \texttt{off}}{% + If \texttt{colvarsTrajFrequency} is defined, the total force + applied on this colvar by Colvars biases are + written to the trajectory under the label + \texttt{fa\_}$<$\texttt{name}$>$''. + For extended Lagrangian colvars, this force is actually applied to the + extended degree of freedom rather than the geometric colvar itself. + The physical unit for this + force is \cvnamdonly{kcal/mol}\cvvmdonly{kcal/mol}\cvlammpsonly{the unit of energy specified by \texttt{units}} divided by the colvar unit.} -\end{cvcoptions} +\end{itemize} -\textbf{Tip: stopping the rotation of a protein.} To stop the -rotation of an elongated macromolecule in solution (and use an -anisotropic box to save water molecules), it is possible to define a -colvar with an \texttt{orientation} component, and restrain it throuh -the \texttt{harmonic} bias around the identity rotation, \texttt{(1.0, - 0.0, 0.0, 0.0)}. Only the overall orientation of the macromolecule -is affected, and \emph{not} its internal degrees of freedom. The user -should also take care that the macromolecule is composed by a single -chain, or disable \texttt{wrapAll} otherwise. +\cvsubsec{Extended Lagrangian}{sec:colvar_extended} +The following options enable extended-system +dynamics, where a colvar is coupled to an additional degree of freedom +(fictitious particle) by a harmonic spring. +All biasing and confining forces are then applied to the extended degree +of freedom. The actual'' geometric colvar (function of Cartesian +coordinates) only feels the force from the harmonic spring. +This is particularly useful when combined with an ABF bias (\ref{sec:colvarbias_abf}) +to perform eABF simulations (\ref{sec:eABF}). -\cvsubsubsec{\texttt{orientationAngle}: angle of rotation from reference coordinates.}{sec:cvc_orientationAngle} -The block \texttt{orientationAngle~\{...\}} accepts the same base options as -the component \texttt{orientation}: \texttt{atoms}, \texttt{refPositions}, \texttt{refPositionsFile}\cvnamebasedonly{, \texttt{refPositionsCol} and \texttt{refPositionsColValue}}. -The returned value is the angle of rotation$\theta$between the current and the reference positions. -This angle is expressed in degrees within the range [0$^{\circ}$:180$^{\circ}$]. +\begin{itemize} +\item % + \keydef + {extendedLagrangian}{% + \texttt{colvar}}{% + Add extended degree of freedom}{% + boolean}{% + \texttt{off}}{% + Adds a fictitious particle to be coupled to the colvar by a harmonic + spring. The fictitious mass and the force constant of the coupling + potential are derived from the parameters \texttt{extendedTimeConstant} + and \texttt{extendedFluctuation}, described below. Biasing forces on the + colvar are applied to this fictitious particle, rather than to the + atoms directly. This implements the extended Lagrangian formalism + used in some metadynamics simulations~\cite{Iannuzzi2003}. + \cvnamdonly{The energy associated with the extended degree of freedom is reported + under the MISC title in NAMD's energy output.} + } + +\item % + \key + {extendedFluctuation}{% + \texttt{colvar}}{% + Standard deviation between the colvar and the fictitious + particle (colvar unit)}{% + positive decimal}{% + Defines the spring stiffness for the \texttt{extendedLagrangian} + mode, by setting the typical deviation between the colvar and the extended + degree of freedom due to thermal fluctuation. + The spring force constant is calculated internally as$k_B T / \sigma^2$, + where$\sigma$is the value of \texttt{extendedFluctuation}.} + +\item % + \keydef + {extendedTimeConstant}{% + \texttt{colvar}}{% + Oscillation period of the fictitious particle (fs)}{% + positive decimal}{% + \texttt{200}}{% + Defines the inertial mass of the fictitious particle, by setting the + oscillation period of the harmonic oscillator formed by the fictitious + particle and the spring. The period + should be much larger than the MD time step to ensure accurate integration + of the extended particle's equation of motion. + The fictitious mass is calculated internally as$k_B T (\tau/2 \pi \sigma)^2$, + where$\tau$is the period and$\sigma$is the typical fluctuation (see above).} -\begin{cvcoptions} \item % - \dupkey{atoms}{\texttt{orientationAngle}}{colvar|rmsd|atoms}{\texttt{rmsd} component} -\item % - \dupkey{refPositions}{\texttt{orientationAngle}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} -\item % - \dupkey{refPositionsFile}{\texttt{orientationAngle}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} + \keydef + {extendedTemp}{% + \texttt{colvar}}{% + Temperature for the extended degree of freedom (K)}{% + positive decimal}{% + thermostat temperature}{% + Temperature used for calculating the coupling force constant of the + extended variable (see \texttt{extendedFluctuation}) and, if needed, as a + target temperature for extended Langevin dynamics (see + \texttt{extendedLangevinDamping}). This should normally be left at its + default value.} -\cvnamebasedonly{ -\item % - \dupkey{refPositionsCol}{\texttt{orientationAngle}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} \item % - \dupkey{refPositionsColValue}{\texttt{orientationAngle}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} -} + \keydef + {extendedLangevinDamping}{% + \texttt{colvar}}{% + Damping factor for extended Langevin dynamics + (ps$^{-1}$)}{% + positive decimal}{% + \texttt{1.0}}{% + If this is non-zero, the extended degree of freedom undergoes Langevin dynamics + at temperature \texttt{extendedTemp}. The friction force is minus + \texttt{extendedLangevinDamping} times the velocity. This is useful because + the extended dynamics coordinate may heat up in the transient + non-equilibrium regime of ABF. Use moderate damping values, to limit + viscous friction (potentially slowing down diffusive sampling) and stochastic + noise (increasing the variance of statistical measurements). In + doubt, use the default value.} +\end{itemize} -\end{cvcoptions} +\cvsubsec{Backward-compatibility}{sec:colvar_compatibility} -\cvsubsubsec{\texttt{orientationProj}: cosine of the angle of rotation from reference coordinates.} {sec:cvc_orientationProj} -The block \texttt{orientationProj~\{...\}} accepts the same base options as -the component \texttt{orientation}: \texttt{atoms}, \texttt{refPositions}, \texttt{refPositionsFile}\cvnamebasedonly{, \texttt{refPositionsCol} and \texttt{refPositionsColValue}}. -The returned value is the cosine of the angle of rotation$\theta$between the current and the reference positions. -The range of values is [-1:1]. +\begin{itemize} -\begin{cvcoptions} -\item % - \dupkey{atoms}{\texttt{orientationProj}}{colvar|rmsd|atoms}{\texttt{rmsd} component} -\item % - \dupkey{refPositions}{\texttt{orientationProj}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} -\item % - \dupkey{refPositionsFile}{\texttt{orientationProj}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} -\cvnamebasedonly{ \item % - \dupkey{refPositionsCol}{\texttt{orientationProj}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} -\item % - \dupkey{refPositionsColValue}{\texttt{orientationProj}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} + \keydef + {subtractAppliedForce}{% + \texttt{colvar}}{% + Do not include biasing forces in the total force for this colvar}{% + boolean}{% + \texttt{off}}{% + If the colvar supports total force calculation (see \ref{sec:cvc_sys_forces}), all forces applied to this colvar by biases will be removed from the total force. + This keyword allows to recover some of the system force'' calculation available in the Colvars module before version 2016-08-10. + Please note that removal of \emph{all} other external forces (including biasing forces applied to a different colvar) is \emph{no longer supported}, due to changes in the underlying simulation engines (primarily NAMD). + This option may be useful when continuing a previous simulation where the removal of external/applied forces is essential. + \emph{For all new simulations, the use of this option is not recommended.} } -\end{cvcoptions} - -\cvsubsubsec{\texttt{spinAngle}: angle of rotation around a given axis.}{sec:cvc_spinAngle} -The complete rotation described by \texttt{orientation} can optionally be decomposed into two sub-rotations: one is a \emph{spin}'' rotation around \textbf{e}, and the other a \emph{tilt}'' rotation around an axis orthogonal to \textbf{e}. -The component \texttt{spinAngle} measures the angle of the spin'' sub-rotation around \textbf{e}. +\end{itemize} -\begin{cvcoptions} -\item % - \dupkey{atoms}{\texttt{spinAngle}}{colvar|rmsd|atoms}{\texttt{rmsd} component} -\item % - \dupkey{refPositions}{\texttt{spinAngle}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} -\item % - \dupkey{refPositionsFile}{\texttt{spinAngle}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} -\cvnamebasedonly{ -\item % - \dupkey{refPositionsCol}{\texttt{spinAngle}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} -\item % - \dupkey{refPositionsColValue}{\texttt{spinAngle}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} -} -\item % - \labelkey{colvar|spinAngle|axis} - \keydef - {axis}{% - \texttt{tilt}}{% - Special rotation axis (\AA{})}{% - \texttt{(x, y, z)} triplet}{% - \texttt{(0.0, 0.0, 1.0)}}{% - The three components of this vector define (when normalized) the special rotation axis used to calculate the \texttt{tilt} and \texttt{spinAngle} components.} -\end{cvcoptions} -The component \texttt{spinAngle} returns an angle (in degrees) within the periodic interval$[-180:180]$. -\textbf{Note:} the value of \texttt{spinAngle} is a continuous function almost everywhere, with the exception of configurations with the corresponding tilt'' angle equal to 180$^\circ$(i.e.~the \texttt{tilt} component is equal to$-1$): in those cases, \texttt{spinAngle} is undefined. If such configurations are expected, consider defining a \texttt{tilt} colvar using the same axis \textbf{e}, and restraining it with a lower wall away from$-1$. +\cvsubsec{Statistical analysis}{sec:colvar_acf} +When the global keyword \texttt{analysis} is defined in the +configuration file, run-time calculations of statistical properties for +individual colvars can be performed. At the moment, several types of +time correlation functions, running averages and running standard +deviations are available. -\cvsubsubsec{\texttt{tilt}: cosine of the rotation orthogonal to a given axis.}{sec:cvc_tilt} -The component \texttt{tilt} measures the cosine of the angle of the tilt'' sub-rotation, which combined with the spin'' sub-rotation provides the complete rotation of a group of atoms. -The cosine of the tilt angle rather than the tilt angle itself is implemented, because the latter is unevenly distributed even for an isotropic system: consider as an analogy the angle$\theta$in the spherical coordinate system. -The component \texttt{tilt} relies on the same options as \texttt{spinAngle}, including the definition of the axis \textbf{e}. -The values of \texttt{tilt} are real numbers in the interval$[-1:1]$: the value$1$represents an orientation fully parallel to \textbf{e} (tilt angle = 0$^\circ$), and the value$-1$represents an anti-parallel orientation. +\begin{itemize} -\begin{cvcoptions} -\item % - \dupkey{atoms}{\texttt{tilt}}{colvar|rmsd|atoms}{\texttt{rmsd} component} -\item % - \dupkey{refPositions}{\texttt{tilt}}{colvar|rmsd|refPositions}{\texttt{rmsd} component} -\item % - \dupkey{refPositionsFile}{\texttt{tilt}}{colvar|rmsd|refPositionsFile}{\texttt{rmsd} component} -\cvnamebasedonly{ -\item % - \dupkey{refPositionsCol}{\texttt{tilt}}{colvar|rmsd|refPositionsCol}{\texttt{rmsd} component} -\item % - \dupkey{refPositionsColValue}{\texttt{tilt}}{colvar|rmsd|refPositionsColValue}{\texttt{rmsd} component} -} \item % - \dupkey{axis}{\texttt{tilt}}{colvar|spinAngle|axis}{\texttt{spinAngle} component} -\end{cvcoptions} - -\cvnamebasedonly{ - -\cvsubsubsec{\texttt{alpha}:$\alpha$-helix content of a protein segment.}{sec:cvc_alpha} -The block \texttt{alpha~\{...\}} defines the -parameters to calculate the helical content of a segment of protein -residues. The$\alpha$-helical content across the$N+1$residues -$N_{0}$to$N_{0}+N$is calculated by the formula: -\begin{eqnarray} - \label{eq:colvars_alpha} - { - \alpha\left( - \mathrm{C}_{\alpha}^{(N_{0})}, - \mathrm{O}^{(N_{0})}, - \mathrm{C}_{\alpha}^{(N_{0}+1)}, - \mathrm{O}^{(N_{0}+1)}, - \ldots - \mathrm{N}^{(N_{0}+5)}, - \mathrm{C}_{\alpha}^{(N_{0}+5)}, - \mathrm{O}^{(N_{0}+5)}, - \ldots - \mathrm{N}^{(N_{0}+N)}, - \mathrm{C}_{\alpha}^{(N_{0}+N)} - \right) - } \; = \; \; \; \; \\ \; \; \; \; { - \nonumber - \frac{1}{2(N-2)} - \sum_{n=N_{0}}^{N_{0}+N-2} - \mathrm{angf}\left( - \mathrm{C}_{\alpha}^{(n)}, - \mathrm{C}_{\alpha}^{(n+1)}, - \mathrm{C}_{\alpha}^{(n+2)}\right) - } \; + \; { - \frac{1}{2(N-4)} - \sum_{n=N_{0}}^{N_{0}+N-4} - \mathrm{hbf}\left( - \mathrm{O}^{(n)}, - \mathrm{N}^{(n+4)}\right) \mathrm{,} - } \\ -\end{eqnarray} -where the score function for the$\mathrm{C}_{\alpha} -
-\mathrm{C}_{\alpha} - \mathrm{C}_{\alpha}$angle is defined as: - - \label{eq:colvars_alpha_Calpha} - { - \mathrm{angf}\left( - \mathrm{C}_{\alpha}^{(n)}, - \mathrm{C}_{\alpha}^{(n+1)}, - \mathrm{C}_{\alpha}^{(n+2)}\right) - } \; = \; { - \frac{1 - \left(\theta( - \mathrm{C}_{\alpha}^{(n)}, - \mathrm{C}_{\alpha}^{(n+1)}, - \mathrm{C}_{\alpha}^{(n+2)}) - - \theta_{0}\right)^{2} / - \left(\Delta\theta_{\mathrm{tol}}\right)^{2}}{ - 1 - \left(\theta( - \mathrm{C}_{\alpha}^{(n)}, - \mathrm{C}_{\alpha}^{(n+1)}, - \mathrm{C}_{\alpha}^{(n+2)}) - - \theta_{0}\right)^{4} / - \left(\Delta\theta_{\mathrm{tol}}\right)^{4}} \mathrm{,} - } - -and the score function for the$\mathrm{O}^{(n)} \leftrightarrow
-\mathrm{N}^{(n+4)}$hydrogen bond is defined through a \texttt{hBond} -colvar component on the same atoms. - -\begin{cvcoptions} + \keydef + {corrFunc}{% + \texttt{colvar}}{% + Calculate a time correlation function?}{% + boolean}{% + \texttt{off}}{% + Whether or not a time correlaction function should be calculated + for this colvar.} \item % - \labelkey{colvar|alpha|residueRange} \key - {residueRange}{% - \texttt{alpha}}{% - Potential$\alpha$-helical residues}{% - $<$Initial residue number$>$-$<$Final residue number$>$''}{% - This option specifies the range of residues on which this - component should be defined. The Colvars module looks for the - atoms within these residues named \texttt{CA}'', \texttt{N}'' - and \texttt{O}'', and raises an error if any of those atoms is - not found.} + {corrFuncWithColvar}{% + \texttt{colvar}}{% + Colvar name for the correlation function}{% + string}{% + By default, the auto-correlation function (ACF) of this colvar, +$\xi_{i}$, is calculated. When this option is specified, the + correlation function is calculated instead with another colvar, +$\xi_{j}$, which must be of the same type (scalar, vector, or + quaternion) as$\xi_{i}$.} + +\item% + \keydef + {corrFuncType}{% + \texttt{colvar}}{% + Type of the correlation function}{% + \texttt{velocity}, \texttt{coordinate} or + \texttt{coordinate\_p2}}{% + \texttt{velocity}}{% + With \texttt{coordinate} or \texttt{velocity}, the correlation + function$C_{i,j}(t)$~=$\left\langle \Pi\left(\xi_{i}(t_{0}),
+        \xi_{j}(t_{0}+t)\right) \right\rangle$is calculated between + the variables$\xi_{i}$and$\xi_{j}$, or their velocities. +$\Pi(\xi_{i}, \xi_{j})$is the scalar product when calculated + between scalar or vector values, whereas for quaternions it is the + cosine between the two corresponding rotation axes. With + \texttt{coordinate\_p2}, the second order Legendre polynomial, +$(3\cos(\theta)^{2}-1)/2$, is used instead of the cosine.} \item % - \labelkey{colvar|alpha|psfSegID} - \key - {psfSegID}{% - \texttt{alpha}}{% - PSF segment identifier}{% - string (max 4 characters)}{% - This option sets the PSF segment identifier for the residues - specified in \texttt{residueRange}. This option is only required - when PSF topologies are used.} + \keydef + {corrFuncNormalize}{% + \texttt{colvar}}{% + Normalize the time correlation function?}{% + boolean}{% + \texttt{on}}{% + If enabled, the value of the correlation function at$t$~= 0 + is normalized to 1; otherwise, it equals to$\left\langle
+      O\left(\xi_{i}, \xi_{j}\right) \right\rangle$.} +\item % + \keydef + {corrFuncLength}{% + \texttt{colvar}}{% + Length of the time correlation function}{% + positive integer}{% + \texttt{1000}}{% + Length (in number of points) of the time correlation function.} \item % - \labelkey{colvar|alpha|hBondCoeff} \keydef - {hBondCoeff}{% - \texttt{alpha}}{% - Coefficient for the hydrogen bond term}{% - positive between 0 and 1}{% - 0.5}{% - This number specifies the contribution to the total value from the - hydrogen bond terms. 0 disables the hydrogen bond terms, 1 - disables the angle terms.} + {corrFuncStride}{% + \texttt{colvar}}{% + Stride of the time correlation function}{% + positive integer}{% + \texttt{1}}{% + Number of steps between two values of the time correlation function.} \item % - \labelkey{colvar|alpha|angleRef} \keydef - {angleRef}{% - \texttt{alpha}}{% - Reference$\mathrm{C}_{\alpha} -
-    \mathrm{C}_{\alpha} - \mathrm{C}_{\alpha}$angle}{% - positive decimal}{% - 88$^{\circ}$}{% - This option sets the reference angle used in the score function - (\ref{eq:colvars_alpha_Calpha}).} + {corrFuncOffset}{% + \texttt{colvar}}{% + Offset of the time correlation function}{% + positive integer}{% + \texttt{0}}{% + The starting time (in number of steps) of the time correlation + function (default:$t$~= 0). \textbf{Note:} \emph{the value at$t$~= 0 is always + used for the normalization}.} \item % - \labelkey{colvar|alpha|angleTol} \keydef - {angleTol}{% - \texttt{alpha}}{% - Tolerance in the$\mathrm{C}_{\alpha} -
-    \mathrm{C}_{\alpha} - \mathrm{C}_{\alpha}$angle}{% - positive decimal}{% - 15$^{\circ}$}{% - This option sets the angle tolerance used in the score function - (\ref{eq:colvars_alpha_Calpha}).} + {corrFuncOutputFile}{% + \texttt{colvar}}{% + Output file for the time correlation function}{% + UNIX filename}{% + \texttt{$<$name$>$.corrfunc.dat}}{% + The time correlation function is saved in this file.} \item % - \labelkey{colvar|alpha|hBondCutoff} \keydef - {hBondCutoff}{% - \texttt{alpha}}{% - Hydrogen bond cutoff}{% - positive decimal}{% - 3.3~\AA{}}{% - Equivalent to the \texttt{cutoff} option in the \texttt{hBond} - component.} + {runAve}{% + \texttt{colvar}}{% + Calculate the running average and standard deviation}{% + boolean}{% + \texttt{off}}{% + Whether or not the running average and standard deviation should + be calculated for this colvar.} \item % - \labelkey{colvar|alpha|hBondExpNumer} \keydef - {hBondExpNumer}{% - \texttt{alpha}}{% - Hydrogen bond numerator exponent}{% + {runAveLength}{% + \texttt{colvar}}{% + Length of the running average window}{% positive integer}{% - 6}{% - Equivalent to the \texttt{expNumer} option in the \texttt{hBond} - component.} + \texttt{1000}}{% + Length (in number of points) of the running average window.} \item % - \labelkey{colvar|alpha|hBondExpDenom} \keydef - {hBondExpDenom}{% - \texttt{alpha}}{% - Hydrogen bond denominator exponent}{% + {runAveStride}{% + \texttt{colvar}}{% + Stride of the running average window values}{% positive integer}{% - 8}{% - Equivalent to the \texttt{expDenom} option in the \texttt{hBond} - component.} - -\end{cvcoptions} - -This component returns positive values, always comprised between 0 -(lowest$\alpha$-helical score) and 1 (highest$\alpha$-helical -score). + \texttt{1}}{% + Number of steps between two values within the running average window.} +\item % + \keydef + {runAveOutputFile}{% + \texttt{colvar}}{% + Output file for the running average and standard deviation}{% + UNIX filename}{% + \texttt{$<$name$>$.runave.dat}}{% + The running average and standard deviation are saved in this file.} -\cvsubsubsec{\texttt{dihedralPC}: protein dihedral pricipal component}{sec:cvc_dihedralPC} -The block \texttt{dihedralPC~\{...\}} defines the -parameters to calculate the projection of backbone dihedral angles within -a protein segment onto a \emph{dihedral principal component}, following -the formalism of dihedral principal component analysis (dPCA) proposed by -Mu et al.\cite{Mu2005} and documented in detail by Altis et -al.\cite{Altis2007}. -Given a peptide or protein segment of$N$residues, each with Ramachandran -angles$\phi_i$and$\psi_i$, dPCA rests on a variance/covariance analysis -of the$4(N-1)$variables$\cos(\psi_1), \sin(\psi_1), \cos(\phi_2), \sin(\phi_2)
-\cdots \cos(\phi_N), \sin(\phi_N)$. Note that angles$\phi_1$and$\psi_N$-have little impact on chain conformation, and are therefore discarded, -following the implementation of dPCA in the analysis software Carma.\cite{Glykos2006} +\end{itemize} -For a given principal component (eigenvector) of coefficients -$(k_i)_{1 \leq i \leq 4(N-1)}$, -the projection of the current backbone conformation is: - -\xi = \sum_{n=1}^{N-1} k_{4n-3} \cos(\psi_n) + k_{4n-2} \sin (\psi_n) -+ k_{4n-1} \cos (\phi_{n+1}) + k_{4n} \sin(\phi_{n+1}) - -\texttt{dihedralPC} expects the same parameters as the \texttt{alpha} -component for defining the relevant residues (\texttt{residueRange} -and \texttt{psfSegID}) in addition to the following: +\cvsec{Selecting atoms}{sec:colvar_atom_groups} -\begin{cvcoptions} +To define collective variables, atoms are usually selected as groups. Each group is defined using an identifier that is unique in the context of the specific colvar component (e.g.~for a distance component, the two groups are \texttt{group1} and \texttt{group2}). +The identifier is followed by a brace-delimited block containing selection keywords and other parameters, including an optional \texttt{name}: -\item % - \dupkey{residueRange}{\texttt{dihedralPC}}{colvar|alpha|residueRange}{\texttt{alpha} component} +\begin{itemize} +\item \key + {name}{% + atom group}{% + Unique name for the atom group}{% + string}{% + This parameter defines a unique name for this atom group, which can be referred to + in the definition of other atom groups (including in other colvars) by invoking + \texttt{atomsOfGroup} as a selection keyword.} +\end{itemize} -\item % - \dupkey{psfSegID}{\texttt{dihedralPC}}{colvar|alpha|psfSegID}{\texttt{alpha} component} -\item % - \key - {vectorFile}{% - \texttt{dihedralPC}}{% - File containing dihedral PCA eigenvector(s)}{% - file name}{% - A text file containing the coefficients of dihedral PCA eigenvectors on the - cosine and sine coordinates. The vectors should be arranged in columns, - as in the files output by Carma.\cite{Glykos2006}} +\cvsubsec{Atom selection keywords}{sec:colvar_atom_groups_sel} -\item % - \key - {vectorNumber}{% - \texttt{dihedralPC}}{% - File containing dihedralPCA eigenvector(s)}{% - positive integer}{% - Number of the eigenvector to be used for this component.} -\end{cvcoptions} +Selection keywords may be used individually or in combination with each other, and each can be repeated any number of times. +Selection is incremental: each keyword adds the corresponding atoms to the selection, so that different sets of atoms can be combined. +However, atoms included by multiple keywords are only counted once. +Below is an example configuration for an atom group called \texttt{atoms}''. +\textbf{Note: }\emph{this is an unusually varied combination of selection keywords, demonstrating how they can be combined together: most simulations only use one of them.}\\ -} % end of \cvnamebasedonly +{% +% verbatim can't appear within commands +\noindent\ttfamily atoms \{\\ +\\ +\-~~\# add atoms 1 and 3 to this group (note: the first atom in the system is 1)\\ +\-~~atomNumbers \{ \\ +\-~~~~1 3\\ +\-~~\}\\ +\\ +\-~~\# add atoms starting from 20 up to and including 50\\ +\-~~atomNumbersRange 20-50\\ +} +\cvnamebasedonly{{% +\noindent\ttfamily\\ +\-~~\# add all the atoms with occupancy 2 in the file atoms.pdb\\ +\-~~atomsFile atoms.pdb\\ +\-~~atomsCol O\\ +\-~~atomsColValue 2.0\\ +\\ +\-~~\# add all the C-alphas within residues 11 to 20 of segments "PR1" and "PR2"\\ +\-~~psfSegID PR1 PR2\\ +\-~~atomNameResidueRange CA 11-20\\ +\-~~atomNameResidueRange CA 11-20\\ +}} +{\noindent\ttfamily\\ +\-~~\# add index group (requires a .ndx file to be provided globally)\\ +\-~~indexGroup Water\\ +\}\\} +The resulting selection includes atoms 1 and 3, those between 20 and 50,\cvnamebasedonly{ the$\mathrm{C}_{\alpha}$atoms between residues 11 and 20 of the two segments \texttt{PR1} and \texttt{PR2},} and those in the index group called Water''. +The indices of this group are read from the file provided by the global keyword \refkey{indexFile}{Colvars-global|indexFile}. -\cvsubsec{Configuration keywords shared by all components}{sec:cvc_common} +\cvvmdonly{In the current version, the Colvars module does not manipulate VMD atom selections directly: however, these can be converted to atom groups within the Colvars configuration string, using selection keywords such as \texttt{atomNumbers}.} +The complete list of selection keywords available in \MDENGINE{} is: -The following options can be used for any of the above colvar components in order to obtain a polynomial combination\cvscriptonly{ or any user-supplied function provided by \refkey{scriptedFunction}{sec:cvc_superp}}. \begin{itemize} -\item % - \keydef - {name}{% - any component}{% - Name of this component}{% - string}{% - type of component + numeric id}{% - The name is an unique case-sensitive string which allows the - Colvars module to identify this component. This is useful, for example, - when combining multiple components via a \texttt{scriptedFunction}. - \cvnamdonly{It also defines the variable name representing the component's value in a \refkey{customFunction}{colvar|customFunction} expression.}} \item % - \keydef - {scalable}{% - any component}{% - Attempt to calculate this component in parallel?}{% - boolean}{% - \texttt{on}, if available}{% - If set to \texttt{on} (default), the Colvars module will attempt to calculate this component in parallel to reduce overhead. - Whether this option is available depends on the type of component: currently supported are \texttt{distance}, \texttt{distanceZ}, \texttt{distanceXY}, \texttt{distanceVec}, \texttt{distanceDir}, \texttt{angle} and \texttt{dihedral}. - This flag influences computational cost, but does not affect numerical results: therefore, it should only be turned off for debugging or testing purposes. + \key + {atomNumbers}{% + atom group}{% + List of atom numbers}{% + space-separated list of positive integers}{% + This option adds to the group all the atoms whose numbers are in + the list. \emph{The number of the first atom in the system is 1: to convert from a VMD selection, use atomselect get serial''.} } -\end{itemize} - - -\cvsubsec{Advanced usage and special considerations}{sec:cvc_advanced} - -\cvsubsubsec{Periodic components.}{sec:cvc_periodic} -The following components returns -real numbers that lie in a periodic interval: -\begin{itemize} -\item \texttt{dihedral}: torsional angle between four groups; -\item \texttt{spinAngle}: angle of rotation around a predefined axis - in the best-fit from a set of reference coordinates. -\end{itemize} -In certain conditions, \texttt{distanceZ} can also be periodic, namely -when periodic boundary conditions (PBCs) are defined in the simulation -and \texttt{distanceZ}'s axis is parallel to a unit cell vector. -In addition, a custom \cvscriptonly{or scripted} scalar colvar may be periodic -depending on its user-defined expression. It will only be treated as such by -the Colvars module if the period is specified using the \texttt{period} keyword, -while \texttt{wrapAround} is optional. +\item % + \key + {indexGroup}{% + atom group}{% + Name of index group to be used (GROMACS format)}{% + string}{% + If the name of an index file has been provided by \texttt{indexFile}, this option allows to select one index group from that file: the atoms from that index group will be used to define the current group.} -The following keywords can be used within periodic components, and within -custom \cvscriptonly{or scripted} colvars (\ref{sec:colvar_custom_function} -\cvscriptonly{, \ref{sec:colvar_scripted}}). +\item % + \key + {atomsOfGroup}{% + atom group}{% + Name of group defined previously}{% + string}{% + Refers to a group defined previously using its user-defined \texttt{name}. + This adds all atoms of that named group to the current group.} -\begin{itemize} \item % - \keydef - {period}{% - \texttt{distanceZ}, custom colvars}{% - Period of the component}{% - positive decimal}{% - 0.0}{% - Setting this number enables the treatment of \texttt{distanceZ} as - a periodic component: by default, \texttt{distanceZ} is not - considered periodic. The keyword is supported, but irrelevant - within \texttt{dihedral} or \texttt{spinAngle}, because their - period is always 360~degrees.} + \key + {atomNumbersRange}{% + atom group}{% + Atoms within a number range}{% +$<$Starting number$>$-$<$Ending number$>$}{% + This option includes in the group all atoms whose numbers are within the range specified. \emph{The number of the first atom in the system is 1.} + } +\cvnamebasedonly{ \item % - \keydef - {wrapAround}{% - \texttt{distanceZ}, \texttt{dihedral}, \texttt{spinAngle}, custom colvars}{% - Center of the wrapping interval for periodic variables}{% - decimal}{% - 0.0}{% - By default, values of the periodic components are centered around zero, ranging from$-P/2$to$P/2$, where$P$is the period. - Setting this number centers the interval around this value. - This can be useful for convenience of output, or to set the walls for a \texttt{harmonicWalls} in an order that would not otherwise be allowed.} -\end{itemize} + \key + {atomNameResidueRange}{% + atom group}{% + Named atoms within a range of residue numbers}{% +$<$Atom name$><$Starting residue$>$-$<$Ending residue$>$}{% + This option adds to the group all the atoms with the provided + name, within residues in the given range.} -Internally, all differences between two values of a periodic colvar -follow the minimum image convention: they are calculated based on -the two periodic images that are closest to each other. +\item % + \key + {psfSegID}{% + atom group}{% + PSF segment identifier}{% + space-separated list of strings (max 4 characters)}{% + This option sets the PSF segment identifier for + \texttt{atomNameResidueRange}. Multiple values may be provided, + which correspond to multiple instances of + \texttt{atomNameResidueRange}, in order of their occurrence. + This option is only necessary if a PSF topology file is used.} -\emph{Note: linear or polynomial combinations of periodic components (see \ref{sec:cvc_superp}) may become meaningless when components cross the periodic boundary. Use such combinations carefully: estimate the range of possible values of each component in a given simulation, and make use of \texttt{wrapAround} to limit this problem whenever possible.} +\item % + \key + {atomsFile}{% + atom group}{% + PDB file name for atom selection}{% + UNIX filename}{% + This option selects atoms from the PDB file provided and adds them + to the group according to numerical flags in the column + \texttt{atomsCol}. \textbf{Note:} \emph{the sequence of atoms in the PDB file + provided must match that in the system's topology}.} +\item % + \key + {atomsCol}{% + atom group}{% + PDB column to use for atom selection flags}{% + \texttt{O}, \texttt{B}, \texttt{X}, \texttt{Y}, or \texttt{Z}}{% + This option specifies which PDB column in \texttt{atomsFile} is used to determine which atoms are to be included in the group. + } -\cvsubsubsec{Non-scalar components.}{sec:cvc_non_scalar} +\item % + \key + {atomsColValue}{% + atom group}{% + Atom selection flag in the PDB column}{% + positive decimal}{% + If defined, this value in \texttt{atomsCol} identifies atoms in \texttt{atomsFile} that are included in the group. + If undefined, all atoms with a non-zero value in \texttt{atomsCol} are included.} +} -When one of the following components are used, the defined colvar returns a value that is not a scalar number: -\begin{itemize} -\item \texttt{distanceVec}: 3-dimensional vector of the distance - between two groups; -\item \texttt{distanceDir}: 3-dimensional unit vector of the distance - between two groups; -\item \texttt{orientation}: 4-dimensional unit quaternion representing - the best-fit rotation from a set of reference coordinates. -\end{itemize} -The distance between two 3-dimensional unit vectors is computed as the -angle between them. The distance between two quaternions is computed -as the angle between the two 4-dimensional unit vectors: because the -orientation represented by$\mathsf{q}$is the same as the one -represented by$-\mathsf{q}$, distances between two quaternions are -computed considering the closest of the two symmetric images. +\item % + \key + {dummyAtom}{% + atom group}{% + Dummy atom position (\AA{})}{% + \texttt{(x, y, z)} triplet}{% + Instead of selecting any atom, this option makes the group a virtual particle at a fixed position in space. This is useful e.g.~to replace a group's center of geometry with a user-defined position.} -Non-scalar components carry the following restrictions: -\begin{itemize} -\item Calculation of total forces (\texttt{outputTotalForce} option) - is currently not implemented. -\item Each colvar can only contain one non-scalar component. -\item Binning on a grid (\texttt{abf}, \texttt{histogram} and - \texttt{metadynamics} with \texttt{useGrids} enabled) is currently - not implemented for colvars based on such components. \end{itemize} -\emph{Note: while these restrictions apply to individual colvars based - on non-scalar components, no limit is set to the number of scalar - colvars. To compute multi-dimensional histograms and PMFs, use sets - of scalar colvars of arbitrary size.} - - -\cvsubsubsec{Calculating total forces.}{sec:cvc_sys_forces} -In addition to the restrictions due to the type of value computed (scalar or non-scalar), -a final restriction can arise when calculating total force -(\texttt{outputTotalForce} option or application of a \texttt{abf} -bias). total forces are available currently only for the following -components: \texttt{distance}, \texttt{distanceZ}, -\texttt{distanceXY}, \texttt{angle}, \texttt{dihedral}, \texttt{rmsd}, -\texttt{eigenvector} and \texttt{gyration}. +\cvsubsec{Moving frame of reference.}{sec:colvar_atom_groups_ref_frame} +The following options define an automatic calculation of an optimal translation (\texttt{centerReference}) or optimal rotation (\texttt{rotateReference}), that superimposes the positions of this group to a provided set of reference coordinates. +This can allow, for example, to effectively remove from certain colvars the effects of molecular tumbling and of diffusion. +Given the set of atomic positions$\mathbf{x}_{i}$, the colvar$\xi$can be defined on a set of roto-translated positions$\mathbf{x}_{i}' = R(\mathbf{x}_{i} - \mathbf{x}^{\mathrm{C}}) + \mathbf{x}^{\mathrm{ref}}$. +$\mathbf{x}^{\mathrm{C}}$is the geometric center of the$\mathbf{x}_{i}$,$R$is the optimal rotation matrix to the reference positions and$\mathbf{x}^{\mathrm{ref}}$is the geometric center of the reference positions. +Components that are defined based on pairwise distances are naturally invariant under global roto-translations. +Other components are instead affected by global rotations or translations: however, they can be made invariant if they are expressed in the frame of reference of a chosen group of atoms, using the \texttt{centerReference} and \texttt{rotateReference} options. +Finally, a few components are defined by convention using a roto-translated frame (e.g. the minimal RMSD): for these components, \texttt{centerReference} and \texttt{rotateReference} are enabled by default. +In typical applications, the default settings result in the expected behavior. -\cvsubsec{Linear and polynomial combinations of components}{sec:cvc_superp} +\paragraph*{Warning on rotating frames of reference and periodic boundary conditions.} +\texttt{rotateReference} affects coordinates that depend on minimum-image distances in periodic boundary conditions (PBC). +After rotation of the coordinates, the periodic cell vectors become irrelevant: the rotated system is effectively non-periodic. +A safe way to handle this is to ensure that the relevant inter-group distance vectors remain smaller than the half-size of the periodic cell. +If this is not desirable, one should avoid the rotating frame of reference, and apply orientational restraints to the reference group instead, in order to keep the orientation of the reference group consistent with the orientation of the periodic cell. -To extend the set of possible definitions of colvars$\xi(\mathbf{r})$, multiple components -$q_i(\mathbf{r})$can be summed with the formula: - - \label{eq:colvar_combination} - \xi(\mathbf{r}) = \sum_i c_i [q_i(\mathbf{r})]^{n_i} - -where each component appears with a unique coefficient$c_i$(1.0 by -default) the positive integer exponent$n_i(1 by default). +\paragraph*{Warning on rotating frames of reference and ABF.} +Note that \texttt{centerReference} and \texttt{rotateReference} may affect the Jacobian derivative of colvar components in a way that is not taken into account by default. +Be careful when using these options in ABF simulations or when using total force values. -Any set of components can be combined within a colvar, provided that -they return the same type of values (scalar, unit vector, vector, or -quaternion). By default, the colvar is the sum of its components. -Linear or polynomial combinations (following -equation~(\ref{eq:colvar_combination})) can be obtained by setting the -following parameters, which are common to all components: \begin{itemize} + \item % \keydef - {componentCoeff}{% - any component}{% - Coefficient of this component in the colvar}{% - decimal}{% - \texttt{1.0}}{% - Defines the coefficient by which this component is multiplied - (after being raised to \texttt{componentExp}) before being added - to the sum.} + {centerReference}{% + atom group}{% + Implicitly remove translations for this group}{% + boolean}{% + \texttt{off}}{% + If this option is \texttt{on}, the center of geometry of the group will be aligned with that of the reference positions provided by \cvnamebasedonly{either} \texttt{refPositions} or \texttt{refPositionsFile}. + Colvar components will only have access to the aligned positions. +\textbf{Note}: unless otherwise specified, \texttt{rmsd} and \texttt{eigenvector} set this option to \texttt{on} \emph{by default}. +} \item % \keydef - {componentExp}{% - any component}{% - Exponent of this component in the colvar}{% - integer}{% - \texttt{1}}{% - Defines the power at which the value of this component is raised - before being added to the sum. When this exponent is - different than 1 (non-linear sum), total forces and the Jacobian - force are not available, making the colvar unsuitable for ABF calculations.} -\end{itemize} + {rotateReference}{% + atom group}{% + Implicitly remove rotations for this group}{% + boolean}{% + \texttt{off}}{% + If this option is \texttt{on}, the coordinates of this group will be optimally superimposed to the reference positions provided by \cvnamebasedonly{either} \texttt{refPositions} or \texttt{refPositionsFile}. + The rotation will be performed around the center of geometry if \texttt{centerReference} is \texttt{on}, or around the origin otherwise. + The algorithm used is the same employed by the \texttt{orientation} colvar component~\cite{Coutsias2004}. + Forces applied to the atoms of this group will also be implicitly rotated back to the original frame. + \textbf{Note}: unless otherwise specified, \texttt{rmsd} and \texttt{eigenvector} set this option to \texttt{on} \emph{by default}. +} -\textbf{Example:} To define the \emph{average} of a colvar across -different parts of the system, simply define within the same colvar -block a series of components of the same type (applied to different -atom groups), and assign to each component a \texttt{componentCoeff} -of1/N$. +\item % + \labelkey{atom-group|refPositions} + \key + {refPositions}{% + atom group}{% + Reference positions for fitting (\AA)}{% + space-separated list of \texttt{(x, y, z)} triplets}{% + \label{key:colvars:atom_group:refPositions} + This option provides a list of reference coordinates for \texttt{centerReference} and/or \texttt{rotateReference}, and is mutually exclusive with \texttt{refPositionsFile}. + If only \texttt{centerReference} is \texttt{on}, the list may contain a single (x, y, z) triplet; if also \texttt{rotateReference} is \texttt{on}, the list should be as long as the atom group, and \emph{its order must match the order in which atoms were defined}. +} + +\item % + \labelkey{atom-group|refPositionsFile} + \key + {refPositionsFile}{% + atom group}{% + File containing the reference positions for fitting}{% + UNIX filename}{% + \label{key:colvars:atom_group:refPositionsFile} + This option provides a list of reference coordinates for \texttt{centerReference} and/or \texttt{rotateReference}, and is mutually exclusive with \texttt{refPositions}. + The acceptable file format is XYZ, which is read in double precision\cvnamebasedonly{, or PDB; \emph{the latter is discouraged if the precision of the reference coordinates is a concern}}. + Atomic positions are read differently depending on the following scenarios: + \textbf{(i)} the file contains exactly as many records as the atoms in the group: all positions are read in sequence; + \textbf{(ii)} (most common case) the file contains coordinates for the entire system: only the positions corresponding to the numeric indices of the atom group are read\cvnamebasedonly{; + \textbf{(iii)} if the file is a PDB file and \texttt{refPositionsCol} is specified, positions are read according to the value of the column \texttt{refPositionsCol} (which may be the same as \texttt{atomsCol})}. + In each case, atoms are read from the file \emph{in order of increasing number}. +} + + +\cvnamebasedonly{ +\item % + \key + {refPositionsCol}{% + atom group}{% + PDB column containing atom flags}{% + \texttt{O}, \texttt{B}, \texttt{X}, \texttt{Y}, or \texttt{Z}}{% + Like \texttt{atomsCol} for \texttt{atomsFile}, indicates which column to use to identify the atoms in \texttt{refPositionsFile} (if this is a PDB file).} +\item % + \key + {refPositionsColValue}{% + atom group}{% + Atom selection flag in the PDB column}{% + positive decimal}{% + Analogous to \texttt{atomsColValue}, but applied to \texttt{refPositionsCol}.} +} -\cvleptononly{ -\cvsubsec{Colvars as custom functions of components}{sec:colvar_custom_function} +\item % + \labelkey{atom-group|fittingGroup} % + \keydef + {fittingGroup}{% + atom group}{% + Use an alternate set of atoms to define the roto-translation}{% + Block \texttt{fittingGroup \{ ... \}}}{% + This group itself}{% + If either \texttt{centerReference} or \texttt{rotateReference} is defined, this keyword defines an alternate atom group to calculate the optimal roto-translation. + Use this option to define a continuous rotation if the structure of the group involved changes significantly (a typical symptom would be the message Warning: discontinuous rotation!''). -Collective variables may be defined by specifying a custom function as an analytical -expression such as \texttt{cos(x) + y\^{}2}. -The expression is parsed by the Lepton expression parser (written by Peter Eastman), -which produces efficient evaluation routines for the function itself as well as its derivatives. -The expression may use the collective variable components as variables, refered to as their \texttt{name} string. -Scalar elements of vector components may be accessed by appending a 1-based index to their \texttt{name}. -When implementing generic functions of Cartesian coordinates rather -than functions of existing components, the \texttt{cartesian} component -may be particularly useful. -A scalar-valued custom variable may be manually defined as periodic by providing -the keyword \texttt{period}, and the optional keyword \texttt{wrapAround}, with the -same meaning as in periodic components (see \ref{sec:cvc_periodic} for details). -A vector variable may be defined by specifying the \texttt{customFunction} parameter several times: each expression defines one scalar element of the vector colvar. -This is illustrated in the example below. +\cvnamebasedonly{ + The following example illustrates the syntax of \texttt{fittingGroup}: a group called \texttt{atoms}'' is defined, including 8 C$_{\alpha}$atoms of a protein of 100 residues. + An optimal roto-translation is calculated automatically by fitting the C$_{\alpha}$trace of the rest of the protein onto the coordinates provided by a PDB file.} -\bigskip {% -% verbatim can't appear within commands -\noindent\ttfamily colvar \{\\ -\-~~name custom\\ +\noindent\ttfamily +\# Example: defining a group "atoms", with its coordinates expressed \\ +\# on a roto-translated frame of reference defined by a second group\\ +atoms \{\\ \\ -\-~~\# A 2-dimensional vector function of a scalar x and a 3-vector r\\ -\-~~customFunction cos(x) * (r1 + r2 + r3)\\ -\-~~customFunction sqrt(r1 * r2)\\ +\-~~psfSegID PROT\\ +\-~~atomNameResidueRange CA 41-48\\ \\ -\-~~distance \{\\ -\-~~~~name x\\ -\-~~~~group1 \{ atomNumbers 1 \}\\ -\-~~~~group2 \{ atomNumbers 50 \}\\ -\-~~\}\\ -\-~~distanceVec \{\\ -\-~~~~name r\\ -\-~~~~group1 \{ atomNumbers 10 11 12 \}\\ -\-~~~~group2 \{ atomNumbers 20 21 22 \}\\ -\-~~\}\\ -\}} +\-~~centerReference yes\\ +\-~~rotateReference yes\\ +\-~~fittingGroup \{\\ +\-~~~~\# define the frame by fitting the rest of the protein\\ +\-~~~~psfSegID PROT PROT\\ +\-~~~~atomNameResidueRange CA 1-40\\ +\-~~~~atomNameResidueRange CA 49-100\\ +\-~~\} \\ +\-~~refPositionsFile all.pdb \# can be the entire system\\ +\}\\} +} +\end{itemize} +The following two options have default values appropriate for the vast majority of applications, and are only provided to support rare, special cases. \begin{itemize} + \item % - \labelkey{colvar|customFunction} - \key - {customFunction}{% - \texttt{colvar}}{% - Compute colvar as a custom function of its components}{% - string}{% - Defines the colvar as a scalar expression of its colvar components. Multiple mentions can - be used to define a vector variable (as in the example above).} + \keydef + {enableFitGradients}{% + atom group}{% + Include the roto-translational contribution to colvar gradients}{% + boolean}{% + \texttt{on}}{% + When either \texttt{centerReference} or \texttt{rotateReference} is on, + the gradients of some colvars include terms proportional to +$\partial{}R/\partial\mathbf{x}_{i}$(rotational gradients) and +$\partial\mathbf{x}^{\mathrm{C}}/\partial\mathbf{x}_{i}$(translational gradients). + By default, these terms are calculated and included in the total gradients; + if this option is set to \texttt{off}, they are neglected. + In the case of a minimum RMSD component, this flag is automatically disabled + because the contributions of those derivatives to the gradients cancel out. +} \item % \keydef - {customFunctionType}{% - \texttt{colvar}}{% - Type of value returned by the scripted colvar}{% - string}{% - \texttt{scalar}}{% - With this flag, the user may specify whether the - colvar is a scalar or one of the following vector types: \texttt{vector3} - (a 3D vector), \texttt{unit\_vector3} (a normalized 3D vector), or - \texttt{unit\_quaternion} (a normalized quaternion), or \texttt{vector}. - Note that the scalar and vector cases are not necessary, as they are detected automatically.} + {enableForces}{% + atom group}{% + Apply forces from this colvar to this group}{% + boolean}{% + \texttt{on}}{% + If this option is \texttt{off}, no forces are applied the atoms in the group. + Other forces are not affected (i.e. those + from the MD engine, from other colvars, and other external forces). + For dummy atoms, this option is \texttt{off} by default. + } + \end{itemize} + + +\cvsubsec{Treatment of periodic boundary conditions.}{sec:colvar_atom_groups_wrapping} + +\cvnamdonly{ + In simulations with periodic boundary conditions, NAMD maintains + the coordinates of all the atoms within a molecule contiguous to + each other (i.e.~there are no spurious jumps'' in the molecular + bonds). The Colvars module relies on this when calculating a group's + center of geometry, but this condition may fail if the group spans + different molecules. In that case, writing the NAMD output and restart files + using \texttt{wrapAll} or \texttt{wrapWater} could produce wrong results + when a simulation run is continued from a previous one. + The user should then determine, according to which + type of colvars are being calculated, whether \texttt{wrapAll} or + \texttt{wrapWater} can be enabled. + + In general, internal coordinate wrapping by NAMD does not affect the calculation of colvars if each atom group satisfies one or more of the following: } +\cvlammpsonly{ +In simulations with periodic boundary conditions, many of the implemented colvar components rely on the fact that each position within a group of atoms is at the nearest periodic image from the center of geometry of the group itself. +However, due to the internal wrapping of individual atomic positions done by LAMMPS, this assumption is broken if the group straddles one of the unit cell's boundaries. +For this reason, within the Colvars module all coordinates are unwrapped by default to avoid discontinuities (see \texttt{unwrap} keyword in \ref{sec:colvars_mdengine_parameters}). +The user should determine whether maintaining the default value of \texttt{unwrap}, depending on the specifics of each system. +In general, internal coordinate wrapping by LAMMPS does not affect the calculation of colvars if each atom group satisfies one or more of the following: +} +\cvvmdonly{ + When periodic boundary conditions are defined, the Colvars module requires that the coordinates of each molecular fragment are contiguous, without jumps'' when a fragment is partially wrapped near a periodic boundary. + The Colvars module relies on this assumption when calculating a group's center of geometry, but the condition may fail if the group spans different molecules. + In general, coordinate wrapping does not affect the calculation of colvars if each atom group satisfies one or more of the following: +} -\cvscriptonly{ -\cvsubsec{Colvars as scripted functions of components}{sec:colvar_scripted} -When scripting is supported\cvnamdonly{ (default in NAMD)}\cvvmdonly{ (default in VMD)}, -a colvar may be defined as a scripted function of its components, -rather than a linear or polynomial combination. -When implementing generic functions of Cartesian coordinates rather -than functions of existing components, the \texttt{cartesian} component -may be particularly useful. -A scalar-valued scripted variable may be manually defined as periodic by providing -the keyword \texttt{period}, and the optional keyword \texttt{wrapAround}, with the -same meaning as in periodic components (see \ref{sec:cvc_periodic} for details). +\begin{enumerate} + \item[\emph{i)}] it is composed by only one atom; + \item[\emph{ii)}] it is used by a colvar component which does not make use of its center of geometry, but only of pairwise distances (\texttt{distanceInv}, \texttt{coordNum}, \texttt{hBond}, \texttt{alpha}, \texttt{dihedralPC}); + \item[\emph{iii)}] it is used by a colvar component that ignores the ill-defined Cartesian components of its center of mass (such as the$x$and$y$components of a membrane's center of mass modeled with \texttt{distanceZ})\cvnamdonly{; + \item[\emph{iv)}] it has all of its atoms within the same molecular fragment% +}. +\end{enumerate} +\cvvmdonly{If none of these conditions are met, wrapping may affect the calculation of collective variables: a possible solution is to use \texttt{pbc wrap} or \texttt{pbc unwrap} prior to processing a trajectory with the Colvars module.} -An example of elaborate scripted colvar is given in example 10, in the -form of path-based collective variables as defined by Branduardi et al\cite{Branduardi2007} -(\ref{sec:pathcv}). +\cvsubsec{Performance of a Colvars calculation based on group size.}{sec:colvar_atom_groups_scaling} + +In simulations performed with message-passing programs (such as NAMD or LAMMPS), the calculation of energy and forces is distributed (i.e., parallelized) across multiple nodes, as well as over the processor cores of each node. +When Colvars is enabled, certain atomic coordinates are collected on a single node, where the calculation of collective variables and of their biases is executed. +This means that for simulations over large numbers of nodes, a Colvars calculation may produce a significant overhead, coming from the costs of transmitting atomic coordinates to one node and of processing them. +\cvnamdonly{The latency-tolerant design and dynamic load balancing of NAMD may alleviate both factors, but a noticeable performance impact may be observed.} + +Performance can be improved in multiple ways: \begin{itemize} - \item % - \labelkey{colvar|scriptedFunction} - \key - {scriptedFunction}{% - \texttt{colvar}}{% - Compute colvar as a scripted function of its components}{% - string}{% - If this option is specified, the colvar will be computed as a - scripted function of the values of its components. - To that effect, the user should define two Tcl procedures: - \texttt{calc\_$<$scriptedFunction$>$} and \texttt{calc\_$<$scriptedFunction$>$\_gradient}, - both accepting as many parameters as the colvar has components. - Values of the components will be passed to those procedures in the - order defined by their sorted \texttt{name} strings. Note that if all - components are of the same type, their default names are sorted in the - order in which they are defined, so that names need only be specified - for combinations of components of different types. - \texttt{calc\_$<$scriptedFunction$>$} should return one value of - type$<$scriptedFunctionType$>$, corresponding to the colvar value. - \texttt{calc\_$<$scriptedFunction$>$\_gradient} should return a Tcl list - containing the derivatives of the function with respect to each - component. - If both the function and some of the components are vectors, the gradient - is really a Jacobian matrix that should be passed as a linear vector in - row-major order, i.e. for a function$f_i(x_j)$:$\nabla_x f_1 \nabla_x f_2 \cdots$. - } +\item The calculation of variables, components and biases can be distributed over the processor cores of the node where the Colvars module is executed. + Currently, an equal weight is assigned to each colvar, or to each component of those colvars that include more than one component. + The performance of simulations that use many colvars or components is improved automatically. + For simulations that use a single large colvar, it may be advisable to partition it in multiple components, which will be then distributed across the available cores. + \cvnamdonly{In NAMD, this feature is enabled in all binaries compiled using SMP builds of Charm++ with the CkLoop extension.} + \cvlammpsonly{In LAMMPS, this feature is supported automatically when LAMMPS is compiled with OpenMP support.} + If printed, the message SMP parallelism is available.'' indicates the availability of the option\cvvmdonly{ (will be supported in a future relase of VMD)}. + If available, the option is turned on by default, but may be disabled using the keyword \refkey{smp}{Colvars-global|smp} if required for debugging. -\item% - \keydef - {scriptedFunctionType}{% - \texttt{colvar}}{% - Type of value returned by the scripted colvar}{% - string}{% - \texttt{scalar}}{% - If a colvar is defined as a scripted function, its type is not constrained by - the types of its components. With this flag, the user may specify whether the - colvar is a scalar or one of the following vector types: \texttt{vector3} - (a 3D vector), \texttt{unit\_vector3} (a normalized 3D vector), or - \texttt{unit\_quaternion} (a normalized quaternion), or \texttt{vector} - (a vector whose size is specified by \texttt{scriptedFunctionVectorSize}). - Non-scalar values should be passed as space-separated lists.} +\cvnamdonly{ + % Use the following command to identify them: + % grep -B10 'provide(f_cvc_com_based' * |grep '\:\:'|grep '(std::string const &conf)' +\item NAMD also offers a parallelized calculation of the centers of mass of groups of atoms. + This option is on by default for all components that are simple functions of centers of mass, and is controlled by the keyword \refkey{scalable}{sec:cvc_common}. + When supported, the message Will enable scalable calculation for group \ldots'' is printed for each group. +} - \item % - \key - {scriptedFunctionVectorSize}{% - \texttt{colvar}}{% - Dimension of the vector value of a scripted colvar}{% - positive integer}{% - This parameter is only valid when \texttt{scriptedFunctionType} is - set to \texttt{vector}. It defines the vector length of the colvar value - returned by the function.} +\item As a general rule, the size of atom groups should be kept relatively small (up to a few thousands of atoms, depending on the size of the entire system in comparison). +To gain an estimate of the computational cost of a large colvar, one can use a test calculation of the same colvar in VMD (hint: use the \texttt{time} Tcl command to measure the cost of running \texttt{cv update}). \end{itemize} -} % \cvscriptonly \cvsec{Biasing and analysis methods}{sec:colvarbias} @@ -3200,7 +3356,7 @@ group used in a \texttt{distance} component. ABF depends on parameters from collective variables to define the grid on which free energy gradients are computed. In the direction of each colvar, the grid ranges from \texttt{lowerBoundary} to \texttt{upperBoundary}, and the bin width (grid spacing) -is set by the \texttt{width} parameter (see~\ref{sec:colvar_general}). +is set by the \refkey{width}{colvar|width} parameter. The following specific parameters can be set in the ABF configuration block: \begin{itemize} @@ -3657,7 +3813,7 @@ The parameters$\delta{}t$and$\sigma_{\xi}$specified by the optional keywords Relative width of a Gaussian hill with respect to the colvar's width}{% positive decimal}{%$\sqrt{2\pi}/2$}{% - The Gaussian width along each colvar,$2\sigma_{\xi_{i}}$, is set as the product between this number and the colvar's parameter$w_{i}$given by \texttt{width} (see \ref{sec:colvar_general}); such product is printed in the standard output of \MDENGINE{}. + The Gaussian width along each colvar,$2\sigma_{\xi_{i}}$, is set as the product between this number and the colvar's parameter$w_{i}$given by \refkey{width}{colvar|width}; such product is printed in the standard output of \MDENGINE{}. The default value of this number gives a Gaussian hill function whose volume is equal to the product of$\prod_{i}w_{i}$, the volume of one \texttt{histogram} bin (see \ref{sec:colvarbias_histogram}), and$W$. \textbf{Tip:} \emph{use this property to estimate the fraction of colvar space covered by the Gaussian bias within a given simulation time.} When \texttt{useGrids} is on, the default value also gives acceptable discretization errors~\cite{Fiorin2013}: for smoother visualization, this parameter may be increased and the \texttt{width}$w_{i}$decreased in the same proportion. @@ -3693,7 +3849,7 @@ The following two options allow to control this behavior and to visually track s \end{itemize} -\textbf{Note:} when Gaussian hills are deposited near \texttt{lowerBoundary} or \texttt{upperBoundary} (see \ref{sec:colvar_general}) and interpolating grids are used (default behavior), their truncation can give rise to accumulating errors. +\textbf{Note:} when Gaussian hills are deposited near the \refkey{lowerBoundary}{colvar|lowerBoundary} or \refkey{upperBoundary}{colvar|upperBoundary} and interpolating grids are used (default behavior), their truncation can give rise to accumulating errors. In these cases, as a measure of fault-tolerance all Gaussian hills near the boundaries are included in the output state file, and are recalculated analytically whenever the colvar falls outside the grid's boundaries. (Such measure protects the accuracy of the calculation, and can only be disabled by \texttt{hardLowerBoundary} or \texttt{hardUpperBoundary}.) To avoid gradual loss of performance and growth of the state file, either one of the following solutions is recommended: @@ -3953,7 +4109,7 @@ Note that this differs from harmonic bond and angle potentials in common force fields, where the factor of one half is typically omitted, resulting in a non-standard definition of the force constant. -The formula above includes the characteristic length scale$w_{\xi}$of the colvar$\xi$(keyword \texttt{width}, see \ref{sec:colvar_general}) to allow the definition of a multi-dimensional restraint with a unified force constant: +The formula above includes the characteristic length scale$w_{\xi}$of the colvar$\xi\$ (keyword \refkey{lowerBoundary}{colvar|lowerBoundary}) to allow the definition of a multi-dimensional restraint with a unified force constant:

\label{eq:colvarbias_harmonic_multi}
V(\xi_{1}, \ldots, \xi_{M}) = \frac{1}{2} k \sum_{i=1}^{M} \left(\frac{\xi_{i} - \xi_0}{w_{\xi}}\right)^2
@@ -4623,14 +4779,16 @@ The list of options is as follows:

\cvscriptonly{%
-\cvsubsec{Scripted biases}{sec:scripted_biases}
+\cvsubsec{Defining scripted biases}{sec:scripted_biases}

Rather than using the biasing methods described above, it is possible to apply biases
provided at run time as a Tcl script\cvnamdonly{, in the spirit of \texttt{TclForces}}.
\cvvmdonly{This option, also available in NAMD, can be useful to test a new algorithm to be used in a MD simulation.}

\begin{itemize}
+
\item %
+  \labelkey{Colvars-global|scriptedColvarForces}
\keydef
{scriptedColvarForces}{%
global}{%
@@ -4641,39 +4799,39 @@ provided at run time as a Tcl script\cvnamdonly{, in the spirit of \texttt{TclFo
accepting one parameter should be defined by the user. It is executed
at each timestep, with the current step number as parameter, between the
calculation of colvars and the application of bias forces. This procedure
-    may use\cvnamdonly{ the scripting interface (see \ref{sec:scripting})}\cvvmdonly{ the \texttt{cv} command} to access
+    may use\cvnamdonly{ the scripting interface (see \ref{sec:cv_scripting})}\cvvmdonly{ the \texttt{cv} command} to access
the values of colvars and apply forces on them, effectively defining custom
collective variable biases.}
\end{itemize}
-}

-\cvnamdonly{\cvsec{Colvars scripting}{sec:scripting}

-This interface is particularly useful to implement custom biases as scripted colvar forces.
-See the \texttt{scriptedColvarForces} option in \ref{sec:scripted_biases}.
-Note that scripting commands may not be used directly in the NAMD configuration file before the first
-\texttt{run} or \texttt{minimize} statement. They may be used either within the callback procedures
-(e.g. \texttt{calc\_colvar\_forces}) or in the NAMD config file after a \texttt{run} or
-\texttt{minimize} statement.
+\cvsubsec{Performance of scripted biases}{sec:cv_scripting_opt}
+
+If concurrent computation over multiple threads is available (this is indicated by the message `SMP parallelism is available.'' printed at initialization time), it is useful to take advantage of the scripting interface to combine many components, all computed in parallel, into a single variable.

-\cvnamdugonly{\input{ug_colvars-cv.tex}}
-\cvrefmanonly{\input{colvars-cv.tex}}
+The default SMP schedule is the following:
+\begin{enumerate}
+\item distribute the computation of all components across available threads;
+\item on a single thread, collect the results of multi-component variables using polynomial combinations (see \ref{sec:cvc_superp})\cvleptononly{, or custom functions (see \ref{sec:colvar_custom_function})}\cvscriptonly{, or scripted functions (see \ref{sec:colvar_scripted})};
+\item distribute the computation of all biases across available threads;
+\item compute on a single thread any scripted biases implemented via the keyword \refkey{scriptedColvarForces}.
+\item communicate on a single thread forces to \MDENGINE{}.
+\end{enumerate}

-The following configuration options can modify the behavior of the scripting interface for optimization purposes:
+The following options allow to fine-tune this schedule:
\begin{itemize}
\item %
\keydef
{scriptingAfterBiases}{%
global}{%
-    Scripted colvar forces need updated biases}{%
+    Scripted colvar forces need updated biases?}{%
boolean}{%
\texttt{on}}{%
-    This flag specifies that the \texttt{calc\_colvar\_forces} procedure, when defined,
-    is executed only after all biases have been updated.
+    This flag specifies that the \texttt{calc\_colvar\_forces} procedure (last step in the list above) is executed only after all biases have been updated (next-to-last step)
For example, this allows using the energy of a restraint bias, or the force applied on a colvar,
to calculate additional scripted forces, such as boundary constraints.
-    When this flag is set to \texttt{off}, it is assumed that only the values of their colvars
-    (but not their energy or forces) will be used by \texttt{calc\_colvar\_forces}:
+    When this flag is set to \texttt{off}, it is assumed that only the values of the variables
+    (but not the energy of the biases or applied forces) will be used by \texttt{calc\_colvar\_forces}:
this can be used to schedule the calculation of scripted forces and biases concurrently
to increase performance.}
\end{itemize}
@@ -4684,16 +4842,12 @@ The following configuration options can modify the behavior of the scripting int
\cvsec{Syntax changes from older versions}{sec:colvars_config_changes}

The following is a list of syntax changes in Colvars since its first release.
-Many of the older keywords are still recognized by the current code: the primary purpose of this list is to make you aware of improvements that affect old configuration files.
+Many of the older keywords are still recognized by the current code, thanks to specific compatibility code.
+\emph{This is not a list of new features:} its primary purpose is to make you aware of those improvements that affect the use of old configuration files with new versions of the code.

-The syntax of the latest version of Colvars can always be accessed at:\\
-\cvnamdonly{\cvurl{https://colvars.github.io/colvars-refman-namd/colvars-refman-namd.html}}
-\cvlammpsonly{\cvurl{https://colvars.github.io/colvars-refman-lammps/colvars-refman-lammps.html}}
-\cvvmdonly{\cvurl{https://colvars.github.io/colvars-refman-vmd/colvars-refman-vmd.html}}
-
-If you are using any of the NAMD and VMD tutorials:\\
+\noindent\textbf{Note:} if you are using any of the NAMD and VMD tutorials:\\
\cvurl{https://www.ks.uiuc.edu/Training/Tutorials/}\\
-please be aware that \emph{several of these tutorials are not actively maintained}: for those cases, the following list will help you reconcile any inconsistencies.
+please be aware that \emph{several of these tutorials are not actively maintained}: for those cases, this list will help you reconcile any inconsistencies.

\begin{itemize}

@@ -4714,3 +4868,8 @@ please be aware that \emph{several of these tutorials are not actively maintaine
When using a metadynamics bias, it is recommended to set the two walls \emph{within} \refkey{lowerBoundary}{colvar|lowerBoundary} and \refkey{upperBoundary}{colvar|upperBoundary}.  This guarantees that the tails of each Gaussian hill are accounted in the region between the grid boundaries and the wall potentials.  See also \refkey{expandBoundaries}{colvar|expandBoundaries} for an automatic definition of the PMF grid boundaries.

\end{itemize}
+
+\noindent{}Up-to-date documentation can always be accessed at:\\
+\cvnamdonly{\cvurl{https://colvars.github.io/colvars-refman-namd/colvars-refman-namd.html}}
+\cvlammpsonly{\cvurl{https://colvars.github.io/colvars-refman-lammps/colvars-refman-lammps.html}}
+\cvvmdonly{\cvurl{https://colvars.github.io/colvars-refman-vmd/colvars-refman-vmd.html}}