2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 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.
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.
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.
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.
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.
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.
36 /*! \libinternal \file
38 * \brief SIMD operations corresponding to Gromacs rvec datatypes.
40 * \author Erik Lindahl <erik.lindahl@scilifelab.se>
43 * \ingroup module_simd
46 #ifndef GMX_SIMD_VECTOR_OPERATIONS_H
47 #define GMX_SIMD_VECTOR_OPERATIONS_H
49 #include "gromacs/simd/simd.h"
54 /*! \addtogroup module_simd */
57 #if (defined GMX_SIMD_HAVE_FLOAT) || (defined DOXYGEN)
58 /*! \brief SIMD float inner product of multiple float vectors.
60 * For normal usage you should always call the real-precision \ref gmx_simd_iprod_r.
62 * \param ax X components of first vectors
63 * \param ay Y components of first vectors
64 * \param az Z components of first vectors
65 * \param bx X components of second vectors
66 * \param by Y components of second vectors
67 * \param bz Z components of second vectors
69 * \return Element i will be res[i] = ax[i]*bx[i]+ay[i]*by[i]+az[i]*bz[i].
71 * \note The SIMD part is that we calculate many scalar products in one call.
73 static gmx_inline gmx_simd_float_t gmx_simdcall
74 gmx_simd_iprod_f(gmx_simd_float_t ax, gmx_simd_float_t ay, gmx_simd_float_t az,
75 gmx_simd_float_t bx, gmx_simd_float_t by, gmx_simd_float_t bz)
79 ret = gmx_simd_mul_f(ax, bx);
80 ret = gmx_simd_fmadd_f(ay, by, ret);
81 ret = gmx_simd_fmadd_f(az, bz, ret);
86 /*! \brief SIMD float norm squared of multiple vectors.
88 * For normal usage you should always call the real-precision \ref gmx_simd_norm2_r.
90 * \param ax X components of vectors
91 * \param ay Y components of vectors
92 * \param az Z components of vectors
94 * \return Element i will be res[i] = ax[i]*ax[i]+ay[i]*ay[i]+az[i]*az[i].
96 * \note This corresponds to the scalar product of the vector with itself, but
97 * the compiler might be able to optimize it better with identical vectors.
99 static gmx_inline gmx_simd_float_t gmx_simdcall
100 gmx_simd_norm2_f(gmx_simd_float_t ax, gmx_simd_float_t ay, gmx_simd_float_t az)
102 gmx_simd_float_t ret;
104 ret = gmx_simd_mul_f(ax, ax);
105 ret = gmx_simd_fmadd_f(ay, ay, ret);
106 ret = gmx_simd_fmadd_f(az, az, ret);
111 /*! \brief Calculating r^2 is the same as evaluating the norm of dx*dx.
113 * For details, see \ref gmx_simd_norm2_f.
115 #define gmx_simd_calc_rsq_f gmx_simd_norm2_f
117 /*! \brief SIMD float cross-product of multiple vectors.
119 * For normal usage you should always call the real-precision \ref gmx_simd_cprod_r.
121 * \param ax X components of first vectors
122 * \param ay Y components of first vectors
123 * \param az Z components of first vectors
124 * \param bx X components of second vectors
125 * \param by Y components of second vectors
126 * \param bz Z components of second vectors
127 * \param[out] cx X components of cross product vectors
128 * \param[out] cy Y components of cross product vectors
129 * \param[out] cz Z components of cross product vectors
133 * This calculates C = A x B, where the cross denotes the cross product.
134 * The arguments x/y/z denotes the different components, and each element
135 * corresponds to a separate vector.
137 static gmx_inline void gmx_simdcall
138 gmx_simd_cprod_f(gmx_simd_float_t ax, gmx_simd_float_t ay, gmx_simd_float_t az,
139 gmx_simd_float_t bx, gmx_simd_float_t by, gmx_simd_float_t bz,
140 gmx_simd_float_t *cx, gmx_simd_float_t *cy, gmx_simd_float_t *cz)
142 *cx = gmx_simd_mul_f(ay, bz);
143 *cx = gmx_simd_fnmadd_f(az, by, *cx);
145 *cy = gmx_simd_mul_f(az, bx);
146 *cy = gmx_simd_fnmadd_f(ax, bz, *cy);
148 *cz = gmx_simd_mul_f(ax, by);
149 *cz = gmx_simd_fnmadd_f(ay, bx, *cz);
151 #endif /* GMX_SIMD_HAVE_FLOAT */
153 #if (defined GMX_SIMD_HAVE_DOUBLE) || (defined DOXYGEN)
154 /*! \brief SIMD double inner product of multiple double vectors.
156 * \copydetails gmx_simd_iprod_f
158 static gmx_inline gmx_simd_double_t gmx_simdcall
159 gmx_simd_iprod_d(gmx_simd_double_t ax, gmx_simd_double_t ay, gmx_simd_double_t az,
160 gmx_simd_double_t bx, gmx_simd_double_t by, gmx_simd_double_t bz)
162 gmx_simd_double_t ret;
164 ret = gmx_simd_mul_d(ax, bx);
165 ret = gmx_simd_fmadd_d(ay, by, ret);
166 ret = gmx_simd_fmadd_d(az, bz, ret);
171 /*! \brief SIMD double norm squared of multiple vectors.
173 * \copydetails gmx_simd_norm2_f
175 static gmx_inline gmx_simd_double_t gmx_simdcall
176 gmx_simd_norm2_d(gmx_simd_double_t ax, gmx_simd_double_t ay, gmx_simd_double_t az)
178 gmx_simd_double_t ret;
180 ret = gmx_simd_mul_d(ax, ax);
181 ret = gmx_simd_fmadd_d(ay, ay, ret);
182 ret = gmx_simd_fmadd_d(az, az, ret);
187 /*! \brief Calculating r^2 is the same as evaluating the norm of dx*dx.
189 * For details, see \ref gmx_simd_norm2_d.
191 #define gmx_simd_calc_rsq_d gmx_simd_norm2_d
193 /*! \brief SIMD double cross-product of multiple vectors.
195 * \copydetails gmx_simd_cprod_f
197 static gmx_inline void gmx_simdcall
198 gmx_simd_cprod_d(gmx_simd_double_t ax, gmx_simd_double_t ay, gmx_simd_double_t az,
199 gmx_simd_double_t bx, gmx_simd_double_t by, gmx_simd_double_t bz,
200 gmx_simd_double_t *cx, gmx_simd_double_t *cy, gmx_simd_double_t *cz)
202 *cx = gmx_simd_mul_d(ay, bz);
203 *cx = gmx_simd_fnmadd_d(az, by, *cx);
205 *cy = gmx_simd_mul_d(az, bx);
206 *cy = gmx_simd_fnmadd_d(ax, bz, *cy);
208 *cz = gmx_simd_mul_d(ax, by);
209 *cz = gmx_simd_fnmadd_d(ay, bx, *cz);
211 #endif /* GMX_SIMD_HAVE_DOUBLE */
214 #if (defined GMX_SIMD4_HAVE_FLOAT) || (defined DOXYGEN)
215 /*! \brief SIMD4 float inner product of four float vectors.
217 * \copydetails gmx_simd_norm2_f
219 static gmx_inline gmx_simd4_float_t gmx_simdcall
220 gmx_simd4_norm2_f(gmx_simd4_float_t ax, gmx_simd4_float_t ay, gmx_simd4_float_t az)
222 gmx_simd4_float_t ret;
224 ret = gmx_simd4_mul_f(ax, ax);
225 ret = gmx_simd4_fmadd_f(ay, ay, ret);
226 ret = gmx_simd4_fmadd_f(az, az, ret);
231 /*! \brief Calculating r^2 is the same as evaluating the norm of dx*dx.
233 * For details, see \ref gmx_simd4_norm2_f
235 #define gmx_simd4_calc_rsq_f gmx_simd4_norm2_f
237 #endif /* GMX_SIMD4_HAVE_FLOAT */
239 #if (defined GMX_SIMD4_HAVE_DOUBLE) || (defined DOXYGEN)
240 /*! \brief SIMD4 double norm squared of multiple vectors.
242 * \copydetails gmx_simd_norm2_f
244 static gmx_inline gmx_simd4_double_t gmx_simdcall
245 gmx_simd4_norm2_d(gmx_simd4_double_t ax, gmx_simd4_double_t ay, gmx_simd4_double_t az)
247 gmx_simd4_double_t ret;
249 ret = gmx_simd4_mul_d(ax, ax);
250 ret = gmx_simd4_fmadd_d(ay, ay, ret);
251 ret = gmx_simd4_fmadd_d(az, az, ret);
256 /*! \brief Calculating r^2 is the same as evaluating the norm of dx*dx.
258 * For details, see \ref gmx_simd4_norm2_d.
260 #define gmx_simd4_calc_rsq_d gmx_simd4_norm2_d
262 #endif /* GMX_SIMD4_HAVE_DOUBLE */
266 /* Documented for the single branch below */
267 # define gmx_simd_iprod_r gmx_simd_iprod_d
268 # define gmx_simd_norm2_r gmx_simd_norm2_d
269 # define gmx_simd_calc_rsq_r gmx_simd_calc_rsq_d
270 # define gmx_simd_cprod_r gmx_simd_cprod_d
271 # define gmx_simd4_norm2_r gmx_simd4_norm2_d
272 # define gmx_simd4_calc_rsq_r gmx_simd4_calc_rsq_d
273 #else /* GMX_DOUBLE */
275 /*! \brief SIMD real inner product of multiple real vectors.
277 * This will call \ref gmx_simd_iprod_d if GMX_DOUBLE is defined, otherwise
278 * \ref gmx_simd_iprod_f.
280 * \copydetails gmx_simd_iprod_f
282 # define gmx_simd_iprod_r gmx_simd_iprod_f
284 /*! \brief SIMD real norm squared of multiple real vectors.
286 * This will call \ref gmx_simd_norm2_d if GMX_DOUBLE is defined, otherwise
287 * \ref gmx_simd_norm2_f.
289 * \copydetails gmx_simd_norm2_f
291 # define gmx_simd_norm2_r gmx_simd_norm2_f
293 /*! \brief Calculating r^2 is the same as evaluating the norm of dx*dx.
295 * This will call \ref gmx_simd_calc_rsq_d if GMX_DOUBLE is defined, otherwise
296 * \ref gmx_simd_calc_rsq_f.
298 * \copydetails gmx_simd_calc_rsq_f
300 # define gmx_simd_calc_rsq_r gmx_simd_calc_rsq_f
302 /*! \brief SIMD real cross-product of multiple real vectors.
304 * This will call \ref gmx_simd_cprod_d if GMX_DOUBLE is defined, otherwise
305 * \ref gmx_simd_cprod_f.
307 * \copydetails gmx_simd_cprod_f
309 # define gmx_simd_cprod_r gmx_simd_cprod_f
311 /*! \brief SIMD4 real norm squared of multiple vectors.
313 * This will call \ref gmx_simd4_norm2_d if GMX_DOUBLE is defined, otherwise
314 * \ref gmx_simd4_norm2_f.
316 * \copydetails gmx_simd4_norm2_f
318 # define gmx_simd4_norm2_r gmx_simd4_norm2_f
320 /*! \brief Calculating r^2 is the same as evaluating the norm of dx*dx.
322 * This will call \ref gmx_simd4_calc_rsq_d if GMX_DOUBLE is defined, otherwise
323 * \ref gmx_simd4_calc_rsq_f.
325 * \copydetails gmx_simd4_calc_rsq_f
327 # define gmx_simd4_calc_rsq_r gmx_simd4_calc_rsq_f
329 #endif /* GMX_DOUBLE */
334 #endif /* GMX_SIMD_VECTOR_OPERATIONS_H */