2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
37 /*! \libinternal \file
39 * \brief This file contains declarations necessary for low-level
40 * functions for computing energies and forces for bonded interactions.
42 * \author Mark Abraham <mark.j.abraham@gmail.com>
44 * \ingroup module_listed_forces
47 #ifndef GMX_LISTED_FORCES_BONDED_H
48 #define GMX_LISTED_FORCES_BONDED_H
50 #include "gromacs/math/vectypes.h"
51 #include "gromacs/topology/ifunc.h"
52 #include "gromacs/utility/basedefinitions.h"
61 /*! \brief Calculate bond-angle. No PBC is taken into account (use mol-shift) */
62 real bond_angle(const rvec xi, const rvec xj, const rvec xk,
63 const struct t_pbc *pbc,
64 rvec r_ij, rvec r_kj, real *costh,
65 int *t1, int *t2); /* out */
67 /*! \brief Calculate dihedral-angle. No PBC is taken into account (use mol-shift) */
68 real dih_angle(const rvec xi, const rvec xj, const rvec xk, const rvec xl,
69 const struct t_pbc *pbc,
70 rvec r_ij, rvec r_kj, rvec r_kl, rvec m, rvec n, /* out */
71 int *t1, int *t2, int *t3);
73 /*! \brief Do an update of the forces for dihedral potentials */
74 void do_dih_fup(int i, int j, int k, int l, real ddphi,
75 rvec r_ij, rvec r_kj, rvec r_kl,
76 rvec m, rvec n, rvec4 f[], rvec fshift[],
77 const struct t_pbc *pbc, const struct t_graph *g,
78 const rvec *x, int t1, int t2, int t3);
80 /*! \brief Make a dihedral fall in the range (-pi,pi) */
81 void make_dp_periodic(real *dp);
83 /*! \brief Compute CMAP dihedral energies and forces */
86 const t_iatom forceatoms[], const t_iparams forceparams[],
87 const gmx_cmap_t *cmap_grid,
88 const rvec x[], rvec4 f[], rvec fshift[],
89 const struct t_pbc *pbc, const struct t_graph *g,
90 real gmx_unused lambda, real gmx_unused *dvdlambda,
91 const t_mdatoms gmx_unused *md, t_fcdata gmx_unused *fcd,
92 int gmx_unused *global_atom_index);
95 /*************************************************************************
97 * Bonded force functions
99 *************************************************************************/
100 real bonds(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
101 const rvec x[], rvec4 f[], rvec fshift[],
102 const t_pbc *pbc, const t_graph *g,
103 real lambda, real *dvdlambda,
104 const t_mdatoms *md, t_fcdata *fcd,
105 int *global_atom_index);
107 real g96bonds(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
108 const rvec x[], rvec4 f[], rvec fshift[],
109 const t_pbc *pbc, const t_graph *g,
110 real lambda, real *dvdlambda,
111 const t_mdatoms *md, t_fcdata *fcd,
112 int *global_atom_index);
114 real morse_bonds(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
115 const rvec x[], rvec4 f[], rvec fshift[],
116 const t_pbc *pbc, const t_graph *g,
117 real lambda, real *dvdlambda,
118 const t_mdatoms *md, t_fcdata *fcd,
119 int *global_atom_index);
121 real cubic_bonds(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
122 const rvec x[], rvec4 f[], rvec fshift[],
123 const t_pbc *pbc, const t_graph *g,
124 real lambda, real *dvdlambda,
125 const t_mdatoms *md, t_fcdata *fcd,
126 int *global_atom_index);
128 real FENE_bonds(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
129 const rvec x[], rvec4 f[], rvec fshift[],
130 const t_pbc *pbc, const t_graph *g,
131 real lambda, real *dvdlambda,
132 const t_mdatoms *md, t_fcdata *fcd,
133 int *global_atom_index);
135 real restraint_bonds(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
136 const rvec x[], rvec4 f[], rvec fshift[],
137 const t_pbc *pbc, const t_graph *g,
138 real lambda, real *dvdlambda,
139 const t_mdatoms *md, t_fcdata *fcd,
140 int *global_atom_index);
142 real angles(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
143 const rvec x[], rvec4 f[], rvec fshift[],
144 const t_pbc *pbc, const t_graph *g,
145 real lambda, real *dvdlambda,
146 const t_mdatoms *md, t_fcdata *fcd,
147 int *global_atom_index);
149 real g96angles(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
150 const rvec x[], rvec4 f[], rvec fshift[],
151 const t_pbc *pbc, const t_graph *g,
152 real lambda, real *dvdlambda,
153 const t_mdatoms *md, t_fcdata *fcd,
154 int *global_atom_index);
156 real cross_bond_bond(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
157 const rvec x[], rvec4 f[], rvec fshift[],
158 const t_pbc *pbc, const t_graph *g,
159 real lambda, real *dvdlambda,
160 const t_mdatoms *md, t_fcdata *fcd,
161 int *global_atom_index);
163 real cross_bond_angle(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
164 const rvec x[], rvec4 f[], rvec fshift[],
165 const t_pbc *pbc, const t_graph *g,
166 real lambda, real *dvdlambda,
167 const t_mdatoms *md, t_fcdata *fcd,
168 int *global_atom_index);
170 real urey_bradley(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
171 const rvec x[], rvec4 f[], rvec fshift[],
172 const t_pbc *pbc, const t_graph *g,
173 real lambda, real *dvdlambda,
174 const t_mdatoms *md, t_fcdata *fcd,
175 int *global_atom_index);
177 real quartic_angles(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
178 const rvec x[], rvec4 f[], rvec fshift[],
179 const t_pbc *pbc, const t_graph *g,
180 real lambda, real *dvdlambda,
181 const t_mdatoms *md, t_fcdata *fcd,
182 int *global_atom_index);
184 real linear_angles(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
185 const rvec x[], rvec4 f[], rvec fshift[],
186 const t_pbc *pbc, const t_graph *g,
187 real lambda, real *dvdlambda,
188 const t_mdatoms *md, t_fcdata *fcd,
189 int *global_atom_index);
191 real restrangles(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
192 const rvec x[], rvec4 f[], rvec fshift[],
193 const t_pbc *pbc, const t_graph *g,
194 real lambda, real *dvdlambda,
195 const t_mdatoms *md, t_fcdata *fcd,
196 int *global_atom_index);
198 real pdihs(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
199 const rvec x[], rvec4 f[], rvec fshift[],
200 const t_pbc *pbc, const t_graph *g,
201 real lambda, real *dvdlambda,
202 const t_mdatoms *md, t_fcdata *fcd,
203 int *global_atom_index);
205 real idihs(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
206 const rvec x[], rvec4 f[], rvec fshift[],
207 const t_pbc *pbc, const t_graph *g,
208 real lambda, real *dvdlambda,
209 const t_mdatoms *md, t_fcdata *fcd,
210 int *global_atom_index);
212 real rbdihs(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
213 const rvec x[], rvec4 f[], rvec fshift[],
214 const t_pbc *pbc, const t_graph *g,
215 real lambda, real *dvdlambda,
216 const t_mdatoms *md, t_fcdata *fcd,
217 int *global_atom_index);
219 real restrdihs(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
220 const rvec x[], rvec4 f[], rvec fshift[],
221 const t_pbc *pbc, const t_graph *g,
222 real lambda, real *dvdlambda,
223 const t_mdatoms *md, t_fcdata *fcd,
224 int *global_atom_index);
226 real cbtdihs(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
227 const rvec x[], rvec4 f[], rvec fshift[],
228 const t_pbc *pbc, const t_graph *g,
229 real lambda, real *dvdlambda,
230 const t_mdatoms *md, t_fcdata *fcd,
231 int *global_atom_index);
233 real tab_bonds(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
234 const rvec x[], rvec4 f[], rvec fshift[],
235 const t_pbc *pbc, const t_graph *g,
236 real lambda, real *dvdlambda,
237 const t_mdatoms *md, t_fcdata *fcd,
238 int *global_atom_index);
240 real tab_angles(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
241 const rvec x[], rvec4 f[], rvec fshift[],
242 const t_pbc *pbc, const t_graph *g,
243 real lambda, real *dvdlambda,
244 const t_mdatoms *md, t_fcdata *fcd,
245 int *global_atom_index);
247 real tab_dihs(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
248 const rvec x[], rvec4 f[], rvec fshift[],
249 const t_pbc *pbc, const t_graph *g,
250 real lambda, real *dvdlambda,
251 const t_mdatoms *md, t_fcdata *fcd,
252 int *global_atom_index);
254 real polarize(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
255 const rvec x[], rvec4 f[], rvec fshift[],
256 const t_pbc *pbc, const t_graph *g,
257 real lambda, real *dvdlambda,
258 const t_mdatoms *md, t_fcdata *fcd,
259 int *global_atom_index);
261 real anharm_polarize(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
262 const rvec x[], rvec4 f[], rvec fshift[],
263 const t_pbc *pbc, const t_graph *g,
264 real lambda, real *dvdlambda,
265 const t_mdatoms *md, t_fcdata *fcd,
266 int *global_atom_index);
268 real water_pol(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
269 const rvec x[], rvec4 f[], rvec fshift[],
270 const t_pbc *pbc, const t_graph *g,
271 real lambda, real *dvdlambda,
272 const t_mdatoms *md, t_fcdata *fcd,
273 int *global_atom_index);
275 real thole_pol(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
276 const rvec x[], rvec4 f[], rvec fshift[],
277 const t_pbc *pbc, const t_graph *g,
278 real lambda, real *dvdlambda,
279 const t_mdatoms *md, t_fcdata *fcd,
280 int *global_atom_index);
282 real angres(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
283 const rvec x[], rvec4 f[], rvec fshift[],
284 const t_pbc *pbc, const t_graph *g,
285 real lambda, real *dvdlambda,
286 const t_mdatoms *md, t_fcdata *fcd,
287 int *global_atom_index);
289 real angresz(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
290 const rvec x[], rvec4 f[], rvec fshift[],
291 const t_pbc *pbc, const t_graph *g,
292 real lambda, real *dvdlambda,
293 const t_mdatoms *md, t_fcdata *fcd,
294 int *global_atom_index);
296 real dihres(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
297 const rvec x[], rvec4 f[], rvec fshift[],
298 const t_pbc *pbc, const t_graph *g,
299 real lambda, real *dvdlambda,
300 const t_mdatoms *md, t_fcdata *fcd,
301 int *global_atom_index);
303 real unimplemented(int nfa, const t_iatom forceatoms[], const t_iparams ip[],
304 const rvec x[], rvec4 f[], rvec fshift[],
305 const t_pbc *pbc, const t_graph *g,
306 real lambda, real *dvdlambda,
307 const t_mdatoms *md, t_fcdata *fcd,
308 int *global_atom_index);
311 /* As pdihs(), but without calculating energies and shift forces */
313 pdihs_noener(int nbonds,
314 const t_iatom forceatoms[], const t_iparams forceparams[],
315 const rvec x[], rvec4 f[],
316 const struct t_pbc gmx_unused *pbc,
317 const struct t_graph gmx_unused *g,
319 const t_mdatoms gmx_unused *md, t_fcdata gmx_unused *fcd,
320 int gmx_unused *global_atom_index);
322 /* TODO these declarations should be internal to the module */
324 /* As angles(), but using SIMD to calculate many angles at once.
325 * This routines does not calculate energies and shift forces.
328 angles_noener_simd(int nbonds,
329 const t_iatom forceatoms[], const t_iparams forceparams[],
330 const rvec x[], rvec4 f[],
331 const struct t_pbc *pbc,
332 const struct t_graph gmx_unused *g,
333 real gmx_unused lambda,
334 const t_mdatoms gmx_unused *md, t_fcdata gmx_unused *fcd,
335 int gmx_unused *global_atom_index);
337 /* As urey_bradley, but using SIMD to calculate many potentials at once.
338 * This routines does not calculate energies and shift forces.
340 void urey_bradley_noener_simd(int nbonds,
341 const t_iatom forceatoms[], const t_iparams forceparams[],
342 const rvec x[], rvec4 f[],
343 const t_pbc *pbc, const t_graph gmx_unused *g,
344 real gmx_unused lambda,
345 const t_mdatoms gmx_unused *md, t_fcdata gmx_unused *fcd,
346 int gmx_unused *global_atom_index);
348 /* As pdihs_noener(), but using SIMD to calculate many dihedrals at once. */
350 pdihs_noener_simd(int nbonds,
351 const t_iatom forceatoms[], const t_iparams forceparams[],
352 const rvec x[], rvec4 f[],
353 const struct t_pbc *pbc,
354 const struct t_graph gmx_unused *g,
355 real gmx_unused lambda,
356 const t_mdatoms gmx_unused *md, t_fcdata gmx_unused *fcd,
357 int gmx_unused *global_atom_index);
359 /* As rbdihs(), when not needing energy or shift force, using SIMD to calculate many dihedrals at once. */
361 rbdihs_noener_simd(int nbonds,
362 const t_iatom forceatoms[], const t_iparams forceparams[],
363 const rvec x[], rvec4 f[],
364 const struct t_pbc *pbc,
365 const struct t_graph gmx_unused *g,
366 real gmx_unused lambda,
367 const t_mdatoms gmx_unused *md, t_fcdata gmx_unused *fcd,
368 int gmx_unused *global_atom_index);