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

Comparing trunk/src/brains/SimCreator.cpp (file contents):
Revision 1025 by gezelter, Wed Aug 30 20:33:44 2006 UTC vs.
Revision 1790 by gezelter, Thu Aug 30 17:18:22 2012 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]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40 + * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41   */
42  
43   /**
# Line 55 | Line 56
56   #include "brains/SimCreator.hpp"
57   #include "brains/SimSnapshotManager.hpp"
58   #include "io/DumpReader.hpp"
59 < #include "UseTheForce/ForceFieldFactory.hpp"
59 > #include "brains/ForceField.hpp"
60   #include "utils/simError.h"
61   #include "utils/StringUtils.hpp"
62   #include "math/SeqRandNumGen.hpp"
# Line 75 | Line 76
76   #include "antlr/NoViableAltForCharException.hpp"
77   #include "antlr/NoViableAltException.hpp"
78  
79 + #include "types/DirectionalAdapter.hpp"
80 + #include "types/MultipoleAdapter.hpp"
81 + #include "types/EAMAdapter.hpp"
82 + #include "types/SuttonChenAdapter.hpp"
83 + #include "types/PolarizableAdapter.hpp"
84 + #include "types/FixedChargeAdapter.hpp"
85 + #include "types/FluctuatingChargeAdapter.hpp"
86 +
87   #ifdef IS_MPI
88 + #include "mpi.h"
89   #include "math/ParallelRandNumGen.hpp"
90   #endif
91  
92 < namespace oopse {
92 > namespace OpenMD {
93    
94 <  Globals* SimCreator::parseFile(std::istream& rawMetaDataStream, const std::string& filename, int startOfMetaDataBlock ){
94 >  Globals* SimCreator::parseFile(std::istream& rawMetaDataStream, const std::string& filename, int mdFileVersion, int startOfMetaDataBlock ){
95      Globals* simParams = NULL;
96      try {
97  
# Line 92 | Line 102 | namespace oopse {
102        const int masterNode = 0;
103        int commStatus;
104        if (worldRank == masterNode) {
105 < #endif
106 <                
105 >        commStatus = MPI_Bcast(&mdFileVersion, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
106 > #endif                
107          SimplePreprocessor preprocessor;
108          preprocessor.preprocess(rawMetaDataStream, filename, startOfMetaDataBlock, ppStream);
109                  
# Line 106 | Line 116 | namespace oopse {
116              
117                  
118        } else {
119 +
120 +        commStatus = MPI_Bcast(&mdFileVersion, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
121 +
122          //get stream size
123          commStatus = MPI_Bcast(&streamSize, 1, MPI_LONG, masterNode, MPI_COMM_WORLD);  
124 <                
124 >
125          char* buf = new char[streamSize];
126          assert(buf);
127                  
# Line 116 | Line 129 | namespace oopse {
129          commStatus = MPI_Bcast(buf, streamSize, MPI_CHAR, masterNode, MPI_COMM_WORLD);
130                  
131          ppStream.str(buf);
132 <        delete buf;
132 >        delete [] buf;
133  
134        }
135   #endif            
# Line 145 | Line 158 | namespace oopse {
158        treeParser.initializeASTFactory(factory);
159        treeParser.setASTFactory(&factory);
160        simParams = treeParser.walkTree(parser.getAST());
148
161      }
162  
163        
# Line 215 | Line 227 | namespace oopse {
227        painCave.isFatal = 1;
228        simError();        
229      }
230 <    catch (OOPSEException& e) {
230 >    catch (OpenMDException& e) {
231        sprintf(painCave.errMsg,
232                "%s\n",
233                e.getMessage().c_str());
# Line 230 | Line 242 | namespace oopse {
242        simError();
243      }
244  
245 +    simParams->setMDfileVersion(mdFileVersion);
246      return simParams;
247    }
248    
249    SimInfo*  SimCreator::createSim(const std::string & mdFileName,
250                                    bool loadInitCoords) {
251 <
251 >    
252      const int bufferSize = 65535;
253      char buffer[bufferSize];
254      int lineNo = 0;
# Line 243 | Line 256 | namespace oopse {
256      int metaDataBlockStart = -1;
257      int metaDataBlockEnd = -1;
258      int i;
259 <    int mdOffset;
259 >    streamoff mdOffset(0);
260 >    int mdFileVersion;
261  
262 +
263   #ifdef IS_MPI            
264      const int masterNode = 0;
265      if (worldRank == masterNode) {
266   #endif
267  
268 <      std::ifstream mdFile_(mdFileName.c_str());
268 >      std::ifstream mdFile_;
269 >      mdFile_.open(mdFileName.c_str(), ifstream::in | ifstream::binary);
270        
271        if (mdFile_.fail()) {
272          sprintf(painCave.errMsg,
# Line 263 | Line 279 | namespace oopse {
279        mdFile_.getline(buffer, bufferSize);
280        ++lineNo;
281        std::string line = trimLeftCopy(buffer);
282 <      i = CaseInsensitiveFind(line, "<OOPSE");
283 <      if (i == string::npos) {
282 >      i = CaseInsensitiveFind(line, "<OpenMD");
283 >      if (static_cast<size_t>(i) == string::npos) {
284 >        // try the older file strings to see if that works:
285 >        i = CaseInsensitiveFind(line, "<OOPSE");
286 >      }
287 >      
288 >      if (static_cast<size_t>(i) == string::npos) {
289 >        // still no luck!
290          sprintf(painCave.errMsg,
291 <                "SimCreator: File: %s is not an OOPSE file!\n",
291 >                "SimCreator: File: %s is not a valid OpenMD file!\n",
292                  mdFileName.c_str());
293          painCave.isFatal = 1;
294          simError();
295        }
296 +      
297 +      // found the correct opening string, now try to get the file
298 +      // format version number.
299  
300 +      StringTokenizer tokenizer(line, "=<> \t\n\r");
301 +      std::string fileType = tokenizer.nextToken();
302 +      toUpper(fileType);
303 +
304 +      mdFileVersion = 0;
305 +
306 +      if (fileType == "OPENMD") {
307 +        while (tokenizer.hasMoreTokens()) {
308 +          std::string token(tokenizer.nextToken());
309 +          toUpper(token);
310 +          if (token == "VERSION") {
311 +            mdFileVersion = tokenizer.nextTokenAsInt();
312 +            break;
313 +          }
314 +        }
315 +      }
316 +            
317        //scan through the input stream and find MetaData tag        
318        while(mdFile_.getline(buffer, bufferSize)) {
319          ++lineNo;
# Line 327 | Line 369 | namespace oopse {
369      std::stringstream rawMetaDataStream(mdRawData);
370  
371      //parse meta-data file
372 <    Globals* simParams = parseFile(rawMetaDataStream, mdFileName, metaDataBlockStart+1);
372 >    Globals* simParams = parseFile(rawMetaDataStream, mdFileName, mdFileVersion,
373 >                                   metaDataBlockStart + 1);
374      
375      //create the force field
376 <    ForceField * ff = ForceFieldFactory::getInstance()
377 <      ->createForceField(simParams->getForceField());
335 <    
376 >    ForceField * ff = new ForceField(simParams->getForceField());
377 >
378      if (ff == NULL) {
379        sprintf(painCave.errMsg,
380                "ForceField Factory can not create %s force field\n",
# Line 365 | Line 407 | namespace oopse {
407      }
408      
409      ff->parse(forcefieldFileName);
368    ff->setFortranForceOptions();
410      //create SimInfo
411      SimInfo * info = new SimInfo(ff, simParams);
412  
# Line 383 | Line 424 | namespace oopse {
424      //create the molecules
425      createMolecules(info);
426      
427 <    
427 >    //find the storage layout
428 >
429 >    int storageLayout = computeStorageLayout(info);
430 >
431      //allocate memory for DataStorage(circular reference, need to
432      //break it)
433 <    info->setSnapshotManager(new SimSnapshotManager(info));
433 >    info->setSnapshotManager(new SimSnapshotManager(info, storageLayout));
434      
435      //set the global index of atoms, rigidbodies and cutoffgroups
436      //(only need to be set once, the global index will never change
# Line 395 | Line 439 | namespace oopse {
439      //responsibility to LocalIndexManager.
440      setGlobalIndex(info);
441      
442 <    //Although addExcludePairs is called inside SimInfo's addMolecule
442 >    //Although addInteractionPairs is called inside SimInfo's addMolecule
443      //method, at that point atoms don't have the global index yet
444      //(their global index are all initialized to -1).  Therefore we
445 <    //have to call addExcludePairs explicitly here. A way to work
445 >    //have to call addInteractionPairs explicitly here. A way to work
446      //around is that we can determine the beginning global indices of
447      //atoms before they get created.
448      SimInfo::MoleculeIterator mi;
449      Molecule* mol;
450      for (mol= info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
451 <      info->addExcludePairs(mol);
451 >      info->addInteractionPairs(mol);
452      }
453      
454      if (loadInitCoords)
455        loadCoordinates(info, mdFileName);    
412    
456      return info;
457    }
458    
# Line 473 | Line 516 | namespace oopse {
516                "\tthe number of molecules.  This will not result in a \n"
517                "\tusable division of atoms for force decomposition.\n"
518                "\tEither try a smaller number of processors, or run the\n"
519 <              "\tsingle-processor version of OOPSE.\n", nProcessors, nGlobalMols);
519 >              "\tsingle-processor version of OpenMD.\n", nProcessors, nGlobalMols);
520        
521        painCave.isFatal = 1;
522        simError();
# Line 591 | Line 634 | namespace oopse {
634      info->setMolToProcMap(molToProcMap);
635      sprintf(checkPointMsg,
636              "Successfully divided the molecules among the processors.\n");
637 <    MPIcheckPoint();
637 >    errorCheckPoint();
638    }
639    
640   #endif
# Line 608 | Line 651 | namespace oopse {
651   #endif
652          
653          stampId = info->getMoleculeStampId(i);
654 <        Molecule * mol = molCreator.createMolecule(info->getForceField(), info->getMoleculeStamp(stampId),
655 <                                                   stampId, i, info->getLocalIndexManager());
654 >        Molecule * mol = molCreator.createMolecule(info->getForceField(),
655 >                                                   info->getMoleculeStamp(stampId),
656 >                                                   stampId, i,
657 >                                                   info->getLocalIndexManager());
658          
659          info->addMolecule(mol);
660          
# Line 621 | Line 666 | namespace oopse {
666        
667      } //end for(int i=0)  
668    }
669 +    
670 +  int SimCreator::computeStorageLayout(SimInfo* info) {
671 +
672 +    Globals* simParams = info->getSimParams();
673 +    int nRigidBodies = info->getNGlobalRigidBodies();
674 +    set<AtomType*> atomTypes = info->getSimulatedAtomTypes();
675 +    set<AtomType*>::iterator i;
676 +    bool hasDirectionalAtoms = false;
677 +    bool hasFixedCharge = false;
678 +    bool hasMultipoles = false;    
679 +    bool hasPolarizable = false;    
680 +    bool hasFluctuatingCharge = false;    
681 +    bool hasMetallic = false;
682 +    int storageLayout = 0;
683 +    storageLayout |= DataStorage::dslPosition;
684 +    storageLayout |= DataStorage::dslVelocity;
685 +    storageLayout |= DataStorage::dslForce;
686 +
687 +    for (i = atomTypes.begin(); i != atomTypes.end(); ++i) {
688 +
689 +      DirectionalAdapter da = DirectionalAdapter( (*i) );
690 +      MultipoleAdapter ma = MultipoleAdapter( (*i) );
691 +      EAMAdapter ea = EAMAdapter( (*i) );
692 +      SuttonChenAdapter sca = SuttonChenAdapter( (*i) );
693 +      PolarizableAdapter pa = PolarizableAdapter( (*i) );
694 +      FixedChargeAdapter fca = FixedChargeAdapter( (*i) );
695 +      FluctuatingChargeAdapter fqa = FluctuatingChargeAdapter( (*i) );
696 +
697 +      if (da.isDirectional()){
698 +        hasDirectionalAtoms = true;
699 +      }
700 +      if (ma.isMultipole()){
701 +        hasMultipoles = true;
702 +      }
703 +      if (ea.isEAM() || sca.isSuttonChen()){
704 +        hasMetallic = true;
705 +      }
706 +      if ( fca.isFixedCharge() ){
707 +        hasFixedCharge = true;
708 +      }
709 +      if ( fqa.isFluctuatingCharge() ){
710 +        hasFluctuatingCharge = true;
711 +      }
712 +      if ( pa.isPolarizable() ){
713 +        hasPolarizable = true;
714 +      }
715 +    }
716      
717 +    if (nRigidBodies > 0 || hasDirectionalAtoms) {
718 +      storageLayout |= DataStorage::dslAmat;
719 +      if(storageLayout & DataStorage::dslVelocity) {
720 +        storageLayout |= DataStorage::dslAngularMomentum;
721 +      }
722 +      if (storageLayout & DataStorage::dslForce) {
723 +        storageLayout |= DataStorage::dslTorque;
724 +      }
725 +    }
726 +    if (hasMultipoles) {
727 +      storageLayout |= DataStorage::dslElectroFrame;
728 +    }
729 +    if (hasFixedCharge || hasFluctuatingCharge) {
730 +      storageLayout |= DataStorage::dslSkippedCharge;
731 +    }
732 +    if (hasMetallic) {
733 +      storageLayout |= DataStorage::dslDensity;
734 +      storageLayout |= DataStorage::dslFunctional;
735 +      storageLayout |= DataStorage::dslFunctionalDerivative;
736 +    }
737 +    if (hasPolarizable) {
738 +      storageLayout |= DataStorage::dslElectricField;
739 +    }
740 +    if (hasFluctuatingCharge){
741 +      storageLayout |= DataStorage::dslFlucQPosition;
742 +      if(storageLayout & DataStorage::dslVelocity) {
743 +        storageLayout |= DataStorage::dslFlucQVelocity;
744 +      }
745 +      if (storageLayout & DataStorage::dslForce) {
746 +        storageLayout |= DataStorage::dslFlucQForce;
747 +      }
748 +    }
749 +    
750 +    // if the user has asked for them, make sure we've got the memory for the
751 +    // objects defined.
752 +
753 +    if (simParams->getOutputParticlePotential()) {
754 +      storageLayout |= DataStorage::dslParticlePot;
755 +    }
756 +
757 +    if (simParams->havePrintHeatFlux()) {
758 +      if (simParams->getPrintHeatFlux()) {
759 +        storageLayout |= DataStorage::dslParticlePot;
760 +      }
761 +    }
762 +
763 +    if (simParams->getOutputElectricField()) {
764 +      storageLayout |= DataStorage::dslElectricField;
765 +    }
766 +    if (simParams->getOutputFluctuatingCharges()) {
767 +      storageLayout |= DataStorage::dslFlucQPosition;
768 +      storageLayout |= DataStorage::dslFlucQVelocity;
769 +      storageLayout |= DataStorage::dslFlucQForce;
770 +    }
771 +
772 +    return storageLayout;
773 +  }
774 +
775    void SimCreator::setGlobalIndex(SimInfo *info) {
776      SimInfo::MoleculeIterator mi;
777      Molecule::AtomIterator ai;
# Line 636 | Line 786 | namespace oopse {
786      int beginRigidBodyIndex;
787      int beginCutoffGroupIndex;
788      int nGlobalAtoms = info->getNGlobalAtoms();
639
640    /**@todo fixme */
641 #ifndef IS_MPI
789      
790      beginAtomIndex = 0;
644    beginRigidBodyIndex = 0;
645    beginCutoffGroupIndex = 0;
646    
647 #else
648    
649    int nproc;
650    int myNode;
651    
652    myNode = worldRank;
653    MPI_Comm_size(MPI_COMM_WORLD, &nproc);
654    
655    std::vector < int > tmpAtomsInProc(nproc, 0);
656    std::vector < int > tmpRigidBodiesInProc(nproc, 0);
657    std::vector < int > tmpCutoffGroupsInProc(nproc, 0);
658    std::vector < int > NumAtomsInProc(nproc, 0);
659    std::vector < int > NumRigidBodiesInProc(nproc, 0);
660    std::vector < int > NumCutoffGroupsInProc(nproc, 0);
661    
662    tmpAtomsInProc[myNode] = info->getNAtoms();
663    tmpRigidBodiesInProc[myNode] = info->getNRigidBodies();
664    tmpCutoffGroupsInProc[myNode] = info->getNCutoffGroups();
665    
666    //do MPI_ALLREDUCE to exchange the total number of atoms, rigidbodies and cutoff groups
667    MPI_Allreduce(&tmpAtomsInProc[0], &NumAtomsInProc[0], nproc, MPI_INT,
668                  MPI_SUM, MPI_COMM_WORLD);
669    MPI_Allreduce(&tmpRigidBodiesInProc[0], &NumRigidBodiesInProc[0], nproc,
670                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
671    MPI_Allreduce(&tmpCutoffGroupsInProc[0], &NumCutoffGroupsInProc[0], nproc,
672                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
673    
674    beginAtomIndex = 0;
675    beginRigidBodyIndex = 0;
676    beginCutoffGroupIndex = 0;
677    
678    for(int i = 0; i < myNode; i++) {
679      beginAtomIndex += NumAtomsInProc[i];
680      beginRigidBodyIndex += NumRigidBodiesInProc[i];
681      beginCutoffGroupIndex += NumCutoffGroupsInProc[i];
682    }
683    
684 #endif
685    
791      //rigidbody's index begins right after atom's
792 <    beginRigidBodyIndex += info->getNGlobalAtoms();
793 <    
794 <    for(mol = info->beginMolecule(mi); mol != NULL;
795 <        mol = info->nextMolecule(mi)) {
792 >    beginRigidBodyIndex = info->getNGlobalAtoms();
793 >    beginCutoffGroupIndex = 0;
794 >
795 >    for(int i = 0; i < info->getNGlobalMolecules(); i++) {
796        
797 <      //local index(index in DataStorge) of atom is important
798 <      for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
799 <        atom->setGlobalIndex(beginAtomIndex++);
797 > #ifdef IS_MPI      
798 >      if (info->getMolToProc(i) == worldRank) {
799 > #endif        
800 >        // stuff to do if I own this molecule
801 >        mol = info->getMoleculeByGlobalIndex(i);
802 >
803 >        //local index(index in DataStorge) of atom is important
804 >        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
805 >          atom->setGlobalIndex(beginAtomIndex++);
806 >        }
807 >        
808 >        for(rb = mol->beginRigidBody(ri); rb != NULL;
809 >            rb = mol->nextRigidBody(ri)) {
810 >          rb->setGlobalIndex(beginRigidBodyIndex++);
811 >        }
812 >        
813 >        //local index of cutoff group is trivial, it only depends on
814 >        //the order of travesing
815 >        for(cg = mol->beginCutoffGroup(ci); cg != NULL;
816 >            cg = mol->nextCutoffGroup(ci)) {
817 >          cg->setGlobalIndex(beginCutoffGroupIndex++);
818 >        }        
819 >        
820 > #ifdef IS_MPI        
821 >      }  else {
822 >
823 >        // stuff to do if I don't own this molecule
824 >        
825 >        int stampId = info->getMoleculeStampId(i);
826 >        MoleculeStamp* stamp = info->getMoleculeStamp(stampId);
827 >
828 >        beginAtomIndex += stamp->getNAtoms();
829 >        beginRigidBodyIndex += stamp->getNRigidBodies();
830 >        beginCutoffGroupIndex += stamp->getNCutoffGroups() + stamp->getNFreeAtoms();
831        }
832 <      
833 <      for(rb = mol->beginRigidBody(ri); rb != NULL;
834 <          rb = mol->nextRigidBody(ri)) {
835 <        rb->setGlobalIndex(beginRigidBodyIndex++);
700 <      }
701 <      
702 <      //local index of cutoff group is trivial, it only depends on the order of travesing
703 <      for(cg = mol->beginCutoffGroup(ci); cg != NULL;
704 <          cg = mol->nextCutoffGroup(ci)) {
705 <        cg->setGlobalIndex(beginCutoffGroupIndex++);
706 <      }
707 <    }
708 <    
832 > #endif          
833 >
834 >    } //end for(int i=0)  
835 >
836      //fill globalGroupMembership
837      std::vector<int> globalGroupMembership(info->getNGlobalAtoms(), 0);
838      for(mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {        
# Line 717 | Line 844 | namespace oopse {
844          
845        }      
846      }
847 <    
847 >  
848   #ifdef IS_MPI    
849      // Since the globalGroupMembership has been zero filled and we've only
850      // poked values into the atoms we know, we can do an Allreduce
851      // to get the full globalGroupMembership array (We think).
852      // This would be prettier if we could use MPI_IN_PLACE like the MPI-2
853      // docs said we could.
854 <    std::vector<int> tmpGroupMembership(nGlobalAtoms, 0);
854 >    std::vector<int> tmpGroupMembership(info->getNGlobalAtoms(), 0);
855      MPI_Allreduce(&globalGroupMembership[0], &tmpGroupMembership[0], nGlobalAtoms,
856                    MPI_INT, MPI_SUM, MPI_COMM_WORLD);
857      info->setGlobalGroupMembership(tmpGroupMembership);
# Line 736 | Line 863 | namespace oopse {
863      std::vector<int> globalMolMembership(info->getNGlobalAtoms(), 0);
864      
865      for(mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
739      
866        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
867          globalMolMembership[atom->getGlobalIndex()] = mol->getGlobalIndex();
868        }
869      }
870      
871   #ifdef IS_MPI
872 <    std::vector<int> tmpMolMembership(nGlobalAtoms, 0);
872 >    std::vector<int> tmpMolMembership(info->getNGlobalAtoms(), 0);
873      
874      MPI_Allreduce(&globalMolMembership[0], &tmpMolMembership[0], nGlobalAtoms,
875                    MPI_INT, MPI_SUM, MPI_COMM_WORLD);
# Line 769 | Line 895 | namespace oopse {
895      std::vector<int> numIntegrableObjectsPerMol = nIOPerMol;
896   #endif    
897  
898 < std::vector<int> startingIOIndexForMol(info->getNGlobalMolecules());
899 <
900 < int startingIndex = 0;
901 < for (int i = 0; i < info->getNGlobalMolecules(); i++) {
902 <  startingIOIndexForMol[i] = startingIndex;
903 <  startingIndex += numIntegrableObjectsPerMol[i];
904 < }
905 <
906 < std::vector<StuntDouble*> IOIndexToIntegrableObject(info->getNGlobalIntegrableObjects(), (StuntDouble*)NULL);
907 < for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
898 >    std::vector<int> startingIOIndexForMol(info->getNGlobalMolecules());
899 >    
900 >    int startingIndex = 0;
901 >    for (int i = 0; i < info->getNGlobalMolecules(); i++) {
902 >      startingIOIndexForMol[i] = startingIndex;
903 >      startingIndex += numIntegrableObjectsPerMol[i];
904 >    }
905 >    
906 >    std::vector<StuntDouble*> IOIndexToIntegrableObject(info->getNGlobalIntegrableObjects(), (StuntDouble*)NULL);
907 >    for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
908        int myGlobalIndex = mol->getGlobalIndex();
909        int globalIO = startingIOIndexForMol[myGlobalIndex];
910 <      for (StuntDouble* integrableObject = mol->beginIntegrableObject(ioi); integrableObject != NULL;
911 <           integrableObject = mol->nextIntegrableObject(ioi)) {
912 <            integrableObject->setGlobalIntegrableObjectIndex(globalIO);
913 <            IOIndexToIntegrableObject[globalIO] = integrableObject;
914 <            globalIO++;
910 >      for (StuntDouble* sd = mol->beginIntegrableObject(ioi); sd != NULL;
911 >           sd = mol->nextIntegrableObject(ioi)) {
912 >        sd->setGlobalIntegrableObjectIndex(globalIO);
913 >        IOIndexToIntegrableObject[globalIO] = sd;
914 >        globalIO++;
915        }
916      }
917 <
918 <  info->setIOIndexToIntegrableObject(IOIndexToIntegrableObject);
919 <  
917 >      
918 >    info->setIOIndexToIntegrableObject(IOIndexToIntegrableObject);
919 >    
920    }
921    
922    void SimCreator::loadCoordinates(SimInfo* info, const std::string& mdFileName) {
923      Globals* simParams;
924 +
925      simParams = info->getSimParams();
926      
800    
927      DumpReader reader(info, mdFileName);
928      int nframes = reader.getNFrames();
929 <    
929 >
930      if (nframes > 0) {
931        reader.readFrame(nframes - 1);
932      } else {
# Line 811 | Line 937 | int startingIndex = 0;
937        painCave.isFatal = 1;
938        simError();
939      }
814    
940      //copy the current snapshot to previous snapshot
941      info->getSnapshotManager()->advance();
942    }
943    
944 < } //end namespace oopse
944 > } //end namespace OpenMD
945  
946  

Comparing trunk/src/brains/SimCreator.cpp (property svn:keywords):
Revision 1025 by gezelter, Wed Aug 30 20:33:44 2006 UTC vs.
Revision 1790 by gezelter, Thu Aug 30 17:18:22 2012 UTC

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines