From b1f5c9dd69b9a854ab392aa022470de403603160 Mon Sep 17 00:00:00 2001 From: ejjordan Date: Thu, 18 Feb 2021 12:59:01 +0100 Subject: [PATCH] Pass restraint data structures to bonded functions This adds explicitly the t_disresdata and t_oriresdata structures to all bonded force calculation calls, instead of passing them as members of the t_fcdata structure. This is an unfortunate but necessary step in order to be able to separate bondedtable_t from the two restraint data structures. The splitting of t_fcdata into separate table data and restraint data structures is needed to simplify the initialization of the forcerec. --- src/gromacs/gmxana/gmx_disre.cpp | 8 +- src/gromacs/listed_forces/bonded.cpp | 82 ++++++++++++++++++++- src/gromacs/listed_forces/bonded.h | 8 +- src/gromacs/listed_forces/disre.cpp | 72 +++++++++--------- src/gromacs/listed_forces/disre.h | 17 +++-- src/gromacs/listed_forces/listed_forces.cpp | 4 + src/gromacs/listed_forces/listed_forces.h | 4 + src/gromacs/listed_forces/orires.cpp | 32 ++++---- src/gromacs/listed_forces/orires.h | 4 + src/gromacs/listed_forces/tests/bonded.cpp | 2 + 10 files changed, 165 insertions(+), 68 deletions(-) diff --git a/src/gromacs/gmxana/gmx_disre.cpp b/src/gromacs/gmxana/gmx_disre.cpp index 8b8eade002..bd3feef043 100644 --- a/src/gromacs/gmxana/gmx_disre.cpp +++ b/src/gromacs/gmxana/gmx_disre.cpp @@ -4,7 +4,7 @@ * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. * Copyright (c) 2013,2014,2015,2016,2017 The GROMACS development team. - * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -205,10 +205,6 @@ static void check_viol(FILE* log, } } - // Set up t_fcdata, only needed for calling ta_disres() - t_fcdata fcd; - fcd.disres = disresdata; - // Get offset for label index label_old = forceparams[forceatoms[0]].disres.label; for (i = 0; (i < disres.size());) @@ -241,7 +237,7 @@ static void check_viol(FILE* log, dr[clust_id].aver_6[ndr] += disresdata->Rt_6[label]; snew(fshift, SHIFTS); - ta_disres(n, &forceatoms[i], forceparams.data(), x, f, fshift, pbc, lam, &dvdl, nullptr, &fcd, nullptr); + ta_disres(n, &forceatoms[i], forceparams.data(), x, f, fshift, pbc, lam, &dvdl, nullptr, nullptr, disresdata, nullptr, nullptr); sfree(fshift); viol = disresdata->sumviol; diff --git a/src/gromacs/listed_forces/bonded.cpp b/src/gromacs/listed_forces/bonded.cpp index d5b939696b..80ae6741de 100644 --- a/src/gromacs/listed_forces/bonded.cpp +++ b/src/gromacs/listed_forces/bonded.cpp @@ -4,7 +4,7 @@ * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. * Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team. - * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -102,6 +102,8 @@ using BondedFunction = real (*)(int nbonds, real* dvdlambda, const t_mdatoms* md, t_fcdata* fcd, + t_disresdata* disresdata, + t_oriresdata* oriresdata, int* ddgatindex); /*! \brief Mysterious CMAP coefficient matrix */ @@ -256,6 +258,8 @@ real morse_bonds(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { const real one = 1.0; @@ -324,6 +328,8 @@ real cubic_bonds(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { const real three = 3.0; @@ -379,7 +385,9 @@ real FENE_bonds(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, - int* global_atom_index) + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, + int* global_atom_index) { const real half = 0.5; const real one = 1.0; @@ -468,6 +476,8 @@ bonds(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ki, ai, aj, type; @@ -528,6 +538,8 @@ bonds(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { constexpr int nfa1 = 3; @@ -618,6 +630,8 @@ real restraint_bonds(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ki, ai, aj, type; @@ -706,6 +720,8 @@ real polarize(int nbonds, real* dvdlambda, const t_mdatoms* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ki, ai, aj, type; @@ -751,6 +767,8 @@ real anharm_polarize(int nbonds, real* dvdlambda, const t_mdatoms* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ki, ai, aj, type; @@ -805,6 +823,8 @@ real water_pol(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { /* This routine implements anisotropic polarizibility for water, through @@ -948,6 +968,8 @@ real thole_pol(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { /* Interaction between two pairs of particles with opposite charge */ @@ -999,6 +1021,8 @@ angles(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ai, aj, ak, t1, t2, type; @@ -1091,6 +1115,8 @@ angles(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { const int nfa1 = 4; @@ -1249,6 +1275,8 @@ real linear_angles(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, m, ai, aj, ak, t1, t2, type; @@ -1319,6 +1347,8 @@ urey_bradley(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, m, ai, aj, ak, t1, t2, type, ki; @@ -1429,6 +1459,8 @@ urey_bradley(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { constexpr int nfa1 = 4; @@ -1571,6 +1603,8 @@ real quartic_angles(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, j, ai, aj, ak, t1, t2, type; @@ -1936,6 +1970,8 @@ pdihs(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int t1, t2, t3; @@ -1997,6 +2033,8 @@ pdihs(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { const int nfa1 = 5; @@ -2111,6 +2149,8 @@ rbdihs(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { const int nfa1 = 5; @@ -2237,6 +2277,8 @@ real idihs(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, type, ai, aj, ak, al; @@ -2403,6 +2445,8 @@ real angres(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { return low_angres(nbonds, forceatoms, forceparams, x, f, fshift, pbc, lambda, dvdlambda, FALSE); @@ -2420,6 +2464,8 @@ real angresz(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { return low_angres(nbonds, forceatoms, forceparams, x, f, fshift, pbc, lambda, dvdlambda, TRUE); @@ -2437,6 +2483,8 @@ real dihres(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { real vtot = 0; @@ -2530,6 +2578,8 @@ real unimplemented(int gmx_unused nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { gmx_impl("*** you are using a not implemented function"); @@ -2547,6 +2597,8 @@ real restrangles(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, d, ai, aj, ak, type, m; @@ -2649,6 +2701,8 @@ real restrdihs(int nbonds, real gmx_unused* dvlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, d, type, ai, aj, ak, al; @@ -2771,6 +2825,8 @@ real cbtdihs(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int type, ai, aj, ak, al, i, d; @@ -2886,6 +2942,8 @@ rbdihs(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { const real c0 = 0.0, c1 = 1.0, c2 = 2.0, c3 = 3.0, c4 = 4.0, c5 = 5.0; @@ -3036,6 +3094,8 @@ real cmap_dihs(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, n; @@ -3440,6 +3500,8 @@ real g96bonds(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ki, ai, aj, type; @@ -3497,6 +3559,8 @@ real g96angles(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ai, aj, ak, type, m, t1, t2; @@ -3564,6 +3628,8 @@ real cross_bond_bond(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { /* Potential from Lawrence and Skimmer, Chem. Phys. Lett. 372 (2003) @@ -3636,6 +3702,8 @@ real cross_bond_angle(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { /* Potential from Lawrence and Skimmer, Chem. Phys. Lett. 372 (2003) @@ -3772,6 +3840,8 @@ real tab_bonds(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ki, ai, aj, type, table; @@ -3827,6 +3897,8 @@ real tab_angles(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, ai, aj, ak, t1, t2, type, table; @@ -3907,6 +3979,8 @@ real tab_dihs(int nbonds, real* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { int i, type, ai, aj, ak, al, table; @@ -4083,13 +4157,15 @@ real calculateSimpleBond(const int ftype, real* dvdlambda, const t_mdatoms* md, t_fcdata* fcd, + t_disresdata* disresdata, + t_oriresdata* oriresdata, int gmx_unused* global_atom_index, const BondedKernelFlavor bondedKernelFlavor) { const BondedInteractions& bonded = c_bondedInteractionFunctionsPerFlavor[bondedKernelFlavor][ftype]; real v = bonded.function( - numForceatoms, forceatoms, forceparams, x, f, fshift, pbc, lambda, dvdlambda, md, fcd, global_atom_index); + numForceatoms, forceatoms, forceparams, x, f, fshift, pbc, lambda, dvdlambda, md, fcd, disresdata, oriresdata, global_atom_index); return v; } diff --git a/src/gromacs/listed_forces/bonded.h b/src/gromacs/listed_forces/bonded.h index beb95b5528..d3c006145e 100644 --- a/src/gromacs/listed_forces/bonded.h +++ b/src/gromacs/listed_forces/bonded.h @@ -4,7 +4,7 @@ * Copyright (c) 1991-2000, University of Groningen, The Netherlands. * Copyright (c) 2001-2004, The GROMACS development team. * Copyright (c) 2013,2014,2015,2016,2017 by the GROMACS development team. - * Copyright (c) 2018,2019,2020, by the GROMACS development team, led by + * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -59,6 +59,8 @@ struct t_fcdata; struct t_mdatom; struct t_nrnb; struct t_pbc; +struct t_disresdata; +struct t_oriresdata; namespace gmx { @@ -127,6 +129,8 @@ real cmap_dihs(int nbonds, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index); /*! \brief For selecting which flavor of bonded kernel is used for simple bonded types */ @@ -180,6 +184,8 @@ real calculateSimpleBond(int ftype, real* dvdlambda, const t_mdatoms* md, t_fcdata* fcd, + t_disresdata* disresdata, + t_oriresdata* oriresdata, int gmx_unused* global_atom_index, BondedKernelFlavor bondedKernelFlavor); diff --git a/src/gromacs/listed_forces/disre.cpp b/src/gromacs/listed_forces/disre.cpp index d08261a7a8..1095ab94df 100644 --- a/src/gromacs/listed_forces/disre.cpp +++ b/src/gromacs/listed_forces/disre.cpp @@ -415,48 +415,48 @@ void calc_disres_R_6(const t_commrec* cr, dd->sumviol = 0; } -real ta_disres(int nfa, - const t_iatom forceatoms[], - const t_iparams ip[], - const rvec x[], - rvec4 f[], - rvec fshift[], - const t_pbc* pbc, +real ta_disres(int nfa, + const t_iatom* forceatoms, + const t_iparams* ip, + const rvec* x, + rvec4* f, + rvec* fshift, + const t_pbc* pbc, real gmx_unused lambda, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, - t_fcdata* fcd, + t_fcdata gmx_unused* fcd, + t_disresdata* disresdata, + t_oriresdata gmx_unused* oriresdata, int gmx_unused* global_atom_index) { const real seven_three = 7.0 / 3.0; - rvec dx; - real weight_rt_1; - real smooth_fc, Rt, Rtav, rt2, *Rtl_6, *Rt_6, *Rtav_6; - real k0, f_scal = 0, fmax_scal, fk_scal, fij; - real tav_viol, instant_viol, mixed_viol, violtot, vtot; - real tav_viol_Rtav7, instant_viol_Rtav7; - real up1, up2, low; - gmx_bool bConservative, bMixed, bViolation; - t_disresdata* dd; - int dr_weighting; - gmx_bool dr_bMixed; - - dd = fcd->disres; - dr_weighting = dd->dr_weighting; - dr_bMixed = dd->dr_bMixed; - Rtl_6 = dd->Rtl_6; - Rt_6 = dd->Rt_6; - Rtav_6 = dd->Rtav_6; + rvec dx; + real weight_rt_1; + real smooth_fc, Rt, Rtav, rt2, *Rtl_6, *Rt_6, *Rtav_6; + real k0, f_scal = 0, fmax_scal, fk_scal, fij; + real tav_viol, instant_viol, mixed_viol, violtot, vtot; + real tav_viol_Rtav7, instant_viol_Rtav7; + real up1, up2, low; + gmx_bool bConservative, bMixed, bViolation; + int dr_weighting; + gmx_bool dr_bMixed; + + dr_weighting = disresdata->dr_weighting; + dr_bMixed = disresdata->dr_bMixed; + Rtl_6 = disresdata->Rtl_6; + Rt_6 = disresdata->Rt_6; + Rtav_6 = disresdata->Rtav_6; tav_viol = instant_viol = mixed_viol = tav_viol_Rtav7 = instant_viol_Rtav7 = 0; - smooth_fc = dd->dr_fc; - if (dd->dr_tau != 0) + smooth_fc = disresdata->dr_fc; + if (disresdata->dr_tau != 0) { /* scaling factor to smoothly turn on the restraint forces * * when using time averaging */ - smooth_fc *= (1.0 - dd->exp_min_t_tau); + smooth_fc *= (1.0 - disresdata->exp_min_t_tau); } violtot = 0; @@ -464,7 +464,7 @@ real ta_disres(int nfa, /* 'loop' over all atom pairs (pair_nr=fa/3) involved in restraints, * * the total number of atoms pairs is nfa/3 */ - int faOffset = static_cast(forceatoms - dd->forceatomsStart); + int faOffset = static_cast(forceatoms - disresdata->forceatomsStart); for (int fa = 0; fa < nfa; fa += 3) { int type = forceatoms[fa]; @@ -474,7 +474,7 @@ real ta_disres(int nfa, low = ip[type].disres.low; k0 = smooth_fc * ip[type].disres.kfac; - int res = type - dd->type_min; + int res = type - disresdata->type_min; /* save some flops when there is only one pair */ if (ip[type].disres.type != 2) @@ -605,12 +605,14 @@ real ta_disres(int nfa, { if (!dr_bMixed) { - weight_rt_1 *= std::pow(dd->rm3tav[pair], seven_three); + weight_rt_1 *= std::pow(disresdata->rm3tav[pair], seven_three); } else { - weight_rt_1 *= tav_viol_Rtav7 * std::pow(dd->rm3tav[pair], seven_three) - + instant_viol_Rtav7 / (dd->rt[pair] * gmx::power6(dd->rt[pair])); + weight_rt_1 *= + tav_viol_Rtav7 * std::pow(disresdata->rm3tav[pair], seven_three) + + instant_viol_Rtav7 + / (disresdata->rt[pair] * gmx::power6(disresdata->rt[pair])); } } @@ -632,7 +634,7 @@ real ta_disres(int nfa, } #pragma omp atomic - dd->sumviol += violtot; + disresdata->sumviol += violtot; /* Return energy */ return vtot; diff --git a/src/gromacs/listed_forces/disre.h b/src/gromacs/listed_forces/disre.h index cbe81d778f..756cce341c 100644 --- a/src/gromacs/listed_forces/disre.h +++ b/src/gromacs/listed_forces/disre.h @@ -56,6 +56,7 @@ struct gmx_multisim_t; class history_t; struct t_commrec; struct t_disresdata; +struct t_oriresdata; struct t_fcdata; struct t_inputrec; struct t_pbc; @@ -107,17 +108,19 @@ void calc_disres_R_6(const t_commrec* cr, //! Calculates the distance restraint forces, return the potential. real ta_disres(int nfa, - const t_iatom forceatoms[], - const t_iparams ip[], - const rvec x[], - rvec4 f[], - rvec fshift[], + const t_iatom* forceatoms, + const t_iparams* ip, + const rvec* x, + rvec4* f, + rvec* fshift, const t_pbc* pbc, real lambda, real* dvdlambda, const t_mdatoms* md, - t_fcdata* fcdata, - int* global_atom_index); + t_fcdata gmx_unused* fcd, + t_disresdata* disresdata, + t_oriresdata gmx_unused* oriresdata, + int* global_atom_index); //! Copies the new time averages that have been calculated in calc_disres_R_6. void update_disres_history(const t_disresdata& disresdata, history_t* hist); diff --git a/src/gromacs/listed_forces/listed_forces.cpp b/src/gromacs/listed_forces/listed_forces.cpp index d855b2e096..10e54fb419 100644 --- a/src/gromacs/listed_forces/listed_forces.cpp +++ b/src/gromacs/listed_forces/listed_forces.cpp @@ -470,6 +470,8 @@ real calc_one_bond(int thread, &(dvdl[efptFTYPE]), md, fcd, + nullptr, + nullptr, global_atom_index); } else @@ -486,6 +488,8 @@ real calc_one_bond(int thread, &(dvdl[efptFTYPE]), md, fcd, + fcd->disres, + fcd->orires, global_atom_index, flavor); } diff --git a/src/gromacs/listed_forces/listed_forces.h b/src/gromacs/listed_forces/listed_forces.h index 4bed406827..d82a293b61 100644 --- a/src/gromacs/listed_forces/listed_forces.h +++ b/src/gromacs/listed_forces/listed_forces.h @@ -94,6 +94,8 @@ struct t_lambda; struct t_mdatoms; struct t_nrnb; class t_state; +struct t_disresdata; +struct t_oriresdata; namespace gmx { @@ -117,6 +119,8 @@ using BondedFunction = real (*)(int nbonds, real* dvdlambda, const t_mdatoms* md, t_fcdata* fcd, + t_disresdata* disresdata, + t_oriresdata* oriresdata, int* ddgatindex); //! Getter for finding a callable CPU function to compute an \c ftype interaction. diff --git a/src/gromacs/listed_forces/orires.cpp b/src/gromacs/listed_forces/orires.cpp index 47e3eb3714..ff9a7b5cec 100644 --- a/src/gromacs/listed_forces/orires.cpp +++ b/src/gromacs/listed_forces/orires.cpp @@ -658,28 +658,28 @@ real orires(int nfa, real gmx_unused lambda, real gmx_unused* dvdlambda, const t_mdatoms gmx_unused* md, - t_fcdata* fcd, + t_fcdata gmx_unused* fcd, + t_disresdata gmx_unused* disresdata, + t_oriresdata* oriresdata, int gmx_unused* global_atom_index) { - int ex, power, ki = CENTRAL; - real r2, invr, invr2, fc, smooth_fc, dev, devins, pfac; - rvec r, Sr, fij; - real vtot; - const t_oriresdata* od; - gmx_bool bTAV; + int ex, power, ki = CENTRAL; + real r2, invr, invr2, fc, smooth_fc, dev, devins, pfac; + rvec r, Sr, fij; + real vtot; + gmx_bool bTAV; vtot = 0; - od = fcd->orires; - if (od->fc != 0) + if (oriresdata->fc != 0) { - bTAV = (od->edt != 0); + bTAV = (oriresdata->edt != 0); - smooth_fc = od->fc; + smooth_fc = oriresdata->fc; if (bTAV) { /* Smoothly switch on the restraining when time averaging is used */ - smooth_fc *= (1.0 - od->exp_min_t_tau); + smooth_fc *= (1.0 - oriresdata->exp_min_t_tau); } for (int fa = 0; fa < nfa; fa += 3) @@ -687,7 +687,7 @@ real orires(int nfa, const int type = forceatoms[fa]; const int ai = forceatoms[fa + 1]; const int aj = forceatoms[fa + 2]; - const int restraintIndex = type - od->typeMin; + const int restraintIndex = type - oriresdata->typeMin; if (pbc) { ki = pbc_dx_aiuc(pbc, x[ai], x[aj], r); @@ -702,7 +702,7 @@ real orires(int nfa, ex = ip[type].orires.ex; power = ip[type].orires.power; fc = smooth_fc * ip[type].orires.kfac; - dev = od->otav[restraintIndex] - ip[type].orires.obs; + dev = oriresdata->otav[restraintIndex] - ip[type].orires.obs; /* NOTE: * there is no real potential when time averaging is applied @@ -712,7 +712,7 @@ real orires(int nfa, if (bTAV) { /* Calculate the force as the sqrt of tav times instantaneous */ - devins = od->oins[restraintIndex] - ip[type].orires.obs; + devins = oriresdata->oins[restraintIndex] - ip[type].orires.obs; if (dev * devins <= 0) { dev = 0; @@ -732,7 +732,7 @@ real orires(int nfa, { pfac *= invr; } - mvmul(od->S[ex], r, Sr); + mvmul(oriresdata->S[ex], r, Sr); for (int i = 0; i < DIM; i++) { fij[i] = -pfac * dev * (4 * Sr[i] - 2 * (2 + power) * invr2 * iprod(Sr, r) * r[i]); diff --git a/src/gromacs/listed_forces/orires.h b/src/gromacs/listed_forces/orires.h index 0ad20a07af..7c269b0217 100644 --- a/src/gromacs/listed_forces/orires.h +++ b/src/gromacs/listed_forces/orires.h @@ -56,6 +56,8 @@ struct t_inputrec; struct t_pbc; struct t_commrec; struct t_oriresdata; +struct t_disresdata; +struct t_fcdata; class t_state; namespace gmx @@ -119,6 +121,8 @@ real orires(int nfa, real* dvdlambda, const t_mdatoms* md, t_fcdata* fcd, + t_disresdata* disresdata, + t_oriresdata* oriresdata, int* global_atom_index); //! Copies the new time averages that have been calculated in calc_orires_dev. diff --git a/src/gromacs/listed_forces/tests/bonded.cpp b/src/gromacs/listed_forces/tests/bonded.cpp index 554213953e..2564a77f3c 100644 --- a/src/gromacs/listed_forces/tests/bonded.cpp +++ b/src/gromacs/listed_forces/tests/bonded.cpp @@ -606,6 +606,8 @@ protected: &output.dvdlambda, &mdatoms, /* struct t_fcdata * */ nullptr, + nullptr, + nullptr, ddgatindex.data(), flavor); // Internal consistency test of both test input -- 2.22.0