From 142d115d31f9c3589188ebd64ca31bb836351f25 Mon Sep 17 00:00:00 2001
From: sjplimp <sjplimp@f3b2605a-c512-4ea7-a41b-209d697bcdaa>
Date: Fri, 16 Feb 2007 14:50:22 +0000
Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@311
 f3b2605a-c512-4ea7-a41b-209d697bcdaa

---
 src/min_cg.cpp | 22 ++++++++++++++--------
 1 file changed, 14 insertions(+), 8 deletions(-)

diff --git a/src/min_cg.cpp b/src/min_cg.cpp
index 44abd303e7..f21a98825b 100644
--- a/src/min_cg.cpp
+++ b/src/min_cg.cpp
@@ -570,8 +570,8 @@ int MinCG::linemin_scan(int n, double *x, double *dir, double eng,
 /* ----------------------------------------------------------------------
    linemin: use secant approximation to estimate parabola minimum at each step
    should converge more quickly/accurately than "scan", but may be less robust
-   how to prevent func evals at new x bigger than maxdist?
    initial secant from two points: 0 and sigma0 = mindist
+   prevents successvive func evals further apart in x than maxdist
 ------------------------------------------------------------------------- */
 
 int MinCG::linemin_secant(int n, double *x, double *dir, double eng,
@@ -579,20 +579,18 @@ int MinCG::linemin_secant(int n, double *x, double *dir, double eng,
 			  double &alpha, int &nfunc)
 {
   int i,iter;
-  double eta,eta_prev,sigma0,alphadelta,fme,fmax,dsq,e0,tmp;
+  double eta,eta_prev,sigma0,sigmamax,alphadelta,fme,fmax,dsq,e0,tmp;
   double *f;
   double epssq = SECANT_EPS * SECANT_EPS;
 
   // stopping criterion for secant iterations
-  // change this
 
   fme = 0.0;
   for (i = 0; i < n; i++) fme += dir[i]*dir[i];
   MPI_Allreduce(&fme,&dsq,1,MPI_DOUBLE,MPI_SUM,world);
 
   // sigma0 = smallest allowed step of mindist
-  // eval func at sigma0
-  // test if minstep is already uphill
+  // sigmamax = largest allowed step (in single iteration) of maxdist
 
   fme = 0.0;
   for (i = 0; i < n; i++) fme = MAX(fme,fabs(dir[i]));
@@ -600,6 +598,10 @@ int MinCG::linemin_secant(int n, double *x, double *dir, double eng,
   if (fmax == 0.0) return 1;
 
   sigma0 = mindist/fmax;
+  sigmamax = maxdist/fmax;
+
+  // eval func at sigma0
+  // test if minstep is already uphill
 
   e0 = eng;
   for (i = 0; i < n; i++) x[i] += sigma0*dir[i];
@@ -613,14 +615,14 @@ int MinCG::linemin_secant(int n, double *x, double *dir, double eng,
     return 1;
   }
 
+  // secant iterations
+  // alphadelta = new increment to move, alpha = accumulated move
+
   f = atom->f[0];
   tmp = 0.0;
   for (i = 0; i < n; i++) tmp -= f[i]*dir[i];
   MPI_Allreduce(&tmp,&eta_prev,1,MPI_DOUBLE,MPI_SUM,world);
 
-  // secant iterations
-  // alphadelta = new increment to move, alpha = accumulated move
-
   alpha = sigma0;
   alphadelta = -sigma0;
 
@@ -638,6 +640,10 @@ int MinCG::linemin_secant(int n, double *x, double *dir, double eng,
     alphadelta *= eta / (eta_prev - eta);
     eta_prev = eta;
     if (alphadelta*alphadelta*dsq <= epssq) break;
+    if (fabs(alphadelta) > sigmamax) {
+      if (alphadelta > 0.0) alphadelta = sigmamax;
+      else alphadelta = -sigmamax;
+    }
   }
 
   // if exited loop on first iteration, func eval was at alpha = 0.0