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

Comparing branches/development/src/brains/SimCreator.cpp (file contents):
Revision 1593 by gezelter, Fri Jul 15 21:35:14 2011 UTC vs.
Revision 1808 by gezelter, Mon Oct 22 20:42:10 2012 UTC

# Line 36 | Line 36
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).                        
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   /**
44   * @file SimCreator.cpp
45   * @author tlin
46   * @date 11/03/2004
46 – * @time 13:51am
47   * @version 1.0
48   */
49   #include <exception>
# Line 55 | Line 55
55   #include "brains/SimCreator.hpp"
56   #include "brains/SimSnapshotManager.hpp"
57   #include "io/DumpReader.hpp"
58 < #include "UseTheForce/ForceFieldFactory.hpp"
58 > #include "brains/ForceField.hpp"
59   #include "utils/simError.h"
60   #include "utils/StringUtils.hpp"
61   #include "math/SeqRandNumGen.hpp"
# Line 75 | Line 75
75   #include "antlr/NoViableAltForCharException.hpp"
76   #include "antlr/NoViableAltException.hpp"
77  
78 + #include "types/DirectionalAdapter.hpp"
79 + #include "types/MultipoleAdapter.hpp"
80 + #include "types/EAMAdapter.hpp"
81 + #include "types/SuttonChenAdapter.hpp"
82 + #include "types/PolarizableAdapter.hpp"
83 + #include "types/FixedChargeAdapter.hpp"
84 + #include "types/FluctuatingChargeAdapter.hpp"
85 +
86   #ifdef IS_MPI
87 + #include "mpi.h"
88   #include "math/ParallelRandNumGen.hpp"
89   #endif
90  
91   namespace OpenMD {
92    
93 <  Globals* SimCreator::parseFile(std::istream& rawMetaDataStream, const std::string& filename, int startOfMetaDataBlock ){
93 >  Globals* SimCreator::parseFile(std::istream& rawMetaDataStream, const std::string& filename, int mdFileVersion, int startOfMetaDataBlock ){
94      Globals* simParams = NULL;
95      try {
96  
# Line 92 | Line 101 | namespace OpenMD {
101        const int masterNode = 0;
102        int commStatus;
103        if (worldRank == masterNode) {
104 < #endif
105 <                
104 >        commStatus = MPI_Bcast(&mdFileVersion, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
105 > #endif                
106          SimplePreprocessor preprocessor;
107          preprocessor.preprocess(rawMetaDataStream, filename, startOfMetaDataBlock, ppStream);
108                  
# Line 106 | Line 115 | namespace OpenMD {
115              
116                  
117        } else {
118 +
119 +        commStatus = MPI_Bcast(&mdFileVersion, 1, MPI_INT, masterNode, MPI_COMM_WORLD);
120 +
121          //get stream size
122          commStatus = MPI_Bcast(&streamSize, 1, MPI_LONG, masterNode, MPI_COMM_WORLD);  
123  
# Line 229 | Line 241 | namespace OpenMD {
241        simError();
242      }
243  
244 +    simParams->setMDfileVersion(mdFileVersion);
245      return simParams;
246    }
247    
# Line 242 | Line 255 | namespace OpenMD {
255      int metaDataBlockStart = -1;
256      int metaDataBlockEnd = -1;
257      int i;
258 <    int mdOffset;
258 >    streamoff mdOffset;
259 >    int mdFileVersion;
260  
261 +
262   #ifdef IS_MPI            
263      const int masterNode = 0;
264      if (worldRank == masterNode) {
265   #endif
266  
267 <      std::ifstream mdFile_(mdFileName.c_str());
267 >      std::ifstream mdFile_;
268 >      mdFile_.open(mdFileName.c_str(), ifstream::in | ifstream::binary);
269        
270        if (mdFile_.fail()) {
271          sprintf(painCave.errMsg,
# Line 276 | Line 292 | namespace OpenMD {
292          painCave.isFatal = 1;
293          simError();
294        }
295 +      
296 +      // found the correct opening string, now try to get the file
297 +      // format version number.
298  
299 +      StringTokenizer tokenizer(line, "=<> \t\n\r");
300 +      std::string fileType = tokenizer.nextToken();
301 +      toUpper(fileType);
302 +
303 +      mdFileVersion = 0;
304 +
305 +      if (fileType == "OPENMD") {
306 +        while (tokenizer.hasMoreTokens()) {
307 +          std::string token(tokenizer.nextToken());
308 +          toUpper(token);
309 +          if (token == "VERSION") {
310 +            mdFileVersion = tokenizer.nextTokenAsInt();
311 +            break;
312 +          }
313 +        }
314 +      }
315 +            
316        //scan through the input stream and find MetaData tag        
317        while(mdFile_.getline(buffer, bufferSize)) {
318          ++lineNo;
# Line 332 | Line 368 | namespace OpenMD {
368      std::stringstream rawMetaDataStream(mdRawData);
369  
370      //parse meta-data file
371 <    Globals* simParams = parseFile(rawMetaDataStream, mdFileName, metaDataBlockStart+1);
371 >    Globals* simParams = parseFile(rawMetaDataStream, mdFileName, mdFileVersion,
372 >                                   metaDataBlockStart + 1);
373      
374      //create the force field
375 <    ForceField * ff = ForceFieldFactory::getInstance()->createForceField(simParams->getForceField());
375 >    ForceField * ff = new ForceField(simParams->getForceField());
376  
377      if (ff == NULL) {
378        sprintf(painCave.errMsg,
# Line 386 | Line 423 | namespace OpenMD {
423      //create the molecules
424      createMolecules(info);
425      
426 +    //find the storage layout
427 +
428 +    int storageLayout = computeStorageLayout(info);
429 +
430      //allocate memory for DataStorage(circular reference, need to
431      //break it)
432 <    info->setSnapshotManager(new SimSnapshotManager(info));
432 >    info->setSnapshotManager(new SimSnapshotManager(info, storageLayout));
433      
434      //set the global index of atoms, rigidbodies and cutoffgroups
435      //(only need to be set once, the global index will never change
# Line 465 | Line 506 | namespace OpenMD {
506      int nGlobalMols = info->getNGlobalMolecules();
507      std::vector<int> molToProcMap(nGlobalMols, -1); // default to an error condition:
508      
509 <    MPI_Comm_size(MPI_COMM_WORLD, &nProcessors);
509 >    nProcessors = MPI::COMM_WORLD.Get_size();
510      
511      if (nProcessors > nGlobalMols) {
512        sprintf(painCave.errMsg,
# Line 503 | Line 544 | namespace OpenMD {
544        nTarget = (int)(precast + 0.5);
545        
546        for(i = 0; i < nGlobalMols; i++) {
547 +
548          done = 0;
549          loops = 0;
550          
# Line 527 | Line 569 | namespace OpenMD {
569            // and be done with it.
570            
571            if (loops > 100) {
572 +
573              sprintf(painCave.errMsg,
574 <                    "I've tried 100 times to assign molecule %d to a "
575 <                    " processor, but can't find a good spot.\n"
576 <                    "I'm assigning it at random to processor %d.\n",
577 <                    i, which_proc);
578 <            
574 >                    "There have been 100 attempts to assign molecule %d to an\n"
575 >                    "\tunderworked processor, but there's no good place to\n"
576 >                    "\tleave it.  OpenMD is assigning it at random to processor %d.\n",
577 >                    i, which_proc);
578 >          
579              painCave.isFatal = 0;
580 +            painCave.severity = OPENMD_INFO;
581              simError();
582              
583              molToProcMap[i] = which_proc;
# Line 578 | Line 622 | namespace OpenMD {
622        }
623        
624        delete myRandom;
625 <      
625 >
626        // Spray out this nonsense to all other processors:
627 <      
584 <      MPI_Bcast(&molToProcMap[0], nGlobalMols, MPI_INT, 0, MPI_COMM_WORLD);
627 >      MPI::COMM_WORLD.Bcast(&molToProcMap[0], nGlobalMols, MPI::INT, 0);
628      } else {
629        
630        // Listen to your marching orders from processor 0:
631 <      
632 <      MPI_Bcast(&molToProcMap[0], nGlobalMols, MPI_INT, 0, MPI_COMM_WORLD);
631 >      MPI::COMM_WORLD.Bcast(&molToProcMap[0], nGlobalMols, MPI::INT, 0);
632 >
633      }
634      
635      info->setMolToProcMap(molToProcMap);
# Line 609 | Line 652 | namespace OpenMD {
652   #endif
653          
654          stampId = info->getMoleculeStampId(i);
655 <        Molecule * mol = molCreator.createMolecule(info->getForceField(), info->getMoleculeStamp(stampId),
656 <                                                   stampId, i, info->getLocalIndexManager());
655 >        Molecule * mol = molCreator.createMolecule(info->getForceField(),
656 >                                                   info->getMoleculeStamp(stampId),
657 >                                                   stampId, i,
658 >                                                   info->getLocalIndexManager());
659          
660          info->addMolecule(mol);
661          
# Line 622 | Line 667 | namespace OpenMD {
667        
668      } //end for(int i=0)  
669    }
670 +    
671 +  int SimCreator::computeStorageLayout(SimInfo* info) {
672 +
673 +    Globals* simParams = info->getSimParams();
674 +    int nRigidBodies = info->getNGlobalRigidBodies();
675 +    set<AtomType*> atomTypes = info->getSimulatedAtomTypes();
676 +    set<AtomType*>::iterator i;
677 +    bool hasDirectionalAtoms = false;
678 +    bool hasFixedCharge = false;
679 +    bool hasDipoles = false;    
680 +    bool hasQuadrupoles = false;    
681 +    bool hasPolarizable = false;    
682 +    bool hasFluctuatingCharge = false;    
683 +    bool hasMetallic = false;
684 +    int storageLayout = 0;
685 +    storageLayout |= DataStorage::dslPosition;
686 +    storageLayout |= DataStorage::dslVelocity;
687 +    storageLayout |= DataStorage::dslForce;
688 +
689 +    for (i = atomTypes.begin(); i != atomTypes.end(); ++i) {
690 +
691 +      DirectionalAdapter da = DirectionalAdapter( (*i) );
692 +      MultipoleAdapter ma = MultipoleAdapter( (*i) );
693 +      EAMAdapter ea = EAMAdapter( (*i) );
694 +      SuttonChenAdapter sca = SuttonChenAdapter( (*i) );
695 +      PolarizableAdapter pa = PolarizableAdapter( (*i) );
696 +      FixedChargeAdapter fca = FixedChargeAdapter( (*i) );
697 +      FluctuatingChargeAdapter fqa = FluctuatingChargeAdapter( (*i) );
698 +
699 +      if (da.isDirectional()){
700 +        hasDirectionalAtoms = true;
701 +      }
702 +      if (ma.isDipole()){
703 +        hasDipoles = true;
704 +      }
705 +      if (ma.isQuadrupole()){
706 +        hasQuadrupoles = true;
707 +      }
708 +      if (ea.isEAM() || sca.isSuttonChen()){
709 +        hasMetallic = true;
710 +      }
711 +      if ( fca.isFixedCharge() ){
712 +        hasFixedCharge = true;
713 +      }
714 +      if ( fqa.isFluctuatingCharge() ){
715 +        hasFluctuatingCharge = true;
716 +      }
717 +      if ( pa.isPolarizable() ){
718 +        hasPolarizable = true;
719 +      }
720 +    }
721      
722 +    if (nRigidBodies > 0 || hasDirectionalAtoms) {
723 +      storageLayout |= DataStorage::dslAmat;
724 +      if(storageLayout & DataStorage::dslVelocity) {
725 +        storageLayout |= DataStorage::dslAngularMomentum;
726 +      }
727 +      if (storageLayout & DataStorage::dslForce) {
728 +        storageLayout |= DataStorage::dslTorque;
729 +      }
730 +    }
731 +    if (hasDipoles) {
732 +      storageLayout |= DataStorage::dslDipole;
733 +    }
734 +    if (hasQuadrupoles) {
735 +      storageLayout |= DataStorage::dslQuadrupole;
736 +    }
737 +    if (hasFixedCharge || hasFluctuatingCharge) {
738 +      storageLayout |= DataStorage::dslSkippedCharge;
739 +    }
740 +    if (hasMetallic) {
741 +      storageLayout |= DataStorage::dslDensity;
742 +      storageLayout |= DataStorage::dslFunctional;
743 +      storageLayout |= DataStorage::dslFunctionalDerivative;
744 +    }
745 +    if (hasPolarizable) {
746 +      storageLayout |= DataStorage::dslElectricField;
747 +    }
748 +    if (hasFluctuatingCharge){
749 +      storageLayout |= DataStorage::dslFlucQPosition;
750 +      if(storageLayout & DataStorage::dslVelocity) {
751 +        storageLayout |= DataStorage::dslFlucQVelocity;
752 +      }
753 +      if (storageLayout & DataStorage::dslForce) {
754 +        storageLayout |= DataStorage::dslFlucQForce;
755 +      }
756 +    }
757 +    
758 +    // if the user has asked for them, make sure we've got the memory for the
759 +    // objects defined.
760 +
761 +    if (simParams->getOutputParticlePotential()) {
762 +      storageLayout |= DataStorage::dslParticlePot;
763 +    }
764 +
765 +    if (simParams->havePrintHeatFlux()) {
766 +      if (simParams->getPrintHeatFlux()) {
767 +        storageLayout |= DataStorage::dslParticlePot;
768 +      }
769 +    }
770 +
771 +    if (simParams->getOutputElectricField()) {
772 +      storageLayout |= DataStorage::dslElectricField;
773 +    }
774 +
775 +    if (simParams->getOutputFluctuatingCharges()) {
776 +      storageLayout |= DataStorage::dslFlucQPosition;
777 +      storageLayout |= DataStorage::dslFlucQVelocity;
778 +      storageLayout |= DataStorage::dslFlucQForce;
779 +    }
780 +
781 +    return storageLayout;
782 +  }
783 +
784    void SimCreator::setGlobalIndex(SimInfo *info) {
785      SimInfo::MoleculeIterator mi;
786      Molecule::AtomIterator ai;
# Line 637 | Line 795 | namespace OpenMD {
795      int beginRigidBodyIndex;
796      int beginCutoffGroupIndex;
797      int nGlobalAtoms = info->getNGlobalAtoms();
798 <
641 <    /**@todo fixme */
642 < #ifndef IS_MPI
798 >    int nGlobalRigidBodies = info->getNGlobalRigidBodies();
799      
800      beginAtomIndex = 0;
645 –    beginRigidBodyIndex = 0;
646 –    beginCutoffGroupIndex = 0;
647 –    
648 – #else
649 –    
650 –    int nproc;
651 –    int myNode;
652 –    
653 –    myNode = worldRank;
654 –    MPI_Comm_size(MPI_COMM_WORLD, &nproc);
655 –    
656 –    std::vector < int > tmpAtomsInProc(nproc, 0);
657 –    std::vector < int > tmpRigidBodiesInProc(nproc, 0);
658 –    std::vector < int > tmpCutoffGroupsInProc(nproc, 0);
659 –    std::vector < int > NumAtomsInProc(nproc, 0);
660 –    std::vector < int > NumRigidBodiesInProc(nproc, 0);
661 –    std::vector < int > NumCutoffGroupsInProc(nproc, 0);
662 –    
663 –    tmpAtomsInProc[myNode] = info->getNAtoms();
664 –    tmpRigidBodiesInProc[myNode] = info->getNRigidBodies();
665 –    tmpCutoffGroupsInProc[myNode] = info->getNCutoffGroups();
666 –    
667 –    //do MPI_ALLREDUCE to exchange the total number of atoms, rigidbodies and cutoff groups
668 –    MPI_Allreduce(&tmpAtomsInProc[0], &NumAtomsInProc[0], nproc, MPI_INT,
669 –                  MPI_SUM, MPI_COMM_WORLD);
670 –    MPI_Allreduce(&tmpRigidBodiesInProc[0], &NumRigidBodiesInProc[0], nproc,
671 –                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
672 –    MPI_Allreduce(&tmpCutoffGroupsInProc[0], &NumCutoffGroupsInProc[0], nproc,
673 –                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
674 –    
675 –    beginAtomIndex = 0;
676 –    beginRigidBodyIndex = 0;
677 –    beginCutoffGroupIndex = 0;
678 –    
679 –    for(int i = 0; i < myNode; i++) {
680 –      beginAtomIndex += NumAtomsInProc[i];
681 –      beginRigidBodyIndex += NumRigidBodiesInProc[i];
682 –      beginCutoffGroupIndex += NumCutoffGroupsInProc[i];
683 –    }
684 –    
685 – #endif
686 –    
801      //rigidbody's index begins right after atom's
802 <    beginRigidBodyIndex += info->getNGlobalAtoms();
803 <    
804 <    for(mol = info->beginMolecule(mi); mol != NULL;
805 <        mol = info->nextMolecule(mi)) {
802 >    beginRigidBodyIndex = info->getNGlobalAtoms();
803 >    beginCutoffGroupIndex = 0;
804 >
805 >    for(int i = 0; i < info->getNGlobalMolecules(); i++) {
806        
807 <      //local index(index in DataStorge) of atom is important
808 <      for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
809 <        atom->setGlobalIndex(beginAtomIndex++);
807 > #ifdef IS_MPI      
808 >      if (info->getMolToProc(i) == worldRank) {
809 > #endif        
810 >        // stuff to do if I own this molecule
811 >        mol = info->getMoleculeByGlobalIndex(i);
812 >
813 >        //local index(index in DataStorge) of atom is important
814 >        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
815 >          atom->setGlobalIndex(beginAtomIndex++);
816 >        }
817 >        
818 >        for(rb = mol->beginRigidBody(ri); rb != NULL;
819 >            rb = mol->nextRigidBody(ri)) {
820 >          rb->setGlobalIndex(beginRigidBodyIndex++);
821 >        }
822 >        
823 >        //local index of cutoff group is trivial, it only depends on
824 >        //the order of travesing
825 >        for(cg = mol->beginCutoffGroup(ci); cg != NULL;
826 >            cg = mol->nextCutoffGroup(ci)) {
827 >          cg->setGlobalIndex(beginCutoffGroupIndex++);
828 >        }        
829 >        
830 > #ifdef IS_MPI        
831 >      }  else {
832 >
833 >        // stuff to do if I don't own this molecule
834 >        
835 >        int stampId = info->getMoleculeStampId(i);
836 >        MoleculeStamp* stamp = info->getMoleculeStamp(stampId);
837 >
838 >        beginAtomIndex += stamp->getNAtoms();
839 >        beginRigidBodyIndex += stamp->getNRigidBodies();
840 >        beginCutoffGroupIndex += stamp->getNCutoffGroups() + stamp->getNFreeAtoms();
841        }
842 <      
843 <      for(rb = mol->beginRigidBody(ri); rb != NULL;
844 <          rb = mol->nextRigidBody(ri)) {
845 <        rb->setGlobalIndex(beginRigidBodyIndex++);
701 <      }
702 <      
703 <      //local index of cutoff group is trivial, it only depends on the order of travesing
704 <      for(cg = mol->beginCutoffGroup(ci); cg != NULL;
705 <          cg = mol->nextCutoffGroup(ci)) {
706 <        cg->setGlobalIndex(beginCutoffGroupIndex++);
707 <      }
708 <    }
709 <    
842 > #endif          
843 >
844 >    } //end for(int i=0)  
845 >
846      //fill globalGroupMembership
847 <    std::vector<int> globalGroupMembership(info->getNGlobalAtoms(), -1);
847 >    std::vector<int> globalGroupMembership(info->getNGlobalAtoms(), 0);
848      for(mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {        
849        for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci)) {
850          
# Line 726 | Line 862 | namespace OpenMD {
862      // This would be prettier if we could use MPI_IN_PLACE like the MPI-2
863      // docs said we could.
864      std::vector<int> tmpGroupMembership(info->getNGlobalAtoms(), 0);
865 <    MPI_Allreduce(&globalGroupMembership[0], &tmpGroupMembership[0], nGlobalAtoms,
866 <                  MPI_INT, MPI_MAX, MPI_COMM_WORLD);
865 >    MPI::COMM_WORLD.Allreduce(&globalGroupMembership[0],
866 >                              &tmpGroupMembership[0], nGlobalAtoms,
867 >                              MPI::INT, MPI::SUM);
868      info->setGlobalGroupMembership(tmpGroupMembership);
732 –
733 –    cerr << "ggm:\n";
734 –    for (int i = 0; i < tmpGroupMembership.size(); i++)
735 –      cerr << "i = " << i << "\t ggm(i) = " << tmpGroupMembership[i] << "\n";
736 –
869   #else
870      info->setGlobalGroupMembership(globalGroupMembership);
871   #endif
872      
873      //fill molMembership
874 <    std::vector<int> globalMolMembership(info->getNGlobalAtoms(), 0);
874 >    std::vector<int> globalMolMembership(info->getNGlobalAtoms() +
875 >                                         info->getNGlobalRigidBodies(), 0);
876      
877 <    for(mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
877 >    for(mol = info->beginMolecule(mi); mol != NULL;
878 >        mol = info->nextMolecule(mi)) {
879        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
880          globalMolMembership[atom->getGlobalIndex()] = mol->getGlobalIndex();
881        }
882 +      for (rb = mol->beginRigidBody(ri); rb != NULL;
883 +           rb = mol->nextRigidBody(ri)) {
884 +        globalMolMembership[rb->getGlobalIndex()] = mol->getGlobalIndex();
885 +      }
886      }
887      
888   #ifdef IS_MPI
889 <    std::vector<int> tmpMolMembership(info->getNGlobalAtoms(), 0);
889 >    std::vector<int> tmpMolMembership(info->getNGlobalAtoms() +
890 >                                      info->getNGlobalRigidBodies(), 0);
891 >    MPI::COMM_WORLD.Allreduce(&globalMolMembership[0], &tmpMolMembership[0],
892 >                              nGlobalAtoms + nGlobalRigidBodies,
893 >                              MPI::INT, MPI::SUM);
894      
753 –    MPI_Allreduce(&globalMolMembership[0], &tmpMolMembership[0], nGlobalAtoms,
754 –                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
755 –    
895      info->setGlobalMolMembership(tmpMolMembership);
896   #else
897      info->setGlobalMolMembership(globalMolMembership);
# Line 762 | Line 901 | namespace OpenMD {
901      // here the molecules are listed by their global indices.
902  
903      std::vector<int> nIOPerMol(info->getNGlobalMolecules(), 0);
904 <    for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
904 >    for (mol = info->beginMolecule(mi); mol != NULL;
905 >         mol = info->nextMolecule(mi)) {
906        nIOPerMol[mol->getGlobalIndex()] = mol->getNIntegrableObjects();      
907      }
908      
909   #ifdef IS_MPI
910      std::vector<int> numIntegrableObjectsPerMol(info->getNGlobalMolecules(), 0);
911 <    MPI_Allreduce(&nIOPerMol[0], &numIntegrableObjectsPerMol[0],
912 <                  info->getNGlobalMolecules(), MPI_INT, MPI_SUM, MPI_COMM_WORLD);
911 >    MPI::COMM_WORLD.Allreduce(&nIOPerMol[0], &numIntegrableObjectsPerMol[0],
912 >                              info->getNGlobalMolecules(), MPI::INT, MPI::SUM);
913   #else
914      std::vector<int> numIntegrableObjectsPerMol = nIOPerMol;
915   #endif    
# Line 783 | Line 923 | namespace OpenMD {
923      }
924      
925      std::vector<StuntDouble*> IOIndexToIntegrableObject(info->getNGlobalIntegrableObjects(), (StuntDouble*)NULL);
926 <    for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
926 >    for (mol = info->beginMolecule(mi); mol != NULL;
927 >         mol = info->nextMolecule(mi)) {
928        int myGlobalIndex = mol->getGlobalIndex();
929        int globalIO = startingIOIndexForMol[myGlobalIndex];
930 <      for (StuntDouble* integrableObject = mol->beginIntegrableObject(ioi); integrableObject != NULL;
931 <           integrableObject = mol->nextIntegrableObject(ioi)) {
932 <        integrableObject->setGlobalIntegrableObjectIndex(globalIO);
933 <        IOIndexToIntegrableObject[globalIO] = integrableObject;
930 >      for (StuntDouble* sd = mol->beginIntegrableObject(ioi); sd != NULL;
931 >           sd = mol->nextIntegrableObject(ioi)) {
932 >        sd->setGlobalIntegrableObjectIndex(globalIO);
933 >        IOIndexToIntegrableObject[globalIO] = sd;
934          globalIO++;
935        }
936      }
937 <    
937 >      
938      info->setIOIndexToIntegrableObject(IOIndexToIntegrableObject);
939      
940    }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines