Added new documentation for code in symfunc module

src/symfunc/AngularTetra.cpp

@@ -28,7 +28,27 @@
 /*
 Calculate the angular tetra CV
 
-\par Examples
+This shortcut calculates a [symmetry function](https://www.plumed-tutorials.org/lessons/23/001/data/SymmetryFunction.html). The particular function that is being 
+evaluated for the coordination sphere here is as follows:
+
+$$
+s_i = 1 - \frac{3}{8}\sum_{j=2}^4 \sum_{k=1}^{j-1} \left( \cos(\theta_{jik} + \frac{1}{2} \right)^2
+$$
+
+In this expression the 4 atoms in the sums over $j$ and $k$ are the four atoms that are nearest to atom $i$.  $\theta_{jik}$ is the angle between the vectors connecting 
+atoms $i$ and $j$ and the vector connecting atoms $i$ and $k$.  The CV is large if the four atoms nearest atom $i$ are arranged on the vertices of a regular tetrahedron 
+and small otherwise.  The following example shows how you can use this action to measure the degree of tetrahedral order in a system.
+
+```plumed
+# Calculate a vector that contains 64 values for the symmetry function.
+# Sum the elements of the vector and calculate the mean value on the atoms from this sum.
+acv: TETRA_ANGULAR SPECIES=1-64 SUM MEAN
+# Print out the positions of the 64 atoms for which the symmetry function was calculated
+# to an xyz file along with the values of the symmetry function
+DUMPATOMS ATOMS=1-64 ARG=acv FILE=mcolv.xyz 
+# Print out the average value of the symmetry function and the sum of all the symmetry functions
+PRINT ARG=acv_sum,acv_mean FILE=colvar 
+```
 
 */
 //+ENDPLUMEDOC

src/symfunc/AtomicSMAC.cpp

@@ -30,7 +30,30 @@
 /*
 Calculate the atomic smac CV
 
-\par Examples
+This shortcut offers another example of a [symmetry function](https://www.plumed-tutorials.org/lessons/23/001/data/SymmetryFunction.html).  
+This action was inspired by [SMAC](SMAC.md) and the distribution of angles in the first coordination sphere around an atom to determine if the 
+environment around the atom is ordered. For atom $i$ the symmetry function is calculated using:
+
+$$
+s_i = [ 1 - \gamma(c_i) ] \frac{ \sum_j \sum{k \ne j} \sigma(r_{ij})\sigma(r_{ik}) \sum_n G\left( \frac{ \theta_{jik} - \phi_n}{b_n}\right) }{\sum{k \ne j} \sigma(r_{ij})\sigma(r_{ik})} \qquad \textrm{where} \qquad c_i = \sum_j \sigma(r_{ij}) 
+$$ 
+
+In this expression $r_{ij}$ is the distance between atom $i$ and atom $j$ and $\sigma$ is a switching function that acts upon this distance.   $c_i$ is thus the number of atoms that are within 
+a certain cutoff of atom $i$ and $\gamma$ is another switching function that acts upon this quantity.  This switching function ensures that the symmetry function is zero for atoms that are 
+regions where the density is low.  $\theta_{jik}$ is the angle between the vector connecting atoms $i$ and $j$ and the vector connecting atoms $i$ and $k$.  This angle is the argument for the 
+set of Gaussian kernel functions, $G$, that are centered on $\phi_n$ and that have bandwidths of $b_n$.  The function above is thus determining if the angles between the bonds in the first coordination 
+sphere around atom $i$ are similar to the $\phi_n$ values that have been specified by the user or not.
+
+The following example demonstrates how this symmetry function can be used in practise.
+
+```plumed
+smac: ATOMIC_SMAC SPECIES=1-64 KERNEL1={GAUSSIAN CENTER=pi/2 SIGMA=1.0} KERNEL2={GAUSSIAN CENTER=pi/4 SIGMA=1.0} SWITCH_COORD={EXP R_0=4.0} SWITCH={RATIONAL R_0=2.0 D_0=2.0} SUM
+PRINT ARG=smac.* FILE=colvar 
+``` 
+
+The input above would calculate 64 instances of $s_i$ using the formula above.  In each of these two Gaussian Kernels are used in the sum over $n$.  The parameters for the switching function 
+$\sigma$ are specified using the SWITCH keyword, while the parameters for $\gamma$ are specified using SWITCH_COORD.  As you can see if you expand the shortcut in the input above, the 64 values 
+for $s_i$ are stored in a vector.  All the elements of this vector are then added together to produce the single quantity that is output in the colvar file. 
 
 */
 //+ENDPLUMEDOC
@@ -48,12 +71,10 @@ PLUMED_REGISTER_ACTION(AtomicSMAC,"ATOMIC_SMAC")
 
 void AtomicSMAC::registerKeywords(Keywords& keys) {
   ActionShortcut::registerKeywords( keys );
-  keys.add("optional","SPECIES","");
-  keys.add("optional","SPECIESA","");
-  keys.add("optional","SPECIESB","");
-  keys.add("optional","SWITCH","This keyword is used if you want to employ an alternative to the continuous swiching function defined above. "
-           "The following provides information on the \\ref switchingfunction that are available. "
-           "When this keyword is present you no longer need the NN, MM, D_0 and R_0 keywords.");
+  keys.add("atoms-3","SPECIES","the list of atoms for which the symmetry function is being calculated and the atoms that can be in the environments");
+  keys.add("atoms-4","SPECIESA","the list of atoms for which the symmetry function is being calculated.  This keyword must be used in conjunction with SPECIESB, which specifies the atoms that are in the environment.");
+  keys.add("atoms-4","SPECIESB","the list of atoms that can be in the environments of each of the atoms for which the symmetry function is being calculated.  This keyword must be used in conjunction with SPECIESA, which specifies the atoms for which the symmetry function is being calculated.");
+  keys.add("optional","SWITCH","the switching function that it used in the construction of the contact matrix");
   keys.add("numbered","KERNEL","The kernels used in the function of the angle");
   keys.add("optional","SWITCH_COORD","This keyword is used to define the coordination switching function.");
   keys.reset_style("KERNEL","optional");

src/symfunc/CoordShellVectorFunction.cpp

@@ -37,7 +37,57 @@ namespace symfunc {
 /*
 Calculate an arbitrary function of all the bond vectors in the first coordination sphere of an atom
 
-\par Examples
+This shortcut allows you to calculate the sum for an arbitrary function of the bond vectors that connect an atom to each of its neighbours.
+In other words, this action allows you to compute the following:
+
+$$
+s_i = \sum_{i \ne j} \sigma(r_{ij}) f(x_{ij}, y_{ij}, z_{ij}, r_{ij}) )
+$$
+
+In this expression, $x_{ij}, y_{ij}, z_{ij}$ are the components of the vector connecting atoms $i$ and $j$ and $r_{ij}$ is the magnitude of this vector.
+$\sigma(r_{ij})$ is then a switching function that ensures that the aritrary function $f$ is only evaluated for if atom $j$ is within a certain cutoff distance  
+of atom $i$.
+
+Below you can see a simple example that shows how this action can be used in practise.
+
+```plumed
+cv: COORDINATION_SHELL_FUNCTION SPECIES=1-64 SWITCH={RATIONAL D_0=3.0 R_0=1.5} FUNCTION=((x+y+z)/r)^3+((x-y-z)/r)^3+((-x+y-z)/r)^3+((-x-y+z)/r)^3
+PRINT ARG=cv FILE=colvar
+```
+
+The above input calculates 64 $s_i$ values - one $s_i$ values for each of the atoms specified using the SPECIES keyword.  These 64 numbers are then output to a file.
+The function of x, y, z and r to be evaluated is specified using the FUNCTION keyword.  Obviously, if your function does not depend on all four of these variables
+they can be excluded from your function.
+
+As discussed in [this paper](https://journals.aps.org/prb/abstract/10.1103/PhysRevB.92.180102) it is sometimes useful to rotate the bond vectors before computing the
+arbitrary function $f$ in the above expression.  To perform such rotations you use the PHI, THETA and PSI keywords.  The $x_{ij}, y_{ij}$ and $z_{ij}$ values that enter $f$ in the 
+expression above are then calculated as:
+
+$$
+\left( 
+\begin{matrix}
+x_{ij} \\
+y_{ij} \\
+z_{ij} 
+\end{matrix}
+\right) = 
+\left(
+\begin{matrix}
+\cos(\psi)\cos(\phi) - \cos(\theta)\sin(\phi)\sin(\psi) & \cos(\psi)*\sin(\phi)+\cos(\theta)*\cos(\phi)*\sin(\psi) & \sin(\psi)*sin(\theta) \\
+-\sin(\psi)*\cos(\phi)-\cos(\theta)*\sin(\phi)*\cos(\psi) & -\sin(\psi)*\sin(\phi)+\cos(\theta)*\cos(\phi)*std::cos(\psi), & \cos(\psi)*\sin(\theta) \\
+\sin(\theta)*\sin(\phi) & \sin(\theta)*\cos(\phi) & \cos(\theta)
+\end{matrix}
+\left( 
+\begin{matrix}
+x_{ij}' \\
+y_{ij}' \\
+z_{ij}' 
+\end{matrix}
+\right)
+$$
+
+$x_{ij}', y_{ij}'$ and $z_{ij}'$ in this expression are the bond vectors that are calculated in the lab frame.  The matrix in the above expression is thus just a 
+[rotation matrix](https://en.wikipedia.org/wiki/Rotation_matrix) that converts the lab frame to some frame of interest.
 
 
 */
@@ -47,8 +97,30 @@ Calculate an arbitrary function of all the bond vectors in the first coordinatio
 /*
 Calculate an arbitrary function of all the bond vectors in the first coordination sphere of an atom and take an average
 
-\par Examples
+This shortcut allows you to calculate the average for an arbitrary function of the bond vectors that connect an atom to each of its neighbours.
+In other words, this action allows you to compute the following:
+    
+$$
+s_i = \frac{\sum_{i \ne j} \sigma(r_{ij}) f(x_{ij}, y_{ij}, z_{ij}, r_{ij}) )}{\sum_{i \ne j} \sigma(r_{ij})}
+$$
 
+In this expression, $x_{ij}, y_{ij}, z_{ij}$ are the components of the vector connecting atoms $i$ and $j$ and $r_{ij}$ is the magnitude of this vector.
+$\sigma(r_{ij})$ is then a switching function that ensures that the aritrary function $f$ is only evaluated for if atom $j$ is within a certain cutoff distance  
+of atom $i$.
+
+Below you can see a simple example that shows how this action can be used in practise.
+
+```plumed
+cv: COORDINATION_SHELL_AVERAGE SPECIES=1-64 SWITCH={RATIONAL D_0=3.0 R_0=1.5} FUNCTION=((x+y+z)/r)^3+((x-y-z)/r)^3+((-x+y-z)/r)^3+((-x-y+z)/r)^3
+PRINT ARG=cv FILE=colvar
+```
+
+The above input calculates 64 $s_i$ values - one $s_i$ values for each of the atoms specified using the SPECIES keyword.  These 64 numbers are then output to a file.
+The function of x, y, z and r to be evaluated is specified using the FUNCTION keyword.  Obviously, if your function does not depend on all four of these variables
+they can be excluded from your function.
+
+Notice that you can you can rotate the bond vectors before computing the
+arbitrary function $f$ in the above expression as is discussed in the documentation for [COORDINATION_SHELL_FUNCTION](COORDINATION_SHELL_FUNCTION.md)
 
 */
 //+ENDPLUMEDOC
@@ -57,35 +129,38 @@ Calculate an arbitrary function of all the bond vectors in the first coordinatio
 /*
 Calculate whether or not the coordination spheres of atoms are arranged as they would be in a simple cubic structure.
 
-We can measure how similar the environment around atom \f$i\f$ is to a simple cubic structure is by evaluating
+This shortcut is an example of a [COORDINATION_SHELL_AVERAGE](COORDINATION_SHELL_AVERAGE.md), 
+which we can use to measure how similar the environment around atom $i$ is to a simple cubic structure.  We perform this comparison by evaluating
 the following quantity:
 
-\f[
+$$
 s_i = \frac{ \sum_{i \ne j} \sigma(r_{ij}) \left[ \frac{ x_{ij}^4 + y_{ij}^4 + z_{ij}^4 }{r_{ij}^4} \right] }{ \sum_{i \ne j} \sigma(r_{ij}) }
-\f]
+$$
 
-In this expression \f$x_{ij}\f$, \f$y_{ij}\f$ and \f$z_{ij}\f$ are the \f$x\f$, \f$y\f$ and \f$z\f$ components of the vector connecting atom \f$i\f$ to
-atom \f$j\f$ and \f$r_{ij}\f$ is the magnitude of this vector.  \f$\sigma(r_{ij})\f$ is a \ref switchingfunction that acts on the distance between atom \f$i\f$ and atom \f$j\f$ and its inclusion in the numerator and the denominator of the above expression as well as the fact that we are summing
+In this expression $x_{ij}$, $y_{ij}$ and $z_{ij}$ are the $x$, $y$ and $z$ components of the vector connecting atom $i$ to
+atom $j$ and $r_{ij}$ is the magnitude of this vector.  $\sigma(r_{ij})$ is a switching function that acts on the distance between atom $i$ and atom $j$ and its inclusion in the numerator and the denominator of the above expression as well as the fact that we are summing
 over all of the other atoms in the system ensures that we are calculating an average
-of the function of \f$x_{ij}\f$, \f$y_{ij}\f$ and \f$z_{ij}\f$ for the atoms in the first coordination sphere around atom \f$i\f$.
+of the function of $x_{ij}$, $y_{ij}$ and $z_{ij}$ for the atoms in the first coordination sphere around atom $i$.
 This quantity is once again a multicolvar so you can compute it for multiple atoms using a single PLUMED action and then compute
-the average value for the atoms in your system, the number of atoms that have an \f$s_i\f$ value that is more that some target and
+the average value for the atoms in your system, the number of atoms that have an $s_i$ value that is more that some target and
 so on.  Notice also that you can rotate the reference frame if you are using a non-standard unit cell.
 
-
-\par Examples
-
 The following input tells plumed to calculate the simple cubic parameter for the atoms 1-100 with themselves.
 The mean value is then calculated.
-\plumedfile
+
+```plumed
 SIMPLECUBIC SPECIES=1-100 R_0=1.0 MEAN
-\endplumedfile
+```
 
 The following input tells plumed to look at the ways atoms 1-100 are within 3.0 are arranged about atoms
 from 101-110.  The number of simple cubic parameters that are greater than 0.8 is then output
-\plumedfile
+
+```plumed
 SIMPLECUBIC SPECIESA=101-110 SPECIESB=1-100 R_0=3.0 MORE_THAN={RATIONAL R_0=0.8 NN=6 MM=12 D_0=0}
-\endplumedfile
+```
+
+Notice that you can you can rotate the bond vectors before computing the
+function in the above expression as is discussed in the documentation for [COORDINATION_SHELL_FUNCTION](COORDINATION_SHELL_FUNCTION.md)
 
 */
 //+ENDPLUMEDOC
@@ -94,40 +169,42 @@ SIMPLECUBIC SPECIESA=101-110 SPECIESB=1-100 R_0=3.0 MORE_THAN={RATIONAL R_0=0.8
 /*
 Calculate the degree to which the environment about ions has a tetrahedral order.
 
-We can measure the degree to which the atoms in the first coordination shell around any atom, \f$i\f$ is
-is arranged like a tetrahedron using the following function.
+This shortcut is an example of a [COORDINATION_SHELL_AVERAGE](COORDINATION_SHELL_AVERAGE.md),
+which we can use to measure the degree to which the atoms in the first coordination shell around any atom, $i$ is
+is arranged like a tetrahedron.  We perform this comparison by evaluating the following function.
 
-\f[
+$$
  s(i) = \frac{1}{\sum_j \sigma( r_{ij} )} \sum_j \sigma( r_{ij} )\left[ \frac{(x_{ij} + y_{ij} + z_{ij})^3}{r_{ij}^3} +
                                                                         \frac{(x_{ij} - y_{ij} - z_{ij})^3}{r_{ij}^3} +
                                                                         \frac{(-x_{ij} + y_{ij} - z_{ij})^3}{r_{ij}^3} +
                                                                         \frac{(-x_{ij} - y_{ij} + z_{ij})^3}{r_{ij}^3} \right]
-\f]
+$$
 
-Here \f$r_{ij}\f$ is the magnitude of the vector connecting atom \f$i\f$ to atom \f$j\f$ and \f$x_{ij}\f$, \f$y_{ij}\f$ and \f$z_{ij}\f$
-are its three components.  The function  \f$\sigma( r_{ij} )\f$ is a \ref switchingfunction that acts on the distance between
-atoms \f$i\f$ and \f$j\f$.  The parameters of this function should be set so that the function is equal to one
-when atom \f$j\f$ is in the first coordination sphere of atom \f$i\f$ and is zero otherwise.
-
-\par Examples
+Here $r_{ij}$ is the magnitude of the vector connecting atom $i$ to atom $j$ and $x_{ij}$, $y_{ij}$ and $z_{ij}$
+are its three components.  The function  $\sigma( r_{ij} )$ is a switching function that acts on the distance between
+atoms $i$ and $j$.  The parameters of this function should be set so that the function is equal to one
+when atom $j$ is in the first coordination sphere of atom $i$ and is zero otherwise.
 
 The following command calculates the average value of the TETRAHEDRAL parameter for a set of 64 atoms all of the same type
 and outputs this quantity to a file called colvar.
 
-\plumedfile
+```plumed
 tt: TETRAHEDRAL SPECIES=1-64 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN
 PRINT ARG=tt.mean FILE=colvar
-\endplumedfile
+```
 
 The following command calculates the number of TETRAHEDRAL parameters that are greater than 0.8 in a set of 10 atoms.
 In this calculation it is assumed that there are two atom types A and B and that the first coordination sphere of the
 10 atoms of type A contains atoms of type B.  The formula above is thus calculated for ten different A atoms and within
-it the sum over \f$j\f$ runs over 40 atoms of type B that could be in the first coordination sphere.
+it the sum over $j$ runs over 40 atoms of type B that could be in the first coordination sphere.
 
-\plumedfile
+```plumed
 tt: TETRAHEDRAL SPECIESA=1-10 SPECIESB=11-40 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MORE_THAN={RATIONAL R_0=0.8}
 PRINT ARG=tt.* FILE=colvar
-\endplumedfile
+```
+
+Notice that you can you can rotate the bond vectors before computing the
+function in the above expression as is discussed in the documentation for [COORDINATION_SHELL_FUNCTION](COORDINATION_SHELL_FUNCTION.md)
 
 */
 //+ENDPLUMEDOC
@@ -155,6 +232,8 @@ void CoordShellVectorFunction::registerKeywords( Keywords& keys ) {
   keys.setValueDescription("vector","the symmetry function for each of the specified atoms");
   keys.needsAction("CONTACT_MATRIX"); keys.needsAction("FCCUBIC_FUNC"); keys.needsAction("CUSTOM");
   keys.needsAction("ONES"); keys.needsAction("MATRIX_VECTOR_PRODUCT");
+  keys.addDOI("10.1103/PhysRevB.81.125416"); keys.addDOI("10.1103/PhysRevB.92.180102");
+  keys.addDOI("10.1063/1.4997180"); keys.addDOI("10.1063/1.5134461");
 }
 
 CoordShellVectorFunction::CoordShellVectorFunction(const ActionOptions& ao):