diff --git a/CHANGES/v2.4.md b/CHANGES/v2.4.md
index da6025dd115718e7c37f2c9887cdcd65e3e226d4..399ff9d990af2e02f2e1bb4bd5c7c0625a0c4d52 100644
--- a/CHANGES/v2.4.md
+++ b/CHANGES/v2.4.md
@@ -207,8 +207,7 @@ For users:
   - Fixed some performances regression issue with OpenMP
   - Updated NAMD patches to version 2.12 and 2.13. Old patches have been removed.
   - GROMACS patch for gromacs-2018.4.
-  - Fixed a bug in CS2BACKBONE when using more than 2 chains
-  - Fixed a threadsafety issue using forces on histograms
+  - Fixed a threadsafety issue using forces on \ref HISTOGRAM 
   - Fixed error message suggesting wrong actions (see \issue{421}).
 
 For developers:
diff --git a/CHANGES/v2.5.md b/CHANGES/v2.5.md
index b2ad6fc562c9a544e33b6b2f6c5e772a5122cdde..5c63c2230ff85fd6e7f0a54fa4461dbb33932e7b 100644
--- a/CHANGES/v2.5.md
+++ b/CHANGES/v2.5.md
@@ -49,10 +49,20 @@ Changes from version 2.4 which are relevant for users:
   - \ref completion (used to generate command line completion scripts).
   - \ref pdbrenumber (see \issue{371}).
 
+- New modules:
+  - A new PIV module has been included, contributed by Silvio Pipolo and Fabio Pietrucci.
+    This module implements the following collective variable:
+    - \ref PIV
+  - A new LOGMFD module has been included, contributed by Tetsuya Morishita.
+    This module implements the following bias:
+    - \ref LOGMFD
+
 - Changes in the ISDB module
   - \ref CS2BACKBONE is now mpi parallelized in particular with DOSCORE and CAMSHIFT
   - \ref SAXS has an additional implementation based on Bessel functions that can be faster for large systems (new keyword BESSEL)
   - \ref SAXS keyword SCEXP has been renamed into SCALEINT
+  - \ref SAXS includes the MARTINI bead structure factors for Proteins and Nucleic Acids
+  - \ref SAXS includes a GPU implementation based on arrayfire (need to be linked at compile time) that can be activated with GPU
   - \ref METAINFERENCE and all related methods has a new keyword REGRES_ZERO to scale data using a linear scale fit
   - \ref CALIBER new bias to perform Maximum Caliber replica-averaged restrained simulations 
 
@@ -144,6 +154,3 @@ Changes from version 2.4 which are relevant for developers:
 - Using CXX compiler to link the main program.
 - plumed can be compiled with arrayfire to enable for gpu code. \ref SAXS collective variable is available as part of the isdb module to provide an example of a gpu implementation for a CV
 
-Fixes done after beta (this list will be merged to the one above):
-
-
diff --git a/src/bias/MetaD.cpp b/src/bias/MetaD.cpp
index ce8199a6fc57edba6df762d7d40849dea827e259..c0ef7f1ae76823de71d831ef033949ccb3f87eb1 100644
--- a/src/bias/MetaD.cpp
+++ b/src/bias/MetaD.cpp
@@ -135,16 +135,16 @@ for replica exchange methods to be computed correctly.
 Multiple walkers  \cite multiplewalkers can also be used. See below the examples.
 
 
-The \\f$c(t)\\f$ reweighting factor can also be calculated on the fly using the equations
+The \f$c(t)\f$ reweighting factor can also be calculated on the fly using the equations
 presented in \cite Tiwary_jp504920s.
-The expression used to calculate \\f$c(t)\\f$ follows directly from Eq. 3 in \cite Tiwary_jp504920s,
-where \\f$F(\vec{s})=-\gamma/(\gamma-1) V(\vec{s})\\f$.
+The expression used to calculate \f$c(t)\f$ follows directly from Eq. 3 in \cite Tiwary_jp504920s,
+where \f$F(\vec{s})=-\gamma/(\gamma-1) V(\vec{s})\f$.
 This gives smoother results than equivalent Eqs. 13 and Eqs. 14 in that paper.
-The \\f$c(t)\\f$ is given by the rct component while the bias
-normalized by \\f$c(t)\\f$ is given by the rbias component (rbias=bias-rct) which can be used
+The \f$c(t)\f$ is given by the rct component while the bias
+normalized by \f$c(t)\f$ is given by the rbias component (rbias=bias-rct) which can be used
 to obtain a reweighted histogram.
-The calculation of \\f$c(t)\\f$ is enabled by using the keyword CALC_RCT.
-By default \\f$c(t)\\f$ is updated every time the bias changes, but if this slows down the simulation
+The calculation of \f$c(t)\f$ is enabled by using the keyword CALC_RCT.
+By default \f$c(t)\f$ is updated every time the bias changes, but if this slows down the simulation
 the keyword RCT_USTRIDE can be set to a value higher than 1.
 This option requires that a grid is used.
 
@@ -233,7 +233,7 @@ one update and the other. Since version 2.2.5, hills files are automatically
 flushed every WALKERS_RSTRIDE steps.
 
 \par
-The \\f$c(t)\\f$ reweighting factor can be calculated on the fly using the equations
+The \f$c(t)\f$ reweighting factor can be calculated on the fly using the equations
 presented in \cite Tiwary_jp504920s as described above.
 This is enabled by using the keyword CALC_RCT,
 and can be done only if the bias is defined on a grid.
@@ -248,9 +248,9 @@ METAD ...
 \endplumedfile
 Here we have asked that the calculation is performed every 10 hills deposition by using
 RCT_USTRIDE keyword. If this keyword is not given, the calculation will
-by default be performed every time the bias changes. The \\f$c(t)\\f$ reweighting factor will be given
+by default be performed every time the bias changes. The \f$c(t)\f$ reweighting factor will be given
 in the rct component while the instantaneous value of the bias potential
-normalized using the \\f$c(t)\\f$ reweighting factor is given in the rbias component
+normalized using the \f$c(t)\f$ reweighting factor is given in the rbias component
 [rbias=bias-rct] which can be used to obtain a reweighted histogram or
 free energy surface using the \ref HISTOGRAM analysis.
 
diff --git a/src/isdb/SAXS.cpp b/src/isdb/SAXS.cpp
index 4712736b13bb08c91334a9c220b8ef0ac3b38976..b568ceb34d9fc3fe042161184d05179032a05c71 100644
--- a/src/isdb/SAXS.cpp
+++ b/src/isdb/SAXS.cpp
@@ -354,8 +354,10 @@ SAXS::SAXS(const ActionOptions&ao):
     if(bessel&&i>0&&q_list[i]<q_list[i-1]) plumed_merror("With BESSEL the Q values should be ordered from the smallest to the largest");
   }
   log<<"  Bibliography ";
-  log<<plumed.cite("Jussupow, et al. (in preparation)");
-  if(martini)   log<<plumed.cite("Niebling, Björling, Westenhoff, J Appl Crystallogr 47, 1190–1198 (2014).");
+  if(martini) {
+    log<<plumed.cite("Niebling, Björling, Westenhoff, J Appl Crystallogr 47, 1190–1198 (2014).");
+    log<<plumed.cite("Paissoni, Jussupow, Camilloni, bioRxiv 498147.");
+  }
   if(atomistic) {
     log<<plumed.cite("Fraser, MacRae, Suzuki, J. Appl. Crystallogr., 11, 693–694 (1978).");
     log<<plumed.cite("Brown, Fox, Maslen, O'Keefe, Willis, International Tables for Crystallography C, 554–595 (International Union of Crystallography, 2006).");
diff --git a/src/logmfd/.gitignore b/src/logmfd/.gitignore
new file mode 100644
index 0000000000000000000000000000000000000000..4888552338441b5d99a97ec7a2cb2ce448e9e4e6
--- /dev/null
+++ b/src/logmfd/.gitignore
@@ -0,0 +1,9 @@
+/*
+# in this directory, only accept source, Makefile and README
+!/.gitignore
+!/*.c
+!/*.cpp
+!/*.h
+!/Makefile
+!/README
+!/module.type
diff --git a/src/piv/PIV.cpp b/src/piv/PIV.cpp
index 8aaa1a96f76d7a93204dde6e2ef0023255f017dc..732f61efee3fdcd5346f42f70ce6f0e50660d009 100644
--- a/src/piv/PIV.cpp
+++ b/src/piv/PIV.cpp
@@ -36,7 +36,7 @@ namespace PLMD
 namespace piv
 {
 
-//+PLUMEDOC COLVAR PIV
+//+PLUMEDOC PIVMOD_COLVAR PIV
 /*
 Calculates the PIV-distance.
 
diff --git a/user-doc/Analysis.md b/user-doc/Analysis.md
index 623d12f6443194baa919b8e19df278a160a75fb0..b3b602218fcdef4470c05b8226bd077f708c566f 100644
--- a/user-doc/Analysis.md
+++ b/user-doc/Analysis.md
@@ -107,11 +107,11 @@ within colvars and functions.  One place where this is very useful is when you a
 not you have implemented the derivatives of a new collective variables correctly.  So for example if
 we wanted to do such a test on the distance CV we would employ an input file something like this:
 
-\verbatim
+\plumedfile
 d1: DISTANCE ATOMS=1,2
 d1n: DISTANCE ATOMS=1,2 NUMERICAL_DERIVATIVES
 DUMPDERIVATIVES ARG=d1,d1n FILE=derivatives
-\endverbatim
+\endplumedfile
 
 The first of these two distance commands calculates the analytical derivatives of the distance
 while the second calculates these derivatives numerically.  Obviously, if your CV is implemented
@@ -170,12 +170,12 @@ that are available in PLUMED are as follows
 
 In general most of these landmark selection algorithms must be used in tandem with a \ref dissimilaritym "dissimilarity matrix" object as as follows:
 
-\verbatim
+\plumedfile
 data: COLLECT_FRAMES ARG=d1 STRIDE=1
 ss1: EUCLIDEAN_DISSIMILARITIES USE_OUTPUT_DATA_FROM=data 
 ll2: LANDMARK_SELECT_FPS USE_OUTPUT_DATA_FROM=ss1 NLANDMARKS=300
 OUTPUT_COLVAR_FILE USE_OUTPUT_DATA_FROM=ll2 FILE=mylandmarks
-\endverbatim
+\endplumedfile
 
 When landmark selection is performed in this way a weight is ascribed to each of the landmark configurations.  This weight is
 calculated by summing the weights of all the trajectory frames in each of the landmarks Voronoi polyhedron 
@@ -207,22 +207,22 @@ the following <a href="https://www.youtube.com/watch?v=ofC2qz0_9_A&feature=youtu
 
 Within PLUMED running an input to run a dimensionality reduction algorithm can be as simple as:
 
-\verbatim
+\plumedfile
 data: COLLECT_FRAMES STRIDE=1 ARG=d1
 ss1: EUCLIDEAN_DISSIMILARITIES USE_OUTPUT_DATA_FROM=data 
 mds: CLASSICAL_MDS USE_OUTPUT_DATA_FROM=ss1 NLOW_DIM=2
-\endverbatim
+\endplumedfile
 
 Where we have to use the \ref EUCLIDEAN_DISSIMILARITIES action here in order to calculate the matrix of dissimilarities between trajectory frames.
 We can even throw some landmark selection into this procedure and perform
 
-\verbatim
+\plumedfile
 data: COLLECT_FRAMES STRIDE=1 ARG=d1
 matrix: EUCLIDEAN_DISSIMILARITIES USE_OUTPUT_DATA_FROM=data
 ll2: LANDMARK_SELECT_FPS USE_OUTPUT_DATA_FROM=matrix NLANDMARKS=300
 mds: CLASSICAL_MDS USE_OUTPUT_DATA_FROM=ll2 NLOW_DIM=2
 osample: PROJECT_ALL_ANALYSIS_DATA USE_OUTPUT_DATA_FROM=matrix PROJECTION=smap
-\endverbatim
+\endplumedfile
 
 Notice here that the final command allows us to calculate the projections of all the non-landmark points that were collected by the action with
 label matrix.
diff --git a/user-doc/PIVMOD.md b/user-doc/PIVMOD.md
new file mode 100644
index 0000000000000000000000000000000000000000..6ec666f580745c62522c453d5479dc6676f1c84e
--- /dev/null
+++ b/user-doc/PIVMOD.md
@@ -0,0 +1,26 @@
+\page PIVMOD PIV collective variable
+
+<!-- 
+description: To be completed
+authors: To be completed
+reference: To be completed
+-->
+
+## Overview
+
+To be completed
+
+## Installation 
+This module is not installed by default. Add '\-\-enable-modules=piv' to your './configure' command when building PLUMED to enable these features.
+
+## Usage
+Currently, all features of the PIV module are included in a single PIV collective variable: \ref PIV
+
+## Module Contents
+- \subpage PIVMODColvar
+
+\page PIVMODColvar CVs Documentation
+
+The following list contains descriptions of biases developed for the PLUMED-PIV module. They can be used in combination with other actions outside of the PIV module.
+
+@PIVMOD_COLVAR@
diff --git a/user-doc/Performances.md b/user-doc/Performances.md
index b4dad4c1c01ee5f51fa3a2c62fd2d76ecde773b0..d62c18d19f51286930e3246ace3959ac8567c5cf 100644
--- a/user-doc/Performances.md
+++ b/user-doc/Performances.md
@@ -271,12 +271,12 @@ You are done!
 
 In some case using a custom expression is almost as fast as using a hard-coded
 function. For instance, with an input like this one:
-\verbatim
+\plumedfile
 ...
 c: COORDINATION GROUPA=1-108 GROUPB=1-108 R_0=1
 dfast: COORDINATION GROUPA=1-108 GROUPB=1-108 SWITCH={CUSTOM FUNC=1/(1+x2^3) R_0=1}
 ...
-\endverbatim
+\endplumedfile
 I (GB) obtained the following timings (on a Macbook laptop):
 \verbatim
 ...