alg067.txt

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

> restart;
> # CROUT FACTORIZATION FOR TRIDIAGONAL LINEAR SYSTEMS ALGORITHM 6.7
> #
> # To solve the n x n linear system
> #
> # E1:  A[1,1] X[1] + A[1,2] X[2]                  = A[1,n+1]
> # E2:  A[2,1] X[1] + A[2,2] X[2] + A[2,3] X[3]    = A[2,n+1]
> # :
> # .
> # E(n):          A[n,n-1] X[n-1] + A[n,n] X[n]    = A[n,n+1]
> #
> # INPUT:   the dimension n; the entries of A.
> #
> # OUTPUT:  the solution X(1), ..., X(N).
> alg067 := proc() local AA, OK, NAME, INP, N, I, A, B, NN, C, BB, Z, X, II, FLAG, OUP;
> printf(`This is Crout Method for tridiagonal linear systems.\n`);
> printf(`The array will be input from a text file in the order:\n`);
> printf(`diagonal entries, lower sub-diagonal entries,\n`);
> printf(`upper sub-diagonal entries, inhomogeneous term.\n\n`);
> printf(`Place as many entries as desired on each line, but separate `);
> printf(`entries \n`);
> printf(`with at least one blank.\n\n\n`);
> printf(`Has the input file been created? - enter Y or N.\n`);
> AA := scanf(`%c`)[1];
> OK := FALSE;
> if AA = "Y" or AA = "y" then
> printf(`Input the file name in the form - drive:\\name.ext\n`);
> printf(`for example:   A:\\DATA.DTA\n`);
> NAME := scanf(`%s`)[1];
> INP := fopen(NAME,READ,TEXT);
> OK := FALSE;
> while OK = FALSE do
> printf(`Input the number of equations - an integer.\n`);
> N := scanf(`%d`)[1];
> if N > 0 then
> # A(I,I) is stored in A(I), 1 <= I <= n 
> for I from 1 to N do
> A[I-1] := fscanf(INP, `%f`)[1];
> od;
> # the lower sub-diagonal A(I,I-1) is stored
> #in B(I), 2 <= I <= n 
> for I from 2 to N do
> B[I-1] := fscanf(INP, `%f`)[1];
> od;
> # the upper sub-diagonal A(I,I+1) is stored
> #in C(I), 1 <= I <= n-1 
> NN := N-1;
> for I from 1 to NN do
> C[I-1] := fscanf(INP, `%f`)[1];
> od;
> # A(I,N+1) is stored in BB(I), 1 <= I <= n 
> for I from 1 to N do
> BB[I-1] := fscanf(INP, `%f`)[1];
> od;
> OK := TRUE;
> fclose(INP);
> else
> printf(`The number must be a positive integer.\n`);
> fi;
> od;
> else
> printf(`The program will end so the input file can be created.\n`)
> fi;
> if OK = TRUE then
> # Steps 1 - 3 set up and solve Lz = 0
> # Step 1
> # The entries of U overwrite C and the entries of L overwrite A
> C[0] := C[0]/A[0];
> Z[0] := BB[0]/A[0];
> # Step 2
> for I from 2 to NN do
> A[I-1] := A[I-1]-B[I-1]*C[I-2];
> C[I-1] := C[I-1]/A[I-1];
> Z[I-1] := (BB[I-1]-B[I-1]*Z[I-2])/A[I-1];
> od;
> # Step 3
> A[N-1] := A[N-1]-B[N-1]*C[N-2];
> Z[N-1] := (BB[N-1]-B[N-1]*Z[N-2])/A[N-1];
> # Step 4
> # Steps 4 and 5 solve Ux = z
> X[N-1] := Z[N-1];
> # Step 5
> for II from 1 to NN do
> I := NN-II+1;
> X[I-1] := Z[I-1]-C[I-1]*X[I];
> od;
> # Step 6
> 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, `CROUT METHOD FOR TRIDIAGONAL LINEAR SYSTEMS\n\n`);
> fprintf(OUP, `The solution is\n`);
> for I from 1 to N do
> fprintf(OUP, `  %12.8f`, X[I-1]);
> od;
> fprintf(OUP, `\n`);
> if OUP <> default then
> fclose(OUP):
> printf(`Output file %s created successfully`,NAME);
> fi;
> fi;
> RETURN(0);
> end;
> alg067();