ComputeGlobal Class Reference

#include <ComputeGlobal.h>

Inheritance diagram for ComputeGlobal:

Public Member Functions
	ComputeGlobal (ComputeID, ComputeMgr *)
virtual	~ComputeGlobal ()
void	doWork ()
void	recvResults (ComputeGlobalResultsMsg *)
void	saveTotalForces (HomePatch *)
int	getForceSendActive () const
Public Member Functions inherited from ComputeHomePatches
	ComputeHomePatches (ComputeID c)
virtual	~ComputeHomePatches ()
virtual void	initialize ()
virtual void	atomUpdate ()
Flags *	getFlags (void)
Public Member Functions inherited from Compute
	Compute (ComputeID)
int	type ()
virtual	~Compute ()
void	setNumPatches (int n)
int	getNumPatches ()
virtual void	patchReady (PatchID, int doneMigration, int seq)
virtual int	noWork ()
virtual void	finishPatch (int)
int	sequence (void)
int	priority (void)
int	getGBISPhase (void)
virtual void	gbisP2PatchReady (PatchID, int seq)
virtual void	gbisP3PatchReady (PatchID, int seq)

Additional Inherited Members
Public Attributes inherited from Compute
const ComputeID	cid
LDObjHandle	ldObjHandle
LocalWorkMsg *const	localWorkMsg
Protected Member Functions inherited from Compute
void	enqueueWork ()
Protected Attributes inherited from ComputeHomePatches
int	useAvgPositions
int	hasPatchZero
ComputeHomePatchList	patchList
PatchMap *	patchMap
Protected Attributes inherited from Compute
int	computeType
int	basePriority
int	gbisPhase
int	gbisPhasePriority [3]

Detailed Description

Definition at line 35 of file ComputeGlobal.h.

Constructor & Destructor Documentation

ComputeGlobal()

ComputeGlobal::ComputeGlobal	(	ComputeID	c,
		ComputeMgr *	m
	)

Definition at line 38 of file ComputeGlobal.C.

References SimParameters::colvarsOn, SimParameters::CUDASOAintegrateMode, DebugM, SimParameters::FMAOn, SimParameters::fullDirectOn, SimParameters::GBISOn, SimParameters::GBISserOn, SimParameters::IMDon, PatchMap::numPatches(), PatchMap::Object(), Node::Object(), ReductionMgr::Object(), SimParameters::PMEOn, REDUCTIONS_BASIC, REDUCTIONS_GPURESIDENT, ResizeArray< Elem >::resize(), Node::simParameters, SimParameters::tclForcesOn, and ReductionMgr::willSubmit().

  : ComputeHomePatches(c)
{
  DebugM(3,"Constructing client\n");
  aid.resize(0);
  gdef.resize(0);
  comm = m;
  firsttime = 1;
  isRequested = 0;
  isRequestedAllocSize = 0;
  endRequested = 0;
  numGroupsRequested = 0;
  SimParameters *sp = Node::Object()->simParameters;
  dofull = (sp->GBISserOn || sp->GBISOn || sp->fullDirectOn || sp->FMAOn || sp->PMEOn);
  forceSendEnabled = 0;
  if ( sp->tclForcesOn ) forceSendEnabled = 1;
  if ( sp->colvarsOn ) forceSendEnabled = 1;
  if ( sp->IMDon ) forceSendEnabled = 1;
  forceSendActive = 0;
  fid.resize(0);
  totalForce.resize(0);
  gfcount = 0;
  groupTotalForce.resize(0);
  if(Node::Object()->simParameters->CUDASOAintegrate) {
    reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_GPURESIDENT);
  } else {
    reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC);
  }
  int numPatches = PatchMap::Object()->numPatches();
  forcePtrs = new Force*[numPatches];
  atomPtrs = new FullAtom*[numPatches];
  for ( int i = 0; i < numPatches; ++i ) { forcePtrs[i] = 0; atomPtrs[i] = 0; }

  if (sp->CUDASOAintegrateMode) {
    // Allocate memory for numPatches to access SOA data
    mass_soa = new float*[numPatches];
    pos_soa_x = new double*[numPatches];
    pos_soa_y = new double*[numPatches];
    pos_soa_z = new double*[numPatches];
    force_soa_x = new double*[numPatches];
    force_soa_y = new double*[numPatches];
    force_soa_z = new double*[numPatches];
    transform_soa_i = new int*[numPatches];
    transform_soa_j = new int*[numPatches];
    transform_soa_k = new int*[numPatches];
    for ( int i = 0; i < numPatches; ++i ) { 
      mass_soa[i] = NULL;
      pos_soa_x[i] = NULL; 
      pos_soa_y[i] = NULL; 
      pos_soa_z[i] = NULL;
      force_soa_x[i] = NULL;
      force_soa_y[i] = NULL;
      force_soa_z[i] = NULL;
      transform_soa_i[i] = NULL; 
      transform_soa_j[i] = NULL; 
      transform_soa_k[i] = NULL;
    }
  } else {
    mass_soa = NULL;
    pos_soa_x = NULL;
    pos_soa_y = NULL;
    pos_soa_z = NULL;
    force_soa_x = NULL;
    force_soa_y = NULL;
    force_soa_z = NULL;
    transform_soa_i = NULL;
    transform_soa_j = NULL;
    transform_soa_k = NULL;
  }
  gridForcesPtrs = new ForceList **[numPatches];
  numGridObjects = numActiveGridObjects = 0;
  for ( int i = 0; i < numPatches; ++i ) {
    forcePtrs[i] = NULL; atomPtrs[i] = NULL;
    gridForcesPtrs[i] = NULL;
  }
}

~ComputeGlobal()

ComputeGlobal::~ComputeGlobal ( )

virtual

Definition at line 115 of file ComputeGlobal.C.

{
  delete[] isRequested;
  delete[] forcePtrs;
  deleteGridObjects();
  delete[] gridForcesPtrs;
  delete[] atomPtrs;
  delete reduction;

  if(mass_soa) delete [] mass_soa;
  if(pos_soa_x) delete [] pos_soa_x;
  if(pos_soa_y) delete [] pos_soa_y;
  if(pos_soa_z) delete [] pos_soa_z;
  if(force_soa_x) delete [] force_soa_x;
  if(force_soa_y) delete [] force_soa_y;
  if(force_soa_z) delete [] force_soa_z;
  if(transform_soa_i) delete [] transform_soa_i;
  if(transform_soa_j) delete [] transform_soa_j;
  if(transform_soa_k) delete [] transform_soa_k;
}

Member Function Documentation

doWork()

void ComputeGlobal::doWork ( void  )

virtual

Reimplemented from Compute.

Definition at line 519 of file ComputeGlobal.C.

References ResizeArray< Elem >::add(), ResizeArrayIter< T >::begin(), ComputeGlobalDataMsg::count, SimParameters::CUDASOAintegrate, DebugM, ComputeMgr::enableComputeGlobalResults(), ResizeArrayIter< T >::end(), endi(), SimParameters::globalMasterFrequency, ComputeHomePatches::hasPatchZero, ComputeGlobalDataMsg::lat, Node::Object(), ComputeGlobalDataMsg::patchcount, ComputeHomePatches::patchList, ComputeMgr::recvComputeGlobalResults(), ComputeMgr::sendComputeGlobalData(), Node::simParameters, and ComputeGlobalDataMsg::step.

{
  DebugM(2,"doWork thread " << CthGetToken(CthSelf())->serialNo << "\n");

  SimParameters *sp = Node::Object()->simParameters;
  ResizeArrayIter<PatchElem> ap(patchList);
  FullAtom **t = atomPtrs;
  int step = patchList[0].p->flags.step;
  if((step % sp->globalMasterFrequency) ==0)
    {
      DebugM(3,"doWork for step " << step <<"\n"<<endi);
      // if(sp->CUDASOAintegrateOn) {
      //   hasPatchZero = 0;
      // }

      for (ap = ap.begin(); ap != ap.end(); ap++) {
        CompAtom *x = (*ap).positionBox->open();
        t[(*ap).patchID] = (*ap).p->getAtomList().begin();

        if (sp->CUDASOAintegrate) {
          // Assigne the pointer to SOA data structure
          PatchID pId = (*ap).patchID;
          mass_soa[pId] = (*ap).p->patchDataSOA.mass;
          pos_soa_x[pId] = (*ap).p->patchDataSOA.pos_x;
          pos_soa_y[pId] = (*ap).p->patchDataSOA.pos_y;
          pos_soa_z[pId] = (*ap).p->patchDataSOA.pos_z;
          transform_soa_i[pId] = (*ap).p->patchDataSOA.transform_i;
          transform_soa_j[pId] = (*ap).p->patchDataSOA.transform_j;
          transform_soa_k[pId] = (*ap).p->patchDataSOA.transform_k;
          // if(sp->CUDASOAintegrateOn && (pId == 0)) {
          //   hasPatchZero = 1;
          // }
        }
      }

      if(!firsttime) {
        //    CkPrintf("*** Start NoFirstTime on PE %d \n", CkMyPe());
        sendData();
        //    CkPrintf("*** End NoFirstTime on PE %d \n", CkMyPe());
      } else {
        //    CkPrintf("*** Start FirstTime on PE %d \n", CkMyPe());
        if ( hasPatchZero ) {
          ComputeGlobalDataMsg *msg = new ComputeGlobalDataMsg;
          msg->lat.add(patchList[0].p->lattice);
          msg->step = -1;
          msg->count = 1;
          msg->patchcount = 0;
          //      CkPrintf("***DoWork calling sendComputeGlobalData PE %d \n", CkMyPe());
          comm->sendComputeGlobalData(msg);
        }
#ifdef NODEGROUP_FORCE_REGISTER
        else if (sp->CUDASOAintegrate) {

          //      CkPrintf("***DoWork FirstTime barrier 1 on PE %d \n", CkMyPe());
          comm->stowSuspendULT();
          //      CmiNodeBarrier();
          //      CkPrintf("***DoWork FirstTime barrier 2 on PE %d \n", CkMyPe());
          comm->stowSuspendULT();
          //      CkPrintf("***DoWork out of barrier 2 on PE %d \n", CkMyPe());
          //      CmiNodeBarrier();
          ComputeGlobalResultsMsg* resultsMsg = CkpvAccess(ComputeGlobalResultsMsg_instance);
          //      CkPrintf("*** ComputeGlobal::doWork PE (%d) calling recvComputeGlobalResults in doWork at step: %d \n",CkMyPe(), patchList[0].p->flags.step);
          comm->recvComputeGlobalResults(resultsMsg);
        }
#endif // NODEGROUP_FORCE_REGISTER
        firsttime = 0;
        //    CkPrintf("*** ComputeGlobal::doWork PE (%d) calling enableComputeGlobalResults in doWork at step: %d \n",CkMyPe(), patchList[0].p->flags.step);
        comm->enableComputeGlobalResults();

        //    CkPrintf("*** End FirstTime on PE %d \n", CkMyPe());
      }
    }
  else
    {
      DebugM(2,"skipping step "<< step <<"\n"<<endi);
      /* TODO to support CPU side MTS we need to do something to avoid hang           some distillation from sendData(); and the reductions
      ADD_VECTOR_OBJECT(reduction,REDUCTION_EXT_FORCE_NORMAL,extForce);
      ADD_TENSOR_OBJECT(reduction,REDUCTION_VIRIAL_NORMAL,extVirial);
      reduction->submit();
      and as yet undetermined message handling
      */
    }
  DebugM(2,"done with doWork\n");
}

getForceSendActive()

int ComputeGlobal::getForceSendActive ( ) const

inline

Definition at line 46 of file ComputeGlobal.h.

{return forceSendActive;}

recvResults()

void ComputeGlobal::recvResults

(

ComputeGlobalResultsMsg *

msg

)

Definition at line 268 of file ComputeGlobal.C.

References ADD_TENSOR_OBJECT, ADD_VECTOR_OBJECT, ComputeGlobalResultsMsg::aid, ResizeArray< Elem >::begin(), ResizeArrayIter< T >::begin(), SimParameters::CUDASOAintegrate, DebugM, ResizeArray< Elem >::end(), ResizeArrayIter< T >::end(), ComputeGlobalResultsMsg::f, ComputeGlobalResultsMsg::gforce, SimParameters::globalMasterFrequency, SimParameters::globalMasterScaleByFrequency, Transform::i, LocalID::index, Transform::j, Transform::k, AtomMap::localID(), FullAtom::mass, NAMD_bug(), ComputeGlobalResultsMsg::newaid, ComputeGlobalResultsMsg::newgdef, ComputeGlobalResultsMsg::newgridobjid, Results::normal, notUsed, AtomMap::Object(), Node::Object(), outer(), ComputeHomePatches::patchList, LocalID::pid, CompAtom::position, ComputeGlobalResultsMsg::reconfig, ComputeGlobalResultsMsg::resendCoordinates, ResizeArray< Elem >::resize(), Lattice::reverse_transform(), Node::simParameters, ResizeArray< Elem >::size(), SubmitReduction::submit(), ComputeGlobalResultsMsg::totalforces, FullAtom::transform, Vector::x, Vector::y, and Vector::z.

Referenced by ComputeMgr::recvComputeGlobalResults().

                                                            {
  DebugM(3,"Receiving results (" << msg->aid.size() << " forces, "
         << msg->newgdef.size() << " new group atoms) on client thread " << CthGetToken(CthSelf())->serialNo <<" msg->resendCoordinates " << msg->resendCoordinates << " msg->totalforces " << msg->totalforces<< "\n");

  forceSendActive = msg->totalforces;
  if ( forceSendActive && ! forceSendEnabled ) NAMD_bug("ComputeGlobal::recvResults forceSendActive without forceSendEnabled");

  // set the forces only if we aren't going to resend the data
  int setForces = !msg->resendCoordinates;
  SimParameters *sp = Node::Object()->simParameters;
  if(setForces) { // we are requested to 
    // Store forces to patches
    AtomMap *atomMap = AtomMap::Object();
    const Lattice & lattice = patchList[0].p->lattice;
    ResizeArrayIter<PatchElem> ap(patchList);
    Force **f = forcePtrs;
    FullAtom **t = atomPtrs;
    Force extForce = 0.;
    Tensor extVirial;

    for (ap = ap.begin(); ap != ap.end(); ap++) {
      (*ap).r = (*ap).forceBox->open();
      f[(*ap).patchID] = (*ap).r->f[Results::normal];
      t[(*ap).patchID] = (*ap).p->getAtomList().begin();

      if (sp->CUDASOAintegrate) {
        // Assigne the pointer to SOA data structure
        PatchID pId = (*ap).patchID;
        mass_soa[pId] = (*ap).p->patchDataSOA.mass;
        force_soa_x[pId] = (*ap).p->patchDataSOA.f_global_x;
        force_soa_y[pId] = (*ap).p->patchDataSOA.f_global_y;
        force_soa_z[pId] = (*ap).p->patchDataSOA.f_global_z;
        transform_soa_i[pId] = (*ap).p->patchDataSOA.transform_i;
        transform_soa_j[pId] = (*ap).p->patchDataSOA.transform_j;
        transform_soa_k[pId] = (*ap).p->patchDataSOA.transform_k;
      }
    }


    AtomIDList::iterator a = msg->aid.begin();
    AtomIDList::iterator a_e = msg->aid.end();
    ForceList::iterator f2 = msg->f.begin();
    if(sp->globalMasterScaleByFrequency)
      for ( ; a != a_e; ++a, ++f2 ) {
        Force f_atom;
        f_atom = (*f2);
        f_atom.x*=(float) sp->globalMasterFrequency;
        f_atom.y*=(float) sp->globalMasterFrequency;
        f_atom.z*=(float) sp->globalMasterFrequency;
      }
    if (sp->CUDASOAintegrate) {
      LocalID localID;
      PatchID lpid;
      int lidx;
      Position x_orig, x_atom;
      Transform trans;
      Force f_atom;
      for ( ; a != a_e; ++a, ++f2 ) {
        DebugM(1,"processing atom "<<(*a)<<", F="<<(*f2)<<"...\n");
        /* XXX if (*a) is out of bounds here we get a segfault */
        localID = atomMap->localID(*a);
        lpid = localID.pid;
        lidx = localID.index;
        if ( lpid == notUsed || ! f[lpid] ) continue;
        f_atom = (*f2);
        // printf("NAMD3-recv: atom %d, Before Force (%8.6f, %8.6f, %8.6f) \n", 
        //   *a, force_soa_x[lpid][lidx], force_soa_y[lpid][lidx], force_soa_z[lpid][lidx]);
        //      printf("NAMD3-recv: atom %d, Added Force (%8.6f, %8.6f, %8.6f) \n", *a, f_atom.x, f_atom.y, f_atom.z);
        force_soa_x[lpid][lidx] += f_atom.x;
        force_soa_y[lpid][lidx] += f_atom.y;
        force_soa_z[lpid][lidx] += f_atom.z;
#ifndef USE_GLOBALMASTER_VIRIAL_KERNEL
        x_orig.x = pos_soa_x[lpid][lidx];
        x_orig.y = pos_soa_y[lpid][lidx];
        x_orig.z = pos_soa_z[lpid][lidx];
        trans.i = transform_soa_i[lpid][lidx];
        trans.j = transform_soa_j[lpid][lidx];
        trans.k = transform_soa_k[lpid][lidx];
        x_atom = lattice.reverse_transform(x_orig,trans);
        extForce += f_atom;
        extVirial += outer(f_atom,x_atom);
#endif
      }
    } else {
      for ( ; a != a_e; ++a, ++f2 ) {
        DebugM(1,"processing atom "<<(*a)<<", F="<<(*f2)<<"...\n");
        /* XXX if (*a) is out of bounds here we get a segfault */
        LocalID localID = atomMap->localID(*a);
        if ( localID.pid == notUsed || ! f[localID.pid] ) continue;
        Force f_atom = (*f2);
        // printf("NAMD3-recv: atom %d, Before Force (%8.6f, %8.6f, %8.6f) \n", 
        //   *a, f[localID.pid][localID.index].x, f[localID.pid][localID.index].y, f[localID.pid][localID.index].z);
        // printf("NAMD3-recv: atom %d, Added Force (%8.6f, %8.6f, %8.6f) \n", *a, f_atom.x, f_atom.y, f_atom.z);
        f[localID.pid][localID.index] += f_atom;
        FullAtom &atom = t[localID.pid][localID.index];
        Position x_orig = atom.position;
        Transform trans = atom.transform;
        Position x_atom = lattice.reverse_transform(x_orig,trans);
        extForce += f_atom;
        extVirial += outer(f_atom,x_atom);
      }
    }
    DebugM(1,"done with the loop\n");

    // calculate forces for atoms in groups
    AtomIDList::iterator g_i, g_e;
    g_i = gdef.begin(); g_e = gdef.end();
    ForceList::iterator gf_i = msg->gforce.begin();
    //iout << iDEBUG << "recvResults\n" << endi;
    if (sp->CUDASOAintegrate) {
      LocalID localID;
      PatchID lpid;
      int lidx;
      Position x_orig, x_atom;
      Transform trans;
      Force f_atom;
      for ( ; g_i != g_e; ++g_i, ++gf_i ) {
        //iout << iDEBUG << *gf_i << '\n' << endi;
        Vector accel = (*gf_i);
        for ( ; *g_i != -1; ++g_i ) {
          //iout << iDEBUG << *g_i << '\n' << endi;
          localID = atomMap->localID(*g_i);
          lpid = localID.pid;
          lidx = localID.index;
          if ( lpid == notUsed || ! f[lpid] ) continue;
          f_atom = accel * mass_soa[lpid][lidx];
#if 0     
          if (*g_i < 20) {
            CkPrintf("NAMD3-recv: group %d, Before Force (%8.6f, %8.6f, %8.6f) \n", 
              *g_i, force_soa_x[lpid][lidx], force_soa_y[lpid][lidx], force_soa_z[lpid][lidx]);
            CkPrintf("NAMD3-recv: group %d, Added Force (%8.6f, %8.6f, %8.6f) \n", *g_i, f_atom.x, f_atom.y, f_atom.z);
          }
#endif    
          force_soa_x[lpid][lidx] += f_atom.x;
          force_soa_y[lpid][lidx] += f_atom.y;
          force_soa_z[lpid][lidx] += f_atom.z;
#ifndef USE_GLOBALMASTER_VIRIAL_KERNEL
          x_orig.x = pos_soa_x[lpid][lidx];
          x_orig.y = pos_soa_y[lpid][lidx];
          x_orig.z = pos_soa_z[lpid][lidx];
          trans.i = transform_soa_i[lpid][lidx];
          trans.j = transform_soa_j[lpid][lidx];
          trans.k = transform_soa_k[lpid][lidx];
          x_atom = lattice.reverse_transform(x_orig,trans);
          extForce += f_atom;
          extVirial += outer(f_atom,x_atom);
#endif
        }
      }
    } else {
      for ( ; g_i != g_e; ++g_i, ++gf_i ) {
        //iout << iDEBUG << *gf_i << '\n' << endi;
        Vector accel = (*gf_i);
        for ( ; *g_i != -1; ++g_i ) {
          //iout << iDEBUG << *g_i << '\n' << endi;
          LocalID localID = atomMap->localID(*g_i);
          if ( localID.pid == notUsed || ! f[localID.pid] ) continue;
          FullAtom &atom = t[localID.pid][localID.index];
          Force f_atom = accel * atom.mass;
#if 0
          if (*g_i < 20) {
            CkPrintf("NAMD2-recv: group %d, Before Force (%8.6f, %8.6f, %8.6f) \n", 
              *g_i, f[localID.pid][localID.index].x, f[localID.pid][localID.index].y, f[localID.pid][localID.index].z);
            CkPrintf("NAMD2-recv: group %d, Added Force (%8.6f, %8.6f, %8.6f) \n", *g_i, f_atom.x, f_atom.y, f_atom.z);
          }
#endif    
          f[localID.pid][localID.index] += f_atom;

          Position x_orig = atom.position;
          Transform trans = atom.transform;
          Position x_atom = lattice.reverse_transform(x_orig,trans);
          extForce += f_atom;
          extVirial += outer(f_atom,x_atom);
        }
      }
    }
    DebugM(1,"done with the groups\n");

    if (numActiveGridObjects > 0) {
      applyGridObjectForces(msg, &extForce, &extVirial);
    }
    //    printf("Finish receiving at step: %d ####################################################\n", 
    //      patchList[0].p->flags.step);

    ADD_VECTOR_OBJECT(reduction,REDUCTION_EXT_FORCE_NORMAL,extForce);
    ADD_TENSOR_OBJECT(reduction,REDUCTION_VIRIAL_NORMAL,extVirial);
    reduction->submit();
  }
  // done setting the forces, close boxes below

  // Get reconfiguration if present
  if ( msg->reconfig ) {
    DebugM(3,"Reconfiguring\n");
    configure(msg->newaid, msg->newgdef, msg->newgridobjid);
  }

  // send another round of data if requested

  if(msg->resendCoordinates) {
    DebugM(3,"Sending requested data right away\n");
    //    CkPrintf("*** Resending data on PE %d \n", CkMyPe());
    sendData();
  }

  groupTotalForce.resize(numGroupsRequested);
  for ( int i=0; i<numGroupsRequested; ++i ) groupTotalForce[i] = 0;
  DebugM(3,"resized\n");
  if(setForces) {
    DebugM(3,"setting forces\n");    
    ResizeArrayIter<PatchElem> ap(patchList);
    Force **f = forcePtrs;
    FullAtom **t = atomPtrs;
    for (ap = ap.begin(); ap != ap.end(); ap++) {
      CompAtom *x;
      PatchID pId = (*ap).patchID;
      if (!sp->CUDASOAintegrate) {
        (*ap).positionBox->close(&x);
        (*ap).forceBox->close(&((*ap).r));
         DebugM(1,"closing boxes\n");
      }
      f[pId] = 0;
      t[pId] = 0;
      if (sp->CUDASOAintegrate) {
        // XXX Possibly code below is needed by SOAintegrate mode
        mass_soa[pId] = NULL;
        pos_soa_x[pId] = NULL;
        pos_soa_y[pId] = NULL;
        pos_soa_z[pId] = NULL;
        force_soa_x[pId] = NULL;
        force_soa_y[pId] = NULL;
        force_soa_z[pId] = NULL;
        transform_soa_i[pId] = NULL;
        transform_soa_j[pId] = NULL;
        transform_soa_k[pId] = NULL;
        DebugM(2,"nulling ptrs\n");
      }
    }
    DebugM(3,"done setting forces\n");    
  }

  #ifdef NODEGROUP_FORCE_REGISTER
  if (!sp->CUDASOAintegrate) {
    // CUDASOAintegrate handles this on PE 0 in sendComputeGlobalResults
    delete msg;
  }
  #else
  delete msg;
  #endif
  DebugM(3,"Done processing results\n");
}

saveTotalForces()

void ComputeGlobal::saveTotalForces

(

HomePatch *

homePatch

)

Definition at line 968 of file ComputeGlobal.C.

References SimParameters::accelMDDebugOn, SimParameters::accelMDdihe, SimParameters::accelMDOn, ResizeArray< Elem >::add(), Results::amdf, ResizeArray< Elem >::begin(), SimParameters::CUDASOAintegrate, ResizeArray< Elem >::end(), Patch::f, PatchDataSOA::f_normal_x, PatchDataSOA::f_normal_y, PatchDataSOA::f_normal_z, PatchDataSOA::f_saved_nbond_x, PatchDataSOA::f_saved_nbond_y, PatchDataSOA::f_saved_nbond_z, PatchDataSOA::f_saved_slow_x, PatchDataSOA::f_saved_slow_y, PatchDataSOA::f_saved_slow_z, intpair::first, SimParameters::fixedAtomsOn, NAMD_bug(), Results::nbond, Results::normal, Patch::numAtoms, Node::Object(), Patch::patchID, ComputeHomePatches::patchList, intpair::second, Node::simParameters, Results::slow, Vector::x, Vector::y, and Vector::z.

Referenced by Sequencer::integrate(), Sequencer::integrate_SOA(), and Sequencer::minimize().

{
  if ( ! forceSendEnabled ) NAMD_bug("ComputeGlobal::saveTotalForces called unexpectedly");
  if ( ! forceSendActive ) return;

  SimParameters *simParms = Node::Object()->simParameters;
  if ( simParms->accelMDOn && simParms->accelMDDebugOn && simParms->accelMDdihe ) {
    int num=homePatch->numAtoms;
    FullAtomList &atoms = homePatch->atom;
    ForceList &af=homePatch->f[Results::amdf];

    for (int i=0; i<num; ++i) {
      int index = atoms[i].id;
      if (index < endRequested && isRequested[index] & 1) {
        fid.add(index);
        totalForce.add(af[i]);
      }
    }
    return;
  }

  //  printf("Start saving force at step: %d ####################################################\n", 
  //        patchList[0].p->flags.step);
  int fixedAtomsOn = simParms->fixedAtomsOn;
  int num=homePatch->numAtoms;
  FullAtomList &atoms = homePatch->atom;
  ForceList &f1=homePatch->f[Results::normal], &f2=homePatch->f_saved[Results::nbond],
            &f3=homePatch->f_saved[Results::slow];
  
  double *f1_soa_x = homePatch->patchDataSOA.f_normal_x;
  double *f1_soa_y = homePatch->patchDataSOA.f_normal_y;
  double *f1_soa_z = homePatch->patchDataSOA.f_normal_z;
  double *f2_soa_x = homePatch->patchDataSOA.f_saved_nbond_x;
  double *f2_soa_y = homePatch->patchDataSOA.f_saved_nbond_y;
  double *f2_soa_z = homePatch->patchDataSOA.f_saved_nbond_z;
  double *f3_soa_x = homePatch->patchDataSOA.f_saved_slow_x;
  double *f3_soa_y = homePatch->patchDataSOA.f_saved_slow_y;
  double *f3_soa_z = homePatch->patchDataSOA.f_saved_slow_z;
  int hasSOA = (simParms->CUDASOAintegrate);
  Force f_sum;
  double f_sum_x, f_sum_y, f_sum_z;
  
  #if 0
  for (int i=0; i<num; ++i) {
    int index = atoms[i].id;
    if (index < 20) {
      if (hasSOA) {
        CkPrintf("ForceSaved: atom %d, ForceN  (%8.6f, %8.6f, %8.6f) \n", index, f1_soa_x[i], f1_soa_y[i], f1_soa_z[i]);
        CkPrintf("            atom %d, ForceNB (%8.6f, %8.6f, %8.6f) \n", index, f2_soa_x[i], f2_soa_y[i], f2_soa_z[i]);
        CkPrintf("            atom %d, ForceSL (%8.6f, %8.6f, %8.6f) \n", index, f3_soa_x[i], f3_soa_y[i], f3_soa_z[i]);
      } else {
        CkPrintf("ForceSaved: atom %d, ForceN  (%8.6f, %8.6f, %8.6f) \n", index, f1[i].x, f1[i].y, f1[i].z);
        CkPrintf("            atom %d, ForceNB (%8.6f, %8.6f, %8.6f) \n", index, f2[i].x, f2[i].y, f2[i].z);
        // not memory safe to access slow forces all the time like this
        //        CkPrintf("            atom %d, ForceSL (%8.6f, %8.6f, %8.6f) \n", index, f3[i].x, f3[i].y, f3[i].z);
      }
    }
  }

  printf("PE, PId (%d, %d) Stop saving at step: %d ####################################################\n", 
    CkMyPe(), homePatch->patchID, patchList[0].p->flags.step);
  #endif
  if ( ! forceSendActive ) return;
  for (int i=0; i<num; ++i) {
    int index = atoms[i].id;
    char reqflag;
    if (index < endRequested && (reqflag = isRequested[index])) {
      if (hasSOA) {
        f_sum_x = f1_soa_x[i] + f2_soa_x[i];
        f_sum_y = f1_soa_y[i] + f2_soa_y[i];
        f_sum_z = f1_soa_z[i] + f2_soa_z[i];
        if (dofull) {
          f_sum_x += f3_soa_x[i];
          f_sum_y += f3_soa_y[i];
          f_sum_z += f3_soa_z[i];
        }
        f_sum.x = f_sum_x;
        f_sum.y = f_sum_y;
        f_sum.z = f_sum_z;
      } else {
        f_sum = f1[i]+f2[i];
        if (dofull)
          f_sum += f3[i];
      }

      if ( fixedAtomsOn && atoms[i].atomFixed )
        f_sum = 0.;
      
      if ( reqflag  & 1 ) {  // individual atom
        fid.add(index);
        totalForce.add(f_sum);
      }
      if ( reqflag  & 2 ) {  // part of group
        intpair *gpend = gpair.end();
        intpair *gpi = std::lower_bound(gpair.begin(),gpend,intpair(index,0));
        if ( gpi == gpend || gpi->first != index )
          NAMD_bug("ComputeGlobal::saveTotalForces gpair corrupted.");
        do {
          ++gfcount;
          groupTotalForce[gpi->second] += f_sum;
        } while ( ++gpi != gpend && gpi->first == index );
      }
    }
  }
}

---

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

• ComputeGlobal.h
• ComputeGlobal.C