--- trunk/src/brains/SimInfo.cpp 2012/08/31 21:16:10 1793 +++ trunk/src/brains/SimInfo.cpp 2014/04/17 19:07:31 1987 @@ -35,7 +35,7 @@ * * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). - * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008). + * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008). * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). */ @@ -47,6 +47,9 @@ * @version 1.0 */ +#ifdef IS_MPI +#include +#endif #include #include #include @@ -61,9 +64,6 @@ #include "io/ForceFieldOptions.hpp" #include "brains/ForceField.hpp" #include "nonbonded/SwitchingFunction.hpp" -#ifdef IS_MPI -#include -#endif using namespace std; namespace OpenMD { @@ -72,11 +72,14 @@ namespace OpenMD { forceField_(ff), simParams_(simParams), ndf_(0), fdf_local(0), ndfRaw_(0), ndfTrans_(0), nZconstraint_(0), nGlobalMols_(0), nGlobalAtoms_(0), nGlobalCutoffGroups_(0), - nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), nGlobalFluctuatingCharges_(0), - nAtoms_(0), nBonds_(0), nBends_(0), nTorsions_(0), nInversions_(0), + nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), + nGlobalFluctuatingCharges_(0), nGlobalBonds_(0), nGlobalBends_(0), + nGlobalTorsions_(0), nGlobalInversions_(0), nGlobalConstraints_(0), + nAtoms_(0), nBonds_(0), nBends_(0), nTorsions_(0), nInversions_(0), nRigidBodies_(0), nIntegrableObjects_(0), nCutoffGroups_(0), - nConstraints_(0), nFluctuatingCharges_(0), sman_(NULL), topologyDone_(false), - calcBoxDipole_(false), useAtomicVirial_(true) { + nConstraints_(0), nFluctuatingCharges_(0), sman_(NULL), + topologyDone_(false), calcBoxDipole_(false), useAtomicVirial_(true), + hasNGlobalConstraints_(false) { MoleculeStamp* molStamp; int nMolWithSameStamp; @@ -91,12 +94,23 @@ namespace OpenMD { for (vector::iterator i = components.begin(); i !=components.end(); ++i) { molStamp = (*i)->getMoleculeStamp(); + if ( (*i)->haveRegion() ) { + molStamp->setRegion( (*i)->getRegion() ); + } else { + // set the region to a disallowed value: + molStamp->setRegion( -1 ); + } + nMolWithSameStamp = (*i)->getNMol(); addMoleculeStamp(molStamp, nMolWithSameStamp); //calculate atoms in molecules - nGlobalAtoms_ += molStamp->getNAtoms() *nMolWithSameStamp; + nGlobalAtoms_ += molStamp->getNAtoms() * nMolWithSameStamp; + nGlobalBonds_ += molStamp->getNBonds() * nMolWithSameStamp; + nGlobalBends_ += molStamp->getNBends() * nMolWithSameStamp; + nGlobalTorsions_ += molStamp->getNTorsions() * nMolWithSameStamp; + nGlobalInversions_ += molStamp->getNInversions() * nMolWithSameStamp; //calculate atoms in cutoff groups int nAtomsInGroups = 0; @@ -268,8 +282,9 @@ namespace OpenMD { ndf_local -= nConstraints_; #ifdef IS_MPI - MPI_Allreduce(&ndf_local,&ndf_,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); - MPI_Allreduce(&nfq_local,&nGlobalFluctuatingCharges_,1, MPI_INT, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&ndf_local, &ndf_, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&nfq_local, &nGlobalFluctuatingCharges_, 1, + MPI_INT, MPI_SUM, MPI_COMM_WORLD); #else ndf_ = ndf_local; nGlobalFluctuatingCharges_ = nfq_local; @@ -283,7 +298,7 @@ namespace OpenMD { int SimInfo::getFdf() { #ifdef IS_MPI - MPI_Allreduce(&fdf_local,&fdf_,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&fdf_local, &fdf_, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD); #else fdf_ = fdf_local; #endif @@ -339,7 +354,7 @@ namespace OpenMD { } #ifdef IS_MPI - MPI_Allreduce(&ndfRaw_local,&ndfRaw_,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&ndfRaw_local, &ndfRaw_, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD); #else ndfRaw_ = ndfRaw_local; #endif @@ -350,15 +365,14 @@ namespace OpenMD { ndfTrans_local = 3 * nIntegrableObjects_ - nConstraints_; - #ifdef IS_MPI - MPI_Allreduce(&ndfTrans_local,&ndfTrans_,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&ndfTrans_local, &ndfTrans_, 1, MPI_INT, MPI_SUM, + MPI_COMM_WORLD); #else ndfTrans_ = ndfTrans_local; #endif ndfTrans_ = ndfTrans_ - 3 - nZconstraint_; - } void SimInfo::addInteractionPairs(Molecule* mol) { @@ -409,6 +423,7 @@ namespace OpenMD { atomGroups.insert(map >::value_type(sd->getGlobalIndex(), oneAtomSet)); } } + for (bond= mol->beginBond(bondIter); bond != NULL; bond = mol->nextBond(bondIter)) { @@ -694,6 +709,7 @@ namespace OpenMD { */ void SimInfo::update() { setupSimVariables(); + calcNConstraints(); calcNdf(); calcNdfRaw(); calcNdfTrans(); @@ -733,7 +749,8 @@ namespace OpenMD { // count_local holds the number of found types on this processor int count_local = foundTypes.size(); - int nproc = MPI::COMM_WORLD.Get_size(); + int nproc; + MPI_Comm_size( MPI_COMM_WORLD, &nproc); // we need arrays to hold the counts and displacement vectors for // all processors @@ -741,8 +758,8 @@ namespace OpenMD { vector disps(nproc, 0); // fill the counts array - MPI::COMM_WORLD.Allgather(&count_local, 1, MPI::INT, &counts[0], - 1, MPI::INT); + MPI_Allgather(&count_local, 1, MPI_INT, &counts[0], + 1, MPI_INT, MPI_COMM_WORLD); // use the processor counts to compute the displacement array disps[0] = 0; @@ -756,9 +773,9 @@ namespace OpenMD { vector ftGlobal(totalCount); // now spray out the foundTypes to all the other processors: - MPI::COMM_WORLD.Allgatherv(&foundTypes[0], count_local, MPI::INT, - &ftGlobal[0], &counts[0], &disps[0], - MPI::INT); + MPI_Allgatherv(&foundTypes[0], count_local, MPI_INT, + &ftGlobal[0], &counts[0], &disps[0], + MPI_INT, MPI_COMM_WORLD); vector::iterator j; @@ -780,6 +797,22 @@ namespace OpenMD { return atomTypes; } + + int getGlobalCountOfType(AtomType* atype) { + /* + set atypes = getSimulatedAtomTypes(); + map counts_; + + for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { + for(atom = mol->beginAtom(ai); atom != NULL; + atom = mol->nextAtom(ai)) { + atom->getAtomType(); + } + } + */ + return 0; + } + void SimInfo::setupSimVariables() { useAtomicVirial_ = simParams_->getUseAtomicVirial(); // we only call setAccumulateBoxDipole if the accumulateBoxDipole @@ -806,22 +839,23 @@ namespace OpenMD { } #ifdef IS_MPI - bool temp; + int temp; + temp = usesDirectional; - MPI::COMM_WORLD.Allreduce(&temp, &usesDirectionalAtoms_, 1, MPI::BOOL, - MPI::LOR); - - temp = usesMetallic; - MPI::COMM_WORLD.Allreduce(&temp, &usesMetallicAtoms_, 1, MPI::BOOL, - MPI::LOR); + MPI_Allreduce(MPI_IN_PLACE, &temp, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); + usesDirectionalAtoms_ = (temp == 0) ? false : true; + temp = usesMetallic; + MPI_Allreduce(MPI_IN_PLACE, &temp, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); + usesMetallicAtoms_ = (temp == 0) ? false : true; + temp = usesElectrostatic; - MPI::COMM_WORLD.Allreduce(&temp, &usesElectrostaticAtoms_, 1, MPI::BOOL, - MPI::LOR); + MPI_Allreduce(MPI_IN_PLACE, &temp, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); + usesElectrostaticAtoms_ = (temp == 0) ? false : true; temp = usesFluctuatingCharges; - MPI::COMM_WORLD.Allreduce(&temp, &usesFluctuatingCharges_, 1, MPI::BOOL, - MPI::LOR); + MPI_Allreduce(MPI_IN_PLACE, &temp, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); + usesFluctuatingCharges_ = (temp == 0) ? false : true; #else usesDirectionalAtoms_ = usesDirectional; @@ -913,23 +947,21 @@ namespace OpenMD { } } - // Build the identArray_ + // Build the identArray_ and regions_ identArray_.clear(); - identArray_.reserve(getNAtoms()); - for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { + identArray_.reserve(getNAtoms()); + regions_.clear(); + regions_.reserve(getNAtoms()); + + for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { + int reg = mol->getRegion(); for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { identArray_.push_back(atom->getIdent()); + regions_.push_back(reg); } } - - //scan topology - - int* excludeList = excludedInteractions_.getPairList(); - int* oneTwoList = oneTwoInteractions_.getPairList(); - int* oneThreeList = oneThreeInteractions_.getPairList(); - int* oneFourList = oneFourInteractions_.getPairList(); - + topologyDone_ = true; } @@ -964,33 +996,55 @@ namespace OpenMD { delete sman_; sman_ = sman; - Molecule* mol; - RigidBody* rb; - Atom* atom; - CutoffGroup* cg; SimInfo::MoleculeIterator mi; + Molecule::AtomIterator ai; Molecule::RigidBodyIterator rbIter; - Molecule::AtomIterator atomIter; Molecule::CutoffGroupIterator cgIter; + Molecule::BondIterator bondIter; + Molecule::BendIterator bendIter; + Molecule::TorsionIterator torsionIter; + Molecule::InversionIterator inversionIter; + Molecule* mol; + Atom* atom; + RigidBody* rb; + CutoffGroup* cg; + Bond* bond; + Bend* bend; + Torsion* torsion; + Inversion* inversion; + for (mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { - for (atom = mol->beginAtom(atomIter); atom != NULL; - atom = mol->nextAtom(atomIter)) { + for (atom = mol->beginAtom(ai); atom != NULL; + atom = mol->nextAtom(ai)) { atom->setSnapshotManager(sman_); - } - + } for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) { rb->setSnapshotManager(sman_); } - for (cg = mol->beginCutoffGroup(cgIter); cg != NULL; cg = mol->nextCutoffGroup(cgIter)) { cg->setSnapshotManager(sman_); } - } - + for (bond = mol->beginBond(bondIter); bond != NULL; + bond = mol->nextBond(bondIter)) { + bond->setSnapshotManager(sman_); + } + for (bend = mol->beginBend(bendIter); bend != NULL; + bend = mol->nextBend(bendIter)) { + bend->setSnapshotManager(sman_); + } + for (torsion = mol->beginTorsion(torsionIter); torsion != NULL; + torsion = mol->nextTorsion(torsionIter)) { + torsion->setSnapshotManager(sman_); + } + for (inversion = mol->beginInversion(inversionIter); inversion != NULL; + inversion = mol->nextInversion(inversionIter)) { + inversion->setSnapshotManager(sman_); + } + } } @@ -1001,7 +1055,7 @@ namespace OpenMD { StuntDouble* SimInfo::getIOIndexToIntegrableObject(int index) { - if (index >= IOIndexToIntegrableObject.size()) { + if (index >= int(IOIndexToIntegrableObject.size())) { sprintf(painCave.errMsg, "SimInfo::getIOIndexToIntegrableObject Error: Integrable Object\n" "\tindex exceeds number of known objects!\n"); @@ -1016,15 +1070,13 @@ namespace OpenMD { IOIndexToIntegrableObject= v; } - int SimInfo::getNGlobalConstraints() { - int nGlobalConstraints; + void SimInfo::calcNConstraints() { #ifdef IS_MPI - MPI_Allreduce(&nConstraints_, &nGlobalConstraints, 1, MPI_INT, MPI_SUM, - MPI_COMM_WORLD); + MPI_Allreduce(&nConstraints_, &nGlobalConstraints_, 1, + MPI_INT, MPI_SUM, MPI_COMM_WORLD); #else - nGlobalConstraints = nConstraints_; + nGlobalConstraints_ = nConstraints_; #endif - return nGlobalConstraints; } }//end namespace OpenMD