ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/brains/ForceManager.cpp
(Generate patch)

Comparing branches/development/src/brains/ForceManager.cpp (file contents):
Revision 1767 by gezelter, Fri Jul 6 22:01:58 2012 UTC vs.
Revision 1847 by gezelter, Thu Jan 31 17:57:07 2013 UTC

# Line 44 | Line 44
44   * @file ForceManager.cpp
45   * @author tlin
46   * @date 11/09/2004
47 – * @time 10:39am
47   * @version 1.0
48   */
49  
# Line 58 | Line 57
57   #include "primitives/Torsion.hpp"
58   #include "primitives/Inversion.hpp"
59   #include "nonbonded/NonBondedInteraction.hpp"
60 + #include "perturbations/ElectricField.hpp"
61   #include "parallel/ForceMatrixDecomposition.hpp"
62  
63   #include <cstdio>
# Line 369 | Line 369 | namespace OpenMD {
369      }
370      switcher_->setSwitchType(sft_);
371      switcher_->setSwitch(rSwitch_, rCut_);
372 –    interactionMan_->setSwitchingRadius(rSwitch_);
372    }
373  
374  
# Line 393 | Line 392 | namespace OpenMD {
392        doParticlePot_ = info_->getSimParams()->getOutputParticlePotential();
393        doHeatFlux_ = info_->getSimParams()->getPrintHeatFlux();
394        if (doHeatFlux_) doParticlePot_ = true;
395 +
396 +      doElectricField_ = info_->getSimParams()->getOutputElectricField();
397    
398      }
399  
# Line 423 | Line 424 | namespace OpenMD {
424      electrostaticScale_[2] = fopts.getelectrostatic13scale();
425      electrostaticScale_[3] = fopts.getelectrostatic14scale();    
426      
427 +    if (info_->getSimParams()->haveElectricField()) {
428 +      ElectricField* eField = new ElectricField(info_);
429 +      perturbations_.push_back(eField);
430 +    }
431 +
432      fDecomp_->distributeInitialData();
433  
434      initialized_ = true;
# Line 681 | Line 687 | namespace OpenMD {
687      potVec longRangePotential(0.0);
688      potVec workPot(0.0);
689      potVec exPot(0.0);
690 +    Vector3d eField1(0.0);
691 +    Vector3d eField2(0.0);
692      vector<int>::iterator ia, jb;
693  
694      int loopStart, loopEnd;
# Line 694 | Line 702 | namespace OpenMD {
702      idat.vpair = &vpair;
703      idat.dVdFQ1 = &dVdFQ1;
704      idat.dVdFQ2 = &dVdFQ2;
705 +    idat.eField1 = &eField1;
706 +    idat.eField2 = &eField2;  
707      idat.f1 = &f1;
708      idat.sw = &sw;
709      idat.shiftedPot = (cutoffMethod_ == SHIFTED_POTENTIAL) ? true : false;
710      idat.shiftedForce = (cutoffMethod_ == SHIFTED_FORCE) ? true : false;
711      idat.doParticlePot = doParticlePot_;
712 +    idat.doElectricField = doElectricField_;
713      sdat.doParticlePot = doParticlePot_;
714      
715      loopEnd = PAIR_LOOP;
# Line 734 | Line 745 | namespace OpenMD {
745            if (iLoop == PAIR_LOOP) {
746              vij = 0.0;
747              fij = V3Zero;
748 +            eField1 = V3Zero;
749 +            eField2 = V3Zero;
750            }
751            
752            in_switching_region = switcher_->getSwitch(rgrpsq, sw, dswdr,
# Line 808 | Line 821 | namespace OpenMD {
821                fij += fg;
822  
823                if (atomListRow.size() == 1 && atomListColumn.size() == 1) {
824 <                stressTensor -= outProduct( *(idat.d), fg);
825 <                if (doHeatFlux_)
826 <                  fDecomp_->addToHeatFlux(*(idat.d) * dot(fg, vel2));
827 <                
824 >                if (!fDecomp_->skipAtomPair(atomListRow[0],
825 >                                            atomListColumn[0],
826 >                                            cg1, cg2)) {
827 >                  stressTensor -= outProduct( *(idat.d), fg);
828 >                  if (doHeatFlux_)
829 >                    fDecomp_->addToHeatFlux(*(idat.d) * dot(fg, vel2));
830 >                }                
831                }
832            
833                for (ia = atomListRow.begin();
# Line 896 | Line 912 | namespace OpenMD {
912        *(fDecomp_->getPairwisePotential());
913  
914      curSnapshot->setLongRangePotential(longRangePotential);
915 +
916 +    // collects single-atom information
917 +    fDecomp_->collectSelfData();
918 +
919 +    longRangePotential = *(fDecomp_->getEmbeddingPotential()) +
920 +      *(fDecomp_->getPairwisePotential());
921 +
922 +    curSnapshot->setLongRangePotential(longRangePotential);
923      
924      curSnapshot->setExcludedPotentials(*(fDecomp_->getExcludedSelfPotential()) +
925                                           *(fDecomp_->getExcludedPotential()));
# Line 904 | Line 928 | namespace OpenMD {
928  
929    
930    void ForceManager::postCalculation() {
931 +
932 +    vector<Perturbation*>::iterator pi;
933 +    for (pi = perturbations_.begin(); pi != perturbations_.end(); ++pi) {
934 +      (*pi)->applyPerturbation();
935 +    }
936 +
937      SimInfo::MoleculeIterator mi;
938      Molecule* mol;
939      Molecule::RigidBodyIterator rbIter;
# Line 928 | Line 958 | namespace OpenMD {
958      curSnapshot->setStressTensor(stressTensor);
959      
960    }
931 –
961   } //end namespace OpenMD

Diff Legend

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