Changeset b03d7d

Timestamp:

May 23, 2017, 8:58:02 PM (8 years ago)

Author:

Frederik Heber <frederik.heber@…>

Branches:

ForceAnnealing_tocheck

Children:

175543

Parents:

9190e6

Message:

tempcommit: Merge with 9190e6a

Location:

src

Files:

• Actions/MoleculeAction/ForceAnnealingAction.cpp (modified) (1 diff)
• Actions/MoleculeAction/ForceAnnealingAction.def (modified) (1 diff)
• Dynamics/ForceAnnealing.hpp (modified) (8 diffs)

src/Actions/MoleculeAction/ForceAnnealingAction.cpp

@@ -93 +93 @@
       set,
       true,
-      params.steps.get());
+      params.steps.get(),
+      params.MaxDistance.get());
   // parse forces into next step
   if (!params.forcesfile.get().string().empty()) {

src/Actions/MoleculeAction/ForceAnnealingAction.def

@@ -16 +16 @@
 // ValueStorage by the token "Z" -> first column: int, Z, "Z"
 // "undefine" if no parameters are required, use (NOPARAM_DEFAULT) for each (undefined) default value
-#define paramtypes (boost::filesystem::path)(unsigned int)(bool)
-#define paramtokens ("forces-file")("steps")("output-every-step")
-#define paramdescriptions ("file containing")("fixed number of optimization steps to be performed")("whether WorldTime should be increased and output written after every step, useful if optimization might hang")
-#define paramdefaults (PARAM_DEFAULT(""))(NOPARAM_DEFAULT)(PARAM_DEFAULT("0"))
-#define paramreferences (forcesfile)(steps)(DoOutput)
+#define paramtypes (boost::filesystem::path)(unsigned int)(bool)(int)
+#define paramtokens ("forces-file")("steps")("output-every-step")("max-distance")
+#define paramdescriptions ("file containing")("fixed number of optimization steps to be performed")("whether WorldTime should be increased and output written after every step, useful if optimization might hang")("maximum distance to which bond graph is taken into account")
+#define paramdefaults (PARAM_DEFAULT(""))(NOPARAM_DEFAULT)(PARAM_DEFAULT("0"))(PARAM_DEFAULT(0))
+#define paramreferences (forcesfile)(steps)(DoOutput)(MaxDistance)
 #define paramvalids \
 (DummyValidator< boost::filesystem::path >()) \
 (NotZeroValidator< unsigned int >()) \
-(DummyValidator<bool>())
+(DummyValidator<bool>()) \
+(DummyValidator< int >())
 
 #define statetypes (std::vector<AtomicInfo>)(std::vector<AtomicInfo>)

src/Dynamics/ForceAnnealing.hpp

@@ -20 +20 @@
 #include "CodePatterns/Log.hpp"
 #include "CodePatterns/Verbose.hpp"
+#include "Descriptors/AtomIdDescriptor.hpp"
 #include "Dynamics/AtomicForceManipulator.hpp"
 #include "Fragmentation/ForceMatrix.hpp"
@@ -38 +39 @@
  * step width (we cannot perform a line search, as we have no control over the
  * calculation of the forces).
+ *
+ * However, we do use the bond graph, i.e. if a single atom needs to be shifted
+ * to the left, then the whole molecule left of it is shifted, too. This is
+ * controlled by the \a max_distance parameter.
  */
 template <class T>
@@ -51 +56 @@
    * \param _IsAngstroem whether length units are in angstroem or bohr radii
    * \param _maxSteps number of optimization steps to perform
+   * \param _max_distance up to this bond order is bond graph taken into account.
    */
   ForceAnnealing(
       AtomSetMixin<T> &_atoms,
       bool _IsAngstroem,
-      const size_t _maxSteps) :
+      const size_t _maxSteps,
+      const int _max_distance) :
     AtomicForceManipulator<T>(_atoms, 1., _IsAngstroem),
-    maxSteps(_maxSteps)
+    maxSteps(_maxSteps),
+    max_distance(_max_distance),
+    damping_factor(0.5)
   {}
   /** Destructor of class ForceAnnealing.
@@ -130 +139 @@
       LOG(3, "DEBUG: Update would be " << PositionUpdate);
 
-      // add update to central atom
-      const atomId_t atomid = (*iter)->getId();
-      if (GatheredUpdates.count(atomid)) {
-        GatheredUpdates[atomid] += PositionUpdate;
-      } else
-        GatheredUpdates.insert( std::make_pair(atomid, PositionUpdate) );
+//      // add update to central atom
+//      const atomId_t atomid = (*iter)->getId();
+//      if (GatheredUpdates.count(atomid)) {
+//        GatheredUpdates[atomid] += PositionUpdate;
+//      } else
+//        GatheredUpdates.insert( std::make_pair(atomid, PositionUpdate) );
 
       // We assume that a force is local, i.e. a bond is too short yet and hence
@@ -145 +154 @@
 
       /// get all nodes from bonds in the direction of the current force
-      const size_t max_distance = 4;
-      const double damping_factor = 0.9;
 
       // remove edges facing in the wrong direction
@@ -184 +191 @@
             "ForceAnnealing() - could not find distance to an atom.");
         const double factor = pow(damping_factor, diter->second);
+        LOG(3, "DEBUG: Update for atom #" << *setiter << " will be "
+            << factor << "*" << PositionUpdate);
         if (GatheredUpdates.count((*setiter))) {
           GatheredUpdates[(*setiter)] += factor*PositionUpdate;
@@ -199 +208 @@
           maxComponents[i] = currentGradient[i];
 
-      // are we in initial step? Then don't check against velocity
-      if ((currentStep > 1) && (!(*iter)->getAtomicVelocity().IsZero()))
-        // update with currentDelta tells us how the current gradient relates to
-        // the last one: If it has become larger, reduce currentDelta
-        if ((PositionUpdate.ScalarProduct((*iter)->getAtomicVelocity()) < 0)
-            && (currentDeltat > MinimumDeltat)) {
-          currentDeltat = .5*currentDeltat;
-          LOG(2, "DEBUG: Upgrade in other direction: " << PositionUpdate.NormSquared()
-              << " > " << (*iter)->getAtomicVelocity().NormSquared()
-              << ", decreasing deltat: " << currentDeltat);
-          PositionUpdate = currentDeltat * currentGradient;
-        }
-
-      // finally set new values
-      (*iter)->setPosition(currentPosition + PositionUpdate);
-      (*iter)->setAtomicVelocity(PositionUpdate);
-      //std::cout << "Id of atom is " << (*iter)->getId() << std::endl;
-//        (*iter)->VelocityVerletUpdateU((*iter)->getId(), NextStep-1, Deltat, IsAngstroem);
-
       // reset force vector for next step except on final one
       if (currentStep != maxSteps)
         (*iter)->setAtomicForce(zeroVec);
+    }
+    // apply the gathered updates
+    for (std::map<atomId_t, Vector>::const_iterator iter = GatheredUpdates.begin();
+        iter != GatheredUpdates.end(); ++iter) {
+      const atomId_t &atomid = iter->first;
+      const Vector &update = iter->second;
+      atom* const walker = World::getInstance().getAtom(AtomById(atomid));
+      ASSERT( walker != NULL,
+          "ForceAnnealing() - walker with id "+toString(atomid)+" has suddenly disappeared.");
+      walker->setPosition( walker->getPosition() + update );
+      walker->setAtomicVelocity(update);
+      LOG(3, "DEBUG: Applying update " << update << " to atom " << *walker);
     }
 
@@ -255 +257 @@
   //!> minimum deltat for internal while loop (adaptive step width)
   static double MinimumDeltat;
+  //!> contains the maximum bond graph distance up to which shifts of a single atom are spread
+  const int max_distance;
+  //!> the shifted is dampened by this factor with the power of the bond graph distance to the shift causing atom
+  const double damping_factor;
 };