Correct CUDA kernel energy flag
authorBerk Hess <hess@kth.se>
Mon, 2 Feb 2015 14:11:11 +0000 (15:11 +0100)
committerGerrit Code Review <gerrit@gerrit.gromacs.org>
Wed, 4 Feb 2015 09:57:37 +0000 (10:57 +0100)
The CUDA kernels calculated energies based on the GMX_FORCE_VIRIAL
flag. This did not cause errors, since (currently) GMX_FORCE_ENERGY
is always set when the virial flag is set. But using the latter flag
gives a small performance improvement when using pressure coupling.

Change-Id: If874e651058dc06c464f0fa810b17ba83146c9a3

src/gromacs/mdlib/nbnxn_cuda/nbnxn_cuda.cu

index d02e9034104480e733b0d4bdadd3827ba00eccb2..29af9eb52f937c73df36f5b037f7bee99a8cab72 100644 (file)
@@ -1,7 +1,7 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2012,2013,2014, by the GROMACS development team, led by
+ * Copyright (c) 2012,2013,2014,2015, 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.
@@ -295,7 +295,7 @@ void nbnxn_cuda_launch_kernel(nbnxn_cuda_ptr_t        cu_nb,
     cu_timers_t         *t       = cu_nb->timers;
     cudaStream_t         stream  = cu_nb->stream[iloc];
 
-    bool                 bCalcEner   = flags & GMX_FORCE_VIRIAL;
+    bool                 bCalcEner   = flags & GMX_FORCE_ENERGY;
     bool                 bCalcFshift = flags & GMX_FORCE_VIRIAL;
     bool                 bDoTime     = cu_nb->bDoTime;
 
@@ -422,7 +422,7 @@ void nbnxn_cuda_launch_cpyback(nbnxn_cuda_ptr_t        cu_nb,
     bool             bDoTime = cu_nb->bDoTime;
     cudaStream_t     stream  = cu_nb->stream[iloc];
 
-    bool             bCalcEner   = flags & GMX_FORCE_VIRIAL;
+    bool             bCalcEner   = flags & GMX_FORCE_ENERGY;
     bool             bCalcFshift = flags & GMX_FORCE_VIRIAL;
 
     /* don't launch copy-back if there was no work to do */
@@ -578,7 +578,7 @@ void nbnxn_cuda_wait_gpu(nbnxn_cuda_ptr_t cu_nb,
     wallclock_gpu_t *timings  = cu_nb->timings;
     nb_staging       nbst     = cu_nb->nbst;
 
-    bool             bCalcEner   = flags & GMX_FORCE_VIRIAL;
+    bool             bCalcEner   = flags & GMX_FORCE_ENERGY;
     bool             bCalcFshift = flags & GMX_FORCE_VIRIAL;
 
     /* turn energy calculation always on/off (for debugging/testing only) */