Template sc-gapsys functions on computeForces.
authorSebastian Kehl <sebastian.kehl@mpcdf.mpg.de>
Fri, 1 Oct 2021 12:22:21 +0000 (14:22 +0200)
committerMagnus Lundborg <magnus.lundborg@scilifelab.se>
Tue, 5 Oct 2021 11:38:46 +0000 (11:38 +0000)
src/gromacs/gmxlib/nonbonded/nb_free_energy.cpp
src/gromacs/gmxlib/nonbonded/nb_softcore.h

index 8e0ce26e43b4fc4160132d25cd40502528ca14b0..d441bd199a0762acdd889996ad0f4cc3d0c71964 100644 (file)
@@ -808,18 +808,18 @@ static void nb_free_energy_kernel(const t_nblist&
 
                                 if constexpr (softcoreType == KernelSoftcoreType::Gapsys)
                                 {
-                                    ewaldQuadraticPotential(qq[i],
-                                                            facel,
-                                                            rC,
-                                                            rCutoffCoul,
-                                                            LFC[i],
-                                                            DLF[i],
-                                                            gapsysScaleLinpointCoulEff,
-                                                            sh_ewald,
-                                                            &fScalC[i],
-                                                            &vCoul[i],
-                                                            &dvdlCoul,
-                                                            computeElecInteraction);
+                                    ewaldQuadraticPotential<computeForces>(qq[i],
+                                                                           facel,
+                                                                           rC,
+                                                                           rCutoffCoul,
+                                                                           LFC[i],
+                                                                           DLF[i],
+                                                                           gapsysScaleLinpointCoulEff,
+                                                                           sh_ewald,
+                                                                           &fScalC[i],
+                                                                           &vCoul[i],
+                                                                           &dvdlCoul,
+                                                                           computeElecInteraction);
                                 }
                             }
                             else
@@ -832,19 +832,20 @@ static void nb_free_energy_kernel(const t_nblist&
 
                                 if constexpr (softcoreType == KernelSoftcoreType::Gapsys)
                                 {
-                                    reactionFieldQuadraticPotential(qq[i],
-                                                                    facel,
-                                                                    rC,
-                                                                    rCutoffCoul,
-                                                                    LFC[i],
-                                                                    DLF[i],
-                                                                    gapsysScaleLinpointCoulEff,
-                                                                    krf,
-                                                                    crf,
-                                                                    &fScalC[i],
-                                                                    &vCoul[i],
-                                                                    &dvdlCoul,
-                                                                    computeElecInteraction);
+                                    reactionFieldQuadraticPotential<computeForces>(
+                                            qq[i],
+                                            facel,
+                                            rC,
+                                            rCutoffCoul,
+                                            LFC[i],
+                                            DLF[i],
+                                            gapsysScaleLinpointCoulEff,
+                                            krf,
+                                            crf,
+                                            &fScalC[i],
+                                            &vCoul[i],
+                                            &dvdlCoul,
+                                            computeElecInteraction);
                                 }
                             }
 
@@ -899,20 +900,20 @@ static void nb_free_energy_kernel(const t_nblist&
 
                             if constexpr (softcoreType == KernelSoftcoreType::Gapsys)
                             {
-                                lennardJonesQuadraticPotential(c6[i],
-                                                               c12[i],
-                                                               r,
-                                                               rSq,
-                                                               LFV[i],
-                                                               DLF[i],
-                                                               gapsysSigma6VdWEff[i],
-                                                               gapsysScaleLinpointVdWEff,
-                                                               repulsionShift,
-                                                               dispersionShift,
-                                                               &fScalV[i],
-                                                               &vVdw[i],
-                                                               &dvdlVdw,
-                                                               computeVdwInteraction);
+                                lennardJonesQuadraticPotential<computeForces>(c6[i],
+                                                                              c12[i],
+                                                                              r,
+                                                                              rSq,
+                                                                              LFV[i],
+                                                                              DLF[i],
+                                                                              gapsysSigma6VdWEff[i],
+                                                                              gapsysScaleLinpointVdWEff,
+                                                                              repulsionShift,
+                                                                              dispersionShift,
+                                                                              &fScalV[i],
+                                                                              &vVdw[i],
+                                                                              &dvdlVdw,
+                                                                              computeVdwInteraction);
                             }
 
                             if constexpr (vdwInteractionTypeIsEwald)
index a2615ede4a7aedd78199d58b71bf94aa4de1bbda..dabbaa28534a0026ab965c312dbc03690d7a31b7 100644 (file)
@@ -68,7 +68,7 @@ static inline void quadraticApproximationCoulomb(const RealType qq,
 }
 
 /* reaction-field linearized electrostatics */
-template<class RealType, class BoolType>
+template<bool computeForces, class RealType, class BoolType>
 static inline void reactionFieldQuadraticPotential(const RealType qq,
                                                    const real     facel,
                                                    const RealType r,
@@ -118,7 +118,10 @@ static inline void reactionFieldQuadraticPotential(const RealType qq,
             potentialQuad = potentialQuad + qq * (krf * r * r - potentialShift);
 
             // update
-            *force     = gmx::blend(*force, forceQuad, computeValues);
+            if constexpr (computeForces)
+            {
+                *force     = gmx::blend(*force, forceQuad, computeValues);
+            }
             *potential = gmx::blend(*potential, potentialQuad, computeValues);
             *dvdl      = *dvdl + gmx::selectByMask(dvdlQuad, computeValues && withinCutoff);
         }
@@ -126,7 +129,7 @@ static inline void reactionFieldQuadraticPotential(const RealType qq,
 }
 
 /* ewald linearized electrostatics */
-template<class RealType, class BoolType>
+template<bool computeForces, class RealType, class BoolType>
 static inline void ewaldQuadraticPotential(const RealType qq,
                                            const real     facel,
                                            const RealType r,
@@ -174,7 +177,10 @@ static inline void ewaldQuadraticPotential(const RealType qq,
             potentialQuad = potentialQuad - qq * potentialShift;
 
             // update
-            *force     = gmx::blend(*force, forceQuad, computeValues);
+            if constexpr (computeForces)
+            {
+                *force     = gmx::blend(*force, forceQuad, computeValues);
+            }
             *potential = gmx::blend(*potential, potentialQuad, computeValues);
             *dvdl      = *dvdl + gmx::selectByMask(dvdlQuad, computeValues && withinCutoff);
         }
@@ -182,7 +188,7 @@ static inline void ewaldQuadraticPotential(const RealType qq,
 }
 
 /* cutoff LJ with quadratic appximation of lj-potential */
-template<class RealType, class BoolType>
+template<bool computeForces, class RealType, class BoolType>
 static inline void lennardJonesQuadraticPotential(const RealType c6,
                                                   const RealType c12,
                                                   const RealType r,
@@ -252,7 +258,10 @@ static inline void lennardJonesQuadraticPotential(const RealType c6,
             potentialQuad = potentialQuad
                             + gmx::selectByMask(((c12s * repulsionShift) - (c6s * dispersionShift)),
                                                 computeValues);
-            *force     = gmx::blend(*force, forceQuad, computeValues);
+            if constexpr (computeForces)
+            {
+                *force     = gmx::blend(*force, forceQuad, computeValues);
+            }
             *potential = gmx::blend(*potential, potentialQuad, computeValues);
             *dvdl      = *dvdl + gmx::selectByMask(dvdlQuad, computeValues);
         }