aboutsummaryrefslogtreecommitdiff
path: root/src/GRHydro_HLLEM.F90
diff options
context:
space:
mode:
Diffstat (limited to 'src/GRHydro_HLLEM.F90')
-rw-r--r--src/GRHydro_HLLEM.F9073
1 files changed, 37 insertions, 36 deletions
diff --git a/src/GRHydro_HLLEM.F90 b/src/GRHydro_HLLEM.F90
index 6f1791e..995dd67 100644
--- a/src/GRHydro_HLLEM.F90
+++ b/src/GRHydro_HLLEM.F90
@@ -226,10 +226,10 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS)
vxtp,vytp,vztp,pressstarp,bxlowp,bylowp,bzlowp,ab0p,wp, &
avg_det,avg_alp,avg_beta)
if(clean_divergence.ne.0) then
- f1(6)=f1(6) + avg_alp*sdet*uxxh*psidcp - cons_p(6)*avg_betax
- f1(7)=f1(7) + avg_alp*sdet*uxyh*psidcp - cons_p(6)*avg_betay
- f1(8)=f1(8) + avg_alp*sdet*uxzh*psidcp - cons_p(6)*avg_betaz
- psidcf = avg_alp*cons_p(6)/sdet-psidcp*avg_betax
+ f1(6)=f1(6) + 1.0*sdet*uxxh*psidcp - cons_p(6)*avg_betax/avg_alp
+ f1(7)=f1(7) + 1.0*sdet*uxyh*psidcp - cons_p(6)*avg_betay/avg_alp
+ f1(8)=f1(8) + 1.0*sdet*uxzh*psidcp - cons_p(6)*avg_betaz/avg_alp
+ psidcf = cons_p(6)/sdet-psidcp*avg_betax/avg_alp
endif
else if (flux_direction == 2) then
@@ -239,10 +239,10 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS)
vytp,vztp,vxtp,pressstarp,bylowp,bzlowp,bxlowp,ab0p,wp, &
avg_det,avg_alp,avg_beta)
if(clean_divergence.ne.0) then
- f1(6)=f1(6) + avg_alp*sdet*uxyh*psidcp - cons_p(7)*avg_betax
- f1(7)=f1(7) + avg_alp*sdet*uyyh*psidcp - cons_p(7)*avg_betay
- f1(8)=f1(8) + avg_alp*sdet*uyzh*psidcp - cons_p(7)*avg_betaz
- psidcf = avg_alp*cons_p(7)/sdet-psidcp*avg_betay
+ f1(6)=f1(6) + 1.0*sdet*uxyh*psidcp - cons_p(7)*avg_betax/avg_alp
+ f1(7)=f1(7) + 1.0*sdet*uyyh*psidcp - cons_p(7)*avg_betay/avg_alp
+ f1(8)=f1(8) + 1.0*sdet*uyzh*psidcp - cons_p(7)*avg_betaz/avg_alp
+ psidcf = cons_p(7)/sdet-psidcp*avg_betay/avg_alp
endif
else if (flux_direction == 3) then
@@ -252,10 +252,10 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS)
vztp,vxtp,vytp,pressstarp,bzlowp,bxlowp,bylowp,ab0p,wp, &
avg_det,avg_alp,avg_beta)
if(clean_divergence.ne.0) then
- f1(6)=f1(6) + avg_alp*sdet*uxzh*psidcp - cons_p(8)*avg_betax
- f1(7)=f1(7) + avg_alp*sdet*uyzh*psidcp - cons_p(8)*avg_betay
- f1(8)=f1(8) + avg_alp*sdet*uzzh*psidcp - cons_p(8)*avg_betaz
- psidcf = avg_alp*cons_p(8)/sdet-psidcp*avg_betaz
+ f1(6)=f1(6) + 1.0*sdet*uxzh*psidcp - cons_p(8)*avg_betax/avg_alp
+ f1(7)=f1(7) + 1.0*sdet*uyzh*psidcp - cons_p(8)*avg_betay/avg_alp
+ f1(8)=f1(8) + 1.0*sdet*uzzh*psidcp - cons_p(8)*avg_betaz/avg_alp
+ psidcf = cons_p(8)/sdet-psidcp*avg_betaz/avg_alp
endif
else
@@ -315,14 +315,14 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS)
avg_det,avg_alp,avg_beta)
if(clean_divergence.ne.0) then
- fminus(6)=fminus(6) + avg_alp*sdet*uxxh*psidcm - cons_m(6)*avg_betax
- fminus(7)=fminus(7) + avg_alp*sdet*uxyh*psidcm - cons_m(6)*avg_betay
- fminus(8)=fminus(8) + avg_alp*sdet*uxzh*psidcm - cons_m(6)*avg_betaz
- fplus(6)=fplus(6) + avg_alp*sdet*uxxh*psidcp - cons_p(6)*avg_betax
- fplus(7)=fplus(7) + avg_alp*sdet*uxyh*psidcp - cons_p(6)*avg_betay
- fplus(8)=fplus(8) + avg_alp*sdet*uxzh*psidcp - cons_p(6)*avg_betaz
- psidcfp = avg_alp*cons_p(6)/sdet-avg_betax*psidcp
- psidcfm = avg_alp*cons_m(6)/sdet-avg_betax*psidcm
+ fminus(6)=fminus(6) + 1.0*sdet*uxxh*psidcm - cons_m(6)*avg_betax/avg_alp
+ fminus(7)=fminus(7) + 1.0*sdet*uxyh*psidcm - cons_m(6)*avg_betay/avg_alp
+ fminus(8)=fminus(8) + 1.0*sdet*uxzh*psidcm - cons_m(6)*avg_betaz/avg_alp
+ fplus(6)=fplus(6) + 1.0*sdet*uxxh*psidcp - cons_p(6)*avg_betax/avg_alp
+ fplus(7)=fplus(7) + 1.0*sdet*uxyh*psidcp - cons_p(6)*avg_betay/avg_alp
+ fplus(8)=fplus(8) + 1.0*sdet*uxzh*psidcp - cons_p(6)*avg_betaz/avg_alp
+ psidcfp = cons_p(6)/sdet-avg_betax*psidcp/avg_alp
+ psidcfm = cons_m(6)/sdet-avg_betax*psidcm/avg_alp
endif
else if (flux_direction == 2) then
@@ -349,14 +349,14 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS)
avg_det,avg_alp,avg_beta)
if(clean_divergence.ne.0) then
- fminus(6)=fminus(6) + avg_alp*sdet*uxyh*psidcm - cons_m(7)*avg_betax
- fminus(7)=fminus(7) + avg_alp*sdet*uyyh*psidcm - cons_m(7)*avg_betay
- fminus(8)=fminus(8) + avg_alp*sdet*uyzh*psidcm - cons_m(7)*avg_betaz
- fplus(6)=fplus(6) + avg_alp*sdet*uxyh*psidcp - cons_p(7)*avg_betax
- fplus(7)=fplus(7) + avg_alp*sdet*uyyh*psidcp - cons_p(7)*avg_betay
- fplus(8)=fplus(8) + avg_alp*sdet*uyzh*psidcp - cons_p(7)*avg_betaz
- psidcfp = avg_alp*cons_p(7)/sdet-avg_betay*psidcp
- psidcfm = avg_alp*cons_m(7)/sdet-avg_betay*psidcm
+ fminus(6)=fminus(6) + 1.0*sdet*uxyh*psidcm - cons_m(7)*avg_betax/avg_alp
+ fminus(7)=fminus(7) + 1.0*sdet*uyyh*psidcm - cons_m(7)*avg_betay/avg_alp
+ fminus(8)=fminus(8) + 1.0*sdet*uyzh*psidcm - cons_m(7)*avg_betaz/avg_alp
+ fplus(6)=fplus(6) + 1.0*sdet*uxyh*psidcp - cons_p(7)*avg_betax/avg_alp
+ fplus(7)=fplus(7) + 1.0*sdet*uyyh*psidcp - cons_p(7)*avg_betay/avg_alp
+ fplus(8)=fplus(8) + 1.0*sdet*uyzh*psidcp - cons_p(7)*avg_betaz/avg_alp
+ psidcfp = cons_p(7)/sdet-avg_betay*psidcp/avg_alp
+ psidcfm = cons_m(7)/sdet-avg_betay*psidcm/avg_alp
endif
else if (flux_direction == 3) then
@@ -383,14 +383,14 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS)
avg_det,avg_alp,avg_beta)
if(clean_divergence.ne.0) then
- fminus(6)=fminus(6) + avg_alp*sdet*uxzh*psidcm - cons_m(8)*avg_betax
- fminus(7)=fminus(7) + avg_alp*sdet*uyzh*psidcm - cons_m(8)*avg_betay
- fminus(8)=fminus(8) + avg_alp*sdet*uzzh*psidcm - cons_m(8)*avg_betaz
- fplus(6)=fplus(6) + avg_alp*sdet*uxzh*psidcp - cons_p(8)*avg_betax
- fplus(7)=fplus(7) + avg_alp*sdet*uyzh*psidcp - cons_p(8)*avg_betay
- fplus(8)=fplus(8) + avg_alp*sdet*uzzh*psidcp - cons_p(8)*avg_betaz
- psidcfp = avg_alp*cons_p(8)/sdet-avg_betaz*psidcp
- psidcfm = avg_alp*cons_m(8)/sdet-avg_betaz*psidcm
+ fminus(6)=fminus(6) + 1.0*sdet*uxzh*psidcm - cons_m(8)*avg_betax/avg_alp
+ fminus(7)=fminus(7) + 1.0*sdet*uyzh*psidcm - cons_m(8)*avg_betay/avg_alp
+ fminus(8)=fminus(8) + 1.0*sdet*uzzh*psidcm - cons_m(8)*avg_betaz/avg_alp
+ fplus(6)=fplus(6) + 1.0*sdet*uxzh*psidcp - cons_p(8)*avg_betax/avg_alp
+ fplus(7)=fplus(7) + 1.0*sdet*uyzh*psidcp - cons_p(8)*avg_betay/avg_alp
+ fplus(8)=fplus(8) + 1.0*sdet*uzzh*psidcp - cons_p(8)*avg_betaz/avg_alp
+ psidcfp = cons_p(8)/sdet-avg_betaz*psidcp/avg_alp
+ psidcfm = cons_m(8)/sdet-avg_betaz*psidcm/avg_alp
endif
else
@@ -428,6 +428,7 @@ subroutine GRHydro_HLLEM(CCTK_ARGUMENTS)
end do
if(clean_divergence.ne.0) then
+
psidcdiff = psidcm - psidcp
select case(whichpsidcspeed)
case(0)
generated by cgit v1.2.3 (git 2.39.1) at 2024-07-06 02:39:57 +0000