! D02SAF Example Program Text
! Mark 28.3 Release. NAG Copyright 2022.
Module d02safe_mod
! D02SAF 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 :: bc, constr, eqn, fcn, prsol, range
! .. Parameters ..
Real (Kind=nag_wp), Parameter :: alpha = 0.032_nag_wp
Real (Kind=nag_wp), Parameter :: beta = 0.02_nag_wp
Real (Kind=nag_wp), Parameter :: xend = 5.0_nag_wp
Integer, Parameter, Public :: iset = 1, m = 4, n = 3, nin = 5, &
nout = 6
Contains
Subroutine eqn(e,q,p,m)
! .. Scalar Arguments ..
Integer, Intent (In) :: m, q
! .. Array Arguments ..
Real (Kind=nag_wp), Intent (Out) :: e(q)
Real (Kind=nag_wp), Intent (In) :: p(m)
! .. Executable Statements ..
e(1) = 0.02_nag_wp - p(4) - 1.0E-5_nag_wp*p(3)
Return
End Subroutine eqn
Subroutine fcn(x,y,f,n,p,m,i)
! .. Scalar Arguments ..
Real (Kind=nag_wp), Intent (In) :: x
Integer, Intent (In) :: i, m, n
! .. Array Arguments ..
Real (Kind=nag_wp), Intent (Out) :: f(n)
Real (Kind=nag_wp), Intent (In) :: p(m), y(n)
! .. Intrinsic Procedures ..
Intrinsic :: cos, tan
! .. Executable Statements ..
f(1) = tan(y(3))
If (i==1) Then
f(2) = -alpha*tan(y(3))/y(2) - beta*y(2)/cos(y(3))
f(3) = -alpha/y(2)**2
Else
f(2) = -p(2)*tan(y(3))/y(2) - p(4)*y(2)/cos(y(3))
f(3) = -p(2)/y(2)**2
End If
Return
End Subroutine fcn
Subroutine bc(g1,g2,p,m,n)
! .. Scalar Arguments ..
Integer, Intent (In) :: m, n
! .. Array Arguments ..
Real (Kind=nag_wp), Intent (Out) :: g1(n), g2(n)
Real (Kind=nag_wp), Intent (In) :: p(m)
! .. Executable Statements ..
g1(1) = 0.0_nag_wp
g1(2) = 0.5_nag_wp
g1(3) = p(1)
g2(1) = 0.0_nag_wp
g2(2) = 0.45_nag_wp
g2(3) = -1.2_nag_wp
Return
End Subroutine bc
Subroutine range(x,npoint,p,m)
! .. Scalar Arguments ..
Integer, Intent (In) :: m, npoint
! .. Array Arguments ..
Real (Kind=nag_wp), Intent (In) :: p(m)
Real (Kind=nag_wp), Intent (Out) :: x(npoint)
! .. Executable Statements ..
x(1) = 0.0_nag_wp
x(2) = p(3)
x(3) = xend
Return
End Subroutine range
Subroutine prsol(z,y,n)
! .. Scalar Arguments ..
Real (Kind=nag_wp), Intent (Inout) :: z
Integer, Intent (In) :: n
! .. Array Arguments ..
Real (Kind=nag_wp), Intent (In) :: y(n)
! .. Local Scalars ..
Integer :: i
! .. Intrinsic Procedures ..
Intrinsic :: abs
! .. Executable Statements ..
If (z==0.0E0_nag_wp) Then
Write (nout,*)
Write (nout,*) ' Z Y(1) Y(2) Y(3)'
End If
Write (nout,99999) z, (y(i),i=1,n)
z = z + 0.5_nag_wp
If (abs(z-xend)<0.25_nag_wp) Then
z = xend
End If
Return
99999 Format (1X,F9.3,3F10.4)
End Subroutine prsol
Function constr(p,m)
! .. Function Return Value ..
Logical :: constr
! .. Scalar Arguments ..
Integer, Intent (In) :: m
! .. Array Arguments ..
Real (Kind=nag_wp), Intent (In) :: p(m)
! .. Intrinsic Procedures ..
Intrinsic :: any
! .. Executable Statements ..
If (any(p(1:m)<0.0_nag_wp) .Or. p(3)>5.0_nag_wp) Then
constr = .False.
Else
constr = .True.
End If
Return
End Function constr
End Module d02safe_mod
Program d02safe
! D02SAF Example Main Program
! .. Use Statements ..
Use d02safe_mod, Only: bc, constr, eqn, fcn, iset, m, n, nin, nout, &
prsol, range
Use nag_library, Only: d02saf, d02sas, nag_wp, x04abf
! .. Implicit None Statement ..
Implicit None
! .. Local Scalars ..
Real (Kind=nag_wp) :: ymax
Integer :: i, icount, ifail, ldswp, ldw, n1, &
npoint, outchn, sdw
! .. Local Arrays ..
Real (Kind=nag_wp), Allocatable :: dp(:), e(:), p(:), pe(:), pf(:), &
swp(:,:), w(:,:)
! .. Intrinsic Procedures ..
Intrinsic :: max
! .. Executable Statements ..
Write (nout,*) 'D02SAF Example Program Results'
! Skip heading in data file
Read (nin,*)
Read (nin,*) npoint
n1 = n
sdw = 3*m + 23
ldswp = npoint
ldw = max(m,n)
Allocate (dp(m),e(n),p(m),pe(m),pf(m),swp(ldswp,6),w(ldw,sdw))
outchn = nout
Read (nin,*) icount
Read (nin,*) ymax
Read (nin,*) pe(1:m)
Read (nin,*) pf(1:m)
Read (nin,*) dp(1:m)
Read (nin,*) e(1:n)
Call x04abf(iset,outchn)
swp(1:npoint-1,1:3) = 0.0_nag_wp
Read (nin,*) p(1:m)
! * To obtain monitoring information, replace the name d02sas by d02hbx
! in the next statement and USE nag_library : d02hbx
! ifail: behaviour on error exit
! =0 for hard exit, =1 for quiet-soft, =-1 for noisy-soft
ifail = 1
Call d02saf(p,m,n,n1,pe,pf,e,dp,npoint,swp,ldswp,icount,range,bc,fcn, &
eqn,constr,ymax,d02sas,prsol,w,ldw,sdw,ifail)
If (ifail/=0) Then
Write (nout,99999) ifail
End If
If (ifail>=4 .And. ifail<=12) Then
Write (nout,99998) 'SWP(NPOINT,1) = ', swp(npoint,1)
If (ifail<=6) Then
Write (nout,99998) 'SWP(1,5) = ', swp(1,5)
Write (nout,*) ' i W(i,1) '
Write (nout,99997)(i,w(i,1),i=1,n)
End If
End If
99999 Format (1X,/,1X,' ** D02SAF returned with IFAIL = ',I5)
99998 Format (1X,A,F10.4)
99997 Format (1X,I4,1X,E10.3)
End Program d02safe