--- branches/development/src/brains/SimInfo.cpp 2011/04/29 17:25:12 1553 +++ branches/development/src/brains/SimInfo.cpp 2011/06/21 06:34:35 1586 @@ -71,7 +71,7 @@ namespace OpenMD { nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), nAtoms_(0), nBonds_(0), nBends_(0), nTorsions_(0), nInversions_(0), nRigidBodies_(0), nIntegrableObjects_(0), nCutoffGroups_(0), - nConstraints_(0), sman_(NULL), fortranInitialized_(false), + nConstraints_(0), sman_(NULL), topologyDone_(false), calcBoxDipole_(false), useAtomicVirial_(true) { MoleculeStamp* molStamp; @@ -125,13 +125,8 @@ namespace OpenMD { //equal to the total number of atoms minus number of atoms belong to //cutoff group defined in meta-data file plus the number of cutoff //groups defined in meta-data file - std::cerr << "nGA = " << nGlobalAtoms_ << "\n"; - std::cerr << "nCA = " << nCutoffAtoms << "\n"; - std::cerr << "nG = " << nGroups << "\n"; nGlobalCutoffGroups_ = nGlobalAtoms_ - nCutoffAtoms + nGroups; - - std::cerr << "nGCG = " << nGlobalCutoffGroups_ << "\n"; //every free atom (atom does not belong to rigid bodies) is an //integrable object therefore the total number of integrable objects @@ -273,6 +268,25 @@ namespace OpenMD { fdf_ = fdf_local; #endif return fdf_; + } + + unsigned int SimInfo::getNLocalCutoffGroups(){ + int nLocalCutoffAtoms = 0; + Molecule* mol; + MoleculeIterator mi; + CutoffGroup* cg; + Molecule::CutoffGroupIterator ci; + + for (mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { + + for (cg = mol->beginCutoffGroup(ci); cg != NULL; + cg = mol->nextCutoffGroup(ci)) { + nLocalCutoffAtoms += cg->getNumAtom(); + + } + } + + return nAtoms_ - nLocalCutoffAtoms + nCutoffGroups_; } void SimInfo::calcNdfRaw() { @@ -745,7 +759,7 @@ namespace OpenMD { if ( simParams_->getAccumulateBoxDipole() ) { calcBoxDipole_ = true; } - + set::iterator i; set atomTypes; atomTypes = getSimulatedAtomTypes(); @@ -758,18 +772,28 @@ namespace OpenMD { usesMetallic |= (*i)->isMetal(); usesDirectional |= (*i)->isDirectional(); } - + #ifdef IS_MPI int temp; temp = usesDirectional; MPI_Allreduce(&temp, &usesDirectionalAtoms_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); - + temp = usesMetallic; MPI_Allreduce(&temp, &usesMetallicAtoms_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); - + temp = usesElectrostatic; MPI_Allreduce(&temp, &usesElectrostaticAtoms_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); +#else + + usesDirectionalAtoms_ = usesDirectional; + usesMetallicAtoms_ = usesMetallic; + usesElectrostaticAtoms_ = usesElectrostatic; + #endif + + requiresPrepair_ = usesMetallicAtoms_ ? true : false; + requiresSkipCorrection_ = usesElectrostaticAtoms_ ? true : false; + requiresSelfCorrection_ = usesElectrostaticAtoms_ ? true : false; } @@ -812,20 +836,10 @@ namespace OpenMD { } - void SimInfo::setupFortran() { - int isError; + void SimInfo::prepareTopology() { int nExclude, nOneTwo, nOneThree, nOneFour; - vector fortranGlobalGroupMembership; - - isError = 0; - //globalGroupMembership_ is filled by SimCreator - for (int i = 0; i < nGlobalAtoms_; i++) { - fortranGlobalGroupMembership.push_back(globalGroupMembership_[i] + 1); - } - //calculate mass ratio of cutoff group - vector mfact; SimInfo::MoleculeIterator mi; Molecule* mol; Molecule::CutoffGroupIterator ci; @@ -834,19 +848,29 @@ namespace OpenMD { Atom* atom; RealType totalMass; - //to avoid memory reallocation, reserve enough space for mfact - mfact.reserve(getNCutoffGroups()); + /** + * The mass factor is the relative mass of an atom to the total + * mass of the cutoff group it belongs to. By default, all atoms + * are their own cutoff groups, and therefore have mass factors of + * 1. We need some special handling for massless atoms, which + * will be treated as carrying the entire mass of the cutoff + * group. + */ + massFactors_.clear(); + massFactors_.resize(getNAtoms(), 1.0); + cerr << "mfs in si = " << massFactors_.size() << "\n"; for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { - for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci)) { + for (cg = mol->beginCutoffGroup(ci); cg != NULL; + cg = mol->nextCutoffGroup(ci)) { totalMass = cg->getMass(); for(atom = cg->beginAtom(ai); atom != NULL; atom = cg->nextAtom(ai)) { // Check for massless groups - set mfact to 1 if true - if (totalMass != 0) - mfact.push_back(atom->getMass()/totalMass); + if (totalMass != 0) + massFactors_[atom->getLocalIndex()] = atom->getMass()/totalMass; else - mfact.push_back( 1.0 ); + massFactors_[atom->getLocalIndex()] = 1.0; } } } @@ -860,15 +884,8 @@ namespace OpenMD { identArray_.push_back(atom->getIdent()); } } - - //fill molMembershipArray - //molMembershipArray is filled by SimCreator - vector molMembershipArray(nGlobalAtoms_); - for (int i = 0; i < nGlobalAtoms_; i++) { - molMembershipArray[i] = globalMolMembership_[i] + 1; - } - //setup fortran simulation + //scan topology nExclude = excludedInteractions_.getSize(); nOneTwo = oneTwoInteractions_.getSize(); @@ -880,72 +897,7 @@ namespace OpenMD { int* oneThreeList = oneThreeInteractions_.getPairList(); int* oneFourList = oneFourInteractions_.getPairList(); - //setFortranSim( &fInfo_, &nGlobalAtoms_, &nAtoms_, &identArray_[0], - // &nExclude, excludeList, - // &nOneTwo, oneTwoList, - // &nOneThree, oneThreeList, - // &nOneFour, oneFourList, - // &molMembershipArray[0], &mfact[0], &nCutoffGroups_, - // &fortranGlobalGroupMembership[0], &isError); - - // if( isError ){ - // - // sprintf( painCave.errMsg, - // "There was an error setting the simulation information in fortran.\n" ); - // painCave.isFatal = 1; - // painCave.severity = OPENMD_ERROR; - // simError(); - //} - - - // sprintf( checkPointMsg, - // "succesfully sent the simulation information to fortran.\n"); - - // errorCheckPoint(); - - // Setup number of neighbors in neighbor list if present - //if (simParams_->haveNeighborListNeighbors()) { - // int nlistNeighbors = simParams_->getNeighborListNeighbors(); - // setNeighbors(&nlistNeighbors); - //} - -#ifdef IS_MPI - // mpiSimData parallelData; - - //fill up mpiSimData struct - // parallelData.nMolGlobal = getNGlobalMolecules(); - // parallelData.nMolLocal = getNMolecules(); - // parallelData.nAtomsGlobal = getNGlobalAtoms(); - // parallelData.nAtomsLocal = getNAtoms(); - // parallelData.nGroupsGlobal = getNGlobalCutoffGroups(); - // parallelData.nGroupsLocal = getNCutoffGroups(); - // parallelData.myNode = worldRank; - // MPI_Comm_size(MPI_COMM_WORLD, &(parallelData.nProcessors)); - - //pass mpiSimData struct and index arrays to fortran - //setFsimParallel(¶llelData, &(parallelData.nAtomsLocal), - // &localToGlobalAtomIndex[0], &(parallelData.nGroupsLocal), - // &localToGlobalCutoffGroupIndex[0], &isError); - - // if (isError) { - // sprintf(painCave.errMsg, - // "mpiRefresh errror: fortran didn't like something we gave it.\n"); - // painCave.isFatal = 1; - // simError(); - // } - - // sprintf(checkPointMsg, " mpiRefresh successful.\n"); - // errorCheckPoint(); -#endif - - // initFortranFF(&isError); - // if (isError) { - // sprintf(painCave.errMsg, - // "initFortranFF errror: fortran didn't like something we gave it.\n"); - // painCave.isFatal = 1; - // simError(); - // } - // fortranInitialized_ = true; + topologyDone_ = true; } void SimInfo::addProperty(GenericData* genData) {