Gromacs  2024.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Classes | Typedefs | Functions
minimize.cpp File Reference
#include "gmxpre.h"
#include "config.h"
#include <cmath>
#include <cstring>
#include <ctime>
#include <algorithm>
#include <limits>
#include <vector>
#include "gromacs/commandline/filenm.h"
#include "gromacs/domdec/collect.h"
#include "gromacs/domdec/dlbtiming.h"
#include "gromacs/domdec/domdec.h"
#include "gromacs/domdec/domdec_struct.h"
#include "gromacs/domdec/mdsetup.h"
#include "gromacs/domdec/partition.h"
#include "gromacs/ewald/pme_pp.h"
#include "gromacs/fileio/confio.h"
#include "gromacs/fileio/mtxio.h"
#include "gromacs/gmxlib/network.h"
#include "gromacs/gmxlib/nrnb.h"
#include "gromacs/imd/imd.h"
#include "gromacs/linearalgebra/sparsematrix.h"
#include "gromacs/listed_forces/listed_forces.h"
#include "gromacs/listed_forces/listed_forces_gpu.h"
#include "gromacs/math/functions.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/constr.h"
#include "gromacs/mdlib/coupling.h"
#include "gromacs/mdlib/ebin.h"
#include "gromacs/mdlib/enerdata_utils.h"
#include "gromacs/mdlib/energyoutput.h"
#include "gromacs/mdlib/force.h"
#include "gromacs/mdlib/force_flags.h"
#include "gromacs/mdlib/forcerec.h"
#include "gromacs/mdlib/gmx_omp_nthreads.h"
#include "gromacs/mdlib/md_support.h"
#include "gromacs/mdlib/mdatoms.h"
#include "gromacs/mdlib/stat.h"
#include "gromacs/mdlib/tgroup.h"
#include "gromacs/mdlib/trajectory_writing.h"
#include "gromacs/mdlib/update.h"
#include "gromacs/mdlib/vsite.h"
#include "gromacs/mdrunutility/handlerestart.h"
#include "gromacs/mdrunutility/printtime.h"
#include "gromacs/mdtypes/checkpointdata.h"
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/forcebuffers.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/inputrec.h"
#include "gromacs/mdtypes/interaction_const.h"
#include "gromacs/mdtypes/md_enums.h"
#include "gromacs/mdtypes/mdatom.h"
#include "gromacs/mdtypes/mdrunoptions.h"
#include "gromacs/mdtypes/multipletimestepping.h"
#include "gromacs/mdtypes/observablesreducer.h"
#include "gromacs/mdtypes/state.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/taskassignment/include/gromacs/taskassignment/decidesimulationworkload.h"
#include "gromacs/timing/wallcycle.h"
#include "gromacs/timing/walltime_accounting.h"
#include "gromacs/topology/mtop_util.h"
#include "gromacs/topology/topology.h"
#include "gromacs/utility/cstringutil.h"
#include "gromacs/utility/exceptions.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/logger.h"
#include "gromacs/utility/smalloc.h"
#include "legacysimulator.h"
#include "shellfc.h"
+ Include dependency graph for minimize.cpp:

Description

This file defines integrators for energy minimization.

Author
Berk Hess hess@.nosp@m.kth..nosp@m.se
Erik Lindahl erik@.nosp@m.kth..nosp@m.se

Classes

struct  em_state
 Utility structure for manipulating states during EM. More...
 
class  anonymous_namespace{minimize.cpp}::EnergyEvaluator
 Class to handle the work of setting and doing an energy evaluation. More...
 

Typedefs

typedef struct em_state em_state_t
 Utility structure for manipulating states during EM.
 

Functions

static void print_em_start (FILE *fplog, const t_commrec *cr, gmx_walltime_accounting_t walltime_accounting, gmx_wallcycle *wcycle, const char *name)
 Print the EM starting conditions.
 
static void em_time_end (gmx_walltime_accounting_t walltime_accounting, gmx_wallcycle *wcycle)
 Stop counting time for EM.
 
static void sp_header (FILE *out, const char *minimizer, real ftol, int nsteps)
 Printing a log file and console header.
 
static void warn_step (FILE *fp, real ftol, real fmax, gmx_bool bLastStep, gmx_bool bConstrain)
 Print warning message.
 
static void print_converged (FILE *fp, const char *alg, real ftol, int64_t count, gmx_bool bDone, int64_t nsteps, const em_state_t *ems, double sqrtNumAtoms)
 Print message about convergence of the EM.
 
static void get_f_norm_max (const t_commrec *cr, const t_grpopts *opts, t_mdatoms *mdatoms, gmx::ArrayRef< const gmx::RVec > f, real *fnorm, real *fmax, int *a_fmax)
 Compute the norm and max of the force array in parallel.
 
static void get_state_f_norm_max (const t_commrec *cr, const t_grpopts *opts, t_mdatoms *mdatoms, em_state_t *ems)
 Compute the norm of the force.
 
static void init_em (FILE *fplog, const gmx::MDLogger &mdlog, const char *title, const t_commrec *cr, const t_inputrec *ir, const MDModulesNotifiers &mdModulesNotifiers, gmx::ImdSession *imdSession, pull_t *pull_work, t_state *state_global, const gmx_mtop_t &top_global, em_state_t *ems, gmx_localtop_t *top, t_nrnb *nrnb, t_forcerec *fr, gmx::MDAtoms *mdAtoms, gmx_global_stat_t *gstat, VirtualSitesHandler *vsite, gmx::Constraints *constr, gmx_shellfc_t **shellfc)
 Initialize the energy minimization.
 
static void finish_em (const t_commrec *cr, gmx_mdoutf_t outf, gmx_walltime_accounting_t walltime_accounting, gmx_wallcycle *wcycle)
 Finalize the minimization.
 
static void swap_em_state (em_state_t **ems1, em_state_t **ems2)
 Swap two different EM states during minimization.
 
static void write_em_traj (FILE *fplog, const t_commrec *cr, gmx_mdoutf_t outf, gmx_bool bX, gmx_bool bF, const char *confout, const gmx_mtop_t &top_global, const t_inputrec *ir, int64_t step, em_state_t *state, t_state *state_global, ObservablesHistory *observablesHistory)
 Save the EM trajectory.
 
static bool do_em_step (const t_commrec *cr, const t_inputrec *ir, t_mdatoms *md, em_state_t *ems1, real a, gmx::ArrayRefWithPadding< const gmx::RVec > force, em_state_t *ems2, gmx::Constraints *constr, int64_t count)
 Do one minimization step.
 
static void em_dd_partition_system (FILE *fplog, const gmx::MDLogger &mdlog, int step, const t_commrec *cr, const gmx_mtop_t &top_global, const t_inputrec *ir, const MDModulesNotifiers &mdModulesNotifiers, gmx::ImdSession *imdSession, pull_t *pull_work, em_state_t *ems, gmx_localtop_t *top, gmx::MDAtoms *mdAtoms, t_forcerec *fr, VirtualSitesHandler *vsite, gmx::Constraints *constr, t_nrnb *nrnb, gmx_wallcycle *wcycle)
 Prepare EM for using domain decomposition parallellization.
 
void anonymous_namespace{minimize.cpp}::setCoordinates (std::vector< RVec > *coords, ArrayRef< const RVec > refCoords)
 Copy coordinates, OpenMP parallelized, from refCoords to coords.
 
real anonymous_namespace{minimize.cpp}::maxCoordinateDifference (ArrayRef< const RVec > coords1, ArrayRef< const RVec > coords2, MPI_Comm mpiCommMyGroup)
 Returns the maximum difference an atom moved between two coordinate sets, over all ranks.
 
static double reorder_partsum (const t_commrec *cr, const t_grpopts *opts, const gmx_mtop_t &top_global, const em_state_t *s_min, const em_state_t *s_b)
 Parallel utility summing energies and forces.
 
static real pr_beta (const t_commrec *cr, const t_grpopts *opts, t_mdatoms *mdatoms, const gmx_mtop_t &top_global, const em_state_t *s_min, const em_state_t *s_b)
 Print some stuff, like beta, whatever that means.