160 likes | 278 Views
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 ?
E N D
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;