! E02DHF Example Program Text ! Mark 24 Release. NAG Copyright 2012. Module e02dhfe_mod ! E02DHF Example Program Module: ! Parameters and User-defined Routines ! .. Use Statements .. Use nag_library, Only: nag_wp ! .. Implicit None Statement .. Implicit None ! .. Parameters .. Integer, Parameter :: nin = 5, nout = 6 Contains Subroutine print_spline(ngx,gridx,ngy,gridy,z,zder) ! Print spline function and spline derivative evaluation ! .. Use Statements .. Use nag_library, Only: x04cbf ! .. Parameters .. Integer, Parameter :: indent = 0, ncols = 80 Character (1), Parameter :: chlabel = 'C', diag = 'N', & matrix = 'G' Character (4), Parameter :: form = 'F8.3' ! .. Scalar Arguments .. Integer, Intent (In) :: ngx, ngy ! .. Array Arguments .. Real (Kind=nag_wp), Intent (In) :: gridx(ngx), gridy(ngy), & z(ngx*ngy), zder(ngx*ngy) ! .. Local Scalars .. Integer :: i, ifail Character (48) :: title ! .. Local Arrays .. Character (10), Allocatable :: clabs(:), rlabs(:) ! .. Executable Statements .. ! Allocate for row and column label Allocate (clabs(ngx),rlabs(ngy)) ! Generate column and row labels to print the results with. Do i = 1, ngx Write (clabs(i),99999) gridx(i) End Do Do i = 1, ngy Write (rlabs(i),99999) gridy(i) End Do ! Print the spline evaluations. title = 'Spline evaluated on X-Y grid (X across, Y down):' Write (nout,*) Flush (nout) ifail = 0 Call x04cbf(matrix,diag,ngy,ngx,z,ngy,form,title,chlabel,rlabs, & chlabel,clabs,ncols,indent,ifail) ! Print the spline derivative evaluations. title = 'Spline derivative evaluated on X-Y grid:' Write (nout,*) Flush (nout) ifail = 0 Call x04cbf(matrix,diag,ngy,ngx,zder,ngy,form,title(1:40),chlabel, & rlabs,chlabel,clabs,ncols,indent,ifail) Deallocate (clabs,rlabs) 99999 Format (F5.2) End Subroutine print_spline End Module e02dhfe_mod Program e02dhfe ! E02DHF Example Main Program ! .. Use Statements .. Use nag_library, Only: e02dcf, e02dhf, nag_wp Use e02dhfe_mod, Only: nin, nout, print_spline ! .. Implicit None Statement .. Implicit None ! .. Local Scalars .. Real (Kind=nag_wp) :: delta, fp, s, xhi, xlo, yhi, ylo Integer :: i, ifail, liwrk, lwrk, mx, my, & nc, ngx, ngy, nux, nuy, nx, & nxest, ny, nyest Character (1) :: start ! .. Local Arrays .. Real (Kind=nag_wp), Allocatable :: c(:), f(:), gridx(:), gridy(:), & lamda(:), mu(:), wrk(:), x(:), & y(:), z(:), zder(:) Integer, Allocatable :: iwrk(:) ! .. Intrinsic Procedures .. Intrinsic :: max, real ! .. Executable Statements .. Write (nout,*) 'E02DHF Example Program Results' ! Skip heading in data file Read (nin,*) ! Input the number of X, Y co-ordinates MX, MY. Read (nin,*) mx, my nxest = mx + 4 nyest = my + 4 nc = (nxest-4)*(nyest-4) ! Allocations for spline fit Allocate (lamda(nxest),mu(nyest),c(nc)) ! Allocations for e02dcf only lwrk = 4*(mx+my) + 11*(nxest+nyest) + nxest*my + max(my,nxest) + 54 liwrk = 3 + mx + my + nxest + nyest Allocate (x(mx),y(my),f(mx*my),wrk(lwrk),iwrk(liwrk)) Read (nin,*) x(1:mx) Read (nin,*) y(1:my) ! Input the MX*MY function values F at grid points nad smoothing factor. Read (nin,*) f(1:mx*my) Read (nin,*) s ! Determine the spline approximation. start = 'C' ifail = 0 Call e02dcf(start,mx,x,my,y,f,s,nxest,nyest,nx,lamda,ny,mu,c,fp,wrk, & lwrk,iwrk,liwrk,ifail) Deallocate (x,y,f,wrk,iwrk) Write (nout,*) Write (nout,99999) 'Spline fit used smoothing factor S =', s, '.' Write (nout,99998) 'Number of knots in each direction =', nx, ny Write (nout,*) Write (nout,99999) 'Sum of squared residuals =', fp, '.' ! Spline and its derivative to be evaluated on rectangular grid with ! ngx*ngy points on the domain [xlo,xhi]x[ylo,yhi]. Read (nin,*) ngx, xlo, xhi Read (nin,*) ngy, ylo, yhi ! Allocations for e02dhf (spline evaluation). Allocate (gridx(ngx),gridy(ngy),z(ngx*ngy),zder(ngx*ngy)) delta = (xhi-xlo)/real(ngx-1,kind=nag_wp) gridx(1) = xlo Do i = 2, ngx - 1 gridx(i) = gridx(i-1) + delta End Do gridx(ngx) = xhi delta = (yhi-ylo)/real(ngy-1,kind=nag_wp) gridy(1) = ylo Do i = 2, ngy - 1 gridy(i) = gridy(i-1) + delta End Do gridy(ngy) = yhi ! Evaluate spline (nux=nuy=0) nux = 0 nuy = 0 ifail = 0 Call e02dhf(ngx,ngy,nx,ny,gridx,gridy,lamda,mu,c,nux,nuy,z,ifail) ! Evaluate spline partial derivative of order (nux,nuy) Read (nin,*) nux, nuy Write (nout,*) Write (nout,99998) 'Derivative of spline has order nux, nuy =', nux, nuy ifail = 0 Call e02dhf(ngx,ngy,nx,ny,gridx,gridy,lamda,mu,c,nux,nuy,zder,ifail) ! Print tabulated spline and derivative evaluations. Call print_spline(ngx,gridx,ngy,gridy,z,zder) 99999 Format (1X,A,1P,E13.4,A) 99998 Format (1X,A,I5,',',I5,'.') End Program e02dhfe