alg058.txt

Numerical Anaysis 8th Edition Burden and Faires (Maple Source)

> restart;
> # TRAPEZOIDAL WITH NEWTON ITERATION ALGORITHM 5.8
> #
> # TO APPROXIMATE THE SOLUTION OF THE INITIAL VALUE PROBLEM:
> #            Y' = F(T,Y), A <= T <= B, Y(A) = ALPHA,
> # AT (N+1) EQUALLY SPACED NUMBERS IN THE INTERVAL [A,B].
> #
> # INPUT:   ENDPOINTS A,B; INITIAL CONDITION ALPHA; INTEGER N;
> #          TOLERANCE TOL; MAXIMUM NUMBER OF ITERATIONS M AT ANY ONE STEP.
> #
> # OUTPUT:  APPROXIMATION W TO Y AT THE (N+1) VALUES OF T
> #          OR A MESSAGE OF FAILURE.
> alg058 := proc() local F, FYP, OK, A, B, ALPHA, N, TOL, M, FLAG, NAME, OUP, W, T, H, I, XK1, W0, J, IFLAG;
> printf(`This is the Implicit Trapezoidal Method.\n`);
> printf(`Input the function F(t,y) in terms of t and y\n`);
> printf(`For example: y-t^2+1\n`);
> F := scanf(`%a`)[1];
> FYP := diff(F,y);
> F := unapply(F,t,y);
> FYP := unapply(FYP,t,y);
> OK := FALSE;
> while  OK = FALSE do
> printf(`Input left and right endpoints separated by blank.\n`);
> A := scanf(`%f`)[1];
> B := scanf(`%f`)[1];
> if A >= B then
> printf(`Left endpoint must be less than right endpoint.\n`);
> else
> OK := TRUE;
> fi;
> od;
> printf(`Input the initial condition.\n`);
> ALPHA := scanf(`%f`)[1];
> OK := FALSE;
> printf(`Input a positive integer for the number of subintervals.\n`);
> while OK = FALSE do
> N := scanf(`%d`)[1];
> if N <= 0 then
> printf(`Number must be a postiive integer.\n`);
> else
> OK := TRUE;
> fi;
> od;
> OK := FALSE;
> while OK = FALSE do
> printf(`Input tolerance.\n`);
> TOL := scanf(`%f`)[1];
> if TOL <= 0 then
> printf(`Tolerance must be positive.\n`);
> else
> OK := TRUE;
> fi;
> od;
> OK := FALSE;
> while OK = FALSE do
> printf(`Input maximum number of iterations.\n`);
> M := scanf(`%f`)[1];
> if M > 0 then
> OK := TRUE;
> else
> printf(`Number of iterations must be positive.\n`);
> fi;
> od;
> if OK = TRUE then
> printf(`Choice of output method:\n`);
> printf(`1. Output to screen\n`);
> printf(`2. Output to text file\n`);
> printf(`Please enter 1 or 2\n`);
> FLAG := scanf(`%d`)[1];
> if FLAG = 2 then
> printf(`Input the file name in the form - drive:\\name.ext\n`);
> printf(`For example   A:\\OUTPUT.DTA\n`);
> NAME := scanf(`%s`)[1];
> OUP := fopen(NAME,WRITE,TEXT);
> else
> OUP := default;
> fi;
> fprintf(OUP, `IMPLICIT TRAPEZOIDAL METHOD USING NEWTONS METHOD\n\n`);
> fprintf(OUP, `    t           w #iter\n`);
> # Step 1
> W := ALPHA;
> T := A;
> H := (B-A)/N;
> fprintf(OUP, `%5.3f %11.8f   0\n`, T, W);
> I := 1;
> OK := TRUE;
> # Step 2
> while I <= N and OK = TRUE do
> # Step 3
> XK1 := W+0.5*H*F(T, W);
> W0 := XK1;
> J := 1;
> IFLAG := 0;
> # Step 4
> while IFLAG = 0 and OK = TRUE do
> # Step 5
> W := W0-(W0-XK1-0.5*H*F(T+H, W0))/(1-0.5*H*FYP(T+H, W0));
> # Step 6
> if abs(W-W0) < TOL then
> # Step 7
> IFLAG := 1;
> T := A+I*H;
> fprintf(OUP,`%5.3f %11.8f %3d\n`, T, W, J);
> I := I+1;
> else
> J := J+1;
> W0 := W;
> if J > M then
> OK := FALSE;
> fi;
> fi;
> od;
> od;
> if OK = FALSE then
> fprintf(OUP, `Maximum Number of Iterations Exceeded\n`);
> fi;
> # Step 8
> if OUP <> default then
> fclose(OUP):
> printf(`Output file %s created successfully`,NAME);
> fi;
> fi;
> RETURN(0);
> end;
> alg058();