Numerical Algorithms for Function Manipulation: Tabulation, Plotting, Derivatives, and Roots
Chapter 3.6 delves into programming numerical algorithms for mathematical functions. It presents generic procedures for tabulating, plotting, and computing zeros, derivatives, and integrals of given functions. It introduces function types and demonstrates how to define functions, such as sine and polynomial expressions, using structured and simple types. Additionally, the chapter covers Newton's algorithm for finding function roots, discussing its iterative method and state management. This comprehensive guide is essential for anyone interested in numerical analysis and programming techniques.
Numerical Algorithms for Function Manipulation: Tabulation, Plotting, Derivatives, and Roots
E N D
Presentation Transcript
Chapter 3.6 Programming of Numerical Algorithms
Numerical Algorithms • Generic procedures to • tabulate, plot, ... a given function • to find the zeros of a given function • to compute the derivative of a given function • to integrate a given function • How is the given function defined ? • By a function procedure • Passed as a parameter to the generic procedures
TYPES • Simple Types • Ordinal Types • REAL Type • POINTER Type • Structured Types • ARRAY • RECORD • SET • PROCEDURE
PROCEDURE Type TYPE RFunction = PROCEDURE(REAL):REAL; PROCEDURE SIN(x : REAL):REAL; BEGIN RETURN VAL(REAL,Sin(VAL(LONGREAL,x))) END SIN; PROCEDURE F1(x:REAL) : REAL; BEGIN RETURN x - VAL(REAL,Exp(VAL(LONGREAL,1/x))) END F1; PROCEDURE F2(x:REAL) : REAL; BEGIN RETURN x*x*x/3.0 - 3.5 * x*x + 10.0 * x - 5.0; END F2;
Tabulating a function PROCEDURE Tab(F : RFunction; First, Last : REAL; NPoints : CARDINAL); VAR x,Step : REAL; BEGIN IF Last > First THEN x := First; Step := (Last - First) / FLOAT(NPoints-1); REPEAT WrReal(x,15,15); WrStr(" : "); WrReal(F(x),15,15); WrLn; x := x + Step UNTIL x > Last; END (* IF *) END Tab;
Plot ProcedureProcedure Header PROCEDURE PlotRF(F : RFunction; First, Last, Min, Max : REAL; NXPoints, NYPoints : CARDINAL); (* F : Y = F(X) is the function to be *) (* plotted by PlotRF *) (* First : First value of X for the plot *) (* Last : Last value of X for the plot *) (* Min,Max : Extreme Y values of the plot *) (* NXPoints : Number of X points on screen *) (* NYPoints : Number of Y points on screen *)
Plot ProcedureCore Code PROCEDURE PlotRF(F : RFunction; First, Last, Min, Max : REAL; NXPoints, NYPoints : CARDINAL); BEGIN XStep := (Last-First)/FLOAT(NXPoints); YStep := (Max-Min)/FLOAT(NYPoints); IF Init(0,0,NXPoints,NYPoints) THEN FOR x := 0 TO NXPoints DO YValue := (F(XStep*FLOAT(x)+First)-Min)/YStep; y := NYPoints - TRUNC(YValue); Plot(x,y,_clrLIGHTRED) END (* FOR *) END; (* IF Init *) END PlotRF;
Plot ProcedureAdding the X axis FOR x := 0 TO NXPoints DO z := TRUNC((0.0 - Min)/YStep); Plot(x,z,_clrGREEN); ... END (* FOR *)
Plot ProcedureTesting for reasonable ranges BEGIN IF Last > First THEN ... FOR x := 0 TO NXPoints DO YValue := (F(XStep*FLOAT(x)+First)-Min)/YStep; IF (YValue >= 0) AND (YValue <= FLOAT(NYPoints) THEN y := NYPoints - TRUNC(YValue); Plot(x,y,_clrLIGHTRED) END (* IF in plot range test *) END (* FOR *) ... END (* IF Last > First *) END PlotRF;
Plot Procedure PROCEDURE PlotRF(F : RFunction; First, Last, Min, Max : REAL; NXPoints, NYPoints : CARDINAL); VAR x,y,z : CARDINAL; YValue : REAL; XStep, YStep : REAL; BEGIN IF Last > First THEN XStep := (Last-First)/FLOAT(NXPoints); YStep := (Max-Min)/FLOAT(NYPoints); IF Init(0,0,NXPoints,NYPoints) THEN FOR x := 0 TO NXPoints DO z := TRUNC((0.0 - Min)/YStep); Plot(x,z,_clrGREEN); YValue := (F(XStep*FLOAT(x)+First)-Min)/YStep; IF (YValue >= 0) AND (YValue <= FLOAT(NYPoints)) THEN y := NYPoints - TRUNC(YValue); Plot(x,y,_clrLIGHTRED) END (* IF in plot range test *) END (* FOR *) END; (* IF Init *) END (* IF Last > First *) END PlotRF;
Newton's Algorithmfor finding a root of a function x0 x x1 |X1-X| < |X0-X|
Procedure NEWTONParameters PROCEDURE Newton (f :RFunction (*Function *); df:RFunction (*First derivative *); GuessedRoot:REAL (*Initial guess *); VAR Root:REAL (*Computed root *); tol:REAL (*Desired precision *); MaxIter:CARDINAL (*Max.Nbr.iterations*); VAR Iter:CARDINAL (*Act.Nbr.iterations*); VAR State:IterState (*Procedure state *));
Procedure NEWTONState Variable IterState = (iterating,rootfound, toomanyiter,toonearzero,tooflat); Iter := 0; State := iterating; REPEAT IF Iter = MaxIter THEN State := toomanyiter ELSIF ... THEN State := toonearzero ELSIF ... THEN State := tooflat ELSE Iter := Iter + 1; ... IF ... < tol THEN State := rootfound END END UNTIL State # iterating;
Procedure NEWTONtoonearzero - tooflat CONST NearZero = 1.0E-20; ... BEGIN ... REPEAT ... ELSIF ABS(Root) < NearZero THEN State := toonearzero ELSIF ABS(dfx) < NearZero THEN State := tooflat ... UNTIL State # iterating; END Newton;
Procedure NEWTONrootfound Iter := 0; Root := GuessedRoot; fx := f(Root); REPEAT dfx := df(Root); ... Iter := Iter + 1; OldRoot := Root; Root := Root - fx/dfx; fx := f(Root); IF ABS((Root-OldRoot)/OldRoot) < tol THEN State := rootfound END UNTIL State # iterating;
Procedure NEWTON PROCEDURE Newton (f,df:RFunction;GuessedRoot:REAL; VAR Root:REAL;tol:REAL; MaxIter:CARDINAL;VAR Iter:CARDINAL;VAR State:IterState); CONST NearZero = 1.0E-20; VAR fx,dfx,OldRoot : REAL; BEGIN Root := GuessedRoot; fx := f(Root); Iter := 0; State := iterating; REPEAT dfx := df(Root); IF Iter = MaxIter THEN State := toomanyiter ELSIF ABS(Root) < NearZero THEN State := toonearzero ELSIF ABS(dfx) < NearZero THEN State := tooflat ELSE Iter := Iter + 1; OldRoot := Root; Root := Root - fx/dfx; fx := f(Root); IF ABS((Root-OldRoot)/OldRoot) < tol THEN State := rootfound END END UNTIL State # iterating; END Newton;