| 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). | 
| 38 | > | * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008). | 
| 39 |  | * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010). | 
| 40 |  | * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). | 
| 41 |  | */ | 
| 47 |  | * @version 1.0 | 
| 48 |  | */ | 
| 49 |  |  | 
| 50 | + | #ifdef IS_MPI | 
| 51 | + | #include <mpi.h> | 
| 52 | + | #endif | 
| 53 |  | #include <algorithm> | 
| 54 |  | #include <set> | 
| 55 |  | #include <map> | 
| 64 |  | #include "io/ForceFieldOptions.hpp" | 
| 65 |  | #include "brains/ForceField.hpp" | 
| 66 |  | #include "nonbonded/SwitchingFunction.hpp" | 
| 64 | – | #ifdef IS_MPI | 
| 65 | – | #include <mpi.h> | 
| 66 | – | #endif | 
| 67 |  |  | 
| 68 |  | using namespace std; | 
| 69 |  | namespace OpenMD { | 
| 72 |  | forceField_(ff), simParams_(simParams), | 
| 73 |  | ndf_(0), fdf_local(0), ndfRaw_(0), ndfTrans_(0), nZconstraint_(0), | 
| 74 |  | nGlobalMols_(0), nGlobalAtoms_(0), nGlobalCutoffGroups_(0), | 
| 75 | < | nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), nGlobalFluctuatingCharges_(0), | 
| 76 | < | nAtoms_(0), nBonds_(0),  nBends_(0), nTorsions_(0), nInversions_(0), | 
| 77 | < | nRigidBodies_(0), nIntegrableObjects_(0), nCutoffGroups_(0), | 
| 78 | < | nConstraints_(0), nFluctuatingCharges_(0), sman_(NULL), topologyDone_(false), | 
| 75 | > | nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), | 
| 76 | > | nGlobalFluctuatingCharges_(0), nGlobalBonds_(0), nGlobalBends_(0), | 
| 77 | > | nGlobalTorsions_(0), nGlobalInversions_(0), nAtoms_(0), nBonds_(0), | 
| 78 | > | nBends_(0), nTorsions_(0), nInversions_(0), nRigidBodies_(0), | 
| 79 | > | nIntegrableObjects_(0), nCutoffGroups_(0), nConstraints_(0), | 
| 80 | > | nFluctuatingCharges_(0), sman_(NULL), topologyDone_(false), | 
| 81 |  | calcBoxDipole_(false), useAtomicVirial_(true) { | 
| 82 |  |  | 
| 83 |  | MoleculeStamp* molStamp; | 
| 93 |  | for (vector<Component*>::iterator i = components.begin(); | 
| 94 |  | i !=components.end(); ++i) { | 
| 95 |  | molStamp = (*i)->getMoleculeStamp(); | 
| 96 | + | if ( (*i)->haveRegion() ) { | 
| 97 | + | molStamp->setRegion( (*i)->getRegion() ); | 
| 98 | + | } else { | 
| 99 | + | // set the region to a disallowed value: | 
| 100 | + | molStamp->setRegion( -1 ); | 
| 101 | + | } | 
| 102 | + |  | 
| 103 |  | nMolWithSameStamp = (*i)->getNMol(); | 
| 104 |  |  | 
| 105 |  | addMoleculeStamp(molStamp, nMolWithSameStamp); | 
| 106 |  |  | 
| 107 |  | //calculate atoms in molecules | 
| 108 | < | nGlobalAtoms_ += molStamp->getNAtoms() *nMolWithSameStamp; | 
| 108 | > | nGlobalAtoms_ += molStamp->getNAtoms() * nMolWithSameStamp; | 
| 109 | > | nGlobalBonds_ += molStamp->getNBonds() * nMolWithSameStamp; | 
| 110 | > | nGlobalBends_ += molStamp->getNBends() * nMolWithSameStamp; | 
| 111 | > | nGlobalTorsions_ += molStamp->getNTorsions() * nMolWithSameStamp; | 
| 112 | > | nGlobalInversions_ += molStamp->getNInversions() * nMolWithSameStamp; | 
| 113 |  |  | 
| 114 |  | //calculate atoms in cutoff groups | 
| 115 |  | int nAtomsInGroups = 0; | 
| 424 |  | atomGroups.insert(map<int, set<int> >::value_type(sd->getGlobalIndex(), oneAtomSet)); | 
| 425 |  | } | 
| 426 |  | } | 
| 427 | + |  | 
| 428 |  |  | 
| 429 |  | for (bond= mol->beginBond(bondIter); bond != NULL; | 
| 430 |  | bond = mol->nextBond(bondIter)) { | 
| 795 |  |  | 
| 796 |  | return atomTypes; | 
| 797 |  | } | 
| 798 | + |  | 
| 799 | + |  | 
| 800 | + | int getGlobalCountOfType(AtomType* atype) { | 
| 801 | + | /* | 
| 802 | + | set<AtomType*> atypes = getSimulatedAtomTypes(); | 
| 803 | + | map<AtomType*, int> counts_; | 
| 804 |  |  | 
| 805 | + | for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { | 
| 806 | + | for(atom = mol->beginAtom(ai); atom != NULL; | 
| 807 | + | atom = mol->nextAtom(ai)) { | 
| 808 | + | atom->getAtomType(); | 
| 809 | + | } | 
| 810 | + | } | 
| 811 | + | */ | 
| 812 | + | return 0; | 
| 813 | + | } | 
| 814 | + |  | 
| 815 |  | void SimInfo::setupSimVariables() { | 
| 816 |  | useAtomicVirial_ = simParams_->getUseAtomicVirial(); | 
| 817 |  | // we only call setAccumulateBoxDipole if the accumulateBoxDipole | 
| 945 |  | } | 
| 946 |  | } | 
| 947 |  |  | 
| 948 | < | // Build the identArray_ | 
| 948 | > | // Build the identArray_ and regions_ | 
| 949 |  |  | 
| 950 |  | identArray_.clear(); | 
| 951 | < | identArray_.reserve(getNAtoms()); | 
| 952 | < | for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { | 
| 951 | > | identArray_.reserve(getNAtoms()); | 
| 952 | > | regions_.clear(); | 
| 953 | > | regions_.reserve(getNAtoms()); | 
| 954 | > |  | 
| 955 | > | for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { | 
| 956 | > | int reg = mol->getRegion(); | 
| 957 |  | for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { | 
| 958 |  | identArray_.push_back(atom->getIdent()); | 
| 959 | + | regions_.push_back(reg); | 
| 960 |  | } | 
| 961 |  | } | 
| 962 | < |  | 
| 928 | < | //scan topology | 
| 929 | < |  | 
| 930 | < | int* excludeList = excludedInteractions_.getPairList(); | 
| 931 | < | int* oneTwoList = oneTwoInteractions_.getPairList(); | 
| 932 | < | int* oneThreeList = oneThreeInteractions_.getPairList(); | 
| 933 | < | int* oneFourList = oneFourInteractions_.getPairList(); | 
| 934 | < |  | 
| 962 | > |  | 
| 963 |  | topologyDone_ = true; | 
| 964 |  | } | 
| 965 |  |  | 
| 994 |  | delete sman_; | 
| 995 |  | sman_ = sman; | 
| 996 |  |  | 
| 969 | – | Molecule* mol; | 
| 970 | – | RigidBody* rb; | 
| 971 | – | Atom* atom; | 
| 972 | – | CutoffGroup* cg; | 
| 997 |  | SimInfo::MoleculeIterator mi; | 
| 998 | + | Molecule::AtomIterator ai; | 
| 999 |  | Molecule::RigidBodyIterator rbIter; | 
| 975 | – | Molecule::AtomIterator atomIter; | 
| 1000 |  | Molecule::CutoffGroupIterator cgIter; | 
| 1001 | + | Molecule::BondIterator bondIter; | 
| 1002 | + | Molecule::BendIterator bendIter; | 
| 1003 | + | Molecule::TorsionIterator torsionIter; | 
| 1004 | + | Molecule::InversionIterator inversionIter; | 
| 1005 |  |  | 
| 1006 | + | Molecule* mol; | 
| 1007 | + | Atom* atom; | 
| 1008 | + | RigidBody* rb; | 
| 1009 | + | CutoffGroup* cg; | 
| 1010 | + | Bond* bond; | 
| 1011 | + | Bend* bend; | 
| 1012 | + | Torsion* torsion; | 
| 1013 | + | Inversion* inversion; | 
| 1014 | + |  | 
| 1015 |  | for (mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { | 
| 1016 |  |  | 
| 1017 | < | for (atom = mol->beginAtom(atomIter); atom != NULL; | 
| 1018 | < | atom = mol->nextAtom(atomIter)) { | 
| 1017 | > | for (atom = mol->beginAtom(ai); atom != NULL; | 
| 1018 | > | atom = mol->nextAtom(ai)) { | 
| 1019 |  | atom->setSnapshotManager(sman_); | 
| 1020 | < | } | 
| 984 | < |  | 
| 1020 | > | } | 
| 1021 |  | for (rb = mol->beginRigidBody(rbIter); rb != NULL; | 
| 1022 |  | rb = mol->nextRigidBody(rbIter)) { | 
| 1023 |  | rb->setSnapshotManager(sman_); | 
| 1024 |  | } | 
| 989 | – |  | 
| 1025 |  | for (cg = mol->beginCutoffGroup(cgIter); cg != NULL; | 
| 1026 |  | cg = mol->nextCutoffGroup(cgIter)) { | 
| 1027 |  | cg->setSnapshotManager(sman_); | 
| 1028 |  | } | 
| 1029 | < | } | 
| 1030 | < |  | 
| 1029 | > | for (bond = mol->beginBond(bondIter); bond != NULL; | 
| 1030 | > | bond = mol->nextBond(bondIter)) { | 
| 1031 | > | bond->setSnapshotManager(sman_); | 
| 1032 | > | } | 
| 1033 | > | for (bend = mol->beginBend(bendIter); bend != NULL; | 
| 1034 | > | bend = mol->nextBend(bendIter)) { | 
| 1035 | > | bend->setSnapshotManager(sman_); | 
| 1036 | > | } | 
| 1037 | > | for (torsion = mol->beginTorsion(torsionIter); torsion != NULL; | 
| 1038 | > | torsion = mol->nextTorsion(torsionIter)) { | 
| 1039 | > | torsion->setSnapshotManager(sman_); | 
| 1040 | > | } | 
| 1041 | > | for (inversion = mol->beginInversion(inversionIter); inversion != NULL; | 
| 1042 | > | inversion = mol->nextInversion(inversionIter)) { | 
| 1043 | > | inversion->setSnapshotManager(sman_); | 
| 1044 | > | } | 
| 1045 | > | } | 
| 1046 |  | } | 
| 1047 |  |  | 
| 1048 |  |  | 
| 1053 |  |  | 
| 1054 |  |  | 
| 1055 |  | StuntDouble* SimInfo::getIOIndexToIntegrableObject(int index) { | 
| 1056 | < | if (index >= IOIndexToIntegrableObject.size()) { | 
| 1056 | > | if (index >= int(IOIndexToIntegrableObject.size())) { | 
| 1057 |  | sprintf(painCave.errMsg, | 
| 1058 |  | "SimInfo::getIOIndexToIntegrableObject Error: Integrable Object\n" | 
| 1059 |  | "\tindex exceeds number of known objects!\n"); |