diff --git a/patches/gromacs-4.5.3.diff b/patches/gromacs-4.5.3.diff
index f5da3ceb6d8751301890f7d63c7652c4d1ee4268..73f1014d0cac3162566b7d92bdcd48052ab7a696 100644
--- a/patches/gromacs-4.5.3.diff
+++ b/patches/gromacs-4.5.3.diff
@@ -1,3 +1,83 @@
+patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/repl_ex.c" << \EOF_EOF
+--- ./src/kernel/repl_ex.c.preplumed
++++ ./src/kernel/repl_ex.c
+@@ -49,10 +49,16 @@
+ #include "names.h"
+ #include "mvdata.h"
+ #include "domdec.h"
+ #include "partdec.h"
+
++/* PLUMED */
++#include "../../Plumed.h"
++extern int plumedswitch;
++extern plumed plumedmain;
++/* END PLUMED */
++
+ typedef struct gmx_repl_ex {
+ int repl;
+ int nrepl;
+ real temp;
+ int type;
+@@ -538,10 +544,25 @@
+ snew(bEx,re->nrepl);
+ snew(prob,re->nrepl);
+
+ exchange = -1;
+ m = (step / re->nst) % 2;
++
++ if(plumedswitch){
++ int partner=re->repl;
++
++ for(i=1; i<re->nrepl; i++) {
++ if (i % 2 != m) continue;
++ a = re->ind[i-1];
++ b = re->ind[i];
++ if(re->repl==a) partner=b;
++ if(re->repl==b) partner=a;
++ }
++ plumed_cmd(plumedmain,"GREX setPartner",&partner);
++ plumed_cmd(plumedmain,"GREX calculate",NULL);
++ }
++
+ for(i=1; i<re->nrepl; i++) {
+ a = re->ind[i-1];
+ b = re->ind[i];
+ bPrint = (re->repl==a || re->repl==b);
+ if (i % 2 == m) {
+@@ -572,10 +593,19 @@
+ dpV = (betaA*re->pres[a]-betaB*re->pres[b])*(Vol[b]-Vol[a])/PRESFAC;
+ if (bPrint)
+ fprintf(fplog," dpV = %10.3e d = %10.3e",dpV,delta + dpV);
+ delta += dpV;
+ }
++ if(plumedswitch){
++ real ldb,fdb,dplumed;
++ plumed_cmd(plumedmain,"GREX getLocalDeltaBias",&ldb);
++ plumed_cmd(plumedmain,"GREX getForeignDeltaBias",&fdb);
++ dplumed=ldb*betaA+fdb*betaB;
++ delta+=dplumed;
++ if (bPrint)
++ fprintf(fplog," dplumed = %10.3e d = %10.3e",dplumed,delta);
++ }
+ if (bPrint)
+ fprintf(fplog,"\n");
+ if (delta <= 0) {
+ prob[i] = 1;
+ bEx[i] = TRUE;
+@@ -633,10 +663,12 @@
+ gmx_multisim_t *ms;
+ int exchange=-1,shift;
+ gmx_bool bExchanged=FALSE;
+
+ ms = cr->ms;
++
++ if(plumedswitch)plumed_cmd(plumedmain,"GREX prepare",NULL);
+
+ if (MASTER(cr))
+ {
+ exchange = get_replica_exchange(fplog,ms,re,ener,det(state->box),
+ step,time);
+EOF_EOF
patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/Makefile" << \EOF_EOF
--- ./src/kernel/Makefile.preplumed
+++ ./src/kernel/Makefile
@@ -47,7 +127,46 @@ patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/md.c" << \EOF_EOF
#endif
#ifdef GMX_THREADS
#include "tmpi.h"
-@@ -1522,10 +1528,41 @@
+@@ -1375,13 +1381,37 @@
+ for(i=0; (i<gnx); i++) {
+ grpindex[i] = i;
+ }
+ }
+
+- if (repl_ex_nst > 0 && MASTER(cr))
++ if (repl_ex_nst > 0 && MASTER(cr)){
+ repl_ex = init_replica_exchange(fplog,cr->ms,state_global,ir,
+ repl_ex_nst,repl_ex_seed);
++ if(plumedswitch) {
++ plumed_cmd(plumedmain,"GREX setMPIIntercomm",&cr->ms->mpi_comm_masters);
++ if(PAR(cr)){
++ if(DOMAINDECOMP(cr)) {
++ plumed_cmd(plumedmain,"GREX setMPIIntracomm",&cr->dd->mpi_comm_all);
++ }else{
++ plumed_cmd(plumedmain,"GREX setMPIIntracomm",&cr->mpi_comm_mysim);
++ }
++ }
++ plumed_cmd(plumedmain,"GREX init",NULL);
++ }
++ }
++
++ if (repl_ex_nst>0 && plumedswitch && ! MASTER(cr)){
++ if(PAR(cr)){
++ if(DOMAINDECOMP(cr)) {
++ plumed_cmd(plumedmain,"GREX setMPIIntracomm",&cr->dd->mpi_comm_all);
++ }else{
++ plumed_cmd(plumedmain,"GREX setMPIIntracomm",&cr->mpi_comm_mysim);
++ }
++ }
++ plumed_cmd(plumedmain,"GREX init",NULL);
++ }
++
+
+ if (!ir->bContinuation && !bRerunMD)
+ {
+ if (mdatoms->cFREEZE && (state->flags & (1<<estV)))
+ {
+@@ -1522,10 +1552,41 @@
}
}
fprintf(fplog,"\n");
@@ -89,7 +208,7 @@ patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/md.c" << \EOF_EOF
print_date_and_time(fplog,cr->nodeid,"Started mdrun",runtime);
wallcycle_start(wcycle,ewcRUN);
if (fplog)
-@@ -1802,10 +1839,17 @@
+@@ -1802,10 +1863,17 @@
state,&f,mdatoms,top,fr,
vsite,shellfc,constr,
nrnb,wcycle,do_verbose);
@@ -107,7 +226,7 @@ patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/md.c" << \EOF_EOF
if (MASTER(cr) && do_log && !bFFscan)
{
-@@ -1928,16 +1972,37 @@
+@@ -1928,16 +1996,40 @@
* in do_force.
* This is parallellized as well, and does communication too.
* Check comments in sim_util.c
@@ -121,6 +240,7 @@ patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/md.c" << \EOF_EOF
+ plumed_cmd(plumedmain,"setMasses",&mdatoms->massT[mdatoms->start]);
+ plumed_cmd(plumedmain,"setCharges",&mdatoms->chargeA[mdatoms->start]);
+ plumed_cmd(plumedmain,"setBox",&state->box[0][0]);
++
+ plumed_cmd(plumedmain,"prepareCalc",NULL);
+ }
+ /* END PLUMED */
@@ -137,6 +257,8 @@ patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/md.c" << \EOF_EOF
+ plumed_cmd(plumedmain,"setForces",&f[mdatoms->start][0]);
+ plumed_cmd(plumedmain,"setVirial",&force_vir[0][0]);
+ plumed_cmd(plumedmain,"performCalc",NULL);
++ if ((repl_ex_nst > 0) && (step > 0) && !bLastStep &&
++ do_per_step(step,repl_ex_nst)) plumed_cmd(plumedmain,"GREX savePositions",NULL);
+ }
+ /* END PLUMED */
+
diff --git a/src/GREX.cpp b/src/GREX.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..9536ef7b7f1615f5be8c8385fc7a028067ecc33d
--- /dev/null
+++ b/src/GREX.cpp
@@ -0,0 +1,105 @@
+#include "GREX.h"
+#include "sstream"
+#include "PlumedMain.h"
+
+using namespace std;
+using namespace PLMD;
+
+GREX::GREX(PlumedMain&p):
+ initialized(false),
+ plumedMain(p),
+ atoms(p.getAtoms()),
+ partner(-1), // = unset
+ myreplica(-1) // = unset
+{
+ p.setSuffix(".NA");
+}
+
+GREX::~GREX(){
+}
+
+void GREX::cmd(const string&key,void*val){
+ if(false){
+ }else if(key=="initialized"){
+ *static_cast<int*>(val)=initialized;
+ }else if(key=="setMPIIntracomm"){
+ assert(!initialized);
+ intracomm.Set_comm(val);
+ }else if(key=="setMPIIntercomm"){
+ assert(!initialized);
+ intercomm.Set_comm(val);
+ }else if(key=="init"){
+ assert(!initialized);
+ initialized=true;
+ std::string s;
+// note that for PEs!=root this is automatically 0 (comm defaults to MPI_COMM_SELF)
+ myreplica=intercomm.Get_rank();
+ intracomm.Sum(&myreplica,1);
+ Tools::convert(myreplica,s);
+ plumedMain.setSuffix("."+s);
+ }else if(key=="prepare"){
+ assert(initialized);
+ if(intracomm.Get_rank()==0) return;
+ intracomm.Bcast(&partner,1,0);
+ calculate();
+ }else if(key=="setPartner"){
+ assert(initialized);
+ partner=*static_cast<int*>(val);
+ }else if(key=="savePositions"){
+ assert(initialized);
+ savePositions();
+ }else if(key=="calculate"){
+ assert(initialized);
+ if(intracomm.Get_rank()!=0) return;
+ intracomm.Bcast(&partner,1,0);
+ calculate();
+ }else if(key=="getLocalDeltaBias"){
+ assert(initialized);
+ double x=localDeltaBias/(atoms.getMDUnits().energy/atoms.getUnits().energy);
+ atoms.double2MD(x,val);
+ }else if(key=="getForeignDeltaBias"){
+ assert(initialized);
+ double x=foreignDeltaBias/(atoms.getMDUnits().energy/atoms.getUnits().energy);
+ atoms.double2MD(x,val);
+ } else {
+ fprintf(stderr,"+++ PLUMED GREX ERROR\n");
+ fprintf(stderr,"+++ CANNOT INTERPRET CALL TO cmd() ROUTINE WITH ARG %s\n",key.c_str());
+ fprintf(stderr,"+++ There might be a mistake in the MD code\n");
+ fprintf(stderr,"+++ or you may be using an out-dated plumed version\n");
+ assert(0);
+ }
+}
+
+void GREX::savePositions(){
+ plumedMain.prepareDependencies();
+ atoms.shareAll();
+ plumedMain.waitData();
+ ostringstream o;
+ atoms.writeBinary(o);
+ buffer=o.str();
+}
+
+void GREX::calculate(){
+//fprintf(stderr,"CALCULATE %d %d\n",intercomm.Get_rank(),partner);
+ unsigned nn=buffer.size();
+ vector<char> rbuf(nn);
+ localDeltaBias=-plumedMain.getBias();
+ if(intracomm.Get_rank()==0){
+ PlumedCommunicator::Request req=intercomm.Isend(&buffer.c_str()[0],nn,partner,1066);
+ intercomm.Recv(&rbuf[0],rbuf.size(),partner,1066);
+ req.wait();
+ }
+ intracomm.Bcast(&rbuf[0],nn,0);
+ istringstream i(string(&rbuf[0],rbuf.size()));
+ atoms.readBinary(i);
+ plumedMain.prepareDependencies();
+ plumedMain.justCalculate();
+ localDeltaBias+=plumedMain.getBias();
+ if(intracomm.Get_rank()==0){
+ PlumedCommunicator::Request req=intercomm.Isend(&localDeltaBias,1,partner,1067);
+ intercomm.Recv(&foreignDeltaBias,1,partner,1067);
+ req.wait();
+//fprintf(stderr,">>> %d %d %20.12f %20.12f\n",intercomm.Get_rank(),partner,localDeltaBias,foreignDeltaBias);
+ }
+ intracomm.Bcast(&foreignDeltaBias,1,0);
+}
diff --git a/src/GREX.h b/src/GREX.h
new file mode 100644
index 0000000000000000000000000000000000000000..a240029be39272d8434178f88d8dbc1a6dc13506
--- /dev/null
+++ b/src/GREX.h
@@ -0,0 +1,36 @@
+#ifndef __PLUMED_Grex_h
+#define __PLUMED_Grex_h
+
+#include "WithCmd.h"
+#include "PlumedCommunicator.h"
+#include <string>
+
+namespace PLMD{
+
+class PlumedMain;
+class Atoms;
+
+class GREX:
+ public WithCmd
+{
+ bool initialized;
+ PlumedCommunicator intracomm;
+ PlumedCommunicator intercomm;
+ PlumedMain& plumedMain;
+ Atoms& atoms;
+ int partner;
+ double localDeltaBias;
+ double foreignDeltaBias;
+ std::string buffer;
+ int myreplica;
+public:
+ GREX(PlumedMain&);
+ ~GREX();
+ void cmd(const std::string&key,void*val=NULL);
+ void calculate();
+ void savePositions();
+};
+
+}
+
+#endif
diff --git a/src/PlumedMain.cpp b/src/PlumedMain.cpp
index 2fa4350c190fe1cd20af5666635c062f4a934177..c944c753cf199e4031667e75cfc95ff6cc2b6f03 100644
--- a/src/PlumedMain.cpp
+++ b/src/PlumedMain.cpp
@@ -14,6 +14,8 @@
#include "ActionRegister.h"
+#include "GREX.h"
+
using namespace PLMD;
using namespace std;
@@ -265,6 +267,7 @@ void PlumedMain::cmd(const std::string & word,void*val){
if(actionRegister().check(words[1])) check=1;
*(static_cast<int*>(val))=check;
} else if(nw==2 && words[0]=="GREX"){
+ if(!grex) grex=new GREX(*this);
assert(grex);
grex->cmd(words[1],val);
} else{