ComputeGridForce Class Reference

#include <ComputeGridForce.h>

Inheritance diagram for ComputeGridForce:
ComputeHomePatch Compute

List of all members.

Public Member Functions

 ComputeGridForce (ComputeID c, PatchID pid)
virtual ~ComputeGridForce ()
void doForce (FullAtom *p, Results *r)

Public Attributes

SubmitReductionreduction

Protected Member Functions

template<class T >
void do_calc (T *grid, int gridnum, FullAtom *p, int numAtoms, Molecule *mol, Force *forces, BigReal &energy, Force &extForce, Tensor &extVirial)

Detailed Description

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

Definition at line 17 of file ComputeGridForce.h.


Constructor & Destructor Documentation

ComputeGridForce::ComputeGridForce ( ComputeID  c,
PatchID  pid 
)

Definition at line 23 of file ComputeGridForce.C.

References ReductionMgr::Object(), reduction, REDUCTIONS_BASIC, and ReductionMgr::willSubmit().

00024     : ComputeHomePatch(c,pid)
00025 {
00026 
00027     reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC);
00028 
00029 }

ComputeGridForce::~ComputeGridForce (  )  [virtual]

Definition at line 33 of file ComputeGridForce.C.

References reduction.

00034 {
00035     delete reduction;
00036 }


Member Function Documentation

template<class T >
void ComputeGridForce::do_calc ( T *  grid,
int  gridnum,
FullAtom p,
int  numAtoms,
Molecule mol,
Force forces,
BigReal energy,
Force extForce,
Tensor extVirial 
) [inline, protected]

Definition at line 39 of file ComputeGridForce.C.

References DebugM, endi(), Patch::flags, Molecule::get_gridfrc_params(), ComputeHomePatch::homePatch, CompAtomExt::id, iout, Molecule::is_atom_gridforced(), iWARN(), Transform::j, Transform::k, Patch::lattice, outer(), CompAtom::position, Lattice::reverse_transform(), Flags::step, FullAtom::transform, Vector::x, Vector::y, and Vector::z.

Referenced by doForce().

00040 {
00041     Real scale;                 // Scaling factor
00042     Charge charge;              // Charge
00043     Vector dV;
00044     float V;
00045     
00046     Vector gfScale = grid->get_scale();
00047     DebugM(3, "doCalc()\n" << endi);
00048 
00049     //  Loop through and check each atom
00050     for (int i = 0; i < numAtoms; i++) {
00051         if (mol->is_atom_gridforced(p[i].id, gridnum))
00052         {
00053             DebugM(1, "Atom " << p[i].id << " is gridforced\n" << endi);
00054             
00055             mol->get_gridfrc_params(scale, charge, p[i].id, gridnum);
00056             
00057             // Wrap coordinates using grid center
00058             Position pos = grid->wrap_position(p[i].position, homePatch->lattice);
00059             DebugM(1, "pos = " << pos << "\n" << endi);
00060             
00061             // Here's where the action happens
00062             int err = grid->compute_VdV(pos, V, dV);
00063             
00064             if (err) {
00065                 DebugM(2, "V = 0\n" << endi);
00066                 DebugM(2, "dV = 0 0 0\n" << endi);
00067                 continue;  // This means the current atom is outside the potential
00068             }
00069             
00070             //Force force = scale * Tensor::diagonal(gfScale) * (-charge * dV);
00071             Force force = -charge * scale * Vector(gfScale.x * dV.x, gfScale.y * dV.y, gfScale.z * dV.z);
00072             
00073 #ifdef DEBUGM
00074             DebugM(2, "scale = " << scale << " gfScale = " << gfScale << " charge = " << charge << "\n" << endi);
00075             
00076             DebugM(2, "V = " << V << "\n" << endi);
00077             DebugM(2, "dV = " << dV << "\n" << endi);
00078             DebugM(2, "grid = " << gridnum << " force = " << force << " pos = " << pos << " V = " << V << " dV = " << dV << " step = " << homePatch->flags.step << " index = " << p[i].id << "\n" << endi);
00079             
00080             DebugM(1, "transform = " << (int)p[i].transform.i << " "
00081                    << (int)p[i].transform.j << " " << (int)p[i].transform.k << "\n" << endi);
00082             
00083             if (V != V) {
00084                 iout << iWARN << "V is NaN!\natomid = " << p[i].id << " loc = " << p[i].position << " V = " << V << "\n" << endi;
00085             }
00086 #endif
00087             
00088             forces[i] += force;
00089             extForce += force;
00090             Position vpos = homePatch->lattice.reverse_transform(p[i].position, p[i].transform);
00091             
00092             //energy -= force * (vpos - homePatch->lattice.origin());
00093             if (gfScale.x == gfScale.y && gfScale.x == gfScale.z)
00094             {
00095                 // only makes sense when scaling is isotropic
00096                 energy += scale * gfScale.x * (charge * V);
00097                 
00098                 // add something when we're off the grid? I'm thinking no
00099             }
00100             extVirial += outer(force,vpos);
00101         }
00102     }
00103     DebugM(3, "doCalc() done\n" << endi);
00104 }

void ComputeGridForce::doForce ( FullAtom p,
Results r 
) [virtual]

Implements ComputeHomePatch.

Definition at line 106 of file ComputeGridForce.C.

References ADD_TENSOR_OBJECT, ADD_VECTOR_OBJECT, SimParameters::berendsenPressureOn, DebugM, do_calc(), endi(), Results::f, GridforceGrid::fits_lattice(), Patch::flags, forces, GridforceGrid::get_center(), GridforceGrid::get_checksize(), GridforceGrid::get_e(), GridforceGrid::get_grid_type(), Molecule::get_gridfrc_grid(), GridforceGrid::get_scale(), Patch::getNumAtoms(), GF_OVERLAPCHECK_FREQ, GridforceGrid::GridforceGridTypeFull, GridforceGrid::GridforceGridTypeLite, ComputeHomePatch::homePatch, SubmitReduction::item(), SimParameters::langevinPistonOn, Patch::lattice, Node::molecule, NAMD_bug(), NAMD_die(), Results::normal, ComputeHomePatch::numAtoms, Molecule::numGridforceGrids, Node::Object(), reduction, REDUCTION_MISC_ENERGY, Node::simParameters, simParams, Flags::step, SubmitReduction::submit(), Vector::x, Vector::y, and Vector::z.

00107 {
00108     SimParameters *simParams = Node::Object()->simParameters;
00109     Molecule *mol = Node::Object()->molecule;
00110     
00111     Force *forces = r->f[Results::normal];
00112     BigReal energy = 0;
00113     Force extForce = 0.;
00114     Tensor extVirial;
00115     
00116     int numAtoms = homePatch->getNumAtoms();
00117 
00118     if ( mol->numGridforceGrids < 1 ) NAMD_bug("No grids loaded in ComputeGridForce::doForce()");
00119 
00120     DebugM(3, "doForce()\n" << endi);
00121     
00122     for (int gridnum = 0; gridnum < mol->numGridforceGrids; gridnum++) {
00123         GridforceGrid *grid = mol->get_gridfrc_grid(gridnum);
00124 
00125         const Vector gfScale = grid->get_scale();
00126         if ((gfScale.x == 0.0) && (gfScale.y == 0.0) && (gfScale.z == 0.0)) {
00127           DebugM(3, "Skipping grid index " << gridnum << "\n" << endi);
00128           continue;
00129         }
00130         
00131         if (homePatch->flags.step % GF_OVERLAPCHECK_FREQ == 0) {
00132             // only check on node 0 and every GF_OVERLAPCHECK_FREQ steps
00133           if (simParams->langevinPistonOn || simParams->berendsenPressureOn) {
00134                 // check for grid overlap if pressure control is on
00135                 // not needed without pressure control, since the check is also performed on startup
00136       if (!grid->fits_lattice(homePatch->lattice)) {
00137         char errmsg[512];
00138         if (grid->get_checksize()) {
00139           sprintf(errmsg, "Warning: Periodic cell basis too small for Gridforce grid %d.  Set gridforcechecksize off in configuration file to ignore.\n", gridnum);
00140           NAMD_die(errmsg);      
00141         }
00142       }
00143          }
00144         }
00145         
00146         if (homePatch->flags.step % 100 == 1) {
00147             Position center = grid->get_center();
00148             DebugM(3, "center = " << center << "\n" << endi);
00149             DebugM(3, "e = " << grid->get_e() << "\n" << endi);
00150         }
00151         
00152         if (grid->get_grid_type() == GridforceGrid::GridforceGridTypeFull) {
00153             GridforceFullMainGrid *g = (GridforceFullMainGrid *)grid;
00154             do_calc(g, gridnum, p, numAtoms, mol, forces, energy, extForce, extVirial);
00155         } else if (grid->get_grid_type() == GridforceGrid::GridforceGridTypeLite) {
00156             GridforceLiteGrid *g = (GridforceLiteGrid *)grid;
00157             do_calc(g, gridnum, p, numAtoms, mol, forces, energy, extForce, extVirial);
00158         }
00159     }
00160     reduction->item(REDUCTION_MISC_ENERGY) += energy;
00161     ADD_VECTOR_OBJECT(reduction,REDUCTION_EXT_FORCE_NORMAL,extForce);
00162     ADD_TENSOR_OBJECT(reduction,REDUCTION_VIRIAL_NORMAL,extVirial);
00163     reduction->submit();
00164     DebugM(3, "doForce() done\n" << endi);
00165 }


Member Data Documentation

Definition at line 28 of file ComputeGridForce.h.

Referenced by ComputeGridForce(), doForce(), and ~ComputeGridForce().


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

Generated on 4 Jun 2020 for NAMD by  doxygen 1.6.1