diff options
author | diener <diener@f69c4107-0314-4c4f-9ad4-17e986b73f4a> | 2006-04-13 22:07:59 +0000 |
---|---|---|
committer | diener <diener@f69c4107-0314-4c4f-9ad4-17e986b73f4a> | 2006-04-13 22:07:59 +0000 |
commit | 690f82fec96d844bc221ba6d6d0c1c5d76bc1552 (patch) | |
tree | 9032e7e2baa1d3de41b9248add47aa71a2422d4f /src/Derivatives_2_1.F90 | |
parent | eb16e69915d6e5c73b3e9fe7b4133f6aaf6631a7 (diff) |
Modify the interface to allow for spatial dependence on the upwinding
direction for the upwinding derivatives. The directions are specified in a
real 3d array with the sign choosing the upwinding direction.
git-svn-id: https://svn.cct.lsu.edu/repos/numrel/LSUThorns/SummationByParts/trunk@64 f69c4107-0314-4c4f-9ad4-17e986b73f4a
Diffstat (limited to 'src/Derivatives_2_1.F90')
-rw-r--r-- | src/Derivatives_2_1.F90 | 197 |
1 files changed, 112 insertions, 85 deletions
diff --git a/src/Derivatives_2_1.F90 b/src/Derivatives_2_1.F90 index ed97d9b..4a0adc1 100644 --- a/src/Derivatives_2_1.F90 +++ b/src/Derivatives_2_1.F90 @@ -120,13 +120,11 @@ subroutine up_deriv_gf_2_1 ( var, ni, nj, nk, dir, bb, gsize, delta, up, dvar ) CCTK_INT, intent(IN) :: bb(2) CCTK_INT, intent(IN) :: gsize CCTK_REAL, intent(IN) :: delta - CCTK_INT, intent(IN) :: up + CCTK_REAL, dimension(ni,nj,nk), intent(IN) :: up CCTK_REAL, dimension(ni,nj,nk), intent(OUT) :: dvar CCTK_REAL, dimension(-1:1), save :: a1, a2 - CCTK_REAL, dimension(-1:1) :: a CCTK_REAL, dimension(3,2), save :: q1, q2 - CCTK_REAL, dimension(3,2) :: q CCTK_REAL :: idel CCTK_INT :: il, ir, jl, jr, kl, kr @@ -161,50 +159,55 @@ subroutine up_deriv_gf_2_1 ( var, ni, nj, nk, dir, bb, gsize, delta, up, dvar ) idel = 1.0_wp / delta - upwinding: select case (up) - case (-1) upwinding - a = a1 - case (1) upwinding - a = a2 - end select upwinding - direction: select case (dir) case (0) direction if ( bb(1) == 0 ) then il = 1 + gsize else - upwinding_xl: select case (up) - case (-1) upwinding_xl - q = q1 - case (1) upwinding_xl - q = q2 - end select upwinding_xl - - dvar(1,:,:) = ( q(1,1) * var(1,:,:) + q(2,1) * var(2,:,:) ) * idel - dvar(2,:,:) = ( q(1,2) * var(1,:,:) + q(2,2) * var(2,:,:) + & - q(3,2) * var(3,:,:) ) * idel + where ( up(1,:,:) < zero ) + dvar(1,:,:) = ( q1(1,1) * var(1,:,:) + q1(2,1) * var(2,:,:) ) * idel + elsewhere + dvar(1,:,:) = ( q2(1,1) * var(1,:,:) + q2(2,1) * var(2,:,:) ) * idel + end where + where ( up(2,:,:) < zero ) + dvar(2,:,:) = ( q1(1,2) * var(1,:,:) + q1(2,2) * var(2,:,:) + & + q1(3,2) * var(3,:,:) ) * idel + elsewhere + dvar(2,:,:) = ( q2(1,2) * var(1,:,:) + q2(2,2) * var(2,:,:) + & + q2(3,2) * var(3,:,:) ) * idel + end where il = 3 end if if ( bb(2) == 0 ) then ir = ni - gsize else - upwinding_xr: select case (up) - case (-1) upwinding_xr - q = q2 - case (1) upwinding_xr - q = q1 - end select upwinding_xr - - dvar(ni-1,:,:) = - ( q(1,2) * var(ni,:,:) + & - q(2,2) * var(ni-1,:,:) + & - q(3,2) * var(ni-2,:,:) ) * idel - dvar(ni,:,:) = - ( q(1,1) * var(ni,:,:) + & - q(2,1) * var(ni-1,:,:) ) * idel + where ( up(ni-1,:,:) < zero ) + dvar(ni-1,:,:) = - ( q2(1,2) * var(ni,:,:) + & + q2(2,2) * var(ni-1,:,:) + & + q2(3,2) * var(ni-2,:,:) ) * idel + elsewhere + dvar(ni-1,:,:) = - ( q1(1,2) * var(ni,:,:) + & + q1(2,2) * var(ni-1,:,:) + & + q1(3,2) * var(ni-2,:,:) ) * idel + end where + where ( up(ni,:,:) < zero ) + dvar(ni,:,:) = - ( q2(1,1) * var(ni,:,:) + & + q2(2,1) * var(ni-1,:,:) ) * idel + elsewhere + dvar(ni,:,:) = - ( q1(1,1) * var(ni,:,:) + & + q1(2,1) * var(ni-1,:,:) ) * idel + end where ir = ni - 2 end if - dvar(il:ir,:,:) = ( a(-1) * var(il-1:ir-1,:,:) + & - a(0) * var(il:ir,:,:) + & - a(1) * var(il+1:ir+1,:,:) ) * idel + where ( up(il:ir,:,:) < zero ) + dvar(il:ir,:,:) = ( a1(-1) * var(il-1:ir-1,:,:) + & + a1(0) * var(il:ir,:,:) + & + a1(1) * var(il+1:ir+1,:,:) ) * idel + elsewhere + dvar(il:ir,:,:) = ( a2(-1) * var(il-1:ir-1,:,:) + & + a2(0) * var(il:ir,:,:) + & + a2(1) * var(il+1:ir+1,:,:) ) * idel + end where case (1) direction if ( zero_derivs_y /= 0 ) then dvar = zero @@ -212,38 +215,50 @@ subroutine up_deriv_gf_2_1 ( var, ni, nj, nk, dir, bb, gsize, delta, up, dvar ) if ( bb(1) == 0 ) then jl = 1 + gsize else - upwinding_yl: select case (up) - case (-1) upwinding_yl - q = q1 - case (1) upwinding_yl - q = q2 - end select upwinding_yl - - dvar(:,1,:) = ( q(1,1) * var(:,1,:) + q(2,1) * var(:,2,:) ) * idel - dvar(:,2,:) = ( q(1,2) * var(:,1,:) + q(2,2) * var(:,2,:) + & - q(3,2) * var(:,3,:) ) * idel + where ( up(:,1,:) < zero ) + dvar(:,1,:) = ( q1(1,1) * var(:,1,:) + q1(2,1) * var(:,2,:) ) * idel + elsewhere + dvar(:,1,:) = ( q2(1,1) * var(:,1,:) + q2(2,1) * var(:,2,:) ) * idel + end where + where ( up(:,2,:) < zero ) + dvar(:,2,:) = ( q1(1,2) * var(:,1,:) + q1(2,2) * var(:,2,:) + & + q1(3,2) * var(:,3,:) ) * idel + elsewhere + dvar(:,2,:) = ( q2(1,2) * var(:,1,:) + q2(2,2) * var(:,2,:) + & + q2(3,2) * var(:,3,:) ) * idel + end where jl = 3 end if if ( bb(2) == 0 ) then jr = nj - gsize else - upwinding_yr: select case (up) - case (-1) upwinding_yr - q = q2 - case (1) upwinding_yr - q = q1 - end select upwinding_yr - - dvar(:,nj-1,:) = - ( q(1,2) * var(:,nj,:) + & - q(2,2) * var(:,nj-1,:) + & - q(3,2) * var(:,nj-2,:) ) * idel - dvar(:,nj,:) = - ( q(1,1) * var(:,nj,:) + & - q(2,1) * var(:,nj-1,:) ) * idel + where ( up(:,nj-1,:) < zero ) + dvar(:,nj-1,:) = - ( q2(1,2) * var(:,nj,:) + & + q2(2,2) * var(:,nj-1,:) + & + q2(3,2) * var(:,nj-2,:) ) * idel + elsewhere + dvar(:,nj-1,:) = - ( q1(1,2) * var(:,nj,:) + & + q1(2,2) * var(:,nj-1,:) + & + q1(3,2) * var(:,nj-2,:) ) * idel + end where + where ( up(:,nj,:) < zero ) + dvar(:,nj,:) = - ( q2(1,1) * var(:,nj,:) + & + q2(2,1) * var(:,nj-1,:) ) * idel + elsewhere + dvar(:,nj,:) = - ( q1(1,1) * var(:,nj,:) + & + q1(2,1) * var(:,nj-1,:) ) * idel + end where jr = nj - 2 end if - dvar(:,jl:jr,:) = ( a(-1) * var(:,jl-1:jr-1,:) + & - a(0) * var(:,jl:jr,:) + & - a(1) * var(:,jl+1:jr+1,:) ) * idel + where ( up(:,jl:jr,:) < zero ) + dvar(:,jl:jr,:) = ( a1(-1) * var(:,jl-1:jr-1,:) + & + a1(0) * var(:,jl:jr,:) + & + a1(1) * var(:,jl+1:jr+1,:) ) * idel + elsewhere + dvar(:,jl:jr,:) = ( a2(-1) * var(:,jl-1:jr-1,:) + & + a2(0) * var(:,jl:jr,:) + & + a2(1) * var(:,jl+1:jr+1,:) ) * idel + end where end if case (2) direction if ( zero_derivs_z /= 0 ) then @@ -252,38 +267,50 @@ subroutine up_deriv_gf_2_1 ( var, ni, nj, nk, dir, bb, gsize, delta, up, dvar ) if ( bb(1) == 0 ) then kl = 1 + gsize else - upwinding_zl: select case (up) - case (-1) upwinding_zl - q = q1 - case (1) upwinding_zl - q = q2 - end select upwinding_zl - - dvar(:,:,1) = ( q(1,1) * var(:,:,1) + q(2,1) * var(:,:,2) ) * idel - dvar(:,:,2) = ( q(1,2) * var(:,:,1) + q(2,2) * var(:,:,2) + & - q(3,2) * var(:,:,3) ) * idel + where ( up(:,:,1) < zero ) + dvar(:,:,1) = ( q1(1,1) * var(:,:,1) + q1(2,1) * var(:,:,2) ) * idel + elsewhere + dvar(:,:,1) = ( q2(1,1) * var(:,:,1) + q2(2,1) * var(:,:,2) ) * idel + end where + where ( up(:,:,2) < zero ) + dvar(:,:,2) = ( q1(1,2) * var(:,:,1) + q1(2,2) * var(:,:,2) + & + q1(3,2) * var(:,:,3) ) * idel + elsewhere + dvar(:,:,2) = ( q2(1,2) * var(:,:,1) + q2(2,2) * var(:,:,2) + & + q2(3,2) * var(:,:,3) ) * idel + end where kl = 3 end if if ( bb(2) == 0 ) then kr = nk - gsize else - upwinding_zr: select case (up) - case (-1) upwinding_zr - q = q2 - case (1) upwinding_zr - q = q1 - end select upwinding_zr - - dvar(:,:,nk-1) = - ( q(1,2) * var(:,:,nk) + & - q(2,2) * var(:,:,nk-1) + & - q(3,2) * var(:,:,nk-2) ) * idel - dvar(:,:,nk) = - ( q(1,1) * var(:,:,nk) + & - q(2,1) * var(:,:,nk-1) ) * idel + where ( up(:,:,nk-1) < zero ) + dvar(:,:,nk-1) = - ( q2(1,2) * var(:,:,nk) + & + q2(2,2) * var(:,:,nk-1) + & + q2(3,2) * var(:,:,nk-2) ) * idel + elsewhere + dvar(:,:,nk-1) = - ( q1(1,2) * var(:,:,nk) + & + q1(2,2) * var(:,:,nk-1) + & + q1(3,2) * var(:,:,nk-2) ) * idel + end where + where ( up(:,:,nk) < zero ) + dvar(:,:,nk) = - ( q2(1,1) * var(:,:,nk) + & + q2(2,1) * var(:,:,nk-1) ) * idel + elsewhere + dvar(:,:,nk) = - ( q1(1,1) * var(:,:,nk) + & + q1(2,1) * var(:,:,nk-1) ) * idel + end where kr = nk - 2 end if - dvar(:,:,kl:kr) = ( a(-1) * var(:,:,kl-1:kr-1) + & - a(0) * var(:,:,kl:kr) + & - a(1) * var(:,:,kl+1:kr+1) ) * idel + where ( up(:,:,kl:kr) < zero ) + dvar(:,:,kl:kr) = ( a1(-1) * var(:,:,kl-1:kr-1) + & + a1(0) * var(:,:,kl:kr) + & + a1(1) * var(:,:,kl+1:kr+1) ) * idel + elsewhere + dvar(:,:,kl:kr) = ( a2(-1) * var(:,:,kl-1:kr-1) + & + a2(0) * var(:,:,kl:kr) + & + a2(1) * var(:,:,kl+1:kr+1) ) * idel + end where end if end select direction end subroutine up_deriv_gf_2_1 |