> 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();