restart:# CLAMPED CUBIC SPLINE ALGORITHM 3.5## To construct the cubic spline interpolant S for the function f,# defined at the numbers x(0) < x(1) < ... < x(n), satisfying# S'(x(0)) = f'(x(0)) and S'(x(n)) = f'(x(n)):## INPUT: n: x(0), x(1), ..., x(n): either generate A(I) = f(x(I))# for i = 0, 1, ..., n or input A(I) for I = 0, 1, ..., n:# FP0 = f'(x(0)): FPN = f'(x(n)).## OUTPUT: A(J), B(J), C(J), D(J) for J = 0, 1, ..., n - 1.## NOTE: S(x) = A(J) + B(J) * ( x - x(J) ) + C(J) * ( x - x(J) )**2 +# D(J) * ( x - x(J) )**3 for x(J) <= x < x(J + 1)print(`This is Clamped Cubic Spline Interpolation.\134n`):OK := FALSE:while OK = FALSE doprint(`Choice of input method: `):print(`1. Input entry by entry from keyboard `):print(`2. Input data from a text file `):print(`3. Generate data using a function F with nodes entered `):print(`from keyboard `):print(`4. Generate data using a function F with nodes from `):print(`a text file `):print(`Choose 1, 2, 3, or 4 please `):FLAG := scanf(`%d`)[1]:print(`Entry = `):print(FLAG):if FLAG >= 1 and FLAG <= 4 thenOK := TRUE:fi:od:if FLAG = 1 thenOK := FALSE:while OK = FALSE doprint(`Input n `):N := scanf(`%d`)[1]:print(`n = `):print(N):if N > 0 thenOK := TRUE:for I1 from 0 to N doprint(`Input X(i) and F(X(i)) for i = `):print(I1):print(`separated by a space `):X[I1] := scanf(`%f`)[1]:A[I1] := scanf(`%f`)[1]:print(`Entries are `):print(X[I1],A[I1]):od:else print(`Number must be a positive integer `):fi:od:fi:if FLAG = 2 thenprint(`Has a text file been created with the data in two columns ? `):print(`Enter Y or N `):AA := scanf(` %c`)[1]:print(`Entry = `):print(AA):if AA = "Y" or AA = "y" thenprint(`Input the file name in the form - `):print(`drive:\134\134name.ext `):print(`For example: A:\134\134DATA.DTA `):NAME := scanf(`%s`)[1]:print(`Entry = `):print(NAME):INP := fopen(NAME,READ,TEXT):OK := FALSE:while OK = FALSE doprint(`Input n `):N := scanf(`%d`)[1]:print(`n = `):print(N):if N > 0 thenfor I1 from 0 to N doX[I1] := fscanf(INP, `%f`)[1]:A[I1] := fscanf(INP, `%f`)[1]:od:fclose(INP):OK := TRUE:else print(`Number must be a positive integer `):fi:od:elseprint(`Please create the input file in two column `):print(`form with the `):print(`X values and F(X) values in the corresponding columns. `):print(`The program will end so the input file can be created. `):OK := FALSE:fi:fi:if FLAG = 3 thenprint(`Input the function F(x) in terms of x `):print(`For example: cos(x) `):F := scanf(`%a`)[1]:print(`F(x) = `):print(F):F := unapply(F,x):OK := FALSE:while OK = FALSE doprint(`Input n `):N := scanf(`%d`)[1]:print(`n = `):print(N):if N > 0 thenfor I1 from 0 to N doprint(`Input X(i) for i = `):print(I1):X[I1] := scanf(`%f`)[1]:print(`Entry = `):print(X[I1]):A[I1] := F(X[I1]):od:OK := TRUE:elseprint(`Number must be a positive integer `):fi:od:fi:if FLAG = 4 thenprint(`Has the text file with X-values been created? `):print(`Enter Y or N `):AA := scanf(` %c`)[1]:print(`Entry = `):print(AA):if AA = "Y" or AA = "y" thenprint(`Input the file name in the form - `):print(`drive:\134\134name.ext `):print(`For example: A:\134\134DATA.DTA `):NAME := scanf(`%s`)[1]:print(`Entry = `):print(NAME):INP := fopen(NAME,READ,TEXT):print(`Input the function F(x) in terms of x `):print(`For example: cos(x) `):F := scanf(`%a`)[1]:print(`F(x) = `):print(F):F := unapply(F,x):OK := FALSE:while OK = FALSE doprint(`Input n `):N := scanf(`%d`)[1]:print(`n = `):print(N):if N > 0 thenOK := TRUE:for I1 from 0 to N doX[I1] := fscanf(INP, `%f`)[1]:A[I1] := F(X[I1]):od:fclose(INP):else print(`Number must be a positive integer `):fi:od:else print(`The program will end so the input file can be created`):OK := FALSE:fi:fi:if OK = TRUE thenprint(`Enter F'(X(0)) and F'(X(N)) separated by a space `):FP0 := scanf(`%f`)[1]:FPN := scanf(`%f`)[1]:print(`Entries are `):print(FP0,FPN):fi:if OK = TRUE thenM := N - 1:# Step 1for I1 from 0 to M doH[I1] := X[I1+1] - X[I1]:od:# Step 2# use XA instead of ALPHAXA[0] := 3.0 * (A[1] - A[0]) / H[0] - 3.0 * FP0:XA[N] := 3.0 * FPN - 3.0 * (A[N] - A[N-1]) / H[N-1]:# Step 3for I1 from 1 to M doXA[I1] := 3.0*(A[I1+1]*H[I1-1]-A[I1]*(X[I1+1]-X[I1-1])+A[I1-1]*H[I1])/(H[I1]*H[I1-1]):od:# Step 4# STEPS 4, 5, 6 and part of 7 solve the tridiagonal system using # Algorithm 6.7.# use XL, XU, XZ in place of L, MU, Z resp.XL[0] := 2.0 * H[0]:XU[0] := 0.5:XZ[0] := XA[0] / XL[0]:# Step 5for I1 from 1 to M doXL[I1] := 2.0 * (X[I1+1] - X[I1-1]) - H[I1-1] * XU[I1-1]:XU[I1] := H[I1] / XL[I1]:XZ[I1] := (XA[I1] - H[I1-1] * XZ[I1-1]) / XL[I1]:od:# Step 6XL[N] := H[N-1] * (2.0 - XU[N-1]):XZ[N] := (XA[N] - H[N-1] * XZ[N-1]) / XL[N]:C[N] := XZ[N]:# Step 7for I1 from 1 to N doJ := N - I1:C[J] := XZ[J] - XU[J] * C[J+1]:B[J] := (A[J+1] - A[J]) / H[J] - H[J] * (C[J+1] + 2.0 * C[J]) / 3.0:D1[J] := (C[J+1] - C[J]) / (3.0 * H[J]):od:# Step 8print(`Select output destination `):print(`1. Screen `):print(`2. Text file `):print(`Enter 1 or 2 `):FLAG := scanf(`%d`)[1]:print(`Entry = `):print(FLAG):if FLAG = 2 thenprint(`Input the file name in the form - drive:\134\134name.ext `):print(`For example: A:\134\134OUTPUT.DTA `):NAME := scanf(`%s`)[1]:print(`Entry = `):print(NAME):OUP := fopen(NAME,WRITE,TEXT):else OUP := default:fi:fprintf(OUP, `CLAMPED CUBIC SPLINE INTERPOLATION\134n\134n`):fprintf(OUP, `The numbers X(0), ..., X(N) are:\134n`):for I1 from 0 to N dofprintf(OUP, ` %11.8f`, X[I1]):od:fprintf(OUP, `\134n\134nThe coefficients of the spline on the subintervals `): fprintf(OUP, `are:\134n`):fprintf(OUP, `for I = 0, ..., N-1\134n`):fprintf(OUP, ` A(I) B(I) C(I) D(I)\134n`):for I1 from 0 to M dofprintf(OUP, ` %13.8f %13.8f %13.8f %13.8f \134n`, A[I1], B[I1], C[I1], D1[I1]):od:if OUP <> default thenfclose(OUP):print(`Output file `,NAME,` created successfully`):fi:fi:JSFH