! D02UCF Example Program Text ! Mark 23 Release. NAG Copyright 2011. MODULE d02ucfe_mod ! D02UCF Example Program Module: ! Parameters and User-defined Routines ! .. Use Statements .. USE nag_library, ONLY : nag_wp ! .. Implicit None Statement .. IMPLICIT NONE ! .. Parameters .. REAL (KIND=nag_wp), PARAMETER :: a = 0.0_nag_wp REAL (KIND=nag_wp), PARAMETER :: b = 1.5_nag_wp REAL (KIND=nag_wp), PARAMETER :: one = 1.0_nag_wp REAL (KIND=nag_wp), PARAMETER :: zero = 0.0_nag_wp INTEGER, PARAMETER :: m = 1, nin = 5, nout = 6 LOGICAL, PARAMETER :: reqerr = .FALSE. CONTAINS FUNCTION exact(x) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Function Return Value .. REAL (KIND=nag_wp) :: exact ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: x ! .. Intrinsic Functions .. INTRINSIC exp ! .. Executable Statements .. exact = exp(-x-one) + one RETURN END FUNCTION exact FUNCTION deriv(x) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Function Return Value .. REAL (KIND=nag_wp) :: deriv ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: x ! .. Intrinsic Functions .. INTRINSIC exp ! ..Executable Statements .. deriv = -exp(-x-one) RETURN END FUNCTION deriv SUBROUTINE bndary(m,y,bmat,bvec) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. INTEGER, INTENT (IN) :: m ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: bmat(m,m+1), bvec(m), y(m) ! .. Executable Statements .. ! Boundary condition of left side of domain y(1) = a ! Set up Dirichlet condition using exact solution at x=a bmat(1:m,1:m+1) = zero bmat(1,1) = one bvec(1) = exact(a) RETURN END SUBROUTINE bndary SUBROUTINE pdedef(m,f) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. INTEGER, INTENT (IN) :: m ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: f(m+1) ! .. Executable Statements .. f(1) = one f(2) = one RETURN END SUBROUTINE pdedef END MODULE d02ucfe_mod PROGRAM d02ucfe ! D02UCF Example Main Program ! .. Use Statements .. USE nag_library, ONLY : d02uaf, d02ubf, d02ucf, d02uef, x02ajf USE d02ucfe_mod, ONLY : a, b, bndary, deriv, exact, m, nag_wp, nin, & nout, one, pdedef, reqerr, zero ! .. Implicit None Statement .. IMPLICIT NONE ! .. Local Scalars .. REAL (KIND=nag_wp) :: resid, teneps, uerr, uxerr INTEGER :: i, ifail, n ! .. Local Arrays .. REAL (KIND=nag_wp) :: bmat(m,m+1), bvec(1:m), f(m+1), & y(m) REAL (KIND=nag_wp), ALLOCATABLE :: c(:), f0(:), u(:), uc(:,:), & ux(:), x(:) ! .. Intrinsic Functions .. INTRINSIC abs, int, max ! .. Executable Statements .. WRITE (nout,*) ' D02UCF Example Program Results ' WRITE (nout,*) READ (nin,*) READ (nin,*) n ALLOCATE (f0(n+1),c(n+1),u(n+1),ux(n+1),uc(n+1,m+1),x(n+1)) ! Set up problem boundary conditions and definition CALL bndary(m,y,bmat,bvec) CALL pdedef(m,f) ! Set up solution grid ifail = 0 CALL d02ucf(n,a,b,x,ifail) ! Set up problem right hand sides for grid and transform f0(1:n+1) = one CALL d02uaf(n,f0,c,ifail) ! Solve in coefficient space ifail = 0 CALL d02uef(n,a,b,m,c,bmat,y,bvec,f,uc,resid,ifail) ! Transform solution and derivative back to real space ifail = 0 CALL d02ubf(n,a,b,0,uc(1,1),u,ifail) ifail = 0 CALL d02ubf(n,a,b,1,uc(1,2),ux,ifail) ! Print solution WRITE (nout,*) ' Numerical solution U and derivative Ux' WRITE (nout,*) WRITE (nout,99999) WRITE (nout,99998) (x(i),u(i),ux(i),i=1,n+1) IF (reqerr) THEN uerr = zero uxerr = zero DO i = 1, n + 1 uerr = max(uerr,abs(u(i)-exact(x(i)))) uxerr = max(uxerr,abs(ux(i)-deriv(x(i)))) END DO teneps = 10.0_nag_wp*x02ajf() WRITE (nout,99997) 10*(int(uerr/teneps)+1) WRITE (nout,99996) 10*(int(uxerr/teneps)+1) END IF 99999 FORMAT (1X,T8,'X',T18,'U',T28,'Ux') 99998 FORMAT (1X,3F10.4) 99997 FORMAT (//1X,'U is within a multiple ',I8,' of machine precision.') 99996 FORMAT (1X,'Ux is within a multiple ',I8,' of machine precision.') END PROGRAM d02ucfe