1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
4 * This source code is part of
8 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
9 * Copyright (c) 2001-2009, The GROMACS Development Team
11 * Gromacs is a library for molecular simulation and trajectory analysis,
12 * written by Erik Lindahl, David van der Spoel, Berk Hess, and others - for
13 * a full list of developers and information, check out http://www.gromacs.org
15 * This program is free software; you can redistribute it and/or modify it under
16 * the terms of the GNU Lesser General Public License as published by the Free
17 * Software Foundation; either version 2 of the License, or (at your option) any
19 * As a special exception, you may use this file as part of a free software
20 * library without restriction. Specifically, if other files instantiate
21 * templates or use macros or inline functions from this file, or you compile
22 * this file and link it with other files to produce an executable, this
23 * file does not by itself cause the resulting executable to be covered by
24 * the GNU Lesser General Public License.
26 * In plain-speak: do not worry about classes/macros/templates either - only
27 * changes to the library have to be LGPL, not an application linking with it.
29 * To help fund GROMACS development, we humbly ask that you cite
30 * the papers people have written on it - you can find them on the website!
36 #include "nbnxn_kernel_common.h"
39 clear_f_all(const nbnxn_atomdata_t *nbat, real *f)
43 for (i = 0; i < nbat->natoms*nbat->fstride; i++)
50 clear_f_flagged(const nbnxn_atomdata_t *nbat, int output_index, real *f)
52 const nbnxn_buffer_flags_t *flags;
56 flags = &nbat->buffer_flags;
58 our_flag = (1U << output_index);
60 for (b = 0; b < flags->nflag; b++)
62 if (flags->flag[b] & our_flag)
64 a0 = b*NBNXN_BUFFERFLAG_SIZE;
65 a1 = a0 + NBNXN_BUFFERFLAG_SIZE;
66 for (i = a0*nbat->fstride; i < a1*nbat->fstride; i++)
75 clear_f(const nbnxn_atomdata_t *nbat, int output_index, real *f)
77 if (nbat->bUseBufferFlags)
79 clear_f_flagged(nbat, output_index, f);
88 clear_fshift(real *fshift)
92 for (i = 0; i < SHIFTS*DIM; i++)
99 reduce_energies_over_lists(const nbnxn_atomdata_t *nbat,
107 for (nb = 0; nb < nlist; nb++)
109 for (i = 0; i < nbat->nenergrp; i++)
111 /* Reduce the diagonal terms */
112 ind = i*nbat->nenergrp + i;
113 Vvdw[ind] += nbat->out[nb].Vvdw[ind];
114 Vc[ind] += nbat->out[nb].Vc[ind];
116 /* Reduce the off-diagonal terms */
117 for (j = i+1; j < nbat->nenergrp; j++)
119 /* The output should contain only one off-diagonal part */
120 ind = i*nbat->nenergrp + j;
121 indr = j*nbat->nenergrp + i;
122 Vvdw[ind] += nbat->out[nb].Vvdw[ind] + nbat->out[nb].Vvdw[indr];
123 Vc[ind] += nbat->out[nb].Vc[ind] + nbat->out[nb].Vc[indr];