! D03PSF Example Program Text ! Mark 23 Release. NAG Copyright 2011. MODULE d03psfe_mod ! D03PSF 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 :: half = 0.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 :: itrace = 0, ncode = 0, nin = 5, & nout = 6, npde = 1, nxfix = 0, & nxi = 0 CONTAINS SUBROUTINE exact(t,u,x,npts) ! Exact solution (for comparison and b.c. purposes) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t INTEGER, INTENT (IN) :: npts ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: u(1,npts) REAL (KIND=nag_wp), INTENT (IN) :: x(*) ! .. Local Scalars .. REAL (KIND=nag_wp) :: del, psi, rm, rn, s INTEGER :: i ! .. Executable Statements .. s = 0.1_nag_wp del = 0.01_nag_wp rm = -one/del rn = one + s/del DO i = 1, npts psi = x(i) - t IF (psi(del+s)) THEN u(1,i) = zero ELSE u(1,i) = rm*psi + rn END IF END DO RETURN END SUBROUTINE exact SUBROUTINE uvin1(npde,npts,nxi,x,xi,u,ncode,v) ! .. Use Statements .. USE nag_library, ONLY : x01aaf ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. INTEGER, INTENT (IN) :: ncode, npde, npts, nxi ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: u(npde,npts), v(ncode) REAL (KIND=nag_wp), INTENT (IN) :: x(npts), xi(nxi) ! .. Local Scalars .. REAL (KIND=nag_wp) :: pi, tmp INTEGER :: i ! .. Intrinsic Functions .. INTRINSIC sin ! .. Executable Statements .. tmp = zero pi = x01aaf(tmp) DO i = 1, npts IF (x(i)>0.2_nag_wp .AND. x(i)<=0.4_nag_wp) THEN tmp = pi*(5.0_nag_wp*x(i)-one) u(1,i) = sin(tmp) ELSE u(1,i) = zero END IF END DO RETURN END SUBROUTINE uvin1 SUBROUTINE pdef1(npde,t,x,u,ux,ncode,v,vdot,p,c,d,s,ires) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t, x INTEGER, INTENT (INOUT) :: ires INTEGER, INTENT (IN) :: ncode, npde ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: c(npde), d(npde), & p(npde,npde), s(npde) REAL (KIND=nag_wp), INTENT (IN) :: u(npde), ux(npde), v(ncode), & vdot(ncode) ! .. Executable Statements .. p(1,1) = one c(1) = 0.002_nag_wp d(1) = ux(1) s(1) = zero RETURN END SUBROUTINE pdef1 SUBROUTINE bndry1(npde,npts,t,x,u,ncode,v,vdot,ibnd,g,ires) ! Zero solution at both boundaries ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t INTEGER, INTENT (IN) :: ibnd, ncode, npde, npts INTEGER, INTENT (INOUT) :: ires ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: g(npde) REAL (KIND=nag_wp), INTENT (IN) :: u(npde,npts), v(ncode), & vdot(ncode), x(npts) ! .. Executable Statements .. IF (ibnd==0) THEN g(1) = u(1,1) ELSE g(1) = u(1,npts) END IF RETURN END SUBROUTINE bndry1 SUBROUTINE monit1(t,npts,npde,x,u,fmon) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t INTEGER, INTENT (IN) :: npde, npts ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: fmon(npts) REAL (KIND=nag_wp), INTENT (IN) :: u(npde,npts), x(npts) ! .. Local Scalars .. REAL (KIND=nag_wp) :: h1, h2, h3 INTEGER :: i ! .. Intrinsic Functions .. INTRINSIC abs ! .. Executable Statements .. DO i = 2, npts - 1 h1 = x(i) - x(i-1) h2 = x(i+1) - x(i) h3 = half*(x(i+1)-x(i-1)) ! Second derivatives .. fmon(i) = abs(((u(1,i+1)-u(1,i))/h2-(u(1,i)-u(1,i-1))/h1)/h3) END DO fmon(1) = fmon(2) fmon(npts) = fmon(npts-1) RETURN END SUBROUTINE monit1 SUBROUTINE nmflx1(npde,t,x,ncode,v,uleft,uright,flux,ires) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t, x INTEGER, INTENT (INOUT) :: ires INTEGER, INTENT (IN) :: ncode, npde ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: flux(npde) REAL (KIND=nag_wp), INTENT (IN) :: uleft(npde), uright(npde), & v(ncode) ! .. Executable Statements .. flux(1) = uleft(1) RETURN END SUBROUTINE nmflx1 SUBROUTINE uvin2(npde,npts,nxi,x,xi,u,ncode,v) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. INTEGER, INTENT (IN) :: ncode, npde, npts, nxi ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: u(npde,npts), v(ncode) REAL (KIND=nag_wp), INTENT (IN) :: x(npts), xi(nxi) ! .. Local Scalars .. REAL (KIND=nag_wp) :: t ! .. Executable Statements .. t = zero CALL exact(t,u,x,npts) RETURN END SUBROUTINE uvin2 SUBROUTINE pdef2(npde,t,x,u,ux,ncode,v,vdot,p,c,d,s,ires) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t, x INTEGER, INTENT (INOUT) :: ires INTEGER, INTENT (IN) :: ncode, npde ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: c(npde), d(npde), & p(npde,npde), s(npde) REAL (KIND=nag_wp), INTENT (IN) :: u(npde), ux(npde), v(ncode), & vdot(ncode) ! .. Executable Statements .. p(1,1) = one c(1) = zero d(1) = zero s(1) = -100.0_nag_wp*u(1)*(u(1)-one)*(u(1)-half) RETURN END SUBROUTINE pdef2 SUBROUTINE bndry2(npde,npts,t,x,u,ncode,v,vdot,ibnd,g,ires) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t INTEGER, INTENT (IN) :: ibnd, ncode, npde, npts INTEGER, INTENT (INOUT) :: ires ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: g(npde) REAL (KIND=nag_wp), INTENT (IN) :: u(npde,npts), v(ncode), & vdot(ncode), x(npts) ! .. Local Arrays .. REAL (KIND=nag_wp) :: ue(1,1) ! .. Executable Statements .. ! Solution known to be constant at both boundaries IF (ibnd==0) THEN CALL exact(t,ue,x(1),1) g(1) = ue(1,1) - u(1,1) ELSE CALL exact(t,ue,x(npts),1) g(1) = ue(1,1) - u(1,npts) END IF RETURN END SUBROUTINE bndry2 SUBROUTINE nmflx2(npde,t,x,ncode,v,uleft,uright,flux,ires) ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t, x INTEGER, INTENT (INOUT) :: ires INTEGER, INTENT (IN) :: ncode, npde ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: flux(npde) REAL (KIND=nag_wp), INTENT (IN) :: uleft(npde), uright(npde), & v(ncode) ! .. Executable Statements .. flux(1) = uleft(1) RETURN END SUBROUTINE nmflx2 SUBROUTINE monit2(t,npts,npde,x,u,fmon) ! .. Use Statements .. USE nag_library, ONLY : x01aaf ! .. Implicit None Statement .. IMPLICIT NONE ! .. Scalar Arguments .. REAL (KIND=nag_wp), INTENT (IN) :: t INTEGER, INTENT (IN) :: npde, npts ! .. Array Arguments .. REAL (KIND=nag_wp), INTENT (OUT) :: fmon(npts) REAL (KIND=nag_wp), INTENT (IN) :: u(npde,npts), x(npts) ! .. Local Scalars .. REAL (KIND=nag_wp) :: h1, pi, ux, uxmax, xl, xleft, & xmax, xr, xright, xx REAL (KIND=nag_wp), SAVE :: xa = zero INTEGER :: i INTEGER, SAVE :: icount = 0 ! .. Intrinsic Functions .. INTRINSIC abs, cos ! .. Executable Statements .. xx = zero pi = x01aaf(xx) ! Locate shock .. uxmax = zero xmax = zero DO i = 2, npts - 1 h1 = x(i) - x(i-1) ux = abs((u(1,i)-u(1,i-1))/h1) IF (ux>uxmax) THEN uxmax = ux xmax = x(i) END IF END DO ! Assign width (on first call only) .. IF (icount==0) THEN icount = 1 xleft = xmax - x(1) xright = x(npts) - xmax IF (xleft>xright) THEN xa = xright ELSE xa = xleft END IF END IF xl = xmax - xa xr = xmax + xa ! Assign monitor function .. DO i = 1, npts IF (x(i)>xl .AND. x(i)