Bug Summary

File:gromacs/mdlib/nbnxn_kernels/nbnxn_kernel_ref_inner.h
Location:line 203, column 21
Description:Value stored to 'VLJ' is never read

Annotated Source Code

1/*
2 * This file is part of the GROMACS molecular simulation package.
3 *
4 * Copyright (c) 2012,2013,2014, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
8 *
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
13 *
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
18 *
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
23 *
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
31 *
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
34 */
35
36/* When calculating RF or Ewald interactions we calculate the electrostatic
37 * forces and energies on excluded atom pairs here in the non-bonded loops.
38 */
39#if defined CHECK_EXCLS && (defined CALC_COULOMB || defined LJ_EWALD)
40#define EXCL_FORCES
41#endif
42
43{
44 int cj;
45#ifdef ENERGY_GROUPS
46 int egp_cj;
47#endif
48 int i;
49
50 cj = l_cj[cjind].cj;
51
52#ifdef ENERGY_GROUPS
53 egp_cj = nbat->energrp[cj];
54#endif
55 for (i = 0; i < UNROLLI; i++)
56 {
57 int ai;
58 int type_i_off;
59 int j;
60
61 ai = ci*UNROLLI + i;
62
63 type_i_off = type[ai]*ntype2;
64
65 for (j = 0; j < UNROLLJ; j++)
66 {
67 int aj;
68 real dx, dy, dz;
69 real rsq, rinv;
70 real rinvsq, rinvsix;
71 real c6, c12;
72 real FrLJ6 = 0, FrLJ12 = 0, frLJ = 0, VLJ = 0;
73#if defined LJ_FORCE_SWITCH || defined LJ_POT_SWITCH
74 real r, rsw;
75#endif
76
77#ifdef CALC_COULOMB
78 real qq;
79 real fcoul;
80#ifdef CALC_COUL_TAB
81 real rs, frac;
82 int ri;
83 real fexcl;
84#endif
85#ifdef CALC_ENERGIES
86 real vcoul;
87#endif
88#endif
89 real fscal;
90 real fx, fy, fz;
91
92 /* A multiply mask used to zero an interaction
93 * when either the distance cutoff is exceeded, or
94 * (if appropriate) the i and j indices are
95 * unsuitable for this kind of inner loop. */
96 real skipmask;
97
98#ifdef CHECK_EXCLS
99 /* A multiply mask used to zero an interaction
100 * when that interaction should be excluded
101 * (e.g. because of bonding). */
102 int interact;
103
104 interact = ((l_cj[cjind].excl>>(i*UNROLLI + j)) & 1);
105#ifndef EXCL_FORCES
106 skipmask = interact;
107#else
108 skipmask = !(cj == ci_sh && j <= i);
109#endif
110#else
111#define interact 1.0
112 skipmask = 1.0;
113#endif
114
115 aj = cj*UNROLLJ + j;
116
117 dx = xi[i*XI_STRIDE+XX0] - x[aj*X_STRIDE+XX0];
118 dy = xi[i*XI_STRIDE+YY1] - x[aj*X_STRIDE+YY1];
119 dz = xi[i*XI_STRIDE+ZZ2] - x[aj*X_STRIDE+ZZ2];
120
121 rsq = dx*dx + dy*dy + dz*dz;
122
123 /* Prepare to enforce the cut-off. */
124 skipmask = (rsq >= rcut2) ? 0 : skipmask;
125 /* 9 flops for r^2 + cut-off check */
126
127#ifdef CHECK_EXCLS
128 /* Excluded atoms are allowed to be on top of each other.
129 * To avoid overflow of rinv, rinvsq and rinvsix
130 * we add a small number to rsq for excluded pairs only.
131 */
132 rsq += (1 - interact)*NBNXN_AVOID_SING_R2_INC1.0e-12f;
133#endif
134
135#ifdef COUNT_PAIRS
136 npair++;
137#endif
138
139 rinv = gmx_invsqrt(rsq)gmx_software_invsqrt(rsq);
140 /* 5 flops for invsqrt */
141
142 /* Partially enforce the cut-off (and perhaps
143 * exclusions) to avoid possible overflow of
144 * rinvsix when computing LJ, and/or overflowing
145 * the Coulomb table during lookup. */
146 rinv = rinv * skipmask;
147
148 rinvsq = rinv*rinv;
149
150#ifdef HALF_LJ
151 if (i < UNROLLI/2)
152#endif
153 {
154 c6 = nbfp[type_i_off+type[aj]*2 ];
155 c12 = nbfp[type_i_off+type[aj]*2+1];
156
157#if defined LJ_CUT || defined LJ_FORCE_SWITCH || defined LJ_POT_SWITCH
158 rinvsix = interact*rinvsq*rinvsq*rinvsq;
159 FrLJ6 = c6*rinvsix;
160 FrLJ12 = c12*rinvsix*rinvsix;
161 frLJ = FrLJ12 - FrLJ6;
162 /* 7 flops for r^-2 + LJ force */
163#if defined CALC_ENERGIES || defined LJ_POT_SWITCH
164 VLJ = (FrLJ12 + c12*ic->repulsion_shift.cpot)/12 -
165 (FrLJ6 + c6*ic->dispersion_shift.cpot)/6;
166 /* 7 flops for LJ energy */
167#endif
168#endif
169
170#if defined LJ_FORCE_SWITCH || defined LJ_POT_SWITCH
171 /* Force or potential switching from ic->rvdw_switch */
172 r = rsq*rinv;
173 rsw = r - ic->rvdw_switch;
174 rsw = (rsw >= 0.0 ? rsw : 0.0);
175#endif
176#ifdef LJ_FORCE_SWITCH
177 frLJ +=
178 -c6*(ic->dispersion_shift.c2 + ic->dispersion_shift.c3*rsw)*rsw*rsw*r
179 + c12*(ic->repulsion_shift.c2 + ic->repulsion_shift.c3*rsw)*rsw*rsw*r;
180#if defined CALC_ENERGIES
181 VLJ +=
182 -c6*(-ic->dispersion_shift.c2/3 - ic->dispersion_shift.c3/4*rsw)*rsw*rsw*rsw
183 + c12*(-ic->repulsion_shift.c2/3 - ic->repulsion_shift.c3/4*rsw)*rsw*rsw*rsw;
184#endif
185#endif
186
187#if defined CALC_ENERGIES || defined LJ_POT_SWITCH
188 /* Masking should be done after force switching,
189 * but before potential switching.
190 */
191 /* Need to zero the interaction if there should be exclusion. */
192 VLJ = VLJ * interact;
193#endif
194
195#ifdef LJ_POT_SWITCH
196 {
197 real sw, dsw;
198
199 sw = 1.0 + (swV3 + (swV4+ swV5*rsw)*rsw)*rsw*rsw*rsw;
200 dsw = (swF2 + (swF3 + swF4*rsw)*rsw)*rsw*rsw;
201
202 frLJ = frLJ*sw - r*VLJ*dsw;
203 VLJ *= sw;
Value stored to 'VLJ' is never read
204 }
205#endif
206
207#ifdef LJ_EWALD
208 {
209 real c6grid, rinvsix_nm, cr2, expmcr2, poly, sh_mask;
210
211#ifdef LJ_EWALD_COMB_GEOM
212 c6grid = ljc[type[ai]*2]*ljc[type[aj]*2];
213#elif defined LJ_EWALD_COMB_LB
214 {
215 real sigma, sigma2, epsilon;
216
217 /* These sigma and epsilon are scaled to give 6*C6 */
218 sigma = ljc[type[ai]*2] + ljc[type[aj]*2];
219 epsilon = ljc[type[ai]*2+1]*ljc[type[aj]*2+1];
220
221 sigma2 = sigma*sigma;
222 c6grid = epsilon*sigma2*sigma2*sigma2;
223 }
224#else
225#error "No LJ Ewald combination rule defined"
226#endif
227
228#ifdef CHECK_EXCLS
229 /* Recalculate rinvsix without exclusion mask */
230 rinvsix_nm = rinvsq*rinvsq*rinvsq;
231#else
232 rinvsix_nm = rinvsix;
233#endif
234 cr2 = lje_coeff2*rsq;
235#ifdef GMX_DOUBLE
236 expmcr2 = exp(-cr2);
237#else
238 expmcr2 = expf(-cr2);
239#endif
240 poly = 1 + cr2 + 0.5*cr2*cr2;
241
242 /* Subtract the grid force from the total LJ force */
243 frLJ += c6grid*(rinvsix_nm - expmcr2*(rinvsix_nm*poly + lje_coeff6_6));
244#ifdef CALC_ENERGIES
245 /* Shift should only be applied to real LJ pairs */
246 sh_mask = lje_vc*interact;
247
248 VLJ += c6grid/6*(rinvsix_nm*(1 - expmcr2*poly) + sh_mask);
249#endif
250 }
251#endif /* LJ_EWALD */
252
253#ifdef VDW_CUTOFF_CHECK
254 /* Mask for VdW cut-off shorter than Coulomb cut-off */
255 {
256 real skipmask_rvdw;
257
258 skipmask_rvdw = (rsq < rvdw2);
259 frLJ *= skipmask_rvdw;
260#ifdef CALC_ENERGIES
261 VLJ *= skipmask_rvdw;
262#endif
263 }
264#else
265#if defined CALC_ENERGIES
266 /* Need to zero the interaction if r >= rcut */
267 VLJ = VLJ * skipmask;
268 /* 1 more flop for LJ energy */
269#endif
270#endif /* VDW_CUTOFF_CHECK */
271
272
273#ifdef CALC_ENERGIES
274#ifdef ENERGY_GROUPS
275 Vvdw[egp_sh_i[i]+((egp_cj>>(nbat->neg_2log*j)) & egp_mask)] += VLJ;
276#else
277 Vvdw_ci += VLJ;
278 /* 1 flop for LJ energy addition */
279#endif
280#endif
281 }
282
283#ifdef CALC_COULOMB
284 /* Enforce the cut-off and perhaps exclusions. In
285 * those cases, rinv is zero because of skipmask,
286 * but fcoul and vcoul will later be non-zero (in
287 * both RF and table cases) because of the
288 * contributions that do not depend on rinv. These
289 * contributions cannot be allowed to accumulate
290 * to the force and potential, and the easiest way
291 * to do this is to zero the charges in
292 * advance. */
293 qq = skipmask * qi[i] * q[aj];
294
295#ifdef CALC_COUL_RF
296 fcoul = qq*(interact*rinv*rinvsq - k_rf2);
297 /* 4 flops for RF force */
298#ifdef CALC_ENERGIES
299 vcoul = qq*(interact*rinv + k_rf*rsq - c_rf);
300 /* 4 flops for RF energy */
301#endif
302#endif
303
304#ifdef CALC_COUL_TAB
305 rs = rsq*rinv*ic->tabq_scale;
306 ri = (int)rs;
307 frac = rs - ri;
308#ifndef GMX_DOUBLE
309 /* fexcl = F_i + frac * (F_(i+1)-F_i) */
310 fexcl = tab_coul_FDV0[ri*4] + frac*tab_coul_FDV0[ri*4+1];
311#else
312 /* fexcl = (1-frac) * F_i + frac * F_(i+1) */
313 fexcl = (1 - frac)*tab_coul_F[ri] + frac*tab_coul_F[ri+1];
314#endif
315 fcoul = interact*rinvsq - fexcl;
316 /* 7 flops for float 1/r-table force */
317#ifdef CALC_ENERGIES
318#ifndef GMX_DOUBLE
319 vcoul = qq*(interact*(rinv - ic->sh_ewald)
320 -(tab_coul_FDV0[ri*4+2]
321 -halfsp*frac*(tab_coul_FDV0[ri*4] + fexcl)));
322 /* 7 flops for float 1/r-table energy (8 with excls) */
323#else
324 vcoul = qq*(interact*(rinv - ic->sh_ewald)
325 -(tab_coul_V[ri]
326 -halfsp*frac*(tab_coul_F[ri] + fexcl)));
327#endif
328#endif
329 fcoul *= qq*rinv;
330#endif
331
332#ifdef CALC_ENERGIES
333#ifdef ENERGY_GROUPS
334 Vc[egp_sh_i[i]+((egp_cj>>(nbat->neg_2log*j)) & egp_mask)] += vcoul;
335#else
336 Vc_ci += vcoul;
337 /* 1 flop for Coulomb energy addition */
338#endif
339#endif
340#endif
341
342#ifdef CALC_COULOMB
343#ifdef HALF_LJ
344 if (i < UNROLLI/2)
345#endif
346 {
347 fscal = frLJ*rinvsq + fcoul;
348 /* 2 flops for scalar LJ+Coulomb force */
349 }
350#ifdef HALF_LJ
351 else
352 {
353 fscal = fcoul;
354 }
355#endif
356#else
357 fscal = frLJ*rinvsq;
358#endif
359 fx = fscal*dx;
360 fy = fscal*dy;
361 fz = fscal*dz;
362
363 /* Increment i-atom force */
364 fi[i*FI_STRIDE+XX0] += fx;
365 fi[i*FI_STRIDE+YY1] += fy;
366 fi[i*FI_STRIDE+ZZ2] += fz;
367 /* Decrement j-atom force */
368 f[aj*F_STRIDE+XX0] -= fx;
369 f[aj*F_STRIDE+YY1] -= fy;
370 f[aj*F_STRIDE+ZZ2] -= fz;
371 /* 9 flops for force addition */
372 }
373 }
374}
375
376#undef interact
377#undef EXCL_FORCES