STORMM Source Documentation
|
Information need for bonded calculations. Templating is used to serve either of two levels of precision: single (which, on almost any machine and HPC accelerator card, will be fp32) or double (fp64). Parameters for bonds, the highest frequency terms in any system by a factor of roughly 5, will always be stored in double precision, as will parameters for bond angles, the next highest frequency group of terms. The precision of other parameters varies according to the template parameter. More...
#include <atomgraph_abstracts.h>
Public Member Functions | |
ValenceKit (int natom_in, int nbond_in, int nangl_in, int ndihe_in, int nbond_param_in, int nangl_param_in, int ndihe_param_in, int ninfr14_in, int nattn14_param_in, int nubrd_in, int ncimp_in, int ncmap_in, int nubrd_param_in, int ncimp_param_in, int ncmap_surf_in, const T *bond_keq_in, const T *bond_leq_in, const T *angl_keq_in, const T *angl_theta_in, const T *dihe_amp_in, const T *dihe_freq_in, const T *dihe_phi_in, const T *attn14_elec_in, const T *attn14_vdw_in, const int *bond_i_atoms_in, const int *bond_j_atoms_in, const int *bond_param_idx_in, const char4 *bond_modifiers_in, const uint *bond_relevance_in, const int *angl_i_atoms_in, const int *angl_j_atoms_in, const int *angl_k_atoms_in, const int *angl_param_idx_in, const char4 *angl_modifiers_in, const uint *angl_relevance_in, const int *dihe_i_atoms_in, const int *dihe_j_atoms_in, const int *dihe_k_atoms_in, const int *dihe_l_atoms_in, const int *dihe_param_idx_in, const int *dihe14_param_idx_in, const char4 *dihe_modifiers_in, const int *infr14_i_atoms_in, const int *infr14_l_atoms_in, const int *infr14_param_idx_in, const int *ubrd_i_atoms_in, const int *ubrd_k_atoms_in, const int *ubrd_param_idx_in, const int *cimp_i_atoms_in, const int *cimp_j_atoms_in, const int *cimp_k_atoms_in, const int *cimp_l_atoms_in, const int *cimp_param_idx_in, const int *cmap_i_atoms_in, const int *cmap_j_atoms_in, const int *cmap_k_atoms_in, const int *cmap_l_atoms_in, const int *cmap_m_atoms_in, const int *cmap_dim_in, const int *cmap_surf_bounds_in, const int *cmap_patch_bounds_in, const int *cmap_surf_idx_in, const T *ubrd_keq_in, const T *ubrd_leq_in, const T *cimp_keq_in, const T *cimp_phi_in, const T *cmap_surf_in, const T *cmap_dphi_in, const T *cmap_dpsi_in, const T *cmap_dphi_dpsi_in, const T *cmap_patches_in, const int *bond_asgn_atoms_in, const int *bond_asgn_index_in, const int *bond_asgn_terms_in, const int *bond_asgn_bounds_in, const int *angl_asgn_atoms_in, const int *angl_asgn_index_in, const int *angl_asgn_terms_in, const int *angl_asgn_bounds_in, const int *dihe_asgn_atoms_in, const int *dihe_asgn_index_in, const int *dihe_asgn_terms_in, const int *dihe_asgn_bounds_in, const int *ubrd_asgn_atoms_in, const int *ubrd_asgn_index_in, const int *ubrd_asgn_terms_in, const int *ubrd_asgn_bounds_in, const int *cimp_asgn_atoms_in, const int *cimp_asgn_index_in, const int *cimp_asgn_terms_in, const int *cimp_asgn_bounds_in, const int *cmap_asgn_atoms_in, const int *cmap_asgn_index_in, const int *cmap_asgn_terms_in, const int *cmap_asgn_bounds_in) | |
The constructor takes a Fortran77-worthy list of arguments. It would not be any cleaner to make this object a friend of or nested struct within AtomGraph, due to its templated nature. | |
ValenceKit (const ValenceKit &original)=default | |
Take the default copy and move constructors. The move assignment operator will get implicitly deleted as this is just a collection of constants. | |
ValenceKit (ValenceKit &&other)=default | |
Public Attributes | |
const int | natom |
const int | nbond |
Number of bonds in the system. | |
const int | nangl |
Number of bond angles in the system. | |
const int | ndihe |
const int | nbond_param |
Number of unique bond parameters. | |
const int | nangl_param |
Number of unqiue bond angle parameters. | |
const int | ndihe_param |
Number of unique cosine dihedral parameters. | |
const int | ninfr14 |
The number of inferred 1:4 exclusions. | |
const int | nattn14_param |
Number of 1:4 attenuation factor pairs (vdW, electrostatic) | |
const int | nubrd |
Number of Urey-Bradley 1:3 harmonic angle interactions. | |
const int | ncimp |
Number of CHARMM harmonic improper interactions. | |
const int | ncmap |
Number of CHARMM CMAP surface spline-based interactions. | |
const int | nubrd_param |
Number of unique Urey-Bradley parameters. | |
const int | ncimp_param |
Number of unique CHARMM harmonic improper parameters. | |
const int | ncmap_surf |
Number of unique CMAP surfaces. | |
const T * | bond_keq |
Equilibrium stiffness constants for all unique bonds. | |
const T * | bond_leq |
Equilibrium lengths for all unique bonds. | |
const T * | angl_keq |
Equilibrium stiffness constants for all unique bond angles. | |
const T * | angl_theta |
Equilibrium lengths for all unique bond angles. | |
const T * | dihe_amp |
Amplitudes for all unique cosine-based dihedrals. | |
const T * | dihe_freq |
Periodicities for cosine-based dihedral / torsion terms. | |
const T * | dihe_phi |
Phase angles for cosine-based dihedral / torsion terms. | |
const T * | attn14_elec |
const T * | attn14_vdw |
const int * | bond_i_atoms |
Array of first atoms in each bond. | |
const int * | bond_j_atoms |
Array of second atoms in each bond. | |
const int * | bond_param_idx |
Parameter indices for each bond in the system. | |
const char4 * | bond_modifiers |
Modifying details of each bond stretching term. | |
const uint * | bond_relevance |
const int * | angl_i_atoms |
Array of first atoms in each bond angle. | |
const int * | angl_j_atoms |
Array of second atoms in each bond angle. | |
const int * | angl_k_atoms |
Array of third atoms in each bond angle. | |
const int * | angl_param_idx |
Parameter indices for each bond angle in the system. | |
const char4 * | angl_modifiers |
Modifying details of each angle bending term. | |
const uint * | angl_relevance |
const int * | dihe_i_atoms |
Array of first atoms in each cosine-based dihedral. | |
const int * | dihe_j_atoms |
Array of second atoms in each cosine-based dihedral. | |
const int * | dihe_k_atoms |
const int * | dihe_l_atoms |
Array of fourth atoms in each cosine-based dihedral. | |
const int * | dihe_param_idx |
Parameter indices for each cosine-based dihedral in the system. | |
const int * | dihe14_param_idx |
const char4 * | dihe_modifiers |
const int * | infr14_i_atoms |
const int * | infr14_l_atoms |
Inferred 1:4 interactions' J atoms. | |
const int * | infr14_param_idx |
const int * | ubrd_i_atoms |
First atoms in each Urey-Bradley harmonic angle interaction. | |
const int * | ubrd_k_atoms |
Second atoms in each Urey-Bradley harmonic angle interaction. | |
const int * | ubrd_param_idx |
Parameter indices for each Urey-Bradley interaction. | |
const int * | cimp_i_atoms |
First atoms in each CHARMM harmonic improper dihedral. | |
const int * | cimp_j_atoms |
Second atoms in each CHARMM harmonic improper dihedral. | |
const int * | cimp_k_atoms |
Third (center) atoms in each CHARMM harmonic improper dihedral. | |
const int * | cimp_l_atoms |
Fourth atoms in each CHARMM harmonic improper dihedral. | |
const int * | cimp_param_idx |
Parameter indices for each CHARMM harmonic improper dihedral. | |
const int * | cmap_i_atoms |
First atoms of the first dihedral in each CHARMM CMAP term. | |
const int * | cmap_j_atoms |
const int * | cmap_k_atoms |
const int * | cmap_l_atoms |
const int * | cmap_m_atoms |
Fourth atoms of the second dihedral in each CHARMM CMAP term. | |
const int * | cmap_dim |
Dimension of each CMAP surface (each surface is a square grid) | |
const int * | cmap_surf_bounds |
Bounds array for individual CMAP surfaces. | |
const int * | cmap_patch_bounds |
const int * | cmap_surf_idx |
Parameter (surface) indices for CHARMM CMAP interactions. | |
const T * | ubrd_keq |
Array of unique Urey-Bradley 1:3 harmonic stiffnesses. | |
const T * | ubrd_leq |
Array of unique Urey-Bradley 1:3 harmonic equilibrium lengths. | |
const T * | cimp_keq |
Array of unique CHARMM harmonic improper dihedral stiffnesses. | |
const T * | cimp_phi |
Array of unique CHARMM harmonic improper dihedral phase angles. | |
const T * | cmap_surf |
Array of CHARMM CMAP surface values. | |
const T * | cmap_dphi |
Array of CHARMM CMAP surface "phi" derivative values. | |
const T * | cmap_dpsi |
Array of CHARMM CMAP surface "psi" derivative values. | |
const T * | cmap_dphi_dpsi |
Array of CHARMM CMAP surface cross derivative values. | |
const T * | cmap_patches |
Array of interlaced CHARMM CMAP patch matrices. | |
const int * | bond_asgn_atoms |
Other atoms for bond terms (one atom per term) | |
const int * | bond_asgn_index |
Parameter indices for bond terms. | |
const int * | bond_asgn_terms |
const int * | bond_asgn_bounds |
Bounds for each atom's assigned bond list. | |
const int * | angl_asgn_atoms |
Other atoms for bond angle terms (two atoms per term) | |
const int * | angl_asgn_index |
Parameter indices for bond angle terms. | |
const int * | angl_asgn_terms |
const int * | angl_asgn_bounds |
Bounds for each atom's assigned bond angle list. | |
const int * | dihe_asgn_atoms |
Other atoms for dihedral terms (three atoms per term) | |
const int * | dihe_asgn_index |
Parameter indices for dihedral terms. | |
const int * | dihe_asgn_terms |
const int * | dihe_asgn_bounds |
Bounds for each atom's assigned dihedral list. | |
const int * | ubrd_asgn_atoms |
Other atoms for Urey-Bradley terms (one atom per term) | |
const int * | ubrd_asgn_index |
Parameter indices for Urey-Bradley terms. | |
const int * | ubrd_asgn_terms |
Urey-Bradley term indices assigned to each atom. | |
const int * | ubrd_asgn_bounds |
Bounds for each atom's assigned Urey-Bradley list. | |
const int * | cimp_asgn_atoms |
Other atoms for CHARMM impropers (three atoms per term) | |
const int * | cimp_asgn_index |
Parameter indices for CHARMM improper dihedral terms. | |
const int * | cimp_asgn_terms |
CHARMM improper term indices assigned to each atom. | |
const int * | cimp_asgn_bounds |
Bounds for each atom's assigned CHARMM improper list. | |
const int * | cmap_asgn_atoms |
Other atoms for CMAP terms (four atoms per term) | |
const int * | cmap_asgn_index |
Parameter indices for CMAP terms. | |
const int * | cmap_asgn_terms |
CMAP term indices assigned to each atom. | |
const int * | cmap_asgn_bounds |
Bounds for each atom's assigned CMAP list. | |
Information need for bonded calculations. Templating is used to serve either of two levels of precision: single (which, on almost any machine and HPC accelerator card, will be fp32) or double (fp64). Parameters for bonds, the highest frequency terms in any system by a factor of roughly 5, will always be stored in double precision, as will parameters for bond angles, the next highest frequency group of terms. The precision of other parameters varies according to the template parameter.
|
explicit |
The constructor takes a Fortran77-worthy list of arguments. It would not be any cleaner to make this object a friend of or nested struct within AtomGraph, due to its templated nature.
See the descriptions of eponymous member variables (minus _in) within the ValenceKit object for descriptions of each parameter.
const int* stormm::topology::ValenceKit< T >::angl_asgn_terms |
Angle term indices assigned to each atom (i.e. the 9th atom is assigned angle terms with indices 46, 48, and 50)
const uint* stormm::topology::ValenceKit< T >::angl_relevance |
Array detailing the relevance of each bond angle interaction, whether it should be evaluated given directives on the rigidity of bonds to hydrogen atoms or virtual sites
const T* stormm::topology::ValenceKit< T >::attn14_elec |
Electrostatic 1:4 attenuation parameters, including a zero for no 1:4 interaction. The arrays dihe14_param_idx and infr14_param_idx both index into this array.
const T* stormm::topology::ValenceKit< T >::attn14_vdw |
van-der Waals 1:4 attenuation parameters, including a zero for no 1:4 interaction. The arrays dihe14_param_idx and infr14_param_idx both index into this array.
const int* stormm::topology::ValenceKit< T >::bond_asgn_terms |
Bond term indices assigned to each atom (i.e. the 5th atom is assigned bond terms with indices 18, 20, and 21)
const uint* stormm::topology::ValenceKit< T >::bond_relevance |
Array detailing the relevance of each bonded interaction, whether it should be evaluated given directives on the rigidity of bonds to hydrogen atoms or virtual sites
const int* stormm::topology::ValenceKit< T >::cmap_j_atoms |
Second atoms of the first dihedral, first atoms of the second dihedral, in each CHARMM CMAP term
const int* stormm::topology::ValenceKit< T >::cmap_k_atoms |
Third atoms of the first dihedral, second atoms of the second dihedral, in each CHARMM CMAP term
const int* stormm::topology::ValenceKit< T >::cmap_l_atoms |
Fourth atoms of the first dihedral, third atoms of the second dihedral, in each CHARMM CMAP term
const int* stormm::topology::ValenceKit< T >::cmap_patch_bounds |
Bounds array for individual CMAP surfaces composed of patch matrices (this array is basically 16 times cmap_surf_bounds)
const int* stormm::topology::ValenceKit< T >::dihe14_param_idx |
Parameter indices for 1:4 attenuated interactions handled by each dihedral through their I and L atoms
const int* stormm::topology::ValenceKit< T >::dihe_asgn_terms |
Dihedral term indices assigned to each atom (i.e. the 4th atom is assigned dihedral terms with indices 27, 28, and 29)
const int* stormm::topology::ValenceKit< T >::dihe_k_atoms |
Array of third atoms in each cosine-based dihedral (these are the center atoms if the dihedral describes a cosine-based improper)
const char4* stormm::topology::ValenceKit< T >::dihe_modifiers |
Modifying details of each dihedral term, including whether it is a proper or improper dihedral with or without a 1:4 attenuated interaction
const int* stormm::topology::ValenceKit< T >::infr14_i_atoms |
Inferred 1:4 interactions' I atoms. These cannot be handled directly by any dihedrals and so must be treated as separate terms.
const int* stormm::topology::ValenceKit< T >::infr14_param_idx |
Indices into the arrays of 1:4 attenuation parameters for the set of 1:4 interactions with inferred parameters
const int stormm::topology::ValenceKit< T >::natom |
The number of atoms in the system (needed for the overall length of bond_asgn_abounds and similar arrays towards the end, if not just for general convenience)
const int stormm::topology::ValenceKit< T >::ndihe |
Number of dihedrals or torsions (includes cosine-based improper torsions) in the system