3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Gromacs Runs On Most of All Computer Systems
50 int glatnr(int *global_atom_index,int i);
51 /* Returns the global topology atom number belonging to local atom index i.
52 * This function is intended for writing ascii output
53 * and returns atom numbers starting at 1.
54 * When global_atom_index=NULL returns i+1.
57 void calc_bonds(FILE *fplog,const gmx_multisim_t *ms,
59 rvec x[],history_t *hist,
60 rvec f[],t_forcerec *fr,
61 const t_pbc *pbc,const t_graph *g,
62 gmx_enerdata_t *enerd,t_nrnb *nrnb,real *lambda,
64 t_fcdata *fcd,int *ddgatindex,
65 t_atomtypes *atype, gmx_genborn_t *born,
67 gmx_bool bPrintSepPot,gmx_large_int_t step);
69 * The function calc_bonds() calculates all bonded force interactions.
70 * The "bonds" are specified as follows:
72 * the total number of bonded interactions.
74 * specifies which atoms are involved in a bond of a certain
75 * type, see also struct t_idef.
76 * t_functype *functype
77 * defines for every bonded force type what type of function to
78 * use, see also struct t_idef.
79 * t_iparams *forceparams
80 * defines the parameters for every bond type, see also struct
83 * total potential energy split up over the function types.
85 * global atom number indices, should be NULL when not using DD.
86 * gmx_bool bPrintSepPot
87 * if TRUE print local potential and dVdlambda for each bonded type.
89 * used with bPrintSepPot
91 * the total potential energy (sum over epot).
94 void calc_bonds_lambda(FILE *fplog,
98 const t_pbc *pbc,const t_graph *g,
99 gmx_enerdata_t *enerd,t_nrnb *nrnb,
102 t_fcdata *fcd,int *global_atom_index);
103 /* As calc_bonds, but only determines the potential energy
104 * for the perturbed interactions.
105 * The shift forces in fr are not affected.
108 real posres(int nbonds,
109 const t_iatom forceatoms[],const t_iparams forceparams[],
110 const rvec x[],rvec f[],rvec vir_diag,
112 real lambda,real *dvdlambda,
113 int refcoord_scaling,int ePBC,rvec comA,rvec comB);
114 /* Position restraints require a different pbc treatment from other bondeds */
116 real bond_angle(const rvec xi,const rvec xj,const rvec xk,
118 rvec r_ij,rvec r_kj,real *costh,
119 int *t1,int *t2); /* out */
120 /* Calculate bond-angle. No PBC is taken into account (use mol-shift) */
122 real dih_angle(const rvec xi,const rvec xj,const rvec xk,const rvec xl,
124 rvec r_ij,rvec r_kj,rvec r_kl,rvec m,rvec n, /* out */
126 int *t1,int *t2,int *t3);
127 /* Calculate dihedral-angle. No PBC is taken into account (use mol-shift) */
129 void do_dih_fup(int i,int j,int k,int l,real ddphi,
130 rvec r_ij,rvec r_kj,rvec r_kl,
131 rvec m,rvec n,rvec f[],rvec fshift[],
132 const t_pbc *pbc,const t_graph *g,
133 const rvec *x,int t1,int t2,int t3);
134 /* Do an update of the forces for dihedral potentials */
136 void make_dp_periodic(real *dp);
137 /* make a dihedral fall in the range (-pi,pi) */
139 /*************************************************************************
141 * Bonded force functions
143 *************************************************************************/
144 t_ifunc bonds,g96bonds,morse_bonds,cubic_bonds,FENE_bonds,restraint_bonds;
145 t_ifunc angles,g96angles,cross_bond_bond,cross_bond_angle,urey_bradley,quartic_angles,linear_angles;
146 t_ifunc pdihs,idihs,rbdihs;
147 t_ifunc tab_bonds,tab_angles,tab_dihs;
148 t_ifunc polarize,anharm_polarize,water_pol,thole_pol,angres,angresz,dihres,unimplemented;
151 /* Initialize the setup for the bonded force buffer reduction
152 * over threads. This should be called each time the bonded setup
153 * changes; i.e. at start-up without domain decomposition and at DD.
155 void init_bonded_thread_force_reduction(t_forcerec *fr,
162 #endif /* _bondf_h */