--- branches/development/src/brains/SimInfo.cpp 2011/09/13 22:05:04 1627 +++ branches/development/src/brains/SimInfo.cpp 2012/06/05 18:07:08 1744 @@ -36,7 +36,8 @@ * [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). - * [4] Vardeman & Gezelter, in progress (2009). + * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). + * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). */ /** @@ -58,7 +59,7 @@ #include "utils/simError.h" #include "selection/SelectionManager.hpp" #include "io/ForceFieldOptions.hpp" -#include "UseTheForce/ForceField.hpp" +#include "brains/ForceField.hpp" #include "nonbonded/SwitchingFunction.hpp" #ifdef IS_MPI #include @@ -71,10 +72,10 @@ 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), + nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), nGlobalFluctuatingCharges_(0), nAtoms_(0), nBonds_(0), nBends_(0), nTorsions_(0), nInversions_(0), nRigidBodies_(0), nIntegrableObjects_(0), nCutoffGroups_(0), - nConstraints_(0), sman_(NULL), topologyDone_(false), + nConstraints_(0), nFluctuatingCharges_(0), sman_(NULL), topologyDone_(false), calcBoxDipole_(false), useAtomicVirial_(true) { MoleculeStamp* molStamp; @@ -224,13 +225,17 @@ namespace OpenMD { void SimInfo::calcNdf() { - int ndf_local; + int ndf_local, nfq_local; MoleculeIterator i; vector::iterator j; + vector::iterator k; + Molecule* mol; StuntDouble* integrableObject; + Atom* atom; ndf_local = 0; + nfq_local = 0; for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { for (integrableObject = mol->beginIntegrableObject(j); integrableObject != NULL; @@ -245,17 +250,27 @@ namespace OpenMD { ndf_local += 3; } } - } + for (atom = mol->beginFluctuatingCharge(k); atom != NULL; + atom = mol->nextFluctuatingCharge(k)) { + if (atom->isFluctuatingCharge()) { + nfq_local++; + } + } } + ndfLocal_ = ndf_local; + cerr << "ndfLocal_ = " << ndfLocal_ << "\n"; + // n_constraints is local, so subtract them on each processor 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); #else ndf_ = ndf_local; + nGlobalFluctuatingCharges_ = nfq_local; #endif // nZconstraints_ is global, as are the 3 COM translations for the @@ -779,11 +794,13 @@ namespace OpenMD { int usesElectrostatic = 0; int usesMetallic = 0; int usesDirectional = 0; + int usesFluctuatingCharges = 0; //loop over all of the atom types for (i = atomTypes.begin(); i != atomTypes.end(); ++i) { usesElectrostatic |= (*i)->isElectrostatic(); usesMetallic |= (*i)->isMetal(); usesDirectional |= (*i)->isDirectional(); + usesFluctuatingCharges |= (*i)->isFluctuatingCharge(); } #ifdef IS_MPI @@ -796,11 +813,15 @@ namespace OpenMD { temp = usesElectrostatic; MPI_Allreduce(&temp, &usesElectrostaticAtoms_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); + + temp = usesFluctuatingCharges; + MPI_Allreduce(&temp, &usesFluctuatingCharges_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); #else usesDirectionalAtoms_ = usesDirectional; usesMetallicAtoms_ = usesMetallic; usesElectrostaticAtoms_ = usesElectrostatic; + usesFluctuatingCharges_ = usesFluctuatingCharges; #endif @@ -1193,7 +1214,7 @@ namespace OpenMD { det = intTensor.determinant(); sysconstants = geomCnst/(RealType)nGlobalIntegrableObjects_; - volume = 4.0/3.0*NumericConstant::PI*pow(sysconstants,3.0/2.0)*sqrt(det); + volume = 4.0/3.0*NumericConstant::PI*pow(sysconstants,geomCnst)*sqrt(det); return; } @@ -1209,7 +1230,7 @@ namespace OpenMD { detI = intTensor.determinant(); sysconstants = geomCnst/(RealType)nGlobalIntegrableObjects_; - volume = 4.0/3.0*NumericConstant::PI*pow(sysconstants,3.0/2.0)*sqrt(detI); + volume = 4.0/3.0*NumericConstant::PI*pow(sysconstants,geomCnst)*sqrt(detI); return; } /*