From: Kevin Boyd Date: Thu, 17 Jan 2019 03:54:11 +0000 (-0500) Subject: Remove references to ecutsGROUP X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=6eb4c2c6afb1873ed53b038efea66966117ddcab;p=alexxy%2Fgromacs.git Remove references to ecutsGROUP Changed warning that group scheme was deprecated to fatal error that group scheme no longer exists Modified warning when no cutoff scheme is specified Change-Id: If971f88003af46c9ccedb98bb2d33987848509ed --- diff --git a/src/gromacs/domdec/domdec_topology.cpp b/src/gromacs/domdec/domdec_topology.cpp index c43971003e..e91353175b 100644 --- a/src/gromacs/domdec/domdec_topology.cpp +++ b/src/gromacs/domdec/domdec_topology.cpp @@ -754,9 +754,7 @@ void dd_make_reverse_top(FILE *fplog, * kernels and only exclusions inside the cut-off lead to exclusion * forces. Since each atom pair is treated at most once in the non-bonded * kernels, it doesn't matter if the exclusions for the same atom pair - * appear multiple times in the exclusion list. In contrast, the, old, - * group cut-off scheme loops over a list of exclusions, so there each - * excluded pair should appear exactly once. + * appear multiple times in the exclusion list. */ rt->bExclRequired = (ir->cutoff_scheme == ecutsGROUP && inputrecExclForces(ir)); diff --git a/src/gromacs/ewald/pme.cpp b/src/gromacs/ewald/pme.cpp index 4ba4204237..1a9cc81276 100644 --- a/src/gromacs/ewald/pme.cpp +++ b/src/gromacs/ewald/pme.cpp @@ -202,10 +202,6 @@ bool pme_gpu_supports_input(const t_inputrec &ir, const gmx_mtop_t &mtop, std::s { errorReasons.emplace_back("Lennard-Jones PME"); } - if (ir.cutoff_scheme == ecutsGROUP) - { - errorReasons.emplace_back("group cutoff scheme"); - } if (!EI_DYNAMICS(ir.eI)) { errorReasons.emplace_back("not a dynamical integrator"); diff --git a/src/gromacs/ewald/pme_load_balancing.cpp b/src/gromacs/ewald/pme_load_balancing.cpp index c188be278d..500e351336 100644 --- a/src/gromacs/ewald/pme_load_balancing.cpp +++ b/src/gromacs/ewald/pme_load_balancing.cpp @@ -182,7 +182,6 @@ void pme_loadbal_init(pme_load_balancing_t **pme_lb_p, gmx_bool bUseGPU, gmx_bool *bPrinting) { - GMX_RELEASE_ASSERT(ir.cutoff_scheme != ecutsGROUP, "PME tuning is not supported with cutoff-scheme=group (because it contains bugs)"); pme_load_balancing_t *pme_lb; real spm, sp; diff --git a/src/gromacs/gmxpreprocess/grompp.cpp b/src/gromacs/gmxpreprocess/grompp.cpp index 1c8b99d808..553cd9f232 100644 --- a/src/gromacs/gmxpreprocess/grompp.cpp +++ b/src/gromacs/gmxpreprocess/grompp.cpp @@ -768,15 +768,6 @@ new_status(const char *topfile, const char *topppfile, const char *confin, warning(wi, buf); } - /* If using the group scheme, make sure charge groups are made whole to avoid errors - * in calculating charge group size later on - */ - if (ir->cutoff_scheme == ecutsGROUP && ir->ePBC != epbcNONE) - { - // Need temporary rvec for coordinates - do_pbc_first_mtop(nullptr, ir->ePBC, state->box, sys, state->x.rvec_array()); - } - /* Do more checks, mostly related to constraints */ if (bVerbose) { diff --git a/src/gromacs/gmxpreprocess/readir.cpp b/src/gromacs/gmxpreprocess/readir.cpp index 60401e169f..2bdcc5a809 100644 --- a/src/gromacs/gmxpreprocess/readir.cpp +++ b/src/gromacs/gmxpreprocess/readir.cpp @@ -291,26 +291,10 @@ void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, if (ir->cutoff_scheme == ecutsGROUP) { - warning_note(wi, - "The group cutoff scheme is deprecated since GROMACS 5.0 and will be removed in a future " - "release when all interaction forms are supported for the verlet scheme. The verlet " - "scheme already scales better, and it is compatible with GPUs and other accelerators."); - - if (ir->rlist > 0 && ir->rlist < ir->rcoulomb) - { - gmx_fatal(FARGS, "rcoulomb must not be greater than rlist (twin-range schemes are not supported)"); - } - if (ir->rlist > 0 && ir->rlist < ir->rvdw) - { - gmx_fatal(FARGS, "rvdw must not be greater than rlist (twin-range schemes are not supported)"); - } - - if (ir->rlist == 0 && ir->ePBC != epbcNONE) - { - warning_error(wi, "Can not have an infinite cut-off with PBC"); - } + gmx_fatal(FARGS, + "The group cutoff scheme has been removed since GROMACS 2020. " + "Please use the Verlet cutoff scheme."); } - if (ir->cutoff_scheme == ecutsVERLET) { real rc_max; @@ -1116,15 +1100,6 @@ void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, CHECK(ir->rcoulomb_switch >= ir->rcoulomb); } } - else if (ir->coulombtype == eelCUT || EEL_RF(ir->coulombtype)) - { - if (ir->cutoff_scheme == ecutsGROUP && ir->coulomb_modifier == eintmodNONE) - { - sprintf(err_buf, "With coulombtype = %s, rcoulomb should be >= rlist unless you use a potential modifier", - eel_names[ir->coulombtype]); - CHECK(ir->rlist > ir->rcoulomb); - } - } if (ir->coulombtype == eelSWITCH || ir->coulombtype == eelSHIFT) { @@ -1177,17 +1152,6 @@ void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, eel_names[ir->coulombtype]); CHECK(ir->rcoulomb > ir->rlist); } - else if (ir->cutoff_scheme == ecutsGROUP && ir->coulomb_modifier == eintmodNONE) - { - if (ir->coulombtype == eelPME || ir->coulombtype == eelP3M_AD) - { - sprintf(err_buf, - "With coulombtype = %s (without modifier), rcoulomb must be equal to rlist.\n" - "For optimal energy conservation,consider using\n" - "a potential modifier.", eel_names[ir->coulombtype]); - CHECK(ir->rcoulomb != ir->rlist); - } - } } if (EEL_PME(ir->coulombtype) || EVDW_PME(ir->vdwtype)) @@ -1250,14 +1214,6 @@ void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, warning_note(wi, warn_buf); } } - else if (ir->vdwtype == evdwCUT || ir->vdwtype == evdwPME) - { - if (ir->cutoff_scheme == ecutsGROUP && ir->vdw_modifier == eintmodNONE) - { - sprintf(err_buf, "With vdwtype = %s, rvdw must be >= rlist unless you use a potential modifier", evdw_names[ir->vdwtype]); - CHECK(ir->rlist > ir->rvdw); - } - } if (ir->vdwtype == evdwPME) { @@ -1271,29 +1227,6 @@ void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, } } - if (ir->cutoff_scheme == ecutsGROUP) - { - if (((ir->coulomb_modifier != eintmodNONE && ir->rcoulomb == ir->rlist) || - (ir->vdw_modifier != eintmodNONE && ir->rvdw == ir->rlist))) - { - warning_note(wi, "With exact cut-offs, rlist should be " - "larger than rcoulomb and rvdw, so that there " - "is a buffer region for particle motion " - "between neighborsearch steps"); - } - - if (ir_coulomb_is_zero_at_cutoff(ir) && ir->rlist <= ir->rcoulomb) - { - sprintf(warn_buf, "For energy conservation with switch/shift potentials, rlist should be 0.1 to 0.3 nm larger than rcoulomb."); - warning_note(wi, warn_buf); - } - if (ir_vdw_switched(ir) && (ir->rlist <= ir->rvdw)) - { - sprintf(warn_buf, "For energy conservation with switch/shift potentials, rlist should be 0.1 to 0.3 nm larger than rvdw."); - warning_note(wi, warn_buf); - } - } - if (ir->vdwtype == evdwUSER && ir->eDispCorr != edispcNO) { warning_note(wi, "You have selected user tables with dispersion correction, the dispersion will be corrected to -C6/r^6 beyond rvdw_switch (the tabulated interaction between rvdw_switch and rvdw will not be double counted). Make sure that you really want dispersion correction to -C6/r^6."); @@ -1310,23 +1243,6 @@ void check_ir(const char *mdparin, t_inputrec *ir, t_gromppopts *opts, warning(wi, "Using L-BFGS with nbfgscorr<=0 just gets you steepest descent."); } - /* ENERGY CONSERVATION */ - if (ir_NVE(ir) && ir->cutoff_scheme == ecutsGROUP) - { - if (!ir_vdw_might_be_zero_at_cutoff(ir) && ir->rvdw > 0 && ir->vdw_modifier == eintmodNONE) - { - sprintf(warn_buf, "You are using a cut-off for VdW interactions with NVE, for good energy conservation use vdwtype = %s (possibly with DispCorr)", - evdw_names[evdwSHIFT]); - warning_note(wi, warn_buf); - } - if (!ir_coulomb_might_be_zero_at_cutoff(ir) && ir->rcoulomb > 0) - { - sprintf(warn_buf, "You are using a cut-off for electrostatics with NVE, for good energy conservation use coulombtype = %s or %s", - eel_names[eelPMESWITCH], eel_names[eelRF_ZERO]); - warning_note(wi, warn_buf); - } - } - /* IMPLICIT SOLVENT */ if (ir->coulombtype == eelGB_NOTUSED) { @@ -1794,11 +1710,10 @@ void get_ir(const char *mdparin, const char *mdparout, { sprintf(warn_buf, "%s did not specify a value for the .mdp option " - "\"cutoff-scheme\". Probably it was first intended for use " - "with GROMACS before 4.6. In 4.6, the Verlet scheme was " - "introduced, but the group scheme was still the default. " - "The default is now the Verlet scheme, so you will observe " - "different behaviour.", mdparin); + "\"cutoff-scheme\". As of GROMACS 2020, the Verlet scheme " + "is the only cutoff scheme supported. This may affect your " + "simulation if you are using an old mdp file that assumes use " + "of the (removed) group cutoff scheme.", mdparin); warning_note(wi, warn_buf); } diff --git a/src/gromacs/listed_forces/gpubonded_impl.cpp b/src/gromacs/listed_forces/gpubonded_impl.cpp index b7b20fdd52..328b225fa5 100644 --- a/src/gromacs/listed_forces/gpubonded_impl.cpp +++ b/src/gromacs/listed_forces/gpubonded_impl.cpp @@ -143,10 +143,6 @@ bool inputSupportsGpuBondeds(const t_inputrec &ir, { errorReasons.emplace_back("No supported bonded interactions are present"); } - if (ir.cutoff_scheme == ecutsGROUP) - { - errorReasons.emplace_back("group cutoff scheme"); - } if (!EI_DYNAMICS(ir.eI)) { errorReasons.emplace_back("not a dynamical integrator"); diff --git a/src/gromacs/mdlib/forcerec.cpp b/src/gromacs/mdlib/forcerec.cpp index 3c625533fb..0814dac43e 100644 --- a/src/gromacs/mdlib/forcerec.cpp +++ b/src/gromacs/mdlib/forcerec.cpp @@ -996,78 +996,6 @@ static real calcBuckinghamBMax(FILE *fplog, const gmx_mtop_t *mtop) return bham_b_max; } -static void make_nbf_tables(FILE *fp, - const interaction_const_t *ic, real rtab, - const char *tabfn, char *eg1, char *eg2, - t_nblists *nbl) -{ - char buf[STRLEN]; - int i, j; - - if (tabfn == nullptr) - { - if (debug) - { - fprintf(debug, "No table file name passed, can not read table, can not do non-bonded interactions\n"); - } - return; - } - - sprintf(buf, "%s", tabfn); - if (eg1 && eg2) - { - /* Append the two energy group names */ - sprintf(buf + strlen(tabfn) - strlen(ftp2ext(efXVG)) - 1, "_%s_%s.%s", - eg1, eg2, ftp2ext(efXVG)); - } - nbl->table_elec_vdw = make_tables(fp, ic, buf, rtab, 0); - /* Copy the contents of the table to separate coulomb and LJ tables too, - * to improve cache performance. - */ - /* For performance reasons we want - * the table data to be aligned to 16-byte. - */ - nbl->table_elec = - new t_forcetable(GMX_TABLE_INTERACTION_ELEC, - nbl->table_elec_vdw->format); - nbl->table_elec->r = nbl->table_elec_vdw->r; - nbl->table_elec->n = nbl->table_elec_vdw->n; - nbl->table_elec->scale = nbl->table_elec_vdw->scale; - nbl->table_elec->formatsize = nbl->table_elec_vdw->formatsize; - nbl->table_elec->ninteractions = 1; - nbl->table_elec->stride = nbl->table_elec->formatsize * nbl->table_elec->ninteractions; - snew_aligned(nbl->table_elec->data, nbl->table_elec->stride*(nbl->table_elec->n+1), 32); - - nbl->table_vdw = - new t_forcetable(GMX_TABLE_INTERACTION_VDWREP_VDWDISP, - nbl->table_elec_vdw->format); - nbl->table_vdw->r = nbl->table_elec_vdw->r; - nbl->table_vdw->n = nbl->table_elec_vdw->n; - nbl->table_vdw->scale = nbl->table_elec_vdw->scale; - nbl->table_vdw->formatsize = nbl->table_elec_vdw->formatsize; - nbl->table_vdw->ninteractions = 2; - nbl->table_vdw->stride = nbl->table_vdw->formatsize * nbl->table_vdw->ninteractions; - snew_aligned(nbl->table_vdw->data, nbl->table_vdw->stride*(nbl->table_vdw->n+1), 32); - - /* NOTE: Using a single i-loop here leads to mix-up of data in table_vdw - * with (at least) gcc 6.2, 6.3 and 6.4 when compiled with -O3 and AVX - */ - for (i = 0; i <= nbl->table_elec_vdw->n; i++) - { - for (j = 0; j < 4; j++) - { - nbl->table_elec->data[4*i+j] = nbl->table_elec_vdw->data[12*i+j]; - } - } - for (i = 0; i <= nbl->table_elec_vdw->n; i++) - { - for (j = 0; j < 8; j++) - { - nbl->table_vdw->data[8*i+j] = nbl->table_elec_vdw->data[12*i+4+j]; - } - } -} - /*!\brief If there's bonded interactions of type \c ftype1 or \c * ftype2 present in the topology, build an array of the number of * interactions present for each bonded interaction index found in the @@ -1316,26 +1244,6 @@ static void initVdwEwaldParameters(FILE *fp, const t_inputrec *ir, } } -gmx_bool uses_simple_tables(int cutoff_scheme, - const nonbonded_verlet_t *nbv) -{ - gmx_bool bUsesSimpleTables = TRUE; - - switch (cutoff_scheme) - { - case ecutsGROUP: - bUsesSimpleTables = TRUE; - break; - case ecutsVERLET: - GMX_RELEASE_ASSERT(nullptr != nbv, "A non-bonded verlet object is required with the Verlet cutoff-scheme"); - bUsesSimpleTables = nbv->pairlistIsSimple(); - break; - default: - gmx_incons("unimplemented"); - } - return bUsesSimpleTables; -} - static void init_ewald_f_table(interaction_const_t *ic, real rtab) { @@ -1540,12 +1448,6 @@ init_interaction_const(FILE *fp, ic->rcoulomb, 0, 0, nullptr, &ic->k_rf, &ic->c_rf); } - - if (ir->cutoff_scheme == ecutsGROUP && ic->eeltype == eelRF_ZERO) - { - /* grompp should have done this, but this scheme is obsolete */ - ic->coulomb_modifier = eintmodEXACTCUTOFF; - } } else { @@ -1616,16 +1518,12 @@ void init_forcerec(FILE *fp, gmx_bool bNoSolvOpt, real print_force) { - int m, negp_pp, negptable, egi, egj; real rtab; char *env; double dbl; const t_block *cgs; gmx_bool bGenericKernelOnly; - gmx_bool needGroupSchemeTables, bSomeNormalNbListsAreInUse; gmx_bool bFEP_NonBonded; - int egp_flags; - /* By default we turn SIMD kernels on, but it might be turned off further down... */ fr->use_simd_kernels = TRUE; @@ -1763,21 +1661,6 @@ void init_forcerec(FILE *fp, fr->bGrid = (ir->ns_type == ensGRID); fr->ePBC = ir->ePBC; - if (fr->cutoff_scheme == ecutsGROUP) - { - const char *note = "NOTE: This file uses the deprecated 'group' cutoff_scheme. This will be\n" - "removed in a future release when 'verlet' supports all interaction forms.\n"; - - if (MASTER(cr)) - { - fprintf(stderr, "\n%s\n", note); - } - if (fp != nullptr) - { - fprintf(fp, "\n%s\n", note); - } - } - /* Determine if we will do PBC for distances in bonded interactions */ if (fr->ePBC == epbcNONE) { @@ -1798,8 +1681,7 @@ void init_forcerec(FILE *fp, * With intermolecular interactions we need PBC for calculating * distances between atoms in different molecules. */ - if ((fr->cutoff_scheme == ecutsGROUP || bSHAKE) && - !mtop->bIntermolecularInteractions) + if (bSHAKE && !mtop->bIntermolecularInteractions) { fr->bMolPBC = ir->bPeriodicMols; @@ -1931,80 +1813,6 @@ void init_forcerec(FILE *fp, } fr->nbkernel_vdw_modifier = ic->vdw_modifier; - if (ir->cutoff_scheme == ecutsGROUP) - { - fr->bvdwtab = ((ic->vdwtype != evdwCUT || !gmx_within_tol(ic->reppow, 12.0, 10*GMX_DOUBLE_EPS)) - && !EVDW_PME(ic->vdwtype)); - /* We have special kernels for standard Ewald and PME, but the pme-switch ones are tabulated above */ - fr->bcoultab = !(ic->eeltype == eelCUT || - ic->eeltype == eelEWALD || - ic->eeltype == eelPME || - ic->eeltype == eelP3M_AD || - ic->eeltype == eelRF || - ic->eeltype == eelRF_ZERO); - - /* If the user absolutely wants different switch/shift settings for coul/vdw, it is likely - * going to be faster to tabulate the interaction than calling the generic kernel. - * However, if generic kernels have been requested we keep things analytically. - */ - if (fr->nbkernel_elec_modifier == eintmodPOTSWITCH && - fr->nbkernel_vdw_modifier == eintmodPOTSWITCH && - !bGenericKernelOnly) - { - if ((ic->rcoulomb_switch != ic->rvdw_switch) || (ic->rcoulomb != ic->rvdw)) - { - fr->bcoultab = TRUE; - /* Once we tabulate electrostatics, we can use the switch function for LJ, - * which would otherwise need two tables. - */ - } - } - else if ((fr->nbkernel_elec_modifier == eintmodPOTSHIFT && fr->nbkernel_vdw_modifier == eintmodPOTSHIFT) || - ((fr->nbkernel_elec_interaction == GMX_NBKERNEL_ELEC_REACTIONFIELD && - fr->nbkernel_elec_modifier == eintmodEXACTCUTOFF && - (fr->nbkernel_vdw_modifier == eintmodPOTSWITCH || fr->nbkernel_vdw_modifier == eintmodPOTSHIFT)))) - { - if ((ic->rcoulomb != ic->rvdw) && (!bGenericKernelOnly)) - { - fr->bcoultab = TRUE; - } - } - - if (fr->nbkernel_elec_modifier == eintmodFORCESWITCH) - { - fr->bcoultab = TRUE; - } - if (fr->nbkernel_vdw_modifier == eintmodFORCESWITCH) - { - fr->bvdwtab = TRUE; - } - - if (getenv("GMX_REQUIRE_TABLES")) - { - fr->bvdwtab = TRUE; - fr->bcoultab = TRUE; - } - - if (fp) - { - fprintf(fp, "Table routines are used for coulomb: %s\n", - gmx::boolToString(fr->bcoultab)); - fprintf(fp, "Table routines are used for vdw: %s\n", - gmx::boolToString(fr->bvdwtab)); - } - - if (fr->bvdwtab) - { - fr->nbkernel_vdw_interaction = GMX_NBKERNEL_VDW_CUBICSPLINETABLE; - fr->nbkernel_vdw_modifier = eintmodNONE; - } - if (fr->bcoultab) - { - fr->nbkernel_elec_interaction = GMX_NBKERNEL_ELEC_CUBICSPLINETABLE; - fr->nbkernel_elec_modifier = eintmodNONE; - } - } - if (ir->cutoff_scheme == ecutsVERLET) { if (!gmx_within_tol(ic->reppow, 12.0, 10*GMX_DOUBLE_EPS)) @@ -2116,41 +1924,6 @@ void init_forcerec(FILE *fp, gmx_fatal(FARGS, "Implict solvation is no longer supported."); } - /* Construct tables for the group scheme. A little unnecessary to - * make both vdw and coul tables sometimes, but what the - * heck. Note that both cutoff schemes construct Ewald tables in - * init_interaction_const_tables. */ - needGroupSchemeTables = (ir->cutoff_scheme == ecutsGROUP && - (fr->bcoultab || fr->bvdwtab)); - - negp_pp = ir->opts.ngener - ir->nwall; - negptable = 0; - if (!needGroupSchemeTables) - { - bSomeNormalNbListsAreInUse = TRUE; - } - else - { - bSomeNormalNbListsAreInUse = FALSE; - for (egi = 0; egi < negp_pp; egi++) - { - for (egj = egi; egj < negp_pp; egj++) - { - egp_flags = ir->opts.egp_flags[GID(egi, egj, ir->opts.ngener)]; - if (!(egp_flags & EGP_EXCL)) - { - if (egp_flags & EGP_TABLE) - { - negptable++; - } - else - { - bSomeNormalNbListsAreInUse = TRUE; - } - } - } - } - } /* This code automatically gives table length tabext without cut-off's, * in that case grompp should already have checked that we do not need @@ -2158,41 +1931,6 @@ void init_forcerec(FILE *fp, */ rtab = ir->rlist + ir->tabext; - if (needGroupSchemeTables) - { - /* make tables for ordinary interactions */ - if (bSomeNormalNbListsAreInUse) - { - make_nbf_tables(fp, ic, rtab, tabfn, nullptr, nullptr, nullptr); - m = 1; - } - else - { - m = 0; - } - if (negptable > 0) - { - /* Read the special tables for certain energy group pairs */ - gmx::ArrayRef nm_ind = mtop->groups.groups[SimulationAtomGroupType::EnergyOutput]; - for (egi = 0; egi < negp_pp; egi++) - { - for (egj = egi; egj < negp_pp; egj++) - { - egp_flags = ir->opts.egp_flags[GID(egi, egj, ir->opts.ngener)]; - if ((egp_flags & EGP_TABLE) && !(egp_flags & EGP_EXCL)) - { - /* Read the table file with the two energy groups names appended */ - make_nbf_tables(fp, ic, rtab, tabfn, - *mtop->groups.groupNames[nm_ind[egi]], - *mtop->groups.groupNames[nm_ind[egj]], - nullptr); - m++; - } - } - } - } - } - /* We want to use unmodified tables for 1-4 coulombic * interactions, so we must in general have an extra set of * tables. */ diff --git a/src/gromacs/mdlib/forcerec.h b/src/gromacs/mdlib/forcerec.h index 51796adf29..52023acd88 100644 --- a/src/gromacs/mdlib/forcerec.h +++ b/src/gromacs/mdlib/forcerec.h @@ -149,12 +149,6 @@ void forcerec_set_excl_load(t_forcerec *fr, */ void update_forcerec(t_forcerec *fr, matrix box); -gmx_bool uses_simple_tables(int cutoff_scheme, - const nonbonded_verlet_t *nbv); -/* Returns whether simple tables (i.e. not for use with GPUs) are used - * with the type of kernel indicated. - */ - void free_gpu_resources(t_forcerec *fr, const gmx::PhysicalNodeCommunicator &physicalNodeCommunicator); diff --git a/src/gromacs/mdlib/perf_est.cpp b/src/gromacs/mdlib/perf_est.cpp index 8c5bc33de8..a66b119536 100644 --- a/src/gromacs/mdlib/perf_est.cpp +++ b/src/gromacs/mdlib/perf_est.cpp @@ -63,18 +63,6 @@ * that the numbers need to be adjusted. */ -/* Cost of a pair interaction in the "group" cut-off scheme */ -static const double c_group_fq = 18.0; -static const double c_group_qlj_cut = 18.0; -static const double c_group_qlj_tab = 24.0; -static const double c_group_lj_cut = 12.0; -static const double c_group_lj_tab = 21.0; -/* Cost of 1 water with one Q/LJ atom */ -static const double c_group_qljw_cut = 24.0; -static const double c_group_qljw_tab = 27.0; -/* Cost of 1 water with one Q atom or with 1/3 water (LJ negligible) */ -static const double c_group_qw = 21.0; - /* Cost of a pair interaction in the "Verlet" cut-off scheme, QEXP is Ewald */ static const double c_nbnxn_lj = 2.5; static const double c_nbnxn_qrf_lj = 2.9; @@ -271,129 +259,6 @@ void count_bonded_distances(const gmx_mtop_t *mtop, const t_inputrec *ir, } } -static void pp_group_load(const gmx_mtop_t *mtop, const t_inputrec *ir, - const matrix box, - int *nq_tot, int *nlj_tot, - double *cost_pp, - gmx_bool *bChargePerturbed, gmx_bool *bTypePerturbed) -{ - int atnr, cg, a0, ncqlj, ncq, nclj; - gmx_bool bBHAM, bLJcut, bWater, bQ, bLJ; - int nw, nqlj, nq, nlj; - double fq, fqlj, flj, fqljw, fqw; - - bBHAM = (mtop->ffparams.functype[0] == F_BHAM); - - bLJcut = ((ir->vdwtype == evdwCUT) && !bBHAM); - - /* Computational cost of bonded, non-bonded and PME calculations. - * This will be machine dependent. - * The numbers here are accurate for Intel Core2 and AMD Athlon 64 - * in single precision. In double precision PME mesh is slightly cheaper, - * although not so much that the numbers need to be adjusted. - */ - fq = c_group_fq; - fqlj = (bLJcut ? c_group_qlj_cut : c_group_qlj_tab); - flj = (bLJcut ? c_group_lj_cut : c_group_lj_tab); - /* Cost of 1 water with one Q/LJ atom */ - fqljw = (bLJcut ? c_group_qljw_cut : c_group_qljw_tab); - /* Cost of 1 water with one Q atom or with 1/3 water (LJ negligible) */ - fqw = c_group_qw; - - gmx::ArrayRef iparams = mtop->ffparams.iparams; - atnr = mtop->ffparams.atnr; - nw = 0; - nqlj = 0; - nq = 0; - nlj = 0; - *bChargePerturbed = FALSE; - *bTypePerturbed = FALSE; - for (const gmx_molblock_t &molb : mtop->molblock) - { - const gmx_moltype_t *molt = &mtop->moltype[molb.type]; - const t_atom *atom = molt->atoms.atom; - int a = 0; - for (cg = 0; cg < molt->cgs.nr; cg++) - { - bWater = !bBHAM; - ncqlj = 0; - ncq = 0; - nclj = 0; - a0 = a; - while (a < molt->cgs.index[cg+1]) - { - bQ = (atom[a].q != 0 || atom[a].qB != 0); - bLJ = (iparams[(atnr+1)*atom[a].type].lj.c6 != 0 || - iparams[(atnr+1)*atom[a].type].lj.c12 != 0); - if (atom[a].q != atom[a].qB) - { - *bChargePerturbed = TRUE; - } - if (atom[a].type != atom[a].typeB) - { - *bTypePerturbed = TRUE; - } - /* This if this atom fits into water optimization */ - if (!((a == a0 && bQ && bLJ) || - (a == a0+1 && bQ && !bLJ) || - (a == a0+2 && bQ && !bLJ && atom[a].q == atom[a-1].q) || - (a == a0+3 && !bQ && bLJ))) - { - bWater = FALSE; - } - if (bQ && bLJ) - { - ncqlj++; - } - else - { - if (bQ) - { - ncq++; - } - if (bLJ) - { - nclj++; - } - } - a++; - } - if (bWater) - { - nw += molb.nmol; - } - else - { - nqlj += molb.nmol*ncqlj; - nq += molb.nmol*ncq; - nlj += molb.nmol*nclj; - } - } - } - - *nq_tot = nq + nqlj + nw*3; - *nlj_tot = nlj + nqlj + nw; - - if (debug) - { - fprintf(debug, "nw %d nqlj %d nq %d nlj %d\n", nw, nqlj, nq, nlj); - } - - /* For the PP non-bonded cost it is (unrealistically) assumed - * that all atoms are distributed homogeneously in space. - * Factor 3 is used because a water molecule has 3 atoms - * (and TIP4P effectively has 3 interactions with (water) atoms)). - */ - *cost_pp = 0.5*(fqljw*nw*nqlj + - fqw *nw*(3*nw + nq) + - fqlj *nqlj*nqlj + - fq *nq*(3*nw + nqlj + nq) + - flj *nlj*(nw + nqlj + nlj)) - *4/3*M_PI*ir->rlist*ir->rlist*ir->rlist/det(box); - - *cost_pp *= simd_cycle_factor(bHaveSIMD); -} - static void pp_verlet_load(const gmx_mtop_t *mtop, const t_inputrec *ir, const matrix box, int *nq_tot, int *nlj_tot, @@ -526,18 +391,9 @@ float pme_load_estimate(const gmx_mtop_t *mtop, const t_inputrec *ir, cost_bond = c_bond*(ndistance_c *simd_cycle_factor(FALSE) + ndistance_simd*simd_cycle_factor(bHaveSIMD)); - if (ir->cutoff_scheme == ecutsGROUP) - { - pp_group_load(mtop, ir, box, - &nq_tot, &nlj_tot, &cost_pp, - &bChargePerturbed, &bTypePerturbed); - } - else - { - pp_verlet_load(mtop, ir, box, - &nq_tot, &nlj_tot, &cost_pp, - &bChargePerturbed, &bTypePerturbed); - } + pp_verlet_load(mtop, ir, box, + &nq_tot, &nlj_tot, &cost_pp, + &bChargePerturbed, &bTypePerturbed); cost_redist = 0; cost_spread = 0;