ComputeHomeTuples.h

Go to the documentation of this file.
00001 
00007 #ifndef COMPUTEHOMETUPLES_H
00008 #define COMPUTEHOMETUPLES_H
00009 
00010 #ifdef USE_HOMETUPLES
00011 #include <vector>
00012 #endif
00013 #include "NamdTypes.h"
00014 #include "common.h"
00015 #include "structures.h"
00016 #include "Compute.h"
00017 #include "HomePatch.h"
00018 
00019 #include "Box.h"
00020 #include "OwnerBox.h"
00021 #include "UniqueSet.h"
00022 
00023 #include "Node.h"
00024 #include "SimParameters.h"
00025 #include "PatchMap.inl"
00026 #include "AtomMap.h"
00027 #include "ComputeHomeTuples.h"
00028 #include "PatchMgr.h"
00029 #include "ProxyMgr.h"
00030 #include "HomePatchList.h"
00031 #include "Molecule.h"
00032 #include "Parameters.h"
00033 #include "ReductionMgr.h"
00034 #include "UniqueSet.h"
00035 #include "UniqueSetIter.h"
00036 #include "Priorities.h"
00037 #include "LdbCoordinator.h"
00038 
00039 class TuplePatchElem {
00040   public:
00041     PatchID patchID;
00042     Patch *p;
00043     Box<Patch,CompAtom> *positionBox;
00044     Box<Patch,CompAtom> *avgPositionBox;
00045     Box<Patch,Results> *forceBox;
00046     CompAtom *x;
00047     CompAtomExt *xExt;
00048     CompAtom *x_avg;
00049     Results *r;
00050     Force *f;
00051     Force *af;
00052 
00053     int hash() const { return patchID; }
00054 
00055   TuplePatchElem(PatchID pid = -1) {
00056     patchID = pid;
00057     p = NULL;
00058     positionBox = NULL;
00059     avgPositionBox = NULL;
00060     forceBox = NULL;
00061     x = NULL;
00062     xExt = NULL;
00063     x_avg = NULL;
00064     r = NULL;
00065     f = NULL;
00066     af = NULL;
00067   }
00068 
00069   TuplePatchElem(Patch *p_param, Compute *cid) {
00070     patchID = p_param->getPatchID();
00071     p = p_param;
00072     positionBox = p_param->registerPositionPickup(cid);
00073     avgPositionBox = p_param->registerAvgPositionPickup(cid);
00074     forceBox = p_param->registerForceDeposit(cid);
00075     x = NULL;
00076     xExt = NULL;
00077     x_avg = NULL;
00078     r = NULL;
00079     f = NULL;
00080     af = NULL;
00081   }
00082     
00083   ~TuplePatchElem() {};
00084 
00085   int operator==(const TuplePatchElem &elem) const {
00086     return (elem.patchID == patchID);
00087   }
00088 
00089   int operator<(const TuplePatchElem &elem) const {
00090     return (patchID < elem.patchID);
00091   }
00092 };
00093 
00094 typedef UniqueSet<TuplePatchElem> TuplePatchList;
00095 typedef UniqueSetIter<TuplePatchElem> TuplePatchListIter;
00096 
00097 class AtomMap;
00098 class ReductionMgr;
00099 
00100 #ifdef MEM_OPT_VERSION
00101 template <class T> struct ElemTraits {
00102   typedef AtomSignature signature;
00103   static signature* get_sig_pointer(Molecule *mol) { return mol->atomSigPool; }
00104   static int get_sig_id(const CompAtomExt &a) { return a.sigId; }
00105 };
00106 
00107 template <> struct ElemTraits <ExclElem> {
00108   typedef ExclusionSignature signature;
00109   static signature* get_sig_pointer(Molecule *mol) { return mol->exclSigPool; }
00110   static int get_sig_id(const CompAtomExt &a) { return a.exclId; }
00111 };
00112 #endif
00113 
00114 #ifdef USE_HOMETUPLES
00115 //
00116 // Simple base class for HomeTuples and SelfTuples that stores the type of the tuple
00117 //
00118 class Tuples {
00119 private:
00120   int type;
00121 protected:
00122   Tuples(int type) : type(type) {}
00123 public:
00124   // Tuple types
00125   enum {BOND=0, ANGLE, DIHEDRAL, IMPROPER, EXCLUSION, CROSSTERM, NUM_TUPLE_TYPES};
00126 
00127   int getType() {return type;}
00128   virtual void submitTupleCount(SubmitReduction *reduction, int tupleCount)=0;
00129   // virtual void copyTupleData(void* tupleData)=0;
00130   virtual int getNumTuples()=0;
00131   virtual void* getTupleList()=0;
00132   virtual void loadTuples(TuplePatchList& tuplePatchList, const char* isBasePatch, AtomMap *atomMap,
00133       const std::vector<int>& pids = std::vector<int>())=0;
00134 };
00135 
00136 //
00137 // HomeTuples class. These are created and stored in ComputeBondedCUDA::registerCompute()
00138 // e.g.: new HomeTuples<BondElem, Bond, BondValue>(BOND)
00139 //
00140 template <class T, class S, class P> class HomeTuples : public Tuples {
00141   protected:
00142     std::vector<T> tupleList;
00143 
00144   public:
00145 
00146     HomeTuples(int type=-1) : Tuples(type) {}
00147 
00148 #if __cplusplus < 201103L
00149 #define final
00150 #endif
00151 
00152     virtual void* getTupleList() final {
00153       return (void*)tupleList.data();
00154     }
00155 
00156     virtual void submitTupleCount(SubmitReduction *reduction, int tupleCount) final {
00157       reduction->item(T::reductionChecksumLabel) += (BigReal)tupleCount;
00158     }
00159 
00160     // virtual void copyTupleData(void* tupleData) final {
00161       // for (int i=0;i < tupleList.size();i++) {
00162       //   tupleData[i] = 
00163       // }
00164       // T::loadTupleData(tupleData);
00165     // }
00166 
00167     virtual int getNumTuples() final {
00168       return tupleList.size();
00169     }
00170 
00171     virtual void loadTuples(TuplePatchList& tuplePatchList, const char* isBasePatch, AtomMap *atomMap,
00172       const std::vector<int>& pids = std::vector<int>()) {
00173 
00174       if (isBasePatch == NULL) {
00175         NAMD_bug("NULL isBasePatch detected in HomeTuples::loadTuples()");
00176       }
00177 
00178       int numTuples;
00179 
00180 #ifdef MEM_OPT_VERSION
00181       typename ElemTraits<T>::signature *allSigs;      
00182 #else
00183       int32 **tuplesByAtom;
00184       /* const (need to propagate const) */ S *tupleStructs;
00185 #endif
00186       
00187       const P *tupleValues;
00188       Node *node = Node::Object();
00189       PatchMap *patchMap = PatchMap::Object();
00190       // AtomMap *atomMap = AtomMap::Object();
00191 
00192 #ifdef MEM_OPT_VERSION
00193       allSigs = ElemTraits<T>::get_sig_pointer(node->molecule);
00194 #else      
00195       T::getMoleculePointers(node->molecule,
00196         &numTuples, &tuplesByAtom, &tupleStructs);      
00197 #endif
00198       
00199       T::getParameterPointers(node->parameters, &tupleValues);
00200 
00201       tupleList.clear();
00202 
00203       LocalID aid[T::size];
00204 
00205       const int lesOn = node->simParameters->lesOn;
00206       Real invLesFactor = lesOn ? 
00207                           1.0/node->simParameters->lesFactor :
00208                           1.0;
00209 
00210       // cycle through each patch and gather all tuples
00211       TuplePatchListIter ai(tuplePatchList);
00212       if (pids.size() == 0) ai = ai.begin();
00213 
00214       int numPid = (pids.size() == 0) ? tuplePatchList.size() : pids.size();
00215 
00216       for (int ipid=0;ipid < numPid;ipid++) {
00217         // Patch *patch;
00218         int numAtoms;
00219         CompAtomExt *atomExt;
00220         // Take next patch
00221         if (pids.size() == 0) {
00222           Patch* patch = (*ai).p;
00223           numAtoms = patch->getNumAtoms();
00224           atomExt = (*ai).xExt;
00225           ai++;
00226         } else {
00227           TuplePatchElem *tpe = tuplePatchList.find(TuplePatchElem(pids[ipid]));
00228           Patch* patch = tpe->p;
00229           numAtoms = patch->getNumAtoms();
00230           atomExt = tpe->xExt;          
00231         }
00232    
00233         // cycle through each atom in the patch and load up tuples
00234         for (int j=0; j < numAtoms; j++)
00235         {
00236           /* cycle through each tuple */
00237 #ifdef MEM_OPT_VERSION
00238           typename ElemTraits<T>::signature *thisAtomSig =
00239                    &allSigs[ElemTraits<T>::get_sig_id(atomExt[j])];
00240           TupleSignature *allTuples;
00241           T::getTupleInfo(thisAtomSig, &numTuples, &allTuples);
00242           for(int k=0; k<numTuples; k++) {
00243             T t(atomExt[j].id, &allTuples[k], tupleValues);
00244 #else
00245           /* get list of all tuples for the atom */
00246           int32 *curTuple = tuplesByAtom[atomExt[j].id];
00247           for( ; *curTuple != -1; ++curTuple) {
00248             T t(&tupleStructs[*curTuple],tupleValues);
00249 #endif            
00250             register int i;
00251             aid[0] = atomMap->localID(t.atomID[0]);
00252             int homepatch = aid[0].pid;
00253             int samepatch = 1;
00254             int has_les = lesOn && node->molecule->get_fep_type(t.atomID[0]);
00255             for (i=1; i < T::size; i++) {
00256               aid[i] = atomMap->localID(t.atomID[i]);
00257               samepatch = samepatch && ( homepatch == aid[i].pid );
00258               has_les |= lesOn && node->molecule->get_fep_type(t.atomID[i]);
00259             }
00260             if ( samepatch ) continue;
00261             t.scale = has_les ? invLesFactor : 1;
00262             for (i=1; i < T::size; i++) {
00263               homepatch = patchMap->downstream(homepatch,aid[i].pid);
00264             }
00265             if ( homepatch != notUsed && isBasePatch[homepatch] ) {
00266               TuplePatchElem *p;
00267               for (i=0; i < T::size; i++) {
00268                 t.p[i] = p = tuplePatchList.find(TuplePatchElem(aid[i].pid));
00269                 if ( ! p ) {
00270 #ifdef MEM_OPT_VERSION
00271                   iout << iWARN << "Tuple with atoms ";
00272 #else
00273                   iout << iWARN << "Tuple " << *curTuple << " with atoms ";
00274 #endif
00275                   int erri;
00276                   for( erri = 0; erri < T::size; erri++ ) {
00277                     iout << t.atomID[erri] << "(" <<  aid[erri].pid << ") ";
00278                   }
00279                   iout << "missing patch " << aid[i].pid << "\n" << endi;
00280                   break;
00281                 }
00282                 t.localIndex[i] = aid[i].index;
00283               }
00284               if ( ! p ) continue;
00285 #ifdef MEM_OPT_VERSION
00286               //avoid adding Tuples whose atoms are all fixed
00287               if(node->simParameters->fixedAtomsOn && !node->simParameters->fixedAtomsForces) {
00288                 int allfixed = 1;
00289                 for(i=0; i<T::size; i++){
00290                   CompAtomExt *one = &(t.p[i]->xExt[aid[i].index]);
00291                   allfixed = allfixed & one->atomFixed;
00292                 }
00293                 if(!allfixed) tupleList.push_back(t);
00294               }else{
00295                 tupleList.push_back(t);
00296               }
00297 #else
00298               tupleList.push_back(t);
00299 #endif               
00300             }
00301           }
00302         }
00303       }
00304     }
00305 
00306 };
00307 #endif
00308 
00309 template <class T, class S, class P> class ComputeHomeTuples : public Compute {
00310 
00311   protected:
00312   
00313 #ifndef USE_HOMETUPLES
00314     virtual void loadTuples(void) {
00315       int numTuples;
00316 
00317       #ifdef MEM_OPT_VERSION
00318       typename ElemTraits<T>::signature *allSigs;      
00319       #else
00320       int32 **tuplesByAtom;
00321       /* const (need to propagate const) */ S *tupleStructs;
00322       #endif
00323       
00324       const P *tupleValues;
00325       Node *node = Node::Object();
00326 
00327       #ifdef MEM_OPT_VERSION
00328       allSigs = ElemTraits<T>::get_sig_pointer(node->molecule);
00329       #else      
00330       T::getMoleculePointers(node->molecule,
00331                     &numTuples, &tuplesByAtom, &tupleStructs);      
00332       #endif
00333       
00334       T::getParameterPointers(node->parameters, &tupleValues);
00335 
00336       tupleList.resize(0);
00337 
00338       LocalID aid[T::size];
00339 
00340       const int lesOn = node->simParameters->lesOn;
00341       Real invLesFactor = lesOn ? 
00342                           1.0/node->simParameters->lesFactor :
00343                           1.0;
00344 
00345       // cycle through each patch and gather all tuples
00346       TuplePatchListIter ai(tuplePatchList);
00347     
00348       for ( ai = ai.begin(); ai != ai.end(); ai++ )
00349       {
00350         // CompAtom *atom = (*ai).x;
00351         Patch *patch = (*ai).p;
00352         int numAtoms = patch->getNumAtoms();
00353         CompAtomExt *atomExt = (*ai).xExt; //patch->getCompAtomExtInfo();
00354     
00355         // cycle through each atom in the patch and load up tuples
00356         for (int j=0; j < numAtoms; j++)
00357         {              
00358            /* cycle through each tuple */
00359            #ifdef MEM_OPT_VERSION
00360            typename ElemTraits<T>::signature *thisAtomSig =
00361                    &allSigs[ElemTraits<T>::get_sig_id(atomExt[j])];
00362            TupleSignature *allTuples;
00363            T::getTupleInfo(thisAtomSig, &numTuples, &allTuples);
00364            for(int k=0; k<numTuples; k++) {
00365                T t(atomExt[j].id, &allTuples[k], tupleValues);
00366            #else
00367            /* get list of all tuples for the atom */
00368            int32 *curTuple = tuplesByAtom[atomExt[j].id];
00369            for( ; *curTuple != -1; ++curTuple) {             
00370              T t(&tupleStructs[*curTuple],tupleValues);
00371            #endif            
00372              register int i;
00373              aid[0] = atomMap->localID(t.atomID[0]);
00374              int homepatch = aid[0].pid;
00375              int samepatch = 1;
00376              int has_les = lesOn && node->molecule->get_fep_type(t.atomID[0]);
00377              for (i=1; i < T::size; i++) {
00378                  aid[i] = atomMap->localID(t.atomID[i]);
00379                  samepatch = samepatch && ( homepatch == aid[i].pid );
00380                  has_les |= lesOn && node->molecule->get_fep_type(t.atomID[i]);
00381              }
00382              if ( samepatch ) continue;
00383              t.scale = has_les ? invLesFactor : 1;
00384              for (i=1; i < T::size; i++) {
00385                  homepatch = patchMap->downstream(homepatch,aid[i].pid);
00386              }
00387              if ( homepatch != notUsed && isBasePatch[homepatch] ) {
00388                TuplePatchElem *p;
00389                for (i=0; i < T::size; i++) {
00390                  t.p[i] = p = tuplePatchList.find(TuplePatchElem(aid[i].pid));
00391                  if ( ! p ) {
00392                      #ifdef MEM_OPT_VERSION
00393                      iout << iWARN << "Tuple with atoms ";
00394                      #else
00395                    iout << iWARN << "Tuple " << *curTuple << " with atoms ";
00396                      #endif
00397                    int erri;
00398                    for( erri = 0; erri < T::size; erri++ ) {
00399                      iout << t.atomID[erri] << "(" <<  aid[erri].pid << ") ";
00400                    }
00401                    iout << "missing patch " << aid[i].pid << "\n" << endi;
00402                    break;
00403                  }
00404                  t.localIndex[i] = aid[i].index;
00405                }
00406                if ( ! p ) continue;
00407              #ifdef MEM_OPT_VERSION
00408                //avoid adding Tuples whose atoms are all fixed
00409                if(node->simParameters->fixedAtomsOn &&
00410                   !node->simParameters->fixedAtomsForces) {
00411                  int allfixed = 1;
00412                  for(i=0; i<T::size; i++){
00413                    CompAtomExt *one = &(t.p[i]->xExt[aid[i].index]);
00414                    allfixed = allfixed & one->atomFixed;
00415                  }
00416                  if(!allfixed) tupleList.add(t);
00417                }else{
00418                  tupleList.add(t);
00419                }
00420              #else
00421                tupleList.add(t);
00422              #endif               
00423              }
00424            }
00425         }
00426       }
00427     }
00428 #endif
00429 
00430     int doLoadTuples;
00431   
00432   protected:
00433   
00434 #ifdef USE_HOMETUPLES
00435     HomeTuples<T, S, P>* tuples;
00436     TuplePatchList tuplePatchList;
00437 #else
00438     ResizeArray<T> tupleList;
00439     TuplePatchList tuplePatchList;
00440 #endif
00441 
00442     PatchMap *patchMap;
00443     AtomMap *atomMap;
00444     SubmitReduction *reduction;
00445     int accelMDdoDihe;
00446     SubmitReduction *pressureProfileReduction;
00447     BigReal *pressureProfileData;
00448     int pressureProfileSlabs;
00449     char *isBasePatch;
00450   
00451     ComputeHomeTuples(ComputeID c) : Compute(c) {
00452       patchMap = PatchMap::Object();
00453       atomMap = AtomMap::Object();
00454       reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC);
00455       
00456       SimParameters *params = Node::Object()->simParameters;
00457       accelMDdoDihe=false;
00458       if (params->accelMDOn) {
00459          if (params->accelMDdihe || params->accelMDdual) accelMDdoDihe=true;
00460       }
00461       if (params->pressureProfileOn) {
00462         pressureProfileSlabs = T::pressureProfileSlabs = 
00463           params->pressureProfileSlabs;
00464         int n = T::pressureProfileAtomTypes = params->pressureProfileAtomTypes;
00465         pressureProfileReduction = ReductionMgr::Object()->willSubmit(
00466           REDUCTIONS_PPROF_BONDED, 3*pressureProfileSlabs*((n*(n+1))/2));
00467         int numAtomTypePairs = n*n;
00468         pressureProfileData = new BigReal[3*pressureProfileSlabs*numAtomTypePairs];
00469       } else {
00470         pressureProfileReduction = NULL;
00471         pressureProfileData = NULL;
00472       }
00473       doLoadTuples = false;
00474       isBasePatch = 0;
00475 #ifdef USE_HOMETUPLES
00476       tuples = NULL;
00477 #endif
00478     }
00479 
00480     ComputeHomeTuples(ComputeID c, PatchIDList &pids) : Compute(c) {
00481       patchMap = PatchMap::Object();
00482       atomMap = AtomMap::Object();
00483       reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC);
00484       SimParameters *params = Node::Object()->simParameters;
00485       accelMDdoDihe=false;
00486       if (params->accelMDOn) {
00487          if (params->accelMDdihe || params->accelMDdual) accelMDdoDihe=true;
00488       }
00489       if (params->pressureProfileOn) {
00490         pressureProfileSlabs = T::pressureProfileSlabs = 
00491           params->pressureProfileSlabs;
00492         int n = T::pressureProfileAtomTypes = params->pressureProfileAtomTypes;
00493         pressureProfileReduction = ReductionMgr::Object()->willSubmit(
00494           REDUCTIONS_PPROF_BONDED, 3*pressureProfileSlabs*((n*(n+1))/2));
00495         int numAtomTypePairs = n*n;
00496         pressureProfileData = new BigReal[3*pressureProfileSlabs*numAtomTypePairs];
00497       } else {
00498         pressureProfileReduction = NULL;
00499         pressureProfileData = NULL;
00500       }
00501       doLoadTuples = false;
00502       int nPatches = patchMap->numPatches();
00503       isBasePatch = new char[nPatches];
00504       int i;
00505       for (i=0; i<nPatches; ++i) { isBasePatch[i] = 0; }
00506       for (i=0; i<pids.size(); ++i) { isBasePatch[pids[i]] = 1; }
00507 #ifdef USE_HOMETUPLES
00508       tuples = NULL;
00509 #endif
00510     }
00511 
00512   public:
00513   
00514     virtual ~ComputeHomeTuples() {
00515       delete reduction;
00516       delete [] isBasePatch;
00517       delete pressureProfileReduction;
00518       delete pressureProfileData;
00519 #ifdef USE_HOMETUPLES
00520       if (tuples != NULL) delete tuples;
00521 #endif
00522     }
00523 
00524     //======================================================================
00525     // initialize() - Method is invoked only the first time
00526     // atom maps, patchmaps etc are ready and we are about to start computations
00527     //======================================================================
00528     virtual void initialize(void) {
00529 
00530 #ifdef NAMD_CUDA
00531       ProxyMgr *proxyMgr = ProxyMgr::Object();
00532 #endif
00533 
00534 #ifdef USE_HOMETUPLES
00535       tuples = new HomeTuples<T, S, P>();
00536 #endif
00537 
00538       // Start with empty list
00539       tuplePatchList.clear();
00540     
00541       int nPatches = patchMap->numPatches();
00542       int pid;
00543       for (pid=0; pid<nPatches; ++pid) {
00544         if ( isBasePatch[pid] ) {
00545 #ifdef NAMD_CUDA
00546           proxyMgr->createProxy(pid);
00547 #endif
00548           Patch *patch = patchMap->patch(pid);
00549           tuplePatchList.add(TuplePatchElem(patch, this));
00550         }
00551       }
00552     
00553       // Gather all proxy patches (neighbors, that is)
00554       PatchID neighbors[PatchMap::MaxOneOrTwoAway];
00555     
00556       for (pid=0; pid<nPatches; ++pid) if ( isBasePatch[pid] ) {
00557         int numNeighbors = patchMap->upstreamNeighbors(pid,neighbors);
00558         for ( int i = 0; i < numNeighbors; ++i ) {
00559           if ( ! tuplePatchList.find(TuplePatchElem(neighbors[i])) ) {
00560 #ifdef NAMD_CUDA
00561             proxyMgr->createProxy(neighbors[i]);
00562 #endif
00563             Patch *patch = patchMap->patch(neighbors[i]);
00564             tuplePatchList.add(TuplePatchElem(patch, this));
00565           }
00566         }
00567       }
00568       setNumPatches(tuplePatchList.size());
00569       doLoadTuples = true;
00570 
00571       basePriority = COMPUTE_PROXY_PRIORITY;  // no patch dependence
00572     }
00573 
00574     //======================================================================
00575     // atomUpdate() - Method is invoked after anytime that atoms have been
00576     // changed in patches used by this Compute object.
00577     //======================================================================
00578     void atomUpdate(void) {
00579       doLoadTuples = true;
00580     }
00581 
00582 //-------------------------------------------------------------------
00583 // Routine which is called by enqueued work msg.  It wraps
00584 // actualy Force computation with the apparatus needed
00585 // to get access to atom positions, return forces etc.
00586 //-------------------------------------------------------------------
00587     virtual void doWork(void) {
00588 
00589       LdbCoordinator::Object()->startWork(ldObjHandle);
00590 
00591       // Open Boxes - register that we are using Positions
00592       // and will be depositing Forces.
00593       UniqueSetIter<TuplePatchElem> ap(tuplePatchList);
00594       for (ap = ap.begin(); ap != ap.end(); ap++) {
00595         ap->x = ap->positionBox->open();
00596         ap->xExt = ap->p->getCompAtomExtInfo();
00597         if ( ap->p->flags.doMolly ) ap->x_avg = ap->avgPositionBox->open();
00598         ap->r = ap->forceBox->open();
00599         ap->f = ap->r->f[Results::normal];
00600         if (accelMDdoDihe) ap->af = ap->r->f[Results::amdf]; // for dihedral-only or dual-boost accelMD
00601       } 
00602     
00603       BigReal reductionData[T::reductionDataSize];
00604       int tupleCount = 0;
00605       int numAtomTypes = T::pressureProfileAtomTypes;
00606       int numAtomTypePairs = numAtomTypes*numAtomTypes;
00607     
00608       for ( int i = 0; i < T::reductionDataSize; ++i ) reductionData[i] = 0;
00609       if (pressureProfileData) {
00610         memset(pressureProfileData, 0, 3*pressureProfileSlabs*numAtomTypePairs*sizeof(BigReal));
00611         // Silly variable hiding of the previous iterator
00612         UniqueSetIter<TuplePatchElem> newap(tuplePatchList);
00613         newap = newap.begin();
00614         const Lattice &lattice = newap->p->lattice;
00615         T::pressureProfileThickness = lattice.c().z / pressureProfileSlabs;
00616         T::pressureProfileMin = lattice.origin().z - 0.5*lattice.c().z;
00617       }
00618 
00619       if ( ! Node::Object()->simParameters->commOnly ) {
00620       if ( doLoadTuples ) {
00621 #ifdef USE_HOMETUPLES
00622         tuples->loadTuples(tuplePatchList, isBasePatch, AtomMap::Object());
00623 #else
00624         loadTuples();
00625 #endif
00626         doLoadTuples = false;
00627       }
00628       // take triplet and pass with tuple info to force eval
00629 #ifdef USE_HOMETUPLES
00630       T *al = (T *)tuples->getTupleList();
00631       const int ntuple = tuples->getNumTuples();
00632 #else
00633       T *al = tupleList.begin();
00634       const int ntuple = tupleList.size();
00635 #endif
00636       if ( ntuple ) T::computeForce(al, ntuple, reductionData, pressureProfileData);
00637       tupleCount += ntuple;
00638       }
00639  
00640     LdbCoordinator::Object()->endWork(ldObjHandle);
00641 
00642       T::submitReductionData(reductionData,reduction);
00643       reduction->item(T::reductionChecksumLabel) += (BigReal)tupleCount;
00644       reduction->submit();
00645 
00646       if (pressureProfileReduction) {
00647         // For ease of calculation we stored interactions between types
00648         // i and j in (ni+j).  For efficiency now we coalesce the
00649         // cross interactions so that just i<=j are stored.
00650         const int arraysize = 3*pressureProfileSlabs;
00651         const BigReal *data = pressureProfileData;
00652         for (int i=0; i<numAtomTypes; i++) {
00653           for (int j=0; j<numAtomTypes; j++) {
00654             int ii=i;
00655             int jj=j;
00656             if (ii > jj) { int tmp=ii; ii=jj; jj=tmp; }
00657             const int reductionOffset = 
00658               (ii*numAtomTypes - (ii*(ii+1))/2 + jj)*arraysize;
00659             for (int k=0; k<arraysize; k++) {
00660               pressureProfileReduction->item(reductionOffset+k) += data[k];
00661             }
00662             data += arraysize;
00663           }
00664         }
00665         pressureProfileReduction->submit();
00666       }
00667     
00668       // Close boxes - i.e. signal we are done with Positions and
00669       // AtomProperties and that we are depositing Forces
00670       for (ap = ap.begin(); ap != ap.end(); ap++) {
00671         ap->positionBox->close(&(ap->x));
00672         if ( ap->p->flags.doMolly ) ap->avgPositionBox->close(&(ap->x_avg));
00673         ap->forceBox->close(&(ap->r));
00674       }
00675     }
00676 };
00677 
00678 
00679 #endif
00680 

Generated on Tue Sep 26 01:17:11 2017 for NAMD by  doxygen 1.4.7