diff options
| author | rhaas <rhaas@c83d129a-5a75-4d5a-9c4d-ed3a5855bf45> | 2012-06-20 16:17:32 +0000 |
|---|---|---|
| committer | rhaas <rhaas@c83d129a-5a75-4d5a-9c4d-ed3a5855bf45> | 2012-06-20 16:17:32 +0000 |
| commit | eab01365bc0d62d8aacaa1a45987b9add116f98d (patch) | |
| tree | 0e5eae235fbd88438df1fb3dfc055b05dd28b1dd /src/GRHydro_HLLEM.F90 | |
| parent | d0f862c149c4f6dfab4bed30b5345d553dcd4001 (diff) | |
Change the characteristic speed for psidc.
It was set to +-c when it should be the maximum coming from the constraint
that the normal vector to the characteristic hypersurface of a wave
propagating along an arbitrary x direction be spacelike. This results
into the following bound for the wavespeed:
+- \alpha sqrt(\gamma^{xx}) - \beta^x)
Patch by Bruno Mundim.
git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinEvolve/GRHydro/trunk@358 c83d129a-5a75-4d5a-9c4d-ed3a5855bf45
Diffstat (limited to 'src/GRHydro_HLLEM.F90')
| -rw-r--r-- | src/GRHydro_HLLEM.F90 | 25 |
1 files changed, 20 insertions, 5 deletions
diff --git a/src/GRHydro_HLLEM.F90 b/src/GRHydro_HLLEM.F90 index 1e274eb..3c694ca 100644 --- a/src/GRHydro_HLLEM.F90 +++ b/src/GRHydro_HLLEM.F90 @@ -55,6 +55,7 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS) CCTK_REAL :: pressstarp,pressstarm,velxlowp,velxlowm,velylowp,velylowm,velzlowp,velzlowm CCTK_REAL :: psidcp, psidcm, psidcf, psidcdiff, psidcfp, psidcfm + CCTK_REAL :: charmax_dc, charmin_dc, charpm_dc CCTK_INT :: type_bits, trivial, not_trivial @@ -440,13 +441,27 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS) end if case(1) !!$ Wavespeeds for psidc are +/-c, not Fast Magnetosonic? - psidcf = 0.5d0 * (1.d0 * psidcfp - (-1.d0) * psidcfm + & - 1.d0 * (-1.d0) * psidcdiff) + !!$ psidcf = 0.5d0 * (1.d0 * psidcfp - (-1.d0) * psidcfm + & + !!$ 1.d0 * (-1.d0) * psidcdiff) + + !!$ The fastest speed for psidc comes from the condition + !!$ that the normal vector to the characteristic hypersurface + !!$ be spacelike (Eq. 60 of Anton et al.) + + charmax_dc = sqrt(usendh) - avg_beta/avg_alp + charmin_dc = -1.d0*sqrt(usendh) - avg_beta/avg_alp + charpm_dc = charmax_dc - charmin_dc + + psidcf = (charmax_dc * psidcfp - charmin_dc * psidcfm + & + charmax_dc * charmin_dc * psidcdiff) / charpm_dc + if(decouple_normal_Bfield .ne. 0) then ! same expression for HLLE and LLF !!$ B^i field decouples from the others and has same propagation - !!$ speed as divergence - f1(5+flux_direction) = 0.5d0 * (1.d0 * fplus(5+flux_direction) - (-1.d0) * fminus(5+flux_direction) + & - 1.d0 * (-1.d0) * qdiff(5+flux_direction)) + !!$ speed as divergence -null direction, + !!$ \pm sqrt(g^{xx}} - beta^x/alpha + f1(5+flux_direction) = (charmax_dc * fplus(5+flux_direction) & + - charmin_dc * fminus(5+flux_direction) + & + charmax_dc * charmin_dc * qdiff(5+flux_direction)) / charpm_dc end if case(2) charmax = setcharmax |