! F04QAF Example Program Text
! Mark 30.2 Release. NAG Copyright 2024.
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
! .. Accessibility Statements ..
Private
Public :: aprod
! .. Parameters ..
Integer, Parameter, Public :: iset = 1, liuser = 0, lruser = 0, &
nin = 5, nout = 6
! .. Local Scalars ..
Integer, Public, Save :: 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) Then
i2 = i2 + 1
End If
i3 = i + nrows
If (j==ncols-2) Then
i3 = i3 - 1
End If
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 f04qafe_mod, Only: aprod, iset, liuser, lruser, ncols, nin, nout, &
nrows
Use nag_library, Only: f04qaf, nag_wp, x04abf
! .. 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