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 1711 by gezelter, Sat May 19 02:58:35 2012 UTC vs.
Revision 1803 by gezelter, Wed Oct 3 14:20:07 2012 UTC

# Line 56 | 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 256 | Line 256 | namespace OpenMD {
256      int metaDataBlockStart = -1;
257      int metaDataBlockEnd = -1;
258      int i;
259 <    int mdOffset;
259 >    streamoff mdOffset;
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 371 | Line 373 | namespace OpenMD {
373                                     metaDataBlockStart + 1);
374      
375      //create the force field
376 <    ForceField * ff = ForceFieldFactory::getInstance()->createForceField(simParams->getForceField());
376 >    ForceField * ff = new ForceField(simParams->getForceField());
377  
378      if (ff == NULL) {
379        sprintf(painCave.errMsg,
# Line 426 | Line 428 | namespace OpenMD {
428  
429      int storageLayout = computeStorageLayout(info);
430  
429 –    cerr << "computed Storage Layout = " << storageLayout << "\n";
430 –
431      //allocate memory for DataStorage(circular reference, need to
432      //break it)
433      info->setSnapshotManager(new SimSnapshotManager(info, storageLayout));
# Line 507 | Line 507 | namespace OpenMD {
507      int nGlobalMols = info->getNGlobalMolecules();
508      std::vector<int> molToProcMap(nGlobalMols, -1); // default to an error condition:
509      
510 <    MPI_Comm_size(MPI_COMM_WORLD, &nProcessors);
510 >    nProcessors = MPI::COMM_WORLD.Get_size();
511      
512      if (nProcessors > nGlobalMols) {
513        sprintf(painCave.errMsg,
# Line 545 | Line 545 | namespace OpenMD {
545        nTarget = (int)(precast + 0.5);
546        
547        for(i = 0; i < nGlobalMols; i++) {
548 +
549          done = 0;
550          loops = 0;
551          
# Line 569 | Line 570 | namespace OpenMD {
570            // and be done with it.
571            
572            if (loops > 100) {
573 +
574              sprintf(painCave.errMsg,
575 <                    "I've tried 100 times to assign molecule %d to a "
576 <                    " processor, but can't find a good spot.\n"
577 <                    "I'm assigning it at random to processor %d.\n",
575 >                    "There have been 100 attempts to assign molecule %d to an\n"
576 >                    "\tunderworked processor, but there's no good place to\n"
577 >                    "\tleave it.  OpenMD is assigning it at random to processor %d.\n",
578                      i, which_proc);
579 <            
579 >          
580              painCave.isFatal = 0;
581 +            painCave.severity = OPENMD_INFO;
582              simError();
583              
584              molToProcMap[i] = which_proc;
# Line 620 | Line 623 | namespace OpenMD {
623        }
624        
625        delete myRandom;
626 <      
626 >
627        // Spray out this nonsense to all other processors:
628 <      
626 <      MPI_Bcast(&molToProcMap[0], nGlobalMols, MPI_INT, 0, MPI_COMM_WORLD);
628 >      MPI::COMM_WORLD.Bcast(&molToProcMap[0], nGlobalMols, MPI::INT, 0);
629      } else {
630        
631        // Listen to your marching orders from processor 0:
632 <      
633 <      MPI_Bcast(&molToProcMap[0], nGlobalMols, MPI_INT, 0, MPI_COMM_WORLD);
632 >      MPI::COMM_WORLD.Bcast(&molToProcMap[0], nGlobalMols, MPI::INT, 0);
633 >
634      }
635      
636      info->setMolToProcMap(molToProcMap);
# Line 675 | Line 677 | namespace OpenMD {
677      set<AtomType*>::iterator i;
678      bool hasDirectionalAtoms = false;
679      bool hasFixedCharge = false;
680 <    bool hasMultipoles = false;    
680 >    bool hasDipoles = false;    
681 >    bool hasQuadrupoles = false;    
682      bool hasPolarizable = false;    
683      bool hasFluctuatingCharge = false;    
684      bool hasMetallic = false;
# Line 697 | Line 700 | namespace OpenMD {
700        if (da.isDirectional()){
701          hasDirectionalAtoms = true;
702        }
703 <      if (ma.isMultipole()){
704 <        hasMultipoles = true;
703 >      if (ma.isDipole()){
704 >        hasDipoles = true;
705        }
706 +      if (ma.isQuadrupole()){
707 +        hasQuadrupoles = true;
708 +      }
709        if (ea.isEAM() || sca.isSuttonChen()){
710          hasMetallic = true;
711        }
# Line 723 | Line 729 | namespace OpenMD {
729          storageLayout |= DataStorage::dslTorque;
730        }
731      }
732 <    if (hasMultipoles) {
733 <      storageLayout |= DataStorage::dslElectroFrame;
732 >    if (hasDipoles) {
733 >      storageLayout |= DataStorage::dslDipole;
734      }
735 +    if (hasQuadrupoles) {
736 +      storageLayout |= DataStorage::dslQuadrupole;
737 +    }
738      if (hasFixedCharge || hasFluctuatingCharge) {
739        storageLayout |= DataStorage::dslSkippedCharge;
740      }
# Line 746 | Line 755 | namespace OpenMD {
755          storageLayout |= DataStorage::dslFlucQForce;
756        }
757      }
758 +    
759 +    // if the user has asked for them, make sure we've got the memory for the
760 +    // objects defined.
761  
762      if (simParams->getOutputParticlePotential()) {
763        storageLayout |= DataStorage::dslParticlePot;
764 +    }
765 +
766 +    if (simParams->havePrintHeatFlux()) {
767 +      if (simParams->getPrintHeatFlux()) {
768 +        storageLayout |= DataStorage::dslParticlePot;
769 +      }
770      }
771  
772 +    if (simParams->getOutputElectricField()) {
773 +      storageLayout |= DataStorage::dslElectricField;
774 +    }
775 +
776 +    if (simParams->getOutputFluctuatingCharges()) {
777 +      storageLayout |= DataStorage::dslFlucQPosition;
778 +      storageLayout |= DataStorage::dslFlucQVelocity;
779 +      storageLayout |= DataStorage::dslFlucQForce;
780 +    }
781 +
782      return storageLayout;
783    }
784  
# Line 768 | Line 796 | namespace OpenMD {
796      int beginRigidBodyIndex;
797      int beginCutoffGroupIndex;
798      int nGlobalAtoms = info->getNGlobalAtoms();
799 +    int nGlobalRigidBodies = info->getNGlobalRigidBodies();
800      
801      beginAtomIndex = 0;
802 <    beginRigidBodyIndex = 0;
802 >    //rigidbody's index begins right after atom's
803 >    beginRigidBodyIndex = info->getNGlobalAtoms();
804      beginCutoffGroupIndex = 0;
805  
806      for(int i = 0; i < info->getNGlobalMolecules(); i++) {
# Line 833 | Line 863 | namespace OpenMD {
863      // This would be prettier if we could use MPI_IN_PLACE like the MPI-2
864      // docs said we could.
865      std::vector<int> tmpGroupMembership(info->getNGlobalAtoms(), 0);
866 <    MPI_Allreduce(&globalGroupMembership[0], &tmpGroupMembership[0], nGlobalAtoms,
867 <                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
866 >    MPI::COMM_WORLD.Allreduce(&globalGroupMembership[0],
867 >                              &tmpGroupMembership[0], nGlobalAtoms,
868 >                              MPI::INT, MPI::SUM);
869      info->setGlobalGroupMembership(tmpGroupMembership);
870   #else
871      info->setGlobalGroupMembership(globalGroupMembership);
872   #endif
873      
874      //fill molMembership
875 <    std::vector<int> globalMolMembership(info->getNGlobalAtoms(), 0);
875 >    std::vector<int> globalMolMembership(info->getNGlobalAtoms() +
876 >                                         info->getNGlobalRigidBodies(), 0);
877      
878 <    for(mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
878 >    for(mol = info->beginMolecule(mi); mol != NULL;
879 >        mol = info->nextMolecule(mi)) {
880        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
881          globalMolMembership[atom->getGlobalIndex()] = mol->getGlobalIndex();
882        }
883 +      for (rb = mol->beginRigidBody(ri); rb != NULL;
884 +           rb = mol->nextRigidBody(ri)) {
885 +        globalMolMembership[rb->getGlobalIndex()] = mol->getGlobalIndex();
886 +      }
887      }
888      
889   #ifdef IS_MPI
890 <    std::vector<int> tmpMolMembership(info->getNGlobalAtoms(), 0);
890 >    std::vector<int> tmpMolMembership(info->getNGlobalAtoms() +
891 >                                      info->getNGlobalRigidBodies(), 0);
892 >    MPI::COMM_WORLD.Allreduce(&globalMolMembership[0], &tmpMolMembership[0],
893 >                              nGlobalAtoms + nGlobalRigidBodies,
894 >                              MPI::INT, MPI::SUM);
895      
855 –    MPI_Allreduce(&globalMolMembership[0], &tmpMolMembership[0], nGlobalAtoms,
856 –                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
857 –    
896      info->setGlobalMolMembership(tmpMolMembership);
897   #else
898      info->setGlobalMolMembership(globalMolMembership);
# Line 864 | Line 902 | namespace OpenMD {
902      // here the molecules are listed by their global indices.
903  
904      std::vector<int> nIOPerMol(info->getNGlobalMolecules(), 0);
905 <    for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
905 >    for (mol = info->beginMolecule(mi); mol != NULL;
906 >         mol = info->nextMolecule(mi)) {
907        nIOPerMol[mol->getGlobalIndex()] = mol->getNIntegrableObjects();      
908      }
909      
910   #ifdef IS_MPI
911      std::vector<int> numIntegrableObjectsPerMol(info->getNGlobalMolecules(), 0);
912 <    MPI_Allreduce(&nIOPerMol[0], &numIntegrableObjectsPerMol[0],
913 <                  info->getNGlobalMolecules(), MPI_INT, MPI_SUM, MPI_COMM_WORLD);
912 >    MPI::COMM_WORLD.Allreduce(&nIOPerMol[0], &numIntegrableObjectsPerMol[0],
913 >                              info->getNGlobalMolecules(), MPI::INT, MPI::SUM);
914   #else
915      std::vector<int> numIntegrableObjectsPerMol = nIOPerMol;
916   #endif    
# Line 885 | Line 924 | namespace OpenMD {
924      }
925      
926      std::vector<StuntDouble*> IOIndexToIntegrableObject(info->getNGlobalIntegrableObjects(), (StuntDouble*)NULL);
927 <    for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
927 >    for (mol = info->beginMolecule(mi); mol != NULL;
928 >         mol = info->nextMolecule(mi)) {
929        int myGlobalIndex = mol->getGlobalIndex();
930        int globalIO = startingIOIndexForMol[myGlobalIndex];
931 <      for (StuntDouble* integrableObject = mol->beginIntegrableObject(ioi); integrableObject != NULL;
932 <           integrableObject = mol->nextIntegrableObject(ioi)) {
933 <        integrableObject->setGlobalIntegrableObjectIndex(globalIO);
934 <        IOIndexToIntegrableObject[globalIO] = integrableObject;
931 >      for (StuntDouble* sd = mol->beginIntegrableObject(ioi); sd != NULL;
932 >           sd = mol->nextIntegrableObject(ioi)) {
933 >        sd->setGlobalIntegrableObjectIndex(globalIO);
934 >        IOIndexToIntegrableObject[globalIO] = sd;
935          globalIO++;
936        }
937      }

Diff Legend

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