Introduction This new article describes the exponential curve fitting method implemented in Graphics-Explorer,my equations grapher program. It replaces the old article, which can be found [here]. New is an exerciser program allowing step by step observation of the curve fitting process.The curve fitter calculates the best fitting exponential function given a set of points. This function is y = a.b
^{x} + cparameters.Contents Exerciser description Below is a reduced picture of the exerciser at work
Operation Points can be added by:-
1. selecting 1 of 4 preset point sets (click radio button)
2. mouseclick on the coordinate system (2nd click removes point) 3. typing into the [x,y] points list To watch the process step by step, mark the [testmode] checkbox. Then click the [GO] button to calculate the a,b,c parameters and paint the graph. To clear parameters or points, click the appropriate [X] button or [RESTART]. Curve fitting theory Input to the curve fitter is a set of points [x_{1},y_{1}]....[x_{n},y_{n}]The minimal required number of points is 3. The maximum number of points is 10. The fitter calculates parameters a,b,c such that the curve y = a.b ^{x} + c has the smallest distance to these points.For internal computations y = a.e ^{bx} + c is used (e^{b} is reported for b)The method has a numerical and an analytical part. These are the steps: -
1. obtain an approximation of b
2. using b, obtain an approximation of a 3. using a,b, obtain an approximation of c 4. discard a,b and replace y = ae ^{bx} + c ....by.... y - c = ae^{bx}5. ln(y - c) = ln(a) + bx ...........a linear function 6. call the least square method to find a,b of the best fitting lineapproximation of b array dx holds the x differences between adjacent pointsdx1 = x2 - x1....dx3 = x4 - x3 array dy holds the y differences dy1 = y2 - y1....dy3 = y4 - y3 array cx holds the centers between points cx1 = (x1+x2)/2....cx3 = (x3+x4)/2 array dq holds the differential quotients dq1 = dy1/dx1....dq3 = dy3/dx3 Note: in the Delphi language, cx1 is written as cx[1]...etc. Starting with -
y = ae
^{bx} + c..........differentiationy' = abe ^{bx}So, we may write -
y
_{1}' = abe^{bx1}y _{2}' = abe^{bx2}y _{3}' = abe^{bx3}with above arrays in mind: The center values (cx) are used, because the tangent (differential quotient dq) is most accurate between two x values. For subsequent points, the average value of b is calculated. This is the approximation of b. approximation of a Starting with-
y
_{1} = ae^{bx1} + cy _{2} = ae^{bx2} + c...... -
dy1 = y
_{2} - y_{1} = a(e^{bx2} - e^{bx1})An average value of a is calculated. This is the approximation of a. approximation of c -
y = ae
^{bx} + c .....soc = y - ae ^{bx}-
c = y
_{1} - ae^{bx1}c = y _{2} - ae^{bx2}.... This is the approximation of c. Intermezzo In equation-
y = g
^{x}
g is called the growth factor.Exponential functions (having c=0, the x-axis is the horizontal asymptote), are the result of constant relative growth. In equation -
y = e
^{bx}
^{b}.b is called the relative growth speed because So, if a function has a constant relative growth speed of b, it's growth factor g = e ^{b}where e is the base of the natural logarithm. (e = 2.71828....) Before, we have obtained an approximation of c so we may write ...y - c = ae ^{bx}If the points (x1,y1)...(xn,yn) are on a true exponential function, we will find the same growth factor, comparing subsequent points. In the next part, fine tuning of c, we are going to vary c to obtain the smallest deviation in growth factors. fine tuning c At this point we forget the a and b values and continue with c and the points [x1,y1]...[xn,yn].For a certain value of c we write -
y
_{1}-c = ae^{bx1}y _{2}-c = ae^{bx2}..... Two successive points produce a value of ee. The values of ee are summed in variable AVG, for now AVG = ee1 + ee2 + .... The growthfactors are stored in array growthfactor: growthfactor1 = e ^{ee/dx1}..........By summing the ee values, we actually sum values b.dx1 + b.dx2 + b.dx3... = b(last X - first X) So, b = AVG / (last X - first X) Average growth factor is e ^{b}-
AVG = e
^{AVG/(lastX-firstX)}
Deviation = abs(AVG-growthFactor1) + abs(AVG-growthfactor2) + ...... Deviation is the result of procedure deviation which is called for a certain value of c.least square method The least-square-method is an algorithm to calculate the best fitting polynomial y = c_{0} + c_{1}x + c_{2}x^{2}+ ....given a set of points. With a final value of c, we may write which represents a linear function like y = a + bx. We present the points (ln(y _{i}-c),x_{i}) to the least square procedures to obtain the best a,b values.Resulting a is replaced by e ^{a} to cancel the ln( ) function.The least square method is a nice application of linear algebra. Look [here] for a description. Aborts Thedeviation procedure produces an abort if c is between two y values.In that case y = c intersects the curve. c must be outside the range y1, y2,y3..... -
ee = (y2 - c)/(y1 - c)
abort = ee <= 0 Below is the flowchart for the fine tuning of c The fine tuning ends when deltaC has become very small. Watchdog timer Proceduremorepointsexp controls the calculation of the a,b,c parameters.Floating point operation are between try...except directives.There are 2 cases for an abort: 1.Illegal Floating point operation This situation occurs if the points are too far off to make an exponential function an example being (1,1) , (5,5) , (8,2) 2.C runaway If C is too far out, an even further position may yield a smaller growth factor deviation. So, C "runs away". The watchdog counter w ends this process. Fine tuning of c must be completed in maximal 100 steps, otherwise the process is aborted. This situation occurs in case of sharp bends. The Delphi-7 project There are 4 units and 1 form:-
form1,unit1: paintbox, events, controls
unit2: data and procedures for the curve fitting algorithm unit3: procedures for testmode reporting unit4: least square procedures Unit1 Graphs are painted in bitmapmap, points are painted on paintbox1.The dimension is 800*800 pixels, coordinate (0,0) is at map position (400,400). x domain is -10< x < 10,.........y range is -10< y <10. The scale is 40:1. There are 2 important flags: -
pointsOK : all points in the list are OK and sorted in ascending order.
paramsOK : parameters are valid If paramsOK = false, pressing [GO] calculates a,b,c and then plots the function. if parameters a,b,c do not represent the actual situation, they are painted in red. Please refer to the source code for details. Unit2 Main call is made to procedure makeExpFunction.Points are collected from unit1 and stored in array Epoints[ ] of type Tpnt. type TPnt = record x,y : double; end; procedure deviation calculates the variation of growthfactors while fine tuning c. These values are stored in variables Bvar1,Bvar2. Unit3 Unit2 calls procecures in unit3 to display the steps in the curve fitting process.Information is written in form1 paintbox1. The waitflag : boolean (in unit1) controls the step-by-step action. While the waitflag is true, the program executes the processmessages procedure. Unit4 Has the necessary procedures for theleast square method.Main call is to procedure buildexp12(var pp : Tpoints; maxp : byte), where Tpoints is the array of points and maxp the number of points supplied. This concludes the description. |
|||||||||||||||||||||||||||||||||