ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/flucq/FluctuatingChargeConstraints.cpp
(Generate patch)

Comparing trunk/src/flucq/FluctuatingChargeConstraints.cpp (file contents):
Revision 1912 by gezelter, Fri Jul 19 21:25:45 2013 UTC vs.
Revision 1913 by jmichalk, Wed Jul 24 20:00:51 2013 UTC

# Line 50 | Line 50 | namespace OpenMD {
50   namespace OpenMD {
51  
52    FluctuatingChargeConstraints::FluctuatingChargeConstraints(SimInfo* info) :
53 <    info_(info), constrainRegions_(false), hasFlucQ_(false) {
53 >    info_(info), constrainRegions_(false), hasFlucQ_(false), initialized_(false) {
54      
55 <    if (info_->usesFluctuatingCharges()) {
56 <      if (info_->getNFluctuatingCharges() > 0) {
57 <        hasFlucQ_ = true;
55 >  }
56 >
57 >  void FluctuatingChargeConstraints::initialize(){
58 >    if(info_->usesFluctuatingCharges()){
59 >      if(info_->getNFluctuatingCharges() > 0){
60 >        hasFlucQ_ = true;
61        }
62      }
63 +    initialized_ = true;
64    }
65  
66    void FluctuatingChargeConstraints::setConstrainRegions(bool cr) {
# Line 75 | Line 79 | namespace OpenMD {
79             mol = info_->nextMolecule(i)) {
80          reg = mol->getRegion();
81          if (reg >= 0) regions.insert(reg);
78 –        
82        }
83 +      // resize the keys vector to the largest found value for regions.
84 +      regionKeys_.resize( *(regions.end()) );
85 +      int which = 0;
86        for (std::set<int>::iterator r=regions.begin(); r!=regions.end(); ++r) {
87 <        regionKeys_.push_back( (*r) );
87 >        regionKeys_[ (*r) ] = which;
88 >        which++;
89        }
90        regionForce_.resize( regionKeys_.size() );
91        regionCharges_.resize( regionKeys_.size() );
# Line 87 | Line 94 | namespace OpenMD {
94  
95  
96    void FluctuatingChargeConstraints::applyConstraints() {
97 +    if (!initialized_) initialize();
98      if (!hasFlucQ_) return;
91 –    
99      SimInfo::MoleculeIterator i;
100      Molecule::FluctuatingChargeIterator  j;
101      Molecule* mol;
102      Atom* atom;
103      
104      RealType totalFrc, totalMolFrc, regionFrc, constrainedFrc;
98 –    
105      // accumulate the total system fluctuating charge forces
106      totalFrc = 0.0;
107      if (constrainRegions_) {
# Line 114 | Line 120 | namespace OpenMD {
120          RealType frc = atom->getFlucQFrc();
121          totalFrc += frc;
122          if (constrainRegions_) {
123 <          regionForce_[regionKeys_[region]] += frc;
124 <          regionCharges_[regionKeys_[region]] += 1;
123 >          if (region >= 0) {        
124 >            regionForce_[regionKeys_[region]] += frc;
125 >            regionCharges_[regionKeys_[region]] += 1;
126 >          }
127          }
128        }
129      }
# Line 150 | Line 158 | namespace OpenMD {
158        
159        if (constrainRegions_) {
160          int region = mol->getRegion();
161 <        regionFrc = regionForce_[regionKeys_[region]];
161 >        if (region >= 0)
162 >          regionFrc = regionForce_[regionKeys_[region]];
163 >        else
164 >          regionFrc = 0.0;
165        } else {
166          regionFrc = 0.0;
167        }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines