diff --git a/patches/gromacs-4.6.3-hrex.config b/patches/gromacs-4.6.3-hrex.config deleted file mode 100644 index 6346748e72bdbd219da14fcd9a05d883ad7c9b5a..0000000000000000000000000000000000000000 --- a/patches/gromacs-4.6.3-hrex.config +++ /dev/null @@ -1,29 +0,0 @@ - - -function plumed_preliminary_test(){ -# check if the README contains the word GROMACS and if gromacs has been already configured - grep -q GROMACS README 1>/dev/null 2>/dev/null -} - -plumed_before_patch(){ -cat << EOF -This is an experimental patch to perform Hamiltonian replica exchange with GROMACS -Implementation described in: - Bussi, Hamiltonian replica-exchange in GROMACS: a flexible implementation, - Mol. Phys., DOI: 10.1080/00268976.2013.824126 (2013) - -To use it: -* Prepare separate topologies (topol0.tpr, topol1.tpr, etc) -* Set environment variable PLUMED_HREX (export PLUMED_HREX=1) -* Run a normal replica exchange with gromacs - -Suggested checks: -* Try with several identical force fields and different seed/starting point - Acceptance should be 1.0 - -Warnings: -* Topologies should have the same number of atoms, same masses and same constraint topology -* Choose neighbor list update that divides replex - -EOF -} diff --git a/patches/gromacs-4.6.3-hrex.diff b/patches/gromacs-4.6.3-hrex.diff deleted file mode 100644 index 529123143e0fed5726f78f0025e198ec157425c7..0000000000000000000000000000000000000000 --- a/patches/gromacs-4.6.3-hrex.diff +++ /dev/null @@ -1,585 +0,0 @@ -patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/CMakeLists.txt" << \EOF_EOF ---- ./src/kernel/CMakeLists.txt.preplumed -+++ ./src/kernel/CMakeLists.txt -@@ -31,10 +31,12 @@ - # - # To help us fund GROMACS development, we humbly ask that you cite - # the research papers on the package. Check out http://www.gromacs.org. - # - -+include(${CMAKE_SOURCE_DIR}/Plumed.cmake) -+ - set(GMXPREPROCESS_SOURCES - add_par.c - calc_verletbuf.c - compute_io.c - convparm.c -@@ -121,11 +123,11 @@ - set_target_properties(${PROGRAM} PROPERTIES OUTPUT_NAME "${PROGRAM}${GMX_BINARY_SUFFIX}") - endforeach() - - add_executable(mdrun ${MDRUN_SOURCES} main.c) - gmx_add_man_page(mdrun) --target_link_libraries(mdrun gmxpreprocess md gmx ${OpenMP_LINKER_FLAGS}) -+target_link_libraries(mdrun gmxpreprocess md gmx ${OpenMP_LINKER_FLAGS} ${PLUMED_LOAD}) - set_target_properties(mdrun PROPERTIES OUTPUT_NAME "mdrun${GMX_BINARY_SUFFIX}" COMPILE_FLAGS "${OpenMP_C_FLAGS}") - - if(GMX_OPENMM) - target_link_libraries(mdrun openmm_api_wrapper) - endif() -EOF_EOF -patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/md.c" << \EOF_EOF ---- ./src/kernel/md.c.preplumed -+++ ./src/kernel/md.c -@@ -91,10 +91,16 @@ - #include "membed.h" - #include "types/nlistheuristics.h" - #include "types/iteratedconstraints.h" - #include "nbnxn_cuda_data_mgmt.h" - -+/* PLUMED */ -+#include "../../Plumed.h" -+extern int plumedswitch; -+extern plumed plumedmain; -+/* END PLUMED */ -+ - #ifdef GMX_LIB_MPI - #include <mpi.h> - #endif - #ifdef GMX_THREAD_MPI - #include "tmpi.h" -@@ -235,10 +241,14 @@ - /* PME load balancing data for GPU kernels */ - pme_load_balancing_t pme_loadbal = NULL; - double cycles_pmes; - gmx_bool bPMETuneTry = FALSE, bPMETuneRunning = FALSE; - -+/* PLUMED */ -+ int plumedNeedsEnergy; -+/* END PLUMED */ -+ - #ifdef GMX_FAHCORE - /* Temporary addition for FAHCORE checkpointing */ - int chkpt_ret; - #endif - -@@ -715,10 +725,50 @@ - } - } - fprintf(fplog, "\n"); - } - -+ /* PLUMED */ -+ if(plumedswitch){ -+ if(cr->ms && cr->ms->nsim>1) { -+ if(MASTER(cr)) 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(PAR(cr)){ -+ if(DOMAINDECOMP(cr)) { -+ plumed_cmd(plumedmain,"setMPIComm",&cr->dd->mpi_comm_all); -+ }else{ -+ plumed_cmd(plumedmain,"setMPIComm",&cr->mpi_comm_mysim); -+ } -+ } -+ plumed_cmd(plumedmain,"setNatoms",&top_global->natoms); -+ plumed_cmd(plumedmain,"setMDEngine","gromacs"); -+ plumed_cmd(plumedmain,"setLog",fplog); -+ real real_delta_t; -+ real_delta_t=ir->delta_t; -+ plumed_cmd(plumedmain,"setTimestep",&real_delta_t); -+ plumed_cmd(plumedmain,"init",NULL); -+ -+ if(PAR(cr)){ -+ if(DOMAINDECOMP(cr)) { -+ plumed_cmd(plumedmain,"setAtomsNlocal",&cr->dd->nat_home); -+ plumed_cmd(plumedmain,"setAtomsGatindex",cr->dd->gatindex); -+ }else{ -+ plumed_cmd(plumedmain,"setAtomsNlocal",&mdatoms->homenr); -+ plumed_cmd(plumedmain,"setAtomsContiguous",&mdatoms->start); -+ } -+ } -+ } -+ /* END PLUMED */ -+ - /* Set and write start time */ - runtime_start(runtime); - print_date_and_time(fplog, cr->nodeid, "Started mdrun", runtime); - wallcycle_start(wcycle, ewcRUN); - if (fplog) -@@ -1031,10 +1081,17 @@ - vsite, shellfc, constr, - nrnb, wcycle, - do_verbose && !bPMETuneRunning); - wallcycle_stop(wcycle, ewcDOMDEC); - /* If using an iterative integrator, reallocate space to match the decomposition */ -+ -+ /* PLUMED */ -+ if(plumedswitch){ -+ plumed_cmd(plumedmain,"setAtomsNlocal",&cr->dd->nat_home); -+ plumed_cmd(plumedmain,"setAtomsGatindex",cr->dd->gatindex); -+ } -+ /* END PLUMED */ - } - } - - if (MASTER(cr) && do_log && !bFFscan) - { -@@ -1080,10 +1137,96 @@ - } - } - - GMX_MPE_LOG(ev_timestep2); - -+ gmx_bool bHREX; -+ bHREX= repl_ex_nst > 0 && (step>0) && !bLastStep && do_per_step(step,repl_ex_nst) -+ && getenv("PLUMED_HREX"); -+ -+/* Hamiltonian Replica Exchange */ -+ if(plumedswitch) if(bHREX){ -+ gmx_enerdata_t *hrex_enerd; -+ snew(hrex_enerd,1); -+ init_enerdata(top_global->groups.grps[egcENER].nr,ir->fepvals->n_lambda,hrex_enerd); -+ int repl=-1; -+ int nrepl=-1; -+ if(MASTER(cr)){ -+ repl=replica_exchange_get_repl(repl_ex); -+ nrepl=replica_exchange_get_nrepl(repl_ex); -+ } -+ if (PAR(cr)) { -+ if (DOMAINDECOMP(cr)) -+ dd_collect_state(cr->dd,state,state_global); -+ else -+ pd_collect_state(cr,state_global); -+ } -+ if(MASTER(cr)){ -+ if(repl%2==step/repl_ex_nst%2){ -+ if(repl-1>=0) exchange_state(cr->ms,repl-1,state_global); -+ }else{ -+ if(repl+1<nrepl) exchange_state(cr->ms,repl+1,state_global); -+ } -+ } -+ if(PAR(cr)){ -+ if (DOMAINDECOMP(cr)) { -+ dd_partition_system(fplog,step,cr,TRUE,1, -+ state_global,top_global,ir, -+ state,&f,mdatoms,top,fr,vsite,shellfc,constr, -+ nrnb,wcycle,FALSE); -+ } else { -+ bcast_state(cr,state,FALSE); -+ } -+ } -+ do_force(fplog, cr, ir, step, nrnb, wcycle, top, top_global, groups, -+ state->box, state->x, &state->hist, -+ f, force_vir, mdatoms, hrex_enerd, fcd, -+ state->lambda, graph, -+ fr, vsite, mu_tot, t, outf->fp_field, ed, bBornRadii, -+ GMX_FORCE_STATECHANGED | -+ ((DYNAMIC_BOX(*ir) || bRerunMD) ? GMX_FORCE_DYNAMICBOX : 0) | -+ GMX_FORCE_ALLFORCES | -+ GMX_FORCE_SEPLRF | -+ GMX_FORCE_VIRIAL | -+ GMX_FORCE_ENERGY | -+ GMX_FORCE_DHDL | -+ GMX_FORCE_NS); -+ plumed_cmd(plumedmain,"GREX cacheLocalUSwap",&hrex_enerd->term[F_EPOT]); -+ sfree(hrex_enerd); -+ -+/* exchange back */ -+ if (PAR(cr)) { -+ if (DOMAINDECOMP(cr)) -+ dd_collect_state(cr->dd,state,state_global); -+ else -+ pd_collect_state(cr,state_global); -+ } -+ if(MASTER(cr)){ -+ if(repl%2==step/repl_ex_nst%2){ -+ if(repl-1>=0) exchange_state(cr->ms,repl-1,state_global); -+ }else{ -+ if(repl+1<nrepl) exchange_state(cr->ms,repl+1,state_global); -+ } -+ } -+ if(PAR(cr)){ -+ if (DOMAINDECOMP(cr)) { -+ dd_partition_system(fplog,step,cr,TRUE,1, -+ state_global,top_global,ir, -+ state,&f,mdatoms,top,fr,vsite,shellfc,constr, -+ nrnb,wcycle,FALSE); -+ if(plumedswitch){ -+ plumed_cmd(plumedmain,"setAtomsNlocal",&cr->dd->nat_home); -+ plumed_cmd(plumedmain,"setAtomsGatindex",cr->dd->gatindex); -+ } -+ } else { -+ bcast_state(cr,state,FALSE); -+ } -+ } -+ -+ } -+/* END Hamiltonian Replica Exchange */ -+ - /* We write a checkpoint at this MD step when: - * either at an NS step when we signalled through gs, - * or at the last step (but not when we do not want confout), - * but never at the first step or with rerun. - */ -@@ -1176,16 +1319,43 @@ - /* The coordinates (x) are shifted (to get whole molecules) - * in do_force. - * This is parallellized as well, and does communication too. - * Check comments in sim_util.c - */ -+ -+ /* PLUMED */ -+ if(plumedswitch){ -+ plumedNeedsEnergy=0; -+ long int lstep=step; plumed_cmd(plumedmain,"setStepLong",&lstep); -+ plumed_cmd(plumedmain,"setPositions",&state->x[mdatoms->start][0]); -+ 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); -+ plumed_cmd(plumedmain,"isEnergyNeeded",&plumedNeedsEnergy); -+ } -+ /* END PLUMED */ - do_force(fplog, cr, ir, step, nrnb, wcycle, top, top_global, groups, - state->box, state->x, &state->hist, - f, force_vir, mdatoms, enerd, fcd, - state->lambda, graph, - fr, vsite, mu_tot, t, outf->fp_field, ed, bBornRadii, -- (bNS ? GMX_FORCE_NS : 0) | force_flags); -+ (plumedNeedsEnergy? GMX_FORCE_ENERGY : 0) | (bNS ? GMX_FORCE_NS : 0) | force_flags); -+ /* PLUMED */ -+ if(plumedswitch){ -+ if(plumedNeedsEnergy) plumed_cmd(plumedmain,"setEnergy",&enerd->term[F_EPOT]); -+ 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); -+ -+ if(bHREX) -+ plumed_cmd(plumedmain,"GREX cacheLocalUNow",&enerd->term[F_EPOT]); -+ } -+ /* END PLUMED */ -+ - } - - GMX_BARRIER(cr->mpi_comm_mygroup); - - if (bTCR) -EOF_EOF -patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/mdrun.c" << \EOF_EOF ---- ./src/kernel/mdrun.c.preplumed -+++ ./src/kernel/mdrun.c -@@ -56,10 +56,16 @@ - #endif - - /* afm stuf */ - #include "pull.h" - -+/* PLUMED */ -+#include "../../Plumed.h" -+int plumedswitch; -+plumed plumedmain; -+/* END PLUMED */ -+ - int cmain(int argc, char *argv[]) - { - const char *desc[] = { - "The [TT]mdrun[tt] program is the main computational chemistry engine", - "within GROMACS. Obviously, it performs Molecular Dynamics simulations,", -@@ -399,10 +405,11 @@ - { efLOG, "-rs", "rotslabs", ffOPTWR }, - { efLOG, "-rt", "rottorque", ffOPTWR }, - { efMTX, "-mtx", "nm", ffOPTWR }, - { efNDX, "-dn", "dipole", ffOPTWR }, - { efRND, "-multidir", NULL, ffOPTRDMULT}, -+ { efDAT, "-plumed", "plumed", ffOPTRD }, /* PLUMED */ - { efDAT, "-membed", "membed", ffOPTRD }, - { efTOP, "-mp", "membed", ffOPTRD }, - { efNDX, "-mn", "membed", ffOPTRD } - }; - #define NFILE asize(fnm) -@@ -731,19 +738,50 @@ - } - - ddxyz[XX] = (int)(realddxyz[XX] + 0.5); - ddxyz[YY] = (int)(realddxyz[YY] + 0.5); - ddxyz[ZZ] = (int)(realddxyz[ZZ] + 0.5); -+ /* PLUMED */ -+ plumedswitch=0; -+ if (opt2bSet("-plumed",NFILE,fnm)) plumedswitch=1; -+ if(plumedswitch){ -+ int plumed_is_there=0; -+ int real_precision=sizeof(real); -+ real energyUnits=1.0; -+ real lengthUnits=1.0; -+ real timeUnits=1.0; -+ -+ -+ if(!plumed_installed()){ -+ gmx_fatal(FARGS,"Plumed is not available. Check your PLUMED_KERNEL variable."); -+ } -+ plumedmain=plumed_create(); -+ plumed_cmd(plumedmain,"setRealPrecision",&real_precision); -+ // this is not necessary for gromacs units: -+ plumed_cmd(plumedmain,"setMDEnergyUnits",&energyUnits); -+ plumed_cmd(plumedmain,"setMDLengthUnits",&lengthUnits); -+ plumed_cmd(plumedmain,"setMDTimeUnits",&timeUnits); -+ // -+ plumed_cmd(plumedmain,"setPlumedDat",ftp2fn(efDAT,NFILE,fnm)); -+ plumedswitch=1; -+ } -+ /* END PLUMED */ - - rc = mdrunner(&hw_opt, fplog, cr, NFILE, fnm, oenv, bVerbose, bCompact, - nstglobalcomm, ddxyz, dd_node_order, rdd, rconstr, - dddlb_opt[0], dlb_scale, ddcsx, ddcsy, ddcsz, - nbpu_opt[0], - nsteps, nstepout, resetstep, - nmultisim, repl_ex_nst, repl_ex_nex, repl_ex_seed, - pforce, cpt_period, max_hours, deviceOptions, Flags); - -+ /* PLUMED */ -+ if(plumedswitch){ -+ plumed_finalize(plumedmain); -+ } -+ /* END PLUMED */ -+ - gmx_finalize_par(); - - if (MULTIMASTER(cr)) - { - thanx(stderr); -EOF_EOF -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 -@@ -51,10 +51,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 */ -+ - #define PROBABILITYCUTOFF 100 - /* we don't bother evaluating if events are more rare than exp(-100) = 3.7x10^-44 */ - - enum { - ereTEMP, ereLAMBDA, ereENDSINGLE, ereTL, ereNR -@@ -111,18 +117,20 @@ - - snew(qall, ms->nsim); - qall[re->repl] = q; - gmx_sum_sim(ms->nsim, qall, ms); - -- bDiff = FALSE; -- for (s = 1; s < ms->nsim; s++) -- { -- if (qall[s] != qall[0]) -- { -+ /* PLUMED */ -+ //bDiff = FALSE; -+ //for (s = 1; s < ms->nsim; s++) -+ //{ -+ // if (qall[s] != qall[0]) -+ // { - bDiff = TRUE; -- } -- } -+ // } -+ //} -+ /* END PLUMED */ - - if (bDiff) - { - /* Set the replica exchange type and quantities */ - re->type = ere; -@@ -255,10 +263,14 @@ - for (i = 0; i < re->nrepl; i++) - { - re->ind[i] = i; - } - -+ /* PLUMED */ -+ // plumed2: check if we want alternative patterns (i.e. for bias-exchange metaD) -+ // in those cases replicas can share the same temperature. -+ /* - if (re->type < ereENDSINGLE) - { - - for (i = 0; i < re->nrepl; i++) - { -@@ -275,10 +287,12 @@ - gmx_fatal(FARGS, "Two replicas have identical %ss", erename[re->type]); - } - } - } - } -+ */ -+ /* END PLUMED */ - - /* keep track of all the swaps, starting with the initial placement. */ - snew(re->allswaps, re->nrepl); - for (i = 0; i < re->nrepl; i++) - { -@@ -524,11 +538,11 @@ - } - sfree(buf); - } - } - --static void exchange_state(const gmx_multisim_t *ms, int b, t_state *state) -+void exchange_state(const gmx_multisim_t *ms, int b, t_state *state) - { - /* When t_state changes, this code should be updated. */ - int ngtc, nnhpres; - ngtc = state->ngtc * state->nhchainlength; - nnhpres = state->nnhpres* state->nhchainlength; -@@ -644,11 +658,11 @@ - svmul(fac, state->v[i], state->v[i]); - } - } - } - --static void pd_collect_state(const t_commrec *cr, t_state *state) -+void pd_collect_state(const t_commrec *cr, t_state *state) - { - int shift; - - if (debug) - { -@@ -869,10 +883,17 @@ - } - if (bPrint) - { - fprintf(fplog, "Repl %d <-> %d dE_term = %10.3e (kT)\n", a, b, delta); - } -+ -+/* PLUMED */ -+/* this is necessary because with plumed HREX the energy contribution is -+ already taken into account */ -+ if(getenv("PLUMED_HREX")) delta=0.0; -+/* END PLUMED */ -+ - if (re->bNPT) - { - /* revist the calculation for 5.0. Might be some improvements. */ - dpV = (beta[ap]*re->pres[ap]-beta[bp]*re->pres[bp])*(Vol[b]-Vol[a])/PRESFAC; - if (bPrint) -@@ -1045,19 +1066,57 @@ - } - else - { - /* standard nearest neighbor replica exchange */ - m = (step / re->nst) % 2; -+ /* PLUMED */ -+ if(plumedswitch){ -+ int partner=re->repl; -+ int test=0; -+ plumed_cmd(plumedmain,"getExchangesFlag",&test); -+ if(test>0){ -+ int *list; -+ snew(list,re->nrepl); -+ plumed_cmd(plumedmain,"setNumberOfReplicas",&(re->nrepl)); -+ plumed_cmd(plumedmain,"getExchangesList",list); -+ for(i=0; i<re->nrepl; i++) re->ind[i]=list[i]; -+ sfree(list); -+ } -+ -+ 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); -+ plumed_cmd(plumedmain,"GREX shareAllDeltaBias",NULL); -+ } -+ /* END PLUMED */ - 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) - { - delta = calc_delta(fplog, bPrint, re, a, b, a, b); -+ /* PLUMED */ -+ if(plumedswitch){ -+ real adb,bdb,dplumed; -+ char buf[300]; -+ sprintf(buf,"GREX getDeltaBias %d",a); plumed_cmd(plumedmain,buf,&adb); -+ sprintf(buf,"GREX getDeltaBias %d",b); plumed_cmd(plumedmain,buf,&bdb); -+ dplumed=adb*re->beta[a]+bdb*re->beta[b]; -+ delta+=dplumed; -+ if (bPrint) -+ fprintf(fplog,"dplumed = %10.3e dE_Term = %10.3e (kT)\n",dplumed,delta); -+ } -+ /* END PLUMED */ - if (delta <= 0) - { - /* accepted */ - prob[i] = 1; - bEx[i] = TRUE; -@@ -1306,10 +1365,14 @@ - * exchanges. */ - /* Where each replica ends up after the exchange attempt(s). */ - /* The order in which multiple exchanges will occur. */ - gmx_bool bThisReplicaExchanged = FALSE; - -+ /* PLUMED */ -+ if(plumedswitch)plumed_cmd(plumedmain,"GREX prepare",NULL); -+ /* END PLUMED */ -+ - if (MASTER(cr)) - { - replica_id = re->repl; - test_for_replica_exchange(fplog, cr->ms, re, enerd, det(state_local->box), step, time); - prepare_to_do_exchange(fplog, re->destinations, replica_id, re->nrepl, &maxswap, -@@ -1437,5 +1500,12 @@ - fprintf(fplog, "\n"); - } - /* print the transition matrix */ - print_transition_matrix(fplog, "", re->nrepl, re->nmoves, re->nattempt); - } -+ -+int replica_exchange_get_repl(const gmx_repl_ex_t re){ -+ return re->repl; -+}; -+int replica_exchange_get_nrepl(const gmx_repl_ex_t re){ -+ return re->nrepl; -+}; -EOF_EOF -patch -u -l -b -F 5 --suffix=.preplumed "./src/kernel/repl_ex.h" << \EOF_EOF ---- ./src/kernel/repl_ex.h.preplumed -+++ ./src/kernel/repl_ex.h -@@ -71,6 +71,12 @@ - /* Should only be called on the master nodes */ - - extern void pd_distribute_state(const t_commrec *cr, t_state *state); - /* Distributes the state after exchange for particle decomposition */ - -+extern int replica_exchange_get_repl(const gmx_repl_ex_t re); -+extern int replica_exchange_get_nrepl(const gmx_repl_ex_t re); -+extern void pd_collect_state(const t_commrec *cr, t_state *state); -+extern void exchange_state(const gmx_multisim_t *ms, int b, t_state *state); -+ -+ - #endif /* _repl_ex_h */ -EOF_EOF diff --git a/patches/gromacs-4.6.3.config b/patches/gromacs-4.6.3.config index a5d00fcea22295d52544fcd345e60947ff7553ef..6346748e72bdbd219da14fcd9a05d883ad7c9b5a 100644 --- a/patches/gromacs-4.6.3.config +++ b/patches/gromacs-4.6.3.config @@ -4,3 +4,26 @@ function plumed_preliminary_test(){ # check if the README contains the word GROMACS and if gromacs has been already configured grep -q GROMACS README 1>/dev/null 2>/dev/null } + +plumed_before_patch(){ +cat << EOF +This is an experimental patch to perform Hamiltonian replica exchange with GROMACS +Implementation described in: + Bussi, Hamiltonian replica-exchange in GROMACS: a flexible implementation, + Mol. Phys., DOI: 10.1080/00268976.2013.824126 (2013) + +To use it: +* Prepare separate topologies (topol0.tpr, topol1.tpr, etc) +* Set environment variable PLUMED_HREX (export PLUMED_HREX=1) +* Run a normal replica exchange with gromacs + +Suggested checks: +* Try with several identical force fields and different seed/starting point + Acceptance should be 1.0 + +Warnings: +* Topologies should have the same number of atoms, same masses and same constraint topology +* Choose neighbor list update that divides replex + +EOF +} diff --git a/patches/gromacs-4.6.3.diff/src/kernel/md.c b/patches/gromacs-4.6.3.diff/src/kernel/md.c index 187b2796a09a2cf80c82d9637db77093835ce358..ec5ffd2894102738520455ae2c85b65df937f04e 100644 --- a/patches/gromacs-4.6.3.diff/src/kernel/md.c +++ b/patches/gromacs-4.6.3.diff/src/kernel/md.c @@ -1139,6 +1139,94 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], GMX_MPE_LOG(ev_timestep2); + gmx_bool bHREX; + bHREX= repl_ex_nst > 0 && (step>0) && !bLastStep && do_per_step(step,repl_ex_nst) + && getenv("PLUMED_HREX"); + +/* Hamiltonian Replica Exchange */ + if(plumedswitch) if(bHREX){ + gmx_enerdata_t *hrex_enerd; + snew(hrex_enerd,1); + init_enerdata(top_global->groups.grps[egcENER].nr,ir->fepvals->n_lambda,hrex_enerd); + int repl=-1; + int nrepl=-1; + if(MASTER(cr)){ + repl=replica_exchange_get_repl(repl_ex); + nrepl=replica_exchange_get_nrepl(repl_ex); + } + if (PAR(cr)) { + if (DOMAINDECOMP(cr)) + dd_collect_state(cr->dd,state,state_global); + else + pd_collect_state(cr,state_global); + } + if(MASTER(cr)){ + if(repl%2==step/repl_ex_nst%2){ + if(repl-1>=0) exchange_state(cr->ms,repl-1,state_global); + }else{ + if(repl+1<nrepl) exchange_state(cr->ms,repl+1,state_global); + } + } + if(PAR(cr)){ + if (DOMAINDECOMP(cr)) { + dd_partition_system(fplog,step,cr,TRUE,1, + state_global,top_global,ir, + state,&f,mdatoms,top,fr,vsite,shellfc,constr, + nrnb,wcycle,FALSE); + } else { + bcast_state(cr,state,FALSE); + } + } + do_force(fplog, cr, ir, step, nrnb, wcycle, top, top_global, groups, + state->box, state->x, &state->hist, + f, force_vir, mdatoms, hrex_enerd, fcd, + state->lambda, graph, + fr, vsite, mu_tot, t, outf->fp_field, ed, bBornRadii, + GMX_FORCE_STATECHANGED | + ((DYNAMIC_BOX(*ir) || bRerunMD) ? GMX_FORCE_DYNAMICBOX : 0) | + GMX_FORCE_ALLFORCES | + GMX_FORCE_SEPLRF | + GMX_FORCE_VIRIAL | + GMX_FORCE_ENERGY | + GMX_FORCE_DHDL | + GMX_FORCE_NS); + plumed_cmd(plumedmain,"GREX cacheLocalUSwap",&hrex_enerd->term[F_EPOT]); + sfree(hrex_enerd); + +/* exchange back */ + if (PAR(cr)) { + if (DOMAINDECOMP(cr)) + dd_collect_state(cr->dd,state,state_global); + else + pd_collect_state(cr,state_global); + } + if(MASTER(cr)){ + if(repl%2==step/repl_ex_nst%2){ + if(repl-1>=0) exchange_state(cr->ms,repl-1,state_global); + }else{ + if(repl+1<nrepl) exchange_state(cr->ms,repl+1,state_global); + } + } + if(PAR(cr)){ + if (DOMAINDECOMP(cr)) { + dd_partition_system(fplog,step,cr,TRUE,1, + state_global,top_global,ir, + state,&f,mdatoms,top,fr,vsite,shellfc,constr, + nrnb,wcycle,FALSE); + if(plumedswitch){ + plumed_cmd(plumedmain,"setAtomsNlocal",&cr->dd->nat_home); + plumed_cmd(plumedmain,"setAtomsGatindex",cr->dd->gatindex); + } + } else { + bcast_state(cr,state,FALSE); + } + } + + } +/* END Hamiltonian Replica Exchange */ + + + /* We write a checkpoint at this MD step when: * either at an NS step when we signalled through gs, * or at the last step (but not when we do not want confout), @@ -1262,6 +1350,9 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], 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); + if(bHREX) + plumed_cmd(plumedmain,"GREX cacheLocalUNow",&enerd->term[F_EPOT]); + } /* END PLUMED */ } diff --git a/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.c b/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.c index 2e79c67ad53769d9b8df84b99b66ca2bfa18de30..2f5f9915908aac5785045757800d1bd75e449201 100644 --- a/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.c +++ b/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.c @@ -540,7 +540,7 @@ static void exchange_rvecs(const gmx_multisim_t *ms, int b, rvec *v, int n) } } -static void exchange_state(const gmx_multisim_t *ms, int b, t_state *state) +void exchange_state(const gmx_multisim_t *ms, int b, t_state *state) { /* When t_state changes, this code should be updated. */ int ngtc, nnhpres; @@ -660,7 +660,7 @@ static void scale_velocities(t_state *state, real fac) } } -static void pd_collect_state(const t_commrec *cr, t_state *state) +void pd_collect_state(const t_commrec *cr, t_state *state) { int shift; @@ -885,6 +885,13 @@ static real calc_delta(FILE *fplog, gmx_bool bPrint, struct gmx_repl_ex *re, int { fprintf(fplog, "Repl %d <-> %d dE_term = %10.3e (kT)\n", a, b, delta); } + +/* PLUMED */ +/* this is necessary because with plumed HREX the energy contribution is + already taken into account */ + if(getenv("PLUMED_HREX")) delta=0.0; +/* END PLUMED */ + if (re->bNPT) { /* revist the calculation for 5.0. Might be some improvements. */ @@ -1518,3 +1525,12 @@ void print_replica_exchange_statistics(FILE *fplog, struct gmx_repl_ex *re) /* print the transition matrix */ print_transition_matrix(fplog, "", re->nrepl, re->nmoves, re->nattempt); } + +int replica_exchange_get_repl(const gmx_repl_ex_t re){ + return re->repl; +}; + +int replica_exchange_get_nrepl(const gmx_repl_ex_t re){ + return re->nrepl; +}; + diff --git a/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.h b/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.h new file mode 100644 index 0000000000000000000000000000000000000000..4a4003b960b0dceac64c43448e5875a1a2d9d8e1 --- /dev/null +++ b/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.h @@ -0,0 +1,82 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 1991-2000, University of Groningen, The Netherlands. + * Copyright (c) 2001-2004, The GROMACS development team, + * check out http://www.gromacs.org for more information. + * Copyright (c) 2012,2013, by the GROMACS development team, led by + * David van der Spoel, Berk Hess, Erik Lindahl, and including many + * others, as listed in the AUTHORS file in the top-level source + * directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ + +#ifndef _repl_ex_h +#define _repl_ex_h + +#include "typedefs.h" +#include "types/commrec.h" + +/* Abstract type for replica exchange */ +typedef struct gmx_repl_ex *gmx_repl_ex_t; + +extern gmx_repl_ex_t init_replica_exchange(FILE *fplog, + const gmx_multisim_t *ms, + const t_state *state, + const t_inputrec *ir, + int nst, int nmultiex, int init_seed); +/* Should only be called on the master nodes */ + +extern gmx_bool replica_exchange(FILE *fplog, + const t_commrec *cr, + gmx_repl_ex_t re, + t_state *state, gmx_enerdata_t *enerd, + t_state *state_local, + gmx_large_int_t step, real time); +/* Attempts replica exchange, should be called on all nodes. + * Returns TRUE if this state has been exchanged. + * When running each replica in parallel, + * this routine collects the state on the master node before exchange. + * With particle the state is redistributed over the nodes after exchange. + * With domain decomposition the global state after exchanged in stored + * in state and still needs to be redistributed over the nodes. + */ + +extern void print_replica_exchange_statistics(FILE *fplog, gmx_repl_ex_t re); +/* Should only be called on the master nodes */ + +extern void pd_distribute_state(const t_commrec *cr, t_state *state); +/* Distributes the state after exchange for particle decomposition */ + +extern int replica_exchange_get_repl(const gmx_repl_ex_t re); +extern int replica_exchange_get_nrepl(const gmx_repl_ex_t re); +extern void pd_collect_state(const t_commrec *cr, t_state *state); +extern void exchange_state(const gmx_multisim_t *ms, int b, t_state *state); + + +#endif /* _repl_ex_h */ diff --git a/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.h.preplumed b/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.h.preplumed new file mode 100644 index 0000000000000000000000000000000000000000..97c45c56b61e77c00b1cf111937ae6e958d62c5a --- /dev/null +++ b/patches/gromacs-4.6.3.diff/src/kernel/repl_ex.h.preplumed @@ -0,0 +1,76 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 1991-2000, University of Groningen, The Netherlands. + * Copyright (c) 2001-2004, The GROMACS development team, + * check out http://www.gromacs.org for more information. + * Copyright (c) 2012,2013, by the GROMACS development team, led by + * David van der Spoel, Berk Hess, Erik Lindahl, and including many + * others, as listed in the AUTHORS file in the top-level source + * directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ + +#ifndef _repl_ex_h +#define _repl_ex_h + +#include "typedefs.h" +#include "types/commrec.h" + +/* Abstract type for replica exchange */ +typedef struct gmx_repl_ex *gmx_repl_ex_t; + +extern gmx_repl_ex_t init_replica_exchange(FILE *fplog, + const gmx_multisim_t *ms, + const t_state *state, + const t_inputrec *ir, + int nst, int nmultiex, int init_seed); +/* Should only be called on the master nodes */ + +extern gmx_bool replica_exchange(FILE *fplog, + const t_commrec *cr, + gmx_repl_ex_t re, + t_state *state, gmx_enerdata_t *enerd, + t_state *state_local, + gmx_large_int_t step, real time); +/* Attempts replica exchange, should be called on all nodes. + * Returns TRUE if this state has been exchanged. + * When running each replica in parallel, + * this routine collects the state on the master node before exchange. + * With particle the state is redistributed over the nodes after exchange. + * With domain decomposition the global state after exchanged in stored + * in state and still needs to be redistributed over the nodes. + */ + +extern void print_replica_exchange_statistics(FILE *fplog, gmx_repl_ex_t re); +/* Should only be called on the master nodes */ + +extern void pd_distribute_state(const t_commrec *cr, t_state *state); +/* Distributes the state after exchange for particle decomposition */ + +#endif /* _repl_ex_h */