diff options
Diffstat (limited to 'src/GRHydro_Con2PrimAM.F90')
-rw-r--r-- | src/GRHydro_Con2PrimAM.F90 | 41 |
1 files changed, 35 insertions, 6 deletions
diff --git a/src/GRHydro_Con2PrimAM.F90 b/src/GRHydro_Con2PrimAM.F90 index 6b8d894..7d2103e 100644 --- a/src/GRHydro_Con2PrimAM.F90 +++ b/src/GRHydro_Con2PrimAM.F90 @@ -13,7 +13,6 @@ #include "cctk_Arguments.h" #include "cctk_Functions.h" #include "SpaceMask.h" -#include "GRHydro_InterfacesAM.h" #include "GRHydro_Macros.h" #define ITER_TOL (1.0e-8) @@ -349,6 +348,18 @@ subroutine Conservative2PrimitiveAM(CCTK_ARGUMENTS) !of the comoving B-field, b^{\mu} b_{\mu}. It is overwritten !in this routine, but we may need to find a better notation !avoid future confusions. + if(GRHydro_eos_handle .eq. 1 .or. GRHydro_eos_handle .eq. 2) then + + call GRHydro_Con2PrimM_ptold(GRHydro_eos_handle, local_gam(1), dens(i,j,k), & + scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), tau(i,j,k), & + Bconsx_tmp,Bconsy_tmp,Bconsz_tmp,rho_tmp, & + velx_tmp,vely_tmp,velz_tmp,eps_tmp,press_tmp, & + Bvecx_tmp,Bvecy_tmp,Bvecz_tmp,b2, w_lorentz_tmp,& + g11(i,j,k),g12(i,j,k),g13(i,j,k),g22(i,j,k),g23(i,j,k),g33(i,j,k), & + uxx,uxy,uxz,uyy,uyz,uzz,det, & + epsnegative,GRHydro_C2P_failed(i,j,k)) + + else call GRHydro_Con2PrimM_pt2(GRHydro_eos_handle, keytemp, GRHydro_eos_rf_prec, GRHydro_perc_ptol, local_gam(1), dens(i,j,k), & scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), tau(i,j,k), & Bconsx_tmp, Bconsy_tmp, Bconsz_tmp,xye(1), & @@ -358,9 +369,16 @@ subroutine Conservative2PrimitiveAM(CCTK_ARGUMENTS) g22(i,j,k),g23(i,j,k),g33(i,j,k), & uxx,uxy,uxz,uyy,uyz,uzz,det, & epsnegative,GRHydro_C2P_failed(i,j,k)) + endif if(evolve_entropy.ne.0) then + if(GRHydro_C2P_failed(i,j,k).ne.0) then + !Use previous time step for rho: entropy(i,j,k) = entropycons(i,j,k)/dens(i,j,k)*rho(i,j,k) + else + !Use the current correct value of rho returned by con2prim: + entropy(i,j,k) = entropycons(i,j,k)/dens(i,j,k)*rho_tmp + endif endif if(GRHydro_C2P_failed(i,j,k).ne.0) then @@ -449,6 +467,19 @@ subroutine Conservative2PrimitiveAM(CCTK_ARGUMENTS) Bconsy_tmp = sdet*Bvecy_tmp Bconsz_tmp = sdet*Bvecz_tmp + if(evolve_entropy.ne.0) then + call GRHydro_Con2PrimM_ptee(GRHydro_eos_handle, keytemp, & + GRHydro_eos_rf_prec, local_gam(1), dens(i,j,k), & + scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), tau(i,j,k), & + Bconsx_tmp,Bconsy_tmp,Bconsz_tmp, & + entropycons(i,j,k), xye(1), & + xtemp(1),rho_tmp,velx_tmp,vely_tmp,velz_tmp,& + eps_tmp,press_tmp,Bvecx_tmp,Bvecy_tmp,Bvecz_tmp,b2,& + w_lorentz_tmp,g11(i,j,k),g12(i,j,k),g13(i,j,k),& + g22(i,j,k),g23(i,j,k),g33(i,j,k), & + uxx,uxy,uxz,uyy,uyz,uzz,det, & + epsnegative,GRHydro_C2P_failed(i,j,k)) + else call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, local_pgam, & dens(i,j,k),scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), sc, & Bconsx_tmp, Bconsy_tmp, Bconsz_tmp, rho_tmp,& @@ -457,6 +488,7 @@ subroutine Conservative2PrimitiveAM(CCTK_ARGUMENTS) g11(i,j,k),g12(i,j,k),g13(i,j,k),g22(i,j,k),g23(i,j,k),g33(i,j,k), & uxx,uxy,uxz,uyy,uyz,uzz,det, & epsnegative,GRHydro_C2P_failed(i,j,k)) + end if rho(i,j,k) = rho_tmp press(i,j,k) = press_tmp @@ -682,11 +714,7 @@ subroutine Conservative2PrimitiveAM(CCTK_ARGUMENTS) !!$ tau does need to take into account the existing B-field !!$ with w_lorentz=1, we find tau = sqrtdet*(rho (1+eps+b^2/2)) - dens [Press drops out] - tau(i,j,k) = sdet * (rho(i,j,k)*(1.0+eps(i,j,k)+b2/2.0)) - dens(i,j,k) - - if(tau(i,j,k).le.sdet*b2*0.5d0)then - tau(i,j,k) = GRHydro_tau_min + sdet*b2*0.5d0 - endif + tau(i,j,k) = sdet * (rho(i,j,k)*eps(i,j,k)+b2/2.0) end if @@ -870,6 +898,7 @@ subroutine Conservative2PrimitiveBoundsAM(CCTK_ARGUMENTS) if(evolve_Y_e.ne.0) then Y_e(i,j,k) = Y_e_con(i,j,k) / dens(i,j,k) endif + Bconsx_tmp = sqrt(avg_detl)*Bvecxminus(i,j,k) Bconsy_tmp = sqrt(avg_detl)*Bvecyminus(i,j,k) Bconsz_tmp = sqrt(avg_detl)*Bveczminus(i,j,k) |