diff options
Diffstat (limited to 'src/GRHydro_SourceM.F90')
-rw-r--r-- | src/GRHydro_SourceM.F90 | 80 |
1 files changed, 29 insertions, 51 deletions
diff --git a/src/GRHydro_SourceM.F90 b/src/GRHydro_SourceM.F90 index ba2f383..5245df9 100644 --- a/src/GRHydro_SourceM.F90 +++ b/src/GRHydro_SourceM.F90 @@ -199,62 +199,40 @@ subroutine SourceTermsM(CCTK_ARGUMENTS) localgxy, localgxz, localgyy, localgyz, localgzz) - if (shift_state .ne. 0) then - - shiftx = betax(i,j,k) - shifty = betay(i,j,k) - shiftz = betaz(i,j,k) - - if (local_spatial_order .eq. 2) then - - dx_betax = DIFF_X_2(betax) - dx_betay = DIFF_X_2(betay) - dx_betaz = DIFF_X_2(betaz) - - dy_betax = DIFF_Y_2(betax) - dy_betay = DIFF_Y_2(betay) - dy_betaz = DIFF_Y_2(betaz) - - dz_betax = DIFF_Z_2(betax) - dz_betay = DIFF_Z_2(betay) - dz_betaz = DIFF_Z_2(betaz) - - else - - dx_betax = DIFF_X_4(betax) - dx_betay = DIFF_X_4(betay) - dx_betaz = DIFF_X_4(betaz) - - dy_betax = DIFF_Y_4(betax) - dy_betay = DIFF_Y_4(betay) - dy_betaz = DIFF_Y_4(betaz) + shiftx = betax(i,j,k) + shifty = betay(i,j,k) + shiftz = betaz(i,j,k) + + if (local_spatial_order .eq. 2) then + + dx_betax = DIFF_X_2(betax) + dx_betay = DIFF_X_2(betay) + dx_betaz = DIFF_X_2(betaz) + + dy_betax = DIFF_Y_2(betax) + dy_betay = DIFF_Y_2(betay) + dy_betaz = DIFF_Y_2(betaz) - dz_betax = DIFF_Z_4(betax) - dz_betay = DIFF_Z_4(betay) - dz_betaz = DIFF_Z_4(betaz) + dz_betax = DIFF_Z_2(betax) + dz_betay = DIFF_Z_2(betay) + dz_betaz = DIFF_Z_2(betaz) - end if - else - shiftx = 0.0d0 - shifty = 0.0d0 - shiftz = 0.0d0 - - dx_betax = 0.0d0 - dx_betay = 0.0d0 - dx_betaz = 0.0d0 - - dy_betax = 0.0d0 - dy_betay = 0.0d0 - dy_betaz = 0.0d0 - - dz_betax = 0.0d0 - dz_betay = 0.0d0 - dz_betaz = 0.0d0 + dx_betax = DIFF_X_4(betax) + dx_betay = DIFF_X_4(betay) + dx_betaz = DIFF_X_4(betaz) + + dy_betax = DIFF_Y_4(betax) + dy_betay = DIFF_Y_4(betay) + dy_betaz = DIFF_Y_4(betaz) + + dz_betax = DIFF_Z_4(betax) + dz_betay = DIFF_Z_4(betay) + dz_betaz = DIFF_Z_4(betaz) + + end if - endif - invalp = 1.0d0 / alp(i,j,k) invalp2 = invalp**2 velxshift = velx(i,j,k) - shiftx*invalp |