ComputeConsForce.C

Go to the documentation of this file.
00001 
00002 #include "InfoStream.h"
00003 #include "SimParameters.h"
00004 #include "ComputeConsForce.h"
00005 #include "Molecule.h"
00006 #include "Node.h"
00007 #include "HomePatch.h"
00008 
00010 // compute constant force
00011 ComputeConsForce::ComputeConsForce(ComputeID c, PatchID pid)
00012   : ComputeHomePatch(c,pid)
00013 {
00014   reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC);
00015 }
00016 
00017 ComputeConsForce::~ComputeConsForce()
00018 {
00019   delete reduction;
00020 }
00021 
00022 void ComputeConsForce::doForce(FullAtom* p, Results* r)
00023 { int localID,forceID;
00024   Molecule *molecule = Node::Object()->molecule;
00025   BigReal scaling = Node::Object()->simParameters->consForceScaling;
00026   int32 *index = molecule->consForceIndexes;  // Indexes into the force array
00027   Vector *cf = molecule->consForce;  // Force array
00028   Vector *forces = r->f[Results::normal];
00029   Force extForce = 0.;
00030   Tensor extVirial;
00031 
00032   for (localID=0; localID<numAtoms; ++localID) {
00033     // When the index is -1, it means there's no constant force on this atom
00034     if ((forceID=index[p[localID].id]) != -1) {
00035       Vector sf = scaling * cf[forceID];
00036       forces[localID] += sf;
00037       extForce += sf;
00038       Position vpos = homePatch->lattice.reverse_transform(
00039                 p[localID].position, p[localID].transform );
00040       extVirial += outer(sf,vpos);
00041     }
00042   }
00043 
00044   ADD_VECTOR_OBJECT(reduction,REDUCTION_EXT_FORCE_NORMAL,extForce);
00045   ADD_TENSOR_OBJECT(reduction,REDUCTION_VIRIAL_NORMAL,extVirial);
00046   reduction->submit();
00047 }
00048 
00050 // compute "constant" torque
00051 ComputeConsTorque::ComputeConsTorque(ComputeID c, PatchID pid)
00052   : ComputeHomePatch(c,pid)
00053 {
00054   reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC);
00055 }
00056 
00057 ComputeConsTorque::~ComputeConsTorque()
00058 {
00059   delete reduction;
00060 }
00061 
00062 void ComputeConsTorque::doForce(FullAtom* p, Results* r)
00063 { int localID,torqueID;
00064   Molecule *molecule = Node::Object()->molecule;
00065   SimParameters *simParams = Node::Object()->simParameters;
00066 
00067   int32 *index = molecule->consTorqueIndexes;  // Indexes into the torque array
00068   Vector *forces = r->f[Results::normal];
00069   Force extForce = 0.;
00070   Tensor extVirial;
00071   const BigReal consTorqueGlobVal = simParams->consTorqueGlobVal;
00072   BigReal consTorqueVal;
00073   Vector consTorqueAxis, consTorquePivot;  
00074   Vector atomRadius;
00075   Vector torque;
00076 
00077   for (localID=0; localID<numAtoms; ++localID) {
00078     // When the index is -1, it means there's no constant torque on this atom
00079     if ((torqueID=index[p[localID].id]) != -1) {
00080       // compute the torqueing force and add it to the net force
00081       molecule->get_constorque_params(consTorqueVal, consTorqueAxis, consTorquePivot, p[localID].id);
00082       consTorqueAxis /= consTorqueAxis.length();
00083       atomRadius = p[localID].position - consTorquePivot;
00084       torque = cross(consTorqueAxis, atomRadius) * consTorqueVal * consTorqueGlobVal;
00085       forces[localID] += torque;
00086       extForce += torque;
00087       Position vpos = homePatch->lattice.reverse_transform(
00088                 p[localID].position, p[localID].transform );
00089       extVirial += outer(torque,vpos);
00090     }
00091   }
00092 
00093   ADD_VECTOR_OBJECT(reduction,REDUCTION_EXT_FORCE_NORMAL,extForce);
00094   ADD_TENSOR_OBJECT(reduction,REDUCTION_VIRIAL_NORMAL,extVirial);
00095   reduction->submit();
00096 }

Generated on Mon Nov 20 01:17:10 2017 for NAMD by  doxygen 1.4.7