Removes memory allocation from free energy foreign lambda calcs.

Removes the memory allocation from the free energy foreign
lambda calcs by declaring some new space in the enerd array, and
reinitializing only this space. Also makes the foreign lambda
calcs conditional on the group scheme and the non-bonded flag.

Change-Id: I98ed2f6ff45a4d33f9aa7dbde91364cf56c38eae

diff --git a/include/bondf.h b/include/bondf.h
index 1dce4f32e10a936c860ba036405519569fdae154..d66faf4de95c6dd249ad97d2ae02d9dd31fa39d0 100644 (file)
--- a/include/bondf.h
+++ b/include/bondf.h
@@ -99,7 +99,7 @@ void calc_bonds_lambda(FILE *fplog,
                              rvec x[],
                              t_forcerec *fr,
                              const t_pbc *pbc,const t_graph *g,
-                             gmx_enerdata_t *enerd,t_nrnb *nrnb,
+                  gmx_grppairener_t *grpp, real *epot,t_nrnb *nrnb,
                              real *lambda,
                              const t_mdatoms *md,
                              t_fcdata *fcd,int *global_atom_index);

diff --git a/include/force.h b/include/force.h
index 0aa20283d7e6461141661ec625cf1117fe244548..a91be65d17037aecb49f8ef95868abebd26e38c2 100644 (file)
--- a/include/force.h
+++ b/include/force.h
@@ -201,13 +201,16 @@ void init_enerdata(int ngener,int n_lambda,gmx_enerdata_t *enerd);
 void destroy_enerdata(gmx_enerdata_t *enerd);
 /* Free all memory associated with enerd */
 
+void reset_foreign_enerdata(gmx_enerdata_t *enerd);
+/* Resets only the foreign energy data */
+
 void reset_enerdata(t_grpopts *opts,
                           t_forcerec *fr,gmx_bool bNS,
                           gmx_enerdata_t *enerd,
                           gmx_bool bMaster);
 /* Resets the energy data, if bNS=TRUE also zeros the long-range part */
 
-void sum_epot(t_grpopts *opts,gmx_enerdata_t *enerd);
+void sum_epot(t_grpopts *opts, gmx_grppairener_t *grpp, real *epot);
 /* Locally sum the non-bonded potential energy terms */
 
 GMX_LIBMD_EXPORT

diff --git a/include/types/forcerec.h b/include/types/forcerec.h
index a5b8c3e15e831fdd88e9494ee7c66fb3cdc7050e..a55a601a1800fee6631827832b7e176c33a3f61e 100644 (file)
--- a/include/types/forcerec.h
+++ b/include/types/forcerec.h
@@ -145,6 +145,8 @@ typedef struct {
   int    n_lambda;
   int    fep_state;              /*current fep state -- just for printing */
   double *enerpart_lambda; /* Partial energy for lambda and flambda[] */
+  real foreign_term[F_NRE];    /* alternate array for storing foreign lambda energies */
+  gmx_grppairener_t foreign_grpp;  /* alternate array for storing foreign lambda energies */
 } gmx_enerdata_t;
 /* The idea is that dvdl terms with linear lambda dependence will be added
  * automatically to enerpart_lambda. Terms with non-linear lambda dependence

diff --git a/src/gmxlib/bondfree.c b/src/gmxlib/bondfree.c
index 5058b899c42d1f5b273d7c9e81db5790eff7fbcb..42a0c25679fa8190c2bcacfe0a2ab9e1f1106036 100644 (file)
--- a/src/gmxlib/bondfree.c
+++ b/src/gmxlib/bondfree.c
@@ -3615,7 +3615,7 @@ static real calc_one_bond(FILE *fplog,int thread,
 static real calc_one_bond_foreign(FILE *fplog,int ftype, const t_idef *idef,
                                   rvec x[], rvec f[], t_forcerec *fr,
                                   const t_pbc *pbc,const t_graph *g,
-                                  gmx_enerdata_t *enerd, t_nrnb *nrnb,
+                                  gmx_grppairener_t *grpp, t_nrnb *nrnb,
                                   real *lambda, real *dvdl,
                                   const t_mdatoms *md,t_fcdata *fcd,
                                   int *global_atom_index, gmx_bool bPrintSepPot)
@@ -3670,7 +3670,7 @@ static real calc_one_bond_foreign(FILE *fplog,int ftype, const t_idef *idef,
                                             idef->iparams,
                                             (const rvec*)x,f,fr->fshift,
                                             pbc,g,lambda,dvdl,
-                                            md,fr,&enerd->grpp,global_atom_index);
+                                            md,fr,grpp,global_atom_index);
                 }
                 if (ind != -1)
                 {
@@ -3824,14 +3824,14 @@ void calc_bonds_lambda(FILE *fplog,
                        rvec x[],
                        t_forcerec *fr,
                        const t_pbc *pbc,const t_graph *g,
-                       gmx_enerdata_t *enerd,t_nrnb *nrnb,
+                       gmx_grppairener_t *grpp, real *epot, t_nrnb *nrnb,
                        real *lambda,
                        const t_mdatoms *md,
                        t_fcdata *fcd,
                        int *global_atom_index)
 {
     int    i,ftype,nbonds_np,nbonds,ind,nat;
-    real   v,dr,dr2,*epot;
+    real   v,dr,dr2;
     real   dvdl_dum[efptNR];
     rvec   *f,*fshift_orig;
     const  t_pbc *pbc_null;
@@ -3846,8 +3846,6 @@ void calc_bonds_lambda(FILE *fplog,
         pbc_null = NULL;
     }
 
-    epot = enerd->term;
-
     snew(f,fr->natoms_force);
     /* We want to preserve the fshift array in forcerec */
     fshift_orig = fr->fshift;
@@ -3857,7 +3855,7 @@ void calc_bonds_lambda(FILE *fplog,
     for(ftype=0; (ftype<F_NRE); ftype++) 
     {
         v = calc_one_bond_foreign(fplog,ftype,idef,x, 
-                                  f,fr,pbc_null,g,enerd,nrnb,lambda,dvdl_dum,
+                                  f,fr,pbc_null,g,grpp,nrnb,lambda,dvdl_dum,
                                   md,fcd,global_atom_index,FALSE);
         epot[ftype] += v;
     }

diff --git a/src/mdlib/force.c b/src/mdlib/force.c
index abfab008019b174a833f4e23185554250e267b89..2432b763c4368c31191e67b418d4b712582de5ed 100644 (file)
--- a/src/mdlib/force.c
+++ b/src/mdlib/force.c
@@ -170,7 +170,6 @@ void do_force_lowlevel(FILE       *fplog,   gmx_large_int_t step,
     t_pbc   pbc;
     real    dvdgb;
     char    buf[22];
-    gmx_enerdata_t ed_lam;
     double  clam_i,vlam_i;
     real    dvdl_dum[efptNR], dvdl, dvdl_nb[efptNR], lam_i[efptNR];
     real    dvdlsum;
@@ -278,35 +277,30 @@ void do_force_lowlevel(FILE       *fplog,   gmx_large_int_t step,
         do_nonbonded(cr,fr,x,f,f_longrange,md,excl,
                     &enerd->grpp,box_size,nrnb,
                     lambda,dvdl_nb,-1,-1,donb_flags);
-        wallcycle_sub_stop(wcycle, ewcsNONBONDED);
-    }
 
-    /* If we do foreign lambda and we have soft-core interactions
-     * we have to recalculate the (non-linear) energies contributions.
-     */
-    if (fepvals->n_lambda > 0 && (flags & GMX_FORCE_DHDL) && fepvals->sc_alpha != 0)
-    {
-        wallcycle_sub_start(wcycle, ewcsNONBONDED);
-        init_enerdata(mtop->groups.grps[egcENER].nr,fepvals->n_lambda,&ed_lam);
-
-        for(i=0; i<enerd->n_lambda; i++)
+        /* If we do foreign lambda and we have soft-core interactions
+         * we have to recalculate the (non-linear) energies contributions.
+         */
+        if (fepvals->n_lambda > 0 && (flags & GMX_FORCE_DHDL) && fepvals->sc_alpha != 0)
         {
-            for (j=0;j<efptNR;j++)
+            for(i=0; i<enerd->n_lambda; i++)
             {
-                lam_i[j] = (i==0 ? lambda[j] : fepvals->all_lambda[j][i-1]);
+                for (j=0;j<efptNR;j++)
+                {
+                    lam_i[j] = (i==0 ? lambda[j] : fepvals->all_lambda[j][i-1]);
+                }
+                reset_foreign_enerdata(enerd);
+                do_nonbonded(cr,fr,x,f,f_longrange,md,excl,
+                             &(enerd->foreign_grpp),box_size,nrnb,
+                             lam_i,dvdl_dum,-1,-1,
+                             (donb_flags & ~GMX_NONBONDED_DO_FORCE) | GMX_NONBONDED_DO_FOREIGNLAMBDA);
+                sum_epot(&ir->opts,&(enerd->foreign_grpp),enerd->foreign_term);
+                enerd->enerpart_lambda[i] += enerd->foreign_term[F_EPOT];
             }
-            reset_enerdata(&ir->opts,fr,TRUE,&ed_lam,FALSE);
-            do_nonbonded(cr,fr,x,f,f_longrange,md,excl,
-                         &(ed_lam.grpp), box_size,nrnb,
-                         lam_i,dvdl_dum,-1,-1,
-                         GMX_NONBONDED_DO_FOREIGNLAMBDA | GMX_NONBONDED_DO_SR);
-            sum_epot(&ir->opts,&ed_lam);
-            enerd->enerpart_lambda[i] += ed_lam.term[F_EPOT];
         }
-        destroy_enerdata(&ed_lam);
         wallcycle_sub_stop(wcycle, ewcsNONBONDED);
+        where();
     }
-    where();
 
        /* If we are doing GB, calculate bonded forces and apply corrections
         * to the solvation forces */
@@ -424,21 +418,18 @@ void do_force_lowlevel(FILE       *fplog,   gmx_large_int_t step,
             {
                 gmx_incons("The bonded interactions are not sorted for free energy");
             }
-            init_enerdata(mtop->groups.grps[egcENER].nr,fepvals->n_lambda,&ed_lam);
-
             for(i=0; i<enerd->n_lambda; i++)
             {
-                reset_enerdata(&ir->opts,fr,TRUE,&ed_lam,FALSE);
+                reset_foreign_enerdata(enerd);
                 for (j=0;j<efptNR;j++)
                 {
                     lam_i[j] = (i==0 ? lambda[j] : fepvals->all_lambda[j][i-1]);
                 }
-                calc_bonds_lambda(fplog,idef,x,fr,&pbc,graph,&ed_lam,nrnb,lam_i,md,
+                calc_bonds_lambda(fplog,idef,x,fr,&pbc,graph,&(enerd->foreign_grpp),enerd->foreign_term,nrnb,lam_i,md,
                                   fcd,DOMAINDECOMP(cr) ? cr->dd->gatindex : NULL);
-                sum_epot(&ir->opts,&ed_lam);
-                enerd->enerpart_lambda[i] += ed_lam.term[F_EPOT];
+                sum_epot(&ir->opts,&(enerd->foreign_grpp),enerd->foreign_term);
+                enerd->enerpart_lambda[i] += enerd->foreign_term[F_EPOT];
             }
-            destroy_enerdata(&ed_lam);
         }
         debug_gmx();
         GMX_MPE_LOG(ev_calc_bonds_finish);
@@ -702,6 +693,7 @@ void init_enerdata(int ngener,int n_lambda,gmx_enerdata_t *enerd)
     for(i=0; i<F_NRE; i++)
     {
         enerd->term[i] = 0;
+        enerd->foreign_term[i] = 0;
     }
 
 
@@ -716,9 +708,11 @@ void init_enerdata(int ngener,int n_lambda,gmx_enerdata_t *enerd)
         fprintf(debug,"Creating %d sized group matrix for energies\n",n2);
     }
     enerd->grpp.nener = n2;
+    enerd->foreign_grpp.nener = n2;
     for(i=0; (i<egNR); i++)
     {
         snew(enerd->grpp.ener[i],n2);
+        snew(enerd->foreign_grpp.ener[i],n2);
     }
 
     if (n_lambda)
@@ -741,6 +735,11 @@ void destroy_enerdata(gmx_enerdata_t *enerd)
         sfree(enerd->grpp.ener[i]);
     }
 
+    for(i=0; (i<egNR); i++)
+    {
+        sfree(enerd->foreign_grpp.ener[i]);
+    }
+
     if (enerd->n_lambda)
     {
         sfree(enerd->enerpart_lambda);
@@ -759,15 +758,10 @@ static real sum_v(int n,real v[])
   return t;
 }
 
-void sum_epot(t_grpopts *opts,gmx_enerdata_t *enerd)
+void sum_epot(t_grpopts *opts, gmx_grppairener_t *grpp, real *epot)
 {
-  gmx_grppairener_t *grpp;
-  real *epot;
   int i;
 
-  grpp = &enerd->grpp;
-  epot = enerd->term;
-
   /* Accumulate energies */
   epot[F_COUL_SR]  = sum_v(grpp->nener,grpp->ener[egCOULSR]);
   epot[F_LJ]       = sum_v(grpp->nener,grpp->ener[egLJSR]);
@@ -887,6 +881,28 @@ void sum_dhdl(gmx_enerdata_t *enerd, real *lambda, t_lambda *fepvals)
     }
 }
 
+
+void reset_foreign_enerdata(gmx_enerdata_t *enerd)
+{
+    int  i,j;
+
+    /* First reset all foreign energy components.  Foreign energies always called on
+       neighbor search steps */
+    for(i=0; (i<egNR); i++)
+    {
+        for(j=0; (j<enerd->grpp.nener); j++)
+        {
+            enerd->foreign_grpp.ener[i][j] = 0.0;
+        }
+    }
+
+    /* potential energy components */
+    for(i=0; (i<=F_EPOT); i++)
+    {
+        enerd->foreign_term[i] = 0.0;
+    }
+}
+
 void reset_enerdata(t_grpopts *opts,
                     t_forcerec *fr,gmx_bool bNS,
                     gmx_enerdata_t *enerd,
@@ -931,4 +947,6 @@ void reset_enerdata(t_grpopts *opts,
             enerd->enerpart_lambda[i] = 0.0;
         }
     }
+    /* reset foreign energy data - separate function since we also call it elsewhere */
+    reset_foreign_enerdata(enerd);
 }

diff --git a/src/mdlib/sim_util.c b/src/mdlib/sim_util.c
index 667287a719ccb58c1c8d0b27f161ba47578ea363..8d8cd81f0515b1aaecc9141cc9525d39d39dcb5c 100644 (file)
--- a/src/mdlib/sim_util.c
+++ b/src/mdlib/sim_util.c
@@ -1407,7 +1407,7 @@ void do_force_cutsVERLET(FILE *fplog,t_commrec *cr,
     }
     
     /* Sum the potential energy terms from group contributions */
-    sum_epot(&(inputrec->opts),enerd);
+    sum_epot(&(inputrec->opts),&(enerd->grpp),enerd->term);
 }
 
 void do_force_cutsGROUP(FILE *fplog,t_commrec *cr,
@@ -1899,7 +1899,7 @@ void do_force_cutsGROUP(FILE *fplog,t_commrec *cr,
     }
 
     /* Sum the potential energy terms from group contributions */
-    sum_epot(&(inputrec->opts),enerd);
+    sum_epot(&(inputrec->opts),&(enerd->grpp),enerd->term);
 }
 
 void do_force(FILE *fplog,t_commrec *cr,