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 1715 by gezelter, Tue May 22 21:55:31 2012 UTC vs.
Revision 1808 by gezelter, Mon Oct 22 20:42:10 2012 UTC

# Line 44 | Line 44
44   * @file SimCreator.cpp
45   * @author tlin
46   * @date 11/03/2004
47 – * @time 13:51am
47   * @version 1.0
48   */
49   #include <exception>
# Line 56 | 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 256 | 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 371 | Line 372 | namespace OpenMD {
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 505 | 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 543 | 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 567 | 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",
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 <            
578 >          
579              painCave.isFatal = 0;
580 +            painCave.severity = OPENMD_INFO;
581              simError();
582              
583              molToProcMap[i] = which_proc;
# Line 618 | Line 622 | namespace OpenMD {
622        }
623        
624        delete myRandom;
625 <      
625 >
626        // Spray out this nonsense to all other processors:
627 <      
624 <      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 673 | Line 676 | namespace OpenMD {
676      set<AtomType*>::iterator i;
677      bool hasDirectionalAtoms = false;
678      bool hasFixedCharge = false;
679 <    bool hasMultipoles = false;    
679 >    bool hasDipoles = false;    
680 >    bool hasQuadrupoles = false;    
681      bool hasPolarizable = false;    
682      bool hasFluctuatingCharge = false;    
683      bool hasMetallic = false;
# Line 695 | Line 699 | namespace OpenMD {
699        if (da.isDirectional()){
700          hasDirectionalAtoms = true;
701        }
702 <      if (ma.isMultipole()){
703 <        hasMultipoles = true;
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        }
# Line 721 | Line 728 | namespace OpenMD {
728          storageLayout |= DataStorage::dslTorque;
729        }
730      }
731 <    if (hasMultipoles) {
732 <      storageLayout |= DataStorage::dslElectroFrame;
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      }
# Line 750 | Line 760 | namespace OpenMD {
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;
# Line 777 | Line 795 | namespace OpenMD {
795      int beginRigidBodyIndex;
796      int beginCutoffGroupIndex;
797      int nGlobalAtoms = info->getNGlobalAtoms();
798 +    int nGlobalRigidBodies = info->getNGlobalRigidBodies();
799      
800      beginAtomIndex = 0;
801 <    beginRigidBodyIndex = 0;
801 >    //rigidbody's index begins right after atom's
802 >    beginRigidBodyIndex = info->getNGlobalAtoms();
803      beginCutoffGroupIndex = 0;
804  
805      for(int i = 0; i < info->getNGlobalMolecules(); i++) {
# Line 842 | 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_SUM, MPI_COMM_WORLD);
865 >    MPI::COMM_WORLD.Allreduce(&globalGroupMembership[0],
866 >                              &tmpGroupMembership[0], nGlobalAtoms,
867 >                              MPI::INT, MPI::SUM);
868      info->setGlobalGroupMembership(tmpGroupMembership);
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      
864 –    MPI_Allreduce(&globalMolMembership[0], &tmpMolMembership[0], nGlobalAtoms,
865 –                  MPI_INT, MPI_SUM, MPI_COMM_WORLD);
866 –    
895      info->setGlobalMolMembership(tmpMolMembership);
896   #else
897      info->setGlobalMolMembership(globalMolMembership);
# Line 873 | 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 894 | 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      }

Diff Legend

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