--- trunk/src/brains/SimInfo.cpp 2005/09/16 16:07:39 603 +++ trunk/src/brains/SimInfo.cpp 2005/10/04 19:34:03 645 @@ -53,8 +53,9 @@ #include "math/Vector3.hpp" #include "primitives/Molecule.hpp" #include "UseTheForce/fCutoffPolicy.h" -#include "UseTheForce/fCoulombicCorrection.h" +#include "UseTheForce/DarkSide/fElectrostaticSummationMethod.h" #include "UseTheForce/doForces_interface.h" +#include "UseTheForce/DarkSide/electrostatic_interface.h" #include "UseTheForce/notifyCutoffs_interface.h" #include "utils/MemoryUtils.hpp" #include "utils/simError.h" @@ -82,7 +83,7 @@ namespace oopse { MoleculeStamp* molStamp; int nMolWithSameStamp; int nCutoffAtoms = 0; // number of atoms belong to cutoff groups - int nGroups = 0; //total cutoff groups defined in meta-data file + int nGroups = 0; //total cutoff groups defined in meta-data file CutoffGroupStamp* cgStamp; RigidBodyStamp* rbStamp; int nRigidAtoms = 0; @@ -107,6 +108,7 @@ namespace oopse { } nGroups += nCutoffGroupsInStamp * nMolWithSameStamp; + nCutoffAtoms += nAtomsInGroups * nMolWithSameStamp; //calculate atoms in rigid bodies @@ -123,18 +125,21 @@ namespace oopse { } - //every free atom (atom does not belong to cutoff groups) is a cutoff group - //therefore the total number of cutoff groups in the system is 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 + //every free atom (atom does not belong to cutoff groups) is a cutoff + //group therefore the total number of cutoff groups in the system is + //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 nGlobalCutoffGroups_ = nGlobalAtoms_ - nCutoffAtoms + nGroups; - //every free atom (atom does not belong to rigid bodies) is an integrable object - //therefore the total number of integrable objects in the system is equal to - //the total number of atoms minus number of atoms belong to rigid body defined in meta-data - //file plus the number of rigid bodies defined in meta-data file - nGlobalIntegrableObjects_ = nGlobalAtoms_ - nRigidAtoms + nGlobalRigidBodies_; - + //every free atom (atom does not belong to rigid bodies) is an + //integrable object therefore the total number of integrable objects + //in the system is equal to the total number of atoms minus number of + //atoms belong to rigid body defined in meta-data file plus the number + //of rigid bodies defined in meta-data file + nGlobalIntegrableObjects_ = nGlobalAtoms_ - nRigidAtoms + + nGlobalRigidBodies_; + nGlobalMols_ = molStampIds_.size(); #ifdef IS_MPI @@ -522,6 +527,14 @@ namespace oopse { int useElectrostatics = 0; //usePBC and useRF are from simParams int usePBC = simParams_->getPBC(); + int useRF; + + // set the useRF logical + std::string myMethod = simParams_->getElectrostaticSummationMethod(); + if (myMethod == "REACTION_FIELD") + useRF = 1; + else + useRF = 0; //loop over all of the atom types for (i = atomTypes.begin(); i != atomTypes.end(); ++i) { @@ -584,6 +597,9 @@ namespace oopse { temp = useFLARB; MPI_Allreduce(&temp, &useFLARB, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); + temp = useRF; + MPI_Allreduce(&temp, &useRF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); + #endif fInfo_.SIM_uses_PBC = usePBC; @@ -598,8 +614,9 @@ namespace oopse { fInfo_.SIM_uses_EAM = useEAM; fInfo_.SIM_uses_Shapes = useShape; fInfo_.SIM_uses_FLARB = useFLARB; + fInfo_.SIM_uses_RF = useRF; - if( fInfo_.SIM_uses_Dipoles && fInfo_.SIM_uses_RF) { + if( fInfo_.SIM_uses_Dipoles && myMethod == "REACTION_FIELD") { if (simParams_->haveDielectric()) { fInfo_.dielect = simParams_->getDielectric(); @@ -649,7 +666,11 @@ namespace oopse { totalMass = cg->getMass(); for(atom = cg->beginAtom(ai); atom != NULL; atom = cg->nextAtom(ai)) { - mfact.push_back(atom->getMass()/totalMass); + // Check for massless groups - set mfact to 1 if true + if (totalMass != 0) + mfact.push_back(atom->getMass()/totalMass); + else + mfact.push_back( 1.0 ); } } @@ -856,7 +877,15 @@ namespace oopse { } } } + + + if (simParams_->haveSkinThickness()) { + double skinThickness = simParams_->getSkinThickness(); + } + notifyFortranCutoffs(&rcut_, &rsw_, &rnblist, &cp); + // also send cutoff notification to electrostatics + setElectrostaticCutoffRadius(&rcut_); } void SimInfo::setupElectrostaticSummationMethod( int isError ) { @@ -864,11 +893,14 @@ namespace oopse { int errorOut; int esm = NONE; double alphaVal; + double dielectric; errorOut = isError; + alphaVal = simParams_->getDampingAlpha(); + dielectric = simParams_->getDielectric(); if (simParams_->haveElectrostaticSummationMethod()) { - std::string myCorrection = simParams_->getElectrostaticSummationMethod(); + std::string myMethod = simParams_->getElectrostaticSummationMethod(); if (myMethod == "NONE") { esm = NONE; } else { @@ -876,15 +908,14 @@ namespace oopse { esm = UNDAMPED_WOLF; } else { if (myMethod == "DAMPED_WOLF") { - esm = WOLF; + esm = DAMPED_WOLF; if (!simParams_->haveDampingAlpha()) { //throw error sprintf( painCave.errMsg, - "SimInfo warning: dampingAlpha was not specified in the input file. A default value of %f (1/ang) will be used for the Damped Wolf Method.", simParams_->getDampingAlpha()); + "SimInfo warning: dampingAlpha was not specified in the input file. A default value of %f (1/ang) will be used for the Damped Wolf Method.", alphaVal); painCave.isFatal = 0; simError(); } - alphaVal = simParams_->getDampingAlpha(); } else { if (myMethod == "REACTION_FIELD") { esm = REACTION_FIELD; @@ -899,7 +930,11 @@ namespace oopse { } } } - initFortranFF( &fInfo_.SIM_uses_RF, &esm, &alphaVal, &errorOut ); + // let's pass some summation method variables to fortran + setElectrostaticSummationMethod( &esm ); + setDampedWolfAlpha( &alphaVal ); + setReactionFieldDielectric( &dielectric ); + initFortranFF( &esm, &errorOut ); } void SimInfo::addProperty(GenericData* genData) {