diff options
Diffstat (limited to 'src/Dissipation_4_3_min_err_coeff.F90')
-rw-r--r-- | src/Dissipation_4_3_min_err_coeff.F90 | 25 |
1 files changed, 19 insertions, 6 deletions
diff --git a/src/Dissipation_4_3_min_err_coeff.F90 b/src/Dissipation_4_3_min_err_coeff.F90 index b4b6ba3..dc4178d 100644 --- a/src/Dissipation_4_3_min_err_coeff.F90 +++ b/src/Dissipation_4_3_min_err_coeff.F90 @@ -4,7 +4,8 @@ #include "cctk_Parameters.h" -subroutine dissipation_4_3_opt (var, lsh, gsh, lbnd, bb, gsize, delta, epsdisl, dfl, rhs) +subroutine dissipation_4_3_opt (var, lsh, gsh, lbnd, bb, gsize, & + delta, epsilon, dfl, rhs) implicit none @@ -17,7 +18,7 @@ subroutine dissipation_4_3_opt (var, lsh, gsh, lbnd, bb, gsize, delta, epsdisl, CCTK_INT, dimension(6), intent(in) :: bb CCTK_INT, dimension(3), intent(in) :: gsize CCTK_REAL, dimension(3), intent(in) :: delta - CCTK_REAL, intent(in) :: epsdisl, dfl + CCTK_REAL, intent(in) :: epsilon, dfl CCTK_REAL :: zero = 0.0 integer, parameter :: wp = kind(zero) @@ -28,6 +29,7 @@ subroutine dissipation_4_3_opt (var, lsh, gsh, lbnd, bb, gsize, delta, epsdisl, CCTK_INT :: il, ir, jl, jr, kl, kr ni = lsh(1); nj = lsh(2); nk = lsh(3) + allocate ( a(ni,ni), d(ni,ni), b(ni,ni), h(ni,ni) ) a = zero; d = zero; b = zero; h = zero; @@ -39,7 +41,11 @@ subroutine dissipation_4_3_opt (var, lsh, gsh, lbnd, bb, gsize, delta, epsdisl, a = - transpose ( matmul ( h, matmul ( transpose(d), matmul ( b, d ) ) ) ) - idel = epsdisl / 16.0_wp + if ( scale_with_h > 0 ) then + idel = epsilon / ( 16 * delta(1) ) + else + idel = epsilon / 16 + end if if ( bb(1) == 0 ) then il = 1 + gsize(1) @@ -132,7 +138,11 @@ subroutine dissipation_4_3_opt (var, lsh, gsh, lbnd, bb, gsize, delta, epsdisl, a = - transpose ( matmul ( h, matmul ( transpose(d), matmul ( b, d ) ) ) ) - idel = epsdisl / 16.0_wp + if ( scale_with_h > 0 ) then + idel = epsilon / ( 16 * delta(2) ) + else + idel = epsilon / 16 + end if if ( bb(3) == 0 ) then jl = 1 + gsize(2) @@ -226,8 +236,11 @@ subroutine dissipation_4_3_opt (var, lsh, gsh, lbnd, bb, gsize, delta, epsdisl, a = - transpose ( matmul ( h, matmul ( transpose(d), matmul ( b, d ) ) ) ) - idel = epsdisl / 16.0_wp - + if ( scale_with_h > 0 ) then + idel = epsilon / ( 16 * delta(3) ) + else + idel = epsilon / 16 + end if if ( bb(5) == 0 ) then kl = 1 + gsize(3) |