Merge branch 'release-2019' into master
[alexxy/gromacs.git] / src / gromacs / mdrun / rerun.cpp
index d3be42c82a89bfe5684bbf914079f8060b534c81..fc6d53e60c2759d20ad3bc842a33aa4521d80e4d 100644 (file)
@@ -53,7 +53,6 @@
 
 #include "gromacs/awh/awh.h"
 #include "gromacs/commandline/filenm.h"
-#include "gromacs/compat/make_unique.h"
 #include "gromacs/domdec/collect.h"
 #include "gromacs/domdec/domdec.h"
 #include "gromacs/domdec/domdec_network.h"
 #include "gromacs/domdec/partition.h"
 #include "gromacs/essentialdynamics/edsam.h"
 #include "gromacs/ewald/pme.h"
-#include "gromacs/ewald/pme-load-balancing.h"
+#include "gromacs/ewald/pme_load_balancing.h"
 #include "gromacs/fileio/trxio.h"
 #include "gromacs/gmxlib/network.h"
 #include "gromacs/gmxlib/nrnb.h"
 #include "gromacs/gpu_utils/gpu_utils.h"
 #include "gromacs/imd/imd.h"
-#include "gromacs/listed-forces/manage-threading.h"
+#include "gromacs/listed_forces/manage_threading.h"
 #include "gromacs/math/functions.h"
 #include "gromacs/math/utilities.h"
 #include "gromacs/math/vec.h"
 #include "gromacs/mdlib/compute_io.h"
 #include "gromacs/mdlib/constr.h"
 #include "gromacs/mdlib/ebin.h"
+#include "gromacs/mdlib/energyoutput.h"
 #include "gromacs/mdlib/expanded.h"
 #include "gromacs/mdlib/force.h"
 #include "gromacs/mdlib/force_flags.h"
 #include "gromacs/mdlib/forcerec.h"
 #include "gromacs/mdlib/md_support.h"
 #include "gromacs/mdlib/mdatoms.h"
-#include "gromacs/mdlib/mdebin.h"
 #include "gromacs/mdlib/mdoutf.h"
 #include "gromacs/mdlib/mdrun.h"
 #include "gromacs/mdlib/mdsetup.h"
 #include "gromacs/mdlib/membed.h"
-#include "gromacs/mdlib/nb_verlet.h"
-#include "gromacs/mdlib/nbnxn_gpu_data_mgmt.h"
 #include "gromacs/mdlib/ns.h"
 #include "gromacs/mdlib/resethandler.h"
 #include "gromacs/mdlib/shellfc.h"
 #include "gromacs/mdlib/update.h"
 #include "gromacs/mdlib/vcm.h"
 #include "gromacs/mdlib/vsite.h"
-#include "gromacs/mdtypes/awh-history.h"
-#include "gromacs/mdtypes/awh-params.h"
+#include "gromacs/mdtypes/awh_history.h"
+#include "gromacs/mdtypes/awh_params.h"
 #include "gromacs/mdtypes/commrec.h"
 #include "gromacs/mdtypes/df_history.h"
 #include "gromacs/mdtypes/energyhistory.h"
@@ -195,7 +192,6 @@ void gmx::Integrator::do_rerun()
     // t_inputrec is being replaced by IMdpOptionsProvider, so this
     // will go away eventually.
     t_inputrec             *ir   = inputrec;
-    gmx_mdoutf             *outf = nullptr;
     int64_t                 step, step_rel;
     double                  t, lam0[efptNR];
     bool                    isLastStep               = false;
@@ -205,14 +201,12 @@ void gmx::Integrator::do_rerun()
     t_trxstatus            *status;
     rvec                    mu_tot;
     t_trxframe              rerun_fr;
-    gmx_localtop_t         *top;
-    t_mdebin               *mdebin   = nullptr;
+    gmx_localtop_t          top;
     gmx_enerdata_t         *enerd;
     PaddedVector<gmx::RVec> f {};
     gmx_global_stat_t       gstat;
     t_graph                *graph = nullptr;
     gmx_groups_t           *groups;
-    gmx_ekindata_t         *ekind;
     gmx_shellfc_t          *shellfc;
 
     double                  cycles;
@@ -302,10 +296,11 @@ void gmx::Integrator::do_rerun()
         top_global->intermolecularExclusionGroup = genQmmmIndices(*top_global);
     }
 
-    /* Initial values */
-    init_rerun(fplog, cr, outputProvider, ir, oenv, mdrunOptions,
-               state_global, lam0, nrnb, top_global,
-               nfile, fnm, &outf, &mdebin, wcycle);
+    initialize_lambdas(fplog, *ir, MASTER(cr), &state_global->fep_state, state_global->lambda, lam0);
+    init_nrnb(nrnb);
+    gmx_mdoutf       *outf = init_mdoutf(fplog, nfile, fnm, mdrunOptions, cr, outputProvider, ir, top_global, oenv, wcycle);
+    gmx::EnergyOutput energyOutput;
+    energyOutput.prepare(mdoutf_get_fp_ene(outf), top_global, ir, mdoutf_get_fp_dhdl(outf), true);
 
     /* Energy terms and groups */
     snew(enerd, 1);
@@ -313,7 +308,8 @@ void gmx::Integrator::do_rerun()
                   enerd);
 
     /* Kinetic energy data */
-    snew(ekind, 1);
+    std::unique_ptr<gmx_ekindata_t> eKinData = std::make_unique<gmx_ekindata_t>();
+    gmx_ekindata_t                 *ekind    = eKinData.get();
     init_ekindata(fplog, top_global, &(ir->opts), ekind);
     /* Copy the cos acceleration to the groups struct */
     ekind->cosacc.cos_accel = ir->cos_accel;
@@ -326,7 +322,7 @@ void gmx::Integrator::do_rerun()
                                  ir->nstcalcenergy, DOMAINDECOMP(cr));
 
     {
-        double io = compute_io(ir, top_global->natoms, groups, mdebin->ebin->nener, 1);
+        double io = compute_io(ir, top_global->natoms, groups, energyOutput.numEnergyTerms(), 1);
         if ((io > 2000) && MASTER(cr))
         {
             fprintf(stderr,
@@ -341,16 +337,16 @@ void gmx::Integrator::do_rerun()
 
     if (DOMAINDECOMP(cr))
     {
-        top = dd_init_local_top(top_global);
+        dd_init_local_top(*top_global, &top);
 
-        stateInstance = compat::make_unique<t_state>();
+        stateInstance = std::make_unique<t_state>();
         state         = stateInstance.get();
         dd_init_local_state(cr->dd, state_global, state);
 
         /* Distribute the charge groups over the nodes from the master node */
         dd_partition_system(fplog, mdlog, ir->init_step, cr, TRUE, 1,
-                            state_global, top_global, ir,
-                            state, &f, mdAtoms, top, fr,
+                            state_global, *top_global, ir,
+                            state, &f, mdAtoms, &top, fr,
                             vsite, constr,
                             nrnb, nullptr, FALSE);
         shouldCheckNumberOfBondedInteractions = true;
@@ -365,8 +361,7 @@ void gmx::Integrator::do_rerun()
         /* Copy the pointer to the global state */
         state = state_global;
 
-        snew(top, 1);
-        mdAlgorithmsSetupAtomData(cr, ir, top_global, top, fr,
+        mdAlgorithmsSetupAtomData(cr, ir, *top_global, &top, fr,
                                   &graph, mdAtoms, constr, vsite, shellfc);
     }
 
@@ -395,7 +390,7 @@ void gmx::Integrator::do_rerun()
                         &totalNumberOfBondedInteractions, &bSumEkinhOld, cglo_flags);
     }
     checkNumberOfBondedInteractions(mdlog, cr, totalNumberOfBondedInteractions,
-                                    top_global, top, state,
+                                    top_global, &top, state,
                                     &shouldCheckNumberOfBondedInteractions);
 
     if (MASTER(cr))
@@ -521,7 +516,7 @@ void gmx::Integrator::do_rerun()
                 gmx_fatal(FARGS, "Vsite recalculation with -rerun is not implemented with domain decomposition, "
                           "use a single rank");
             }
-            prepareRerunState(rerun_fr, state_global, constructVsites, vsite, top->idef, ir->delta_t, *fr, graph);
+            prepareRerunState(rerun_fr, state_global, constructVsites, vsite, top.idef, ir->delta_t, *fr, graph);
         }
 
         isLastStep = isLastStep || stopHandler->stoppingAfterCurrentStep(bNS);
@@ -532,8 +527,8 @@ void gmx::Integrator::do_rerun()
             const bool bMasterState = true;
             dd_partition_system(fplog, mdlog, step, cr,
                                 bMasterState, nstglobalcomm,
-                                state_global, top_global, ir,
-                                state, &f, mdAtoms, top, fr,
+                                state_global, *top_global, ir,
+                                state, &f, mdAtoms, &top, fr,
                                 vsite, constr,
                                 nrnb, wcycle,
                                 mdrunOptions.verbose);
@@ -562,7 +557,7 @@ void gmx::Integrator::do_rerun()
             /* Now is the time to relax the shells */
             relax_shell_flexcon(fplog, cr, ms, mdrunOptions.verbose,
                                 enforcedRotation, step,
-                                ir, bNS, force_flags, top,
+                                ir, bNS, force_flags, &top,
                                 constr, enerd, fcd,
                                 state, f.arrayRefWithPadding(), force_vir, mdatoms,
                                 nrnb, wcycle, graph, groups,
@@ -580,7 +575,7 @@ void gmx::Integrator::do_rerun()
             Awh       *awh = nullptr;
             gmx_edsam *ed  = nullptr;
             do_force(fplog, cr, ms, ir, awh, enforcedRotation,
-                     step, nrnb, wcycle, top, groups,
+                     step, nrnb, wcycle, &top, groups,
                      state->box, state->x.arrayRefWithPadding(), &state->hist,
                      f.arrayRefWithPadding(), force_vir, mdatoms, enerd, fcd,
                      state->lambda, graph,
@@ -599,7 +594,7 @@ void gmx::Integrator::do_rerun()
             do_md_trajectory_writing(fplog, cr, nfile, fnm, step, step_rel, t,
                                      ir, state, state_global, observablesHistory,
                                      top_global, fr,
-                                     outf, mdebin, ekind, f,
+                                     outf, energyOutput, ekind, f,
                                      isCheckpointingStep, doRerun, isLastStep,
                                      mdrunOptions.writeConfout,
                                      bSumEkinhOld);
@@ -621,7 +616,7 @@ void gmx::Integrator::do_rerun()
                 shift_self(graph, state->box, as_rvec_array(state->x.data()));
             }
             construct_vsites(vsite, as_rvec_array(state->x.data()), ir->delta_t, as_rvec_array(state->v.data()),
-                             top->idef.iparams, top->idef.il,
+                             top.idef.iparams, top.idef.il,
                              fr->ePBC, fr->bMolPBC, cr, state->box);
 
             if (graph != nullptr)
@@ -647,7 +642,7 @@ void gmx::Integrator::do_rerun()
                             | (shouldCheckNumberOfBondedInteractions ? CGLO_CHECK_NUMBER_OF_BONDED_INTERACTIONS : 0)
                             );
             checkNumberOfBondedInteractions(mdlog, cr, totalNumberOfBondedInteractions,
-                                            top_global, top, state,
+                                            top_global, &top, state,
                                             &shouldCheckNumberOfBondedInteractions);
         }
 
@@ -667,11 +662,11 @@ void gmx::Integrator::do_rerun()
         if (MASTER(cr))
         {
             const bool bCalcEnerStep = true;
-            upd_mdebin(mdebin, doFreeEnergyPerturbation, bCalcEnerStep,
-                       t, mdatoms->tmass, enerd, state,
-                       ir->fepvals, ir->expandedvals, rerun_fr.box,
-                       shake_vir, force_vir, total_vir, pres,
-                       ekind, mu_tot, constr);
+            energyOutput.addDataAtEnergyStep(doFreeEnergyPerturbation, bCalcEnerStep,
+                                             t, mdatoms->tmass, enerd, state,
+                                             ir->fepvals, ir->expandedvals, state->box,
+                                             shake_vir, force_vir, total_vir, pres,
+                                             ekind, mu_tot, constr);
 
             const bool do_ene = true;
             const bool do_log = true;
@@ -679,9 +674,10 @@ void gmx::Integrator::do_rerun()
             const bool do_dr  = ir->nstdisreout != 0;
             const bool do_or  = ir->nstorireout != 0;
 
-            print_ebin(mdoutf_get_fp_ene(outf), do_ene, do_dr, do_or, do_log ? fplog : nullptr,
-                       step, t,
-                       eprNORMAL, mdebin, fcd, groups, &(ir->opts), awh);
+            energyOutput.printStepToEnergyFile(mdoutf_get_fp_ene(outf), do_ene, do_dr, do_or,
+                                               do_log ? fplog : nullptr,
+                                               step, t,
+                                               eprNORMAL, fcd, groups, &(ir->opts), awh);
 
             if (do_per_step(step, ir->nstlog))
             {
@@ -761,7 +757,6 @@ void gmx::Integrator::do_rerun()
         gmx_pme_send_finish(cr);
     }
 
-    done_mdebin(mdebin);
     done_mdoutf(outf);
 
     done_shellfc(fplog, shellfc, step_rel);
@@ -776,5 +771,4 @@ void gmx::Integrator::do_rerun()
 
     destroy_enerdata(enerd);
     sfree(enerd);
-    sfree(top);
 }