diff --git a/gromacs/CMakeLists.txt.Lepton.diff b/gromacs/CMakeLists.txt.Lepton.patch
similarity index 100%
rename from gromacs/CMakeLists.txt.Lepton.diff
rename to gromacs/CMakeLists.txt.Lepton.patch
diff --git a/gromacs/CMakeLists.txt.diff b/gromacs/CMakeLists.txt.patch
similarity index 100%
rename from gromacs/CMakeLists.txt.diff
rename to gromacs/CMakeLists.txt.patch
diff --git a/gromacs/gmxManageColvars.cmake.diff b/gromacs/gmxManageColvars.cmake.patch
similarity index 100%
rename from gromacs/gmxManageColvars.cmake.diff
rename to gromacs/gmxManageColvars.cmake.patch
diff --git a/lammps/COLVARS.cmake.diff b/lammps/COLVARS.cmake.patch
similarity index 100%
rename from lammps/COLVARS.cmake.diff
rename to lammps/COLVARS.cmake.patch
diff --git a/namd/Make.depends.patch b/namd/Make.depends.patch
new file mode 100644
index 000000000..0c02bf721
--- /dev/null
+++ b/namd/Make.depends.patch
@@ -0,0 +1,59 @@
+diff --git a/Make.depends b/Make.depends
+index 994b3914..1d54373f 100644
+--- a/Make.depends
++++ b/Make.depends
+@@ -6398,6 +6398,45 @@ obj/GlobalMasterTMD.o: \
+ 	src/fitrms.h \
+ 	src/Debug.h
+ 	$(CXX) $(CXXFLAGS) $(COPTO)obj/GlobalMasterTMD.o $(COPTC) src/GlobalMasterTMD.C
++obj/GlobalMasterColvars.o: \
++	obj/.exists \
++	src/GlobalMasterColvars.C \
++	src/Controller.h \
++	src/common.h \
++	src/Vector.h \
++	src/Tensor.h \
++	src/Lattice.h \
++	src/NamdTypes.h \
++	src/CudaRecord.h \
++	src/HipDefines.h \
++	src/ResizeArray.h \
++	src/ResizeArrayRaw.h \
++	src/fstream_namd.h \
++	src/DataExchanger.h \
++	src/main.h \
++	inc/DataExchanger.decl.h \
++	src/Pointer.h \
++	src/GlobalMaster.h \
++	src/GlobalMasterColvars.h \
++	src/colvarproxy_namd.h \
++	src/colvarproxy_namd_version.h \
++	src/SimParameters.h \
++	src/imd.h \
++	src/MGridforceParams.h \
++	src/GroupRestraintsParam.h \
++	src/Random.h \
++	src/ConfigList.h \
++	colvars/src/colvarmodule.h \
++	colvars/src/colvarproxy.h \
++	colvars/src/colvarmodule.h \
++	colvars/src/colvartypes.h \
++	colvars/src/colvarproxy_io.h \
++	colvars/src/colvarproxy_replicas.h \
++	colvars/src/colvarproxy_system.h \
++	colvars/src/colvarproxy_tcl.h \
++	colvars/src/colvarproxy_volmaps.h \
++	colvars/src/colvarvalue.h
++	$(CXX) $(CXXFLAGS) $(COPTO)obj/GlobalMasterColvars.o $(COPTC) src/GlobalMasterColvars.C
+ obj/Matrix4Symmetry.o: \
+ 	obj/.exists \
+ 	src/Matrix4Symmetry.C \
+@@ -6637,6 +6676,8 @@ obj/colvarproxy_namd.o: \
+ 	src/MStream.h \
+ 	plugins/include/molfile_plugin.h \
+ 	plugins/include/vmdplugin.h \
++	src/GlobalMasterColvars.h \
++	src/GlobalMaster.h \
+ 	src/GridForceGrid.inl \
+ 	src/PDB.h \
+ 	src/parm.h \
diff --git a/namd/Makefile.patch b/namd/Makefile.patch
index 580c60b5b..89cdf7611 100644
--- a/namd/Makefile.patch
+++ b/namd/Makefile.patch
@@ -1,13 +1,12 @@
 diff --git a/Makefile b/Makefile
-index fb23fcac..541fe8ca 100644
+index 960f7718..85126dff 100644
 --- a/Makefile
 +++ b/Makefile
-@@ -446,7 +446,7 @@ CXXMICFLAGS = $(CXXBASEFLAGS) $(CXXOPTS) $(CXXMICOPTS)
- CXXTHREADFLAGS = $(CXXBASEFLAGS) $(CXXTHREADOPTS)
- CXXSIMPARAMFLAGS = $(CXXBASEFLAGS) $(CXXSIMPARAMOPTS)
- CXXNOALIASFLAGS = $(CXXBASEFLAGS) $(CXXNOALIASOPTS)
--COLVARSCXXFLAGS = $(CXXBASEFLAGS) $(CXXOPTS) $(COPTI)$(LEPTONINCDIR) -DLEPTON -DLEPTON_USE_STATIC_LIBRARIES
-+COLVARSCXXFLAGS = $(CXXBASEFLAGS) $(CXXOPTS) $(COPTI)$(LEPTONINCDIR) -DLEPTON -DLEPTON_USE_STATIC_LIBRARIES $(EXTRACOLVARSFLAGS)
- GXXFLAGS = $(CXXBASEFLAGS) -DNO_STRSTREAM_H
- CFLAGS = $(COPTI)$(SRCDIR) $(TCL) $(COPTS) $(RELEASE) $(EXTRADEFINES) $(TRACEOBJDEF)
- PLUGINGCCFLAGS = $(COPTI)$(PLUGINSRCDIR) $(COPTI)$(PLUGININCDIR) $(COPTD)STATIC_PLUGIN
+@@ -231,6 +231,7 @@ OBJS = \
+ 	$(DSTDIR)/GlobalMasterTcl.o \
+ 	$(DSTDIR)/GlobalMasterSMD.o \
+ 	$(DSTDIR)/GlobalMasterTMD.o \
++	$(DSTDIR)/GlobalMasterColvars.o \
+ 	$(DSTDIR)/Matrix4Symmetry.o \
+ 	$(DSTDIR)/GlobalMasterSymmetry.o \
+ 	$(DSTDIR)/GlobalMasterFreeEnergy.o \
diff --git a/namd/config.patch b/namd/config.patch
index 76618f2e5..e69de29bb 100644
--- a/namd/config.patch
+++ b/namd/config.patch
@@ -1,59 +0,0 @@
-diff --git a/config b/config
-index ddf9f860..4148d87c 100755
---- a/config
-+++ b/config
-@@ -57,6 +57,8 @@ function error_syntax {
-   echo '  --cuda-dlink arch=<arch>,code=<code> (for cuFFT, may be repeated)'
-   echo '  --with-cuda-profiling (enables CUDA profiling with NVTX)'
-   echo '  --with-rocm-profiling (enables ROCm profiling with ROCtracer API)'
-+  echo '  --with-colvars-torch (enables torchann in Colvars)'
-+  echo '  --torch-prefix <directory containing libTorch>'
-   echo ''
-   if [ -n "${PRINT_ARCH_LIST+set}" ]; then
-     ARCH_PAT=''
-@@ -152,6 +154,7 @@ function error_exists {
-   use_mkl=0
-   use_cuda=0
-   use_hip=0
-+  use_colvars_torch=0
-   use_cuda_prof=0
-   use_rocm_prof=0
-   use_memopt=0
-@@ -257,6 +260,19 @@ function error_exists {
-         ARCH_SUFFIX_ARG=$ARCH_SUFFIX_ARG-$1
-       ;;
- 
-+      --with-colvars-torch)
-+        use_colvars_torch=1
-+      ;;
-+
-+      --torch-prefix)
-+        shift
-+        if [ ! -d "$1" ]; then
-+          echo "ERROR: No such directory $1"
-+          error_syntax
-+        fi
-+        TORCH_PREFIX=$1
-+      ;;
-+
-       --with-debug)
-         use_debug=1
-       ;;
-@@ -1015,7 +1031,16 @@ function error_exists {
-   if [ -n "$CC_OPTS" ]; then
-     echo "COPTS = $CC_OPTS" >> Make.config
-   fi
--
-+  if (( "$use_colvars_torch" )); then
-+    echo "TORCHDIR = $TORCH_PREFIX" >> Make.config
-+    echo 'TORCHINCFLAGS = $(COPTI)$(TORCHDIR)/include $(COPTI)$(TORCHDIR)/include/torch/csrc/api/include' >> Make.config
-+    echo 'EXTRACOLVARSFLAGS = -std=c++17 -DCOLVARS_TORCH $(TORCHINCFLAGS)' >> Make.config
-+    if [[ $use_cuda && -f "${TORCH_PREFIX}/lib/libtorch_cuda.so" ]]; then 
-+      echo 'EXTRALINKLIBS = -Wl,-rpath,$(TORCHDIR)/lib -L$(TORCHDIR)/lib -ltorch -ltorch_cpu -lc10 -ltorch_cuda -lc10_cuda' >> Make.config
-+    else
-+      echo 'EXTRALINKLIBS = -Wl,-rpath,$(TORCHDIR)/lib -L$(TORCHDIR)/lib -ltorch -ltorch_cpu -lc10' >> Make.config
-+    fi
-+  fi
- 
-   if (( $use_debug )); then
-     echo 'CXXOPTS = -g' >> Make.config
diff --git a/namd/src/.clang-format b/namd/src/.clang-format
new file mode 100644
index 000000000..66e78c8b5
--- /dev/null
+++ b/namd/src/.clang-format
@@ -0,0 +1,14 @@
+---
+Language: Cpp
+ColumnLimit: 100
+UseTab: Never
+TabWidth: 2
+IndentWidth: 2
+ObjCBlockIndentWidth: 2
+PenaltyBreakAssignment: 2
+AccessModifierOffset: -2
+BreakBeforeBraces: WebKit
+EmptyLineAfterAccessModifier: Leave
+NamespaceIndentation: None
+MaxEmptyLinesToKeep: 2
+...
diff --git a/namd/src/GlobalMasterColvars.C b/namd/src/GlobalMasterColvars.C
new file mode 100644
index 000000000..1789867c1
--- /dev/null
+++ b/namd/src/GlobalMasterColvars.C
@@ -0,0 +1,31 @@
+// -*- c++ -*-
+
+#include "GlobalMaster.h"
+#include "GlobalMasterColvars.h"
+#include "colvarproxy_namd.h"
+
+
+GlobalMasterColvars::GlobalMasterColvars() : proxy(new colvarproxy_namd(this)) {}
+
+GlobalMasterColvars::~GlobalMasterColvars() { reset(); }
+
+void GlobalMasterColvars::calculate() { proxy->calculate(); }
+
+
+void GlobalMasterColvars::reset()
+{
+  // Unrequest all positions, total forces, etc from NAMD
+  modifyRequestedAtoms().clear();
+  modifyForcedAtoms().clear();
+  modifyAppliedForces().clear();
+
+  modifyRequestedGroups().clear();
+  modifyGroupForces().clear();
+
+#if NAMD_VERSION_NUMBER >= 34471681
+  modifyRequestedGridObjects().clear();
+  modifyGridObjForces().clear();
+#endif
+
+  requestTotalForce(false);
+}
diff --git a/namd/src/GlobalMasterColvars.h b/namd/src/GlobalMasterColvars.h
new file mode 100644
index 000000000..5ceaaf13c
--- /dev/null
+++ b/namd/src/GlobalMasterColvars.h
@@ -0,0 +1,155 @@
+#ifndef GLOBALMASTERCOLVARS_H
+#define GLOBALMASTERCOLVARS_H
+
+#include <memory>
+
+#include "GlobalMaster.h"
+
+
+class colvarproxy_namd;
+
+class GlobalMasterColvars : public GlobalMaster {
+public:
+
+  GlobalMasterColvars();
+
+  ~GlobalMasterColvars();
+
+  void reset();
+
+  void calculate() override;
+
+  inline void requestTotalForcePublic(bool yesno = true)
+  {
+    requestTotalForce(yesno);
+  }
+
+  inline AtomIDList &modifyRequestedAtomsPublic()
+  {
+    return modifyRequestedAtoms();
+  }
+
+  inline AtomIDList &modifyForcedAtomsPublic()
+  {
+    return modifyForcedAtoms();
+  }
+
+  inline ForceList &modifyAppliedForcesPublic()
+  {
+    return modifyAppliedForces();
+  }
+
+  inline ResizeArray<AtomIDList> &modifyRequestedGroupsPublic()
+  {
+    return modifyRequestedGroups();
+  }
+
+  inline ForceList &modifyGroupForcesPublic()
+  {
+    return modifyGroupForces();
+  }
+
+  inline IntList &modifyRequestedGridObjectsPublic()
+  {
+    return modifyRequestedGridObjects();
+  }
+
+  inline BigRealList &modifyGridObjForcesPublic()
+  {
+    return modifyGridObjForces();
+  }
+
+  inline AtomIDList::const_iterator getAtomIdBeginPublic()
+  {
+    return getAtomIdBegin();
+  }
+
+  inline AtomIDList::const_iterator getAtomIdEndPublic()
+  {
+    return getAtomIdEnd();
+  }
+
+  inline PositionList::const_iterator getAtomPositionBeginPublic()
+  {
+    return getAtomPositionBegin();
+  }
+
+  inline PositionList::const_iterator getGroupPositionBeginPublic()
+  {
+    return getGroupPositionBegin();
+  }
+
+  inline PositionList::const_iterator getGroupPositionEndPublic()
+  {
+    return getGroupPositionEnd();
+  }
+
+  inline ForceList::const_iterator getGroupTotalForceBeginPublic()
+  {
+    return getGroupTotalForceBegin();
+  }
+
+  inline ForceList::const_iterator getGroupTotalForceEndPublic()
+  {
+    return getGroupTotalForceEnd();
+  }
+
+  inline IntList::const_iterator getGridObjIndexBeginPublic()
+  {
+    return getGridObjIndexBegin();
+  }
+
+  inline IntList::const_iterator getGridObjIndexEndPublic()
+  {
+    return getGridObjIndexEnd();
+  }
+
+  inline BigRealList::const_iterator getGridObjValueBeginPublic()
+  {
+    return getGridObjValueBegin();
+  }
+
+  inline BigRealList::const_iterator getGridObjValueEndPublic()
+  {
+    return getGridObjValueEnd();
+  }
+
+  inline AtomIDList::const_iterator getForceIdBeginPublic()
+  {
+    return getForceIdBegin();
+  }
+
+  inline AtomIDList::const_iterator getForceIdEndPublic()
+  {
+    return getForceIdEnd();
+  }
+
+  inline ForceList::const_iterator getTotalForcePublic()
+  {
+    return getTotalForce();
+  }
+
+  inline void addReductionEnergyPublic(int reductionTag, BigReal energy)
+  {
+    addReductionEnergy(reductionTag, energy);
+  }
+
+  inline Lattice const *get_lattice() const
+  {
+    return lattice;
+  }
+
+protected:
+
+  std::unique_ptr<colvarproxy_namd> proxy;
+};
+
+
+namespace {
+  // Constants used for profiling CkLoop calls
+  constexpr int32_t GLOBAL_MASTER_CKLOOP_CALC_ITEM = 2000;
+  constexpr int32_t GLOBAL_MASTER_CKLOOP_CALC_BIASES = 2001;
+  constexpr int32_t GLOBAL_MASTER_CKLOOP_CALC_SCRIPTED_BIASES = 2002;
+}
+
+#endif
diff --git a/namd/src/colvarproxy_namd.C b/namd/src/colvarproxy_namd.C
index d84d87b9a..d4612686f 100644
--- a/namd/src/colvarproxy_namd.C
+++ b/namd/src/colvarproxy_namd.C
@@ -9,31 +9,34 @@
 
 #include <errno.h>
 
-#include "common.h"
-#include "fstream_namd.h"
+#include "ConfigList.h"
+#include "Controller.h"
+#include "DataExchanger.h"
 #include "Debug.h"
-#include "BackEnd.h"
-#include "InfoStream.h"
-#include "Node.h"
-#include "Molecule.h"
+#include "GlobalMaster.h"
+#include "GlobalMasterColvars.h"
 #include "GridForceGrid.h"
 #include "GridForceGrid.inl"
+#include "InfoStream.h"
+#include "Lattice.h"
+#include "Molecule.h"
+#include "NamdState.h"
+#include "NamdTypes.h"
+#include "Node.h"
 #include "PDB.h"
-#include "PDBData.h"
+#include "Random.h"
 #include "ReductionMgr.h"
+#include "ResizeArray.h"
 #include "ScriptTcl.h"
-#include "NamdState.h"
-#include "Controller.h"
-#include "PatchData.h"
+#include "SimParameters.h"
+#include "Vector.h"
+#include "common.h"
+#include "fstream_namd.h"
 
 #ifdef NAMD_TCL
 #include <tcl.h>
 #endif
 
-// For replica exchange
-#include "converse.h"
-#include "DataExchanger.h"
-
 #include "colvarmodule.h"
 #include "colvar.h"
 #include "colvarbias.h"
@@ -41,10 +44,10 @@
 #include "colvarproxy.h"
 #include "colvarproxy_namd.h"
 #include "colvarproxy_namd_version.h"
-#include "colvarscript.h"
 
 
-colvarproxy_namd::colvarproxy_namd()
+colvarproxy_namd::colvarproxy_namd(GlobalMasterColvars *gm)
+  : globalmaster(gm)
 {
   engine_name_ = "NAMD";
 #if CMK_SMP && USE_CKLOOP
@@ -60,7 +63,7 @@ colvarproxy_namd::colvarproxy_namd()
   }
 #endif
   first_timestep = true;
-  requestTotalForce(total_force_requested);
+  globalmaster->requestTotalForcePublic(total_force_requested);
 
   boltzmann_ = 0.001987191;
 
@@ -83,7 +86,7 @@ colvarproxy_namd::colvarproxy_namd()
   set_integration_timestep(simparams->dt);
   set_time_step_factor(simparams->globalMasterFrequency);
 
-  random = Random(simparams->randomSeed);
+  random.reset(new Random(simparams->randomSeed));
 
   // both fields are taken from data structures already available
   updated_masses_ = updated_charges_ = true;
@@ -220,7 +223,7 @@ int colvarproxy_namd::update_atoms_map(AtomIDList::const_iterator begin,
       // add it here anyway to avoid having to test for array boundaries at each step
       int const index = add_atom_slot(*a_i);
       atoms_map[*a_i] = index;
-      modifyRequestedAtoms().add(*a_i);
+      globalmaster->modifyRequestedAtomsPublic().add(*a_i);
       update_atom_properties(index);
     }
   }
@@ -265,8 +268,8 @@ int colvarproxy_namd::setup()
   }
 
   size_t n_group_atoms = 0;
-  for (int ig = 0; ig < modifyRequestedGroups().size(); ig++) {
-    n_group_atoms += modifyRequestedGroups()[ig].size();
+  for (int ig = 0; ig < globalmaster->modifyRequestedGroupsPublic().size(); ig++) {
+    n_group_atoms += globalmaster->modifyRequestedGroupsPublic()[ig].size();
   }
 
   log("updating group data ("+cvm::to_str(atom_groups_ids.size())+
@@ -274,7 +277,7 @@ int colvarproxy_namd::setup()
       cvm::to_str(n_group_atoms)+" atoms in total).\n");
 
   // Note: groupMassBegin, groupMassEnd may be used here, but they won't work for charges
-  for (int ig = 0; ig < modifyRequestedGroups().size(); ig++) {
+  for (int ig = 0; ig < globalmaster->modifyRequestedGroupsPublic().size(); ig++) {
 
     // update mass and charge
     update_group_properties(ig);
@@ -287,7 +290,7 @@ int colvarproxy_namd::setup()
 #if NAMD_VERSION_NUMBER >= 34471681
   log("updating grid object data ("+cvm::to_str(volmaps_ids.size())+
       " grid objects in total).\n");
-  for (int imap = 0; imap < modifyGridObjForces().size(); imap++) {
+  for (int imap = 0; imap < globalmaster->modifyGridObjForcesPublic().size(); imap++) {
     volmaps_new_colvar_forces[imap] = 0.0;
   }
 #endif
@@ -320,20 +323,7 @@ int colvarproxy_namd::reset()
 
   int error_code = COLVARS_OK;
 
-  // Unrequest all positions, total forces, etc from NAMD
-  modifyRequestedAtoms().clear();
-  modifyForcedAtoms().clear();
-  modifyAppliedForces().clear();
-
-  modifyRequestedGroups().clear();
-  modifyGroupForces().clear();
-
-#if NAMD_VERSION_NUMBER >= 34471681
-  modifyRequestedGridObjects().clear();
-  modifyGridObjForces().clear();
-#endif
-
-  requestTotalForce(false);
+  globalmaster->reset();
 
   atoms_map.clear();
 
@@ -348,6 +338,8 @@ void colvarproxy_namd::calculate()
 {
   errno = 0;
 
+  auto const step = globalmaster->step;
+
   if (first_timestep) {
 
     // First run after the proxy is constructed
@@ -388,6 +380,8 @@ void colvarproxy_namd::calculate()
   previous_NAMD_step = step;
   if (accelMDOn) update_accelMD_info();
 
+  auto *lattice = globalmaster->get_lattice();
+
   {
     Vector const a = lattice->a();
     Vector const b = lattice->b();
@@ -420,17 +414,17 @@ void colvarproxy_namd::calculate()
   // must delete the forces applied at the previous step: we can do
   // that because they have already been used and copied to other
   // memory locations
-  modifyForcedAtoms().clear();
-  modifyAppliedForces().clear();
+  globalmaster->modifyForcedAtomsPublic().clear();
+  globalmaster->modifyAppliedForcesPublic().clear();
 
   // If new atomic positions or forces have been requested by other
   // GlobalMaster objects, add these to the atom map as well
   size_t const n_all_atoms = Node::Object()->molecule->numAtoms;
   if ( (atoms_map.size() != n_all_atoms) ||
-       (int(atoms_ids.size()) < (getAtomIdEnd() - getAtomIdBegin())) ||
-       (int(atoms_ids.size()) < (getForceIdEnd() - getForceIdBegin())) ) {
-    update_atoms_map(getAtomIdBegin(), getAtomIdEnd());
-    update_atoms_map(getForceIdBegin(), getForceIdEnd());
+       (int(atoms_ids.size()) < (globalmaster->getAtomIdEndPublic() - globalmaster->getAtomIdBeginPublic())) ||
+       (int(atoms_ids.size()) < (globalmaster->getForceIdEndPublic() - globalmaster->getForceIdBeginPublic())) ) {
+    update_atoms_map(globalmaster->getAtomIdBeginPublic(), globalmaster->getAtomIdEndPublic());
+    update_atoms_map(globalmaster->getForceIdBeginPublic(), globalmaster->getForceIdEndPublic());
   }
 
   // prepare local arrays
@@ -456,9 +450,9 @@ void colvarproxy_namd::calculate()
       cvm::log("Updating positions arrays.\n");
     }
     size_t n_positions = 0;
-    AtomIDList::const_iterator a_i = getAtomIdBegin();
-    AtomIDList::const_iterator a_e = getAtomIdEnd();
-    PositionList::const_iterator p_i = getAtomPositionBegin();
+    AtomIDList::const_iterator a_i = globalmaster->getAtomIdBeginPublic();
+    AtomIDList::const_iterator a_e = globalmaster->getAtomIdEndPublic();
+    PositionList::const_iterator p_i = globalmaster->getAtomPositionBeginPublic();
 
     for ( ; a_i != a_e; ++a_i, ++p_i ) {
       atoms_positions[atoms_map[*a_i]] = cvm::rvector((*p_i).x, (*p_i).y, (*p_i).z);
@@ -481,9 +475,9 @@ void colvarproxy_namd::calculate()
         cvm::log("Updating total forces arrays.\n");
       }
       size_t n_total_forces = 0;
-      AtomIDList::const_iterator a_i = getForceIdBegin();
-      AtomIDList::const_iterator a_e = getForceIdEnd();
-      ForceList::const_iterator f_i = getTotalForce();
+      AtomIDList::const_iterator a_i = globalmaster->getForceIdBeginPublic();
+      AtomIDList::const_iterator a_e = globalmaster->getForceIdEndPublic();
+      ForceList::const_iterator f_i = globalmaster->getTotalForcePublic();
 
       for ( ; a_i != a_e; ++a_i, ++f_i ) {
         atoms_total_forces[atoms_map[*a_i]] = cvm::rvector((*f_i).x, (*f_i).y, (*f_i).z);
@@ -504,8 +498,8 @@ void colvarproxy_namd::calculate()
       if (cvm::debug()) {
         cvm::log("Updating group total forces arrays.\n");
       }
-      ForceList::const_iterator f_i = getGroupTotalForceBegin();
-      ForceList::const_iterator f_e = getGroupTotalForceEnd();
+      ForceList::const_iterator f_i = globalmaster->getGroupTotalForceBeginPublic();
+      ForceList::const_iterator f_e = globalmaster->getGroupTotalForceEndPublic();
       size_t i = 0;
       if ( (! b_simulation_continuing) &&
            ((f_e - f_i) != ((int) atom_groups_ids.size())) ) {
@@ -529,8 +523,8 @@ void colvarproxy_namd::calculate()
     size_t ig;
     // note: getGroupMassBegin() could be used here, but masses and charges
     // have already been calculated from the last call to setup()
-    PositionList::const_iterator gp_i = getGroupPositionBegin();
-    for (ig = 0; gp_i != getGroupPositionEnd(); gp_i++, ig++) {
+    PositionList::const_iterator gp_i = globalmaster->getGroupPositionBeginPublic();
+    for (ig = 0; gp_i != globalmaster->getGroupPositionEndPublic(); gp_i++, ig++) {
       atom_groups_coms[ig] = cvm::rvector(gp_i->x, gp_i->y, gp_i->z);
     }
   }
@@ -542,9 +536,12 @@ void colvarproxy_namd::calculate()
     }
     // Using a simple nested loop: there probably won't be so many maps that
     // this becomes performance-limiting
-    IntList::const_iterator goi_i = getGridObjIndexBegin();
-    BigRealList::const_iterator gov_i = getGridObjValueBegin();
-    for ( ; gov_i != getGridObjValueEnd(); goi_i++, gov_i++) {
+    IntList::const_iterator goi_i = globalmaster->getGridObjIndexBeginPublic();
+    BigRealList::const_iterator gov_i = globalmaster->getGridObjValueBeginPublic();
+    for ( ; gov_i != globalmaster->getGridObjValueEndPublic(); goi_i++, gov_i++) {
+      if (cvm::debug()) {
+        log("  ID = " + cvm::to_str(*goi_i) + ", value = " + cvm::to_str(*gov_i) + ".\n");
+      }
       for (size_t imap = 0; imap < volmaps_ids.size(); imap++) {
         if (volmaps_ids[imap] == *goi_i) {
           volmaps_values[imap] = *gov_i;
@@ -571,14 +568,14 @@ void colvarproxy_namd::calculate()
   // communicate all forces to the MD integrator
   for (size_t i = 0; i < atoms_ids.size(); i++) {
     cvm::rvector const &f = atoms_new_colvar_forces[i];
-    modifyForcedAtoms().add(atoms_ids[i]);
-    modifyAppliedForces().add(Vector(f.x, f.y, f.z));
+    globalmaster->modifyForcedAtomsPublic().add(atoms_ids[i]);
+    globalmaster->modifyAppliedForcesPublic().add(Vector(f.x, f.y, f.z));
   }
 
   if (atom_groups_new_colvar_forces.size() > 0) {
-    modifyGroupForces().resize(requestedGroups().size());
-    ForceList::iterator gf_i = modifyGroupForces().begin();
-    for (int ig = 0; gf_i != modifyGroupForces().end(); gf_i++, ig++) {
+    globalmaster->modifyGroupForcesPublic().resize(globalmaster->requestedGroups().size());
+    ForceList::iterator gf_i = globalmaster->modifyGroupForcesPublic().begin();
+    for (int ig = 0; gf_i != globalmaster->modifyGroupForcesPublic().end(); gf_i++, ig++) {
       cvm::rvector const &f = atom_groups_new_colvar_forces[ig];
       *gf_i = Vector(f.x, f.y, f.z);
     }
@@ -586,11 +583,11 @@ void colvarproxy_namd::calculate()
 
 #if NAMD_VERSION_NUMBER >= 34471681
   if (volmaps_new_colvar_forces.size() > 0) {
-    modifyGridObjForces().resize(requestedGridObjs().size());
-    modifyGridObjForces().setall(0.0);
-    IntList::const_iterator goi_i = getGridObjIndexBegin();
-    BigRealList::iterator gof_i = modifyGridObjForces().begin();
-    for ( ; goi_i != getGridObjIndexEnd(); goi_i++, gof_i++) {
+    globalmaster->modifyGridObjForcesPublic().resize(globalmaster->requestedGridObjs().size());
+    globalmaster->modifyGridObjForcesPublic().setall(0.0);
+    IntList::const_iterator goi_i = globalmaster->getGridObjIndexBeginPublic();
+    BigRealList::iterator gof_i = globalmaster->modifyGridObjForcesPublic().begin();
+    for ( ; goi_i != globalmaster->getGridObjIndexEndPublic(); goi_i++, gof_i++) {
       for (size_t imap = 0; imap < volmaps_ids.size(); imap++) {
         if (volmaps_ids[imap] == *goi_i) {
           *gof_i = volmaps_new_colvar_forces[imap];
@@ -621,11 +618,19 @@ void colvarproxy_namd::calculate()
   }
 }
 
+
+cvm::real colvarproxy_namd::rand_gaussian() { return random->gaussian(); }
+
+
 void colvarproxy_namd::update_accelMD_info() {
   // This aMD factor is from previous step!
   amd_weight_factor = std::exp(controller->accelMDdV / (target_temperature() * boltzmann()));
 }
 
+cvm::real colvarproxy_namd::get_accelMD_factor() const { return amd_weight_factor; }
+
+bool colvarproxy_namd::accelMD_enabled() const { return accelMDOn; }
+
 
 // Callback functions
 
@@ -674,7 +679,7 @@ void colvarproxy_namd::add_energy(cvm::real energy)
   #if !defined(NAMD_UNIFIED_REDUCTION)
   reduction->item(REDUCTION_MISC_ENERGY) += energy;
   #else
-  addReductionEnergy(REDUCTION_MISC_ENERGY, energy);
+  globalmaster->addReductionEnergyPublic(REDUCTION_MISC_ENERGY, energy);
   #endif
   #endif
 }
@@ -685,7 +690,7 @@ void colvarproxy_namd::request_total_force(bool yesno)
     cvm::log("colvarproxy_namd::request_total_force()\n");
   }
   total_force_requested = yesno;
-  requestTotalForce(total_force_requested);
+  globalmaster->requestTotalForcePublic(total_force_requested);
   if (cvm::debug()) {
     cvm::log("colvarproxy_namd::request_total_force() end\n");
   }
@@ -770,7 +775,7 @@ int colvarproxy_namd::init_atom(int atom_number)
 
   int const index = add_atom_slot(aid);
   atoms_map[aid] = index;
-  modifyRequestedAtoms().add(aid);
+  globalmaster->modifyRequestedAtomsPublic().add(aid);
   update_atom_properties(index);
   return index;
 }
@@ -832,7 +837,7 @@ int colvarproxy_namd::init_atom(cvm::residue_id const &residue,
 
   int const index = add_atom_slot(aid);
   atoms_map[aid] = index;
-  modifyRequestedAtoms().add(aid);
+  globalmaster->modifyRequestedAtomsPublic().add(aid);
   update_atom_properties(index);
   return index;
 }
@@ -868,7 +873,7 @@ cvm::rvector colvarproxy_namd::position_distance(cvm::atom_pos const &pos1,
   Position const p1(pos1.x, pos1.y, pos1.z);
   Position const p2(pos2.x, pos2.y, pos2.z);
   // return p2 - p1
-  Vector const d = this->lattice->delta(p2, p1);
+  Vector const d = globalmaster->get_lattice()->delta(p2, p1);
   return cvm::rvector(d.x, d.y, d.z);
 }
 
@@ -1210,8 +1215,8 @@ int colvarproxy_namd::init_atom_group(std::vector<int> const &atoms_ids)
 
   // compare this new group to those already allocated inside GlobalMaster
   int ig;
-  for (ig = 0; ig < modifyRequestedGroups().size(); ig++) {
-    AtomIDList const &namd_group = modifyRequestedGroups()[ig];
+  for (ig = 0; ig < globalmaster->modifyRequestedGroupsPublic().size(); ig++) {
+    AtomIDList const &namd_group = globalmaster->modifyRequestedGroupsPublic()[ig];
     bool b_match = true;
 
     if (namd_group.size() != ((int) atoms_ids.size())) {
@@ -1238,10 +1243,10 @@ int colvarproxy_namd::init_atom_group(std::vector<int> const &atoms_ids)
 
   // add this group (note: the argument of add_atom_group_slot() is redundant for NAMD, and provided only for consistency)
   size_t const index = add_atom_group_slot(atom_groups_ids.size());
-  modifyRequestedGroups().resize(atom_groups_ids.size());
+  globalmaster->modifyRequestedGroupsPublic().resize(atom_groups_ids.size());
   // the following is done in calculate()
-  // modifyGroupForces().resize(atom_groups_ids.size());
-  AtomIDList &namd_group = modifyRequestedGroups()[index];
+  // globalmaster->modifyGroupForces().resize(atom_groups_ids.size());
+  AtomIDList &namd_group = globalmaster->modifyRequestedGroupsPublic()[index];
   namd_group.resize(atoms_ids.size());
   int const n_all_atoms = Node::Object()->molecule->numAtoms;
   for (size_t ia = 0; ia < atoms_ids.size(); ia++) {
@@ -1261,8 +1266,8 @@ int colvarproxy_namd::init_atom_group(std::vector<int> const &atoms_ids)
 
   if (cvm::debug()) {
     cvm::log("Group has index "+cvm::to_str(index)+"\n");
-    cvm::log("modifyRequestedGroups length = "+cvm::to_str(modifyRequestedGroups().size())+
-        ", modifyGroupForces length = "+cvm::to_str(modifyGroupForces().size())+"\n");
+    cvm::log("modifyRequestedGroups length = "+cvm::to_str(globalmaster->modifyRequestedGroupsPublic().size())+
+        ", modifyGroupForces length = "+cvm::to_str(globalmaster->modifyGroupForcesPublic().size())+"\n");
   }
 
   return index;
@@ -1278,7 +1283,7 @@ void colvarproxy_namd::clear_atom_group(int index)
 
 int colvarproxy_namd::update_group_properties(int index)
 {
-  AtomIDList const &namd_group = modifyRequestedGroups()[index];
+  AtomIDList const &namd_group = globalmaster->modifyRequestedGroupsPublic()[index];
   if (cvm::debug()) {
     cvm::log("Re-calculating total mass and charge for scalable group no. "+cvm::to_str(index+1)+" ("+
              cvm::to_str(namd_group.size())+" atoms).\n");
@@ -1336,7 +1341,7 @@ int colvarproxy_namd::init_volmap_by_id(int volmap_id)
   int index = -1;
   if (error_code == COLVARS_OK) {
     index = add_volmap_slot(volmap_id);
-    modifyRequestedGridObjects().add(volmap_id);
+    globalmaster->modifyRequestedGridObjectsPublic().add(volmap_id);
   }
 
   return (error_code == COLVARS_OK) ? index : -1;
@@ -1378,7 +1383,7 @@ int colvarproxy_namd::init_volmap_by_name(char const *volmap_name)
     }
 
     index = add_volmap_slot(volmap_id);
-    modifyRequestedGridObjects().add(volmap_id);
+    globalmaster->modifyRequestedGridObjectsPublic().add(volmap_id);
   }
 
   return (error_code == COLVARS_OK) ? index : -1;
@@ -1677,7 +1682,7 @@ int colvarproxy_namd::request_alch_energy_freq(int const freq) {
 
 /// Get value of alchemical lambda parameter from back-end
 int colvarproxy_namd::get_alch_lambda(cvm::real* lambda) {
-  *lambda = simparams->getCurrentLambda(step);
+  *lambda = simparams->getCurrentLambda(globalmaster->step);
   return COLVARS_OK;
 }
 
diff --git a/namd/src/colvarproxy_namd.h b/namd/src/colvarproxy_namd.h
index 1de72d511..4e36e0fcf 100644
--- a/namd/src/colvarproxy_namd.h
+++ b/namd/src/colvarproxy_namd.h
@@ -15,29 +15,29 @@
 #define NAMD_VERSION_NUMBER 34471681
 #endif
 
-#include "Vector.h"
-#include "ResizeArray.h"
-#include "NamdTypes.h"
-#include "SimParameters.h"
-#include "Lattice.h"
-#include "GlobalMaster.h"
-#include "Random.h"
-#include "ConfigList.h"
+#include <memory>
+
+#include "colvarproxy_namd_version.h"
 
 #include "colvarmodule.h"
 #include "colvarproxy.h"
 #include "colvarvalue.h"
 
-#define GLOBAL_MASTER_CKLOOP_CALC_ITEM 2000
-#define GLOBAL_MASTER_CKLOOP_CALC_BIASES 2001
-#define GLOBAL_MASTER_CKLOOP_CALC_SCRIPTED_BIASES 2002
 
-/// \brief Communication between colvars and NAMD (implementation of
-/// \link colvarproxy \endlink)
-class colvarproxy_namd : public colvarproxy, public GlobalMaster {
+class Controller;
+class GlobalMasterColvars;
+class Random;
+class SimParameters;
+
+
+/// Communication between colvars and NAMD (implementation of \link colvarproxy \endlink)
+class colvarproxy_namd : public colvarproxy {
 
 protected:
 
+  /// Pointer to the parent GlobalMaster object
+  GlobalMasterColvars *globalmaster = nullptr;
+
   /// \brief Array of atom indices (relative to the colvarproxy arrays),
   /// usedfor faster copy of atomic data
   std::vector<int> atoms_map;
@@ -49,7 +49,7 @@ class colvarproxy_namd : public colvarproxy, public GlobalMaster {
   Controller const *controller;
 
   /// NAMD-style PRNG object
-  Random random;
+  std::unique_ptr<Random> random;
 
   bool first_timestep;
   cvm::step_number previous_NAMD_step;
@@ -71,7 +71,7 @@ class colvarproxy_namd : public colvarproxy, public GlobalMaster {
 
   void init_tcl_pointers() override;
 
-  colvarproxy_namd();
+  colvarproxy_namd(GlobalMasterColvars *gm);
   ~colvarproxy_namd();
 
   int setup() override;
@@ -108,18 +108,11 @@ class colvarproxy_namd : public colvarproxy, public GlobalMaster {
                                    std::vector<const colvarvalue *> const &cvcs,
                                    std::vector<cvm::matrix2d<cvm::real> > &gradient) override;
 
-  cvm::real rand_gaussian() override
-  {
-    return random.gaussian();
-  }
+  cvm::real rand_gaussian() override;
 
-  cvm::real get_accelMD_factor() const override {
-    return amd_weight_factor;
-  }
+  cvm::real get_accelMD_factor() const override;
 
-  bool accelMD_enabled() const override {
-    return accelMDOn;
-  }
+  bool accelMD_enabled() const override;
 
 #if CMK_SMP && USE_CKLOOP
   colvarproxy::smp_mode_t get_smp_mode() const override;
diff --git a/update-colvars-code.sh b/update-colvars-code.sh
index ddbf9af86..a333e26e4 100755
--- a/update-colvars-code.sh
+++ b/update-colvars-code.sh
@@ -45,12 +45,16 @@ then
   shift
 fi
 
+# Assume forward patching
+PATCH_OPTS="-p1 --forward -s"
+
 # Undocumented flag
 reverse=0
 if [ $1 = "-R" ]
 then
   echo "Reverse: updating git tree from downstream tree"
   reverse=1
+  PATCH_OPTS="-p1 --reverse -s"
   shift
 fi
 
@@ -211,11 +215,9 @@ condcopy() {
   then
     a=$2
     b=$1
-    PATCH_OPT="-R"
   else
     a=$1
     b=$2
-    PATCH_OPT=""
   fi
 
   updated_file=0
@@ -261,6 +263,7 @@ checkfile() {
 }
 
 
+
 # Update LAMMPS tree
 if [ ${code} = "LAMMPS" ]
 then
@@ -284,11 +287,9 @@ then
     condcopy "${src}" "${target}/src/COLVARS/${tgt}"
   done
 
-  if [ -f ${source}/lammps/COLVARS.cmake.diff ] ; then
+  if [ -f ${source}/lammps/COLVARS.cmake.patch ] ; then
     # Do not exit if the patch fails - already applied in development branch
-    set +e
-    patch -p1 -N -d ${target} < ${source}/lammps/COLVARS.cmake.diff
-    set -e
+    patch ${PATCH_OPTS} -d ${target} < ${source}/lammps/COLVARS.cmake.patch || true
   fi
 
   downloaded_pdf=0
@@ -346,9 +347,6 @@ then
              "${target}/lepton/Make.depends"
     condcopy "${source}/namd/lepton/Makefile.namd" \
              "${target}/lepton/Makefile.namd"
-    if ! grep -q lepton/Makefile.namd "${target}/Makefile" ; then
-      patch -p1 -N -d ${target} < namd/Makefile.patch
-    fi
   fi
 
   # Copy library files to the "colvars" folder
@@ -367,13 +365,17 @@ then
 
   # Update NAMD interface files
   for src in \
-      ${source}/namd/src/colvarproxy_namd.h \
-      ${source}/namd/src/colvarproxy_namd_version.h \
-      ${source}/namd/src/colvarproxy_namd.C
+      ${source}/namd/src/*.h \
+      ${source}/namd/src/*.C
   do \
     tgt=$(basename ${src})
     condcopy "${src}" "${target}/src/${tgt}"
   done
+  for patch in ${source}/namd/*.patch ${source}/namd/src/*.patch ; do
+    if [ -s ${patch} ] ; then
+      patch ${PATCH_OPTS} -d "${target}" < "${patch}" || true
+    fi
+  done
 
   # Update abf_integrate
   for src in ${source}/colvartools/*h ${source}/colvartools/*cpp
@@ -438,11 +440,9 @@ then
     condcopy "${src}" "${target}/src/${tgt}.C"
   done
 
-  if [ -f ${source}/vmd/src/tcl_commands.C.diff ] ; then
+  if [ -f ${source}/vmd/src/tcl_commands.C.patch ] ; then
     # Do not exit if the patch fails - already applied in development branch
-    set +e
-    patch -p1 -N -d ${target} < ${source}/vmd/src/tcl_commands.C.diff
-    set -e
+    patch ${PATCH_OPTS} -d ${target} < ${source}/vmd/src/tcl_commands.C.patch || true
   fi
 
   # Update replacement text for the Colvars manual
@@ -495,7 +495,6 @@ if [ ${code} = "GROMACS" ]
 then
 
   target_folder=${target}/src/external/colvars
-  patch_opts="-p1 --forward -s"
 
   mkdir -p ${target_folder}
 
@@ -508,21 +507,21 @@ then
   echo ""
 
   # Patch CMake build recipes when applicable
-  if [ -s ${source}/gromacs/gmxManageColvars.cmake.diff ] ; then
-    patch ${patch_opts} -d ${target} < ${source}/gromacs/gmxManageColvars.cmake.diff || true
+  if [ -s ${source}/gromacs/gmxManageColvars.cmake.patch ] ; then
+    patch ${PATCH_OPTS} -d ${target} < ${source}/gromacs/gmxManageColvars.cmake.patch || true
   fi
-  if [ -s ${source}/gromacs/CMakeLists.txt.diff ] ; then
-    patch ${patch_opts} -d ${target} < ${source}/gromacs/CMakeLists.txt.diff || true
+  if [ -s ${source}/gromacs/CMakeLists.txt.patch ] ; then
+    patch ${PATCH_OPTS} -d ${target} < ${source}/gromacs/CMakeLists.txt.patch || true
   fi
-  if [ -s ${source}/gromacs/CMakeLists.txt.diff ] ; then
-    patch ${patch_opts} -d ${target} < ${source}/gromacs/CMakeLists.txt.diff || true
+  if [ -s ${source}/gromacs/CMakeLists.txt.patch ] ; then
+    patch ${PATCH_OPTS} -d ${target} < ${source}/gromacs/CMakeLists.txt.patch || true
   fi
 
   if [ "x${UPDATE_LEPTON}" == "xyes" ] ; then
     echo -n "(note: adding/updating Lepton)"
     copy_lepton ${target}/src/external/ || exit 1
-    if [ -s ${source}/gromacs/CMakeLists.txt.Lepton.diff ] ; then
-      patch ${patch_opts} -d ${target} < ${source}/gromacs/CMakeLists.txt.Lepton.diff || true
+    if [ -s ${source}/gromacs/CMakeLists.txt.Lepton.patch ] ; then
+      patch ${PATCH_OPTS} -d ${target} < ${source}/gromacs/CMakeLists.txt.Lepton.patch || true
     fi
     condcopy ${source}/gromacs/cmake/gmxManageLepton.cmake "${target}/cmake/gmxManageLepton.cmake"
   fi
@@ -532,12 +531,10 @@ then
   # Copy MDModules files to the "src/gromacs/src/applied_forces/colvars" folder
   target_folder=${target}/src/gromacs/applied_forces/colvars
 
-  for patch_file in ${source}/gromacs/src/applied_forces/colvars/*.diff ; do
+  for patch_file in ${source}/gromacs/src/applied_forces/colvars/*.patch ; do
     if [ -s ${patch_file} ] ; then
       # Do not exit if the patch fails - already applied
-      set +e
-      patch -p1 -N -d ${target} < ${patch_file}
-      set -e
+      patch ${PATCH_OPTS} -d ${target} < ${patch_file} || true
     fi
   done
 
@@ -567,7 +564,7 @@ then
 
   # Apply patch for Gromacs files
   if [ -s ${source}/gromacs/gromacs-mdmodules.patch ] ; then
-    patch ${patch_opts} -d ${target} < ${source}/gromacs/gromacs-mdmodules.patch || true
+    patch ${PATCH_OPTS} -d ${target} < ${source}/gromacs/gromacs-mdmodules.patch || true
   fi
   ret_val=$?
   if [ $ret_val -ne 0 ]
diff --git a/vmd/src/tcl_commands.C.diff b/vmd/src/tcl_commands.C.patch
similarity index 100%
rename from vmd/src/tcl_commands.C.diff
rename to vmd/src/tcl_commands.C.patch