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 1627 by gezelter, Tue Sep 13 22:05:04 2011 UTC vs.
Revision 1794 by gezelter, Thu Sep 6 19:44:06 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   /**
# 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"
# Line 247 | 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 362 | 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 413 | Line 424 | namespace OpenMD {
424      //create the molecules
425      createMolecules(info);
426      
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 652 | Line 667 | namespace OpenMD {
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 hasDipoles = false;    
679 +    bool hasQuadrupoles = false;    
680 +    bool hasPolarizable = false;    
681 +    bool hasFluctuatingCharge = false;    
682 +    bool hasMetallic = false;
683 +    int storageLayout = 0;
684 +    storageLayout |= DataStorage::dslPosition;
685 +    storageLayout |= DataStorage::dslVelocity;
686 +    storageLayout |= DataStorage::dslForce;
687 +
688 +    for (i = atomTypes.begin(); i != atomTypes.end(); ++i) {
689 +
690 +      DirectionalAdapter da = DirectionalAdapter( (*i) );
691 +      MultipoleAdapter ma = MultipoleAdapter( (*i) );
692 +      EAMAdapter ea = EAMAdapter( (*i) );
693 +      SuttonChenAdapter sca = SuttonChenAdapter( (*i) );
694 +      PolarizableAdapter pa = PolarizableAdapter( (*i) );
695 +      FixedChargeAdapter fca = FixedChargeAdapter( (*i) );
696 +      FluctuatingChargeAdapter fqa = FluctuatingChargeAdapter( (*i) );
697 +
698 +      if (da.isDirectional()){
699 +        hasDirectionalAtoms = true;
700 +      }
701 +      if (ma.isDipole()){
702 +        hasDipoles = true;
703 +      }
704 +      if (ma.isQuadrupole()){
705 +        hasQuadrupoles = true;
706 +      }
707 +      if (ea.isEAM() || sca.isSuttonChen()){
708 +        hasMetallic = true;
709 +      }
710 +      if ( fca.isFixedCharge() ){
711 +        hasFixedCharge = true;
712 +      }
713 +      if ( fqa.isFluctuatingCharge() ){
714 +        hasFluctuatingCharge = true;
715 +      }
716 +      if ( pa.isPolarizable() ){
717 +        hasPolarizable = true;
718 +      }
719 +    }
720 +    
721 +    if (nRigidBodies > 0 || hasDirectionalAtoms) {
722 +      storageLayout |= DataStorage::dslAmat;
723 +      if(storageLayout & DataStorage::dslVelocity) {
724 +        storageLayout |= DataStorage::dslAngularMomentum;
725 +      }
726 +      if (storageLayout & DataStorage::dslForce) {
727 +        storageLayout |= DataStorage::dslTorque;
728 +      }
729 +    }
730 +    if (hasDipoles) {
731 +      storageLayout |= DataStorage::dslDipole;
732 +    }
733 +    if (hasQuadrupoles) {
734 +      storageLayout |= DataStorage::dslQuadrupole;
735 +    }
736 +    if (hasFixedCharge || hasFluctuatingCharge) {
737 +      storageLayout |= DataStorage::dslSkippedCharge;
738 +    }
739 +    if (hasMetallic) {
740 +      storageLayout |= DataStorage::dslDensity;
741 +      storageLayout |= DataStorage::dslFunctional;
742 +      storageLayout |= DataStorage::dslFunctionalDerivative;
743 +    }
744 +    if (hasPolarizable) {
745 +      storageLayout |= DataStorage::dslElectricField;
746 +    }
747 +    if (hasFluctuatingCharge){
748 +      storageLayout |= DataStorage::dslFlucQPosition;
749 +      if(storageLayout & DataStorage::dslVelocity) {
750 +        storageLayout |= DataStorage::dslFlucQVelocity;
751 +      }
752 +      if (storageLayout & DataStorage::dslForce) {
753 +        storageLayout |= DataStorage::dslFlucQForce;
754 +      }
755 +    }
756 +    
757 +    // if the user has asked for them, make sure we've got the memory for the
758 +    // objects defined.
759 +
760 +    if (simParams->getOutputParticlePotential()) {
761 +      storageLayout |= DataStorage::dslParticlePot;
762 +    }
763 +
764 +    if (simParams->havePrintHeatFlux()) {
765 +      if (simParams->getPrintHeatFlux()) {
766 +        storageLayout |= DataStorage::dslParticlePot;
767 +      }
768 +    }
769 +
770 +    if (simParams->getOutputElectricField()) {
771 +      storageLayout |= DataStorage::dslElectricField;
772 +    }
773 +
774 +    if (simParams->getOutputFluctuatingCharges()) {
775 +      storageLayout |= DataStorage::dslFlucQPosition;
776 +      storageLayout |= DataStorage::dslFlucQVelocity;
777 +      storageLayout |= DataStorage::dslFlucQForce;
778 +    }
779 +
780 +    return storageLayout;
781 +  }
782 +
783    void SimCreator::setGlobalIndex(SimInfo *info) {
784      SimInfo::MoleculeIterator mi;
785      Molecule::AtomIterator ai;
# Line 668 | Line 796 | namespace OpenMD {
796      int nGlobalAtoms = info->getNGlobalAtoms();
797      
798      beginAtomIndex = 0;
799 <    beginRigidBodyIndex = 0;
799 >    //rigidbody's index begins right after atom's
800 >    beginRigidBodyIndex = info->getNGlobalAtoms();
801      beginCutoffGroupIndex = 0;
802  
803      for(int i = 0; i < info->getNGlobalMolecules(); i++) {
# Line 786 | Line 915 | namespace OpenMD {
915      for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
916        int myGlobalIndex = mol->getGlobalIndex();
917        int globalIO = startingIOIndexForMol[myGlobalIndex];
918 <      for (StuntDouble* integrableObject = mol->beginIntegrableObject(ioi); integrableObject != NULL;
919 <           integrableObject = mol->nextIntegrableObject(ioi)) {
920 <        integrableObject->setGlobalIntegrableObjectIndex(globalIO);
921 <        IOIndexToIntegrableObject[globalIO] = integrableObject;
918 >      for (StuntDouble* sd = mol->beginIntegrableObject(ioi); sd != NULL;
919 >           sd = mol->nextIntegrableObject(ioi)) {
920 >        sd->setGlobalIntegrableObjectIndex(globalIO);
921 >        IOIndexToIntegrableObject[globalIO] = sd;
922          globalIO++;
923        }
924      }

Diff Legend

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