TholeElem Class Reference

#include <ComputeThole.h>

List of all members.

Public Types

enum  { size = 4 }
enum  {
  tholeEnergyIndex, tholeEnergyIndex_f, tholeEnergyIndex_ti_1, tholeEnergyIndex_ti_2,
  TENSOR = (virialIndex), reductionDataSize
}
enum  { reductionChecksumLabel = REDUCTION_THOLE_CHECKSUM }

Public Member Functions

int hash () const
 TholeElem ()
 TholeElem (AtomID atom0, const TupleSignature *sig, const TholeValue *v)
 TholeElem (const Thole *a, const TholeValue *v)
 TholeElem (AtomID atom0, AtomID atom1, AtomID atom2, AtomID atom3)
 ~TholeElem ()
int operator== (const TholeElem &a) const
int operator< (const TholeElem &a) const

Static Public Member Functions

static void computeForce (TholeElem *, int, BigReal *, BigReal *)
static void getMoleculePointers (Molecule *, int *, int32 ***, Thole **)
static void getParameterPointers (Parameters *, const TholeValue **)
static void getTupleInfo (AtomSignature *sig, int *count, TupleSignature **t)
static void submitReductionData (BigReal *, SubmitReduction *)

Public Attributes

AtomID atomID [size]
int localIndex [size]
TuplePatchElemp [size]
Real scale
const TholeValuevalue

Static Public Attributes

static int pressureProfileSlabs = 0
static int pressureProfileAtomTypes = 1
static BigReal pressureProfileThickness = 0
static BigReal pressureProfileMin = 0

Detailed Description

Definition at line 18 of file ComputeThole.h.


Member Enumeration Documentation

anonymous enum
Enumerator:
size 

Definition at line 21 of file ComputeThole.h.

00021 { size = 4 };

anonymous enum
Enumerator:
tholeEnergyIndex 
tholeEnergyIndex_f 
tholeEnergyIndex_ti_1 
tholeEnergyIndex_ti_2 
TENSOR 
reductionDataSize 

Definition at line 49 of file ComputeThole.h.

anonymous enum
Enumerator:
reductionChecksumLabel 

Definition at line 51 of file ComputeThole.h.


Constructor & Destructor Documentation

TholeElem::TholeElem (  )  [inline]

Copyright (c) 1995, 1996, 1997, 1998, 1999, 2000 by The Board of Trustees of the University of Illinois. All rights reserved.

Definition at line 12 of file ComputeThole.inl.

00012 { ; }

TholeElem::TholeElem ( AtomID  atom0,
const TupleSignature sig,
const TholeValue v 
) [inline]

Definition at line 14 of file ComputeThole.inl.

References NAMD_die().

00014                                                                                        {
00015     NAMD_die("Can't use Thole with memory optimized version of NAMD.");
00016     // atomID[0] = atom0;
00017     // atomID[1] = atom0 + sig->offset[0];
00018     // atomID[2] = atom0 + sig->offset[1];
00019     // atomID[3] = atom0 + sig->offset[2];
00020     // value = &v[sig->tupleParamType];
00021 }

TholeElem::TholeElem ( const Thole a,
const TholeValue v 
) [inline]

Definition at line 23 of file ComputeThole.inl.

References thole::atom1, thole::atom2, thole::atom3, thole::atom4, atomID, and value.

00024   {
00025     atomID[0] = a->atom1;
00026     atomID[1] = a->atom2;
00027     atomID[2] = a->atom3;
00028     atomID[3] = a->atom4;
00029     value = a;  // expect v to be NULL
00030   }

TholeElem::TholeElem ( AtomID  atom0,
AtomID  atom1,
AtomID  atom2,
AtomID  atom3 
) [inline]

Definition at line 32 of file ComputeThole.inl.

References atomID.

00034   {
00035     // atoms arranged:  HEAVY DRUDE HEAVY DRUDE
00036     if (atom0 > atom2) {  // Swap heavy atoms so lowest is first!
00037       AtomID tmp = atom2; atom2 = atom0; atom0 = tmp; 
00038       tmp = atom1; atom1 = atom3; atom3 = tmp;
00039     }
00040     atomID[0] = atom0;
00041     atomID[1] = atom1;
00042     atomID[2] = atom2;
00043     atomID[3] = atom3;
00044   }

TholeElem::~TholeElem (  )  [inline]

Definition at line 58 of file ComputeThole.h.

00058 {};


Member Function Documentation

void TholeElem::computeForce ( TholeElem tuples,
int  ntuple,
BigReal reduction,
BigReal pressureProfileData 
) [static]

Definition at line 41 of file ComputeThole.C.

References thole::aa, SimParameters::alchDecouple, SimParameters::alchOn, atomID, DebugM, Lattice::delta(), TuplePatchElem::f, Patch::flags, Molecule::get_fep_type(), SimParameters::getCurrentLambda(), SimParameters::getCurrentLambda2(), SimParameters::getElecLambda(), Patch::lattice, localIndex, Node::molecule, Node::Object(), order, TuplePatchElem::p, p, CompAtom::partition, CompAtom::position, pp_clamp(), pp_reduction(), pressureProfileAtomTypes, pressureProfileMin, pressureProfileSlabs, pressureProfileThickness, thole::qq, Vector::rlength(), scale, Node::simParameters, simParams, size, Flags::step, tholeEnergyIndex, tholeEnergyIndex_f, tholeEnergyIndex_ti_1, tholeEnergyIndex_ti_2, value, TuplePatchElem::x, and Vector::z.

00043 {
00044  const Lattice & lattice = tuples[0].p[0]->p->lattice;
00045 
00046  //fepb BKR
00047  SimParameters *const simParams = Node::Object()->simParameters;
00048  const int step = tuples[0].p[0]->p->flags.step;
00049  const BigReal alchLambda = simParams->getCurrentLambda(step);
00050  const BigReal alchLambda2 = simParams->getCurrentLambda2(step);
00051  const BigReal elec_lambda_1 = simParams->getElecLambda(alchLambda);
00052  const BigReal elec_lambda_2 = simParams->getElecLambda(1-alchLambda);
00053  const BigReal elec_lambda_12 = simParams->getElecLambda(alchLambda2);
00054  const BigReal elec_lambda_22 = simParams->getElecLambda(1-alchLambda2);
00055  Molecule *const mol = Node::Object()->molecule;
00056  //fepe
00057 
00058 
00059  for ( int ituple=0; ituple<ntuple; ++ituple ) {
00060   const TholeElem &tup = tuples[ituple];
00061   enum { size = 4 };
00062   const AtomID (&atomID)[size](tup.atomID);
00063   const int    (&localIndex)[size](tup.localIndex);
00064   TuplePatchElem * const(&p)[size](tup.p);
00065   const Real (&scale)(tup.scale);
00066   const TholeValue * const(&value)(tup.value);
00067 
00068   DebugM(3, "::computeForce() localIndex = " << localIndex[0] << " "
00069                << localIndex[1] << " " << localIndex[2] << " "
00070                << localIndex[3] << std::endl);
00071 
00072 #ifdef CALCULATE_THOLE_CORRECTION
00073   const BigReal aa = value->aa;
00074   const BigReal qq = value->qq;
00075 
00076   //  Calculate the vectors between atoms
00077   const Position & rai = p[0]->x[localIndex[0]].position;  // atom i
00078   const Position & rdi = p[1]->x[localIndex[1]].position;  // atom i's Drude
00079   const Position & raj = p[2]->x[localIndex[2]].position;  // atom j
00080   const Position & rdj = p[3]->x[localIndex[3]].position;  // atom j's Drude
00081 
00082   // r_ij = r_i - r_j
00083   Vector raa = lattice.delta(rai,raj);  // shortest vector image:  rai - raj
00084   Vector rad = lattice.delta(rai,rdj);  // shortest vector image:  rai - rdj
00085   Vector rda = lattice.delta(rdi,raj);  // shortest vector image:  rdi - raj
00086   Vector rdd = lattice.delta(rdi,rdj);  // shortest vector image:  rdi - rdj
00087 
00088   // 1/r, r = |r_ij|
00089   BigReal raa_invlen = raa.rlength();  // reciprocal of length
00090   BigReal rad_invlen = rad.rlength();
00091   BigReal rda_invlen = rda.rlength();
00092   BigReal rdd_invlen = rdd.rlength();
00093 
00094   // ar
00095   BigReal auaa = aa / raa_invlen;
00096   BigReal auad = aa / rad_invlen;
00097   BigReal auda = aa / rda_invlen;
00098   BigReal audd = aa / rdd_invlen;
00099 
00100   // exp(-ar)
00101   BigReal expauaa = exp(-auaa);
00102   BigReal expauad = exp(-auad);
00103   BigReal expauda = exp(-auda);
00104   BigReal expaudd = exp(-audd);
00105 
00106   // (1 + ar/2)
00107   BigReal polyauaa = 1 + 0.5*auaa;
00108   BigReal polyauad = 1 + 0.5*auad;
00109   BigReal polyauda = 1 + 0.5*auda;
00110   BigReal polyaudd = 1 + 0.5*audd;
00111 
00112   // U(r) = qq/r (1 - (1 + ar/2) exp(-ar))
00113   BigReal ethole = 0;
00114   ethole += qq * raa_invlen * (1 - polyauaa * expauaa);
00115   ethole += -qq * rad_invlen * (1 - polyauad * expauad);
00116   ethole += -qq * rda_invlen * (1 - polyauda * expauda);
00117   ethole += qq * rdd_invlen * (1 - polyaudd * expaudd);
00118 
00119   polyauaa = 1 + auaa*polyauaa;
00120   polyauad = 1 + auad*polyauad;
00121   polyauda = 1 + auda*polyauda;
00122   polyaudd = 1 + audd*polyaudd;
00123 
00124   BigReal raa_invlen3 = raa_invlen * raa_invlen * raa_invlen;
00125   BigReal rad_invlen3 = rad_invlen * rad_invlen * rad_invlen;
00126   BigReal rda_invlen3 = rda_invlen * rda_invlen * rda_invlen;
00127   BigReal rdd_invlen3 = rdd_invlen * rdd_invlen * rdd_invlen;
00128 
00129   // df = (1/r) (dU/dr)
00130   BigReal dfaa = qq * raa_invlen3 * (polyauaa*expauaa - 1);
00131   BigReal dfad = -qq * rad_invlen3 * (polyauad*expauad - 1);
00132   BigReal dfda = -qq * rda_invlen3 * (polyauda*expauda - 1);
00133   BigReal dfdd = qq * rdd_invlen3 * (polyaudd*expaudd - 1);
00134 
00135   Vector faa = -dfaa * raa;
00136   Vector fad = -dfad * rad;
00137   Vector fda = -dfda * rda;
00138   Vector fdd = -dfdd * rdd;
00139 
00140   //fepb - BKR scaling of alchemical drude terms
00141   //       NB: TI derivative is the _unscaled_ energy.
00142   if ( simParams->alchOn ) {
00143     // THIS ASSUMES THAT AN ATOM AND ITS DRUDE PARTICLE ARE ALWAYS IN THE SAME
00144     // PARTITION!
00145     int typeSum = 0;
00146     typeSum += (mol->get_fep_type(atomID[0]) == 2 ? -1 :\
00147         mol->get_fep_type(atomID[0]));
00148     typeSum += (mol->get_fep_type(atomID[2]) == 2 ? -1 :\
00149         mol->get_fep_type(atomID[2]));
00150     int order = (!simParams->alchDecouple ? 3 : 2);
00151 
00152     if ( 0 < typeSum && typeSum < order ) {
00153       reduction[tholeEnergyIndex_ti_1] += ethole;
00154       reduction[tholeEnergyIndex_f] += elec_lambda_12*ethole;
00155       ethole *= elec_lambda_1;
00156       faa *= elec_lambda_1;
00157       fad *= elec_lambda_1;
00158       fda *= elec_lambda_1;
00159       fdd *= elec_lambda_1; 
00160     } else if ( 0 > typeSum && typeSum > -order ) {
00161       reduction[tholeEnergyIndex_ti_2] += ethole;
00162       reduction[tholeEnergyIndex_f] += elec_lambda_22*ethole;
00163       ethole *= elec_lambda_2;
00164       faa *= elec_lambda_2;
00165       fad *= elec_lambda_2;
00166       fda *= elec_lambda_2;
00167       fdd *= elec_lambda_2;
00168     }
00169   }
00170   //fepe
00171 
00172   p[0]->f[localIndex[0]] += faa + fad;
00173   p[1]->f[localIndex[1]] += fda + fdd;
00174   p[2]->f[localIndex[2]] -= faa + fda;
00175   p[3]->f[localIndex[3]] -= fad + fdd;
00176 
00177   DebugM(3, "::computeForce() -- ending with delta energy " << ethole
00178       << std::endl);
00179   reduction[tholeEnergyIndex] += ethole;
00180 
00181   reduction[virialIndex_XX] += faa.x * raa.x + fad.x * rad.x
00182     + fda.x * rda.x + fdd.x * rdd.x;
00183   reduction[virialIndex_XY] += faa.x * raa.y + fad.x * rad.y
00184     + fda.x * rda.y + fdd.x * rdd.y;
00185   reduction[virialIndex_XZ] += faa.x * raa.z + fad.x * rad.z
00186     + fda.x * rda.z + fdd.x * rdd.z;
00187   reduction[virialIndex_YX] += faa.y * raa.x + fad.y * rad.x
00188     + fda.y * rda.x + fdd.y * rdd.x;
00189   reduction[virialIndex_YY] += faa.y * raa.y + fad.y * rad.y
00190     + fda.y * rda.y + fdd.y * rdd.y;
00191   reduction[virialIndex_YZ] += faa.y * raa.z + fad.y * rad.z
00192     + fda.y * rda.z + fdd.y * rdd.z;
00193   reduction[virialIndex_ZX] += faa.z * raa.x + fad.z * rad.x
00194     + fda.z * rda.x + fdd.z * rdd.x;
00195   reduction[virialIndex_ZY] += faa.z * raa.y + fad.z * rad.y
00196     + fda.z * rda.y + fdd.z * rdd.y;
00197   reduction[virialIndex_ZZ] += faa.z * raa.z + fad.z * rad.z
00198     + fda.z * rda.z + fdd.z * rdd.z;
00199 
00200   if (pressureProfileData) {
00201     BigReal zai = p[0]->x[localIndex[0]].position.z;
00202     BigReal zdi = p[1]->x[localIndex[1]].position.z;
00203     BigReal zaj = p[2]->x[localIndex[2]].position.z;
00204     BigReal zdj = p[3]->x[localIndex[3]].position.z;
00205     int nai = (int)floor((zai-pressureProfileMin)/pressureProfileThickness);
00206     int ndi = (int)floor((zdi-pressureProfileMin)/pressureProfileThickness);
00207     int naj = (int)floor((zaj-pressureProfileMin)/pressureProfileThickness);
00208     int ndj = (int)floor((zdj-pressureProfileMin)/pressureProfileThickness);
00209     pp_clamp(nai, pressureProfileSlabs);
00210     pp_clamp(ndi, pressureProfileSlabs);
00211     pp_clamp(naj, pressureProfileSlabs);
00212     pp_clamp(ndj, pressureProfileSlabs);
00213     int pai = p[0]->x[localIndex[0]].partition;
00214     int pdi = p[1]->x[localIndex[1]].partition;
00215     int paj = p[2]->x[localIndex[2]].partition;
00216     int pdj = p[3]->x[localIndex[3]].partition;
00217     int pn = pressureProfileAtomTypes;
00218     pp_reduction(pressureProfileSlabs, nai, naj,
00219         pai, paj, pn, faa.x * raa.x, faa.y * raa.y, faa.z * raa.z,
00220         pressureProfileData);
00221     pp_reduction(pressureProfileSlabs, nai, ndj,
00222         pai, pdj, pn, fad.x * rad.x, fad.y * rad.y, fad.z * rad.z,
00223         pressureProfileData);
00224     pp_reduction(pressureProfileSlabs, ndi, naj,
00225         pdi, paj, pn, fda.x * rda.x, fda.y * rda.y, fda.z * rda.z,
00226         pressureProfileData);
00227     pp_reduction(pressureProfileSlabs, ndi, ndj,
00228         pdi, pdj, pn, fdd.x * rdd.x, fdd.y * rdd.y, fdd.z * rdd.z,
00229         pressureProfileData);
00230   }
00231 #endif
00232 
00233  }
00234 }

void TholeElem::getMoleculePointers ( Molecule mol,
int *  count,
int32 ***  byatom,
Thole **  structarray 
) [static]

Definition at line 26 of file ComputeThole.C.

References NAMD_die(), and Molecule::numTholes.

00027 {
00028 #ifdef MEM_OPT_VERSION
00029   NAMD_die("Should not be called in TholeElem::getMoleculePointers in memory optimized version!");
00030 #else
00031   *count = mol->numTholes;
00032   *byatom = mol->tholesByAtom;
00033   *structarray = mol->tholes;
00034 #endif
00035 }

void TholeElem::getParameterPointers ( Parameters p,
const TholeValue **  v 
) [static]

Definition at line 37 of file ComputeThole.C.

00037                                                                         {
00038   *v = NULL;  // parameters are stored in the structure
00039 }

static void TholeElem::getTupleInfo ( AtomSignature sig,
int *  count,
TupleSignature **  t 
) [inline, static]

Definition at line 30 of file ComputeThole.h.

References NAMD_die().

00030                                                                                  {
00031         NAMD_die("Can't use Thole with memory optimized version of NAMD.");
00032         // *count = sig->tholeCnt;
00033         // *t = sig->tholeSigs;
00034     }

int TholeElem::hash ( void   )  const [inline]

Definition at line 45 of file ComputeThole.h.

References atomID.

00045                    { 
00046     return 0x7FFFFFFF &((atomID[0]<<24) + (atomID[1]<<16) + (atomID[2]<<8) + atomID[3]);
00047   }

int TholeElem::operator< ( const TholeElem a  )  const [inline]

Definition at line 52 of file ComputeThole.inl.

References atomID.

00053   {
00054     return  (atomID[0] < a.atomID[0] ||
00055             (atomID[0] == a.atomID[0] &&
00056             (atomID[1] < a.atomID[1] ||
00057             (atomID[1] == a.atomID[1] &&
00058             (atomID[2] < a.atomID[2] ||
00059             (atomID[2] == a.atomID[2] &&
00060              atomID[3] < a.atomID[3] 
00061              ))))));
00062   }

int TholeElem::operator== ( const TholeElem a  )  const [inline]

Definition at line 46 of file ComputeThole.inl.

References atomID.

00047   {
00048     return (a.atomID[0] == atomID[0] && a.atomID[1] == atomID[1] &&
00049         a.atomID[2] == atomID[2] && a.atomID[3] == atomID[3]);
00050   }

void TholeElem::submitReductionData ( BigReal data,
SubmitReduction reduction 
) [static]

Member Data Documentation

Definition at line 22 of file ComputeThole.h.

Referenced by computeForce(), hash(), operator<(), operator==(), and TholeElem().

Definition at line 23 of file ComputeThole.h.

Referenced by computeForce().

Definition at line 24 of file ComputeThole.h.

Referenced by computeForce().

Definition at line 38 of file ComputeThole.h.

Referenced by computeForce().

Definition at line 40 of file ComputeThole.h.

Referenced by computeForce().

Definition at line 37 of file ComputeThole.h.

Referenced by computeForce().

Definition at line 39 of file ComputeThole.h.

Referenced by computeForce().

Definition at line 25 of file ComputeThole.h.

Referenced by computeForce().

Definition at line 43 of file ComputeThole.h.

Referenced by computeForce(), and TholeElem().


The documentation for this class was generated from the following files:

Generated on 12 Nov 2019 for NAMD by  doxygen 1.6.1