Subject: Re: Cartesian vs Internal coordinates

From: [email-masked] (Ron Shepard)
Date: Thu, 14 Jan 93 15:25:59 CST

John Upham writes:
>>Dear All,
>>        I have always assumed that given the same starting geometry any
>>program that has a choice as to how one enters it should (in the end)
>>produce the same optimised final geometry. Is this true or false ?
>>Without wishing to embarasses anyone I have at least one program 
>>genearally available where different final geometries are obtained.
>>  
>>Any comments (apart form use internal coordinates) ?
>>  
>>john upham 
>>  
>>John Upham, Dept. of Chemistry, University of Reading, Berks., RG6 2AD, UK.
>>Email: [email-masked] (BITnet), [email-masked] (Janet)
>>Voice:   +44 734 875123 x7441 (day), Fax: +44 734 311610

Mike Frisch replies:
>If all calculations were done to infinite precision and all geometry
>optimizations were continued until the forces were exactly zero, and all
>optimizations used analytic first AND SECOND derivatives then optimizations
>starting from the same structure but using different coordinate systems
>would go to exactly the same place.
    [text deleted]
>In fact, because of the differences in the Hessians, two optimizations with
>different coordinates started at the same point far from any minimum might
>fall down into different wells and to different minima.  (This will not
>be the case if the two optimizations use analytic rather than approximated
>second derivatives at every point -- then every step should be the same
>regardless of coordinate system.)
>
>All of these considerations apply to comparing cartesian with internal
>coordinates and also to comparing two different sets of internal coordinates.
>
>Mike Frisch

I agree with Mike Frisch's comments regarding multiple minima.
Specifically, different optimization algorithms can lead to different
minima on a surface.  This includes different hessian update
procedures in quasi-Newton approaches, different approximations to
derivatives (analytic vs. forward difference vs. central difference
vs. ...etc).  However, I believe that this also applies to different
definitions of a coordinate system, even for analytic (i.e. exact)
derivatives.

Consider minimixation of the 1-dimensional function f(x)=half*x**2.
Newton-Raphson with analytic derivatives gives the update procedure

    Delta_x = -x0

where x0 is the initial guess.  This gives 

    xnew = x0 + Delta_x = 0

So that the minimum at x=0 is always reached in one iteration
regardless of the initial guess.

Now consider the "coordinate transformation"

    y = exp(x)

This is a one-to-one onto mapping, so there are no problems with
multiply-defined angles or such.  Now, define an F(y) such that
F(y)=f(x).  That is, F(y) is the original function expressed in terms
of the new coordinate, y.  The answser using the above transformation
is

    F(y) = half * (ln(y))**2

The original function was minimized at x=0, whereas the new function
is minized at y=1.  That is, f(0)=F(1)=0.

Again using analytic derivatives, the NR method in this new coordinate
system results in the update procedure

    Delta_y = ln(y) / ( ln(y) - 1 )
            = -x0 / ( 1 - x0 )

The last line gives the y update in terms of the original x coordinate
for comparison.  Clearly, this does not converge in 1 iteration,
although it does converge quadratically in the neighborhood of the
minimum as all good NR methods do.  

Note however that in the neighborhood of x=1 (or y=e) there will be
problems because of the denominator in the NR update expression in the
"y" coordinate system.  Specifically, the NR method diverges near x=1
using the procedure in the y "coordinate", whereas it is globally
convergent, and converges in 1 iteration, when everything is done in
the "x" coordinate system.

If all of this is true for a simple 1-dimensional function, it is also
true for more complicated (3N-6)-dimensional energy surfaces.  In
summary, convergence rates and convergence behavior depend on the
coordinate system.

-Ron Shepard
[email-masked]

---

From: [email-masked] (Ron Shepard)
Date: Thu, 14 Jan 93 17:30:04 CST

In the previous posting I wrote :-(
>    Delta_y = ln(y) / ( ln(y) - 1 )

On closer inspection, I find that this should have been :-)

    Delta_y = y*ln(y) / ( ln(y) - 1)
            = -x0*exp(x0) / ( 1 - x0 )

-Ron Shepard
[email-masked]