From 8c0b1386be6c5afbb5238f1dc26a87dbb3d55cce Mon Sep 17 00:00:00 2001
From: Giovanni Bussi <giovanni.bussi@gmail.com>
Date: Mon, 24 Apr 2017 11:12:10 +0200
Subject: [PATCH] Syntax highlight in doc

I replaced
verbatim
endverbatim
with
plumedfile
endplumedfile

in (hopefully) all the places where it was appropriate.
---
 src/adjmat/ClusterDiameter.cpp                |  4 +-
 src/adjmat/ClusterDistribution.cpp            |  4 +-
 src/adjmat/ClusterProperties.cpp              |  4 +-
 src/adjmat/ClusterSize.cpp                    |  4 +-
 src/adjmat/ClusterWithSurface.cpp             |  4 +-
 src/adjmat/ContactAlignedMatrix.cpp           |  4 +-
 src/adjmat/ContactMatrix.cpp                  |  4 +-
 src/adjmat/DFSClustering.cpp                  |  4 +-
 src/adjmat/HbondMatrix.cpp                    |  4 +-
 src/adjmat/MatrixColumnSums.cpp               |  8 +--
 src/adjmat/MatrixRowSums.cpp                  |  8 +--
 src/adjmat/OutputCluster.cpp                  |  4 +-
 src/adjmat/SMACMatrix.cpp                     |  4 +-
 src/adjmat/Sprint.cpp                         |  8 +--
 src/analysis/Average.cpp                      | 12 ++--
 .../ClassicalMultiDimensionalScaling.cpp      |  4 +-
 src/analysis/Commit.cpp                       |  4 +-
 src/analysis/Histogram.cpp                    | 16 ++---
 src/analysis/PCA.cpp                          |  8 +--
 src/bias/ABMD.cpp                             |  5 +-
 src/bias/BiasValue.cpp                        | 11 ++--
 src/bias/ExtendedLagrangian.cpp               | 10 ++-
 src/bias/External.cpp                         |  9 ++-
 src/bias/LWalls.cpp                           |  4 +-
 src/bias/MetaD.cpp                            | 36 +++++-----
 src/bias/Metainference.cpp                    | 10 ++-
 src/bias/MovingRestraint.cpp                  | 12 ++--
 src/bias/PBMetaD.cpp                          | 16 ++---
 src/bias/Restraint.cpp                        |  5 +-
 src/bias/ReweightBias.cpp                     |  4 +-
 src/bias/ReweightMetad.cpp                    |  5 +-
 src/bias/UWalls.cpp                           |  5 +-
 src/colvar/Angle.cpp                          |  5 +-
 src/colvar/CS2Backbone.cpp                    |  8 +--
 src/colvar/Cell.cpp                           |  5 +-
 src/colvar/Constant.cpp                       |  9 ++-
 src/colvar/ContactMap.cpp                     | 15 ++---
 src/colvar/Coordination.cpp                   | 13 ++--
 src/colvar/DHEnergy.cpp                       |  5 +-
 src/colvar/DRMSD.cpp                          |  8 +--
 src/colvar/Dimer.cpp                          | 20 +++---
 src/colvar/Dipole.cpp                         |  5 +-
 src/colvar/Distance.cpp                       | 23 +++----
 src/colvar/ERMSD.cpp                          |  4 +-
 src/colvar/Energy.cpp                         |  5 +-
 src/colvar/Fake.cpp                           |  5 +-
 src/colvar/FretEfficiency.cpp                 | 10 ++-
 src/colvar/Gyration.cpp                       |  5 +-
 src/colvar/Implicit.cpp                       |  5 +-
 src/colvar/Jcoupling.cpp                      |  8 +--
 src/colvar/MultiRMSD.cpp                      |  8 +--
 src/colvar/NOE.cpp                            |  5 +-
 src/colvar/PCARMSD.cpp                        |  4 +-
 src/colvar/PRE.cpp                            |  4 +-
 src/colvar/PathMSD.cpp                        |  4 +-
 src/colvar/Position.cpp                       |  6 +-
 src/colvar/PropertyMap.cpp                    |  4 +-
 src/colvar/Puckering.cpp                      |  4 +-
 src/colvar/RDC.cpp                            |  5 +-
 src/colvar/RMSD.cpp                           |  4 +-
 src/colvar/Template.cpp                       |  4 +-
 src/colvar/Torsion.cpp                        |  8 +--
 src/colvar/Volume.cpp                         |  5 +-
 src/crystallization/Fccubic.cpp               |  4 +-
 src/crystallization/Gradient.cpp              |  8 +--
 .../InterMolecularTorsions.cpp                |  4 +-
 src/crystallization/MoleculeOrientation.cpp   |  4 +-
 src/crystallization/Q3.cpp                    | 24 +++----
 src/crystallization/Q4.cpp                    | 24 +++----
 src/crystallization/Q6.cpp                    | 24 +++----
 src/crystallization/SMAC.cpp                  |  8 +--
 src/crystallization/SimpleCubic.cpp           |  8 +--
 src/crystallization/Tetrahedral.cpp           |  8 +--
 src/eds/EDS.cpp                               | 12 ++--
 src/function/Combine.cpp                      |  8 +--
 src/function/Ensemble.cpp                     |  5 +-
 src/function/FuncPathMSD.cpp                  |  8 +--
 src/function/LocalEnsemble.cpp                |  4 +-
 src/function/Matheval.cpp                     | 16 ++---
 src/function/Piecewise.cpp                    |  5 +-
 src/function/Sort.cpp                         |  5 +-
 src/function/Stats.cpp                        |  4 +-
 src/function/Target.cpp                       |  4 +-
 src/generic/Debug.cpp                         |  4 +-
 src/generic/DumpAtoms.cpp                     | 10 ++-
 src/generic/DumpDerivatives.cpp               |  4 +-
 src/generic/DumpForces.cpp                    |  6 +-
 src/generic/DumpMassCharge.cpp                |  9 ++-
 src/generic/DumpProjections.cpp               |  4 +-
 src/generic/EffectiveEnergyDrift.cpp          |  8 +--
 src/generic/FitToTemplate.cpp                 |  5 +-
 src/generic/Flush.cpp                         |  4 +-
 src/generic/Group.cpp                         | 24 +++----
 src/generic/Include.cpp                       | 12 ++--
 src/generic/Print.cpp                         |  5 +-
 src/generic/RandomExchanges.cpp               | 12 ++--
 src/generic/Read.cpp                          |  4 +-
 src/generic/ResetCell.cpp                     |  5 +-
 src/generic/Time.cpp                          |  5 +-
 src/generic/UpdateIf.cpp                      |  5 +-
 src/generic/WholeMolecules.cpp                | 15 ++---
 src/generic/WrapAround.cpp                    | 20 +++---
 src/gridtools/ConvertToFES.cpp                |  4 +-
 src/gridtools/DumpCube.cpp                    |  4 +-
 src/gridtools/DumpGrid.cpp                    | 16 ++---
 src/gridtools/FindContour.cpp                 |  4 +-
 src/gridtools/FindSphericalContour.cpp        |  4 +-
 src/gridtools/FourierTransform.cpp            |  4 +-
 src/gridtools/InterpolateGrid.cpp             |  4 +-
 src/maketools/plumedcheck                     |  6 ++
 src/manyrestraints/LWalls.cpp                 |  4 +-
 src/manyrestraints/UWalls.cpp                 |  4 +-
 src/mapping/PCAVars.cpp                       |  4 +-
 src/mapping/Path.cpp                          | 16 ++---
 src/mapping/PropertyMap.cpp                   |  4 +-
 src/multicolvar/AlphaBeta.cpp                 | 12 ++--
 src/multicolvar/Angles.cpp                    | 12 ++--
 src/multicolvar/Bridge.cpp                    |  4 +-
 src/multicolvar/CenterOfMultiColvar.cpp       |  8 +--
 src/multicolvar/CoordinationNumbers.cpp       |  8 +--
 src/multicolvar/Density.cpp                   |  4 +-
 src/multicolvar/DihedralCorrelation.cpp       |  8 +--
 src/multicolvar/DistanceFromContour.cpp       |  4 +-
 src/multicolvar/Distances.cpp                 | 24 +++----
 src/multicolvar/DumpMultiColvar.cpp           |  4 +-
 src/multicolvar/FilterBetween.cpp             | 16 ++---
 src/multicolvar/FilterLessThan.cpp            | 12 ++--
 src/multicolvar/FilterMoreThan.cpp            | 16 ++---
 src/multicolvar/InPlaneDistances.cpp          |  4 +-
 src/multicolvar/LocalAverage.cpp              |  8 +--
 src/multicolvar/MultiColvarDensity.cpp        |  8 +--
 src/multicolvar/NumberOfLinks.cpp             |  8 +--
 src/multicolvar/Torsions.cpp                  |  8 +--
 src/multicolvar/VolumeAround.cpp              |  4 +-
 src/multicolvar/VolumeCavity.cpp              |  8 +--
 src/multicolvar/VolumeInCylinder.cpp          |  4 +-
 src/multicolvar/VolumeInSphere.cpp            |  4 +-
 src/multicolvar/VolumeTetrapore.cpp           |  8 +--
 src/multicolvar/XAngle.cpp                    | 12 ++--
 src/multicolvar/XDistances.cpp                | 48 +++++++-------
 src/multicolvar/XYDistances.cpp               | 12 ++--
 src/multicolvar/XYTorsion.cpp                 | 24 +++----
 src/secondarystructure/AlphaRMSD.cpp          |  5 +-
 src/secondarystructure/AntibetaRMSD.cpp       |  5 +-
 src/secondarystructure/ParabetaRMSD.cpp       |  5 +-
 src/setup/Load.cpp                            | 12 ++--
 src/setup/MolInfo.cpp                         | 10 ++-
 src/setup/Restart.cpp                         | 15 ++---
 src/setup/Units.cpp                           |  9 ++-
 src/vatom/COM.cpp                             |  5 +-
 src/vatom/Center.cpp                          |  5 +-
 src/vatom/FixedAtom.cpp                       | 10 ++-
 src/vatom/Ghost.cpp                           |  5 +-
 user-doc/Analysis.txt                         |  8 +--
 user-doc/Colvar.txt                           |  8 +--
 user-doc/Files.txt                            | 10 ++-
 user-doc/Functions.txt                        |  4 +-
 user-doc/Group.txt                            | 16 ++---
 user-doc/Misc.txt                             | 65 +++++++++----------
 user-doc/Performances.txt                     | 12 ++--
 user-doc/Regex.txt                            | 16 ++---
 user-doc/Syntax.txt                           | 14 ++--
 162 files changed, 675 insertions(+), 749 deletions(-)

diff --git a/src/adjmat/ClusterDiameter.cpp b/src/adjmat/ClusterDiameter.cpp
index ec26f2eeb..8bebeab92 100644
--- a/src/adjmat/ClusterDiameter.cpp
+++ b/src/adjmat/ClusterDiameter.cpp
@@ -40,7 +40,7 @@ than 2.0 and if they are within 6.0 nm of each other.  Depth first search cluste
 between every pair of atoms that are within the largest of the clusters found is then calculated and the largest of these distances is output to a file named
 colvar.
 
-\verbatim
+\plumedfile
 # Calculate coordination numbers
 c1: COORDINATIONNUMBER SPECIES=1-512 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0}
 # Select coordination numbers that are more than 2.0
@@ -52,7 +52,7 @@ dfs: DFSCLUSTERING MATRIX=mat
 clust1: CLUSTER_PROPERTIES CLUSTERS=dfs CLUSTER=1
 dia: CLUSTER_DIAMETER CLUSTERS=dfs CLUSTER=1
 PRINT ARG=dia FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/ClusterDistribution.cpp b/src/adjmat/ClusterDistribution.cpp
index 93f12fb04..afd7c956e 100644
--- a/src/adjmat/ClusterDistribution.cpp
+++ b/src/adjmat/ClusterDistribution.cpp
@@ -44,7 +44,7 @@ algorithm on the corresponding graph. The number of componets in this graph that
 As discussed in \cite tribello-clustering this input was used to analyse the formation of a polycrystal of GeTe from amorphous
 GeTe.
 
-\verbatim
+\plumedfile
 q6: Q6 SPECIES=1-32768 SWITCH={GAUSSIAN D_0=5.29 R_0=0.01 D_MAX=5.3} LOWMEM
 lq6: LOCAL_Q6 SPECIES=q6 SWITCH={GAUSSIAN D_0=5.29 R_0=0.01 D_MAX=5.3} LOWMEM
 flq6: MFILTER_MORE DATA=lq6 SWITCH={GAUSSIAN D_0=0.19 R_0=0.01 D_MAX=0.2}
@@ -54,7 +54,7 @@ mat: CONTACT_MATRIX ATOMS=fcc SWITCH={GAUSSIAN D_0=3.59 R_0=0.01 D_MAX=3.6}
 dfs: DFSCLUSTERING MATRIX=mat
 nclust: CLUSTER_DISTRIBUTION CLUSTERS=dfs TRANSFORM={GAUSSIAN D_0=5.99 R_0=0.01 D_MAX=6.0} MORE_THAN={GAUSSIAN D_0=26.99 R_0=0.01 D_MAX=27}
 PRINT ARG=nclust.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/ClusterProperties.cpp b/src/adjmat/ClusterProperties.cpp
index 9a0d38e5b..b1d891fe0 100644
--- a/src/adjmat/ClusterProperties.cpp
+++ b/src/adjmat/ClusterProperties.cpp
@@ -42,13 +42,13 @@ a graph.  This dfs action then finds the largest connected component in this gra
 numbers for the atoms in this largest connected component are then computed and this quantity is output to a colvar
 file.  The way this input can be used is described in detail in \cite tribello-clustering.
 
-\verbatim
+\plumedfile
 lq: COORDINATIONNUMBER SPECIES=1-100 SWITCH={CUBIC D_0=0.45  D_MAX=0.55} LOWMEM
 cm: CONTACT_MATRIX ATOMS=lq  SWITCH={CUBIC D_0=0.45  D_MAX=0.55}
 dfs: DFSCLUSTERING MATRIX=cm
 clust1: CLUSTER_PROPERTIES CLUSTERS=dfs CLUSTER=1 SUM
 PRINT ARG=clust1.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/ClusterSize.cpp b/src/adjmat/ClusterSize.cpp
index b0b22621d..4ae18bb49 100644
--- a/src/adjmat/ClusterSize.cpp
+++ b/src/adjmat/ClusterSize.cpp
@@ -39,7 +39,7 @@ The following input uses PLUMED to calculate a adjacency matrix that connects a
 than 2.0 and if they are within 6.0 nm of each other.  Depth first search clustering is used to find the connected components in this matrix and then
 the number of atoms in the largest cluster is found.  This quantity is then output to a file called colvar
 
-\verbatim
+\plumedfile
 # Calculate coordination numbers
 c1: COORDINATIONNUMBER SPECIES=1-512 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0}
 # Select coordination numbers that are more than 2.0
@@ -51,7 +51,7 @@ dfs: DFSCLUSTERING MATRIX=mat
 clust1: CLUSTER_PROPERTIES CLUSTERS=dfs CLUSTER=1
 nat: CLUSTER_NATOMS CLUSTERS=dfs CLUSTER=1
 PRINT ARG=nat FILE=COLVAR
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/ClusterWithSurface.cpp b/src/adjmat/ClusterWithSurface.cpp
index ab472eee5..1042accd0 100644
--- a/src/adjmat/ClusterWithSurface.cpp
+++ b/src/adjmat/ClusterWithSurface.cpp
@@ -44,7 +44,7 @@ number of atoms with indices that are between 1 and 1996 and that are either in
 atoms within the the second largest cluster are then counted and this number of atoms is output to a file called size.  In addition the indices of the atoms
 that were counted are output to a file called dfs2.dat.
 
-\verbatim
+\plumedfile
 c1: COORDINATIONNUMBER SPECIES=1-1996 SWITCH={CUBIC D_0=0.34 D_MAX=0.38}
 cf: MFILTER_LESS DATA=c1 SWITCH={CUBIC D_0=13 D_MAX=13.5}
 mat: CONTACT_MATRIX ATOMS=cf SWITCH={CUBIC D_0=0.34 D_MAX=0.38}
@@ -53,7 +53,7 @@ clust2a: CLUSTER_WITHSURFACE CLUSTERS=dfs RCUT_SURF=0.3
 size2a: CLUSTER_NATOMS CLUSTERS=clust2a CLUSTER=2
 PRINT ARG=size2a FILE=size FMT=%8.4f
 OUTPUT_CLUSTER CLUSTERS=clust2a CLUSTER=2 FILE=dfs2.dat
-\endverbatim
+\endplumedfile
 
 
 */
diff --git a/src/adjmat/ContactAlignedMatrix.cpp b/src/adjmat/ContactAlignedMatrix.cpp
index e1bdfcdba..d00542aea 100644
--- a/src/adjmat/ContactAlignedMatrix.cpp
+++ b/src/adjmat/ContactAlignedMatrix.cpp
@@ -57,12 +57,12 @@ is greater than 0.5.  The sum of the rows of this matrix are then computed.  The
 many of the molecules that are within the first coordination sphere of molecule \f$i\f$ have an orientation that is similar to that of
 molecule \f$i\f$.  We thus calculate the number of these "coordination numbers" that are greater than 1.0 and output this quantity to a file.
 
-\verbatim
+\plumedfile
 m1: MOLECULES MOL1=1,2 MOL2=3,4 MOL3=5,6 MOL4=7,8
 mat: ALIGNED_MATRIX ATOMS=m1 SWITCH={RATIONAL R_0=0.1} ORIENTATION_SWITCH={RATIONAL R_0=0.1 D_MAX=0.5}
 rr: ROWSUMS MATRIX=mat MORE_THAN={RATIONAL D_0=1.0 R_0=0.1}
 PRINT ARG=rr.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/ContactMatrix.cpp b/src/adjmat/ContactMatrix.cpp
index 1105a1aa2..e1a332a2a 100644
--- a/src/adjmat/ContactMatrix.cpp
+++ b/src/adjmat/ContactMatrix.cpp
@@ -48,11 +48,11 @@ The input shown below calculates a \f$6 \times 6\f$ matrix whose elements are eq
 of each other and which is zero otherwise.  The columns in this matrix are then summed so as to give the coordination number for each atom.
 The final quantity output in the colvar file is thus the average coordination number.
 
-\verbatim
+\plumedfile
 aa: CONTACT_MATRIX ATOMS=1-6 SWITCH={EXP D_0=0.2 R_0=0.1 D_MAX=0.66}
 COLUMNSUMS MATRIX=mat MEAN LABEL=csums
 PRINT ARG=csums.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/DFSClustering.cpp b/src/adjmat/DFSClustering.cpp
index af0dc327e..23bf5c592 100644
--- a/src/adjmat/DFSClustering.cpp
+++ b/src/adjmat/DFSClustering.cpp
@@ -53,13 +53,13 @@ a graph.  This dfs action then finds the largest connected component in this gra
 numbers for the atoms in this largest connected component are then computed and this quantity is output to a colvar
 file.  The way this input can be used is described in detail in \cite tribello-clustering.
 
-\verbatim
+\plumedfile
 lq: COORDINATIONNUMBER SPECIES=1-100 SWITCH={CUBIC D_0=0.45  D_MAX=0.55} LOWMEM
 cm: CONTACT_MATRIX ATOMS=lq  SWITCH={CUBIC D_0=0.45  D_MAX=0.55}
 dfs: DFSCLUSTERING MATRIX=cm
 clust1: CLUSTER_PROPERTIES CLUSTERS=dfs CLUSTER=1 SUM
 PRINT ARG=clust1.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/HbondMatrix.cpp b/src/adjmat/HbondMatrix.cpp
index f75e95754..32e7989de 100644
--- a/src/adjmat/HbondMatrix.cpp
+++ b/src/adjmat/HbondMatrix.cpp
@@ -66,13 +66,13 @@ on the number of hydrogen bonds each of the water molecules donates and accepts.
 five columns of data.  The first four of these columns are a label for the atom and the x, y and z position of the oxygen.  The last column is then
 the number of accepted/donated hydrogen bonds.
 
-\verbatim
+\plumedfile
 mat: HBOND_MATRIX ATOMS=1-192:3 HYDROGENS=2-192:3,3-192:3 SWITCH={RATIONAL R_0=3.20} HSWITCH={RATIONAL R_0=2.30} ASWITCH={RATIONAL R_0=0.167pi} SUM
 rsums: ROWSUMS MATRIX=mat MEAN
 csums: COLUMNSUMS MATRIX=mat MEAN
 DUMPMULTICOLVAR DATA=rsums FILE=donors.xyz
 DUMPMULTICOLVAR DATA=csums FILE=acceptors.x
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/MatrixColumnSums.cpp b/src/adjmat/MatrixColumnSums.cpp
index 773502349..cae24fcec 100644
--- a/src/adjmat/MatrixColumnSums.cpp
+++ b/src/adjmat/MatrixColumnSums.cpp
@@ -43,21 +43,21 @@ tells you whether atoms \f$i\f$ and \f$j\f$ are within 1.0 nm of each other.  Th
 and the average value is computed.  As such the following input provides an alternative method for calculating the coordination numbers
 of atoms 1 to 10.
 
-\verbatim
+\plumedfile
 mat: CONTACT_MATRIX ATOMS=1-10 SWITCH={RATIONAL R_0=1.0}
 rsums: COLUMNSUMS MATRIX=mat MEAN
 PRINT ARG=rsums.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following input demonstrates another way that an average coordination number can be computed.  This input calculates the number of atoms
 with indices between 1 and 5 that are within the first coordination spheres of each of the atoms within indices between 6 and 15.  The average
 coordination number is then calculated from these fifteen coordination numbers and this quantity is output to a file.
 
-\verbatim
+\plumedfile
 mat2: CONTACT_MATRIX ATOMSA=1-5 ATOMSB=6-15 SWITCH={RATIONAL R_0=1.0}
 rsums: COLUMNSUMS MATRIX=mat2 MEAN
 PRINT ARG=rsums.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/MatrixRowSums.cpp b/src/adjmat/MatrixRowSums.cpp
index 038aa89a7..6ef15fb28 100644
--- a/src/adjmat/MatrixRowSums.cpp
+++ b/src/adjmat/MatrixRowSums.cpp
@@ -43,21 +43,21 @@ tells you whether atoms \f$i\f$ and \f$j\f$ are within 1.0 nm of each other.  Th
 and the average value is computed.  As such the following input provides an alternative method for calculating the coordination numbers
 of atoms 1 to 10.
 
-\verbatim
+\plumedfile
 mat: CONTACT_MATRIX ATOMS=1-10 SWITCH={RATIONAL R_0=1.0}
 rsums: ROWSUMS MATRIX=mat MEAN
 PRINT ARG=rsums.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following input demonstrates another way that an average coordination number can be computed.  This input calculates the number of atoms
 with indices between 6 and 15 that are within the first coordination spheres of each of the atoms within indices between 1 and 5.  The average
 coordination number is then calculated from these five coordination numbers and this quantity is output to a file.
 
-\verbatim
+\plumedfile
 mat2: CONTACT_MATRIX ATOMSA=1-5 ATOMSB=6-15 SWITCH={RATIONAL R_0=1.0}
 rsums: ROWSUMS MATRIX=mat2 MEAN
 PRINT ARG=rsums.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/OutputCluster.cpp b/src/adjmat/OutputCluster.cpp
index c2b752f0d..4dd430f3c 100644
--- a/src/adjmat/OutputCluster.cpp
+++ b/src/adjmat/OutputCluster.cpp
@@ -47,13 +47,13 @@ that satisfy this criteria.  The DFS algorithm is then used to find the connecte
 in this matrix and the indices of the atoms in the largest connected component are then output
 to a file.
 
-\verbatim
+\plumedfile
 c1: COORDINATIONNUMBER SPECIES=1-1996 SWITCH={CUBIC D_0=0.34 D_MAX=0.38}
 cf: MFILTER_LESS DATA=c1 SWITCH={CUBIC D_0=13 D_MAX=13.5}
 mat: CONTACT_MATRIX ATOMS=cf SWITCH={CUBIC D_0=0.34 D_MAX=0.38}
 dfs: DFSCLUSTERING MATRIX=mat
 OUTPUT_CLUSTER CLUSTERS=dfs CLUSTER=1 FILE=dfs.dat
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/SMACMatrix.cpp b/src/adjmat/SMACMatrix.cpp
index ac433076b..a4e8288f4 100644
--- a/src/adjmat/SMACMatrix.cpp
+++ b/src/adjmat/SMACMatrix.cpp
@@ -49,7 +49,7 @@ molecules \f$i\f$ and \f$j\f$ are within 6 angstroms of each other and if the to
 of these molecules is close to 0 or \f$\pi\f$.  The various connected components of this matrix are determined using the
 \ref DFSCLUSTERING algorithm and then the size of the largest cluster of connectes molecules is output to a colvar file
 
-\verbatim
+\plumedfile
 UNITS LENGTH=A
 
 MOLECULES ...
@@ -70,7 +70,7 @@ SMAC_MATRIX ...
 dfs1: DFSCLUSTERING MATRIX=smacm
 cc2: CLUSTER_NATOMS CLUSTERS=dfs1 CLUSTER=1
 PRINT ARG=smac.*,cc1.*,cc2 FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/adjmat/Sprint.cpp b/src/adjmat/Sprint.cpp
index f7e0a0da8..3e90100e8 100644
--- a/src/adjmat/Sprint.cpp
+++ b/src/adjmat/Sprint.cpp
@@ -46,17 +46,17 @@ This example input calculates the 7 SPRINT coordinates for a 7 atom cluster of L
 atoms and prints their values to a file.  In this input the SPRINT coordinates are calculated
 in the manner described in ?? so two atoms are adjacent if they are within a cutoff:
 
-\verbatim
+\plumedfile
 DENSITY SPECIES=1-7 LABEL=d1
 CONTACT_MATRIX ATOMS=d1 SWITCH={RATIONAL R_0=0.1} LABEL=mat
 SPRINT MATRIX=mat LABEL=ss
 PRINT ARG=ss.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 This example input calculates the 14 SPRINT coordinates foa a molecule composed of 7 hydrogen and
 7 carbon atoms.  Once again two atoms are adjacent if they are within a cutoff:
 
-\verbatim
+\plumedfile
 DENSITY SPECIES=1-7 LABEL=c
 DENSITY SPECIES=8-14 LABEL=h
 
@@ -71,7 +71,7 @@ CONTACT_MATRIX ...
 SPRINT MATRIX=mat LABEL=ss
 
 PRINT ARG=ss.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/analysis/Average.cpp b/src/analysis/Average.cpp
index c96b2556b..9d0566da8 100644
--- a/src/analysis/Average.cpp
+++ b/src/analysis/Average.cpp
@@ -51,11 +51,11 @@ and output this to a file called COLVAR.  In this example it is assumed that no
 on the system and that the weights, \f$w(t')\f$ in the formulae above can thus all be set equal
 to one.
 
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,2
 d1a: AVERAGE ARG=d1
 PRINT ARG=d1a FILE=colvar STRIDE=100
-\endverbatim
+\endplumedfile
 
 The following example calculates the ensemble average for the torsional angle involving atoms 1, 2, 3 and 4.
 At variance with the previous example this quantity is periodic so the second formula in the above introduction
@@ -65,23 +65,23 @@ forgotten and the process of averaging is begun again.  The quantities output in
 block averages taken over the first 100 frames of the trajectory, the block average over the second 100 frames
 of trajectory and so on.
 
-\verbatim
+\plumedfile
 t1: TORSION ATOMS=1,2,3,4
 t1a: AVERAGE ARG=t1 CLEAR=100
 PRINT ARG=t1a FILE=colvar STRIDE=100
-\endverbatim
+\endplumedfile
 
 This third example incorporates a bias.  Notice that the effect the bias has on the ensemble average is removed by taking
 advantage of the \ref REWEIGHT_BIAS method.  The final ensemble averages output to the file are thus block ensemble averages for the
 unbiased canononical ensemble at a temperature of 300 K.
 
-\verbatim
+\plumedfile
 t1: TORSION ATOMS=1,2,3,4
 RESTRAINT ARG=t1 AT=pi KAPPA=100.
 ww: REWEIGHT_BIAS TEMP=300
 t1a: AVERAGE ARG=t1 LOGWEIGHTS=ww CLEAR=100
 PRINT ARG=t1a FILE=colvar STRIDE=100
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/analysis/ClassicalMultiDimensionalScaling.cpp b/src/analysis/ClassicalMultiDimensionalScaling.cpp
index 93bb25f1d..68ac6aeae 100644
--- a/src/analysis/ClassicalMultiDimensionalScaling.cpp
+++ b/src/analysis/ClassicalMultiDimensionalScaling.cpp
@@ -55,7 +55,7 @@ can be found in the tutorial \ref belfast-3 and in the following <a href="https:
 The following command instructs plumed to construct a classical multidimensional scaling projection of a trajectory.
 The RMSD distance between atoms 1-256 have moved is used to measure the distances in the high-dimensional space.
 
-\verbatim
+\plumedfile
 CLASSICAL_MDS ...
   ATOMS=1-256
   METRIC=OPTIMAL-FAST
@@ -63,7 +63,7 @@ CLASSICAL_MDS ...
   NLOW_DIM=2
   OUTPUT_FILE=rmsd-embed
 ... CLASSICAL_MDS
-\endverbatim
+\endplumedfile
 
 The following section is for people who are interested in how this method works in detail. A solid understanding of this material is
 not necessary to use MDS.
diff --git a/src/analysis/Commit.cpp b/src/analysis/Commit.cpp
index 5cbacf4db..d51d6b9e5 100644
--- a/src/analysis/Commit.cpp
+++ b/src/analysis/Commit.cpp
@@ -37,7 +37,7 @@ The following input monitors two torsional angles during a simulation,
 defines two basins (A and B) as a function of the two torsions and
 stops the simulation when it falls in one of the two. In the log
 file will be shown the latest values for the CVs and the basin reached.
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 COMMITTOR ...
@@ -48,7 +48,7 @@ COMMITTOR ...
   BASIN_LL2=-0.15,-0.20
   BASIN_UL2=-0.25,-0.40
 ... COMMITTOR
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/analysis/Histogram.cpp b/src/analysis/Histogram.cpp
index b8a74c797..16abb4a6d 100644
--- a/src/analysis/Histogram.cpp
+++ b/src/analysis/Histogram.cpp
@@ -84,7 +84,7 @@ Additional material and examples can be also found in the tutorial \ref belfast-
 
 The following input monitors two torsional angles during a simulation
 and outputs a continuos histogram as a function of them at the end of the simulation.
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 HISTOGRAM ...
@@ -97,11 +97,11 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hh FILE=histo
-\endverbatim
+\endplumedfile
 
 The following input monitors two torsional angles during a simulation
 and outputs a discrete histogram as a function of them at the end of the simulation.
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 HISTOGRAM ...
@@ -115,11 +115,11 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hh FILE=histo
-\endverbatim
+\endplumedfile
 
 The following input monitors two torsional angles during a simulation
 and outputs the histogram accumulated thus far every 100000 steps.
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 HISTOGRAM ...
@@ -132,14 +132,14 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hh FILE=histo STRIDE=100000
-\endverbatim
+\endplumedfile
 
 The following input monitors two torsional angles during a simulation
 and outputs a separate histogram for each 100000 steps worth of trajectory.
 Notice how the CLEAR keyword is used here and how it is not used in the
 previous example.
 
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 HISTOGRAM ...
@@ -153,7 +153,7 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hh FILE=histo STRIDE=100000
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/analysis/PCA.cpp b/src/analysis/PCA.cpp
index 5b31a04fa..528c6290b 100644
--- a/src/analysis/PCA.cpp
+++ b/src/analysis/PCA.cpp
@@ -63,9 +63,9 @@ of the first 22 atoms.  The TYPE=OPTIMAL instruction ensures that translational
 The first two principal components will be output to a file called pca-comp.pdb.  Trajectory frames will be collected on every step and the PCA calculation
 will be performed at the end of the simulation.
 
-\verbatim
+\plumedfile
 PCA METRIC=OPTIMAL ATOMS=1-22 STRIDE=1 USE_ALL_DATA NLOW_DIM=2 OFILE=pca-comp.pdb
-\endverbatim
+\endplumedfile
 
 The following input instructs PLUMED to perform a principal component analysis in which the covariance matrix is calculated from chnages in the six distances
 seen in the previous lines.  Notice that here the TYPE=EUCLIDEAN keyword is used to indicate that no alighment has to be done when calculating the various
@@ -75,7 +75,7 @@ PCA analysis will be performed twice.  The REWEIGHT_BIAS keyword in this input t
 when calculating averages and covariances a reweighting should be performed based and each frames' weight in these calculations should be determined based on
 the current value of the instantaneous bias (see \ref REWEIGHT_BIAS).
 
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,2
 d2: DISTANCE ATOMS=1,3
 d3: DISTANCE ATOMS=1,4
@@ -84,7 +84,7 @@ d5: DISTANCE ATOMS=2,4
 d6: DISTANCE ATOMS=3,4
 
 PCA ARG=d1,d2,d3,d4,d5,d6 METRIC=EUCLIDEAN STRIDE=5 RUN=1000 NLOW_DIM=2 REWEIGHT_BIAS OFILE=pca-comp.pdb
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/bias/ABMD.cpp b/src/bias/ABMD.cpp
index 5e8bef344..a76c5341c 100644
--- a/src/bias/ABMD.cpp
+++ b/src/bias/ABMD.cpp
@@ -70,13 +70,12 @@ an additional white noise acting on the minimum position of the bias.
 The following input sets up two biases, one on the distance between atoms 3 and 5
 and another on the distance between atoms 2 and 4. The two target values are defined
 using TO and the two strength using KAPPA. The total energy of the bias is printed.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 ABMD ARG=d1,d2 TO=1.0,1.5 KAPPA=5.0,5.0 LABEL=abmd
 PRINT ARG=abmd.bias,abmd.d1_min,abmd.d2_min
-\endverbatim
-(See also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/bias/BiasValue.cpp b/src/bias/BiasValue.cpp
index 3e77621a5..ff741a6c0 100644
--- a/src/bias/BiasValue.cpp
+++ b/src/bias/BiasValue.cpp
@@ -42,18 +42,17 @@ to some collective variable then using the value of this function directly as a
 The following input tells plumed to use the value of the distance between atoms 3 and 5
 and the value of the distance between atoms 2 and 4 as biases.
 It then tells plumed to print the energy of the restraint
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=3,6 LABEL=d2
 BIASVALUE ARG=d1,d2 LABEL=b
 PRINT ARG=d1,d2,b.d1,b.d2
-\endverbatim
-(See also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 
 Another thing one can do is asking one system to follow
 a circle in sin/cos according a  time dependence
 
-\verbatim
+\plumedfile
 t: TIME
 # this just print cos and sin of time
 cos: MATHEVAL ARG=t VAR=t FUNC=cos(t) PERIODIC=NO
@@ -72,9 +71,7 @@ vv1:  MATHEVAL ARG=mycos,mysin,cos,sin VAR=mc,ms,c,s  FUNC=100*((mc-c)^2+(ms-s)^
 cc: BIASVALUE ARG=vv1
 # some printout
 PRINT ARG=t,cos,sin,d.x,d.y,d.z,mycos,mysin,cc.bias.vv1 STRIDE=1 FILE=colvar FMT=%8.4f
-\endverbatim
-(see also \ref TIME, \ref MATHEVAL, \ref COM, \ref DISTANCE, and \ref PRINT).
-
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/bias/ExtendedLagrangian.cpp b/src/bias/ExtendedLagrangian.cpp
index 85de570dd..ee3d3ada0 100644
--- a/src/bias/ExtendedLagrangian.cpp
+++ b/src/bias/ExtendedLagrangian.cpp
@@ -88,28 +88,26 @@ and many shorter runs.
 
 The following input tells plumed to perform a metadynamics
 with an extended Lagrangian on two torsional angles.
-\verbatim
+\plumedfile
 phi: TORSION ATOMS=5,7,9,15
 psi: TORSION ATOMS=7,9,15,17
 ex: EXTENDED_LAGRANGIAN ARG=phi,psi KAPPA=20,20.0 TAU=0.1,0.1
 METAD ARG=ex.phi_fict,ex.psi_fict PACE=100 SIGMA=0.35,0.35 HEIGHT=0.1
 # monitor the two variables
 PRINT STRIDE=10 ARG=phi,psi,ex.phi_fict,ex.psi_fict FILE=COLVAR
-\endverbatim
-(See also \ref TORSION, \ref METAD, and \ref PRINT).
+\endplumedfile
 
 The following input tells plumed to perform a TAMD (or dAFED)
 calculation on two torsional angles, keeping the two variables
 at a fictitious temperature of 3000K with a Langevin thermostat
 with friction 10
-\verbatim
+\plumedfile
 phi: TORSION ATOMS=5,7,9,15
 psi: TORSION ATOMS=7,9,15,17
 ex: EXTENDED_LAGRANGIAN ARG=phi,psi KAPPA=20,20.0 TAU=0.1,0.1 FRICTION=10,10 TEMP=3000
 # monitor the two variables
 PRINT STRIDE=10 ARG=phi,psi,ex.phi_fict,ex.psi_fict FILE=COLVAR
-\endverbatim
-(See also \ref TORSION and \ref PRINT)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/bias/External.cpp b/src/bias/External.cpp
index 62bc54cf9..7c714efae 100644
--- a/src/bias/External.cpp
+++ b/src/bias/External.cpp
@@ -39,11 +39,10 @@ Calculate a restraint that is defined on a grid that is read during start up
 \par Examples
 The following is an input for a calculation with an external potential that is
 defined in the file bias.dat and that acts on the distance between atoms 3 and 5.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 EXTERNAL ARG=d1 FILE=bias.dat LABEL=external
-\endverbatim
-(See also \ref DISTANCE \ref PRINT).
+\endplumedfile
 
 The header in the file bias.dat should read:
 \verbatim
@@ -61,11 +60,11 @@ with NOSPLINE you do not need to provide derivative information.
 You can also include grids that are a function of more than one collective
 variable.  For instance the following would be the input for an external
 potential acting on two torsional angles:
-\verbatim
+\plumedfile
 TORSION ATOMS=4,5,6,7 LABEL=t1
 TORSION ATOMS=6,7,8,9 LABEL=t2
 EXTERNAL ARG=t1,t2 FILE=bias.dat LABEL=ext
-\endverbatim
+\endplumedfile
 
 The header in the file bias.dat for this calculation would read:
 \verbatim
diff --git a/src/bias/LWalls.cpp b/src/bias/LWalls.cpp
index 33d7f7c04..26bb720f1 100644
--- a/src/bias/LWalls.cpp
+++ b/src/bias/LWalls.cpp
@@ -52,13 +52,13 @@ The following input tells plumed to add both a lower and an upper walls on the d
 between atoms 3 and 5 and the distance between atoms 2 and 4. The lower and upper limits
 are defined at different values. The strength of the walls is the same for the four cases.
 It also tells plumed to print the energy of the walls.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 UPPER_WALLS ARG=d1,d2 AT=1.0,1.5 KAPPA=150.0,150.0 EXP=2,2 EPS=1,1 OFFSET=0,0 LABEL=uwall
 LOWER_WALLS ARG=d1,d2 AT=0.0,1.0 KAPPA=150.0,150.0 EXP=2,2 EPS=1,1 OFFSET=0,0 LABEL=lwall
 PRINT ARG=uwall.bias,lwall.bias
-\endverbatim
+\endplumedfile
 (See also \ref DISTANCE and \ref PRINT).
 
 */
diff --git a/src/bias/MetaD.cpp b/src/bias/MetaD.cpp
index b58618570..5023c38f2 100644
--- a/src/bias/MetaD.cpp
+++ b/src/bias/MetaD.cpp
@@ -163,35 +163,35 @@ The following input is for a standard metadynamics calculation using as
 collective variables the distance between atoms 3 and 5
 and the distance between atoms 2 and 4. The value of the CVs and
 the metadynamics bias potential are written to the COLVAR file every 100 steps.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 METAD ARG=d1,d2 SIGMA=0.2,0.2 HEIGHT=0.3 PACE=500 LABEL=restraint
 PRINT ARG=d1,d2,restraint.bias STRIDE=100  FILE=COLVAR
-\endverbatim
+\endplumedfile
 (See also \ref DISTANCE \ref PRINT).
 
 \par
 If you use adaptive Gaussians, with diffusion scheme where you use
 a Gaussian that should cover the space of 20 timesteps in collective variables.
 Note that in this case the histogram correction is needed when summing up hills.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 METAD ARG=d1,d2 SIGMA=20 HEIGHT=0.3 PACE=500 LABEL=restraint ADAPTIVE=DIFF
 PRINT ARG=d1,d2,restraint.bias STRIDE=100  FILE=COLVAR
-\endverbatim
+\endplumedfile
 
 \par
 If you use adaptive Gaussians, with geometrical scheme where you use
 a Gaussian that should cover the space of 0.05 nm in Cartesian space.
 Note that in this case the histogram correction is needed when summing up hills.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 METAD ARG=d1,d2 SIGMA=0.05 HEIGHT=0.3 PACE=500 LABEL=restraint ADAPTIVE=GEOM
 PRINT ARG=d1,d2,restraint.bias STRIDE=100  FILE=COLVAR
-\endverbatim
+\endplumedfile
 
 \par
 When using adaptive Gaussians you might want to limit how the hills width can change.
@@ -199,7 +199,7 @@ You can use SIGMA_MIN and SIGMA_MAX keywords.
 The sigmas should specified in terms of CV so you should use the CV units.
 Note that if you use a negative number, this means that the limit is not set.
 Note also that in this case the histogram correction is needed when summing up hills.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 METAD ...
@@ -207,7 +207,7 @@ METAD ...
   SIGMA_MIN=0.2,0.1 SIGMA_MAX=0.5,1.0
 ... METAD
 PRINT ARG=d1,d2,restraint.bias STRIDE=100  FILE=COLVAR
-\endverbatim
+\endplumedfile
 
 \par
 Multiple walkers can be also use as in  \cite multiplewalkers
@@ -215,7 +215,7 @@ These are enabled by setting the number of walker used, the id of the
 current walker which interprets the input file, the directory where the
 hills containing files resides, and the frequency to read the other walkers.
 Here is an example
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 METAD ...
    ARG=d1 SIGMA=0.05 HEIGHT=0.3 PACE=500 LABEL=restraint
@@ -224,7 +224,7 @@ METAD ...
    WALKERS_DIR=../
    WALKERS_RSTRIDE=100
 ... METAD
-\endverbatim
+\endplumedfile
 where  WALKERS_N is the total number of walkers, WALKERS_ID is the
 id of the present walker (starting from 0 ) and the WALKERS_DIR is the directory
 where all the walkers are located. WALKERS_RSTRIDE is the number of step between
@@ -237,7 +237,7 @@ presented in \cite Tiwary_jp504920s as described above.
 This is enabled by using the keyword REWEIGHTING_NGRID where the grid used for
 the calculation is set. The number of grid points given in REWEIGHTING_NGRID
 should be equal or larger than the number of grid points given in GRID_BIN.
-\verbatim
+\plumedfile
 METAD ...
  LABEL=metad
  ARG=phi,psi SIGMA=0.20,0.20 HEIGHT=1.20 BIASFACTOR=5 TEMP=300.0 PACE=500
@@ -245,7 +245,7 @@ METAD ...
  REWEIGHTING_NGRID=150,150
  REWEIGHTING_NHILLS=20
 ... METAD
-\endverbatim
+\endplumedfile
 Here we have asked that the calculation is performed every 20 hills by using
 REWEIGHTING_NHILLS keyword. If this keyword is not given the calculation will
 by default be performed every 50 hills. The c(t) reweighting factor will be given
@@ -285,7 +285,7 @@ Alternatively, if you use a BIASFACTOR yout simulation will converge to a free
 energy that is a linear combination of the target free energy and of the intrinsic free energy
 determined by the original force field.
 
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 METAD ...
  LABEL=t1
@@ -295,7 +295,7 @@ METAD ...
 ... METAD
 
 PRINT ARG=d1,t1.bias STRIDE=100 FILE=COLVAR
-\endverbatim
+\endplumedfile
 
 The header in the file dist.dat for this calculation would read:
 
@@ -309,19 +309,19 @@ The header in the file dist.dat for this calculation would read:
 
 Notice that BIASFACTOR can also be chosen as equal to 1. In this case one will perform
 unbiased sampling. Instead of using HEIGHT, one should provide the TAU parameter.
-\verbatim
+\plumedfile
 d: DISTANCE ATOMS=3,5
 METAD ARG=d SIGMA=0.1 TAU=4.0 TEMP=300 PACE=100 BIASFACTOR=1.0
-\endverbatim
+\endplumedfile
 The HILLS file obtained will still work with `plumed sum_hills` so as to plot a free-energy.
 The case where this makes sense is probably that of RECT simulations.
 
 Regarding RECT simulations, you can also use the RECT keyword so as to avoid using multiple input files.
 For instance, a single input file will be
-\verbatim
+\plumedfile
 d: DISTANCE ATOMS=3,5
 METAD ARG=d SIGMA=0.1 TAU=4.0 TEMP=300 PACE=100 RECT=1.0,1.5,2.0,3.0
-\endverbatim
+\endplumedfile
 The number of elements in the RECT array should be equal to the number of replicas.
 
 
diff --git a/src/bias/Metainference.cpp b/src/bias/Metainference.cpp
index e4d3e8b94..3953c436d 100644
--- a/src/bias/Metainference.cpp
+++ b/src/bias/Metainference.cpp
@@ -58,7 +58,7 @@ them and comparing them with a set of experimental values. RDCs are compared wit
 the experimental data but for a multiplication factor SCALE that is also sampled by
 MC on-the-fly
 
-\verbatim
+\plumedfile
 RDC ...
 LABEL=rdc
 SCALE=0.0001
@@ -83,12 +83,12 @@ LABEL=spe
 ... METAINFERENCE
 
 PRINT ARG=spe.bias FILE=BIAS STRIDE=1
-\endverbatim
+\endplumedfile
 
 in the following example instead of using one uncertainty parameter per data point we use
 a single uncertainty value in a long-tailed gaussian to take into account for outliers.
 
-\verbatim
+\plumedfile
 METAINFERENCE ...
 ARG=ardc.*
 NOISETYPE=OUTLIERS
@@ -99,9 +99,7 @@ SIGMA_MEAN=0.001
 TEMP=300
 LABEL=spe
 ... METAINFERENCE
-\endverbatim
-
-(See also \ref RDC and \ref ENSEMBLE).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/bias/MovingRestraint.cpp b/src/bias/MovingRestraint.cpp
index 13cafce97..d5317af1d 100644
--- a/src/bias/MovingRestraint.cpp
+++ b/src/bias/MovingRestraint.cpp
@@ -53,7 +53,7 @@ Additional material and examples can be also found in the tutorial \ref belfast-
 The following input is dragging the distance between atoms 2 and 4
 from 1 to 2 in the first 1000 steps, then back in the next 1000 steps.
 In the following 500 steps the restraint is progressively switched off.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=2,4 LABEL=d
 MOVINGRESTRAINT ...
   ARG=d
@@ -62,12 +62,12 @@ MOVINGRESTRAINT ...
   STEP2=2000 AT2=1.0
   STEP3=2500         KAPPA3=0.0
 ... MOVINGRESTRAINT
-\endverbatim
+\endplumedfile
 The following input is progressively building restraints
 distances between atoms 1 and 5 and between atoms 2 and 4
 in the first 1000 steps. Afterwards, the restraint is kept
 static.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=1,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 MOVINGRESTRAINT ...
@@ -75,10 +75,10 @@ MOVINGRESTRAINT ...
   STEP0=0    AT0=1.0,1.5 KAPPA0=0.0,0.0
   STEP1=1000 AT1=1.0,1.5 KAPPA1=1.0,1.0
 ... MOVINGRESTRAINT
-\endverbatim
+\endplumedfile
 The following input is progressively bringing atoms 1 and 2
 close to each other with an upper wall
-\verbatim
+\plumedfile
 DISTANCE ATOMS=1,2 LABEL=d1
 MOVINGRESTRAINT ...
   ARG=d1
@@ -86,7 +86,7 @@ MOVINGRESTRAINT ...
   STEP0=0    AT0=1.0 KAPPA0=10.0
   STEP1=1000 AT1=0.0
 ... MOVINGRESTRAINT
-\endverbatim
+\endplumedfile
 
 By default the Action is issuing some values which are
 the work on each degree of freedom, the center of the harmonic potential,
diff --git a/src/bias/PBMetaD.cpp b/src/bias/PBMetaD.cpp
index b50557f7d..76c35154b 100644
--- a/src/bias/PBMetaD.cpp
+++ b/src/bias/PBMetaD.cpp
@@ -150,18 +150,18 @@ The following input is for PBMetaD calculation using as
 collective variables the distance between atoms 3 and 5
 and the distance between atoms 2 and 4. The value of the CVs and
 the PBMetaD bias potential are written to the COLVAR file every 100 steps.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 PBMETAD ARG=d1,d2 SIGMA=0.2,0.2 HEIGHT=0.3 PACE=500 LABEL=pb FILE=HILLS_d1,HILLS_d2
 PRINT ARG=d1,d2,pb.bias STRIDE=100 FILE=COLVAR
-\endverbatim
+\endplumedfile
 (See also \ref DISTANCE and \ref PRINT).
 
 \par
 If you use well-tempered metadynamics, you should specify a single biasfactor and initial
 Gaussian height.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 PBMETAD ...
@@ -170,11 +170,11 @@ PACE=500 BIASFACTOR=8 LABEL=pb
 FILE=HILLS_d1,HILLS_d2
 ... PBMETAD
 PRINT ARG=d1,d2,pb.bias STRIDE=100 FILE=COLVAR
-\endverbatim
+\endplumedfile
 
 \par
 The following input enables the MPI version of multiple-walkers.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 PBMETAD ...
@@ -184,7 +184,7 @@ FILE=HILLS_d1,HILLS_d2
 WALKERS_MPI
 ... PBMETAD
 PRINT ARG=d1,d2,pb.bias STRIDE=100 FILE=COLVAR
-\endverbatim
+\endplumedfile
 
 \par
 The disk version of multiple-walkers can be
@@ -192,7 +192,7 @@ enabled by setting the number of walker used, the id of the
 current walker which interprets the input file, the directory where the
 hills containing files resides, and the frequency to read the other walkers.
 Here is an example
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 PBMETAD ...
@@ -205,7 +205,7 @@ WALKERS_DIR=../
 WALKERS_RSTRIDE=100
 ... PBMETAD
 PRINT ARG=d1,d2,pb.bias STRIDE=100 FILE=COLVAR
-\endverbatim
+\endplumedfile
 where  WALKERS_N is the total number of walkers, WALKERS_ID is the
 id of the present walker (starting from 0 ) and the WALKERS_DIR is the directory
 where all the walkers are located. WALKERS_RSTRIDE is the number of step between
diff --git a/src/bias/Restraint.cpp b/src/bias/Restraint.cpp
index 8bb3df202..85208d35b 100644
--- a/src/bias/Restraint.cpp
+++ b/src/bias/Restraint.cpp
@@ -49,13 +49,12 @@ Additional material and examples can be also found in the tutorial \ref belfast-
 The following input tells plumed to restrain the distance between atoms 3 and 5
 and the distance between atoms 2 and 4, at different equilibrium
 values, and to print the energy of the restraint
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 RESTRAINT ARG=d1,d2 AT=1.0,1.5 KAPPA=150.0,150.0 LABEL=restraint
 PRINT ARG=restraint.bias
-\endverbatim
-(See also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/bias/ReweightBias.cpp b/src/bias/ReweightBias.cpp
index d0b2d90a3..60976e644 100644
--- a/src/bias/ReweightBias.cpp
+++ b/src/bias/ReweightBias.cpp
@@ -43,7 +43,7 @@ restraint will have an effect on the region of phase space that will be sampled
 run using this variable.  Consequently, when the histogram as a function of the distance, \f$x\f$, is accumulated,
 we use reweighting into order to discount the effect of the bias from our final histogram.
 
-\verbatim
+\plumedfile
 x: DISTANCE ATOMS=1,2
 RESTRAINT ARG=x SLOPE=1.0 AT=0.0
 as: REWEIGHT_BIAS TEMP=300
@@ -59,7 +59,7 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hB FILE=histoB STRIDE=1 FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 
 */
diff --git a/src/bias/ReweightMetad.cpp b/src/bias/ReweightMetad.cpp
index 766ed40df..42a4855c2 100644
--- a/src/bias/ReweightMetad.cpp
+++ b/src/bias/ReweightMetad.cpp
@@ -39,7 +39,7 @@ we use reweighting into order to discount the effect of the bias from our final
 metadynamics instead.  Notice also that we have to specify how often we would like to calculate the c(t) reweighting
 factor and the grid over which we calculate c(t) in the input to the METAD command.
 
-\verbatim
+\plumedfile
 a: ANGLE ATOMS=1,2,3
 x: DISTANCE ATOMS=1,2
 METAD ARG=x PACE=100 SIGMA=0.1 HEIGHT=1.5 BIASFACTOR=5 GRID_MIN=0 GRID_MAX=10 GRID_BIN=100 REWEIGHTING_NGRID=100 REWEIGHTING_NHILLS=50
@@ -57,8 +57,7 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hB FILE=histoB STRIDE=1 FMT=%8.4f
-
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/bias/UWalls.cpp b/src/bias/UWalls.cpp
index 0f096f11b..bfe4d8929 100644
--- a/src/bias/UWalls.cpp
+++ b/src/bias/UWalls.cpp
@@ -52,14 +52,13 @@ The following input tells plumed to add both a lower and an upper walls on the d
 between atoms 3 and 5 and the distance between atoms 2 and 4. The lower and upper limits
 are defined at different values. The strength of the walls is the same for the four cases.
 It also tells plumed to print the energy of the walls.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=3,5 LABEL=d1
 DISTANCE ATOMS=2,4 LABEL=d2
 UPPER_WALLS ARG=d1,d2 AT=1.0,1.5 KAPPA=150.0,150.0 EXP=2,2 EPS=1,1 OFFSET=0,0 LABEL=uwall
 LOWER_WALLS ARG=d1,d2 AT=0.0,1.0 KAPPA=150.0,150.0 EXP=2,2 EPS=1,1 OFFSET=0,0 LABEL=lwall
 PRINT ARG=uwall.bias,lwall.bias
-\endverbatim
-(See also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/Angle.cpp b/src/colvar/Angle.cpp
index ed57d03ea..285d6fa50 100644
--- a/src/colvar/Angle.cpp
+++ b/src/colvar/Angle.cpp
@@ -67,7 +67,7 @@ This command tells plumed to calculate the angle between the vector connecting a
 the vector connecting atom 2 to atom 3 and to print it on file COLVAR1. At the same time,
 the angle between vector connecting atom 1 to atom 2 and the vector connecting atom 3 to atom 4 is printed
 on file COLVAR2.
-\verbatim
+\plumedfile
 
 a: ANGLE ATOMS=1,2,3
 # equivalently one could state:
@@ -77,8 +77,7 @@ b: ANGLE ATOMS=1,2,3,4
 
 PRINT ARG=a FILE=COLVAR1
 PRINT ARG=b FILE=COLVAR2
-\endverbatim
-(see also \ref PRINT)
+\endplumedfile
 
 
 */
diff --git a/src/colvar/CS2Backbone.cpp b/src/colvar/CS2Backbone.cpp
index 040e0dfad..6c606e4e5 100644
--- a/src/colvar/CS2Backbone.cpp
+++ b/src/colvar/CS2Backbone.cpp
@@ -116,17 +116,17 @@ Additional material and examples can be also found in the tutorial \ref belfast-
 In this first example the chemical shifts are used to calculate a scoring function to be used
 in NMR driven Metadynamics \cite Granata:2013dk :
 
-\verbatim
+\plumedfile
 whole: GROUP ATOMS=2612-2514:-1,961-1:-1,2466-962:-1,2513-2467:-1
 WHOLEMOLECULES ENTITY0=whole
 cs: CS2BACKBONE ATOMS=1-2612 NRES=176 DATA=../data/ TEMPLATE=template.pdb CAMSHIFT NOPBC
 metad: METAD ARG=cs HEIGHT=0.5 SIGMA=0.1 PACE=200 BIASFACTOR=10
 PRINT ARG=cs,metad.bias FILE=COLVAR STRIDE=100
-\endverbatim
+\endplumedfile
 
 In this second example the chemical shifts are used as replica-averaged restrained as in \cite Camilloni:2012je \cite Camilloni:2013hs.
 
-\verbatim
+\plumedfile
 cs: CS2BACKBONE ATOMS=1-174 DATA=data/ NRES=13
 encs: ENSEMBLE ARG=(cs\.hn_.*),(cs\.nh_.*)
 stcs: STATS ARG=encs.* SQDEVSUM PARARG=(cs\.exphn_.*),(cs\.expnh_.*)
@@ -134,7 +134,7 @@ RESTRAINT ARG=stcs.sqdevsum AT=0 KAPPA=0 SLOPE=24
 
 PRINT ARG=(cs\.hn_.*),(cs\.nh_.*) FILE=RESTRAINT STRIDE=100
 
-\endverbatim
+\endplumedfile
 
 (See also \ref WHOLEMOLECULES, \ref STATS, \ref METAD, \ref RESTRAINT and \ref PRINT)
 
diff --git a/src/colvar/Cell.cpp b/src/colvar/Cell.cpp
index 5c59e9f49..082fa237d 100644
--- a/src/colvar/Cell.cpp
+++ b/src/colvar/Cell.cpp
@@ -36,14 +36,13 @@ Calculate the components of the simulation cell
 
 \par Examples
 The following input tells plumed to print the squared modulo of each of the three lattice vectors
-\verbatim
+\plumedfile
 cell: CELL
 aaa:    COMBINE ARG=cell.ax,cell.ay,cell.az POWERS=2,2,2 PERIODIC=NO
 bbb:    COMBINE ARG=cell.bx,cell.by,cell.bz POWERS=2,2,2 PERIODIC=NO
 ccc:    COMBINE ARG=cell.cx,cell.cy,cell.cz POWERS=2,2,2 PERIODIC=NO
 PRINT ARG=aaa,bbb,ccc
-\endverbatim
-(See also \ref COMBINE and \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/Constant.cpp b/src/colvar/Constant.cpp
index 49423413d..247f9ac6d 100644
--- a/src/colvar/Constant.cpp
+++ b/src/colvar/Constant.cpp
@@ -44,21 +44,20 @@ The following input instructs plumed to compute the distance
 between atoms 1 and 2. If this distance is between 1.0 and 2.0, it is
 printed. If it is lower than 1.0 (larger than 2.0), 1.0 (2.0) is printed
 
-\verbatim
+\plumedfile
 cn: CONSTANT VALUES=1.0,2.0
 dis: DISTANCE ATOMS=1,2
 sss: SORT ARG=cn.v_0,dis,cn.v_1
 PRINT ARG=sss.2
-\endverbatim
-(See also \ref DISTANCE, \ref SORT, and \ref PRINT).
+\endplumedfile
 
 In case you want to pass a single value you can use VALUE:
-\verbatim
+\plumedfile
 cn: CONSTANT VALUE=1.0
 dis: DISTANCE ATOMS=1
 sss: SORT ARG=cn,dis
 PRINT ARG=sss.1
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/ContactMap.cpp b/src/colvar/ContactMap.cpp
index 9405f172f..9fe8f82c9 100644
--- a/src/colvar/ContactMap.cpp
+++ b/src/colvar/ContactMap.cpp
@@ -45,17 +45,17 @@ The following example calculates switching functions based on the distances betw
 1 and 2, 3 and 4 and 4 and 5. The values of these three switching functions are then output
 to a file named colvar.
 
-\verbatim
+\plumedfile
 CONTACTMAP ATOMS1=1,2 ATOMS2=3,4 ATOMS3=4,5 ATOMS4=5,6 SWITCH={RATIONAL R_0=1.5} LABEL=f1
 PRINT ARG=f1.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following example calculates the difference of the current contact map with respect
 to a reference provided. In this case REFERENCE is the fraction of contact that is formed
 (i.e. the distance between two atoms transformed with the SWITH), while R_0 is the contact
 distance. WEIGHT gives the relative weight of each contact to the final distance measure.
 
-\verbatim
+\plumedfile
 CONTACTMAP ...
 ATOMS1=1,2 REFERENCE1=0.1 WEIGHT1=0.5
 ATOMS2=3,4 REFERENCE2=0.5 WEIGHT2=1.0
@@ -67,7 +67,7 @@ CMDIST
 ... CONTACTMAP
 
 PRINT ARG=cmap FILE=colvar
-\endverbatim
+\endplumedfile
 
 The next example calculates calculates fraction of native contacts (Q)
 for Trp-cage mini-protein. R_0 is the distance at which the switch function is guaranteed to
@@ -80,7 +80,7 @@ WEIGHT is the 1/(number of contacts) giving equal weight to each contact.
 
 When using native contact Q switch function, please cite \cite best2013
 
-\verbatim
+\plumedfile
 # Full example available in regtest/basic/rt72/
 
 CONTACTMAP ...
@@ -97,9 +97,8 @@ SUM
 ... CONTACTMAP
 
 PRINT ARG=cmap FILE=colvar
-\endverbatim
-
-(See also \ref PRINT and \ref switchingfunction)
+\endplumedfile
+(See also \ref switchingfunction)
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/Coordination.cpp b/src/colvar/Coordination.cpp
index c1aa2da4f..fe0a8b08d 100644
--- a/src/colvar/Coordination.cpp
+++ b/src/colvar/Coordination.cpp
@@ -62,19 +62,19 @@ so that they actually count as "zero".
 \par Examples
 
 The following example instructs plumed to calculate the total coordination number of the atoms in group 1-10 with the atoms in group 20-100.  For atoms 1-10 coordination numbers are calculated that count the number of atoms from the second group that are within 0.3 nm of the central atom.  A neighbour list is used to make this calculation faster, this neighbour list is updated every 100 steps.
-\verbatim
+\plumedfile
 COORDINATION GROUPA=1-10 GROUPB=20-100 R_0=0.3 NLIST NL_CUTOFF=0.5 NL_STRIDE=100
-\endverbatim
+\endplumedfile
 
 The following is a dummy example which should compute the value 0 because the self interaction
 of atom 1 is skipped. Notice that in plumed 2.0 "self interactions" were not skipped, and the
 same calculation should return 1.
-\verbatim
+\plumedfile
 c: COORDINATION GROUPA=1 GROUPB=1 R_0=0.3
 PRINT ARG=c STRIDE=10
-\endverbatim
+\endplumedfile
 
-\verbatim
+\plumedfile
 c1: COORDINATION GROUPA=1-10 GROUPB=1-10 R_0=0.3
 x: COORDINATION GROUPA=1-10 R_0=0.3
 c2: COMBINE ARG=x COEFFICIENTS=2
@@ -82,8 +82,7 @@ c2: COMBINE ARG=x COEFFICIENTS=2
 # since it runs on half of the pairs. Notice that to get the same result you
 # should double it
 PRINT ARG=c1,c2 STRIDE=10
-\endverbatim
-See also \ref PRINT and \ref COMBINE
+\endplumedfile
 
 
 
diff --git a/src/colvar/DHEnergy.cpp b/src/colvar/DHEnergy.cpp
index 23ed21d58..d13aeecdc 100644
--- a/src/colvar/DHEnergy.cpp
+++ b/src/colvar/DHEnergy.cpp
@@ -56,12 +56,11 @@ Notice that if there are common atoms between GROUPA and GROUPB their interactio
 
 
 \par Examples
-\verbatim
+\plumedfile
 # this is printing the electrostatic interaction between two groups of atoms
 dh: DHENERGY GROUPA=1-10 GROUPB=11-20 EPSILON=80.0 I=0.1 TEMP=300.0
 PRINT ARG=dh
-\endverbatim
-(see also \ref PRINT)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/DRMSD.cpp b/src/colvar/DRMSD.cpp
index 5a3152ad4..49efb882a 100644
--- a/src/colvar/DRMSD.cpp
+++ b/src/colvar/DRMSD.cpp
@@ -65,15 +65,15 @@ the positions of the atoms in the reference file and their instantaneous
 position. Only pairs of atoms whose distance in the reference structure is within
 0.1 and 0.8 nm are considered.
 
-\verbatim
+\plumedfile
 DRMSD REFERENCE=file.pdb LOWER_CUTOFF=0.1 UPPER_CUTOFF=0.8
-\endverbatim
+\endplumedfile
 
 The following tells plumed to calculate a DRMSD value for a pair of molecules.
 
-\verbatim
+\plumedfile
 DRMSD REFERENCE=file.pdb LOWER_CUTOFF=0.1 UPPER_CUTOFF=0.8 TYPE=INTER-DRMSD
-\endverbatim
+\endplumedfile
 
 In the input reference file (file.pdb) the atoms in each of the two molecules are separated by a TER
 command as shown below.
diff --git a/src/colvar/Dimer.cpp b/src/colvar/Dimer.cpp
index 827c84365..baead008e 100644
--- a/src/colvar/Dimer.cpp
+++ b/src/colvar/Dimer.cpp
@@ -85,43 +85,43 @@ the temperature of the system.
 
 This line tells Plumed to compute the Dimer interaction energy for every dimer in the system.
 
-\verbatim
+\plumedfile
 dim: DIMER TEMP=300 Q=0.5 ALLATOMS DSIGMA=0.002
-\endverbatim
+\endplumedfile
 
 If the simulation doesn't use virtual sites for the dimers centers of mass,
 Plumed has to know in order to determine correctly the total number of dimers from
 the total number of atoms:
-\verbatim
+\plumedfile
 dim: DIMER TEMP=300 Q=0.5 ALLATOMS DSIGMA=0.002 NOVSITES
-\endverbatim
+\endplumedfile
 
 The NOVSITES flag is not required if one provides the atom serial of each Dimer. This is
 defined as the atom serial of the first bead of the dimer and is thus a number between 1 and N.
 Along with the ATOMS list also the number N of lines describing the first beads has to be given.
 For example, the Dimer interaction energy of dimers 1,5,7 is:
-\verbatim
+\plumedfile
 dim: DIMER TEMP=300 Q=0.5 NATOMS=N ATOMS=1,5,7 DSIGMA=0.002
-\endverbatim
+\endplumedfile
 
 In a Replica Exchange simulation the keyword DSIGMA can be used in two ways:
 if a plumed.n.dat file is provided for each replica, then DSIGMA is passed as a single value,
 like in the previous examples, and each replica will read its own DSIGMA value. If
 a unique plumed.dat is given, DSIGMA has to be a list containing a value for each replica.
 For 4 replicas:
-\verbatim
+\plumedfile
 dim: DIMER TEMP=300 Q=0.5 NATOMS=N ATOMS=1,5,7 DSIGMA=0.002,0.002,0.004,0.01
-\endverbatim
+\endplumedfile
 
 
 \par Usage of the CV
 
 The dimer interaction is not coded in the driver program and has to be inserted
 in the hamiltonian of the system as a linear RESTRAINT (see \ref RESTRAINT):
-\verbatim
+\plumedfile
 dim: DIMER TEMP=300 Q=0.5 ALLATOMS DSIGMA=0.002
 RESTRAINT ARG=dim AT=0 KAPPA=0 SLOPE=1 LABEL=dimforces
-\endverbatim
+\endplumedfile
 
 In a replica exchange, Metadynamics (see \ref METAD) can be used on the Dimer CV to reduce
 the number of replicas. Just keep in mind that METAD SIGMA values should be tuned
diff --git a/src/colvar/Dipole.cpp b/src/colvar/Dipole.cpp
index 09a4a6580..3cffa3553 100644
--- a/src/colvar/Dipole.cpp
+++ b/src/colvar/Dipole.cpp
@@ -37,11 +37,10 @@ Calculate the dipole moment for a group of atoms.
 \par Examples
 The following tells plumed to calculate the dipole of the group of atoms containing
 the atoms from 1-10 and print it every 5 steps
-\verbatim
+\plumedfile
 d: DIPOLE GROUP=1-10
 PRINT FILE=output STRIDE=5 ARG=5
-\endverbatim
-(see also \ref PRINT)
+\endplumedfile
 
 \attention
 If the total charge Q of the group in non zero, then a charge Q/N will be subtracted to every atom,
diff --git a/src/colvar/Distance.cpp b/src/colvar/Distance.cpp
index c3a981771..279961b76 100644
--- a/src/colvar/Distance.cpp
+++ b/src/colvar/Distance.cpp
@@ -49,46 +49,44 @@ better to use SCALED_COMPONENTS.
 
 The following input tells plumed to print the distance between atoms 3 and 5,
 the distance between atoms 2 and 4 and the x component of the distance between atoms 2 and 4.
-\verbatim
+\plumedfile
 d1:  DISTANCE ATOMS=3,5
 d2:  DISTANCE ATOMS=2,4
 d2c: DISTANCE ATOMS=2,4 COMPONENTS
 PRINT ARG=d1,d2,d2c.x
-\endverbatim
-(See also \ref PRINT).
+\endplumedfile
 
 The following input computes the end-to-end distance for a polymer
 of 100 atoms and keeps it at a value around 5.
-\verbatim
+\plumedfile
 WHOLEMOLECULES ENTITY0=1-100
 e2e: DISTANCE ATOMS=1,100 NOPBC
 RESTRAINT ARG=e2e KAPPA=1 AT=5
-\endverbatim
-(See also \ref WHOLEMOLECULES and \ref RESTRAINT).
+\endplumedfile
 
 Notice that NOPBC is used
 to be sure that if the end-to-end distance is larger than half the simulation
 box the distance is compute properly. Also notice that, since many MD
 codes break molecules across cell boundary, it might be necessary to
 use the \ref WHOLEMOLECULES keyword (also notice that it should be
-_before_ distance). The list of atoms provided to WHOLEMOLECULES
+_before_ distance). The list of atoms provided to \ref WHOLEMOLECULES
 here contains all the atoms between 1 and 100. Strictly speaking, this
 is not necessary. If you know for sure that atoms with difference in
 the index say equal to 10 are _not_ going to be farther than half cell
 you can e.g. use
-\verbatim
+\plumedfile
 WHOLEMOLECULES ENTITY0=1,10,20,30,40,50,60,70,80,90,100
 e2e: DISTANCE ATOMS=1,100 NOPBC
 RESTRAINT ARG=e2e KAPPA=1 AT=5
-\endverbatim
-Just be sure that the ordered list provide to WHOLEMOLECULES has the following
+\endplumedfile
+Just be sure that the ordered list provide to \ref WHOLEMOLECULES has the following
 properties:
 - Consecutive atoms should be closer than half-cell throughout the entire simulation.
 - Atoms required later for the distance (e.g. 1 and 100) should be included in the list
 
 The following example shows how to take into account periodicity e.g.
 in z-component of a distance
-\verbatim
+\plumedfile
 # this is a center of mass of a large group
 c: COM ATOMS=1-100
 # this is the distance between atom 101 and the group
@@ -99,8 +97,7 @@ d: DISTANCE ATOMS=c,101 COMPONENTS
 dz: COMBINE ARG=d.z PERIODIC=-10,10
 # metadynamics on dd
 METAD ARG=dz SIGMA=0.1 HEIGHT=0.1 PACE=200
-\endverbatim
-(see also \ref COM, \ref COMBINE, and \ref METAD)
+\endplumedfile
 
 Using SCALED_COMPONENTS this problem should not arise because they are always periodic
 with domain (-0.5,+0.5).
diff --git a/src/colvar/ERMSD.cpp b/src/colvar/ERMSD.cpp
index 3602df086..fb5ff6cad 100644
--- a/src/colvar/ERMSD.cpp
+++ b/src/colvar/ERMSD.cpp
@@ -84,9 +84,9 @@ for a multi molecular system.
 Calculate the eRMSD from reference structure reference.pdb using the default cutoff (2.4). The list of residues involved in the calculation has to be specified. In this example, the eRMSD is calculated
 considering residues 1,2,3,4,5,6.
 
-\verbatim
+\plumedfile
 eRMSD1: ERMSD REFERENCE=reference.pdb ATOMS=@lcs-1,@lcs-2,@lcs-3,@lcs-4,@lcs-5,@lcs-6
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/Energy.cpp b/src/colvar/Energy.cpp
index 6f79eb883..4b3f80736 100644
--- a/src/colvar/Energy.cpp
+++ b/src/colvar/Energy.cpp
@@ -48,11 +48,10 @@ using GROMACS with lambda replica exchange of with plumed-hrex branch.
 
 \par Examples
 The following input instructs plumed to print the energy of the system
-\verbatim
+\plumedfile
 ene: ENERGY
 PRINT ARG=ene
-\endverbatim
-(See also \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/Fake.cpp b/src/colvar/Fake.cpp
index 70bd7f278..20d70aed9 100644
--- a/src/colvar/Fake.cpp
+++ b/src/colvar/Fake.cpp
@@ -39,10 +39,9 @@ and just support input and period definition
 
 \par Examples
 
-\verbatim
+\plumedfile
 FAKE ATOMS=1 PERIODIC=-3.14,3.14   LABEL=d2
-\endverbatim
-(See also \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/FretEfficiency.cpp b/src/colvar/FretEfficiency.cpp
index 02600e105..e6f149452 100644
--- a/src/colvar/FretEfficiency.cpp
+++ b/src/colvar/FretEfficiency.cpp
@@ -51,22 +51,20 @@ boundary conditions. This behavior can be changed with the NOPBC flag.
 The following input tells plumed to print the FRET efficiencies
 calculated as a function of the distance between atoms 3 and 5 and
 the distance between atoms 2 and 4.
-\verbatim
+\plumedfile
 fe1:  FRET ATOMS=3,5 R0=5.5
 fe2:  FRET ATOMS=2,4 R0=5.5
 PRINT ARG=fe1,fe2
-\endverbatim
-(See also \ref PRINT).
+\endplumedfile
 
 The following input computes the FRET efficiency calculated on the
 terminal atoms of a polymer
 of 100 atoms and keeps it at a value around 0.5.
-\verbatim
+\plumedfile
 WHOLEMOLECULES ENTITY0=1-100
 fe: FRET ATOMS=1,100 R0=5.5 NOPBC
 RESTRAINT ARG=fe KAPPA=100 AT=0.5
-\endverbatim
-(See also \ref WHOLEMOLECULES and \ref RESTRAINT).
+\endplumedfile
 
 Notice that NOPBC is used
 to be sure that if the distance is larger than half the simulation
diff --git a/src/colvar/Gyration.cpp b/src/colvar/Gyration.cpp
index 791f6a3ed..3764ad404 100644
--- a/src/colvar/Gyration.cpp
+++ b/src/colvar/Gyration.cpp
@@ -75,11 +75,10 @@ periodic image.
 
 The following input tells plumed to print the radius of gyration of the
 chain containing atoms 10 to 20.
-\verbatim
+\plumedfile
 GYRATION TYPE=RADIUS ATOMS=10-20 LABEL=rg
 PRINT ARG=rg STRIDE=1 FILE=colvar
-\endverbatim
-(See also \ref PRINT)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/Implicit.cpp b/src/colvar/Implicit.cpp
index ad53f243c..b43b3ad60 100644
--- a/src/colvar/Implicit.cpp
+++ b/src/colvar/Implicit.cpp
@@ -58,7 +58,7 @@ where \f$\Delta G^\mathrm{free}_i\f$ is the solvation free energy of the isolate
 The output from this collective variable, the free energy of solvation, can be used with the \ref BIASVALUE keyword to provide implicit solvation to a system. All parameters are designed to be used with a modified CHARMM36 force field. It takes only non-hydrogen atoms as input, these can be conveniently specified using the \ref GROUP action with the NDX_GROUP parameter. To speed up the calculation, IMPLICIT internally uses a neighbourlist with a cutoff dependent on the type of atom (maximum of 1.95 nm). This cutoff can be extended further by using the NL_BUFFER keyword.
 
 \par Examples
-\verbatim
+\plumedfile
 MOLINFO MOLTYPE=protein STRUCTURE=peptide.pdb
 WHOLEMOLECULES ENTITY0=1-111
 
@@ -72,8 +72,7 @@ solv: IMPLICIT ATOMS=protein-h NL_STRIDE=10 NL_BUFFER=0.2
 bias: BIASVALUE ARG=solv
 
 PRINT ARG=solv FILE=SOLV
-\endverbatim
-(see also \ref PRINT, \ref GROUP, \ref MOLINFO, \ref WHOLEMOLECULES)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/Jcoupling.cpp b/src/colvar/Jcoupling.cpp
index d192b8bbd..8d8975f91 100644
--- a/src/colvar/Jcoupling.cpp
+++ b/src/colvar/Jcoupling.cpp
@@ -62,7 +62,7 @@ In the following example we calculate the Ha-N J-coupling from a set of atoms in
 dihedral \f$\psi\f$ angles in the peptide backbone. We also add the experimental datapoints and compute
 the correlation and other measures and finally print the results.
 
-\verbatim
+\plumedfile
 
 MOLINFO MOLTYPE=protein STRUCTURE=peptide.pdb
 WHOLEMOLECULES ENTITY0=1-111
@@ -81,11 +81,7 @@ JCOUPLING ...
 jhanst: STATS ARG=(jhan\.j_.*) PARARG=(jhan\.exp_.*)
 
 PRINT ARG=jhanst.*,jhan.* FILE=COLVAR STRIDE=100
-
-ENDPLUMED
-
-\endverbatim
-(See also \ref PRINT, \ref STATS)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/MultiRMSD.cpp b/src/colvar/MultiRMSD.cpp
index 7656c232b..4cf4ff08b 100644
--- a/src/colvar/MultiRMSD.cpp
+++ b/src/colvar/MultiRMSD.cpp
@@ -99,17 +99,17 @@ The following tells plumed to calculate the RMSD distance between
 the positions of the atoms in the reference file and their instantaneous
 position.  The Kearseley algorithm for each of the domains.
 
-\verbatim
+\plumedfile
 MULTI-RMSD REFERENCE=file.pdb TYPE=MULTI-OPTIMAL
-\endverbatim
+\endplumedfile
 
 The following tells plumed to calculate the RMSD distance btween the positions of
 the atoms in the domains of reference the reference structure and their instantaneous
 positions.  Here distances are calculated using the \ref DRMSD measure.
 
-\verbatim
+\plumedfile
 MULTI-RMSD REFERENCE=file.pdb TYPE=MULTI-DRMSD
-\endverbatim
+\endplumedfile
 
 in this case it is possible to use the following DRMSD options in the pdb file using the REMARK syntax:
 \verbatim
diff --git a/src/colvar/NOE.cpp b/src/colvar/NOE.cpp
index 82d2b4372..b46bd6597 100644
--- a/src/colvar/NOE.cpp
+++ b/src/colvar/NOE.cpp
@@ -52,7 +52,7 @@ In the following examples three noes are defined, the first is calculated based
 of atom 1-2 and 3-2; the second is defined by the distance 5-7 and the third by the distances
 4-15,4-16,8-15,8-16.
 
-\verbatim
+\plumedfile
 NOE ...
 GROUPA1=1,3 GROUPB1=2,2
 GROUPA2=5 GROUPB2=7
@@ -61,8 +61,7 @@ LABEL=noes
 ... NOE
 
 PRINT ARG=noes.* FILE=colvar
-\endverbatim
-(See also \ref PRINT)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/PCARMSD.cpp b/src/colvar/PCARMSD.cpp
index f5d1cf6f2..02861e109 100644
--- a/src/colvar/PCARMSD.cpp
+++ b/src/colvar/PCARMSD.cpp
@@ -57,9 +57,9 @@ Note that beta and occupancy values in the pdb are neglected and all the weights
 
 \par Examples
 
-\verbatim
+\plumedfile
 PCARMSD AVERAGE=file.pdb EIGENVECTORS=eigenvectors.pdb
-\endverbatim
+\endplumedfile
 
 The input is taken so to be compatible with the output you get from g_covar utility of gromacs (suitably adapted to have a pdb input format).
 
diff --git a/src/colvar/PRE.cpp b/src/colvar/PRE.cpp
index 906826b1f..82c94eee2 100644
--- a/src/colvar/PRE.cpp
+++ b/src/colvar/PRE.cpp
@@ -47,7 +47,7 @@ In the following example five PRE intensities are calculated using the distance
 oxigen of the spin label and the backbone hydrogens. Omega is the NMR frequency, RTWO the
 R2 for the hydrogens, INEPT of 8 ms for the experiment and a TAUC of 1.21 ns
 
-\verbatim
+\plumedfile
 PRE ...
 LABEL=HN_pre
 INEPT=8
@@ -63,7 +63,7 @@ GROUPA5=451 RTWO5=0.0086341843
 
 PRINT ARG=HN_pre.* FILE=PRE.dat STRIDE=1
 
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/PathMSD.cpp b/src/colvar/PathMSD.cpp
index cdfd96ee1..17567f2d4 100644
--- a/src/colvar/PathMSD.cpp
+++ b/src/colvar/PathMSD.cpp
@@ -42,10 +42,10 @@ in input ("sss" component) and the distance from them ("zzz" component).
 Here below is a case where you have defined three frames and you want to
 calculate the progress along the path and the distance from it in p1
 
-\verbatim
+\plumedfile
 p1: PATHMSD REFERENCE=file.pdb  LAMBDA=500.0 NEIGH_STRIDE=4 NEIGH_SIZE=8
 PRINT ARG=p1.sss,p1.zzz STRIDE=1 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 note that NEIGH_STRIDE=4 NEIGH_SIZE=8 control the neighborlist parameter (optional but
 recommended for performance) and states that the neighbor list will be calculated every 4
diff --git a/src/colvar/Position.cpp b/src/colvar/Position.cpp
index a381e039e..b66ae5e6c 100644
--- a/src/colvar/Position.cpp
+++ b/src/colvar/Position.cpp
@@ -54,14 +54,12 @@ This can be done e.g. using \ref FIT_TO_TEMPLATE.
 
 \par Examples
 
-\verbatim
+\plumedfile
 # align to a template
 FIT_TO_TEMPLATE REFERENCE=ref.pdb
 p: POSITION ATOM=3
 PRINT ARG=p.x,p.y,p.z
-\endverbatim
-(see also \ref FIT_TO_TEMPLATE)
-
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/PropertyMap.cpp b/src/colvar/PropertyMap.cpp
index a8a21fef2..eb1c27f02 100644
--- a/src/colvar/PropertyMap.cpp
+++ b/src/colvar/PropertyMap.cpp
@@ -46,10 +46,10 @@ where the parameters \f$X_i\f$  and  \f$Y_i\f$ are provided in the input pdb (al
  \f$D_i(x)\f$  is the MSD after optimal alignment calculated on the pdb frames you input (see Kearsley).
 
 \par Examples
-\verbatim
+\plumedfile
 p3: PROPERTYMAP REFERENCE=../../trajectories/path_msd/allv.pdb PROPERTY=X,Y LAMBDA=69087 NEIGH_SIZE=8 NEIGH_STRIDE=4
 PRINT ARG=p3.X,p3.Y,p3.zzz STRIDE=1 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 note that NEIGH_STRIDE=4 NEIGH_SIZE=8 control the neighborlist parameter (optional but
 recommended for performance) and states that the neighbor list will be calculated every 4
diff --git a/src/colvar/Puckering.cpp b/src/colvar/Puckering.cpp
index b90bdfbc4..66edb54c5 100644
--- a/src/colvar/Puckering.cpp
+++ b/src/colvar/Puckering.cpp
@@ -60,11 +60,11 @@ namespace colvar {
  \par Examples
 
  This input tells plumed to print the puckering phase angle of the 3rd nucleotide of a RNA molecule on file COLVAR.
- \verbatim
+ \plumedfile
  MOLINFO STRUCTURE=rna.pdb MOLTYPE=rna
  PUCKERING ATOMS=@sugar-3 LABEL=puck
  PRINT ARG=puck.phs FILE=COLVAR
- \endverbatim
+ \endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/RDC.cpp b/src/colvar/RDC.cpp
index 18db38014..43dc4c1ee 100644
--- a/src/colvar/RDC.cpp
+++ b/src/colvar/RDC.cpp
@@ -94,7 +94,7 @@ respect to a set of experimental data is calculated and restrained. In addition,
 and only for analysis purposes, the same RDCs are calculated using a Single Value
 Decomposition algorithm.
 
-\verbatim
+\plumedfile
 RDC ...
 GYROM=-72.5388
 SCALE=1.0
@@ -124,8 +124,7 @@ LABEL=svd
 
 PRINT ARG=nh.corr,rdce.bias FILE=colvar
 PRINT ARG=svd.* FILE=svd
-\endverbatim
-(See also \ref PRINT, \ref RESTRAINT)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/colvar/RMSD.cpp b/src/colvar/RMSD.cpp
index 8801dfeed..ffdcefeeb 100644
--- a/src/colvar/RMSD.cpp
+++ b/src/colvar/RMSD.cpp
@@ -133,9 +133,9 @@ The following tells plumed to calculate the RMSD distance between
 the positions of the atoms in the reference file and their instantaneous
 position.  The Kearseley algorithm is used so this is done optimally.
 
-\verbatim
+\plumedfile
 RMSD REFERENCE=file.pdb TYPE=OPTIMAL
-\endverbatim
+\endplumedfile
 
 ...
 
diff --git a/src/colvar/Template.cpp b/src/colvar/Template.cpp
index 85cb95d27..6fb3c5845 100644
--- a/src/colvar/Template.cpp
+++ b/src/colvar/Template.cpp
@@ -40,11 +40,11 @@ This file provides a template for if you want to introduce a new CV.
 
 <!---You should put an example of how to use your CV here--->
 
-\verbatim
+\plumedfile
 # This should be a sample input.
 t: TEMPLATE ATOMS=1,2
 PRINT ARG=t STRIDE=100 FILE=COLVAR
-\endverbatim
+\endplumedfile
 <!---You should reference here the other actions used in this example--->
 (see also \ref PRINT)
 
diff --git a/src/colvar/Torsion.cpp b/src/colvar/Torsion.cpp
index 274d0565a..4b697b069 100644
--- a/src/colvar/Torsion.cpp
+++ b/src/colvar/Torsion.cpp
@@ -43,23 +43,23 @@ orthogonal to an axis.
 
 This input tells plumed to print the torsional angle between atoms 1, 2, 3 and 4
 on file COLVAR.
-\verbatim
+\plumedfile
 t: TORSION ATOMS=1,2,3,4
 # this is an alternative, equivalent, definition:
 # t: TORSION VECTOR1=2,1 AXIS=2,3 VECTOR2=3,4
 PRINT ARG=t FILE=COLVAR
-\endverbatim
+\endplumedfile
 
 If you are working with a protein you can specify the special named torsion angles \f$\phi\f$, \f$\psi\f$, \f$\omega\f$ and \f$\chi_1\f$
 by using TORSION in combination with the \ref MOLINFO command.  This can be done by using the following
 syntax.
 
-\verbatim
+\plumedfile
 MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
 t1: TORSION ATOMS=@phi-3
 t2: TORSION ATOMS=@psi-4
 PRINT ARG=t1,t2 FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 Here, \@phi-3 tells plumed that you would like to calculate the \f$\phi\f$ angle in the third residue of the protein.
 Similarly \@psi-4 tells plumed that you want to calculate the \f$\psi\f$ angle of the 4th residue of the protein.
diff --git a/src/colvar/Volume.cpp b/src/colvar/Volume.cpp
index 993dc802a..837f6ce75 100644
--- a/src/colvar/Volume.cpp
+++ b/src/colvar/Volume.cpp
@@ -36,11 +36,10 @@ Calculate the volume of the simulation box.
 
 \par Examples
 The following input tells plumed to print the volume of the system
-\verbatim
+\plumedfile
 vol: VOLUME
 PRINT ARG=vol
-\endverbatim
-(See also \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/Fccubic.cpp b/src/crystallization/Fccubic.cpp
index fc1793ab3..47154f4c6 100644
--- a/src/crystallization/Fccubic.cpp
+++ b/src/crystallization/Fccubic.cpp
@@ -63,10 +63,10 @@ so on.  Notice also that you can rotate the reference frame if you are using a n
 The following input calculates the FCCUBIC parameter for the 64 atoms in the system
 and then calculates and prints the average value for this quantity.
 
-\verbatim
+\plumedfile
 FCCUBIC SPECIES=1-64 SWITCH={RATIONAL D_0=3.0 R_0=1.5} MEAN LABEL=d
 PRINT ARG=d.* FILE=colv
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/Gradient.cpp b/src/crystallization/Gradient.cpp
index 173d9e5ad..e3af25856 100644
--- a/src/crystallization/Gradient.cpp
+++ b/src/crystallization/Gradient.cpp
@@ -39,22 +39,22 @@ The input below calculates the gradient of the density of atoms in the manner
 described in \cite fede-grad in order to detect whether or not atoms are distributed
 uniformly along the x-axis of the simulation cell.
 
-\verbatim
+\plumedfile
 d1: DENSITY SPECIES=1-50
 s1: GRADIENT ORIGIN=1 DATA=d1 DIR=x NBINS=4 SIGMA=1.0
 PRINT ARG=s1 FILE=colvar
-\endverbatim
+\endplumedfile
 
 The input below calculates the coordination numbers of the 50 atoms in the simulation cell.
 The gradient of this quantity is then evaluated in the manner described using the equation above
 to detect whether the average values of the coordination number are uniformly distributed along the
 x-axis of the simulation cell.
 
-\verbatim
+\plumedfile
 d2: COORDINATIONNUMBER SPECIES=1-50 SWITCH={RATIONAL R_0=2.0} MORE_THAN={EXP R_0=4.0}
 s2: GRADIENT ORIGIN=1 DATA=d2 DIR=x NBINS=4 SIGMA=1.0
 PRINT ARG=s2 FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/InterMolecularTorsions.cpp b/src/crystallization/InterMolecularTorsions.cpp
index 3839e094d..58294492f 100644
--- a/src/crystallization/InterMolecularTorsions.cpp
+++ b/src/crystallization/InterMolecularTorsions.cpp
@@ -64,12 +64,12 @@ between the two molecules.  As such the torsional angles between molecules that
 histogram while the torsional angles between molecules that are far apart does not contribute to the histogram.  The histogram is
 averaged over the whole trajectory and output once all the trajectory frames have been read.
 
-\verbatim
+\plumedfile
 m1: MOLECULES MOL1=1,2 MOL2=3,4 MOL3=5,6 MOL4=7,8
 tt_p: INTERMOLECULARTORSIONS MOLS=m1 SWITCH={RATIONAL R_0=0.25 D_0=2.0 D_MAX=3.0}
 htt_p: HISTOGRAM DATA=tt_p GRID_MIN=-pi GRID_MAX=pi BANDWIDTH=0.1 GRID_BIN=200 STRIDE=1
 DUMPGRID GRID=htt_p FILE=myhist.out
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/MoleculeOrientation.cpp b/src/crystallization/MoleculeOrientation.cpp
index 5bade7022..6001f6f09 100644
--- a/src/crystallization/MoleculeOrientation.cpp
+++ b/src/crystallization/MoleculeOrientation.cpp
@@ -39,10 +39,10 @@ The following input tells plumed to calculate the distances between two of the a
 This is done for the same set of atoms four different molecules and the average separation is then
 calculated.
 
-\verbatim
+\plumedfile
 MOLECULES MOL1=1,2 MOL2=3,4 MOL3=5,6 MOL4=7,8 MEAN LABEL=mm
 PRINT ARG=mm.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 
 */
diff --git a/src/crystallization/Q3.cpp b/src/crystallization/Q3.cpp
index 6c3aeedc2..038115a92 100644
--- a/src/crystallization/Q3.cpp
+++ b/src/crystallization/Q3.cpp
@@ -60,28 +60,28 @@ the \f$q_{3}\f$ vectors on adjacent atoms.  More information on these variables
 The following command calculates the average Q3 parameter for the 64 atoms in a box of Lennard Jones and prints this
 quantity to a file called colvar:
 
-\verbatim
+\plumedfile
 Q3 SPECIES=1-64 D_0=1.3 R_0=0.2 MEAN LABEL=q3
 PRINT ARG=q3.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following command calculates the histogram of Q3 parameters for the 64 atoms in a box of Lennard Jones and prints these
 quantities to a file called colvar:
 
-\verbatim
+\plumedfile
 Q3 SPECIES=1-64 D_0=1.3 R_0=0.2 HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=1.0 NBINS=20 SMEAR=0.1} LABEL=q3
 PRINT ARG=q3.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following command could be used to measure the Q3 paramters that describe the arrangement of chlorine ions around the
 sodium atoms in NaCl.  The imagined system here is composed of 64 NaCl formula units and the atoms are arranged in the input
 with the 64 Na\f$^+\f$ ions followed by the 64 Cl\f$-\f$ ions.  Once again the average Q3 paramter is calculated and output to a
 file called colvar
 
-\verbatim
+\plumedfile
 Q3 SPECIESA=1-64 SPECIESB=65-128 D_0=1.3 R_0=0.2 MEAN LABEL=q3
 PRINT ARG=q3.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
@@ -136,30 +136,30 @@ adjacent atoms is correlated.
 The following command calculates the average value of the LOCAL_Q3 parameter for the 64 Lennard Jones atoms in the system under study and prints this
 quantity to a file called colvar.
 
-\verbatim
+\plumedfile
 Q3 SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=q3
 LOCAL_Q3 ARG=q3 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LABEL=lq3
 PRINT ARG=lq3.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following input calculates the distribution of LOCAL_Q3 parameters at any given time and outputs this information to a file.
 
-\verbatim
+\plumedfile
 Q3 SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=q3
 LOCAL_Q3 ARG=q3 SWITCH={RATIONAL D_0=1.3 R_0=0.2} HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=1.0 NBINS=20 SMEAR=0.1} LABEL=lq3
 PRINT ARG=lq3.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following calculates the LOCAL_Q3 parameters for atoms 1-5 only. For each of these atoms comparisons of the geometry of the coordination sphere
 are done with those of all the other atoms in the system.  The final quantity is the average and is outputted to a file
 
-\verbatim
+\plumedfile
 Q3 SPECIESA=1-5 SPECIESB=1-64 D_0=1.3 R_0=0.2 LABEL=q3a
 Q3 SPECIESA=6-64 SPECIESB=1-64 D_0=1.3 R_0=0.2 LABEL=q3b
 
 LOCAL_Q3 ARG=q3a,q3b SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LOWMEM LABEL=w3
 PRINT ARG=w3.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/Q4.cpp b/src/crystallization/Q4.cpp
index fa76b0d10..a4f7d327b 100644
--- a/src/crystallization/Q4.cpp
+++ b/src/crystallization/Q4.cpp
@@ -60,28 +60,28 @@ the \f$q_{4}\f$ vectors on adjacent atoms.  More information on these variables
 The following command calculates the average Q4 parameter for the 64 atoms in a box of Lennard Jones and prints this
 quantity to a file called colvar:
 
-\verbatim
+\plumedfile
 Q4 SPECIES=1-64 D_0=1.3 R_0=0.2 MEAN LABEL=q4
 PRINT ARG=q4.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following command calculates the histogram of Q4 parameters for the 64 atoms in a box of Lennard Jones and prints these
 quantities to a file called colvar:
 
-\verbatim
+\plumedfile
 Q4 SPECIES=1-64 D_0=1.3 R_0=0.2 HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=1.0 NBINS=20 SMEAR=0.1} LABEL=q4
 PRINT ARG=q4.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following command could be used to measure the Q4 paramters that describe the arrangement of chlorine ions around the
 sodium atoms in NaCl.  The imagined system here is composed of 64 NaCl formula units and the atoms are arranged in the input
 with the 64 Na\f$^+\f$ ions followed by the 64 Cl\f$-\f$ ions.  Once again the average Q4 paramter is calculated and output to a
 file called colvar
 
-\verbatim
+\plumedfile
 Q4 SPECIESA=1-64 SPECIESB=65-128 D_0=1.3 R_0=0.2 MEAN LABEL=q4
 PRINT ARG=q4.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
@@ -136,30 +136,30 @@ adjacent atoms is correlated.
 The following command calculates the average value of the LOCAL_Q4 parameter for the 64 Lennard Jones atoms in the system under study and prints this
 quantity to a file called colvar.
 
-\verbatim
+\plumedfile
 Q4 SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=q4
 LOCAL_Q4 ARG=q4 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LABEL=lq4
 PRINT ARG=lq4.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following input calculates the distribution of LOCAL_Q4 parameters at any given time and outputs this information to a file.
 
-\verbatim
+\plumedfile
 Q4 SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=q4
 LOCAL_Q4 ARG=q4 SWITCH={RATIONAL D_0=1.3 R_0=0.2} HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=1.0 NBINS=20 SMEAR=0.1} LABEL=lq4
 PRINT ARG=lq4.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following calculates the LOCAL_Q4 parameters for atoms 1-5 only. For each of these atoms comparisons of the geometry of the coordination sphere
 are done with those of all the other atoms in the system.  The final quantity is the average and is outputted to a file
 
-\verbatim
+\plumedfile
 Q4 SPECIESA=1-5 SPECIESB=1-64 D_0=1.3 R_0=0.2 LABEL=q4a
 Q4 SPECIESA=6-64 SPECIESB=1-64 D_0=1.3 R_0=0.2 LABEL=q4b
 
 LOCAL_Q4 ARG=q4a,q4b SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LOWMEM LABEL=w4
 PRINT ARG=w4.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/Q6.cpp b/src/crystallization/Q6.cpp
index 447b0e1cb..a848174ab 100644
--- a/src/crystallization/Q6.cpp
+++ b/src/crystallization/Q6.cpp
@@ -60,28 +60,28 @@ the \f$q_{6}\f$ vectors on adjacent atoms.  More information on these variables
 The following command calculates the average Q6 parameter for the 64 atoms in a box of Lennard Jones and prints this
 quantity to a file called colvar:
 
-\verbatim
+\plumedfile
 Q6 SPECIES=1-64 D_0=1.3 R_0=0.2 MEAN LABEL=q6
 PRINT ARG=q6.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following command calculates the histogram of Q6 parameters for the 64 atoms in a box of Lennard Jones and prints these
 quantities to a file called colvar:
 
-\verbatim
+\plumedfile
 Q6 SPECIES=1-64 D_0=1.3 R_0=0.2 HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=1.0 NBINS=20 SMEAR=0.1} LABEL=q6
 PRINT ARG=q6.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following command could be used to measure the Q6 paramters that describe the arrangement of chlorine ions around the
 sodium atoms in NaCl.  The imagined system here is composed of 64 NaCl formula units and the atoms are arranged in the input
 with the 64 Na\f$^+\f$ ions followed by the 64 Cl\f$-\f$ ions.  Once again the average Q6 paramter is calculated and output to a
 file called colvar
 
-\verbatim
+\plumedfile
 Q6 SPECIESA=1-64 SPECIESB=65-128 D_0=1.3 R_0=0.2 MEAN LABEL=q6
 PRINT ARG=q6.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
@@ -136,30 +136,30 @@ adjacent atoms is correlated.
 The following command calculates the average value of the LOCAL_Q6 parameter for the 64 Lennard Jones atoms in the system under study and prints this
 quantity to a file called colvar.
 
-\verbatim
+\plumedfile
 Q6 SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=q6
 LOCAL_Q6 ARG=q6 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LABEL=lq6
 PRINT ARG=lq6.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following input calculates the distribution of LOCAL_Q6 parameters at any given time and outputs this information to a file.
 
-\verbatim
+\plumedfile
 Q6 SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=q6
 LOCAL_Q6 ARG=q6 SWITCH={RATIONAL D_0=1.3 R_0=0.2} HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=1.0 NBINS=20 SMEAR=0.1} LABEL=lq6
 PRINT ARG=lq6.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following calculates the LOCAL_Q6 parameters for atoms 1-5 only. For each of these atoms comparisons of the geometry of the coordination sphere
 are done with those of all the other atoms in the system.  The final quantity is the average and is outputted to a file
 
-\verbatim
+\plumedfile
 Q6 SPECIESA=1-5 SPECIESB=1-64 D_0=1.3 R_0=0.2 LABEL=q6a
 Q6 SPECIESA=6-64 SPECIESB=1-64 D_0=1.3 R_0=0.2 LABEL=q6b
 
 LOCAL_Q6 ARG=q4a,q4b SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LOWMEM LABEL=w4
 PRINT ARG=w6.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/SMAC.cpp b/src/crystallization/SMAC.cpp
index ce83df738..3d7085a97 100644
--- a/src/crystallization/SMAC.cpp
+++ b/src/crystallization/SMAC.cpp
@@ -65,7 +65,7 @@ the indices of three atoms for each of the MOL keywords below we are telling PLU
 numbers to determine the orientation of the molecule that will ultimately be used when calculating the \f$\theta_{ij}\f$
 terms in the formula above.  The atom with the third index meanwhile is used when we calculate \f$r_{ij}\f$.
 
-\verbatim
+\plumedfile
 MOLECULES ...
 MOL1=9,10,9
 MOL2=89,90,89
@@ -86,13 +86,13 @@ SMAC ...
 ... SMAC
 
 PRINT ARG=s2.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 This second example works in a way that is very similar to the previous command.  Now, however,
 the orientation of the molecules is determined by finding the plane that contains the positions
 of three atoms.
 
-\verbatim
+\plumedfile
 PLANES ...
 MOL1=9,10,11
 MOL2=89,90,91
@@ -114,7 +114,7 @@ SMAC ...
 ... SMAC
 
 PRINT ARG=s2.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/SimpleCubic.cpp b/src/crystallization/SimpleCubic.cpp
index 4d50cc7d5..2c15767c2 100644
--- a/src/crystallization/SimpleCubic.cpp
+++ b/src/crystallization/SimpleCubic.cpp
@@ -55,15 +55,15 @@ so on.  Notice also that you can rotate the reference frame if you are using a n
 
 The following input tells plumed to calculate the simple cubic parameter for the atoms 1-100 with themselves.
 The mean value is then calculated.
-\verbatim
+\plumedfile
 SIMPLECUBIC SPECIES=1-100 R_0=1.0 MEAN
-\endverbatim
+\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
-\verbatim
+\plumedfile
 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}
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/crystallization/Tetrahedral.cpp b/src/crystallization/Tetrahedral.cpp
index 616cf766b..87c6e33ae 100644
--- a/src/crystallization/Tetrahedral.cpp
+++ b/src/crystallization/Tetrahedral.cpp
@@ -54,20 +54,20 @@ when atom \f$j\f$ is in the first coordination sphere of atom \f$i\f$ and is zer
 The following command calculates the average value of the tetrahedrality parameter for a set of 64 atoms all of the same type
 and outputs this quantity to a file called colvar.
 
-\verbatim
+\plumedfile
 tt: TETRAHEDRAL SPECIES=1-64 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN
 PRINT ARG=tt.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following command calculates the number of tetrahedrality 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.
 
-\verbatim
+\plumedfile
 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
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/eds/EDS.cpp b/src/eds/EDS.cpp
index e5e575abe..e7e0e4331 100644
--- a/src/eds/EDS.cpp
+++ b/src/eds/EDS.cpp
@@ -63,7 +63,7 @@ It is not possible to set the target value of the observable to zero with the de
 
 The following input for a harmonic oscillator of two beads will adaptively find a linear bias to change the mean and variance to the target values. The PRINT line shows how to access the value of the coupling constants.
 
-\verbatim
+\plumedfile
 dist: DISTANCE ATOMS=1,2
 # this is the squared of the distance
 dist2: COMBINE ARG=dist POWERS=2 PERIODIC=NO
@@ -71,19 +71,19 @@ dist2: COMBINE ARG=dist POWERS=2 PERIODIC=NO
 #bias mean and variance
 eds: EDS ARG=dist,dist2 CENTER=2.0,1.0 PERIOD=50000 TEMP=1.0
 PRINT ARG=dist,dist2,eds.dist_coupling,eds.dist2_coupling,eds.bias,eds.force2 FILE=colvars.dat STRIDE=100
-\endverbatim
+\endplumedfile
 
 Rather than trying to find the coupling constants adaptively, can ramp up to a constant value.
-\verbatim
+\plumedfile
 #ramp couplings from 0,0 to -1,1 over 50000 steps
 eds: EDS ARG=dist,dist2 CENTER=2.0,1.0 FIXED=-1,1 RAMP PERIOD=50000 TEMP=1.0
 
 #same as above, except starting at -0.5,0.5 rather than default of 0,0
 eds: EDS ARG=dist,dist2 CENTER=2.0,1.0 FIXED=-1,1 INIT=-0.5,0.5 RAMP PERIOD=50000 TEMP=1.0
-\endverbatim
+\endplumedfile
 
 A restart file can be added to dump information needed to restart/continue simulation using these parameters every STRIDE.
-\verbatim
+\plumedfile
 #add the option to write to a restart file
 eds: EDS ARG=dist,dist2 CENTER=2.0,1.0 PERIOD=50000 TEMP=1.0 OUT_RESTART=restart.dat
 
@@ -98,7 +98,7 @@ eds: EDS ARG=dist,dist2 CENTER=2.0,1.0 PERIOD=50000 TEMP=1.0 IN_RESTART=restart.
 
 #add the option to read in a previous restart file and continue the bias, but use the mean from the previous run as the starting point
 eds: EDS ARG=dist,dist2 CENTER=2.0,1.0 PERIOD=50000 TEMP=1.0 IN_RESTART=restart.dat MEAN
-\endverbatim
+\endplumedfile
 
 
 */
diff --git a/src/function/Combine.cpp b/src/function/Combine.cpp
index e76ed7b8b..62f66e35d 100644
--- a/src/function/Combine.cpp
+++ b/src/function/Combine.cpp
@@ -52,23 +52,23 @@ is periodic.
 The following input tells plumed to print the distance between atoms 3 and 5
 its square (as computed from the x,y,z components) and the distance
 again as computed from the square root of the square.
-\verbatim
+\plumedfile
 DISTANCE LABEL=dist      ATOMS=3,5 COMPONENTS
 COMBINE  LABEL=distance2 ARG=dist.x,dist.y,dist.z POWERS=2,2,2 PERIODIC=NO
 COMBINE  LABEL=distance  ARG=distance2 POWERS=0.5 PERIODIC=NO
 PRINT ARG=distance,distance2
-\endverbatim
+\endplumedfile
 (See also \ref PRINT and \ref DISTANCE).
 
 The following input tells plumed to add a restraint on the
 cube of a dihedral angle. Notice that since the angle has a
 periodic domain
 -pi,pi its cube has a domain -pi**3,pi**3.
-\verbatim
+\plumedfile
 t: TORSION ATOMS=1,3,5,7
 c: COMBINE ARG=t POWERS=3 PERIODIC=-31.0062766802998,31.0062766802998
 RESTRAINT ARG=c KAPPA=10 AT=0
-\endverbatim
+\endplumedfile
 
 
 
diff --git a/src/function/Ensemble.cpp b/src/function/Ensemble.cpp
index 9976c2501..bc1e05fb2 100644
--- a/src/function/Ensemble.cpp
+++ b/src/function/Ensemble.cpp
@@ -38,12 +38,11 @@ Each collective variable is averaged separately and stored in a component labell
 \par Examples
 The following input tells plumed to calculate the distance between atoms 3 and 5
 and the average it over the available replicas.
-\verbatim
+\plumedfile
 dist: DISTANCE ATOMS=3,5
 ens: ENSEMBLE ARG=dist
 PRINT ARG=dist,ens.dist
-\endverbatim
-(See also \ref PRINT and \ref DISTANCE).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/function/FuncPathMSD.cpp b/src/function/FuncPathMSD.cpp
index f4d481d73..25a87e718 100644
--- a/src/function/FuncPathMSD.cpp
+++ b/src/function/FuncPathMSD.cpp
@@ -49,17 +49,17 @@ It is a function of MSD that are obtained by the joint use of MSD variable and S
 Here below is a case where you have defined three frames and you want to
 calculate the progress alng the path and the distance from it in p1
 
-\verbatim
+\plumedfile
 t1: RMSD REFERENCE=frame_1.dat TYPE=OPTIMAL SQUARED
 t2: RMSD REFERENCE=frame_21.dat TYPE=OPTIMAL SQUARED
 t3: RMSD REFERENCE=frame_42.dat TYPE=OPTIMAL SQUARED
 p1: FUNCPATHMSD ARG=t1,t2,t3 LAMBDA=500.0
 PRINT ARG=t1,t2,t3,p1.s,p1.z STRIDE=1 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 In this second example is shown how to define a PATH in the \ref CONTACTMAP space:
 
-\verbatim
+\plumedfile
 CONTACTMAP ...
 ATOMS1=1,2 REFERENCE1=0.1
 ATOMS2=3,4 REFERENCE2=0.5
@@ -92,7 +92,7 @@ CMDIST
 
 p1: FUNCPATHMSD ARG=c1,c2,c3 LAMBDA=500.0
 PRINT ARG=c1,c2,c3,p1.s,p1.z STRIDE=1 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/function/LocalEnsemble.cpp b/src/function/LocalEnsemble.cpp
index 6e3bd96f6..429a5a143 100644
--- a/src/function/LocalEnsemble.cpp
+++ b/src/function/LocalEnsemble.cpp
@@ -40,7 +40,7 @@ The following input tells plumed to calculate the chemical shifts for four
 different proteins in the same simulation box then average them, calcualated
 the sum of the squared deviation with respect to the experiemntal values and
 applies a linear restraint.
-\verbatim
+\plumedfile
 MOLINFO STRUCTURE=data/template.pdb
 
 chaina: GROUP ATOMS=1-1640
@@ -68,7 +68,7 @@ sthn: STATS ARG=(enshn\.csa\.hn_.*) PARARG=(csa\.exphn_.*) SQDEVSUM
 stnh: STATS ARG=(ensnh\.csa\.nh_.*) PARARG=(csa\.expnh_.*) SQDEVSUM
 
 res: RESTRAINT ARG=stca.*,stcb.*,stco.*,sthn.*,stnh.* AT=0.,0.,0.,0.,0. KAPPA=0.,0.,0.,0.,0 SLOPE=16.,16.,12.,24.,0.5
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/function/Matheval.cpp b/src/function/Matheval.cpp
index 59e6e9bf0..d21907e68 100644
--- a/src/function/Matheval.cpp
+++ b/src/function/Matheval.cpp
@@ -54,14 +54,14 @@ PLUMED has been linked to it
 
 The following input tells plumed to perform a metadynamics
 using as a CV the difference between two distances.
-\verbatim
+\plumedfile
 dAB: DISTANCE ARG=10,12
 dAC: DISTANCE ARG=10,15
 diff: MATHEVAL ARG=dAB,dAC FUNC=y-x PERIODIC=NO
 # notice: the previous line could be replaced with the following
 # diff: COMBINE ARG=dAB,dAC COEFFICIENTS=-1,1
 METAD ARG=diff WIDTH=0.1 HEIGHT=0.5 BIASFACTOR=10 PACE=100
-\endverbatim
+\endplumedfile
 (see also \ref DISTANCE, \ref COMBINE, and \ref METAD).
 Notice that forces applied to diff will be correctly propagated
 to atoms 10, 12, and 15.
@@ -75,7 +75,7 @@ The following input tells plumed to print the angle between vectors
 identified by atoms 1,2 and atoms 2,3
 its square (as computed from the x,y,z components) and the distance
 again as computed from the square root of the square.
-\verbatim
+\plumedfile
 DISTANCE LABEL=d1 ATOMS=1,2 COMPONENTS
 DISTANCE LABEL=d2 ATOMS=2,3 COMPONENTS
 MATHEVAL ...
@@ -86,7 +86,7 @@ MATHEVAL ...
   PERIODIC=NO
 ... MATHEVAL
 PRINT ARG=theta
-\endverbatim
+\endplumedfile
 (See also \ref PRINT and \ref DISTANCE).
 
 Notice that the matheval library implements a large number of functions (trigonometric, exp, log, etc).
@@ -96,13 +96,13 @@ a straightforward implementation of if clauses.
 
 For example, imagine that you want to implement a restraint that only acts when a
 distance is larger than 0.5. You can do it with
-\verbatim
+\plumedfile
 d: DISTANCE ATOMS=10,15
 m: MATHEVAL ARG=d FUNC=0.5*step(0.5-x)+x*step(x-0.5) PERIODIC=NO
 # check the function you are applying:
 PRINT ARG=d,n FILE=checkme
 RESTRAINT ARG=d AT=0.5 KAPPA=10.0
-\endverbatim
+\endplumedfile
 (see also \ref DISTANCE, \ref PRINT, and \ref RESTRAINT)
 
 The meaning of the function `0.5*step(0.5-x)+x*step(x-0.5)` is:
@@ -133,7 +133,7 @@ MATHEVAL can be used in combination with \ref DISTANCE to implement variants of
 DISTANCE keyword that were present in PLUMED 1.3 and that allowed to compute
 the distance of a point from a line defined by two other points, or the progression
 along that line.
-\verbatim
+\plumedfile
 # take center of atoms 1 to 10 as reference point 1
 p1: CENTER ATOMS=1-10
 # take center of atoms 11 to 20 as reference point 2
@@ -156,7 +156,7 @@ fromaxis: MATHEVAL ARG=d13,d23,d12,onaxis VAR=x,y,z,o FUNC=(0.5*(y^2+x^2)-o^2-0.
 
 PRINT ARG=onaxis,fromaxis
 
-\endverbatim
+\endplumedfile
 
 Notice that these equations have been used to combine \ref RMSD
 from different snapshots of a protein so as to define
diff --git a/src/function/Piecewise.cpp b/src/function/Piecewise.cpp
index e0f4bbd52..ebd886de0 100644
--- a/src/function/Piecewise.cpp
+++ b/src/function/Piecewise.cpp
@@ -55,15 +55,14 @@ in a vector of values. Each value will be named as the name of the original
 argument with suffix _pfunc.
 
 \par Examples
-\verbatim
+\plumedfile
 dist1: DISTANCE ATOMS=1,10
 dist2: DISTANCE ATOMS=2,11
 
 pw: PIECEWISE POINT0=1,10 POINT1=1,PI POINT2=3,10 ARG=dist1
 ppww: PIECEWISE POINT0=1,10 POINT1=1,PI POINT2=3,10 ARG=dist1,dist2
 PRINT ARG=pw,ppww.dist1_pfunc,ppww.dist2_pfunc
-\endverbatim
-(See also \ref PRINT and \ref DISTANCE).
+\endplumedfile
 
 
 */
diff --git a/src/function/Sort.cpp b/src/function/Sort.cpp
index f1e9b2a71..fc8e4ca23 100644
--- a/src/function/Sort.cpp
+++ b/src/function/Sort.cpp
@@ -43,15 +43,14 @@ labelled <em>label</em>.1, the second lowest will be labelled <em>label</em>.2 a
 \par Examples
 The following input tells plumed to print the distance of the closest and of
 the farthest atoms to atom 1, chosen among atoms from 2 to 5
-\verbatim
+\plumedfile
 d12:  DISTANCE ATOMS=1,2
 d13:  DISTANCE ATOMS=1,3
 d14:  DISTANCE ATOMS=1,4
 d15:  DISTANCE ATOMS=1,5
 sort: SORT ARG=d12,d13,d14,d15
 PRINT ARG=sort.1,sort.4
-\endverbatim
-(See also \ref PRINT and \ref DISTANCE).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/function/Stats.cpp b/src/function/Stats.cpp
index 429442e6d..438d8f98b 100644
--- a/src/function/Stats.cpp
+++ b/src/function/Stats.cpp
@@ -41,13 +41,13 @@ from other actions using PARARG (for example using experimental values from coll
 The following input tells plumed to print the distance between three couple of atoms
 and compare them with three reference distances.
 
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=10,50
 d2: DISTANCE ATOMS=1,100
 d3: DISTANCE ATOMS=45,75
 st: STATS ARG=d1,d2,d3 PARAMETERS=1.5,4.0,2.0
 PRINT ARG=d1,d2,d3,st.*
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/function/Target.cpp b/src/function/Target.cpp
index 5ce912495..1724cf3ff 100644
--- a/src/function/Target.cpp
+++ b/src/function/Target.cpp
@@ -75,12 +75,12 @@ specified in the input.
 The following input calculates the distance between a reference configuration and the instaneous position of the system in the trajectory.
 The position of the reference configuration is specified by providing the values of the distance between atoms 1 and 2 and atoms 3 and 4.
 
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,2
 d2: DISTANCE ATOMS=3,4
 t1: TARGET REFERENCE=myref.pdb TYPE=EUCLIDEAN
 PRINT ARG=t1 FILE=colvar
-\endverbatim
+\endplumedfile
 
 The contents of the file containing the reference structure (myref.pdb) is shown below.  As you can see you must provide information on the
 labels of the CVs that are being used to define the position of the reference configuration in this file together with the values that these
diff --git a/src/generic/Debug.cpp b/src/generic/Debug.cpp
index 742b3b854..36654b0d9 100644
--- a/src/generic/Debug.cpp
+++ b/src/generic/Debug.cpp
@@ -37,12 +37,12 @@ Can be used while debugging or optimizing plumed.
 
 \par Examples
 
-\verbatim
+\plumedfile
 # print detailed (action-by-action) timers at the end of simulation
 DEBUG DETAILED_TIMERS
 # dump every two steps which are the atoms required from the MD code
 DEBUG logRequestedAtoms STRIDE=2
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/DumpAtoms.cpp b/src/generic/DumpAtoms.cpp
index d6ee37748..c144b5f07 100644
--- a/src/generic/DumpAtoms.cpp
+++ b/src/generic/DumpAtoms.cpp
@@ -70,15 +70,14 @@ Notice that gro/xtc/trr files can only contain coordinates in nm.
 The following input instructs plumed to print out the positions of atoms
 1-10 together with the position of the center of mass of atoms 11-20 every
 10 steps to a file called file.xyz.
-\verbatim
+\plumedfile
 COM ATOMS=11-20 LABEL=c1
 DUMPATOMS STRIDE=10 FILE=file.xyz ATOMS=1-10,c1
-\endverbatim
-(see also \ref COM)
+\endplumedfile
 
 The following input is very similar but dumps a .gro (gromacs) file,
 which also contains atom and residue names.
-\verbatim
+\plumedfile
 # this is required to have proper atom names:
 MOLINFO STRUCTURE=reference.pdb
 # if omitted, atoms will have "X" name...
@@ -87,8 +86,7 @@ COM ATOMS=11-20 LABEL=c1
 DUMPATOMS STRIDE=10 FILE=file.gro ATOMS=1-10,c1
 # notice that last atom is a virtual one and will not have
 # a correct name in the resulting gro file
-\endverbatim
-(see also \ref COM and \ref MOLINFO)
+\endplumedfile
 
 
 */
diff --git a/src/generic/DumpDerivatives.cpp b/src/generic/DumpDerivatives.cpp
index c74343aca..38c686d7b 100644
--- a/src/generic/DumpDerivatives.cpp
+++ b/src/generic/DumpDerivatives.cpp
@@ -42,11 +42,11 @@ can be done by outputting the derivatives calculated analytically and numericall
 \par Examples
 The following input instructs plumed to write a file called deriv that contains both the
 analytical and numerical derivatives of the distance between atoms 1 and 2.
-\verbatim
+\plumedfile
 DISTANCE ATOM=1,2 LABEL=distance
 DISTANCE ATOM=1,2 LABEL=distanceN NUMERICAL_DERIVATIVES
 DUMPDERIVATIVES ARG=distance,distanceN STRIDE=1 FILE=deriv
-\endverbatim
+\endplumedfile
 
 (See also \ref DISTANCE)
 
diff --git a/src/generic/DumpForces.cpp b/src/generic/DumpForces.cpp
index f5e204d65..a376de281 100644
--- a/src/generic/DumpForces.cpp
+++ b/src/generic/DumpForces.cpp
@@ -43,12 +43,10 @@ by specifying more than one argument. You can control the buffering of output us
 \par Examples
 The following input instructs plumed to write a file called forces that contains
 the force acting on the distance between atoms 1 and 2.
-\verbatim
+\plumedfile
 DISTANCE ATOM=1,2 LABEL=distance
 DUMPFORCES ARG=distance STRIDE=1 FILE=forces
-\endverbatim
-
-(See also \ref DISTANCE)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/DumpMassCharge.cpp b/src/generic/DumpMassCharge.cpp
index 4153bc838..25224f3b3 100644
--- a/src/generic/DumpMassCharge.cpp
+++ b/src/generic/DumpMassCharge.cpp
@@ -49,14 +49,13 @@ masses for all atoms are written.
 You can add the DUMPMASSCHARGE action at the end of the plumed.dat
 file that you use during an MD simulations:
 
-\verbatim
+\plumedfile
 c1: COM ATOMS=1-10
 c2: COM ATOMS=11-20
 PRINT ARG=c1,c2 FILE=colvar STRIDE=100
 
 DUMPMASSCHARGE FILE=mcfile
-\endverbatim
-(see also \ref COM and \ref PRINT)
+\endplumedfile
 
 In this way, you will be able to use the same masses while processing
 a trajectory from the \ref driver . To do so, you need to
@@ -67,11 +66,11 @@ plumed driver --mc mcfile --plumed plumed.dat --ixyz traj.xyz
 
 With the following input you can dump only the charges for a specific
 group.
-\verbatim
+\plumedfile
 solute_ions: GROUP ATOMS=1-121,200-2012
 DUMPATOMS FILE=traj.gro ATOMS=solute_ions STRIDE=100
 DUMPMASSCHARGE FILE=mcfile ATOMS=solute_ions
-\endverbatim
+\endplumedfile
 Notice however that if you want to process the charges
 with the driver (e.g. reading traj.gro) you have to fix atom
 numbers first, e.g. with the script
diff --git a/src/generic/DumpProjections.cpp b/src/generic/DumpProjections.cpp
index 3f0d40429..c3a0ca750 100644
--- a/src/generic/DumpProjections.cpp
+++ b/src/generic/DumpProjections.cpp
@@ -39,12 +39,12 @@ Dump the derivatives with respect to the input parameters for one or more object
 Compute the distance between two groups and write on a file the
 derivatives of this distance with respect to all the atoms of the two groups
 
-\verbatim
+\plumedfile
 x1: CENTER ATOMS=1-10
 x2: CENTER ATOMS=11-20
 d: DISTANCE ATOMS=x1,x2
 DUMPPROJECTIONS ARG=d FILE=proj STRIDE=20
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/EffectiveEnergyDrift.cpp b/src/generic/EffectiveEnergyDrift.cpp
index 1c3961390..a9ca15e6c 100644
--- a/src/generic/EffectiveEnergyDrift.cpp
+++ b/src/generic/EffectiveEnergyDrift.cpp
@@ -55,18 +55,18 @@ Print the effective energy drift described in Ref \cite Ferrarotti2015
 This is to monitor the effective energy drift for a metadynamics
 simulation on the Debye-Huckel energy. Since this variable is very expensive,
 it could be conveniently computed every second step.
-\verbatim
+\plumedfile
 dh: DHENERGY GROUPA=1-10 GROUPB=11-20 EPSILON=80.0 I=0.1 TEMP=300.0
 METAD ARG=dh HEIGHT=0.5 SIGMA=0.1 PACE=500 STRIDE=2
 EFFECTIVE_ENERGY_DRIFT PRINT_STRIDE=100 FILE=eff
-\endverbatim
+\endplumedfile
 
 This is to monitor if a restraint is too stiff
-\verbatim
+\plumedfile
 d: DISTANCE ATOMS=10,20
 RESTRAINT ARG=d KAPPA=100000 AT=0.6
 EFFECTIVE_ENERGY_DRIFT PRINT_STRIDE=100 FILE=eff
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/FitToTemplate.cpp b/src/generic/FitToTemplate.cpp
index 37e38333f..6fa379e1f 100644
--- a/src/generic/FitToTemplate.cpp
+++ b/src/generic/FitToTemplate.cpp
@@ -78,13 +78,12 @@ this action is performed at every MD step.
 \par Examples
 
 Align the atomic position to a template then print them
-\verbatim
+\plumedfile
 # to see the effect, one could dump the atoms before alignment
 DUMPATOMS FILE=dump-before.xyz ATOMS=1-20
 FIT_TO_TEMPLATE STRIDE=1 REFERENCE=ref.pdb TYPE=SIMPLE
 DUMPATOMS FILE=dump-after.xyz ATOMS=1-20
-\endverbatim
-(see also \ref DUMPATOMS)
+\endplumedfile
 
 
 
diff --git a/src/generic/Flush.cpp b/src/generic/Flush.cpp
index 78c00f87d..71e9a06d3 100644
--- a/src/generic/Flush.cpp
+++ b/src/generic/Flush.cpp
@@ -41,7 +41,7 @@ plumed input file, it will flush all the open files.
 
 \par Examples
 A command like this in the input will instruct plumed to flush all the output files every 100 steps
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,10
 PRINT ARG=d1 STRIDE=5 FILE=colvar1
 
@@ -50,7 +50,7 @@ FLUSH STRIDE=100
 d2: DISTANCE ATOMS=2,11
 # also this print is flushed every 100 steps:
 PRINT ARG=d2 STRIDE=10 FILE=colvar2
-\endverbatim
+\endplumedfile
 (see also \ref DISTANCE and \ref PRINT).
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/Group.cpp b/src/generic/Group.cpp
index 80393a9ce..0f04cf75f 100644
--- a/src/generic/Group.cpp
+++ b/src/generic/Group.cpp
@@ -65,7 +65,7 @@ the \ref INCLUDE command so as to store long group definitions in a separate fil
 
 This command create a group of atoms containing atoms 1,4,7,11 and 14 (labeled 'o'), and another containing
 atoms 2,3,5,6,8,9,12,13 (labeled 'h'):
-\verbatim
+\plumedfile
 o: GROUP ATOMS=1,4,7,11,14
 h: GROUP ATOMS=2,3,5,6,8,9,12,13
 # compute the coordination among the two groups
@@ -75,45 +75,41 @@ c: COORDINATION GROUPA=o GROUPB=h R_0=0.3
 
 # print the coordination on file 'colvar'
 PRINT ARG=c FILE=colvar
-\endverbatim
-(see also \ref COORDINATION and \ref PRINT)
+\endplumedfile
 
 Groups can be conveniently stored in a separate file.
 E.g. one could create a file named 'groups.dat' which reads
-\verbatim
+\plumedfile
 o: GROUP ATOMS=1,4,7,11,14
 h: GROUP ATOMS=2,3,5,6,8,9,12,13
-\endverbatim
+\endplumedfile
 and then include it in the main 'plumed.dat' file
-\verbatim
+\plumedfile
 INCLUDE FILE=groups.dat
 # compute the coordination among the two groups
 c: COORDINATION GROUPA=o GROUPB=h R_0=0.3
 # print the coordination on file 'colvar'
 PRINT ARG=c FILE=colvar
-\endverbatim
-(see also \ref INCLUDE, \ref COORDINATION, and \ref PRINT).
+\endplumedfile
 The groups.dat file could be very long and include lists of thousand atoms without cluttering the main plumed.dat file.
 
 A GROMACS index file can also be imported
-\verbatim
+\plumedfile
 # import group named 'protein' from file index.ndx
 pro: GROUP NDX_FILE=index.ndx NDX_GROUP=protein
 # dump all the atoms of the protein on a trajectory file
 DUMPATOMS ATOMS=pro FILE=traj.gro
-\endverbatim
-(see also \ref DUMPATOMS)
+\endplumedfile
 
 A list can be edited with REMOVE
-\verbatim
+\plumedfile
 # take one atom every three
 ox: GROUP ATOMS=1-90:3
 # take the remaining atoms
 hy: GROUP ATOMS=1-90 REMOVE=ox
 DUMPATOMS ATOMS=ox FILE=ox.gro
 DUMPATOMS ATOMS=hy FILE=hy.gro
-\endverbatim
-(see also \ref DUMPATOMS)
+\endplumedfile
 
 
 */
diff --git a/src/generic/Include.cpp b/src/generic/Include.cpp
index 664c84ba6..30539f0bb 100644
--- a/src/generic/Include.cpp
+++ b/src/generic/Include.cpp
@@ -38,25 +38,25 @@ Useful to split very large plumed.dat files.
 \par Examples
 
 This input
-\verbatim
+\plumedfile
 c1: COM ATOMS=1-100
 c2: COM ATOMS=101-202
 d: DISTANCE ARG=c1,c2
 PRINT ARG=d
-\endverbatim
+\endplumedfile
 
 can be replaced with
-\verbatim
+\plumedfile
 INCLUDE FILE=pippo.dat
 d: DISTANCE ARG=c1,c2
 PRINT ARG=d
-\endverbatim
+\endplumedfile
 
 where the content of file pippo.dat is
-\verbatim
+\plumedfile
 c1: COM ATOMS=1-100
 c2: COM ATOMS=101-202
-\endverbatim
+\endplumedfile
 
 (see also \ref COM, \ref DISTANCE, and \ref PRINT).
 
diff --git a/src/generic/Print.cpp b/src/generic/Print.cpp
index 2a87d1f7b..74245cf99 100644
--- a/src/generic/Print.cpp
+++ b/src/generic/Print.cpp
@@ -40,13 +40,12 @@ to different files.  You can control the buffering of output using the \subpage
 The following input instructs plumed to print the distance between atoms 3 and 5 on a file
 called COLVAR every 10 steps, and the distance and total energy on a file called COLVAR_ALL
 every 1000 steps.
-\verbatim
+\plumedfile
 DISTANCE ATOMS=2,5 LABEL=distance
 ENERGY             LABEL=energy
 PRINT ARG=distance          STRIDE=10   FILE=COLVAR
 PRINT ARG=distance,energy   STRIDE=1000 FILE=COLVAR_ALL
-\endverbatim
-(See also \ref DISTANCE and \ref ENERGY).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/RandomExchanges.cpp b/src/generic/RandomExchanges.cpp
index be36873bc..60d207120 100644
--- a/src/generic/RandomExchanges.cpp
+++ b/src/generic/RandomExchanges.cpp
@@ -44,25 +44,25 @@ metadynamics simulation using a different angle in each replica.
 Exchanges will be randomly tried between replicas 0-1, 0-2 and 1-2
 
 Here is plumed.0.dat
-\verbatim
+\plumedfile
 RANDOM_EXCHANGES
 t: TORSION ATOMS=1,2,3,4
 METAD ARG=t HEIGHT=0.1 PACE=100 SIGMA=0.3
-\endverbatim
+\endplumedfile
 
 Here is plumed.1.dat
-\verbatim
+\plumedfile
 RANDOM_EXCHANGES
 t: TORSION ATOMS=2,3,4,5
 METAD ARG=t HEIGHT=0.1 PACE=100 SIGMA=0.3
-\endverbatim
+\endplumedfile
 
 Here is plumed.2.dat
-\verbatim
+\plumedfile
 RANDOM_EXCHANGES
 t: TORSION ATOMS=3,4,5,6
 METAD ARG=t HEIGHT=0.1 PACE=100 SIGMA=0.3
-\endverbatim
+\endplumedfile
 
 \warning Multi replica simulations are presently only working with gromacs.
 
diff --git a/src/generic/Read.cpp b/src/generic/Read.cpp
index 88c279630..63b8f8fc8 100644
--- a/src/generic/Read.cpp
+++ b/src/generic/Read.cpp
@@ -52,11 +52,11 @@ This input reads in data from a file called input_colvar.data that was generated
 in a calculation that involved PLUMED.  The first command reads in the data from the
 column headed phi1 while the second reads in the data from the column headed phi2.
 
-\verbatim
+\plumedfile
 rphi1:       READ FILE=input_colvar.data  VALUES=phi1
 rphi2:       READ FILE=input_colvar.data  VALUES=phi2
 PRINT ARG=rphi1,rphi2 STRIDE=500  FILE=output_colvar.data
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/ResetCell.cpp b/src/generic/ResetCell.cpp
index fd14c1fe1..5268d52de 100644
--- a/src/generic/ResetCell.cpp
+++ b/src/generic/ResetCell.cpp
@@ -71,14 +71,13 @@ this action is performed at every MD step.
 \par Examples
 
 Reset cell to be triangular after a rototranslational fit
-\verbatim
+\plumedfile
 DUMPATOMS FILE=dump-original.xyz ATOMS=1-20
 FIT_TO_TEMPLATE STRIDE=1 REFERENCE=ref.pdb TYPE=OPTIMAL
 DUMPATOMS FILE=dump-fit.xyz ATOMS=1-20
 RESET_CELL TYPE=TRIANGULAR
 DUMPATOMS FILE=dump-reset.xyz ATOMS=1-20
-\endverbatim
-(see also \ref DUMPATOMS)
+\endplumedfile
 
 
 */
diff --git a/src/generic/Time.cpp b/src/generic/Time.cpp
index 3512cb53d..9a821f0b6 100644
--- a/src/generic/Time.cpp
+++ b/src/generic/Time.cpp
@@ -35,11 +35,10 @@ retrieve the time of the simulation to be used elsewere
 
 \par Examples
 
-\verbatim
+\plumedfile
 TIME            LABEL=t1
 PRINT ARG=t1
-\endverbatim
-(See also \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/UpdateIf.cpp b/src/generic/UpdateIf.cpp
index 2d067fcf3..696e66c20 100644
--- a/src/generic/UpdateIf.cpp
+++ b/src/generic/UpdateIf.cpp
@@ -62,7 +62,7 @@ can lead to unexpected results.
 \par Examples
 The following input instructs plumed dump all the snapshots where an atom is in touch with
 the solute.
-\verbatim
+\plumedfile
 solute: GROUP ATOMS=1-124
 coord: COORDINATION GROUPA=solute GROUPB=500 R_0=0.5
 
@@ -72,8 +72,7 @@ coord: COORDINATION GROUPA=solute GROUPB=500 R_0=0.5
 UPDATE_IF ARG=coord MORE_THAN=0.5
 DUMPATOMS ATOMS=solute,500 FILE=output.xyz
 UPDATE_IF ARG=coord END
-\endverbatim
-(See also \ref GROUP, \ref COORDINATION, and \ref DUMPATOMS)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/WholeMolecules.cpp b/src/generic/WholeMolecules.cpp
index 12e9c8281..e43c0aa07 100644
--- a/src/generic/WholeMolecules.cpp
+++ b/src/generic/WholeMolecules.cpp
@@ -74,30 +74,27 @@ to skip some atoms, provided consecute chosen atoms are close enough.
 This command instructs plumed to reconstruct the molecule containing atoms 1-20
 at every step of the calculation and dump them on a file.
 
-\verbatim
+\plumedfile
 # to see the effect, one could dump the atoms as they were before molecule reconstruction:
 # DUMPATOMS FILE=dump-broken.xyz ATOMS=1-20
 WHOLEMOLECULES ENTITY0=1-20
 DUMPATOMS FILE=dump.xyz ATOMS=1-20
-\endverbatim
-(see also \ref DUMPATOMS)
+\endplumedfile
 
 This command instructs plumed to reconstruct two molecules containing atoms 1-20 and 30-40
 
-\verbatim
+\plumedfile
 WHOLEMOLECULES ENTITY0=1-20 ENTITY1=30-40
 DUMPATOMS FILE=dump.xyz ATOMS=1-20,30-40
-\endverbatim
-(see also \ref DUMPATOMS)
+\endplumedfile
 
 This command instructs plumed to reconstruct the chain of backbone atoms in a
 protein
 
-\verbatim
+\plumedfile
 MOLINFO STRUCTURE=helix.pdb
 WHOLEMOLECULES RESIDUES=all MOLTYPE=protein
-\endverbatim
-(See also \ref MOLINFO)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/generic/WrapAround.cpp b/src/generic/WrapAround.cpp
index 946a060fd..6320ee4bb 100644
--- a/src/generic/WrapAround.cpp
+++ b/src/generic/WrapAround.cpp
@@ -69,20 +69,19 @@ consider the possibility of using the STRIDE keyword here (with great care).
 This command instructs plumed to move all the ions to their periodic image that is as close as possible to
 the rna group.
 
-\verbatim
+\plumedfile
 rna: GROUP ATOMS=1-100
 ions: GROUP ATOMS=101-110
 # first make the rna molecule whole
 WHOLEMOLECULES ENTITY0=rna
 WRAPAROUND ATOMS=ions AROUND=rna
 DUMPATOMS FILE=dump.xyz ATOMS=rna,ions
-\endverbatim
-(see also \ref WHOLEMOLECULES, \ref GROUP and \ref DUMPATOMS)
+\endplumedfile
 
 In case you want to do it during a simulation and you only care about wrapping the ions in
 the `dump.xyz` file, you can use the following:
 
-\verbatim
+\plumedfile
 # add some restraint that do not require molecules to be whole:
 a: TORSION ATOMS=1,2,10,11
 RESTRAINT ARG=a AT=0.0 KAPPA=5
@@ -98,15 +97,14 @@ ions: GROUP ATOMS=101-110
 WHOLEMOLECULES ENTITY0=rna STRIDE=100
 WRAPAROUND ATOMS=ions AROUND=rna STRIDE=100
 DUMPATOMS FILE=dump.xyz ATOMS=rna,ions STRIDE=100
-\endverbatim
-(see also \ref TORSION, \ref GROUP, \ref WHOLEMOLECULES and \ref DUMPATOMS)
+\endplumedfile
 
 Notice that if the biased variable requires a molecule to be whole, you might have to put
 just the \ref WHOLEMOLECULES command before computing that variable and leave the default STRIDE=1.
 
 This command instructs plumed to center all atoms around the center of mass of a solute molecule.
 
-\verbatim
+\plumedfile
 solute: GROUP ATOMS=1-100
 all: GROUP ATOMS=1-1000
 # center of the solute:
@@ -116,8 +114,7 @@ com: COM ATOMS=solute
 # notice that we wrap around a single atom. this should be fast
 WRAPAROUND ATOMS=all AROUND=com
 DUMPATOMS FILE=dump.xyz ATOMS=all
-\endverbatim
-(see also \ref GROUP \ref COM \ref DUMPATOMS)
+\endplumedfile
 
 Notice that whereas \ref WHOLEMOLECULES is designed to make molecules whole,
 \ref WRAPAROUND can easily break molecules. In the last example,
@@ -136,7 +133,7 @@ in the following examples all the water oxygens will be brought
 close to the solute, and all the hydrogens will be kept close
 to their related oxygen.
 
-\verbatim
+\plumedfile
 solute: GROUP ATOMS=1-100
 water: GROUP ATOMS=101-1000
 com: COM ATOMS=solute
@@ -145,8 +142,7 @@ WRAPAROUND ATOMS=solute AROUND=com
 # notice that we wrap around a single atom. this should be fast
 WRAPAROUND ATOMS=water AROUND=com GROUPBY=3
 DUMPATOMS FILE=dump.xyz ATOMS=solute,water
-\endverbatim
-(see also \ref GROUP \ref COM \ref DUMPATOMS)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/gridtools/ConvertToFES.cpp b/src/gridtools/ConvertToFES.cpp
index 4e3bbdf17..661aa2bde 100644
--- a/src/gridtools/ConvertToFES.cpp
+++ b/src/gridtools/ConvertToFES.cpp
@@ -46,12 +46,12 @@ and the HISTOGRAM action.  All the data within this trajectory is used in the co
 HISTOGRAM.  Finally, once all the data has been read in, the histogram is converted to a free energy
 using the formula above and the free energy is output to a file called fes.dat
 
-\verbatim
+\plumedfile
 x: DISTANCE ATOMS=1,2
 hA1: HISTOGRAM ARG=x GRID_MIN=0.0 GRID_MAX=3.0 GRID_BIN=100 BANDWIDTH=0.1
 ff: CONVERT_TO_FES GRID=hA1 TEMP=300
 DUMPGRID GRID=ff FILE=fes.dat
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/gridtools/DumpCube.cpp b/src/gridtools/DumpCube.cpp
index 0fbc47f56..f1df7406e 100644
--- a/src/gridtools/DumpCube.cpp
+++ b/src/gridtools/DumpCube.cpp
@@ -46,14 +46,14 @@ between atoms 1 and 2, the distance between atom 1 and 3 and the angle between t
 all the kernels have been added the resulting histogram is output to a file called histoA1.cube.  This file has the
 Gaussian cube file format.  The histogram can thus be visualized using tools such as VMD.
 
-\verbatim
+\plumedfile
 x1: DISTANCE ATOMS=1,2
 x2: DISTANCE ATOMS=1,3
 x3: ANGLE ATOMS=1,2,3
 
 hA1: HISTOGRAM ARG=x1,x2,x3 GRID_MIN=0.0,0.0,0.0 GRID_MAX=3.0,3.0,3.0 GRID_BIN=10,10,10 BANDWIDTH=1.0,1.0,1.0
 DUMPCUBE GRID=hA1 FILE=histoA1.cube
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/gridtools/DumpGrid.cpp b/src/gridtools/DumpGrid.cpp
index 608742a8b..e3524d563 100644
--- a/src/gridtools/DumpGrid.cpp
+++ b/src/gridtools/DumpGrid.cpp
@@ -67,7 +67,7 @@ for y.  This block is then followed by a blank line again and this pattern conti
 
 The following input monitors two torsional angles during a simulation
 and outputs a continuos histogram as a function of them at the end of the simulation.
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 HISTOGRAM ...
@@ -80,11 +80,11 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hh FILE=histo
-\endverbatim
+\endplumedfile
 
 The following input monitors two torsional angles during a simulation
 and outputs a discrete histogram as a function of them at the end of the simulation.
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 HISTOGRAM ...
@@ -98,11 +98,11 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hh FILE=histo
-\endverbatim
+\endplumedfile
 
 The following input monitors two torsional angles during a simulation
 and outputs the histogram accumulated thus far every 100000 steps.
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 HISTOGRAM ...
@@ -115,14 +115,14 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hh FILE=histo STRIDE=100000
-\endverbatim
+\endplumedfile
 
 The following input monitors two torsional angles during a simulation
 and outputs a separate histogram for each 100000 steps worth of trajectory.
 Notice how the CLEAR keyword is used here and how it is not used in the
 previous example.
 
-\verbatim
+\plumedfile
 TORSION ATOMS=1,2,3,4 LABEL=r1
 TORSION ATOMS=2,3,4,5 LABEL=r2
 HISTOGRAM ...
@@ -136,7 +136,7 @@ HISTOGRAM ...
 ... HISTOGRAM
 
 DUMPGRID GRID=hh FILE=histo STRIDE=100000
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/gridtools/FindContour.cpp b/src/gridtools/FindContour.cpp
index a58e4f650..1bae86d43 100644
--- a/src/gridtools/FindContour.cpp
+++ b/src/gridtools/FindContour.cpp
@@ -67,7 +67,7 @@ solidness at each point in the simulation cell.  The Willard-Chandler dividing s
 at which the value of this phase field is equal to 0.5.  This set of points is output to file called mycontour.dat.  A new contour
 is found on every single step for each frame that is read in.
 
-\verbatim
+\plumedfile
 UNITS NATURAL
 FCCUBIC ...
   SPECIES=1-96000 SWITCH={CUBIC D_0=1.2 D_MAX=1.5}
@@ -84,7 +84,7 @@ MULTICOLVARDENS ...
 ... MULTICOLVARDENS
 
 FIND_CONTOUR GRID=dens CONTOUR=0.5 FILE=mycontour.dat
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/gridtools/FindSphericalContour.cpp b/src/gridtools/FindSphericalContour.cpp
index a2ee2aa86..7d36e818c 100644
--- a/src/gridtools/FindSphericalContour.cpp
+++ b/src/gridtools/FindSphericalContour.cpp
@@ -81,7 +81,7 @@ by exploiting the functionality within \ref CENTER_OF_MULTICOLVAR.  We can then
 cluster and the values of the coordination numbers of these atoms.  The final line in the input then finds the a set of points on the dividing surface that separates
 teh droplet from the surrounding gas.  The value of the phase field on this isocontour is equal to 0.75.
 
-\verbatim
+\plumedfile
 # Calculate coordination numbers
 c1: COORDINATIONNUMBER SPECIES=1-512 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0}
 # Select coordination numbers that are more than 2.0
@@ -98,7 +98,7 @@ cent: CENTER_OF_MULTICOLVAR DATA=trans1
 dens: MULTICOLVARDENS DATA=trans1 ORIGIN=cent DIR=xyz NBINS=30,30,30 BANDWIDTH=2.0,2.0,2.0
 # Find the isocontour around the nucleus
 FIND_SPHERICAL_CONTOUR GRID=dens CONTOUR=0.85 INNER_RADIUS=10.0 OUTER_RADIUS=40.0 FILE=mysurface.xyz UNITS=A PRECISION=4 NPOINTS=100
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/gridtools/FourierTransform.cpp b/src/gridtools/FourierTransform.cpp
index 5a53b390d..87a078336 100644
--- a/src/gridtools/FourierTransform.cpp
+++ b/src/gridtools/FourierTransform.cpp
@@ -57,9 +57,9 @@ The default values of these parameters are: \f$a=1\f$ and \f$b=1\f$.
 \par Examples
 The following example tells Plumed to compute the complex 2D 'backward' Discrete Fourier Transform by taking the data saved on a grid called 'density', and normalizing the output by \f$ \frac{1}{\sqrt{N_x\, N_y}}\f$, where \f$N_x\f$ and \f$N_y\f$ are the number of data on the grid (it can be the case that \f$N_x \neq N_y\f$):
 
-\verbatim
+\plumedfile
 FOURIER_TRANSFORM STRIDE=1 GRID=density FT_TYPE=complex FOURIER_PARAMETERS=0,-1 FILE=fourier.dat
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/gridtools/InterpolateGrid.cpp b/src/gridtools/InterpolateGrid.cpp
index e8f97c0d4..7fbc1d65d 100644
--- a/src/gridtools/InterpolateGrid.cpp
+++ b/src/gridtools/InterpolateGrid.cpp
@@ -38,12 +38,12 @@ distance between atoms 1 and 2 using kernel density estimation.  During the calc
 are evaluated at 100 points on a uniform grid between 0.0 and 3.0.  Prior to outputting this function at the end of the
 simulation this function is interpolated onto a finer grid of 200 points between 0.0 and 3.0.
 
-\verbatim
+\plumedfile
 x: DISTANCE ATOMS=1,2
 hA1: HISTOGRAM ARG=x GRID_MIN=0.0 GRID_MAX=3.0 GRID_BIN=100 BANDWIDTH=0.1
 ii: INTERPOLATE_GRID GRID=hA1 GRID_BIN=200
 DUMPGRID GRID=ii FILE=histo.dat
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/maketools/plumedcheck b/src/maketools/plumedcheck
index 120808dec..5886ec086 100755
--- a/src/maketools/plumedcheck
+++ b/src/maketools/plumedcheck
@@ -354,6 +354,12 @@ BEGINFILE{
     if(in_plumed_doc && (in_plumed_doc in provide_examples) && match($0,"^ *\\\\verbatim *$")){
       provide_verbatim[in_plumed_doc]=1
     }
+# check if, besides the \par Examples text, there is a plumedfile command
+# now it is considered equivalent to a verbatim
+# later on we might make plumedfile compulsory instead of verbatim
+    if(in_plumed_doc && (in_plumed_doc in provide_examples) && match($0,"^ *\\\\plumedfile *$")){
+      provide_verbatim[in_plumed_doc]=1
+    }
 #print match($0,"^ *\\\\verbatim *$"),"X" $0 "X"
 
 # take note of actions or cltools that provide examples
diff --git a/src/manyrestraints/LWalls.cpp b/src/manyrestraints/LWalls.cpp
index 7c0ad569e..c98405d04 100644
--- a/src/manyrestraints/LWalls.cpp
+++ b/src/manyrestraints/LWalls.cpp
@@ -45,10 +45,10 @@ The following set of commands can be used to stop any of the 800 atoms in group
 in the z direction from atom 34137.  This is done by adding a lower wall on the z-distance between all the atoms in
 group A and the position of 34137.
 
-\verbatim
+\plumedfile
 l: ZDISTANCES GROUPA=1-800 GROUPB=34137 NOPBC
 LWALLS DATA=l AT=2.46465 KAPPA=150.0 EXP=2 EPS=1 OFFSET=0 LABEL=lwall
-\endverbatim
+\endplumedfile
 
 
 */
diff --git a/src/manyrestraints/UWalls.cpp b/src/manyrestraints/UWalls.cpp
index f4c75f0d8..2c0f75367 100644
--- a/src/manyrestraints/UWalls.cpp
+++ b/src/manyrestraints/UWalls.cpp
@@ -48,11 +48,11 @@ and the center of mass of the cluster.  These distances are then passed to the U
 a \ref UPPER_WALLS restraint on each of them and thereby prevents each of them from moving very far from the centre
 of mass of the cluster.
 
-\verbatim
+\plumedfile
 COM ATOMS=1-20 LABEL=c1
 DISTANCES GROUPA=c1 GROUPB=1-20 LABEL=d1
 UWALLS DATA=d1 AT=2.5 KAPPA=0.2 LABEL=sr
-\endverbatim
+\endplumedfile
 
 
 */
diff --git a/src/mapping/PCAVars.cpp b/src/mapping/PCAVars.cpp
index d029720ed..20ae00301 100644
--- a/src/mapping/PCAVars.cpp
+++ b/src/mapping/PCAVars.cpp
@@ -76,10 +76,10 @@ frame are removed from these displacements.  The matrix \f$A\f$ and the referenc
 configuration \f$R^{ref}\f$ are specified in the pdb input file reference.pdb and the
 value of all projections (and the residual) are output to a file called colvar2.
 
-\verbatim
+\plumedfile
 PCAVARS REFERENCE=reference.pdb TYPE=OPTIMAL LABEL=pca2
 PRINT ARG=pca2.* FILE=colvar2
-\endverbatim
+\endplumedfile
 
 The reference configurations can be specified using a pdb file.  The first configuration that you provide is the reference configuration,
 which is refered to in the above as \f$X^{ref}\f$ subsequent configurations give the directions of row vectors that are contained in
diff --git a/src/mapping/Path.cpp b/src/mapping/Path.cpp
index 37db777d2..4808fad15 100644
--- a/src/mapping/Path.cpp
+++ b/src/mapping/Path.cpp
@@ -73,21 +73,21 @@ In the example below the path is defined using RMSD distance from frames.
 The reference frames in the path are defined in the pdb file.  In this frame
 each configuration in the path is separated by a line containing just the word END.
 
-\verbatim
+\plumedfile
 p1: PATH REFERENCE=file.pdb TYPE=OPTIMAL LAMBDA=500.0
 PRINT ARG=p1.sss,p1.zzz STRIDE=1 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 In the example below the path is defined using the values of two torsional angles (t1 and t2).
 In addition, the \f$s\f$ and \f$z\f$ are calculated using the geometric expressions described
 above rather than the alegebraic expressions that are used by default.
 
-\verbatim
+\plumedfile
 t1: TORSION ATOMS=5,7,9,15
 t2: TORSION ATOMS=7,9,15,17
 pp: PATH TYPE=EUCLIDEAN REFERENCE=epath.pdb GPATH NOSPATH NOZPATH
 PRINT ARG=pp.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 Notice that the LAMBDA parameter is not required here as we are not calculating \f$s\f$ and \f$s\f$
 using the algebraic formulae defined earlier.  The positions of the frames in the path are defined
@@ -109,19 +109,19 @@ The following input instructs PLUMED to calculate the values of the path collect
 path are defined in the file all.pdb and all distances are measured using the OPTIMAL metric that is discussed in the manual
 page on \ref RMSD.
 
-\verbatim
+\plumedfile
 p2: PATH REFERENCE=all.pdb LAMBDA=69087
 PRINT ARG=p2.spath,p2.zpath STRIDE=1 FILE=colvar
-\endverbatim
+\endplumedfile
 
 If you wish to use collective variable values in the definition of your path you would use an input file with something like this:
 
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,2
 d2: DISTANCE ATOMS=3,4a
 p2: PATH REFERENCE=mypath.pdb LAMBDA=2 TYPE=EUCLIDEAN
 PRINT ARG=p2.spath,p2.zpath STRIDE=1 FILE=colvar
-\endverbatim
+\endplumedfile
 
 The corresponding pdb file containing the  definitions of the frames in the path would then look like this:
 
diff --git a/src/mapping/PropertyMap.cpp b/src/mapping/PropertyMap.cpp
index 9da32b47b..91c96f4dd 100644
--- a/src/mapping/PropertyMap.cpp
+++ b/src/mapping/PropertyMap.cpp
@@ -47,10 +47,10 @@ that these properties take at a set of reference configurations and using the fo
 between the reference configurations and the instantaneous configurations are calculated using the OPTIMAL metric that is
 discussed at length in the manual pages on \ref RMSD.
 
-\verbatim
+\plumedfile
 p2: GPROPERTYMAP REFERENCE=allv.pdb PROPERTY=X,Y LAMBDA=69087
 PRINT ARG=p2.X,p2.Y,p2.zpath STRIDE=1 FILE=colvar
-\endverbatim
+\endplumedfile
 
 The additional input file for this calculation, which contains the reference frames and the values of X and Y at these reference
 points has the following format.
diff --git a/src/multicolvar/AlphaBeta.cpp b/src/multicolvar/AlphaBeta.cpp
index 01379c583..aea79fc2e 100644
--- a/src/multicolvar/AlphaBeta.cpp
+++ b/src/multicolvar/AlphaBeta.cpp
@@ -49,7 +49,7 @@ The \f$\phi_i^{\textrm{Ref}}\f$ values are the user-specified reference values f
 
 The following provides an example of the input for an alpha beta similarity.
 
-\verbatim
+\plumedfile
 ALPHABETA ...
 ATOMS1=168,170,172,188 REFERENCE1=3.14
 ATOMS2=170,172,188,190 REFERENCE2=3.14
@@ -57,11 +57,11 @@ ATOMS3=188,190,192,230 REFERENCE3=3.14
 LABEL=ab
 ... ALPHABETA
 PRINT ARG=ab FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 Because all the reference values are the same we can calculate the same quantity using
 
-\verbatim
+\plumedfile
 ALPHABETA ...
 ATOMS1=168,170,172,188 REFERENCE=3.14
 ATOMS2=170,172,188,190
@@ -69,13 +69,13 @@ ATOMS3=188,190,192,230
 LABEL=ab
 ... ALPHABETA
 PRINT ARG=ab FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 Writing out the atoms involved in all the torsions in this way can be rather tedious. Thankfully if you are working with protein you
 can avoid this by using the \ref MOLINFO command.  PLUMED uses the pdb file that you provide to this command to learn
 about the topology of the protein molecule.  This means that you can specify torsion angles using the following syntax:
 
-\verbatim
+\plumedfile
 MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
 ALPHABETA ...
 ATOMS1=@phi-3 REFERENCE=3.14
@@ -84,7 +84,7 @@ ATOMS3=@phi-4
 LABEL=ab
 ... ALPHABETA
 PRINT ARG=ab FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 Here, \@phi-3 tells plumed that you would like to calculate the \f$\phi\f$ angle in the third residue of the protein.
 Similarly \@psi-4 tells plumed that you want to calculate the \f$\psi\f$ angle of the 4th residue of the protein.
diff --git a/src/multicolvar/Angles.cpp b/src/multicolvar/Angles.cpp
index a42233cbd..cdc062579 100644
--- a/src/multicolvar/Angles.cpp
+++ b/src/multicolvar/Angles.cpp
@@ -58,28 +58,28 @@ an atom / molecule \cite lj-recon.
 The following example instructs plumed to find the average of two angles and to
 print it to a file
 
-\verbatim
+\plumedfile
 ANGLES ATOMS1=1,2,3 ATOMS2=4,5,6 MEAN LABEL=a1
 PRINT ARG=a1.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following example tells plumed to calculate all angles involving
 at least one atom from GROUPA and two atoms from GROUPB in which the distances
 are less than 1.0. The number of angles between \f$\frac{\pi}{4}\f$ and
 \f$\frac{3\pi}{4}\f$ is then output
 
-\verbatim
+\plumedfile
 ANGLES GROUPA=1-10 GROUPB=11-100 BETWEEN={GAUSSIAN LOWER=0.25pi UPPER=0.75pi} SWITCH={GAUSSIAN R_0=1.0} LABEL=a1
 PRINT ARG=a1.between FILE=colvar
-\endverbatim
+\endplumedfile
 
 This final example instructs plumed to calculate all the angles in the first coordination
 spheres of the atoms. A discretized-normalized histogram of the distribution is then output
 
-\verbatim
+\plumedfile
 ANGLES GROUP=1-38 HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=pi NBINS=20} SWITCH={GAUSSIAN R_0=1.0} LABEL=a1
 PRINT ARG=a1.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/Bridge.cpp b/src/multicolvar/Bridge.cpp
index 42198ffa1..53e76158b 100644
--- a/src/multicolvar/Bridge.cpp
+++ b/src/multicolvar/Bridge.cpp
@@ -51,10 +51,10 @@ The following example instructs plumed to calculate the number of water molecule
 that are bridging betweeen atoms 1-10 and atoms 11-20 and to print the value
 to a file
 
-\verbatim
+\plumedfile
 BRIDGE BRIDGING_ATOMS=100-200 GROUPA=1-10 GROUPB=11-20 LABEL=w1
 PRINT ARG=a1.mean FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/CenterOfMultiColvar.cpp b/src/multicolvar/CenterOfMultiColvar.cpp
index f8b5be660..0708bd402 100644
--- a/src/multicolvar/CenterOfMultiColvar.cpp
+++ b/src/multicolvar/CenterOfMultiColvar.cpp
@@ -55,10 +55,10 @@ As you want to calculate the position of the droplets you thus recognise that th
 numbers should have a high weight in the weighted average you are using to calculate the position of the droplet.
 You can thus calculate the position of the droplet using an input like the one shown below:
 
-\verbatim
+\plumedfile
 c1: COORDINATIONNUMBER SPECIES=1-512 SWITCH={EXP D_0=4.0 R_0=0.5}
 cc: CENTER_OF_MULTICOLVAR DATA=c1
-\endverbatim
+\endplumedfile
 
 The first line here calclates the coordination numbers of all the atoms in the system.  The virtual atom then uses the values
 of the coordination numbers calculated by the action labelled c1 when it calculates the Berry Phase average described above.
@@ -66,11 +66,11 @@ of the coordination numbers calculated by the action labelled c1 when it calcula
 
 The above input is fine we can, however, refine this somewhat by making use of a multicolvar transform action as shown below:
 
-\verbatim
+\plumedfile
 c1: COORDINATIONNUMBER SPECIES=1-512 SWITCH={EXP D_0=4.0 R_0=0.5}
 cf: MTRANSFORM_MORE DATA=c1 SWITCH={RATIONAL D_0=2.0 R_0=0.1} LOWMEM
 cc: CENTER_OF_MULTICOLVAR DATA=cf
-\endverbatim
+\endplumedfile
 
 This input once again calculates the coordination numbers of all the atoms in the system.  The middle line then transforms these
 coordinations numbers to numbers between 0 and 1.  Essentially any atom with a coordination number larger than 2.0 is given a weight
diff --git a/src/multicolvar/CoordinationNumbers.cpp b/src/multicolvar/CoordinationNumbers.cpp
index a3823b488..38c26d348 100644
--- a/src/multicolvar/CoordinationNumbers.cpp
+++ b/src/multicolvar/CoordinationNumbers.cpp
@@ -48,16 +48,16 @@ s = \frac{ 1 - \left(\frac{r-d_0}{r_0}\right)^n } { 1 - \left(\frac{r-d_0}{r_0}\
 
 The following input tells plumed to calculate the coordination numbers of atoms 1-100 with themselves.
 The minimum coordination number is then calculated.
-\verbatim
+\plumedfile
 COORDINATIONNUMBER SPECIES=1-100 R_0=1.0 MIN={BETA=0.1}
-\endverbatim
+\endplumedfile
 
 The following input tells plumed to calculate how many atoms from 1-100 are within 3.0 of each of the atoms
 from 101-110.  In the first 101 is the central atom, in the second 102 is the central atom and so on.  The
 number of coordination numbers more than 6 is then computed.
-\verbatim
+\plumedfile
 COORDINATIONNUMBER SPECIESA=101-110 SPECIESB=1-100 R_0=3.0 MORE_THAN={RATIONAL R_0=6.0 NN=6 MM=12 D_0=0}
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/Density.cpp b/src/multicolvar/Density.cpp
index 0351140a3..598e16b30 100644
--- a/src/multicolvar/Density.cpp
+++ b/src/multicolvar/Density.cpp
@@ -40,11 +40,11 @@ the number of atoms in half the box.
 
 The following example calculates the number of atoms in one half of the simulation box.
 
-\verbatim
+\plumedfile
 DENSITY SPECIES=1-100 LABEL=d
 AROUND ARG=d XLOWER=0.0 XUPPER=0.5 LABEL=d1
 PRINT ARG=d1.* FILE=colvar1 FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/DihedralCorrelation.cpp b/src/multicolvar/DihedralCorrelation.cpp
index e75bde668..a1f5ea872 100644
--- a/src/multicolvar/DihedralCorrelation.cpp
+++ b/src/multicolvar/DihedralCorrelation.cpp
@@ -48,14 +48,14 @@ where the \f$\phi_i\f$ and \f$\psi\f$ values and the instantaneous values for th
 
 The following provides an example input for the dihcor action
 
-\verbatim
+\plumedfile
 DIHCOR ...
   ATOMS1=1,2,3,4,5,6,7,8
   ATOMS2=5,6,7,8,9,10,11,12
   LABEL=dih
 ... DIHCOR
 PRINT ARG=dih FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 In the above input we are calculating the correation between the torsion angle involving atoms 1, 2, 3 and 4 and the torsion angle
 involving atoms 5, 6, 7 and 8.	This is then added to the correlation betwene the torsion angle involving atoms 5, 6, 7 and 8 and the
@@ -65,7 +65,7 @@ Writing out the atoms involved in all the torsions in this way can be rather ted
 can avoid this by using the \ref MOLINFO command.  PLUMED uses the pdb file that you provide to this command to learn
 about the topology of the protein molecule.  This means that you can specify torsion angles using the following syntax:
 
-\verbatim
+\plumedfile
 MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
 DIHCOR ...
 ATOMS1=@phi-3,@psi-3
@@ -73,7 +73,7 @@ ATOMS2=@psi-3,@phi-4
 ATOMS4=@phi-4,@psi-4
 ... DIHCOR
 PRINT ARG=dih FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 Here, \@phi-3 tells plumed that you would like to calculate the \f$\phi\f$ angle in the third residue of the protein.
 Similarly \@psi-4 tells plumed that you want to calculate the \f$\psi\f$ angle of the 4th residue of the protein.
diff --git a/src/multicolvar/DistanceFromContour.cpp b/src/multicolvar/DistanceFromContour.cpp
index 214930392..6f99e214a 100644
--- a/src/multicolvar/DistanceFromContour.cpp
+++ b/src/multicolvar/DistanceFromContour.cpp
@@ -58,10 +58,10 @@ In this example atoms 2-100 are assumed to be concentraed along some part of the
 an interface between a liquid/solid and the vapour.  The quantity dc measures the distance between the
 surface at which the density of 2-100 atoms is equal to 0.2 and the position of the test particle atom 1.
 
-\verbatim
+\plumedfile
 dens: DENSITY SPECIES=2-100
 dc: DISTANCE_FROM_CONTOUR DATA=dens ATOM=1 BANDWIDTH=0.5,0.5,0.5 DIR=z CONTOUR=0.2
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/Distances.cpp b/src/multicolvar/Distances.cpp
index d5cb8fde0..579668b88 100644
--- a/src/multicolvar/Distances.cpp
+++ b/src/multicolvar/Distances.cpp
@@ -42,36 +42,36 @@ distances such as the minimum, the number less than a certain quantity and so on
 
 The following input tells plumed to calculate the distances between atoms 3 and 5 and between atoms 1 and 2 and to
 print the minimum for these two distances.
-\verbatim
+\plumedfile
 DISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 
 The following input tells plumed to calculate the distances between atoms 3 and 5 and between atoms 1 and 2
 and then to calculate the number of these distances that are less than 0.1 nm.  The number of distances
 less than 0.1nm is then printed to a file.
-\verbatim
+\plumedfile
 DISTANCES ATOMS1=3,5 ATOMS2=1,2 LABEL=d1 LESS_THAN={RATIONAL R_0=0.1}
 PRINT ARG=d1.lt0.1
-\endverbatim
+\endplumedfile
 (See also \ref PRINT \ref switchingfunction).
 
 The following input tells plumed to calculate all the distances between atoms 1, 2 and 3 (i.e. the distances between atoms
 1 and 2, atoms 1 and 3 and atoms 2 and 3).  The average of these distances is then calculated.
-\verbatim
+\plumedfile
 DISTANCES GROUP=1-3 MEAN LABEL=d1
 PRINT ARG=d1.mean
-\endverbatim
+\endplumedfile
 (See also \ref PRINT)
 
 The following input tells plumed to calculate all the distances between the atoms in GROUPA and the atoms in GROUPB.
 In other words the distances between atoms 1 and 2 and the distance between atoms 1 and 3.  The number of distances
 more than 0.1 is then printed to a file.
-\verbatim
+\plumedfile
 DISTANCES GROUPA=1 GROUPB=2,3 MORE_THAN={RATIONAL R_0=0.1}
 PRINT ARG=d1.gt0.1
-\endverbatim
+\endplumedfile
 (See also \ref PRINT \ref switchingfunction)
 
 
@@ -79,10 +79,10 @@ PRINT ARG=d1.gt0.1
 
 To calculate and print the minimum distance between two groups of atoms you use the following commands
 
-\verbatim
+\plumedfile
 d1: DISTANCES GROUPA=1-10 GROUPB=11-20 MIN={BETA=500.}
 PRINT ARG=d1.min FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 (see \ref DISTANCES and \ref PRINT)
 
 In order to ensure differentiability the minimum is calculated using the following function:
@@ -102,7 +102,7 @@ allow you to calculate multiple functions of a distribution of simple collective
 can calculate the number of distances less than 1.0, the minimum distance, the number of distances more than
 2.0 and the number of distances between 1.0 and 2.0 by using the following command:
 
-\verbatim
+\plumedfile
 DISTANCES ...
  GROUPA=1-10 GROUPB=11-20
  LESS_THAN={RATIONAL R_0=1.0}
@@ -111,7 +111,7 @@ DISTANCES ...
  MIN={BETA=500.}
 ... DISTANCES
 PRINT ARG=d1.lessthan,d1.morethan,d1.between,d1.min FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 (see \ref DISTANCES and \ref PRINT)
 
 A calculation performed this way is fast because the expensive part of the calculation - the calculation of all the distances - is only
diff --git a/src/multicolvar/DumpMultiColvar.cpp b/src/multicolvar/DumpMultiColvar.cpp
index 27029a436..fb08af5b4 100644
--- a/src/multicolvar/DumpMultiColvar.cpp
+++ b/src/multicolvar/DumpMultiColvar.cpp
@@ -49,13 +49,13 @@ In this examples we calculate the distances between the  atoms of the first and
 group and we write them in the file MULTICOLVAR.xyz. For each couple it writes the
 coordinates of their geometric center and their distance.
 
-\verbatim
+\plumedfile
 pos:   GROUP ATOMS=220,221,235,236,247,248,438,439,450,451,534,535
 neg:   GROUP ATOMS=65,68,138,182,185,267,270,291,313,316,489,583,621,711
 DISTANCES GROUPA=pos GROUPB=neg LABEL=slt
 
 DUMPMULTICOLVAR DATA=slt FILE=MULTICOLVAR.xyz
-\endverbatim
+\endplumedfile
 
 (see also \ref DISTANCES)
 
diff --git a/src/multicolvar/FilterBetween.cpp b/src/multicolvar/FilterBetween.cpp
index d6a418711..17ad0be9e 100644
--- a/src/multicolvar/FilterBetween.cpp
+++ b/src/multicolvar/FilterBetween.cpp
@@ -59,22 +59,22 @@ In other words, you are calculating the mean for the transformed colvar.
 The following input gives an example of how a MTRANSFORM_BETWEEN action can be used to duplicate
 functionality that is elsehwere in PLUMED.
 
-\verbatim
+\plumedfile
 DISTANCES ...
  GROUPA=1-10 GROUPB=11-20
  LABEL=d1
 ... DISTANCES
 MTRANSFORM_BETWEEN DATA=d1 LOWER=1.0 UPPER=2.0 SMEAR=0.5
-\endverbatim
+\endplumedfile
 
 In this case you can achieve the same result by using:
 
-\verbatim
+\plumedfile
 DISTANCES ...
  GROUPA=1-10 GROUPB=11-20
  BETWEEN={GAUSSIAN LOWER=1.0 UPPER=2.0}
 ... DISTANCES
-\endverbatim
+\endplumedfile
 (see \ref DISTANCES)
 
 The advantage of MTRANSFORM_BETWEEN comes, however, if you want to use transformed colvars as input
@@ -112,10 +112,10 @@ One is thus calculating the mean for those colvars that are within the range of
 
 The example shown below calculates the mean for those distances that are between 0 and 3 nm in length
 
-\verbatim
+\plumedfile
 DISTANCES GROUPA=1 GROUPB=2-50 MEAN LABEL=d1
 MFILTER_BETWEEN DATA=d1 LOWER=0 UPPER=3.0 SMEAR=0.0001 MEAN LABEL=d4
-\endverbatim
+\endplumedfile
 
 More complicated things can be done by using the label of a filter as input to a new multicolvar as shown
 in the example below.  Here the coordination numbers of all atoms are computed.  The atoms with a coordination
@@ -123,11 +123,11 @@ number between 4 and 6 are then identified using the filter.  This reduced list
 to a second coordination number calculation.  This second coordination number thus measures the number of atoms
 4-6 coordinated atoms each of the 4-6 coordination atoms is bound to.
 
-\verbatim
+\plumedfile
 c1: COORDINATIONNUMBER SPECIES=1-150 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0}
 cf: MFILTER_BETWEEN DATA=c1 LOWER=4 UPPER=6 SMEAR=0.5 LOWMEM
 c2: COORDINATIONNUMBER SPECIES=cf SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0} MORE_THAN={RATIONAL D_0=2.0 R_0=0.1}
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/FilterLessThan.cpp b/src/multicolvar/FilterLessThan.cpp
index 6363bcbb4..e65d5fcd1 100644
--- a/src/multicolvar/FilterLessThan.cpp
+++ b/src/multicolvar/FilterLessThan.cpp
@@ -51,22 +51,22 @@ In other words, you are calculating the mean for the transformed colvar.
 The following input gives an example of how a MTRANSFORM_LESS action can be used to duplicate
 functionality that is elsehwere in PLUMED.
 
-\verbatim
+\plumedfile
 DISTANCES ...
  GROUPA=1-10 GROUPB=11-20
  LABEL=d1
 ... DISTANCES
 MTRANSFORM_LESS DATA=d1 SWITCH={GAUSSIAN D_0=1.5 R_0=0.00001}
-\endverbatim
+\endplumedfile
 
 In this case you can achieve the same result by using:
 
-\verbatim
+\plumedfile
 DISTANCES ...
  GROUPA=1-10 GROUPB=11-20
  LESS_THAN={GAUSSIAN D_0=1.5 R_0=0.00001}
 ... DISTANCES
-\endverbatim
+\endplumedfile
 (see \ref DISTANCES)
 
 The advantage of MTRANSFORM_LESS comes, however, if you want to use transformed colvars as input
@@ -99,10 +99,10 @@ One is thus calculating the mean for those colvars that are less than the target
 
 The example shown below calculates the mean for those distances that less than 1.5 nm in length
 
-\verbatim
+\plumedfile
 DISTANCES GROUPA=1 GROUPB=2-50 MEAN LABEL=d1
 MFILTER_LESS DATA=d1 SWITCH={GAUSSIAN D_0=1.5 R_0=0.00001} MEAN LABEL=d4
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/FilterMoreThan.cpp b/src/multicolvar/FilterMoreThan.cpp
index 961c33d01..50c1dc59b 100644
--- a/src/multicolvar/FilterMoreThan.cpp
+++ b/src/multicolvar/FilterMoreThan.cpp
@@ -51,22 +51,22 @@ In other words, you are calculating the mean for the transformed colvar.
 The following input gives an example of how a MTRANSFORM_MORE action can be used to duplicate
 functionality that is elsehwere in PLUMED.
 
-\verbatim
+\plumedfile
 DISTANCES ...
  GROUPA=1-10 GROUPB=11-20
  LABEL=d1
 ... DISTANCES
 MTRANSFORM_MORE DATA=d1 SWITCH={GAUSSIAN D_0=1.5 R_0=0.00001}
-\endverbatim
+\endplumedfile
 
 In this case you can achieve the same result by using:
 
-\verbatim
+\plumedfile
 DISTANCES ...
  GROUPA=1-10 GROUPB=11-20
  MORE_THAN={GAUSSIAN D_0=1.5 R_0=0.00001}
 ... DISTANCES
-\endverbatim
+\endplumedfile
 (see \ref DISTANCES)
 
 The advantage of MTRANSFORM_MORE comes, however, if you want to use transformed colvars as input
@@ -104,10 +104,10 @@ One is thus calculating the mean for those colvars that are greater than the tar
 
 The example shown below calculates the mean for those distances that greater than 1.5 nm in length
 
-\verbatim
+\plumedfile
 DISTANCES GROUPA=1 GROUPB=2-50 MEAN LABEL=d1
 MFILTER_MORE DATA=d1 SWITCH={GAUSSIAN D_0=1.5 R_0=0.00001} MEAN LABEL=d4
-\endverbatim
+\endplumedfile
 
 More complicated things can be done by using the label of a filter as input to a new multicolvar as shown
 in the example below.  Here the coordination numbers of all atoms are computed.  The atoms with a coordination
@@ -115,11 +115,11 @@ number greater than 2 are then identified using the filter.  This reduced list o
 to a second coordination number calculation.  This second coordination number thus measures the number of
 two-coordinated atoms that each of the two-coordinated atoms is bound to.
 
-\verbatim
+\plumedfile
 1: COORDINATIONNUMBER SPECIES=1-150 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0}
 cf: MFILTER_MORE DATA=c1 SWITCH={RATIONAL D_0=2.0 R_0=0.1} LOWMEM
 c2: COORDINATIONNUMBER SPECIES=cf SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0} MORE_THAN={RATIONAL D_0=2.0 R_0=0.1}
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/InPlaneDistances.cpp b/src/multicolvar/InPlaneDistances.cpp
index d041fa1e0..ef9f4a55f 100644
--- a/src/multicolvar/InPlaneDistances.cpp
+++ b/src/multicolvar/InPlaneDistances.cpp
@@ -53,10 +53,10 @@ Keywords such as MORE_THAN and LESS_THAN can then be used to calculate the numbe
 The following input can be used to calculate the number of atoms that have indices greater than 3 and less than 101 that
 are within a cylinder with a radius of 0.3 nm that has its long axis aligned with the vector connecting atoms 1 and 2.
 
-\verbatim
+\plumedfile
 d1: INPLANEDISTANCES VECTORSTART=1 VECTOREND=2 GROUP=3-100 LESS_THAN={RATIONAL D_0=0.2 R_0=0.1}
 PRINT ARG=d1.lessthan FILE=colvar
-\endverbatim
+\endplumedfile
 
 
 */
diff --git a/src/multicolvar/LocalAverage.cpp b/src/multicolvar/LocalAverage.cpp
index fdc3e09c4..c8bae1d56 100644
--- a/src/multicolvar/LocalAverage.cpp
+++ b/src/multicolvar/LocalAverage.cpp
@@ -59,21 +59,21 @@ and so on.  You can also probe the value of these averaged variables in regions
 This example input calculates the coordination numbers for all the atoms in the system.  These coordination numbers are then averaged over
 spherical regions.  The number of averaged coordination numbers that are greater than 4 is then output to a file.
 
-\verbatim
+\plumedfile
 COORDINATIONNUMBER SPECIES=1-64 D_0=1.3 R_0=0.2 LABEL=d1
 LOCAL_AVERAGE ARG=d1 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MORE_THAN={RATIONAL R_0=4} LABEL=la
 PRINT ARG=la.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 This example input calculates the \f$q_4\f$ (see \ref Q4) vectors for each of the atoms in the system.  These vectors are then averaged
 component by component over a spherical region.  The average value for this quantity is then outputeed to a file.  This calculates the
 quantities that were used in the paper by Lechner and Dellago \cite dellago-q6
 
-\verbatim
+\plumedfile
 Q4 SPECIES=1-64 SWITCH={RATIONAL D_0=1.3 R_0=0.2} LABEL=q4
 LOCAL_AVERAGE ARG=q4 SWITCH={RATIONAL D_0=1.3 R_0=0.2} MEAN LABEL=la
 PRINT ARG=la.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/MultiColvarDensity.cpp b/src/multicolvar/MultiColvarDensity.cpp
index 7ea3dc009..8e63b6e51 100644
--- a/src/multicolvar/MultiColvarDensity.cpp
+++ b/src/multicolvar/MultiColvarDensity.cpp
@@ -58,11 +58,11 @@ The following example shows perhaps the simplest way in which this action can be
 input computes the density of atoms at each point on the grid and ouptuts this quantity to a file.  In
 other words this input instructs plumed to calculate \f$\rho(\mathbf{r}) = \sum_i K(\mathbf{r} - \mathbf{r}_i )\f$
 
-\verbatim
+\plumedfile
 dens: DENSITY SPECIES=1-100
 grid: MULTICOLVARDENS DATA=dens ORIGIN=1 DIR=xyz NBINS=100,100,100 BANDWIDTH=0.05,0.05,0.05 STRIDE=1
 DUMPGRID GRID=grid STRIDE=500 FILE=density
-\endverbatim
+\endplumedfile
 
 In the above example density is added to the grid on every step.  The PRINT_GRID instruction thus tells PLUMED to
 output the average density at each point on the grid every 500 steps of simulation.  Notice that the that grid output
@@ -72,11 +72,11 @@ of data separately you must use the CLEAR flag.
 This second example computes an order parameter (in this case \ref FCCUBIC) and constructs a phase field model
 for this order parameter using the equation above.
 
-\verbatim
+\plumedfile
 fcc: FCCUBIC SPECIES=1-5184 SWITCH={CUBIC D_0=1.2 D_MAX=1.5} ALPHA=27
 dens: MULTICOLVARDENS DATA=fcc ORIGIN=1 DIR=xyz NBINS=14,14,28 BANDWIDTH=1.0,1.0,1.0 STRIDE=1 CLEAR=1
 DUMPCUBE GRID=dens STRIDE=1 FILE=dens.cube
-\endverbatim
+\endplumedfile
 
 In this example the phase field model is computed and output to a file on every step of the simulation.  Furthermore,
 because the CLEAR=1 keyword is set on the MULTICOLVARDENS line each Gaussian cube file output is a phase field
diff --git a/src/multicolvar/NumberOfLinks.cpp b/src/multicolvar/NumberOfLinks.cpp
index b8257c5e4..c8850755b 100644
--- a/src/multicolvar/NumberOfLinks.cpp
+++ b/src/multicolvar/NumberOfLinks.cpp
@@ -43,22 +43,22 @@ similar orientations.  The vectors on individual atoms could be Steinhardt param
 The following calculates how many bonds there are in a system containing 64 atoms and outputs
 this quantity to a file.
 
-\verbatim
+\plumedfile
 DENSITY SPECIES=1-64 LABEL=d1
 NLINKS ARG=d1 SWITCH={RATIONAL D_0=1.3 R_0=0.2} LABEL=dd
 PRINT ARG=dd FILE=colvar
-\endverbatim
+\endplumedfile
 
 The following calculates how many pairs of neighbouring atoms in a system containg 64 atoms have
 similar dispositions for the atoms in their coordination sphere.  This calculation uses the
 dot product of the Q6 vectors on adjacent atoms to measure whether or not two atoms have the same
 ``orientation"
 
-\verbatim
+\plumedfile
 Q6 SPECIES=1-64 SWITCH={RATIONAL D_0=1.3 R_0=0.2} LABEL=q6
 NLINKS ARG=q6 SWITCH={RATIONAL D_0=1.3 R_0=0.2} LABEL=dd
 PRINT ARG=dd FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/Torsions.cpp b/src/multicolvar/Torsions.cpp
index 59d02cbfd..4e3fcfbb2 100644
--- a/src/multicolvar/Torsions.cpp
+++ b/src/multicolvar/Torsions.cpp
@@ -40,7 +40,7 @@ Calculate whether or not a set of torsional angles are within a particular range
 
 The following provides an example of the input for the torsions command
 
-\verbatim
+\plumedfile
 TORSIONS ...
 ATOMS1=168,170,172,188
 ATOMS2=170,172,188,190
@@ -48,13 +48,13 @@ ATOMS3=188,190,192,230
 LABEL=ab
 ... TORSIONS
 PRINT ARG=ab.* FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 Writing out the atoms involved in all the torsions in this way can be rather tedious. Thankfully if you are working with protein you
 can avoid this by using the \ref MOLINFO command.  PLUMED uses the pdb file that you provide to this command to learn
 about the topology of the protein molecule.  This means that you can specify torsion angles using the following syntax:
 
-\verbatim
+\plumedfile
 MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
 TORSIONS ...
 ATOMS1=@phi-3
@@ -63,7 +63,7 @@ ATOMS3=@phi-4
 LABEL=ab
 ... TORSIONS
 PRINT ARG=ab FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 Here, \@phi-3 tells plumed that you would like to calculate the \f$\phi\f$ angle in the third residue of the protein.
 Similarly \@psi-4 tells plumed that you want to calculate the \f$\psi\f$ angle of the 4th residue of the protein.
diff --git a/src/multicolvar/VolumeAround.cpp b/src/multicolvar/VolumeAround.cpp
index bab79d8e0..d44f36d4b 100644
--- a/src/multicolvar/VolumeAround.cpp
+++ b/src/multicolvar/VolumeAround.cpp
@@ -55,11 +55,11 @@ When AROUND is used with the \ref DENSITY action the number of atoms in the spec
 
 The following commands tell plumed to calculate the average coordination number for the atoms
 that have x (in fractional coordinates) within 2.0 nm of the com of mass c1. The final value will be labeled s.mean.
-\verbatim
+\plumedfile
 COM ATOMS=1-100 LABEL=c1
 COORDINATIONNUMBER SPECIES=1-100 R_0=1.0 LABEL=c
 AROUND DATA=c ORIGIN=c1 XLOWER=-2.0 XUPPER=2.0 SIGMA=0.1 MEAN LABEL=s
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/VolumeCavity.cpp b/src/multicolvar/VolumeCavity.cpp
index af14a2be0..70d8254f3 100644
--- a/src/multicolvar/VolumeCavity.cpp
+++ b/src/multicolvar/VolumeCavity.cpp
@@ -86,19 +86,19 @@ described above and the resulting projections determine the \f$u'\f$, \f$v'\f$ a
 The following commands tell plumed to calculate the number of atoms in an ion chanel in a protein.
 The extent of the chanel is calculated from the positions of atoms 1, 4, 5 and 11. The final value will be labeled cav.
 
-\verbatim
+\plumedfile
 d1: DENSITY SPECIES=20-500
 CAVITY DATA=d1 ATOMS=1,4,5,11 SIGMA=0.1 LABEL=cav
-\endverbatim
+\endplumedfile
 
 The following command tells plumed to calculate the coordination numbers (with other water molecules) for the water
 molecules in the protein channel described above.  The average coordination number and the number of coordination
 numbers more than 4 is then calculated.  The values of these two quantities are given the labels cav.mean and cav.morethan
 
-\verbatim
+\plumedfile
 d1: COORDINATIONNUMBER SPECIES=20-500
 CAVITY DATA=d1 ATOMS=1,4,5,11 SIGMA=0.1 MEAN MORE_THAN={RATIONAL R_0=4} LABEL=cav
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/VolumeInCylinder.cpp b/src/multicolvar/VolumeInCylinder.cpp
index eb6bc611e..3e8acfb22 100644
--- a/src/multicolvar/VolumeInCylinder.cpp
+++ b/src/multicolvar/VolumeInCylinder.cpp
@@ -59,11 +59,11 @@ When INCYLINDER is used with the \ref DENSITY action the number of atoms in the
 The input below can be use to calculate the average coordination numbers for those atoms that are within a cylindrical tube
 of radius 1.5 nm that is centered on the position of atom 101 and that has its long axis parallel to the z-axis.
 
-\verbatim
+\plumedfile
 c1: COORDINATIONNUMBER SPECIES=1-100 SWITCH={RATIONAL R_0=0.1}
 d2: INCYLINDER ATOM=101 DATA=d1 DIRECTION=Z RADIUS={TANH R_0=1.5} SIGMA=0.1 LOWER=-0.1 UPPER=0.1 MEAN
 PRINT ARG=d2.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/VolumeInSphere.cpp b/src/multicolvar/VolumeInSphere.cpp
index 632af88df..937e067dd 100644
--- a/src/multicolvar/VolumeInSphere.cpp
+++ b/src/multicolvar/VolumeInSphere.cpp
@@ -58,11 +58,11 @@ When INCYLINDER is used with the \ref DENSITY action the number of atoms in the
 The input below can be use to calculate the average coordination numbers for those atoms that are within a sphere
 of radius 1.5 nm that is centered on the position of atom 101.
 
-\verbatim
+\plumedfile
 c1: COORDINATIONNUMBER SPECIES=1-100 SWITCH={RATIONAL R_0=0.1}
 d2: INSPHERE ATOM=101 DATA=d1 RADIUS={TANH R_0=1.5} SIGMA=0.1 LOWER=-0.1 UPPER=0.1 MEAN
 PRINT ARG=d2.* FILE=colvar
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/VolumeTetrapore.cpp b/src/multicolvar/VolumeTetrapore.cpp
index 98fe915c3..c7dcd66d0 100644
--- a/src/multicolvar/VolumeTetrapore.cpp
+++ b/src/multicolvar/VolumeTetrapore.cpp
@@ -94,19 +94,19 @@ This is in fact the only point of difference between these two actions.
 The following commands tell plumed to calculate the number of atom inside a tetrahedral cavity.  The extent of the tetrahedral
 cavity is calculated from the positions of atoms 1, 4, 5, and 11,  The final value will be labeled cav.
 
-\verbatim
+\plumedfile
 d1: DENSITY SPECIES=20-500
 TETRAHEDRALPORE DATA=d1 ATOMS=1,4,5,11 SIGMA=0.1 LABEL=cav
-\endverbatim
+\endplumedfile
 
 The following command tells plumed to calculate the coordination numbers (with other water molecules) for the water
 molecules in the tetrahedral cavity described above.  The average coordination number and the number of coordination
 numbers more than 4 is then calculated.  The values of these two quantities are given the labels cav.mean and cav.morethan
 
-\verbatim
+\plumedfile
 d1: COORDINATIONNUMBER SPECIES=20-500
 CAVITY DATA=d1 ATOMS=1,4,5,11 SIGMA=0.1 MEAN MORE_THAN={RATIONAL R_0=4} LABEL=cav
-\endverbatim
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/XAngle.cpp b/src/multicolvar/XAngle.cpp
index fbc2ab1d7..dd90c187f 100644
--- a/src/multicolvar/XAngle.cpp
+++ b/src/multicolvar/XAngle.cpp
@@ -41,10 +41,10 @@ Calculate the angles between the vector connecting two atoms and the x axis.
 
 The following input tells plumed to calculate the angles between the x-axis and the vector connecting atom 3 to atom 5 and between the x-axis
 and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then
-\verbatim
+\plumedfile
 XANLGES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
@@ -57,10 +57,10 @@ Calculate the angles between the vector connecting two atoms and the y axis.
 
 The following input tells plumed to calculate the angles between the y-axis and the vector connecting atom 3 to atom 5 and between the y-axis
 and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then
-\verbatim
+\plumedfile
 YANLGES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
@@ -73,10 +73,10 @@ Calculate the angles between the vector connecting two atoms and the z axis.
 
 The following input tells plumed to calculate the angles between the z-axis and the vector connecting atom 3 to atom 5 and between the z-axis
 and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then
-\verbatim
+\plumedfile
 ZANLGES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
diff --git a/src/multicolvar/XDistances.cpp b/src/multicolvar/XDistances.cpp
index 009b04de5..fe29a8b65 100644
--- a/src/multicolvar/XDistances.cpp
+++ b/src/multicolvar/XDistances.cpp
@@ -42,38 +42,38 @@ values such as the minimum, the number less than a certain quantity and so on.
 The following input tells plumed to calculate the x-component of the vector connecting atom 3 to atom 5 and
 the x-component of the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then
 printed
-\verbatim
+\plumedfile
 XDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 
 
 The following input tells plumed to calculate the x-component of the vector connecting atom 3 to atom 5 and
 the x-component of the vector connecting atom 1 to atom 2.  The number of values that are
 less than 0.1nm is then printed to a file.
-\verbatim
+\plumedfile
 XDISTANCES ATOMS1=3,5 ATOMS2=1,2 LABEL=d1 LESS_THAN={RATIONAL R_0=0.1}
 PRINT ARG=d1.lt0.1
-\endverbatim
+\endplumedfile
 (See also \ref PRINT \ref switchingfunction).
 
 The following input tells plumed to calculate the x-components of all the distinct vectors that can be created
 between atoms 1, 2 and 3 (i.e. the vectors between atoms 1 and 2, atoms 1 and 3 and atoms 2 and 3).
 The average of these quantities is then calculated.
-\verbatim
+\plumedfile
 XDISTANCES GROUP=1-3 AVERAGE LABEL=d1
 PRINT ARG=d1.average
-\endverbatim
+\endplumedfile
 (See also \ref PRINT)
 
 The following input tells plumed to calculate all the vectors connecting the the atoms in GROUPA to the atoms in GROUPB.
 In other words the vector between atoms 1 and 2 and the vector between atoms 1 and 3.  The number of values
 more than 0.1 is then printed to a file.
-\verbatim
+\plumedfile
 XDISTANCES GROUPA=1 GROUPB=2,3 MORE_THAN={RATIONAL R_0=0.1}
 PRINT ARG=d1.gt0.1
-\endverbatim
+\endplumedfile
 (See also \ref PRINT \ref switchingfunction)
 */
 //+ENDPLUMEDOC
@@ -89,38 +89,38 @@ values such as the minimum, the number less than a certain quantity and so on.
 The following input tells plumed to calculate the y-component of the vector connecting atom 3 to atom 5 and
 the y-component of the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then
 printed
-\verbatim
+\plumedfile
 YDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 
 
 The following input tells plumed to calculate the y-component of the vector connecting atom 3 to atom 5 and
 the y-component of the vector connecting atom 1 to atom 2.  The number of values that are
 less than 0.1nm is then printed to a file.
-\verbatim
+\plumedfile
 YDISTANCES ATOMS1=3,5 ATOMS2=1,2 LABEL=d1 LESS_THAN={RATIONAL R_0=0.1}
 PRINT ARG=d1.lt0.1
-\endverbatim
+\endplumedfile
 (See also \ref PRINT \ref switchingfunction).
 
 The following input tells plumed to calculate the y-components of all the distinct vectors that can be created
 between atoms 1, 2 and 3 (i.e. the vectors between atoms 1 and 2, atoms 1 and 3 and atoms 2 and 3).
 The average of these quantities is then calculated.
-\verbatim
+\plumedfile
 YDISTANCES GROUP=1-3 AVERAGE LABEL=d1
 PRINT ARG=d1.average
-\endverbatim
+\endplumedfile
 (See also \ref PRINT)
 
 The following input tells plumed to calculate all the vectors connecting the the atoms in GROUPA to the atoms in GROUPB.
 In other words the vector between atoms 1 and 2 and the vector between atoms 1 and 3.  The number of values
 more than 0.1 is then printed to a file.
-\verbatim
+\plumedfile
 YDISTANCES GROUPA=1 GROUPB=2,3 MORE_THAN={RATIONAL R_0=0.1}
 PRINT ARG=d1.gt0.1
-\endverbatim
+\endplumedfile
 (See also \ref PRINT \ref switchingfunction)
 
 */
@@ -137,38 +137,38 @@ values such as the minimum, the number less than a certain quantity and so on.
 The following input tells plumed to calculate the z-component of the vector connecting atom 3 to atom 5 and
 the z-component of the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then
 printed
-\verbatim
+\plumedfile
 ZDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 
 
 The following input tells plumed to calculate the z-component of the vector connecting atom 3 to atom 5 and
 the z-component of the vector connecting atom 1 to atom 2.  The number of values that are
 less than 0.1nm is then printed to a file.
-\verbatim
+\plumedfile
 ZDISTANCES ATOMS1=3,5 ATOMS2=1,2 LABEL=d1 LESS_THAN={RATIONAL R_0=0.1}
 PRINT ARG=d1.lt0.1
-\endverbatim
+\endplumedfile
 (See also \ref PRINT \ref switchingfunction).
 
 The following input tells plumed to calculate the z-components of all the distinct vectors that can be created
 between atoms 1, 2 and 3 (i.e. the vectors between atoms 1 and 2, atoms 1 and 3 and atoms 2 and 3).
 The average of these quantities is then calculated.
-\verbatim
+\plumedfile
 ZDISTANCES GROUP=1-3 AVERAGE LABEL=d1
 PRINT ARG=d1.average
-\endverbatim
+\endplumedfile
 (See also \ref PRINT)
 
 The following input tells plumed to calculate all the vectors connecting the the atoms in GROUPA to the atoms in GROUPB.
 In other words the vector between atoms 1 and 2 and the vector between atoms 1 and 3.  The number of values
 more than 0.1 is then printed to a file.
-\verbatim
+\plumedfile
 ZDISTANCES GROUPA=1 GROUPB=2,3 MORE_THAN={RATIONAL R_0=0.1}
 PRINT ARG=d1.gt0.1
-\endverbatim
+\endplumedfile
 (See also \ref PRINT \ref switchingfunction)
 
 */
diff --git a/src/multicolvar/XYDistances.cpp b/src/multicolvar/XYDistances.cpp
index 767c64a86..c0c8d4aa5 100644
--- a/src/multicolvar/XYDistances.cpp
+++ b/src/multicolvar/XYDistances.cpp
@@ -43,10 +43,10 @@ The following input tells plumed to calculate the projection of the length of th
 to atom 5 projected in the xy-plane and the projection of the length of the vector
 the vector connecting atom 1 to atom 2 in the xy-plane.  The minimum of these two quantities is then
 printed
-\verbatim
+\plumedfile
 XYDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 
 */
@@ -64,10 +64,10 @@ The following input tells plumed to calculate the projection of the length of th
 to atom 5 projected in the xz-plane and the projection of the length of the vector
 the vector connecting atom 1 to atom 2 in the xz-plane.  The minimum of these two quantities is then
 printed
-\verbatim
+\plumedfile
 XZDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 
 */
@@ -85,10 +85,10 @@ The following input tells plumed to calculate the projection of the length of th
 to atom 5 in the yz-plane and the projection of the length of the vector
 the vector connecting atom 1 to atom 2 in the yz-plane.  The minimum of these two quantities is then
 printed
-\verbatim
+\plumedfile
 YZDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 
 */
diff --git a/src/multicolvar/XYTorsion.cpp b/src/multicolvar/XYTorsion.cpp
index fd512799a..4e0804a50 100644
--- a/src/multicolvar/XYTorsion.cpp
+++ b/src/multicolvar/XYTorsion.cpp
@@ -41,10 +41,10 @@ Calculate the torsional angle around the x axis from the positive y direction.
 
 The following input tells plumed to calculate the angle around the x direction between the positive y-axis and the vector connecting atom 3 to atom 5 and
 the angle around the x direction between the positive y axis and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then output
-\verbatim
+\plumedfile
 XYTORSIONS ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
@@ -57,10 +57,10 @@ Calculate the torsional angle around the x axis from the positive z direction.
 
 The following input tells plumed to calculate the angle around the x direction between the positive z-axis and the vector connecting atom 3 to atom 5 and
 the angle around the x direction between the positive z direction and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then output
-\verbatim
+\plumedfile
 XZTORSIONS ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
@@ -73,10 +73,10 @@ Calculate the torsional angle around the y axis from the positive x direction.
 
 The following input tells plumed to calculate the angle around the y direction between the positive x-direction and the vector connecting atom 3 to atom 5 and
 the angle around the y direction between the positive x axis and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then output
-\verbatim
+\plumedfile
 YXTORSIONS ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
@@ -89,10 +89,10 @@ Calculate the torsional angle around the y axis from the positive z direction.
 
 The following input tells plumed to calculate the angle around the y direction between the positive z-direction and the vector connecting atom 3 to atom 5 and
 the angle around the y direction between the positive z direction and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then output
-\verbatim
+\plumedfile
 YZTORSIONS ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
@@ -105,10 +105,10 @@ Calculate the torsional angle around the z axis from the positive x direction.
 
 The following input tells plumed to calculate the angle around the z direction between the positive x-direction and the vector connecting atom 3 to atom 5 and
 the angle around the z direction between the positive x-direction and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then output
-\verbatim
+\plumedfile
 ZXTORSIONS ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
@@ -121,10 +121,10 @@ Calculate the torsional angle around the z axis from the positive y direction.
 
 The following input tells plumed to calculate the angle around the z direction between the positive y-axis and the vector connecting atom 3 to atom 5 and
 the angle around the z direction between the positive y axis and the vector connecting atom 1 to atom 2.  The minimum of these two quantities is then output
-\verbatim
+\plumedfile
 ZYTORSIONS ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1} LABEL=d1
 PRINT ARG=d1.min
-\endverbatim
+\endplumedfile
 (See also \ref PRINT).
 */
 //+ENDPLUMEDOC
diff --git a/src/secondarystructure/AlphaRMSD.cpp b/src/secondarystructure/AlphaRMSD.cpp
index 6c9b49083..ecf10eacb 100644
--- a/src/secondarystructure/AlphaRMSD.cpp
+++ b/src/secondarystructure/AlphaRMSD.cpp
@@ -68,11 +68,10 @@ anthing other than TYPE=DRMSD.  For more details as to how to do this see \ref W
 The following input calculates the number of six residue segments of
 protein that are in an alpha helical configuration.
 
-\verbatim
+\plumedfile
 MOLINFO STRUCTURE=helix.pdb
 ALPHARMSD RESIDUES=all TYPE=DRMSD LESS_THAN={RATIONAL R_0=0.08 NN=8 MM=12} LABEL=a
-\endverbatim
-(see also \ref MOLINFO)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/secondarystructure/AntibetaRMSD.cpp b/src/secondarystructure/AntibetaRMSD.cpp
index 6e60d5445..637e4c0b3 100644
--- a/src/secondarystructure/AntibetaRMSD.cpp
+++ b/src/secondarystructure/AntibetaRMSD.cpp
@@ -70,11 +70,10 @@ anthing other than TYPE=DRMSD.  For more details as to how to do this see \ref W
 The following input calculates the number of six residue segments of
 protein that are in an antiparallel beta sheet configuration.
 
-\verbatim
+\plumedfile
 MOLINFO STRUCTURE=helix.pdb
 ANTIBETARMSD RESIDUES=all TYPE=DRMSD LESS_THAN={RATIONAL R_0=0.08 NN=8 MM=12} LABEL=a
-\endverbatim
-(see also \ref MOLINFO)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/secondarystructure/ParabetaRMSD.cpp b/src/secondarystructure/ParabetaRMSD.cpp
index fe8740864..1a67ec8af 100644
--- a/src/secondarystructure/ParabetaRMSD.cpp
+++ b/src/secondarystructure/ParabetaRMSD.cpp
@@ -70,11 +70,10 @@ anthing other than TYPE=DRMSD.  For more details as to how to do this see \ref W
 The following input calculates the number of six residue segments of
 protein that are in an parallel beta sheet configuration.
 
-\verbatim
+\plumedfile
 MOLINFO STRUCTURE=helix.pdb
 PARABETARMSD RESIDUES=all TYPE=DRMSD LESS_THAN={RATIONAL R_0=0.08 NN=8 MM=12} LABEL=a
-\endverbatim
-(see also \ref MOLINFO)
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/setup/Load.cpp b/src/setup/Load.cpp
index 93f9d39b8..4002ed0f0 100644
--- a/src/setup/Load.cpp
+++ b/src/setup/Load.cpp
@@ -43,9 +43,9 @@ If you have a shared object named extensions.so and want to
 use the functionalities implemented in it within PLUMED you can
 load it with the following syntax
 
-\verbatim
+\plumedfile
 LOAD FILE=extensions.so
-\endverbatim
+\endplumedfile
 
 As a more practical example, imagine that you want to make a
 small change to one collective variable that is already implemented
@@ -61,7 +61,7 @@ with different names. Then you can compile it into a shared object using
 This will generate a file `Distance2.so` (or `Distance2.dylib` on a mac)
 that can be loaded.
 Now you can use your new implementation with the following input
-\verbatim
+\plumedfile
 # load the new library
 LOAD FILE=Distance2.so
 # compute standard distance
@@ -70,11 +70,11 @@ d: DISTANCE ATOMS=1,10
 d2: DISTANCE2 ATOMS=1,10
 # print them on a file
 PRINT ARG=d,d2 FILE=compare-them
-\endverbatim
+\endplumedfile
 
 You can even skip the initial step and directly feed PLUMED
 with the `Distance2.cpp` file: it will be compiled on the fly.
-\verbatim
+\plumedfile
 # load the new definition
 # this is a cpp file so it will be compiled
 LOAD FILE=Distance2.cpp
@@ -84,7 +84,7 @@ d: DISTANCE ATOMS=1,10
 d2: DISTANCE2 ATOMS=1,10
 # print them on a file
 PRINT ARG=d,d2 FILE=compare-them
-\endverbatim
+\endplumedfile
 
 This will allow to make quick tests while developing your own
 variables. Of course, after your implementation is ready you might
diff --git a/src/setup/MolInfo.cpp b/src/setup/MolInfo.cpp
index 03a0012f9..24f2cb689 100644
--- a/src/setup/MolInfo.cpp
+++ b/src/setup/MolInfo.cpp
@@ -111,23 +111,21 @@ interpret terminal residue 1.
 In the following example the MOLINFO command is used to provide the information on which atoms
 are in the backbone of a protein to the ALPHARMSD CV.
 
-\verbatim
+\plumedfile
 MOLINFO STRUCTURE=reference.pdb
 ALPHARMSD RESIDUES=all TYPE=DRMSD LESS_THAN={RATIONAL R_0=0.08 NN=8 MM=12} LABEL=a
-\endverbatim
-(see also \ref ALPHARMSD)
+\endplumedfile
 
 The following example prints the distance corresponding to the hydrogen bonds
 in a GC Watson-Crick pair.
 
-\verbatim
+\plumedfile
 MOLINFO STRUCTURE=reference.pdb
 hb1: DISTANCE ATOMS=@N2-1,@O2-14
 hb2: DISTANCE ATOMS=@N1-1,@N3-14
 hb3: DISTANCE ATOMS=@O6-1,@N4-14
 PRINT ARG=hb1,hb2,hb3
-\endverbatim
-(see also \ref DISTANCE).
+\endplumedfile
 
 
 */
diff --git a/src/setup/Restart.cpp b/src/setup/Restart.cpp
index ffa4f46b1..4cfb03ef0 100644
--- a/src/setup/Restart.cpp
+++ b/src/setup/Restart.cpp
@@ -51,31 +51,28 @@ and \ref PBMETAD and on some analysis action.
 \par Examples
 
 Using the following input:
-\verbatim
+\plumedfile
 d: DISTANCE ATOMS=1,2
 PRINT ARG=d FILE=out
-\endverbatim
-(See also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 a new 'out' file will be created. If an old one is on the way, it will be automatically backed up.
 
 On the other hand, using the following input:
-\verbatim
+\plumedfile
 RESTART
 d: DISTANCE ATOMS=1,2
 PRINT ARG=d FILE=out
-\endverbatim
+\endplumedfile
 the file 'out' will be appended.
-(See also \ref DISTANCE and \ref PRINT).
 
 In the following case, file out1 will be backed up and file out2 will be concatenated
-\verbatim
+\plumedfile
 RESTART
 d1: DISTANCE ATOMS=1,2
 d2: DISTANCE ATOMS=1,2
 PRINT ARG=d1 FILE=out1 RESTART=NO
 PRINT ARG=d2 FILE=out2
-\endverbatim
-(See also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 
 
 
diff --git a/src/setup/Units.cpp b/src/setup/Units.cpp
index 1b4225d24..617e8c0a6 100644
--- a/src/setup/Units.cpp
+++ b/src/setup/Units.cpp
@@ -44,16 +44,15 @@ the units. For example, trajectories written in .gro format (with \ref DUMPATOMS
 are going to be always in nm.
 
 \par Examples
-\verbatim
+\plumedfile
 # this is using nm - kj/mol - fs
 UNITS LENGTH=nm TIME=fs
-\endverbatim
+\endplumedfile
 If a number, x, is found, the new unit is equal to x (default units)
-\verbatim
+\plumedfile
 # this is using nm - kj/mol - fs
 UNITS LENGTH=nm TIME=0.001
-\endverbatim
-
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/vatom/COM.cpp b/src/vatom/COM.cpp
index 8eb60b64d..d0cfd0b1f 100644
--- a/src/vatom/COM.cpp
+++ b/src/vatom/COM.cpp
@@ -54,13 +54,12 @@ periodic image.
 
 The following input instructs plumed to print the distance between the
 center of mass for atoms 1,2,3,4,5,6,7 and that for atoms 15,20:
-\verbatim
+\plumedfile
 c1: COM ATOMS=1-7
 c2: COM ATOMS=15,20
 d1: DISTANCE ATOMS=c1,c2
 PRINT ARG=d1
-\endverbatim
-(See also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/vatom/Center.cpp b/src/vatom/Center.cpp
index 7e4cdd440..64eea7d74 100644
--- a/src/vatom/Center.cpp
+++ b/src/vatom/Center.cpp
@@ -54,7 +54,7 @@ periodic image.
 
 \par Examples
 
-\verbatim
+\plumedfile
 # a point which is on the line connecting atoms 1 and 10, so that its distance
 # from 10 is twice its distance from 1:
 c1: CENTER ATOMS=1,1,10
@@ -67,8 +67,7 @@ c2: CENTER ATOMS=2,3,4,5 MASS
 d1: DISTANCE ATOMS=c1,c2
 
 PRINT ARG=d1
-\endverbatim
-(See also \ref DISTANCE, \ref COM and \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/src/vatom/FixedAtom.cpp b/src/vatom/FixedAtom.cpp
index ccf7eb52c..3baf40379 100644
--- a/src/vatom/FixedAtom.cpp
+++ b/src/vatom/FixedAtom.cpp
@@ -51,24 +51,22 @@ then it is safe to add further fixed atoms without breaking translational invari
 
 The following input instructs plumed to compute the angle between
 distance of atoms 15 and 20 and the z axis and keeping it close to zero.
-\verbatim
+\plumedfile
 a: FIXEDATOM AT=0,0,0
 b: FIXEDATOM AT=0,0,1
 an: ANGLE ATOMS=a,b,15,20
 RESTRAINT ARG=an AT=0.0 KAPPA=100.0
-\endverbatim
-(See also \ref ANGLE and \ref RESTRAINT).
+\endplumedfile
 
 The following input instructs plumed to align a protein on a template
 and then compute the distance of one of its atom from the point
 (10,20,30).
-\verbatim
+\plumedfile
 FIT_TO_TEMPLATE STRIDE=1 REFERENCE=ref.pdb TYPE=SIMPLE
 a: FIXEDATOM AT=10,20,30
 d: DISTANCE ARG=a,20
 PRINT ARG=d FILE=colvar
-\endverbatim
-(See also \ref FIT_TO_TEMPLATE and \ref DISTANCE).
+\endplumedfile
 
 
 */
diff --git a/src/vatom/Ghost.cpp b/src/vatom/Ghost.cpp
index b2c783717..7c57cee4f 100644
--- a/src/vatom/Ghost.cpp
+++ b/src/vatom/Ghost.cpp
@@ -39,13 +39,12 @@ an atom list through the the label for the GHOST action that creates it.
 \par Examples
 The following input instructs plumed to print the distance between the
 ghost atom and the center of mass for atoms 15,20:
-\verbatim
+\plumedfile
 c1: GHOST ATOMS=1,5,10 COORDINATES=10.0,10.0,10.0
 c2: COM ATOMS=15,20
 d1: DISTANCE ATOMS=c1,c2
 PRINT ARG=d1
-\endverbatim
-(See also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 
 */
 //+ENDPLUMEDOC
diff --git a/user-doc/Analysis.txt b/user-doc/Analysis.txt
index a5a0e5572..ef3a00589 100644
--- a/user-doc/Analysis.txt
+++ b/user-doc/Analysis.txt
@@ -45,11 +45,11 @@ As an example the following set of commands instructs PLUMED to calculate the di
 atoms 1 and 2 for every 5th frame in the trajectory and to accumulate a histogram from this data
 which will be output every 100 steps (i.e. when 20 distances have been added to the histogram).
 
-\verbatim
+\plumedfile
 x: DISTANCE ATOMS=1,2
 h: HISTOGRAM ARG=x GRID_MIN=0.0 GRID_MAX=3.0 GRID_BIN=100 BANDWIDTH=0.1 STRIDE=5
 DUMPGRID GRID=h FILE=histo STRIDE=100 
-\endverbatim
+\endplumedfile
 
 It is important to note when using commands such as the above the first frame in the trajectory is assumed 
 to be the initial configuration that was input to the MD code. It is thus ignored.  Furthermore, if you are 
@@ -57,11 +57,11 @@ running with driver and you would like to analyse the whole trajectory (without
 and then print the result you simply call \ref DUMPGRID (or any of the commands above) without a STRIDE 
 keyword as shown in the example below. 
 
-\verbatim
+\plumedfile
 x: DISTANCE ATOMS=1,2
 h: HISTOGRAM ARG=x GRID_MIN=0.0 GRID_MAX=3.0 GRID_BIN=100 BANDWIDTH=0.1 STRIDE=5
 DUMPGRID GRID=h FILE=histo 
-\endverbatim
+\endplumedfile
 
 Please note that even with this calculation the first frame in the trajectory is ignored when computing the 
 histogram.
diff --git a/user-doc/Colvar.txt b/user-doc/Colvar.txt
index 7c6557796..e3e20d2a1 100644
--- a/user-doc/Colvar.txt
+++ b/user-doc/Colvar.txt
@@ -92,9 +92,9 @@ Oftentimes the simplest way to specify the atoms involved is to use multiple ins
 i.e. ATOMS1, ATOMS2, ATOMS3,...  Separate instances of the quantity specified by NAME are then calculated for 
 each of the sets of atoms.  For example if the command issued contains the following:
 
-\verbatim
+\plumedfile
 DISTANCES ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6
-\endverbatim
+\endplumedfile
 
 The distances between atoms 1 and 2, atoms 3 and 4, and atoms 5 and 6 are calculated. Obviously, generating 
 this sort of input is rather tedious so short cuts are also available many of the collective variables. 
@@ -163,11 +163,11 @@ prevent the cluster subliming.  Alternatively, you may wish to insist that a par
 coordination number greater than 2.  You can add these sorts of restraints by employing the following biases, which all act 
 on the set of collective variable values calculated by a multicolvar.  So for example the following set of commands:
 
-\verbatim
+\plumedfile
 COM ATOMS=1-20 LABEL=c1
 DISTANCES GROUPA=c1 GROUPB=1-20 LABEL=d1
 UWALLS DATA=d1 AT=2.5 KAPPA=0.2 LABEL=sr
-\endverbatim
+\endplumedfile
 
 creates the aforementioned set of restraints on the distances between the 20 atoms in a cluster and the center of mass of the cluster.
 
diff --git a/user-doc/Files.txt b/user-doc/Files.txt
index f360af5d2..2ed8f21ee 100644
--- a/user-doc/Files.txt
+++ b/user-doc/Files.txt
@@ -34,11 +34,10 @@ since your disk might fill up quickly with this setting.
 When running with multiple replicas (e.g., with GROMACS, -multi option) PLUMED adds the replica index as a suffix to
 all the files. The following command
 will thus print files named COLVAR.0, COLVAR.1, etc for the different replicas.
-\verbatim
+\plumedfile
 d: DISTANCE ATOMS=1,2
 PRINT ARG=d FILE=COLVAR
-\endverbatim
-(see also \ref DISTANCE and \ref PRINT).
+\endplumedfile
 
 When reading a file, PLUMED will try to add the suffix. If the file is not found, it will fall back to
 the name without suffix. The most important case is the reading of the plumed input file.
@@ -51,12 +50,11 @@ extension. Before PLUMED 2.2, the only recognized suffix was ".gz". Since 2.2, a
 less or equal to five letters is recognized.
 
 This means that using in a multireplica context an input such as
-\verbatim
+\plumedfile
 d: DISTANCE ATOMS=1,2
 PRINT ARG=d FILE=COLVAR.gz
 METAD ARG=d FILE=test.HILLS SIGMA=0.1 HEIGHT=0.1
-\endverbatim
-(see \ref DISTANCE, \ref PRINT, and \ref METAD)
+\endplumedfile
 PLUMED will write files named COLVAR.0.gz, COLVAR.1.gz, test.0.HILLS, test.1.HILLS, etc
 etc. This is useful since the preserved extension makes it easy
 to process the files later.
diff --git a/user-doc/Functions.txt b/user-doc/Functions.txt
index 8bf37c42d..53de40f36 100644
--- a/user-doc/Functions.txt
+++ b/user-doc/Functions.txt
@@ -23,7 +23,7 @@ is the smallest value and `B` is the largest value. In case the answer is that t
 values at the discontinuity are not always the same, then you cannot construct a variable that
 can be biased with PLUMED. Consider the following examples:
 
-\verbatim
+\plumedfile
 t: TORSION ATOMS=1,2,3,4
 # When atoms are moved, t could jump suddenly from -pi to +pi
 
@@ -39,7 +39,7 @@ d: DISTANCE ATOMS=1,10 COMPONENTS
 # make a new variable equal to d.z but with the correct periodicity
 dz: COMBINE ARG=d.z PERIODIC=-10,10
 # here we assumed the system is in a orthorhombic box with z side = 20
-\endverbatim
+\endplumedfile
 
 @FUNCTION@
 
diff --git a/user-doc/Group.txt b/user-doc/Group.txt
index 6365d9f7d..6edc0225b 100644
--- a/user-doc/Group.txt
+++ b/user-doc/Group.txt
@@ -56,12 +56,12 @@ See http://en.wikipedia.org/wiki/Ramachandran_plot
 The following example shows how to use \ref MOLINFO with \ref TORSION to calculate the torsion angles phi and psi for the first and fourth residue
 of the protein:
 
-\verbatim
+\plumedfile
 MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
 t1: TORSION ATOMS=@phi-3
 t2: TORSION ATOMS=@psi-4
 PRINT ARG=t1,t2 FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 \subsection pbc Broken Molecules and PBC 
 
@@ -76,11 +76,11 @@ using the ALIGN_ATOMS keyword.  In PLUMED 2 the same effect can be achieved usin
 
 The following input computes the end-to-end distance for a polymer of 100 atoms and keeps it at a value around 5.
 
-\verbatim
+\plumedfile
 WHOLEMOLECULES ENTITY0=1-100
 e2e: DISTANCE ATOMS=1,100 NOPBC
 RESTRAINT ARG=e2e KAPPA=1 AT=5
-\endverbatim
+\endplumedfile
 
 Notice that NOPBC is used to be sure in \ref DISTANCE that if the end-to-end distance is larger than half the simulation box the distance 
 is compute properly. Also notice that, since many MD codes break molecules across cell boundary, it might be necessary to use the 
@@ -90,12 +90,12 @@ Notice that most expressions are invariant with respect to a change in the order
 but some of them depend on that order. E.g., with \ref WHOLEMOLECULES it could be useful to
 specify atom lists in a reversed order.
 
-\verbatim
+\plumedfile
 # to see the effect, one could dump the atoms as they were before molecule reconstruction:
 # DUMPATOMS FILE=dump-broken.xyz ATOMS=1-20
 WHOLEMOLECULES STRIDE=1 ENTITY0=1-20
 DUMPATOMS FILE=dump.xyz ATOMS=1-20
-\endverbatim
+\endplumedfile
 
 Notice that there are other ways to manipulate the coordinates stored within PLUMED:
 - Using the \subpage FIT_TO_TEMPLATE they can be aligned to a template structure.
@@ -116,10 +116,10 @@ To specify to a colvar that you want to use the position of a virtual atom to ca
 in your system you simply use the label for your virtual atom in place of the usual numerical index.  Virtual
 atoms and normal atoms can be mixed together in the input to colvars as shown below:
 
-\verbatim
+\plumedfile
 COM ATOMS=1,10 LABEL=com1
 DISTANCE ATOMS=11,com1
-\endverbatim 
+\endplumedfile
 
 If you don't want to calculate CVs from the virtual atom.  That is to say you just want to monitor the position of a virtual atom 
 (or any set of atoms) over the course of your trajectory you can do this using \ref DUMPATOMS.
diff --git a/user-doc/Misc.txt b/user-doc/Misc.txt
index 284a2a6e9..4fceb133d 100644
--- a/user-doc/Misc.txt
+++ b/user-doc/Misc.txt
@@ -19,12 +19,12 @@ particular simulation you might find it useful to put comments in your input fil
 comments can be added using a # sign.  On any given line everything after the # sign is ignored so 
 erm... yes add lines of comments or trailing comments to your hearts content as shown below (using Shakespeare is optional):
 
-\verbatim
+\plumedfile
 # This is the distance between two atoms:
 DISTANCE ATOM=1,2 LABEL=d1
 UPPER_WALLS ARG=d1 AT=3.0 KAPPA=3.0 LABEL=Snout # In this same interlude it doth befall.
 # That I, one Snout by name, present a wall.
-\endverbatim
+\endplumedfile
 (see \ref DISTANCE and \ref UPPER_WALLS)
 
 An alternative to including comments in this way is to use line starting ENDPLUMED.  Everything in the PLUMED input after this
@@ -37,14 +37,14 @@ We at PLUMED are aware of this fact and thus have provided a way of doing line c
 easier - aren't we kind?  Well no not really, we have to use this code too.  Anyway, you can do continuations by using the "..." syntax
 as this makes this: 
 
-\verbatim
+\plumedfile
 DISTANCES ATOMS1=1,300 ATOMS2=1,400 ATOMS3=1,500 LABEL=dist
-\endverbatim
+\endplumedfile
 (see \ref DISTANCES)
 
 equivalent to this:
 
-\verbatim
+\plumedfile
 DISTANCES ...
   LABEL=dist
 # we can also insert comments here
@@ -54,11 +54,11 @@ DISTANCES ...
 #empty lines are also allowed
 
 ... DISTANCES
-\endverbatim
+\endplumedfile
 
 Notice that the closing `...` is followed by the word `DISTANCES`. This is optional, but might be
 useful to find more easily which is the matching start of the statement. The following is equally correct
-\verbatim
+\plumedfile
 DISTANCES ...
   LABEL=dist
 # we can also insert comments here
@@ -68,12 +68,12 @@ DISTANCES ...
 #empty lines are also allowed
 
 ...
-\endverbatim
+\endplumedfile
 
 Notice that PLUMED makes a check that the word following the closing `...` is actually identical to
 the first word in the line with the first `...`. If not, it will throw an error.
 Also notice that you might put more than one word in the first line. E.g.
-\verbatim
+\plumedfile
 DISTANCES LABEL=dist ...
 # we can also insert comments here
   ATOMS1=1,300
@@ -81,9 +81,9 @@ DISTANCES LABEL=dist ...
   ATOMS2=1,400 ATOMS3=1,500
 #empty lines are also allowed
 ...
-\endverbatim
+\endplumedfile
 or, equivalently,
-\verbatim
+\plumedfile
 dist: DISTANCES ...
 # we can also insert comments here
   ATOMS1=1,300
@@ -91,7 +91,7 @@ dist: DISTANCES ...
   ATOMS2=1,400 ATOMS3=1,500
 #empty lines are also allowed
 ...  
-\endverbatim
+\endplumedfile
 
 
 \page VimSyntax Using VIM syntax file
@@ -171,11 +171,11 @@ Add to your `.vimrc` file the following commands
 :set modelines=5
 \endverbatim
 Then, at the beginning of your PLUMED input file, put the following comment:
-\verbatim
+\plumedfile
 # vim:ft=plumed
 d: DISTANCE ATOMS=1,2
 RESTRAINT ARG=d AT=0.0 KAPPA=1.0
-\endverbatim
+\endplumedfile
 Now, every time you open this file, you will see it highlighted.
 
 \par Syntax highlighting
@@ -220,27 +220,27 @@ Folded lines can be expanded with `zo` and folded with `zc`. Look at VIM documen
 learn more.
 In case you want to use this feature, we suggest you to put both label
 and action type on the first line of multi-line statements. E.g.
-\verbatim
+\plumedfile
 m: METAD ...
   ARG=d
   HEIGHT=1.0
   SIGMA=0.5
   PACE=100
 ...
-\endverbatim
+\endplumedfile
 will be folded to
 \verbatim
 +--  6 lines: m: METAD ...------------------------------------------------------
 \endverbatim
 and
-\verbatim
+\plumedfile
 METAD LABEL=m ...
   ARG=d
   HEIGHT=1.0
   SIGMA=0.5
   PACE=100
 ...
-\endverbatim
+\endplumedfile
 will be folded to
 \verbatim
 +--  6 lines: METAD LABEL=m ...-------------------------------------------------
@@ -293,9 +293,9 @@ Try to do the following. Enable plumed syntax:
 :set ft=plumed
 \endverbatim
 Then add the following line
-\verbatim
+\plumedfile
 DISTANCE
-\endverbatim
+\endplumedfile
 Now, in normal mode, go with the cursor on the `DISTANCE` line and type
 \verbatim
 :PHelp
@@ -365,38 +365,37 @@ vimrc file.
 
 If, for some reason, you want to spread your PLUMED input over a number of files you can use \subpage INCLUDE as shown below:
 
-\verbatim
+\plumedfile
 INCLUDE FILE=filename
-\endverbatim
+\endplumedfile
 
 So, for example, a single "plumed.dat" file:
 
-\verbatim
+\plumedfile
 DISTANCE ATOMS=0,1 LABEL=dist
 RESTRAINT ARG=dist
-\endverbatim
-(see \ref DISTANCE and \ref RESTRAINT)
+\endplumedfile
 
 could be split up into two files as shown below:
  
-\verbatim
+\plumedfile
 DISTANCE ATOMS=0,1 LABEL=dist
 INCLUDE FILE=toBeIncluded.dat
-\endverbatim
+\endplumedfile
 plus a "toBeIncluded.dat" file
-\verbatim
+\plumedfile
 RESTRAINT ARG=dist
-\endverbatim
+\endplumedfile
 
 However, when you do this it is important to recognise that \ref INCLUDE is a real directive that is only resolved
 after all the \ref comments have been stripped and the \ref ContinuationLines have been unrolled.  This means it
 is not possible to do things like:
 
-\verbatim
+\plumedfile
 # this is wrong:
 DISTANCE INCLUDE FILE=options.dat
 RESTRAINT ARG=dist
-\endverbatim
+\endplumedfile
 
 \page load Loading shared libraries
 
@@ -405,9 +404,9 @@ PLUMED libraries.  Alternatively, if you want to keep your code independent from
 so you can release it independely - we won't be offended), then you can create your own dynamic library.  To use this 
 in conjuction with PLUMED you can then load it at runtime by using the \subpage LOAD keyword as shown below:
 
-\verbatim
+\plumedfile
 LOAD FILE=library.so
-\endverbatim
+\endplumedfile
  
 N.B.  If your system uses a different suffix for dynamic libraries (e.g. macs use .dylib) then PLUMED will try to 
 automatically adjust the suffix accordingly.
diff --git a/user-doc/Performances.txt b/user-doc/Performances.txt
index 57bc2c66a..00fc973ac 100644
--- a/user-doc/Performances.txt
+++ b/user-doc/Performances.txt
@@ -105,12 +105,12 @@ a factor equal to STRIDE.
 This technique can allow your simulation to run faster if you need
 the apply a bias potential on some very expensive collective variable.
 Consider the following input:
-\verbatim
+\plumedfile
 c1: COM ATOMS=1-1000
 c2: COM ATOMS=1001-2000
 d:  DISTANCE ATOMS=c1,c2
 METAD ARG=d HEIGHT=1 SIGMA=0.1 BIASFACTOR=5 PACE=500
-\endverbatim
+\endplumedfile
 This performs a \ref METAD simulation biasing the distance between two
 centers of mass. Since computing these centers requires a lot of atoms
 to be imported from the MD engine, it could slow down significantly the
@@ -118,12 +118,12 @@ simulation. Notice that whereas the bias is changed every PACE=500 steps,
 it is applied every STRIDE step, where STRIDE=1 by default.
 The following input could lead to a significantly faster simulation at the price
 of a negligible systematic error
-\verbatim
+\plumedfile
 c1: COM ATOMS=1-1000
 c2: COM ATOMS=1001-2000
 d:  DISTANCE ATOMS=c1,c2
 METAD ARG=d HEIGHT=1 SIGMA=0.1 BIASFACTOR=5 PACE=500 STRIDE=2
-\endverbatim
+\endplumedfile
 Similarly, the STRIDE keyword can be used with other biases (e.g. \ref RESTRAINT).
 
 The technique is discussed in details here \cite Ferrarotti2015.
@@ -133,9 +133,9 @@ See also \subpage EFFECTIVE_ENERGY_DRIFT.
 
 Whenever you have a multicolvar action such as:
 
-\verbatim
+\plumedfile
 COORDINATIONNUMBER SPECIES=1-100 SWITCH={RATIONAL R_0=1. D_MAX=3.0} MORE_THAN={RATIONAL R_0=6.0 NN=6 MM=12 D_0=0}
-\endverbatim
+\endplumedfile
 
 You will get a collosal speedup by specifying the D_MAX keyword in all switching functions that act on distances.
 D_MAX tells PLUMED that the switching function is strictly zero if the distance is greater than this value.  As a result
diff --git a/user-doc/Regex.txt b/user-doc/Regex.txt
index b24524a0f..93561f632 100644
--- a/user-doc/Regex.txt
+++ b/user-doc/Regex.txt
@@ -14,14 +14,14 @@ Regular expressions are enclosed in round braces and must not contain spaces (th
 names have no spaces indeed, so why use them?).
 
 As an example then command
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,2 COMPONENTS
 PRINT ARG=(d1\.[xy])   STRIDE=100 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 will cause both the d1.x and d1.y components of the DISTANCE action to be printed out in the order that they are created by plumed.
 The "." character must be escaped in order to interpret it as a literal ".". An unescaped dot is a wildcard which is matched by any character,
 So as an example
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,2 COMPONENTS
 dxy: DISTANCE ATOMS=1,3
 
@@ -30,26 +30,26 @@ PRINT ARG=(d1.[xy])   STRIDE=100 FILE=colvar FMT=%8.4f
 
 # while this will match d1.x,d1.y only
 PRINT ARG=(d1\.[xy])   STRIDE=100 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 You can include more than one regular expression by using comma separated regular expressions 
 
-\verbatim
+\plumedfile
 t1: TORSION ATOMS=5,7,9,15
 t2: TORSION ATOMS=7,9,15,17
 d1: DISTANCE ATOMS=7,17 COMPONENTS
 PRINT ARG=(d1\.[xy]),(t[0-9]) STRIDE=100 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 (this selects t1,t2,d1.x and d2.x) Be aware that if you have overlapping selection they will be duplicated so it 
 a better alternative is to use the "or" operator "|". 
 
-\verbatim
+\plumedfile
 t1: TORSION ATOMS=5,7,9,15
 t2: TORSION ATOMS=7,9,15,17
 d1: DISTANCE ATOMS=7,17 COMPONENTS
 PRINT ARG=(d1\.[xy]|t[0-9]) STRIDE=100 FILE=colvar FMT=%8.4f
-\endverbatim
+\endplumedfile
 
 this selects the same set of arguments as the previous example.
 
diff --git a/user-doc/Syntax.txt b/user-doc/Syntax.txt
index 47c73c1f2..277d339e9 100644
--- a/user-doc/Syntax.txt
+++ b/user-doc/Syntax.txt
@@ -25,18 +25,16 @@ is enclosed in curly braces (e.g. ATOMS={1 2 3 4}).  Please note that you can sp
 The most important of these keywords is the label keyword as it is only by using these labels that we can pass data 
 from one action to another.  As an example if you do:
 
-\verbatim
+\plumedfile
 DISTANCE ATOMS=1,2
-\endverbatim
-(see \ref DISTANCE)
+\endplumedfile
 
 Then PLUMED will do nothing other than read in your input file.  In contrast if you do:
 
-\verbatim
+\plumedfile
 DISTANCE ATOMS=1,2 LABEL=d1
 PRINT ARG=d1 FILE=colvar STRIDE=10
-\endverbatim
-(see \ref PRINT)
+\endplumedfile
 
 then PLUMED will print out the value of the distance between atoms 1 and 2 every 10 steps to the file colvar as you have told
 PLUMED to take the value calculated by the action d1 and to print it. You can use any character string to label your actions
@@ -46,10 +44,10 @@ code-generated groups of atoms and to give labels to any Actions for which the u
 Notice that if a word followed by a column is added at the beginning of the line (e.g. pippo:), PLUMED automatically
 removes it and adds an equivalent label (LABEL=pippo).
 Thus, a completely equivalent result can be obtained with the following shortcut:
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,2
 PRINT ARG=d1 FILE=colvar STRIDE=10
-\endverbatim
+\endplumedfile
 
 Also notice that all the actions can be labeled, and that many actions besides normal collective variables can define
 one or more value, which can be then referred using the corresponding label.
-- 
GitLab