! F04QAF Example Program Text ! Mark 24 Release. NAG Copyright 2012. Module f04qafe_mod ! F04QAF Example Program Module: ! Parameters and User-defined Routines ! .. Use Statements .. Use nag_library, Only: nag_wp ! .. Implicit None Statement .. Implicit None ! .. Parameters .. Integer, Parameter :: iset = 1, liuser = 1, & lruser = 1, nin = 5, nout = 6 ! .. Local Scalars .. Integer :: ncols, nrows Contains Subroutine atimes(n,x,y) ! Called by routine aprod. Returns Y = Y + A*X, ! where A is not stored explicitly. ! .. Scalar Arguments .. Integer, Intent (In) :: n ! .. Array Arguments .. Real (Kind=nag_wp), Intent (In) :: x(n) Real (Kind=nag_wp), Intent (Inout) :: y(n) ! .. Local Scalars .. Integer :: i, i1, i2, i3, il, j ! .. Executable Statements .. Do j = 1, nrows - 2 y(j) = y(j) + x(j) - x(j+nrows-1) End Do Do j = 1, ncols - 2 i = j*nrows - 1 y(i) = y(i) + x(i) - x(i+1) i1 = i + 1 il = i1 + nrows - 3 Do i = i1, il i2 = i - nrows If (j==1) i2 = i2 + 1 i3 = i + nrows If (j==ncols-2) i3 = i3 - 1 y(i) = y(i) - x(i2) - x(i-1) + 4.0_nag_wp*x(i) - x(i+1) - x(i3) End Do i = il + 1 y(i) = y(i) - x(i-1) + x(i) End Do Do j = n - nrows + 3, n y(j) = y(j) - x(j-nrows+1) + x(j) End Do Return End Subroutine atimes Subroutine aprod(mode,m,n,x,y,ruser,lruser,iuser,liuser) ! APROD returns ! Y = Y + A*X when MODE = 1 ! X = X + ( A**T )*Y when MODE = 2 ! for a given X and Y. ! .. Scalar Arguments .. Integer, Intent (In) :: liuser, lruser, m, n Integer, Intent (Inout) :: mode ! .. Array Arguments .. Real (Kind=nag_wp), Intent (Inout) :: ruser(lruser), x(n), y(m) Integer, Intent (Inout) :: iuser(liuser) ! .. Local Scalars .. Integer :: j, j1, j2 ! .. Executable Statements .. If (mode/=2) Then Call atimes(n,x,y) Do j = 1, nrows - 2 y(m) = y(m) + x(j) End Do Do j = 1, ncols - 2 y(m) = y(m) + x(j*nrows-1) + x(j*nrows+nrows-2) End Do Do j = m - nrows + 2, n y(m) = y(m) + x(j) End Do Else Call atimes(n,y,x) Do j = 1, nrows - 2 x(j) = x(j) + y(m) End Do Do j = 1, ncols - 2 j1 = j*nrows - 1 j2 = j1 + nrows - 1 x(j1) = x(j1) + y(m) x(j2) = x(j2) + y(m) End Do Do j = m - nrows + 2, n x(j) = x(j) + y(m) End Do End If Return End Subroutine aprod End Module f04qafe_mod Program f04qafe ! F04QAF Example Main Program ! .. Use Statements .. Use nag_library, Only: f04qaf, nag_wp, x04abf Use f04qafe_mod, Only: aprod, iset, liuser, lruser, ncols, nin, nout, & nrows ! .. Implicit None Statement .. Implicit None ! .. Local Scalars .. Real (Kind=nag_wp) :: acond, anorm, arnorm, atol, & btol, c, conlim, damp, h, rnorm, & xnorm Integer :: i1, ifail, inform, itn, itnlim, & k, m, msglvl, n, outchn ! .. Local Arrays .. Real (Kind=nag_wp), Allocatable :: b(:), se(:), work(:,:), x(:) Real (Kind=nag_wp) :: ruser(lruser) Integer :: iuser(liuser) ! .. Executable Statements .. Write (nout,*) 'F04QAF Example Program Results' Write (nout,*) Flush (nout) ! Skip heading in data file Read (nin,*) Read (nin,*) nrows, ncols n = ncols*nrows - 4 m = n + 1 Allocate (b(m),se(n),work(n,2),x(n)) outchn = nout Call x04abf(iset,outchn) h = 0.1_nag_wp ! Initialize rhs and other quantities required by F04QAF. ! Convergence will be sooner if we do not regard A as exact, ! so atol is not set to zero. b(1:n) = 0.0_nag_wp c = -h**2 i1 = nrows Do k = 3, ncols b(i1:(i1+nrows-3)) = c i1 = i1 + nrows End Do b(m) = 1.0_nag_wp/h damp = 0.0_nag_wp atol = 1.0E-5_nag_wp btol = 1.0E-4_nag_wp conlim = 1.0_nag_wp/atol itnlim = 100 ! * Set msglvl to 2 to get output at each iteration * msglvl = 1 ! ifail: behaviour on error exit ! =0 for hard exit, =1 for quiet-soft, =-1 for noisy-soft ifail = 0 Call f04qaf(m,n,b,x,se,aprod,damp,atol,btol,conlim,itnlim,msglvl,itn, & anorm,acond,rnorm,arnorm,xnorm,work,ruser,lruser,iuser,liuser,inform, & ifail) Write (nout,*) Write (nout,*) 'Solution returned by F04QAF' Write (nout,99999) x(1:n) Write (nout,*) Write (nout,99998) 'Norm of the residual = ', rnorm 99999 Format (1X,5F9.3) 99998 Format (1X,A,1P,E12.2) End Program f04qafe