Changeset 1b64b2f

Timestamp:

May 8, 2017, 2:04:19 PM (8 years ago)

Author:

Frederik Heber <frederik.heber@…>

Branches:

ForceAnnealing_goodresults, ForceAnnealing_tocheck

Children:

7ee4b5

Parents:

485eea

git-author:

Frederik Heber <heber@…> (04/06/17 05:09:37)

git-committer:

Frederik Heber <frederik.heber@…> (05/08/17 14:04:19)

Message:

ForceAnnealing now uses step width according to Barzilai-Borwein method.

• this is not as good as what MPQC does but as long as Actions don't have an internal state, we cannot do anything better, e.g. a real line-search.
• removed paramater delta from optimize-structure.
• DOCU: Updated entry in user guide and noted Barzilai-Borwein method.
• TESTS: removed deltat from force annealing regression tests.

Files:

• doc/userguide/userguide.xml (modified) (3 diffs)
• src/Actions/MoleculeAction/ForceAnnealingAction.cpp (modified) (1 diff)
• src/Actions/MoleculeAction/ForceAnnealingAction.def (modified) (1 diff)
• src/Dynamics/ForceAnnealing.hpp (modified) (4 diffs)
• tests/GuiChecks/Molecules/ForceAnnealing/testsuite-molecules-force-annealing.at (modified) (3 diffs)
• tests/regression/Molecules/ForceAnnealing/testsuite-molecules-force-annealing.at (modified) (3 diffs)

doc/userguide/userguide.xml

@@ -2207 +2207 @@
     combines the same Actions as <emphasis role="bold">molecular-dynamics</emphasis> does. However, it
     uses the <emphasis role="bold">force-annealing</emphasis> action instead of <emphasis role="bold">verlet-integration</emphasis>.</para>
+          <remark>Because it is a MacroAction we cannot use any more elaborate 
+          line search method during the optimization. In essence, we minimize the
+          energy function that depends on the nuclei coordinates. The function is
+          non-linear and non-convex. The best choice, given that we have gradient
+          information, would be the Conjugate Gradient method (such as Fletcher-
+          Reeves which work well also for non-linear functions). However, these
+          methods perform a line search along the gradient direction which we 
+          cannot as our Actions do not contain any internal information (apart
+          from the state for Undo/Redo). Therefore, we are restricted to the
+          method of Barzilai-Borwein where no line-search is needed but that also
+          works well for high-dimensional minimization problems.</remark>
           <para>The command below performs a structure optimization of the 
     currently selected atoms (may also be a subset) for up to 100 time
@@ -2219 +2230 @@
           --order 3 \
           --distance 3. \
-          --deltat 0.5 \
           --keep-fixed-CenterOfMass 1 \
           --fragment-executable mpqc \
@@ -2227 +2237 @@
     Otherwise only two time steps would be created: the initial and
     the final one containing the optimized structure.</para>
+        <para>Because of the use of Barzilai-Borwein for computing the 
+        stepwidth  we do not need any 'deltat' parameters. If the computed
+        step width is zero, we use a default step width of 1.</para>
         </section>
         <section xml:id="dynamics.step-world-time">

src/Actions/MoleculeAction/ForceAnnealingAction.cpp

@@ -92 +92 @@
   ForceAnnealing<std::vector<atom *> > optimizer(
       set,
-      params.Deltat.get(),
       true,
       params.steps.get());

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)(double)(unsigned int)(bool)
-#define paramtokens ("forces-file")("deltat")("steps")("output-every-step")
-#define paramdescriptions ("file containing")("time step width")("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(""))(PARAM_DEFAULT(0.1))(NOPARAM_DEFAULT)(PARAM_DEFAULT("0"))
-#define paramreferences (forcesfile)(Deltat)(steps)(DoOutput)
+#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 paramvalids \
 (DummyValidator< boost::filesystem::path >()) \
-(PositiveValidator< double >()) \
 (NotZeroValidator< unsigned int >()) \
 (DummyValidator<bool>())

src/Dynamics/ForceAnnealing.hpp

@@ -32 +32 @@
  *
  * Sadly, we have to use some static instances as so far values cannot be passed
- * between actions. Hence, we need to store the current step and the adaptive
+ * between actions. Hence, we need to store the current step and the adaptive-
  * step width (we cannot perform a linesearch, as we have no control over the
  * calculation of the forces).
@@ -42 +42 @@
   /** Constructor of class ForceAnnealing.
    *
+   * \note We use a fixed delta t of 1.
+   *
    * \param _atoms set of atoms to integrate
    * \param _Deltat time step width in atomic units
@@ -49 +51 @@
   ForceAnnealing(
       AtomSetMixin<T> &_atoms,
-      double _Deltat,
       bool _IsAngstroem,
       const size_t _maxSteps) :
-    AtomicForceManipulator<T>(_atoms, _Deltat, _IsAngstroem),
+    AtomicForceManipulator<T>(_atoms, 1., _IsAngstroem),
     maxSteps(_maxSteps)
   {}
@@ -89 +90 @@
         iter != AtomicForceManipulator<T>::atoms.end(); ++iter) {
       // atom's force vector gives steepest descent direction
+      const Vector oldPosition = (*iter)->getPositionAtStep(NextStep-2 >= 0 ? NextStep - 2 : 0);
       const Vector currentPosition = (*iter)->getPosition();
+      const Vector oldGradient = (*iter)->getAtomicForceAtStep(NextStep-2 >= 0 ? NextStep - 2 : 0);
       const Vector currentGradient = (*iter)->getAtomicForce();
       LOG(4, "DEBUG: Force for atom " << **iter << " is " << currentGradient);
 
-      // artificial update: deltat may be considered as 1/2 s^2 units, mass
-      // is neglected deliberately as this makes all atoms equally fast or
-      // hydrogens slower (and they need to wait for other atoms to arrive at
-      // final position).
-      Vector PositionUpdate = currentDeltat * currentGradient;
+      // we use Barzilai-Borwein update with position reversed to get descent
+      const Vector GradientDifference = (currentGradient - oldGradient);
+      const double stepwidth =
+          fabs((currentPosition - oldPosition).ScalarProduct(GradientDifference))/
+          GradientDifference.NormSquared();
+      Vector PositionUpdate = stepwidth * currentGradient;
+      if (fabs(stepwidth) < 1e-10) {
+        // dont' warn in first step, deltat usage normal
+        if (currentStep != 1)
+          ELOG(1, "INFO: Barzilai-Borwein stepwidth is zero, using deltat " << currentDeltat << " instead.");
+        PositionUpdate = currentDeltat * currentGradient;
+      }
       LOG(3, "DEBUG: Update would be " << PositionUpdate);
 

tests/GuiChecks/Molecules/ForceAnnealing/testsuite-molecules-force-annealing.at

@@ -25 +25 @@
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/pre/test.forces .], 0)
 AT_CHECK([chmod u+w $file], 0)
-AT_CHECK([../../molecuilder --dry-run -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --deltat 0.01 --no-dry-run --store-session session-molecules-force-annealing.py --session-type python], 0, [stdout], [stderr])
+AT_CHECK([../../molecuilder --dry-run -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --no-dry-run --store-session session-molecules-force-annealing.py --session-type python], 0, [stdout], [stderr])
 AT_CHECK([grep -v "Command.*DryRun" session-molecules-force-annealing.py >session-molecules-force-annealing_new.py], 0, [ignore], [ignore])
 AT_CHECK([../../molecuilderguitest session-molecules-force-annealing_new.py], 0, [stdout], [stderr])
@@ -40 +40 @@
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/pre/test.forces .], 0)
 AT_CHECK([chmod u+w $file], 0)
-AT_CHECK([../../molecuilder --dry-run -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --deltat 0.01 --undo --no-dry-run --store-session session-molecules-force-annealing.py --session-type python], 0, [stdout], [stderr])
+AT_CHECK([../../molecuilder --dry-run -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --undo --no-dry-run --store-session session-molecules-force-annealing.py --session-type python], 0, [stdout], [stderr])
 AT_CHECK([grep -v "Command.*DryRun" session-molecules-force-annealing.py >session-molecules-force-annealing_new.py], 0, [ignore], [ignore])
 AT_CHECK([../../molecuilderguitest session-molecules-force-annealing_new.py], 0, [stdout], [stderr])
@@ -55 +55 @@
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/pre/test.forces .], 0)
 AT_CHECK([chmod u+w $file], 0)
-AT_CHECK([../../molecuilder --dry-run -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --deltat 0.01 --undo --redo --no-dry-run --store-session session-molecules-force-annealing.py --session-type python], 0, [stdout], [stderr])
+AT_CHECK([../../molecuilder --dry-run -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --undo --redo --no-dry-run --store-session session-molecules-force-annealing.py --session-type python], 0, [stdout], [stderr])
 AT_CHECK([grep -v "Command.*DryRun" session-molecules-force-annealing.py >session-molecules-force-annealing_new.py], 0, [ignore], [ignore])
 AT_CHECK([../../molecuilderguitest session-molecules-force-annealing_new.py], 0, [stdout], [stderr])

tests/regression/Molecules/ForceAnnealing/testsuite-molecules-force-annealing.at

@@ -25 +25 @@
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/pre/test.forces .], 0)
 AT_CHECK([chmod u+w $file], 0)
-AT_CHECK([../../molecuilder -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --deltat 0.01], 0, [stdout], [stderr])
+AT_CHECK([../../molecuilder -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1], 0, [stdout], [stderr])
 AT_CHECK([diff $file ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/post/test.conf], 0, [ignore], [ignore])
 
@@ -38 +38 @@
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/pre/test.forces .], 0)
 AT_CHECK([chmod u+w $file], 0)
-AT_CHECK([../../molecuilder -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --deltat 0.01 --undo], 0, [stdout], [stderr])
+AT_CHECK([../../molecuilder -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --undo], 0, [stdout], [stderr])
 AT_CHECK([diff $file ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/post/test-undo.conf], 0, [ignore], [ignore])
 
@@ -51 +51 @@
 AT_CHECK([/bin/cp -f ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/pre/test.forces .], 0)
 AT_CHECK([chmod u+w $file], 0)
-AT_CHECK([../../molecuilder -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --deltat 0.01 --undo --redo], 0, [stdout], [stderr])
+AT_CHECK([../../molecuilder -i $file --select-all-atoms --force-annealing --forces-file test.forces --steps 1 --undo --redo], 0, [stdout], [stderr])
 AT_CHECK([diff $file ${abs_top_srcdir}/tests/regression/Molecules/ForceAnnealing/post/test.conf], 0, [ignore], [ignore])