Main Page | Class List | File List | Class Members | File Members

compute.h File Reference

Lower-level routines for force computation. More...

Go to the source code of this file.

Functions

int force_compute_scaled_coords (Force *, const MD_Dvec pos[], const int32 sel[], int32 sel_len)
 Compute scaled coordinates.
int force_compute_domain_update (Force *, const int32 sel[], int32 sel_len)
 Update nonperiodic domain specification.
int force_compute_bonds (Force *fobj, double *u_bond, MD_Dvec f_bond[], double e_bond[], double virial[], const MD_Dvec pos[], const int32 bond_sel[], int32 bond_sel_len)
 Compute contribution of bond springs to force and potential.
int force_compute_angles (Force *fobj, double *u_angle, MD_Dvec f_angle[], double e_angle[], double virial[], const MD_Dvec pos[], const int32 angle_sel[], int32 angle_sel_len)
 Compute contribution of bond angles to force and potential.
int force_compute_dihedrals (Force *fobj, double *u_dihed, MD_Dvec f_dihed[], double e_dihed[], double virial[], const MD_Dvec pos[], const int32 dihed_sel[], int32 dihed_sel_len)
 Compute contribution of dihedrals to force and potential.
int force_compute_impropers (Force *fobj, double *u_impr, MD_Dvec f_impr[], double e_impr[], double virial[], const MD_Dvec pos[], const int32 impr_sel[], int32 impr_sel_len)
 Compute contribution of impropers to force and potential.
double force_compute_bond_interaction (MD_Dvec f_bond[], double virial[], const MD_Dvec pos[], const MD_Bond *bond, const MD_BondPrm *prm)
 Compute an individual bond spring interaction.
double force_compute_angle_interaction (MD_Dvec f_angle[], double virial[], const MD_Dvec pos[], const MD_Angle *angle, const MD_AnglePrm *prm)
 Compute an individual bond angle interaction.
double force_compute_torsion_interaction (MD_Dvec f_tors[], double virial[], const MD_Dvec pos[], const MD_Tors *tors, const MD_TorsPrm *prm)
 Compute an individual torsion angle interaction.
int force_compute_bres_sphere (Force *fobj, double *u_bres, MD_Dvec f_bres[], double e_bres[], const MD_Dvec *pos, const int32 bres_sel[], int32 bres_sel_len)
 Compute spherical boundary restraint.
int force_compute_bres_cylinder (Force *fobj, double *u_bres, MD_Dvec f_bres[], double e_bres[], const MD_Dvec *pos, const int32 *bres_sel, int32 bres_sel_len, int32 bresopts)
 Compute cylindrical boundary restraint.
int force_compute_nbpairs_direct (Force *fobj, double virial[], double *u_elec, MD_Dvec f_elec[], double e_elec[], double e_epot[], int32 is_elec_direct, int32 elec_pair_potential, double *u_vdw, MD_Dvec f_vdw[], double e_vdw[], int32 is_vdw_direct, int32 vdw_pair_potential, const MD_Dvec pos[], const int32 aset_sel[], int32 aset_sel_len, const int32 bset_sel[], int32 bset_sel_len)
 Use direct all pairs algorithm to compute nonbonded interactions.
int force_compute_nbpairs_geomhash (Force *fobj, const MD_Dvec trpos[], const int32 sel[], int32 sel_len)
 Apply geometric hashing of atomic positions into grid cells.
int force_compute_nbpairs_gridcells (Force *fobj, double virial[], double *u_elec, MD_Dvec f_elec[], double e_elec[], double e_epot[], int32 is_elec_gridcells, int32 elec_pair_potential, double *u_vdw, MD_Dvec f_vdw[], double e_vdw[], int32 is_vdw_gridcells, int32 vdw_pair_potential, const MD_Dvec pos[], int32 is_subtracted)
 Compute cutoff nonbonded interactions using grid cells.
int force_compute_nbpairs_subtexcl (Force *fobj, double virial[], double *u_elec, MD_Dvec f_elec[], double e_elec[], double e_epot[], int32 is_elec_subtexcl, int32 elec_pair_potential, double *u_vdw, MD_Dvec f_vdw[], double e_vdw[], int32 is_vdw_subtexcl, int32 vdw_pair_potential, const MD_Dvec pos[], const int32 sel[], int32 sel_len, const int32 mapnb[], int32 map_id)
 Subtract excluded interactions for nonbonded pairwise interactions.
void force_compute_nbpairs_elec_standard (double *u, double *du_r, double r2, double c)
 Compute a standard electrostatic interaction.
void force_compute_nbpairs_elec_shifted (double *u, double *du_r, double r2, double c, double inv_elec_cutoff2)
 Compute a shifted electrostatic interaction.
void force_compute_nbpairs_elec_ewald (double *u, double *du_r, double r2, double c, double ewald_coef, double grad_coef)
void force_compute_nbpairs_vdw_standard (double *u, double *du_r, double r2, double a, double b)
 Compute a standard van der Waals interaction.
void force_compute_nbpairs_vdw_switched (double *u, double *du_r, double r2, double a, double b, double vdw_cutoff2, double switch_dist2, double inv_denom_switch)
 Compute a switched van der Waals interaction.
void force_compute_nbpairs_vdw_buck (double *u, double *du_r, double r2, double a, double b, double c)
 Compute Buckingham potential (replaces van der Waals).
void force_compute_nbpairs_vdw_switchbuck (double *u, double *du_r, double r2, double a, double b, double c, double vdw_cutoff2, double switch_dist2, double inv_denom_switch)
 Compute switched Buckingham potential (replaces van der Waals).
void force_compute_nbpairs_vdw_bucknd (double *u, double *du_r, double r2, double a, double b)
 Compute Buckingham without dispersion term potential (replaces van der Waals).
void force_compute_nbpairs_vdw_switchbucknd (double *u, double *du_r, double r2, double a, double b, double vdw_cutoff2, double switch_dist2, double inv_denom_switch)
 Compute switched Buckingham without dispersion term potential (replaces van der Waals).
void force_compute_nbpairs_vdw_bucksafe (double *u, double *du_r, double r2, double a, double b, double c, double rn2, double an, double bn)
void force_compute_nbpairs_vdw_switchbucksafe (double *u, double *du_r, double r2, double a, double b, double c, double rn2, double an, double bn, double roff2, double ron2, double denom)
int force_compute_nbpairs_isregen_pairlists (Force *fobj, double delta2, const MD_Dvec initpos[], const MD_Dvec pos[], const int32 sel[], int32 sel_len)
int force_compute_nbpairs_regen_pairlists (Force *fobj, MD_Dvec initpos[], const MD_Dvec pos[], const int32 aset_sel[], int32 aset_sel_len, const int32 mapnb[], int32 aset_id, int32 bset_id)
int force_compute_nbpairs_eval_pairlists (Force *fobj, double virial[], double *u_elec, MD_Dvec f_elec[], double e_elec[], double e_epot[], int32 is_elec_gridcells, int32 elec_pair_potential, double *u_vdw, MD_Dvec f_vdw[], double e_vdw[], int32 is_vdw_gridcells, int32 vdw_pair_potential, const MD_Dvec pos[], const int32 sel[], int32 sel_len, const int32 mapnb[], int32 map_id)


Detailed Description

Lower-level routines for force computation.

Author:
David J. Hardy
Date:
2006
The high-level force_compute() routine is the manager that calls the lower-level routines prototyped in this file. The various force computation routines are called depending on the force types indicated in ForceParam_t::forcetypes. The various array and atom/bond selection parameters are determined by the ForceResult_t object passed to force_compute() and the ForceSelect_t object provided to the force_create() constructor.

The medium-level routines:

can all still be viewed as methods of class Force_t, but with all of the array buffer space and data (from ForceResult_t and ForceSelect_t objects) explicitly passed as parameters along with any "control" data originally specified as part of ForceParam_t::forcetypes, etc.

The low-level routines:

compute various types of individual interactions independently of the force and helper classes.

Function Documentation

double force_compute_angle_interaction MD_Dvec  f_angle[],
double  virial[],
const MD_Dvec  pos[],
const MD_Angle *  angle,
const MD_AnglePrm *  prm
 

Compute an individual bond angle interaction.

Parameters:
[in,out] f_angle Receives forces on atoms.
[in,out] virial Accumulates contribution to virial.
[in] pos Atomic position coordinates.
[in] bond Contains the angle.
[in] prm Contains the parameters for this angle.
Compute force and potential for this individual bond angle. The angle contains atom indices for the participating atoms, so the entire f_angle and pos arrays are passed. The prm contains the force constants needed for this bond angle.

Returns:
the potential for this interaction.

int force_compute_angles Force fobj,
double *  u_angle,
MD_Dvec  f_angle[],
double  e_angle[],
double  virial[],
const MD_Dvec  pos[],
const int32  angle_sel[],
int32  angle_sel_len
 

Compute contribution of bond angles to force and potential.

Parameters:
[out] u_angle Receives total potential energy.
[in,out] f_angle Receives forces on atoms.
[out] e_angle Receives interaction potentials.
[in,out] virial Accumulates contribution to virial.
[in] pos Atomic position coordinates.
[in] angle_sel Selects angles from fobj->param->angle array.
[in] angle_sel_len Length of angle selection array.
Compute force and potentials due to bond angles. All parameters must exist: array buffers must all be allocated and be of expected length and data arrays must contain the expected information. The angle_sel array contains integer indices into the fobj->param->angle array. The output arrays are expected to be in a sensible state upon calling (i.e. these arrays are not zeroed). In particular, the forces are accumulated into f_angle.

Returns:
0 for success or FORCE_FAIL on failure.

double force_compute_bond_interaction MD_Dvec  f_bond[],
double  virial[],
const MD_Dvec  pos[],
const MD_Bond *  bond,
const MD_BondPrm *  prm
 

Compute an individual bond spring interaction.

Parameters:
[in,out] f_bond Receives forces on atoms.
[in,out] virial Accumulates contribution to virial.
[in] pos Atomic position coordinates.
[in] bond Contains the bond.
[in] prm Contains the parameters for this bond.
Compute force and potential for this individual spring bond. The bond contains atom indices for the participating atoms, so the entire f_bond and pos arrays are passed. The prm contains the force constants needed for this bond.

Returns:
the potential for this interaction.

int force_compute_bonds Force fobj,
double *  u_bond,
MD_Dvec  f_bond[],
double  e_bond[],
double  virial[],
const MD_Dvec  pos[],
const int32  bond_sel[],
int32  bond_sel_len
 

Compute contribution of bond springs to force and potential.

Parameters:
[out] u_bond Receives total potential energy.
[in,out] f_bond Receives forces on atoms.
[out] e_bond Receives interaction potentials.
[in,out] virial Accumulates contribution to virial.
[in] pos Atomic position coordinates.
[in] bond_sel Selects bonds from fobj->param->bond array.
[in] bond_sel_len Length of bond selection array.
Compute force and potentials due to bond springs. All parameters must exist: array buffers must all be allocated and be of expected length and data arrays must contain the expected information. The bond_sel array contains integer indices into the fobj->param->bond array. The output arrays are expected to be in a sensible state upon calling (i.e. these arrays are not zeroed). In particular, the forces are accumulated into f_bond.

Returns:
0 for success or FORCE_FAIL on failure.

int force_compute_bres_cylinder Force fobj,
double *  u_bres,
MD_Dvec  f_bres[],
double  e_bres[],
const MD_Dvec *  pos,
const int32 *  bres_sel,
int32  bres_sel_len,
int32  bresopts
 

Compute cylindrical boundary restraint.

Parameters:
[out] u_bres Receives total potential energy.
[in,out] f_bres Receives forces on atoms.
[out] e_bres Receives interaction potentials.
[in] pos Atomic position coordinates.
[in] bres_sel Selects atoms from fobj->param->atom.
[in] bres_sel_len Length of atom selection array.
[in] bresopts Gives orientation of cylinder.
Compute force and potentials due to cylindrical boundary restraint. All arguments must exist: array buffers must all be allocated and be of the expected length and data arrays must contain the expected information. The bres_sel array contains integer indices into the fobj->param->atom array. The output arrays are expected to be in a sensible state upon calling (i.e. these arrays are not zeroed). In particular, the forces are accumulated into f_bres. The bresopts parameter gives orientation of cylinder, with value one of FORCE_BRES_X_CYLINDER, FORCE_BRES_Y_CYLINDER, or FORCE_BRES_Z_CYLINDER.

Returns:
0 for success or FORCE_FAIL on failure.

int force_compute_bres_sphere Force fobj,
double *  u_bres,
MD_Dvec  f_bres[],
double  e_bres[],
const MD_Dvec *  pos,
const int32  bres_sel[],
int32  bres_sel_len
 

Compute spherical boundary restraint.

Parameters:
[out] u_bres Receives total potential energy.
[in,out] f_bres Receives forces on atoms.
[out] e_bres Receives interaction potentials.
[in] pos Atomic position coordinates.
[in] bres_sel Selects atoms from fobj->param->atom.
[in] bres_sel_len Length of atom selection array.
Compute force and potentials due to spherical boundary restraint. All arguments must exist: array buffers must all be allocated and be of the expected length and data arrays must contain the expected information. The bres_sel array contains integer indices into the fobj->param->atom array. The output arrays are expected to be in a sensible state upon calling (i.e. these arrays are not zeroed). In particular, the forces are accumulated into f_bres.

Returns:
0 for success or FORCE_FAIL on failure.

int force_compute_dihedrals Force fobj,
double *  u_dihed,
MD_Dvec  f_dihed[],
double  e_dihed[],
double  virial[],
const MD_Dvec  pos[],
const int32  dihed_sel[],
int32  dihed_sel_len
 

Compute contribution of dihedrals to force and potential.

Parameters:
[out] u_dihed Receives total potential energy.
[in,out] f_dihed Receives forces on atoms.
[out] e_dihed Receives interaction potentials.
[in,out] virial Accumulates contribution to virial.
[in] pos Atomic position coordinates.
[in] dihed_sel Selects dihedrals from fobj->param->dihed.
[in] dihed_sel_len Length of dihedral selection array.
Compute force and potentials due to dihedrals. All parameters must exist: array buffers must all be allocated and be of expected length and data arrays must contain the expected information. The dihed_sel array contains integer indices into the fobj->param->dihed array. The output arrays are expected to be in a sensible state upon calling (i.e. these arrays are not zeroed). In particular, the forces are accumulated into f_dihed.

Returns:
0 for success or FORCE_FAIL on failure.

int force_compute_domain_update Force ,
const int32  sel[],
int32  sel_len
 

Update nonperiodic domain specification.

Parameters:
[in] sel Selects atoms from fobj->param->atoms array.
[in] sel_len Length of atom selection array.
This routine is called internally to update the nonperiodic domain specification. The cell basis vectors for a nonperiodic domain are chosen to map most of the atoms into the unit cube (but the choice of these vectors is arbitrary). During the simulation, the atoms might drift away from the initial domain specification, so this routine is called every Force_t::max_steps time intervals. If the center (average position coordinates) has drifted outside of the domain, a new domain specification needs to be determined.

If the user previously called force_setup_lattice() to fix a lattice in space, then this routine chooses the new center to be some multiple of the lattice spacing away from the old center, and does not rescale the domain basis vector lengths, and so the location of the points of the (infinite) lattice is preserved in space.

Returns:
0 if no change is to be made (the average over all positions is not sufficiently far away from the current center) or if a fixed lattice has been defined (in which case the domain center is updated to modify the translation vector, but the transformation to the unit cube remains unchanged); otherwise, the domain center is updated and 1 is returned to indicate that force_setup_domain() needs to be called at the beginning of the next force_compute().

int force_compute_impropers Force fobj,
double *  u_impr,
MD_Dvec  f_impr[],
double  e_impr[],
double  virial[],
const MD_Dvec  pos[],
const int32  impr_sel[],
int32  impr_sel_len
 

Compute contribution of impropers to force and potential.

Parameters:
[out] u_impr Receives total potential energy.
[in,out] f_impr Receives forces on atoms.
[out] e_impr Receives interaction potentials.
[in,out] virial Accumulates contribution to virial.
[in] pos Atomic position coordinates.
[in] impr_sel Selects impropers from fobj->param->impr.
[in] impr_sel_len Length of improper selection array.
Compute force and potentials due to impropers. All parameters must exist: array buffers must all be allocated and be of expected length and data arrays must contain the expected information. The impr_sel array contains integer indices into the fobj->param->impr array. The output arrays are expected to be in a sensible state upon calling (i.e. these arrays are not zeroed). In particular, the forces are accumulated into f_impr.

Returns:
0 for success or FORCE_FAIL on failure.

int force_compute_nbpairs_direct Force fobj,
double  virial[],
double *  u_elec,
MD_Dvec  f_elec[],
double  e_elec[],
double  e_epot[],
int32  is_elec_direct,
int32  elec_pair_potential,
double *  u_vdw,
MD_Dvec  f_vdw[],
double  e_vdw[],
int32  is_vdw_direct,
int32  vdw_pair_potential,
const MD_Dvec  pos[],
const int32  aset_sel[],
int32  aset_sel_len,
const int32  bset_sel[],
int32  bset_sel_len
 

Use direct all pairs algorithm to compute nonbonded interactions.

Parameters:
[in,out] virial Accumulates nonbonded contribution to virial.
[out] u_elec Receives total electrostatic potential energy.
[in,out] f_elec Receives electrostatic forces on atoms.
[in,out] e_elec Receives electrostatic interaction energies.
[in,out] e_epot Receives electrostatic potentials.
[in] is_elec_direct Gives on/off status for electrostatics.
[in] elec_pair_potential Gives electrostatic pairwise potential.
[out] u_vdw Receives total van der Waals potential energy.
[in,out] f_vdw Receives van der Waals forces on atoms.
[in,out] e_vdw Receives van der Waals interaction energies.
[in] is_vdw_direct Gives on/off status for van der Waals.
[in] vdw_pair_potential Gives van der Waals pairwise potential.
[in] pos Atomic position coordinates.
[in] aset_sel Set A selection of atoms.
[in] aset_sel_len Length of set A atom selection array.
[in] bset_sel Set B selection of atoms.
[in] bset_sel_len Length of set B atom selection array.
Compute nonbonded electrostatics and/or van der Waals using simple direct all-pairs algorithm. All arguments must exist: array buffers must all be allocated and be of the expected length and data arrays must contain the expected information. The two selection sets should either point to the same array or to two nonempty arrays containing no overlapping indices.

Set is_elec_direct to nonzero to compute electrostatic interactions. Set is_vdw_direct to nonzero to compute van der Waals interactions. The elec_pair_potential parameter tells specifically what electrostatic pairwise potential to use (e.g. FORCE_ELEC_STANDARD is standard 1/r potential, FORCE_ELEC_SHIFTED is the 1/r potential multiplied by a shifting function for smooth cutoff). Similarly, the vdw_pair_potential parameter tells specifically what van der Waals pairwise potential to use (e.g. FORCE_VDW_STANDARD is the standard potential, FORCE_VDW_SWITCHED is uses a switching function to smoothly truncate the potential at the cutoff).

Returns:
0 for success or FORCE_FAIL on failure.

void force_compute_nbpairs_elec_shifted double *  u,
double *  du_r,
double  r2,
double  c,
double  inv_elec_cutoff2
 

Compute a shifted electrostatic interaction.

Parameters:
[out] u Receives potential energy.
[out] du_r Receives 1/r times derivative of potential energy.
[in] r2 Square of pairwise distance.
[in] c Multiplicative constant.
[in] inv_elec_cutoff2 Inverse of square of cutoff distance.
Computes the shifted electrostatic potential and force scaling. This can be expressed in the form (1/r) s(r), where s(r) is the shifting function that is parameterized by inv_elec_cutoff2. The exact function is documented in the NAMD User Guide, as well as the X-Plor and CHARMM literature.

void force_compute_nbpairs_elec_standard double *  u,
double *  du_r,
double  r2,
double  c
 

Compute a standard electrostatic interaction.

Parameters:
[out] u Receives potential energy.
[out] du_r Receives 1/r times derivative of potential energy.
[in] r2 Square of pairwise distance.
[in] c Multiplicative constant.
Computes the standard c/r potential and the force scaling given by D(c/r)/r.

int force_compute_nbpairs_geomhash Force fobj,
const MD_Dvec  trpos[],
const int32  sel[],
int32  sel_len
 

Apply geometric hashing of atomic positions into grid cells.

Parameters:
[in] trpos Scaled atomic position coordinates.
[in] sel Select atom indices into fobj->param->atom.
[in] sel_len Length of atom selection array.
Atoms are hashed into grid cell linked lists using scaled atomic coordinates. Must be called before calling force_compute_nbpairs_gridcells(). All arguments must exist: data arrays must contain the expected information.

Returns:
0 for success or FORCE_FAIL on failure.

int force_compute_nbpairs_gridcells Force fobj,
double  virial[],
double *  u_elec,
MD_Dvec  f_elec[],
double  e_elec[],
double  e_epot[],
int32  is_elec_gridcells,
int32  elec_pair_potential,
double *  u_vdw,
MD_Dvec  f_vdw[],
double  e_vdw[],
int32  is_vdw_gridcells,
int32  vdw_pair_potential,
const MD_Dvec  pos[],
int32  is_subtracted
 

Compute cutoff nonbonded interactions using grid cells.

Parameters:
[in,out] virial Accumulates nonbonded contribution to virial.
[out] u_elec Receives total electrostatic potential energy.
[in,out] f_elec Receives electrostatic forces on atoms.
[in,out] e_elec Receives electrostatic interaction energies.
[in,out] e_epot Receives electrostatic potentials.
[in] is_elec_gridcells Gives on/off status for electrostatics.
[in] elec_pair_potential Gives electrostatic pairwise potential.
[out] u_vdw Receives total van der Waals potential energy.
[in,out] f_vdw Receives van der Waals forces on atoms.
[in,out] e_vdw Receives van der Waals interaction energies.
[in] is_vdw_gridcells Gives on/off status for van der Waals.
[in] vdw_pair_potential Gives van der Waals pairwise potential.
[in] pos Atomic position coordinates.
[in] is_subtracted Tells to add (zero) or subtract (nonzero) interactions to force and potential arrays.
Compute cutoff nonbonded electrostatics and/or van der Waals using grid cell hashing, with cost linear in the number of atoms. Requires calling force_compute_nbpairs_geomhash() first to fill the grid cells. All arguments must exist: array buffers must all be allocated and be of the expected length and data arrays must contain the expected information. The selection of positions was done during the preceding hashing.

Set is_elec_gridcells to nonzero to compute electrostatic interactions. Set is_vdw_gridcells to nonzero to compute van der Waals interactions. The elec_pair_potential parameter tells specifically what electrostatic pairwise potential to use (e.g. FORCE_ELEC_STANDARD is standard 1/r potential, FORCE_ELEC_SHIFTED is the 1/r potential multiplied by a shifting function for smooth cutoff). Similarly, the vdw_pair_potential parameter tells specifically what van der Waals pairwise potential to use (e.g. FORCE_VDW_STANDARD is the standard potential, FORCE_VDW_SWITCHED is uses a switching function to smoothly truncate the potential at the cutoff).

The is_subtracted parameter indicates whether interactions should be added (value: zero) or subtracted (value: nonzero) to the force and potential arrays. (Note: this is a kludge so that grid cell hashing can be applied conveniently to compute potential differences between disjoint sets of atoms. The idea is to add interactions involving combined sets, then subtract out separate contributions each set individually.)

Returns:
0 for success or FORCE_FAIL on failure.

int force_compute_nbpairs_subtexcl Force fobj,
double  virial[],
double *  u_elec,
MD_Dvec  f_elec[],
double  e_elec[],
double  e_epot[],
int32  is_elec_subtexcl,
int32  elec_pair_potential,
double *  u_vdw,
MD_Dvec  f_vdw[],
double  e_vdw[],
int32  is_vdw_subtexcl,
int32  vdw_pair_potential,
const MD_Dvec  pos[],
const int32  sel[],
int32  sel_len,
const int32  mapnb[],
int32  map_id
 

Subtract excluded interactions for nonbonded pairwise interactions.

Parameters:
[in,out] virial Accumulates nonbonded contribution to virial.
[out] u_elec Receives total electrostatic potential energy.
[in,out] f_elec Receives electrostatic forces on atoms.
[in,out] e_elec Receives electrostatic interaction energies.
[in,out] e_epot Receives electrostatic potentials.
[in] is_elec_subtexcl Gives on/off status for electrostatics.
[in] elec_pair_potential Gives electrostatic pairwise potential.
[out] u_vdw Receives total van der Waals potential energy.
[in,out] f_vdw Receives van der Waals forces on atoms.
[in,out] e_vdw Receives van der Waals interaction energies.
[in] is_vdw_subtexcl Gives on/off status for van der Waals.
[in] vdw_pair_potential Gives van der Waals pairwise potential.
[in] pos Atomic position coordinates.
[in] sel Select atom indices into fobj->param->atom.
[in] sel_len Length of atom selection arrays.
[in] mapnb Per atom map indicating set ownership.
[in] map_id ID of set to match against for exclusion.
Subtract excluded interactions for nonbonded electrostatics and/or van der Waals. All arguments must exist: array buffers must all be allocated and be of the expected length and data arrays must contain the expected information.

The atom selection array sel gives the indices of a set of atoms in which to consider exclusions. The second atom in the excluded interaction pair is found by considering the exclusion list for each atom in sel, and matching these against set ownership given by mapnb. Atoms matching map_id are taken to be the second atom in the exclusion.

Here is how the above is applied in practice: for standard computation of nonbonded interactions involving all atoms, then sel selects all atoms (the identity map, 0..N-1) and mapnb contains (FORCE_SELECT_ASET | FORCE_SELECT_BSET) meaning that every atom belongs to both sets, so that is the value to assign to map_id. When computing potential differences between two disjoint subsets, the excluded interactions should span the two sets, so sel should take on one set, say, ForceSelect_t::aset_sel and map_id should be assigned the other set FORCE_SELECT_BSET.

Set is_elec_subtexcl to nonzero to compute electrostatic interactions. Set is_vdw_subtexcl to nonzero to compute van der Waals interactions. The elec_pair_potential parameter tells specifically what electrostatic pairwise potential to use (e.g. FORCE_ELEC_STANDARD is standard 1/r potential, FORCE_ELEC_SHIFTED is the 1/r potential multiplied by a shifting function for smooth cutoff). Similarly, the vdw_pair_potential parameter tells specifically what van der Waals pairwise potential to use (e.g. FORCE_VDW_STANDARD is the standard potential, FORCE_VDW_SWITCHED is uses a switching function to smoothly truncate the potential at the cutoff).

Returns:
0 for success or FORCE_FAIL on failure.

void force_compute_nbpairs_vdw_buck double *  u,
double *  du_r,
double  r2,
double  a,
double  b,
double  c
 

Compute Buckingham potential (replaces van der Waals).

Parameters:
[out] u Receives potential energy.
[out] du_r Receives 1/r times derivative of potential energy.
[in] r2 Square of pairwise distance.
[in] a constant
[in] b constant
[in] c constant
Compute Buckingham form: a exp(-r/b) - c/r^6

See Flikkema & Bromley, Chem Phys Lett 378 (2003) 622-629, for parameterization for silica.

void force_compute_nbpairs_vdw_bucknd double *  u,
double *  du_r,
double  r2,
double  a,
double  b
 

Compute Buckingham without dispersion term potential (replaces van der Waals).

Parameters:
[out] u Receives potential energy.
[out] du_r Receives 1/r times derivative of potential energy.
[in] r2 Square of pairwise distance.
[in] a constant
[in] b constant
Compute Buckingham without dispersion form: a exp(-r/b)

void force_compute_nbpairs_vdw_standard double *  u,
double *  du_r,
double  r2,
double  a,
double  b
 

Compute a standard van der Waals interaction.

Parameters:
[out] u Receives potential energy.
[out] du_r Receives 1/r times derivative of potential energy.
[in] r2 Square of pairwise distance.
[in] a Multiplicative constant for r^{-12} term.
[in] b Multiplicative constant for r^{-6} term.
Computes the standard van der Waals potential and force scaling.

void force_compute_nbpairs_vdw_switchbuck double *  u,
double *  du_r,
double  r2,
double  a,
double  b,
double  c,
double  vdw_cutoff2,
double  switch_dist2,
double  inv_denom_switch
 

Compute switched Buckingham potential (replaces van der Waals).

Parameters:
[out] u Receives potential energy.
[out] du_r Receives 1/r times derivative of potential energy.
[in] r2 Square of pairwise distance.
[in] a constant
[in] b constant
[in] c constant
[in] vdw_cutoff2 Square of cutoff distance.
[in] switch_dist2 Square of switching distance.
[in] inv_denom_switch Inverse of switching function denominator.
Compute Buckingham form: (a exp(-r/b) - c/r^6) s(r), where s(r) is same switching function as used for van der Waals.

See Flikkema & Bromley, Chem Phys Lett 378 (2003) 622-629, for parameterization for silica.

void force_compute_nbpairs_vdw_switchbucknd double *  u,
double *  du_r,
double  r2,
double  a,
double  b,
double  vdw_cutoff2,
double  switch_dist2,
double  inv_denom_switch
 

Compute switched Buckingham without dispersion term potential (replaces van der Waals).

Parameters:
[out] u Receives potential energy.
[out] du_r Receives 1/r times derivative of potential energy.
[in] r2 Square of pairwise distance.
[in] a constant
[in] b constant
[in] vdw_cutoff2 Square of cutoff distance.
[in] switch_dist2 Square of switching distance.
[in] inv_denom_switch Inverse of switching function denominator.
Compute Buckingham form: a exp(-r/b) s(r), where s(r) is same switching function as used for van der Waals.

void force_compute_nbpairs_vdw_switched double *  u,
double *  du_r,
double  r2,
double  a,
double  b,
double  vdw_cutoff2,
double  switch_dist2,
double  inv_denom_switch
 

Compute a switched van der Waals interaction.

Parameters:
[out] u Receives potential energy.
[out] du_r Receives 1/r times derivative of potential energy.
[in] r2 Square of pairwise distance.
[in] a Multiplicative constant for r^{-12} term.
[in] b Multiplicative constant for r^{-6} term.
[in] vdw_cutoff2 Square of cutoff distance.
[in] switch_dist2 Square of switching distance.
[in] inv_denom_switch Inverse of switching function denominator.
Computes the switched van der Waals potential and force scaling. This is piecewise defined, where it is the standard potential for pairwise distance less than switching distance, or the standard potential times the switching function for pairwise distance between the switching distance and cutoff, and zero beyond the cutoff. The switching function gives C^1 joins and is documented in the NAMD User Guide, as well as the X-Plor and CHARMM literature.

int force_compute_scaled_coords Force ,
const MD_Dvec  pos[],
const int32  sel[],
int32  sel_len
 

Compute scaled coordinates.

Parameters:
[in] pos Atomic position coordinates.
[in] sel Selects atoms from fobj->param->atoms array.
[in] sel_len Length of atom selection array.
Compute the scaled atomic coordinates (Force_t::trpos) and position wrapping array (Force_t::poswrap) for periodic boundaries. These arrays are held internally by the force object. For periodic boundaries, it is assumed that the pos[i]+poswrap[i] result from previous computation is no more than one cell image away. This is ensured by monitoring atom speeds so that an atom travels less than one cell width in a time step.

Returns:
0 for success or FORCE_FAIL on failure.

double force_compute_torsion_interaction MD_Dvec  f_tors[],
double  virial[],
const MD_Dvec  pos[],
const MD_Tors *  tors,
const MD_TorsPrm *  prm
 

Compute an individual torsion angle interaction.

Parameters:
[in,out] f_tors Receives forces on atoms.
[in,out] virial Accumulates contribution to virial.
[in] pos Atomic position coordinates.
[in] tors Contains the torsion.
[in] prm Contains the parameters for this torsion.
Compute force and potential for this individual torsion angle. (This is used for contributions of both dihedrals and impropers.) The tors contains atom indices for the participating atoms, so the entire f_tors and pos arrays are passed. The prm contains the force constants needed for this torsion angle.

Returns:
the potential for this interaction.


Generated on Thu Feb 7 18:11:41 2008 for MDX by  doxygen 1.3.9.1