ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/brains/SimCreator.cpp
(Generate patch)

Comparing:
trunk/src/brains/SimCreator.cpp (file contents), Revision 1287 by gezelter, Wed Sep 10 18:11:32 2008 UTC vs.
branches/development/src/brains/SimCreator.cpp (file contents), Revision 1627 by gezelter, Tue Sep 13 22:05:04 2011 UTC

# Line 6 | Line 6
6   * redistribute this software in source and binary code form, provided
7   * that the following conditions are met:
8   *
9 < * 1. Acknowledgement of the program authors must be made in any
10 < *    publication of scientific results based in part on use of the
11 < *    program.  An acceptable form of acknowledgement is citation of
12 < *    the article in which the program was described (Matthew
13 < *    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14 < *    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15 < *    Parallel Simulation Engine for Molecular Dynamics,"
16 < *    J. Comput. Chem. 26, pp. 252-271 (2005))
17 < *
18 < * 2. Redistributions of source code must retain the above copyright
9 > * 1. Redistributions of source code must retain the above copyright
10   *    notice, this list of conditions and the following disclaimer.
11   *
12 < * 3. Redistributions in binary form must reproduce the above copyright
12 > * 2. Redistributions in binary form must reproduce the above copyright
13   *    notice, this list of conditions and the following disclaimer in the
14   *    documentation and/or other materials provided with the
15   *    distribution.
# Line 37 | Line 28
28   * arising out of the use of or inability to use software, even if the
29   * University of Notre Dame has been advised of the possibility of
30   * such damages.
31 + *
32 + * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your
33 + * research, please cite the appropriate papers when you publish your
34 + * work.  Good starting points are:
35 + *                                                                      
36 + * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).            
37 + * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
38 + * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008).          
39 + * [4]  Vardeman & Gezelter, in progress (2009).                        
40   */
41  
42   /**
# Line 76 | Line 76
76   #include "antlr/NoViableAltException.hpp"
77  
78   #ifdef IS_MPI
79 + #include "mpi.h"
80   #include "math/ParallelRandNumGen.hpp"
81   #endif
82  
83 < namespace oopse {
83 > namespace OpenMD {
84    
85 <  Globals* SimCreator::parseFile(std::istream& rawMetaDataStream, const std::string& filename, int startOfMetaDataBlock ){
85 >  Globals* SimCreator::parseFile(std::istream& rawMetaDataStream, const std::string& filename, int mdFileVersion, int startOfMetaDataBlock ){
86      Globals* simParams = NULL;
87      try {
88  
# Line 92 | Line 93 | namespace oopse {
93        const int masterNode = 0;
94        int commStatus;
95        if (worldRank == masterNode) {
96 < #endif
97 <                
96 >        commStatus = MPI_Bcast(&mdFileVersion, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
97 > #endif                
98          SimplePreprocessor preprocessor;
99          preprocessor.preprocess(rawMetaDataStream, filename, startOfMetaDataBlock, ppStream);
100                  
# Line 106 | Line 107 | namespace oopse {
107              
108                  
109        } else {
110 +
111 +        commStatus = MPI_Bcast(&mdFileVersion, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
112 +
113          //get stream size
114          commStatus = MPI_Bcast(&streamSize, 1, MPI_LONG, masterNode, MPI_COMM_WORLD);  
115 <                
115 >
116          char* buf = new char[streamSize];
117          assert(buf);
118                  
# Line 116 | Line 120 | namespace oopse {
120          commStatus = MPI_Bcast(buf, streamSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);
121                  
122          ppStream.str(buf);
123 <        delete buf;
123 >        delete [] buf;
124  
125        }
126   #endif            
# Line 214 | Line 218 | namespace oopse {
218        painCave.isFatal = 1;
219        simError();        
220      }
221 <    catch (OOPSEException& e) {
221 >    catch (OpenMDException& e) {
222        sprintf(painCave.errMsg,
223                "%s\n",
224                e.getMessage().c_str());
# Line 229 | Line 233 | namespace oopse {
233        simError();
234      }
235  
236 +    simParams->setMDfileVersion(mdFileVersion);
237      return simParams;
238    }
239    
240    SimInfo*  SimCreator::createSim(const std::string & mdFileName,
241                                    bool loadInitCoords) {
242 <
242 >    
243      const int bufferSize = 65535;
244      char buffer[bufferSize];
245      int lineNo = 0;
# Line 243 | Line 248 | namespace oopse {
248      int metaDataBlockEnd = -1;
249      int i;
250      int mdOffset;
251 +    int mdFileVersion;
252  
253   #ifdef IS_MPI            
254      const int masterNode = 0;
# Line 262 | Line 268 | namespace oopse {
268        mdFile_.getline(buffer, bufferSize);
269        ++lineNo;
270        std::string line = trimLeftCopy(buffer);
271 <      i = CaseInsensitiveFind(line, "<OOPSE");
271 >      i = CaseInsensitiveFind(line, "<OpenMD");
272        if (static_cast<size_t>(i) == string::npos) {
273 +        // try the older file strings to see if that works:
274 +        i = CaseInsensitiveFind(line, "<OOPSE");
275 +      }
276 +      
277 +      if (static_cast<size_t>(i) == string::npos) {
278 +        // still no luck!
279          sprintf(painCave.errMsg,
280 <                "SimCreator: File: %s is not an OOPSE file!\n",
280 >                "SimCreator: File: %s is not a valid OpenMD file!\n",
281                  mdFileName.c_str());
282          painCave.isFatal = 1;
283          simError();
284        }
285 +      
286 +      // found the correct opening string, now try to get the file
287 +      // format version number.
288  
289 +      StringTokenizer tokenizer(line, "=<> \t\n\r");
290 +      std::string fileType = tokenizer.nextToken();
291 +      toUpper(fileType);
292 +
293 +      mdFileVersion = 0;
294 +
295 +      if (fileType == "OPENMD") {
296 +        while (tokenizer.hasMoreTokens()) {
297 +          std::string token(tokenizer.nextToken());
298 +          toUpper(token);
299 +          if (token == "VERSION") {
300 +            mdFileVersion = tokenizer.nextTokenAsInt();
301 +            break;
302 +          }
303 +        }
304 +      }
305 +            
306        //scan through the input stream and find MetaData tag        
307        while(mdFile_.getline(buffer, bufferSize)) {
308          ++lineNo;
# Line 326 | Line 358 | namespace oopse {
358      std::stringstream rawMetaDataStream(mdRawData);
359  
360      //parse meta-data file
361 <    Globals* simParams = parseFile(rawMetaDataStream, mdFileName, metaDataBlockStart+1);
361 >    Globals* simParams = parseFile(rawMetaDataStream, mdFileName, mdFileVersion,
362 >                                   metaDataBlockStart + 1);
363      
364      //create the force field
365      ForceField * ff = ForceFieldFactory::getInstance()->createForceField(simParams->getForceField());
# Line 363 | Line 396 | namespace oopse {
396      }
397      
398      ff->parse(forcefieldFileName);
366    ff->setFortranForceOptions();
399      //create SimInfo
400      SimInfo * info = new SimInfo(ff, simParams);
401  
# Line 381 | Line 413 | namespace oopse {
413      //create the molecules
414      createMolecules(info);
415      
384    
416      //allocate memory for DataStorage(circular reference, need to
417      //break it)
418      info->setSnapshotManager(new SimSnapshotManager(info));
# Line 407 | Line 438 | namespace oopse {
438      
439      if (loadInitCoords)
440        loadCoordinates(info, mdFileName);    
410    
441      return info;
442    }
443    
# Line 471 | Line 501 | namespace oopse {
501                "\tthe number of molecules.  This will not result in a \n"
502                "\tusable division of atoms for force decomposition.\n"
503                "\tEither try a smaller number of processors, or run the\n"
504 <              "\tsingle-processor version of OOPSE.\n", nProcessors, nGlobalMols);
504 >              "\tsingle-processor version of OpenMD.\n", nProcessors, nGlobalMols);
505        
506        painCave.isFatal = 1;
507        simError();
# Line 606 | Line 636 | namespace oopse {
636   #endif
637          
638          stampId = info->getMoleculeStampId(i);
639 <        Molecule * mol = molCreator.createMolecule(info->getForceField(), info->getMoleculeStamp(stampId),
640 <                                                   stampId, i, info->getLocalIndexManager());
639 >        Molecule * mol = molCreator.createMolecule(info->getForceField(),
640 >                                                   info->getMoleculeStamp(stampId),
641 >                                                   stampId, i,
642 >                                                   info->getLocalIndexManager());
643          
644          info->addMolecule(mol);
645          
# Line 634 | Line 666 | namespace oopse {
666      int beginRigidBodyIndex;
667      int beginCutoffGroupIndex;
668      int nGlobalAtoms = info->getNGlobalAtoms();
637
638    /**@todo fixme */
639 #ifndef IS_MPI
640    
641    beginAtomIndex = 0;
642    beginRigidBodyIndex = 0;
643    beginCutoffGroupIndex = 0;
644    
645 #else
646    
647    int nproc;
648    int myNode;
649    
650    myNode = worldRank;
651    MPI_Comm_size(MPI_COMM_WORLD, &nproc);
652    
653    std::vector < int > tmpAtomsInProc(nproc, 0);
654    std::vector < int > tmpRigidBodiesInProc(nproc, 0);
655    std::vector < int > tmpCutoffGroupsInProc(nproc, 0);
656    std::vector < int > NumAtomsInProc(nproc, 0);
657    std::vector < int > NumRigidBodiesInProc(nproc, 0);
658    std::vector < int > NumCutoffGroupsInProc(nproc, 0);
659    
660    tmpAtomsInProc[myNode] = info->getNAtoms();
661    tmpRigidBodiesInProc[myNode] = info->getNRigidBodies();
662    tmpCutoffGroupsInProc[myNode] = info->getNCutoffGroups();
663    
664    //do MPI_ALLREDUCE to exchange the total number of atoms, rigidbodies and cutoff groups
665    MPI_Allreduce(&tmpAtomsInProc[0], &NumAtomsInProc[0], nproc, MPI_INT,
666                  MPI_SUM, MPI_COMM_WORLD);
667    MPI_Allreduce(&tmpRigidBodiesInProc[0], &NumRigidBodiesInProc[0], nproc,
668                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
669    MPI_Allreduce(&tmpCutoffGroupsInProc[0], &NumCutoffGroupsInProc[0], nproc,
670                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
669      
670      beginAtomIndex = 0;
671      beginRigidBodyIndex = 0;
672      beginCutoffGroupIndex = 0;
673 <    
674 <    for(int i = 0; i < myNode; i++) {
677 <      beginAtomIndex += NumAtomsInProc[i];
678 <      beginRigidBodyIndex += NumRigidBodiesInProc[i];
679 <      beginCutoffGroupIndex += NumCutoffGroupsInProc[i];
680 <    }
681 <    
682 < #endif
683 <    
684 <    //rigidbody's index begins right after atom's
685 <    beginRigidBodyIndex += info->getNGlobalAtoms();
686 <    
687 <    for(mol = info->beginMolecule(mi); mol != NULL;
688 <        mol = info->nextMolecule(mi)) {
673 >
674 >    for(int i = 0; i < info->getNGlobalMolecules(); i++) {
675        
676 <      //local index(index in DataStorge) of atom is important
677 <      for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
678 <        atom->setGlobalIndex(beginAtomIndex++);
676 > #ifdef IS_MPI      
677 >      if (info->getMolToProc(i) == worldRank) {
678 > #endif        
679 >        // stuff to do if I own this molecule
680 >        mol = info->getMoleculeByGlobalIndex(i);
681 >
682 >        //local index(index in DataStorge) of atom is important
683 >        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
684 >          atom->setGlobalIndex(beginAtomIndex++);
685 >        }
686 >        
687 >        for(rb = mol->beginRigidBody(ri); rb != NULL;
688 >            rb = mol->nextRigidBody(ri)) {
689 >          rb->setGlobalIndex(beginRigidBodyIndex++);
690 >        }
691 >        
692 >        //local index of cutoff group is trivial, it only depends on
693 >        //the order of travesing
694 >        for(cg = mol->beginCutoffGroup(ci); cg != NULL;
695 >            cg = mol->nextCutoffGroup(ci)) {
696 >          cg->setGlobalIndex(beginCutoffGroupIndex++);
697 >        }        
698 >        
699 > #ifdef IS_MPI        
700 >      }  else {
701 >
702 >        // stuff to do if I don't own this molecule
703 >        
704 >        int stampId = info->getMoleculeStampId(i);
705 >        MoleculeStamp* stamp = info->getMoleculeStamp(stampId);
706 >
707 >        beginAtomIndex += stamp->getNAtoms();
708 >        beginRigidBodyIndex += stamp->getNRigidBodies();
709 >        beginCutoffGroupIndex += stamp->getNCutoffGroups() + stamp->getNFreeAtoms();
710        }
711 <      
712 <      for(rb = mol->beginRigidBody(ri); rb != NULL;
713 <          rb = mol->nextRigidBody(ri)) {
714 <        rb->setGlobalIndex(beginRigidBodyIndex++);
698 <      }
699 <      
700 <      //local index of cutoff group is trivial, it only depends on the order of travesing
701 <      for(cg = mol->beginCutoffGroup(ci); cg != NULL;
702 <          cg = mol->nextCutoffGroup(ci)) {
703 <        cg->setGlobalIndex(beginCutoffGroupIndex++);
704 <      }
705 <    }
706 <    
711 > #endif          
712 >
713 >    } //end for(int i=0)  
714 >
715      //fill globalGroupMembership
716      std::vector<int> globalGroupMembership(info->getNGlobalAtoms(), 0);
717      for(mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {        
# Line 715 | Line 723 | namespace oopse {
723          
724        }      
725      }
726 <    
726 >  
727   #ifdef IS_MPI    
728      // Since the globalGroupMembership has been zero filled and we've only
729      // poked values into the atoms we know, we can do an Allreduce
730      // to get the full globalGroupMembership array (We think).
731      // This would be prettier if we could use MPI_IN_PLACE like the MPI-2
732      // docs said we could.
733 <    std::vector<int> tmpGroupMembership(nGlobalAtoms, 0);
733 >    std::vector<int> tmpGroupMembership(info->getNGlobalAtoms(), 0);
734      MPI_Allreduce(&globalGroupMembership[0], &tmpGroupMembership[0], nGlobalAtoms,
735                    MPI_INT, MPI_SUM, MPI_COMM_WORLD);
736      info->setGlobalGroupMembership(tmpGroupMembership);
# Line 734 | Line 742 | namespace oopse {
742      std::vector<int> globalMolMembership(info->getNGlobalAtoms(), 0);
743      
744      for(mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
737      
745        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
746          globalMolMembership[atom->getGlobalIndex()] = mol->getGlobalIndex();
747        }
748      }
749      
750   #ifdef IS_MPI
751 <    std::vector<int> tmpMolMembership(nGlobalAtoms, 0);
751 >    std::vector<int> tmpMolMembership(info->getNGlobalAtoms(), 0);
752      
753      MPI_Allreduce(&globalMolMembership[0], &tmpMolMembership[0], nGlobalAtoms,
754                    MPI_INT, MPI_SUM, MPI_COMM_WORLD);
# Line 767 | Line 774 | namespace oopse {
774      std::vector<int> numIntegrableObjectsPerMol = nIOPerMol;
775   #endif    
776  
777 < std::vector<int> startingIOIndexForMol(info->getNGlobalMolecules());
778 <
779 < int startingIndex = 0;
780 < for (int i = 0; i < info->getNGlobalMolecules(); i++) {
781 <  startingIOIndexForMol[i] = startingIndex;
782 <  startingIndex += numIntegrableObjectsPerMol[i];
783 < }
784 <
785 < std::vector<StuntDouble*> IOIndexToIntegrableObject(info->getNGlobalIntegrableObjects(), (StuntDouble*)NULL);
786 < for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
777 >    std::vector<int> startingIOIndexForMol(info->getNGlobalMolecules());
778 >    
779 >    int startingIndex = 0;
780 >    for (int i = 0; i < info->getNGlobalMolecules(); i++) {
781 >      startingIOIndexForMol[i] = startingIndex;
782 >      startingIndex += numIntegrableObjectsPerMol[i];
783 >    }
784 >    
785 >    std::vector<StuntDouble*> IOIndexToIntegrableObject(info->getNGlobalIntegrableObjects(), (StuntDouble*)NULL);
786 >    for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
787        int myGlobalIndex = mol->getGlobalIndex();
788        int globalIO = startingIOIndexForMol[myGlobalIndex];
789        for (StuntDouble* integrableObject = mol->beginIntegrableObject(ioi); integrableObject != NULL;
790             integrableObject = mol->nextIntegrableObject(ioi)) {
791 <            integrableObject->setGlobalIntegrableObjectIndex(globalIO);
792 <            IOIndexToIntegrableObject[globalIO] = integrableObject;
793 <            globalIO++;
791 >        integrableObject->setGlobalIntegrableObjectIndex(globalIO);
792 >        IOIndexToIntegrableObject[globalIO] = integrableObject;
793 >        globalIO++;
794        }
795      }
796 <
797 <  info->setIOIndexToIntegrableObject(IOIndexToIntegrableObject);
798 <  
796 >      
797 >    info->setIOIndexToIntegrableObject(IOIndexToIntegrableObject);
798 >    
799    }
800    
801    void SimCreator::loadCoordinates(SimInfo* info, const std::string& mdFileName) {
802      Globals* simParams;
803 +
804      simParams = info->getSimParams();
805      
798    
806      DumpReader reader(info, mdFileName);
807      int nframes = reader.getNFrames();
808 <    
808 >
809      if (nframes > 0) {
810        reader.readFrame(nframes - 1);
811      } else {
# Line 809 | Line 816 | int startingIndex = 0;
816        painCave.isFatal = 1;
817        simError();
818      }
812    
819      //copy the current snapshot to previous snapshot
820      info->getSnapshotManager()->advance();
821    }
822    
823 < } //end namespace oopse
823 > } //end namespace OpenMD
824  
825  

Comparing:
trunk/src/brains/SimCreator.cpp (property svn:keywords), Revision 1287 by gezelter, Wed Sep 10 18:11:32 2008 UTC vs.
branches/development/src/brains/SimCreator.cpp (property svn:keywords), Revision 1627 by gezelter, Tue Sep 13 22:05:04 2011 UTC

# Line 0 | Line 1
1 + Author Id Revision Date

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines