Curve Fitting to a Continuous Function A Useful Tool in Theoretical Chemistry M h i o J. C. Varandas Universidade de Coimbra. 3049 Coimbra Codex, Portugal Chemists are generally familiar with discrete leastsquares analysis as shown by the numerous articles in this Journal that discuss the applicability and limitations of discrete linear and nonlinear least-squares fitting techniques (1-5). However, another variant of the least-squares method in which little attention has comparatively been focused is that of curve fitting to a continuous function (5); to our knowledge, no article of this Journal has been devoted to this topic. Yet, there are important areas in the field of theoretical chemistry where it is often desirable to find a simpler representation of a continuous curve or a function that leads to simplifications for a given application in mind. Since textbooks that deal adequately with this subject are often too formidable for quick reference, we plan to discussit in the present paper. In the next section we make a brief presentation of the problem while sample applications of relevance in the field of theoretical chemistry are treated later.
where N is a constant that may, without loss of generality, be chosen as unity. By differentiating eq 5 with respect to the coefficients, one gets
and, by substitution of this result into eq 4, one obtains
This result shows that, under the normalization condition of eq 5, the least-squares method is equivalent to finding the maximum of the weighted overlap integral between F a n d f, which is defined by
Least-Squares Approximation to a Contlnuous Function Let F(x)tC[a,b] be the function to least-squares fit by a where c is a collective variable reprefunction f(x;c)r~?[a,b], senting the least-squares parameters co, CI, . . . ., c,, and i?[a,b] denotes the set of all functions continuous on the closed interval [a,b]. In the spirit of the least-squares method the problem thus reduces to minimize the error function E(c) with respect to t h e n coefficients,
The fact that this extremum of S(c) is a maximum can easily be seen from eq 2 by expanding E(c) as a sum of three terms and noting that two of these terms (those involving the integrals of p and F)are nonnegative and constant; thus, E(c) is minimum when the remaining term of the sum, -2S(c), is most negative, i.e., S(c) is maximum. Fit to a Polynomial Function
Let f(x;c) he a polynomial function, f(x:c) = Z k x k
More generally, one may wish to consider different weights according to the value of x. For this purpose, one introduces a weight function on [a,b] such that w(x) 2 0 for xr[a,b] but w(x) # 0 on any subinterval of [a,b]. Accordingly, eq 1 assumes the form
The least-squares problem therefore reduces to the solution of a system of n equations obtained from
(9)
k=o
Substitution of eq 9 into eq 4 thus leads to the equations
.
which are known as the normal equations. These equations have always a unique solution provided f(x;c)rC[a,b] and a + b. We consider next the explicit formulas for the simplest two-term polynomial representation.
+
cfin (n = 1, 2, . . .) In this case, with n being a known constant, the normal equations are The Two-Parameter Relation f = co
By substitution of eq 2 into eq 3 and applying well-known theory of calculus, one gets
from which one gets In the remaining of this paper we omit the arguments for simplicity. Suppose now that f(x;c) obeys the following normalization condition with respect to the weight function 28
Journal of Chemical Education
Fmlng of an Unnonnalized Gaussian Function to a Slater Function' C"
Fit
=
s
b
C.
S
where
Fit of a Slater Function to a Gaussian Function Molecular electronic structure calculations are usually carried out in t e r n s of analvtic functions throueh which the molecular orbitals are approximated. Two sets i f hasis functions have found common use: the Slater-type functions (SF'S) and the Gaussian-type functions (GF's). If normalized, theSFand GF for thesimplest lsorhital have the form
where r a n d a are the Slater and the Gaussian orhital exponents. resnectivelv: eeneralizations of the S F and GF for orbit& with high& & ~ ~ u lmomenta ar are obtained by introducing polynomials in x , y, and z (components of r) multiplying the same exponential or Gaussian fall-off. The diffuseness or size of the basis functions is thus determined by the magnitude of the orbital exponents, which are positive numbers larger than zero. The two major differences between the SF'S and GF's refer to their behavior a t r = 0 and a t large values of r. At r = 0, the SF'S have a nonzero slope in contrast to the GF's, which have a zero slope. At large r, the GF's have a much more rapid fall-off. As a result, the SF'S give a more realistic representation of the actual wave functions, e.g., a t large which is distances the molecular orbitah decay as exp (-k), typical of the SF'S rather than the GF's. In particular, the exact solution of the 1s orbital for the hydrogen atom is the SF eiven hv ea 15 with t = 1.0. &w to exdlain then the interest of GF's in quantum chemistrv? The answer stems from the fact that the most time-cokuming step in a self-consistent-field (SCF) electronic structure calculation comes from the evaluation of about N4/8 ( N representing the number of electrons) hielectronic integrals of the type
.
where XA. XB, XC, and XD are one-electron orbitals centered on A, B. C, and D, resnectivelv, and the intemation is carried over the complete space of theelectronic coordinates. These integrals are relatively easy to solve in closed form using GF's as distinct from the numerical solution required in the case of SF'S. The reason is that the product of two GF's, each on a different center is, apart from a constant, a GF on a third center. As a result.the four-center molecular integrals reduce to two-center-type ones (6). One thus needs to represent the SF'S in terms of GF's. This may he accomplished through the least-squares method. Consider for example the case of a 1s orbital with f = 1.0. Since the radial part, eq 16, of the 1s Gaussian orhital is required to he normalized (it obeys eq 5 with x = r, w = 4&, and N = 1). the least-squares optimization is equivalent to finding (by numerical methods) the maximum ofthe overlap integral, S(a), between the S F and the GF. It may be shown that S is maximum for' a = 0.2710. The fitting of a GF to a S F can also be used to illustrate a linearized solutiun of this problem that consists of adjusting
In @GF by imposing an appropriate transformation-related weighting (4). The problem is thus to find co and cz in f = In dcF = c, c d with F = In dSF. As shown in the next paragraph, the equivalent linear minimization requires as weight function W = 4ar2(dSF)2. L& he the variable bef&e'transformation and, e.g., Y = In y the transformed one for which the problem becomes linearized. Assuming that Ay and AY are small enough over the entire range of the integration coordinate, one may write AY-dY-dhy-1 (18) 4Y dy dy Y from which it follows that
+
i.e., the weight function for the transformed equation is W = wy2 or, equivalently,
Thus, a least-squares fit in which S b ( ~ y / y ) ~ disx minimized corresponds, within the ahove approximations, to one is minimized. The table summarizes in which J',~W(AY)%X the results of a few linearized least-squares fits. Note that, e.g., for fit V the unnormalized GF can be normalized by multiplication hy the factor (2 X 0.2222/~)3'~/0.2452= 0.9408. Figure 1compares the unnormalized and normalized GF's of fit V with that obtained from the one-parameter (the restriction imnosed bv ea 5 determines the remaining varameter) n o n h e a r le&-squares procedure. Also incLded for comnarison in this fieure is bSF.Note that the intemation limits a-and b can he varied toget an improved representation of the S F over a range of values of r crucial for some specific quantum mechanical application, say the tail of the wave function for calculations of long-range interaction enereies. Of course, a more accurate representation of the S F can be obtained by a sum of GF's. Two situations may again he distinguished. First, the exponents are selected according to a suitable criterion (e.g., as multiples of a prechosen value) and the expansion coefficients then treated as the linear least-squares variables: in this case the ahove equations can still be used by carrying out the transformation exp (-a;r2) x. Thefittedfunction is then normalized aposteriorias in the previous example. Alternatively, one performs a nonlinear least-squares optimization of the GF exponents and expansion coefficients by finding numerically the maximum of the corresponding overlap integral, eq 8. In this case, a more reliable n-Gaussian description (usually denoted STO-nG) +
-
Atomic units are used throughout this paper: a. = au of length = 0.0529177 nm; 4, = au of energy = 4.3597482a J. This property may be of practical value for quasiclassical trajectory ( 11)studies of rare gas atom-rare gas dimer reactions (12)as it offers an easier route (9) for the separation of the vibrational and rotational energies in the products. Volume 67 Number 1 January 1990
29
(r= ---_.
-.
Figure 1. Comparison of he 1s
1)Sieter function wim some Gaussianfits.
Key for symbols: SF: normalized (a= 0.2710) GF obtained by maximization of he overlap inlegal: ..... unnwmalized GF (n = 0.2222) obtained by solvlng he normal equations, eqs 12-14; -.-..OF (n = 0.2222) after nwmaliratlon.
of the SF is obtained for the same number of expansion terms. Examples involving three and more GF's can he found in the especialized literature (6)to which the reader is addressed for further details. Fn ol a Reallstlc Potentlal Curve by a W(10.6) Fundlon In molecular dynamics simulations (7)one is often faced with the dilemma of having to choose between a simple potential function (8),which is not too realistic, and a more realistic model that is unavoidably more time consuming. Perhaps, the most popular simple form for nonhonding interactions is the LennardJones (rn,n) function defined by V'-%,n;$Z)
= CmR-"
- C,,R-"
(21)
where C,, and C, are positive constants. We chose rn = 10 and n = 6 for the fitting of a realistic Krz potential energy curve from the literature (9) (denoted EHFACE) bv ea 21. Figure 2 illustrates the g o n e gets by minimiz
