Subject: constrained opt's

From: Matt Stahl <[email-masked]>
Date: Mon, 31 Jul 1995 09:06:18 -0700 (MST)

Netters,

	A while back I posted a question about the validity of constrained
optimization.  The original post and the summary of responses are included
in this message.  In short, I think constrained opt's are generally
accepted as valid in the context of rotation driving studies.  There are
"better" methods available for reaction path characterization, but as with
basis set/level of theory choices it really depends on how much time and
effort one wants to spend at that step of the research process.
	I would like to thank all of the respondants very much.  Yet again
I have been impressed with the depth of knowledge available in the list.

Sincerely,
matt

__________________________ original post _________________________________

	I am working on a project which involves parameterizing molecular
mechanics torsion potentials.  We used gaussian to generate energy vs. 
torsion data by performing constrained optimizations.  One torsion angle was
held constant while the rest of the coordinates in the molecule were
optimized.  Someone recently made the comment that constrained optimizations
with gaussian are not appropriate for such a study because the 
unconstrained coordinates do not fully optimize.  This suggestion seems 
odd because the literature seems to be repleat with examples of isolating 
the energetics of a molecular motions (distances, angles, and torsions) 
through constrained optimization.  How far off the mark am i in thinking 
that constrained opts. are appropriate for this type of work. 
	The alternative method would be to perform an IRC (internal reaction
coordinate) calculation to characterize the path between the two minima of
interest, but as i understand it IRC calculations are dependant on Linear
Synchronous Transit calculations which ESTIMATE the transition state.  Can
one perform an IRC from a calculated transition state (opt=ts)?  If not, how
much more "valid" is an IRC calculation than constrained optimization if at
any point the IRC has to "guess" the path? 

___________________________ responses ___________________________________

>From [email-masked]
Date: Fri, 14 Jul 1995 09:43:26 +0200
From: Per-Ola Norrby <[email-masked]>
To: Matt Stahl <[email-masked]>
Subject: Re: CCL:validity of constrained calcs.

        Dear Matt,

        If you are interested in saving CPU-time while doing the torsional
parameterization, we have recently published a method that will allow you
to do this in a simple and non-ambiguous manner.  In short, if you do it
correctly, it's actually best to do single point ab initio calculations on
several conformations generated by rigid rotations.  The trick is doing the
MM calculations also as single points on exactly the same confromations,
and adjusting the torsion parameters to fit the energy differences at each
point.  Not only do you save CPU-time, but you decouple possible errors in
the remaining force field (except non-bonded interactions) from the current
parameterization.  The reference is: Norrby et al., J. Comput. Chem., 1995,
16(5), 620.

        Sincerely,

        Per-Ola Norrby

 *  Per-Ola Norrby
 *  The Royal Danish School of Pharmacy, Dept. of Med. Chem.
 *  Universitetsparken 2, DK 2100 Copenhagen, Denmark
 *  tel. +45-35376777-506, +45-35370850    fax +45-35372209
 *  Internet: [email-masked], [email-masked]

_________________________________________________________________________

>From [email-masked]
Date: Fri, 14 Jul 1995 14:55:10 MDT
From: Stefan Fau <[email-masked]>
To: [email-masked]
Subject: Re: CCL:validity of constrained calcs.

Hi Matt,
I didn't hear of problems with unconstrained variables 
that are not fully optimized. G9x determines the 
convergence of the optimization by the forces along 
your geometry variables and the displacement(i.e. the 
change in the value of the variable). You can look at 
two sections in the G9x output: The displacement for 
each variable (frozen or not) is listed directly above 
the line "0    Item   Value   Converged?" and the 
corresponding forces are listed under "Internal 
Coordinate Forces". The units are the same in those two 
places and below the "0   Item   Value   Converged?" 
line. Therefore your partial optimization has converged 
if all _free_ variables have forces and displacements 
below the threshhold values.

The IRC calcs start with second derivatives (calculated 
at the beginning of the IRC calc or read in). They are 
not connected to LST calcs unless you use LST to find 
your TS. You provide the geometry of the transition 
state and you can find the TS however you like. It is 
important to choose a Z-Matrix with enough degrees of 
freedom to allow symmetry lowering on the way to the 
minima. The IRC calcs find a steepest descent path from 
the TS to the minima with the additional restriction 
that the path is built up from arcs of circles that are 
smoothly connected (for more details see: C.Gonzalez and 
H.B.Schlegel, J. Chem. Phys. 90,2154 (1989); J. Chem. 
Phys. 94,5523 (1990) and J. Chem. Phys. 95,5853 (1991)).

Hope it helps, 
                  Stefan

Stefan Fau,               [email-masked]

FB Chemie der Philipps-Universitaet Marburg,
Hans-Meerwein-Str.
D-35032 Marburg

[email-masked]
_________________________________________________________________________

>From [email-masked]
Date: Fri, 14 Jul 95 09:46:51 -0400
From: [email-masked]
To: [email-masked]
Subject: CCL:validity of constrained calcs.

Matt,

  We have have been doing parameterization for a number of years, much
of which was based on the reproduction of ab initio calculated
dihedral surfaces.  These were obtained from the various version of
gaussian by constraining one degree of freedom and minimizing the
remainder.  I am not aware of any problems with the optimization of
the remaining degrees of freedom.  If you do think that there are
problems you can perform a freqency calculation to determine the
number of negative frequencies; ideally you should only obtain one for
the constrained degree of freedom.  

good luck

Alex MacKerell, [email-masked]
School of Pharmacy
University of Maryland at Baltimore
20 North Pine Street
Baltimore, MD  21201
410-706-7442

_________________________________________________________________________

>From [email-masked]
Date: Mon, 17 Jul 95 10:36:22 -0400
From: Joel Polowin <[email-masked]>
To: Matt Stahl <[email-masked]>
Cc: [email-masked]
Subject: Re: CCL:validity of constrained calcs.

> Date: Thu, 13 Jul 1995 13:24:12 -0700 (MST)
> From: Matt Stahl <[email-masked]>
> Subject: CCL:validity of constrained calcs.
> 
> 	I am working on a project which involves parameterizing molecular
> mechanics torsion potentials.  We used gaussian to generate energy vs. 
> torsion data by performing constrained optimizations.  One torsion angle was
> held constant while the rest of the coordinates in the molecule were
> optimized.  Someone recently made the comment that constrained optimizations
> with gaussian are not appropriate for such a study because the 
> unconstrained coordinates do not fully optimize.  This suggestion seems 
> odd because the literature seems to be repleat with examples of isolating 
> the energetics of a molecular motions (distances, angles, and torsions) 
> through constrained optimization.  How far off the mark am i in thinking 
> that constrained opts. are appropriate for this type of work. 

If you add a constraint -- say, a torsion angle -- then none of those atoms
have their normal freedom to move; they and their attachments won't be
fully optimized.  Through non-bonding interactions and propagated bonding
interactions, everything else will be slightly distorted as well.  If this
is what you're talking about, it's a problem inherent in the procedure.  In 
my opinion, it is reasonable to use these results as long as one keeps that 
problem in mind; in general I think that the structural distortions are
not serious.

In principle, I suppose it would be possible to get the same results in
other ways.  Most of my work has involved incrementing a constrained
torsion angle to get the barrier to rotation of part of a system with
respect to the rest.  Instead of doing that, I imagine that one could
set up a molecular-dynamics scheme where part of the system was stationary
and the rest was started with initial velocities applied to the atoms
corresponding to the desired rotation.  For sufficiently high velocities,
the barrier would be overcome.  It would be a matter of repeating the
calculations several times with different starting velocities to find out
how much "rotational energy" would be necessary to overcome the barrier.
I'm not sure how effective such calculations would be compared to the
usual procedures.

Joel
[email-masked]

_________________________________________________________________________

>From [email-masked]
Date: Mon, 17 Jul 1995 16:43:33 -0400
From: "E. Lewars" <[email-masked]>
To: [email-masked]
Subject: CONSTRAINTS

I think the use of constraints for your purpose is valid.  The optimization
will give the other param's optimized *for the particular value of the 
constraint*, which is what you want, I think.  In other words, you will
get a slice through the potential energy hypersurface corresponding to
the constrained value of one torsional.  Of course the other param's won't be fully optimized, but you do not want them to be---you want the values they 
assume when the molecue is held at some particular dihedral angle (torsional)
9or you want the energy of this contrained, partially optimized geometry).
Errol Lewars
===

>From [email-masked]
Date: Tue, 18 Jul 1995 10:01:16 +0200 (MST)
From: Thorsten Koch <[email-masked]>
To: [email-masked]
Subject: Re: CCL:validity of constrained calcs.

> > Date: Thu, 13 Jul 1995 13:24:12 -0700 (MST)
> > From: Matt Stahl <[email-masked]>
> > Subject: CCL:validity of constrained calcs.
> > 
> >  	The alternative method would be to perform an IRC (internal reaction
> > coordinate) calculation to characterize the path between the two minima of
> > interest, but as i understand it IRC calculations are dependant on Linear
> > Synchronous Transit calculations which ESTIMATE the transition state.  Can
> > one perform an IRC from a calculated transition state (opt=ts)?  If not, how
> > much more "valid" is an IRC calculation than constrained optimization if at
> > any point the IRC has to "guess" the path? 
> > 
As far as I can see, the LST takes the maximum of the linear path between 
reactants and products as a guess for a transition state. This has to be 
optimized using opt=ts. If you get a proper transition state from that 
you can use it for a IRC calculation.
Note that in G94 there is a better method to find transition structures. 
Its invoked by opt=qst2 and input requires the route card, title and 
molecule specifications of the product and title and molecule 
specifications of the reactand. Of course, corresponding atoms must be in 
the same order in the z-matrix.

Thorsten

/-----------------------------------------------------------------\
|                        Thorsten Koch                            |
| Institut fuer physikalische Chemie II der Universitaet zu Koeln |
|                    [email-masked]                       |
|                   Tel. +49 [0]221 470 4816                      |
\-----------------------------------------------------------------/

_________________________________________________________________________

>From [email-masked]
Date: 19 Jul 1995 14:07:58 -0800
From: Richard Bone <[email-masked]>
To: Matt Stahl <[email-masked]>
Subject: Re: CCL-validity of constrai

        Reply to:   RE>CCL:validity of constrained calcs

Matt

A couple of comments for you:  you raise a pretty thorny issue of what
constitutes the correct 'reaction path'.  

The general problem with what you have done is that you have assumed
that the coordinate system is the same throughout the process.  This 
will not necessarily be the case, i.e., with full relaxation, the proper
modes of the system maybe such that the internal coordinate you are
interested in does not correspond with anything particularly physical
at all points on the path.  On the other hand, what you are doing does
serve to illustrate particular points effectively because it is easily
understandable.   Often, at intermediate points on a reaction path, the
system is in low symmetry.  This allows any motion to be expressed as
an l.c. of lots of low-symmetry motions.  Put another way, the motion is
not 'pure', it is contaminated by other motions of the same (low) symm-
etry.  This makes analysis fairly tricky.  

If all you are interested in is the start and end-points, then it probably
does not matter too much how you get between them.  An IRC (strictly 
only valid in mass-weighted cartesians) allows you to  roll down hill 
>from a (known) T.S.  to wherever the path takes you.  How you walk
'uphill' (i.e., do you start by displacing along a normal coordinate?) is 
still a problem.  

The whole business of factoring out various degrees of freedom in order
to  facilitate study of a particular motion introduces all sorts of 
artifacts and it depends on how much they bother you and the type of
motion involved as to what you do about it.  

Just because the literature is full of anything of a particular nature,
doesn't mean it's the best/correct way to go !

What you probably have to do is obtain the (internal) coordinate system at
each new point on the path and use the 'mode' which has greatest overlap
with the previous mode that you displaced along to make your next step.
This is of course, time-consuming and tricky, which is part of the reason
why most people don't bother.  

The subject is a minefield;  watch where you tread!

Richard Bone

__________________________________________________________

Richard G. A. Bone, PhD.
Computational Chemist
Terrapin Technologies, Inc.
South San Francisco
USA

E-mail  [email-masked]