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
54 extern int glatnr(int *global_atom_index
,int i
);
55 /* Returns the global topology atom number belonging to local atom index i.
56 * This function is intended for writing ascii output
57 * and returns atom numbers starting at 1.
58 * When global_atom_index=NULL returns i+1.
61 extern void calc_bonds(FILE *fplog
,const gmx_multisim_t
*ms
,
63 rvec x
[],history_t
*hist
,
64 rvec f
[],t_forcerec
*fr
,
65 const t_pbc
*pbc
,const t_graph
*g
,
66 gmx_enerdata_t
*enerd
,t_nrnb
*nrnb
,real lambda
,
68 t_fcdata
*fcd
,int *ddgatindex
,
69 t_atomtypes
*atype
, gmx_genborn_t
*born
,gmx_cmap_t
*cmap
,
70 bool bPrintSepPot
,gmx_step_t step
);
72 * The function calc_bonds() calculates all bonded force interactions.
73 * The "bonds" are specified as follows:
75 * the total number of bonded interactions.
77 * specifies which atoms are involved in a bond of a certain
78 * type, see also struct t_idef.
79 * t_functype *functype
80 * defines for every bonded force type what type of function to
81 * use, see also struct t_idef.
82 * t_iparams *forceparams
83 * defines the parameters for every bond type, see also struct
86 * total potential energy split up over the function types.
88 * global atom number indices, should be NULL when not using DD.
90 * if TRUE print local potential and dVdlambda for each bonded type.
92 * used with bPrintSepPot
94 * the total potential energy (sum over epot).
97 extern void calc_bonds_lambda(FILE *fplog
,
101 const t_pbc
*pbc
,const t_graph
*g
,
102 gmx_enerdata_t
*enerd
,t_nrnb
*nrnb
,
105 t_fcdata
*fcd
,int *global_atom_index
);
106 /* As calc_bonds, but only determines the potential energy
107 * for the perturbed interactions.
108 * The shift forces in fr are not affected.
111 extern real
posres(int nbonds
,
112 const t_iatom forceatoms
[],const t_iparams forceparams
[],
113 const rvec x
[],rvec f
[],rvec vir_diag
,
115 real lambda
,real
*dvdlambda
,
116 int refcoord_scaling
,int ePBC
,rvec comA
,rvec comB
);
117 /* Position restraints require a different pbc treatment from other bondeds */
119 extern real
bond_angle(const rvec xi
,const rvec xj
,const rvec xk
,
121 rvec r_ij
,rvec r_kj
,real
*costh
,
122 int *t1
,int *t2
); /* out */
123 /* Calculate bond-angle. No PBC is taken into account (use mol-shift) */
125 extern real
dih_angle(const rvec xi
,const rvec xj
,const rvec xk
,const rvec xl
,
127 rvec r_ij
,rvec r_kj
,rvec r_kl
,rvec m
,rvec n
, /* out */
128 real
*cos_phi
,real
*sign
,
129 int *t1
,int *t2
,int *t3
);
130 /* Calculate dihedral-angle. No PBC is taken into account (use mol-shift) */
132 extern void do_dih_fup(int i
,int j
,int k
,int l
,real ddphi
,
133 rvec r_ij
,rvec r_kj
,rvec r_kl
,
134 rvec m
,rvec n
,rvec f
[],rvec fshift
[],
135 const t_pbc
*pbc
,const t_graph
*g
,
136 const rvec
*x
,int t1
,int t2
,int t3
);
137 /* Do an update of the forces for dihedral potentials */
139 /*************************************************************************
141 * Bonded force functions
143 *************************************************************************/
144 extern t_ifunc bonds
,g96bonds
,morse_bonds
,cubic_bonds
,FENE_bonds
;
145 extern t_ifunc angles
,g96angles
,cross_bond_bond
,cross_bond_angle
,urey_bradley
,quartic_angles
;
146 extern t_ifunc pdihs
,idihs
,rbdihs
;
147 extern t_ifunc tab_bonds
,tab_angles
,tab_dihs
;
148 extern t_ifunc polarize
,water_pol
,thole_pol
,angres
,angresz
,unimplemented
;
154 #endif /* _bondf_h */