CH-02 (523169), страница 2
Текст из файла (страница 2)
When the third argument is set equalto 4, we then have the case of exact curve-fit of five points by a fourth-orderpolynomial. Readers are referred to the program ExactFit for more discussions.Also, it is of interest to know whether one may select an arbitrary set offunctions and linearly combine them for least-squares fit, instead of the unbrokenset of polynomial terms X0, X1, X2, …, XN.
Program LeastSqG to be presentedin the next section will discuss such generalized least-squares curve-fit. But beforewe do that, let us first look into a situation where program LeastLQ1 can beapplied for a given set of data after some mathematical transformations areemployed to modify the data.Transformed Least-Squares Curve-FitThere are occasions when we know in advance that a given set of data supposedto fall on a curve described by exponential equations of the type:© 2001 by CRC Press LLCy = b1e b2x(5)y = c1xe c2x(6)orTo determine the coefficients b1 and b2, or, c1 and c2 based on the least-squarescriterion, Equation 5 or 6 need to be first transformed into a linear form. To do so,let us first consider Equation 5 and take natural logarithm of both sides.
It gives: n y = n b1 + n e b2x = n b1 + b 2 x(7)If we introduce new variable z, and new coefficients a1 and a2 such that:z = n y,a1 = n b1,anda 2 = b2(8,9,10)Then Equation 7 becomes:z = a1 + a 2 x(11)Hence, if we need to determine the coefficients b1 and b2 for Equation 2.8 basedon N pairs of (xi,yi), for i = 1 to N, values and on the least-squares criterion, wesimply generate N z values according to Equation 2.11 and then use the N (xi,zi)values as input for the program LeastSq1 and expect the program to calculate a1 anda2. Equations 9 and 10 suggest that b2 is to have the value of a2 while b1 should beequal to e raised to the a1 power, or, EXP(a1) with EXP being the exponential functionavailable in the FORTRAN or QuickBASIC library.Equation 6 can be treated in a similar manner by taking logarithms of both sidesto obtain: n y = n c1 + n x + n e c2x = n c1 + n x + c2 xor n y− n x = ny= n c1 + c2 xx(12)If we introduce new variable z, and new coefficients a1 and a2 such that:z = n(y x),a1 = n c1,anda 2 = c2(13,14,15)Then, Equation 12 becomes Equation 11 and a1 and a2 can be obtained by theprogram LeastSq1 using the data set of (xi,yi/xi) values.As a numerical example, consider the case of a set of nine stress-versus-strain( vs.
) data collected from a stretching test of a long bar: (.265,1025), (.4,1400),(.5,1710), (.7,2080), (.95,2425), (1.36,2760), (2.08,3005), (2.45,2850), and© 2001 by CRC Press LLC(2.94,2675) where the units for the strains and stresses are in microinch x 102 andlb/in2, respectively. When program LeastSq1 is applied for the modified data of(,), the coefficients for the linear fit are a1 = 15.288 and a2 = –537.71. Consequently, according to Equations 13 and 14, and realizing that x and y are now and, respectively, we have arrived at = 4.3615 × 106e–537.71. The derived curve andthe given points are plotted in Figure 2 which shows the curve passes the origin asit should.FIGURE 2.
The curve passes the origin as it should.© 2001 by CRC Press LLC2.4 PROGRAM LEASTSQG — GENERALIZEDLEAST-SQUARES CURVEFITProgram LeastSqG is designed for curve-fitting of a set of given data by a linearcombination of a set of selected functions based on the least-squares criterion.2Let us consider N points whose coordinates are (Xk,Yk) for k = 1 to N and letthe M selected function be f1(X) to fM(X) and the equation determined by the leastsquares curve-fit be:MY( X) = a1f1 ( X) + a 2 f2 ( X) + … + a M fM ( X) =∑ a f (X)j j(1)j=1The least-squares curve-fit for finding the coefficients c1–M is to minimize theerrors between the computed Y values based on Equation 1 and the given Y valuesat all Xk’s for k = 1 to N. Let us denote yk to be the value of Y calculated at Xkusing Equation 1 and S be the sum of the errors denoted as ek for k = 1 to N.
Sincethe yk could either be greater or less than Yk, these errors ek’s may cancel from eachother if they are added together. We therefore consider the sum of the squares ofthese errors and write:NS = e12 + e 22 + … + e 2N =∑e2k(2)k =1where for k = 1,2,…,Ne k = y k − Yk = M∑j=1a jfj ( X k ) − Yk(3)It is obvious that since Xk and Yk are constants, the sum S of the errors is afunction of a1 to M.
To find a particular set of values of a1 to M such that S reaches aminimum, we may therefore base on calculus3 and utilize the conditions ∂S/∂ai = 0for i = 1 to M. From Equation 2, we can have the partial derivatives of S with respectto ai’s, for i = 1 to M, as:N ∂e∂e ∂e∂e∂S= 2 e1 1 + e 2 2 + … + e N N = 2ek k∂a i∂a i∂a i∂a i ∂a1k =1∑From Equation 3, we note that ∂ek/∂ai = fi(Xk). Consequently, the M extremumconditions, ∂S/∂ai = 0 for i = 1 to M, lead to M algebraic equations for solving thecoefficients a1 to M. That is, for i = 1 to M:M Nfi ( X k )fj ( X k )a j = k =1∑∑j=1© 2001 by CRC Press LLCN∑ f (X )Yik =1kk(4)If we express Equation 4 in matrix notation, it has the simple form:[C]{A} = {R}(5)where [C] is a MxM square coefficient matrix, and {A} and {R} are Mx1 columnmatrices (vectors).
{A} contains the unknown coefficients a1 to M needed in Equation 1.If we denote the elements in [C] and {R} as cij and ri, respectively, Equation 5indicates that these elements are to be calculated using the formulas, for i = 1,2,…,M:Ncij =∑ f (X )f (X )ikjk(6)k =1andNri =∑ f (X )Yikk(7)k =1The above derivation appears to be too mathematical; a few examples of actualcurve-fit will clarify the procedure involved. As a first example, consider the caseof selecting two (M = 2) functions f1(X) = 1 and f2(X) = X to fit three given points(N = 3), (X1,Y1) = (1,1), (X2,Y2) = (2.6,2), and (X3,Y3) = (2.8,2). Equations 6 and 7then provide the following:andHence, the system of two linear algebraic equations for finding a1 and a2 for thestraight-line equation is:36.4© 2001 by CRC Press LLC6.4 a1 5 =15.6 a 2 11.8The solution can be obtained by application of Cramer’s Rule to be a1 = 0.42466and a2 = 0.58219.
More examples will be given after we discuss how computerprograms can be written to compute [C] and {R} and then solve for {A}.Program LeastSqG provided in both FORTRAN and QuickBASIC versionsis developed using the above two equations. It can be readily applied for calculatingthe coefficients c1 to M. Both QucikBASIC and FORTRAN versions are listed andsample applications are presented below.QUICKBASIC VERSION© 2001 by CRC Press LLCSample ApplicationsWhen four functions are selected as those listed in SUB FS, an interactiveapplication of the program LeastSqG QuickBASIC version using the input dataentered through keyboard has resulted in a screen display of:If three sinusoidal functions, sin(x/20), sin(3x/20), and sin(5x/20) wereselected and replacing those listed in SUB FS, another interactive application of theprogram LeastSqG QuickBASIC version is shown below.© 2001 by CRC Press LLCFORTRAN VERSION© 2001 by CRC Press LLCSample ApplicationBy selecting the four functions listed in Subroutine FS, an interactive applicationof the program LeastSqG using the input data given below has resulted in a screendisplay of:MATLAB APPLICATIONA LeastSqG.m file can be created and added to MATLAB m files which willtake N sets of X and Y points and fitted by a linear combination of M selectedfunctions in the least-squares sense.
The selected functions can be specified inanother m file called FS.m (using the same name as in the FORTRAN and QuickBASIC versions). First, let us look at a version of LeastSqG.m:Notice that the coefficients C’s is obtained by solving [A]{C} = {R} as in thetext. For MATLAB, a simple matrix multiplication of the inverse of [A] to and onthe left of the vector {R}. Complete execution of LeastSqG.m will be indicated bya display of the matrix [A], vector {R}, and the solution vector {C}. The expressionfeval(funfcn,i,X(k)) in the above program is to evaluate the ith function at X(k)defined in a function file to be specified when LeastSqG.m is applied which is tobe illustrated next.© 2001 by CRC Press LLCConsider the case of given 5 (X,Y) points (N = 5) which are (1.4,2.25), (3.2,15),(4.8,26.25), (8,33), and (10,35).
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
