Hasan,

I am delighted that someone else is willing to take a close look at these

expressions. They originate with the idea that in the linesearch, the

objective function is approximated by a quadratic function in alpha, for

sufficiently small alpha. This idea is the same as that used in the

well-known secant method. However, since the change in the objective

function (difference of two finite numbers) is not known as accurately as

the gradient (which is close to zero), all the expressions are written in

terms of gradients. In this way, we can converge the LAMMPS forces much

closer to zero.

So much for the math. In practice, many of the details were arrived in a

highly Edisonian manner. The final code does seem to work in most practical

cases, although it may be less robust that the backtrack linesearch in

situations where the forces do not exactly represent dE/dx, which can occur

for some potentials with discontinuities

Let's see if I can answer your detailed questions:

You are right. The expression for relerr is wrong. It should be:

de_ideal = -BACKTRACK_SLOPE*alpha*fdothall;

if (de_ideal >= -IDEAL_TOL) {

relerr = 0.0;

} else {

relerr = fabs(1.0+0.5*(alpha-alphaprev)*(fh+fhprev)/(ecurrent-engprev));

}

The derivation of this is given below. For whatever reason, both version

seem to work equally well. We are going to make further improvements to the

minimizer soon and will include that change.

The motivation behind IDEAL_TOL is that if the ideal energy change is too

small, then we do not expect to be able to measure the real energy reliably,

due to machine precision. IDEAL_TOL defines what we think is too small.

Finally, if you are trying to relax a crystal with fix box/relax aniso, with

different target pressures in x, y, and z, it is not unusual to encounter

difficulties. There are several reason for this, but it relates to the fact

that this case, unlike the isotropic case, can not be posed as a

minimization problem, unless the gradients include second-derivatives of

energy. We have some ideas for improvement, but have not tried them out yet.

Aidan

*** relerr derivation ***

We know E,Eprev,fh,fhprev. The Taylor series about alpha_prev truncated at

the quadratic term is:

E = Eprev - del_alpha*fhprev + (1/2)del_alpha^2*Hprev

and

fh = -fhprev + del_alpha*Hprev

where del_alpha = alpha-alpha_prev

The only unknown is Hprev, but we have two equations. So we can choose one

of E,Eprev,fh,fhprev as the second unknown, and compare the solution with

the known value. E seems like the best choice, giving:

Esolve = Eprev - del_alpha*(f+fprev)/2

relerr follows from that.

relerr = (Esolve-E)/(E-Eprev)

= -0.5*del_alpha*(f+fprev)/(E-Eprev) - 1