aboutsummaryrefslogtreecommitdiff
path: root/src/GRHydro_Con2PrimM.F90
diff options
context:
space:
mode:
authorbmundim <bmundim@c83d129a-5a75-4d5a-9c4d-ed3a5855bf45>2010-12-25 06:53:00 +0000
committerbmundim <bmundim@c83d129a-5a75-4d5a-9c4d-ed3a5855bf45>2010-12-25 06:53:00 +0000
commit4c09e9bedce52b7c867b47dbe83b23d6802b344c (patch)
tree5e15fa210d276ce1c83df8519fd28bcb39539f34 /src/GRHydro_Con2PrimM.F90
parentb0598b903ba8f3bb5ed89c29e0807cf4cb108110 (diff)
RIT MHD dev:
Fix to GRHydro_RegisterVarsM.cc and allows MHD to work without Y_E bing evolved. Add C2P polytype interface. Change schedule.ccl to represent a more detailed version of choosing what gridfucntions get sync'd in Boundaries depending on what is being evolved. Un-hardwire fixed values of gamma for all MHD Con2Prim routines. git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinEvolve/GRHydro/trunk@199 c83d129a-5a75-4d5a-9c4d-ed3a5855bf45
Diffstat (limited to 'src/GRHydro_Con2PrimM.F90')
-rw-r--r--src/GRHydro_Con2PrimM.F90254
1 files changed, 184 insertions, 70 deletions
diff --git a/src/GRHydro_Con2PrimM.F90 b/src/GRHydro_Con2PrimM.F90
index bbb803d..984569a 100644
--- a/src/GRHydro_Con2PrimM.F90
+++ b/src/GRHydro_Con2PrimM.F90
@@ -58,11 +58,11 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS)
CCTK_INT :: type_bits, atmosphere
CCTK_INT :: type2_bits
- CCTK_REAL :: local_min_tracer
+ CCTK_REAL :: local_min_tracer, local_gam, local_pgam,local_K,sc
! begin EOS Omni vars
integer :: n,keytemp,anyerr,keyerr(1)
- real*8 :: xpress,xtemp,xye,xeps
+ real*8 :: xpress,xtemp,xye,xeps,xrho
n=1;keytemp=0;anyerr=0;keyerr(1)=0
xpress=0.0d0;xtemp=0.0d0;xye=0.0d0;xeps=0.0d0
! end EOS Omni vars
@@ -87,6 +87,11 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS)
call EOS_Omni_EpsFromPress(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
GRHydro_rho_min,xeps,xtemp,xye,pmin,epsmin,keyerr,anyerr)
+ call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ 1.d0,1.d0,xtemp,xye,local_gam,keyerr,anyerr)
+ local_gam = local_gam+1.d0
+
+
!$OMP PARALLEL DO PRIVATE(i,j,itracer,&
!$OMP uxx, uxy, uxz, uyy, uyz, uzz, det, epsnegative)
do k = 1, nz
@@ -154,7 +159,7 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS)
end if
if(evolve_temper.eq.0) then
- call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, dens(i,j,k), &
+ call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam, dens(i,j,k), &
scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), tau(i,j,k), &
rho(i,j,k),&
vel(i,j,k,1),vel(i,j,k,2),vel(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k),&
@@ -219,14 +224,24 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS)
call CCTK_WARN(GRHydro_NaN_verbose+2, 'Specific internal energy just went below 0, trying polytype.')
!$OMP END CRITICAL
-!!$ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, dens(i,j,k), &
-!!$ scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), tau(i,j,k), &
-!!$ rho(i,j,k),&
-!!$ vel(i,j,k,1),vel(i,j,k,2),vel(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k),&
-!!$ gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),gyz(i,j,k),gzz(i,j,k), &
-!!$ uxx,uxy,uxz,uyy,uyz,uzz,det, &
-!!$ Bvec(i,j,k,1), Bvec(i,j,k,2),Bvec(i,j,k,3),b2,&
-!!$ epsnegative,GRHydro_C2P_failed(i,j,k))
+ xrho=10.0d0
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ 1.d0,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_K = xpress
+
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ xrho,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_pgam=log(xpress/local_K)/log(xrho)
+ sc = local_K*dens(i,j,k)
+
+ 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, &
+ rho(i,j,k),&
+ vel(i,j,k,1),vel(i,j,k,2),vel(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k),&
+ gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),gyz(i,j,k),gzz(i,j,k), &
+ uxx,uxy,uxz,uyy,uyz,uzz,det, &
+ Bvec(i,j,k,1), Bvec(i,j,k,2),Bvec(i,j,k,3),b2,&
+ epsnegative,GRHydro_C2P_failed(i,j,k))
#endif
@@ -276,7 +291,7 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS)
uxxr, uxyr, uxzr, uyyr, uyzr, uzzr, pmin, epsmin
CCTK_REAL :: gxxl,gxyl,gxzl,gyyl,gyzl,gzzl,avg_detl,&
gxxr,gxyr,gxzr,gyyr,gyzr,gzzr,avg_detr
- CCTK_REAL :: b2minus, b2plus
+ CCTK_REAL :: b2minus, b2plus, local_gam, local_pgam,local_K,sc
CCTK_INT :: epsnegative
character(len=100) warnline
@@ -285,9 +300,9 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS)
CCTK_REAL :: local_min_tracer
-! begin EOS Omni vars
+! begin EOS Omni vars
integer :: n,keytemp,anyerr,keyerr(1)
- real*8 :: xpress,xtemp,xye,xeps
+ real*8 :: xpress,xtemp,xye,xeps,xrho
n=1;keytemp=0;anyerr=0;keyerr(1)=0
xpress=0.0d0;xeps=0.0d0;xtemp=0.0d0;xye=0.0d0
! end EOS Omni vars
@@ -298,7 +313,11 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS)
call EOS_Omni_EpsFromPress(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
GRHydro_rho_min,epsmin,xtemp,xye,pmin,epsmin,keyerr,anyerr)
-
+
+ call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ 1.d0,1.0d0,xtemp,xye,local_gam,keyerr,anyerr)
+ local_gam=local_gam+1.0
+
call SpaceMask_GetTypeBits(type_bits, "Hydro_Atmosphere")
call SpaceMask_GetStateBits(atmosphere, "Hydro_Atmosphere", "in_atmosphere")
type2_bits = -1
@@ -363,7 +382,7 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS)
endif
- call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, densminus(i,j,k), &
+ call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam,densminus(i,j,k), &
sxminus(i,j,k),syminus(i,j,k),szminus(i,j,k), tauminus(i,j,k), &
rhominus(i,j,k),&
velxminus(i,j,k),velyminus(i,j,k),velzminus(i,j,k),epsminus(i,j,k),pressminus(i,j,k),w_lorentzminus(i,j,k),&
@@ -376,14 +395,25 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS)
!$OMP CRITICAL
call CCTK_WARN(GRHydro_NaN_verbose+2, 'Specific internal energy just went below 0, trying polytype!')
!$OMP END CRITICAL
-!!$ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, densminus(i,j,k), &
-!!$ sxminus(i,j,k),syminus(i,j,k),szminus(i,j,k), tauminus(i,j,k), &
-!!$ rhominus(i,j,k),&
-!!$ velxminus(i,j,k),velyminus(i,j,k),velzminus(i,j,k),epsminus(i,j,k),pressminus(i,j,k),w_lorentzminus(i,j,k),&
-!!$ gxxl,gxyl,gxzl,gyyl,gyzl,gzzl, &
-!!$ uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl, &
-!!$ Bvecxminus(i,j,k), Bvecyminus(i,j,k),Bveczminus(i,j,k),b2minus,&
-!!$ epsnegative,GRHydro_C2P_failed(i,j,k))
+
+ xrho=10.0d0
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ 1.d0,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_K = xpress
+
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ xrho,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_pgam=log(xpress/local_K)/log(xrho)
+ sc = local_K*densminus(i,j,k)
+
+ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, local_pgam, densminus(i,j,k), &
+ sxminus(i,j,k),syminus(i,j,k),szminus(i,j,k), sc, &
+ rhominus(i,j,k),&
+ velxminus(i,j,k),velyminus(i,j,k),velzminus(i,j,k),epsminus(i,j,k),pressminus(i,j,k),w_lorentzminus(i,j,k),&
+ gxxl,gxyl,gxzl,gyyl,gyzl,gzzl, &
+ uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl, &
+ Bvecxminus(i,j,k), Bvecyminus(i,j,k),Bveczminus(i,j,k),b2minus,&
+ epsnegative,GRHydro_C2P_failed(i,j,k))
end if
if (epsminus(i,j,k) .lt. 0.0d0) then
@@ -408,7 +438,7 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS)
endif
epsnegative = 0
- call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, densplus(i,j,k), &
+ call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam, densplus(i,j,k), &
sxplus(i,j,k),syplus(i,j,k),szplus(i,j,k), tauplus(i,j,k),&
rhoplus(i,j,k),&
velxplus(i,j,k),velyplus(i,j,k),velzplus(i,j,k),epsplus(i,j,k),pressplus(i,j,k),w_lorentzplus(i,j,k),&
@@ -421,14 +451,25 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS)
!$OMP CRITICAL
call CCTK_WARN(GRHydro_NaN_verbose+2, 'Specific internal energy just went below 0, trying polytype!!')
!$OMP END CRITICAL
-!!$ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, densplus(i,j,k), &
-!!$ sxplus(i,j,k),syplus(i,j,k),szplus(i,j,k), tauplus(i,j,k),&
-!!$ rhoplus(i,j,k),&
-!!$ velxplus(i,j,k),velyplus(i,j,k),velzplus(i,j,k),epsplus(i,j,k),pressplus(i,j,k),w_lorentzplus(i,j,k),&
-!!$ gxxr,gxyr,gxzr,gyyr,gyzr,gzzr, &
-!!$ uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr, &
-!!$ Bvecxplus(i,j,k), Bvecyplus(i,j,k),Bveczplus(i,j,k),b2plus,&
-!!$ epsnegative,GRHydro_C2P_failed(i,j,k))
+
+ xrho=10.0d0
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ 1.d0,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_K = xpress
+
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ xrho,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_pgam=log(xpress/local_K)/log(xrho)
+ sc = local_K*densplus(i,j,k)
+
+ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, local_pgam, densplus(i,j,k), &
+ sxplus(i,j,k),syplus(i,j,k),szplus(i,j,k), sc,&
+ rhoplus(i,j,k),&
+ velxplus(i,j,k),velyplus(i,j,k),velzplus(i,j,k),epsplus(i,j,k),pressplus(i,j,k),w_lorentzplus(i,j,k),&
+ gxxr,gxyr,gxzr,gyyr,gyzr,gzzr, &
+ uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr, &
+ Bvecxplus(i,j,k), Bvecyplus(i,j,k),Bveczplus(i,j,k),b2plus,&
+ epsnegative,GRHydro_C2P_failed(i,j,k))
end if
if (epsplus(i,j,k) .lt. 0.0d0) then
@@ -492,9 +533,15 @@ subroutine Conservative2PrimitivePolytypeM(CCTK_ARGUMENTS)
CCTK_INT :: epsnegative
- CCTK_REAL :: local_min_tracer
+ CCTK_REAL :: local_min_tracer, local_gam, local_pgam,local_K, sc
! character(len=400) :: warnline
+! begin EOS Omni vars
+ integer :: n,keytemp,anyerr,keyerr(1)
+ real*8 :: xpress,xtemp,xye,xeps,xrho
+ n=1;keytemp=0;anyerr=0;keyerr(1)=0
+ xpress=0.0d0;xtemp=0.0d0;xye=0.0d0;xeps=0.0d0
+! end EOS Omni vars
call SpaceMask_GetTypeBits(type_bits, "Hydro_Atmosphere")
call SpaceMask_GetStateBits(atmosphere, "Hydro_Atmosphere", "in_atmosphere")
@@ -547,15 +594,24 @@ subroutine Conservative2PrimitivePolytypeM(CCTK_ARGUMENTS)
Y_e(i,j,k) = Y_e_con(i,j,k) / dens(i,j,k)
endif
-
-!!$ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, dens(i,j,k), &
-!!$ scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), tau(i,j,k),&
-!!$ rho(i,j,k),&
-!!$ vel(i,j,k,1),vel(i,j,k,2),vel(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k),&
-!!$ gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),gyz(i,j,k),gzz(i,j,k), &
-!!$ uxx,uxy,uxz,uyy,uyz,uzz,det, &
-!!$ Bvec(i,j,k,1), Bvec(i,j,k,2),Bvec(i,j,k,3),b2,&
-!!$ epsnegative,GRHydro_C2P_failed(i,j,k))
+ xrho=10.0d0
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ 1.d0,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_K = xpress
+
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ xrho,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_pgam=log(xpress/local_K)/log(xrho)
+ sc = local_K*dens(i,j,k)
+
+ 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,&
+ rho(i,j,k),&
+ vel(i,j,k,1),vel(i,j,k,2),vel(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k),&
+ gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),gyz(i,j,k),gzz(i,j,k), &
+ uxx,uxy,uxz,uyy,uyz,uzz,det, &
+ Bvec(i,j,k,1), Bvec(i,j,k,2),Bvec(i,j,k,3),b2,&
+ epsnegative,GRHydro_C2P_failed(i,j,k))
end do
end do
@@ -603,6 +659,15 @@ subroutine Con2PrimBoundsPolytypeM(CCTK_ARGUMENTS)
CCTK_INT :: type_bits, atmosphere
CCTK_INT :: type2_bits, epsnegative
+
+ CCTK_REAL :: local_gam, local_pgam,local_K,scplus,scminus
+
+! begin EOS Omni vars
+ integer :: n,keytemp,anyerr,keyerr(1)
+ real*8 :: xpress,xtemp,xye,xeps,xrho
+ n=1;keytemp=0;anyerr=0;keyerr(1)=0
+ xpress=0.0d0;xtemp=0.0d0;xye=0.0d0;xeps=0.0d0
+! end EOS Omni vars
call SpaceMask_GetTypeBits(type_bits, "Hydro_Atmosphere")
call SpaceMask_GetStateBits(atmosphere, "Hydro_Atmosphere", "in_atmosphere")
@@ -640,22 +705,34 @@ subroutine Con2PrimBoundsPolytypeM(CCTK_ARGUMENTS)
call UpperMetric(uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr,&
gxxr,gxyr,gxzr,gyyr,gyzr,gzzr)
-!!$ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, densminus(i,j,k), &
-!!$ sxminus(i,j,k),syminus(i,j,k),szminus(i,j,k), tauminus(i,j,k),&
-!!$ rhominus(i,j,k),&
-!!$ velxminus(i,j,k),velyminus(i,j,k),velzminus(i,j,k),epsminus(i,j,k),pressminus(i,j,k),w_lorentzminus(i,j,k),&
-!!$ gxxl,gxyl,gxzl,gyyl,gyzl,gzzl, &
-!!$ uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl, &
-!!$ Bvecxminus(i,j,k), Bvecyminus(i,j,k),Bveczminus(i,j,k),b2minus,&
-!!$ epsnegative,GRHydro_C2P_failed(i,j,k))
-!!$ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, densplus(i,j,k), &
-!!$ sxplus(i,j,k),syplus(i,j,k),szplus(i,j,k), tauplus(i,j,k),&
-!!$ rhoplus(i,j,k),&
-!!$ velxplus(i,j,k),velyplus(i,j,k),velzplus(i,j,k),epsplus(i,j,k),pressplus(i,j,k),w_lorentzplus(i,j,k),&
-!!$ gxxr,gxyr,gxzr,gyyr,gyzr,gzzr, &
-!!$ uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr, &
-!!$ Bvecxplus(i,j,k), Bvecyplus(i,j,k),Bveczplus(i,j,k),b2plus,&
-!!$ epsnegative,GRHydro_C2P_failed(i,j,k))
+ xrho=10.0d0
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ 1.d0,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_K = xpress
+
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ xrho,1.0d0,xtemp,xye,xpress,keyerr,anyerr)
+ local_pgam=log(xpress/local_K)/log(xrho)
+ scminus = local_K*densminus(i,j,k)
+ scplus = local_K*densplus(i,j,k)
+
+ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, local_pgam,densminus(i,j,k), &
+ sxminus(i,j,k),syminus(i,j,k),szminus(i,j,k), scminus,&
+ rhominus(i,j,k),&
+ velxminus(i,j,k),velyminus(i,j,k),velzminus(i,j,k),epsminus(i,j,k),pressminus(i,j,k),w_lorentzminus(i,j,k),&
+ gxxl,gxyl,gxzl,gyyl,gyzl,gzzl, &
+ uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl, &
+ Bvecxminus(i,j,k), Bvecyminus(i,j,k),Bveczminus(i,j,k),b2minus,&
+ epsnegative,GRHydro_C2P_failed(i,j,k))
+
+ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, local_pgam,densplus(i,j,k), &
+ sxplus(i,j,k),syplus(i,j,k),szplus(i,j,k), scplus,&
+ rhoplus(i,j,k),&
+ velxplus(i,j,k),velyplus(i,j,k),velzplus(i,j,k),epsplus(i,j,k),pressplus(i,j,k),w_lorentzplus(i,j,k),&
+ gxxr,gxyr,gxzr,gyyr,gyzr,gzzr, &
+ uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr, &
+ Bvecxplus(i,j,k), Bvecyplus(i,j,k),Bveczplus(i,j,k),b2plus,&
+ epsnegative,GRHydro_C2P_failed(i,j,k))
end do
end do
end do
@@ -682,6 +759,7 @@ end subroutine Con2PrimBoundsPolytypeM
subroutine Conservative2PrimitiveGeneralM(CCTK_ARGUMENTS)
USE GRHydro_Scalars
+ USE Con2PrimM_fortran_interfaces
implicit none
@@ -699,7 +777,7 @@ subroutine Conservative2PrimitiveGeneralM(CCTK_ARGUMENTS)
CCTK_INT :: count
CCTK_REAL :: b2,det
- CCTK_REAL :: uxx,uxy,uxz,uyy,uyz,uzz
+ CCTK_REAL :: uxx,uxy,uxz,uyy,uyz,uzz, local_gam, local_pgam
CCTK_INT :: type_bits, atmosphere
CCTK_INT :: type2_bits,epsnegative
@@ -707,6 +785,24 @@ subroutine Conservative2PrimitiveGeneralM(CCTK_ARGUMENTS)
CCTK_REAL :: local_min_tracer
CCTK_INT, dimension(3) :: loc, loc2
+
+!!$ begin EOS Omni vars
+ integer :: n,keytemp,anyerr,keyerr(1)
+ real*8 :: xpress(1),xtemp(1),xye(1),xdens(1),xeps(1)
+ n=1
+ keytemp=0
+ anyerr=0
+ keyerr(1)=0
+ xpress(1)=0.0d0
+ xtemp(1)=0.0d0
+ xye(1)=0.0d0
+!!$ end EOS Omni vars
+
+ xdens(1)=1.d0
+ xeps(1)=1.d0
+ call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ xdens,xeps,xtemp,xye,xpress,keyerr,anyerr)
+ local_gam=xpress(1)+1.d0
call SpaceMask_GetTypeBits(type_bits, "Hydro_Atmosphere")
call SpaceMask_GetStateBits(atmosphere, "Hydro_Atmosphere", "in_atmosphere")
@@ -776,7 +872,7 @@ subroutine Conservative2PrimitiveGeneralM(CCTK_ARGUMENTS)
gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),&
gyz(i,j,k),gzz(i,j,k))
- call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, dens(i,j,k), &
+ call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam,dens(i,j,k), &
scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), tau(i,j,k),&
rho(i,j,k),&
vel(i,j,k,1),vel(i,j,k,2),vel(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k),&
@@ -828,7 +924,14 @@ subroutine Con2PrimPolytypeGeneralM(CCTK_ARGUMENTS)
CCTK_INT :: atmosphere
CCTK_INT :: type2_bits,epsnegative
- CCTK_REAL :: local_min_tracer
+ CCTK_REAL :: local_min_tracer, local_pgam, local_k, sc
+
+! begin EOS Omni vars
+ integer :: n,keytemp,anyerr,keyerr(1)
+ real*8 :: xpress(1),xtemp(1),xye(1),xeps(1),xrho(1)
+ n=1;keytemp=0;anyerr=0;keyerr(1)=0
+ xpress(1)=0.0d0;xtemp(1)=0.0d0;xye(1)=0.0d0;xeps(1)=0.0d0
+! end EOS Omni vars
call SpaceMask_GetTypeBits(type_bits, "Hydro_Atmosphere")
call SpaceMask_GetStateBits(atmosphere, "Hydro_Atmosphere", "in_atmosphere")
@@ -894,15 +997,26 @@ subroutine Con2PrimPolytypeGeneralM(CCTK_ARGUMENTS)
call UpperMetric(uxx,uxy,uxz,uyy,uyz,uzz,det,&
gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),&
gyz(i,j,k),gzz(i,j,k))
+
+ xrho(1)=1.d0
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ xrho,xeps,xtemp,xye,xpress,keyerr,anyerr)
+ local_K = xpress(1)
+
+ xrho(1)=10.d0
+ call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,&
+ xrho,xeps,xtemp,xye,xpress,keyerr,anyerr)
+ local_pgam=log(xpress(1)/local_K)/log(xrho(1))
+ sc = local_K*dens(i,j,k)
-!!$ call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, dens(i,j,k), &
-!!$ scon(i,j,k,1),scon(i,j,k,2),scon(i,j,k,3), tau(i,j,k),&
-!!$ rho(i,j,k),&
-!!$ vel(i,j,k,1),vel(i,j,k,2),vel(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k),&
-!!$ gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),gyz(i,j,k),gzz(i,j,k), &
-!!$ uxx,uxy,uxz,uyy,uyz,uzz,det, &
-!!$ Bvec(i,j,k,1), Bvec(i,j,k,2),Bvec(i,j,k,3),b2,&
-!!$ epsnegative,GRHydro_C2P_failed(i,j,k))
+ 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,&
+ rho(i,j,k),&
+ vel(i,j,k,1),vel(i,j,k,2),vel(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k),&
+ gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),gyz(i,j,k),gzz(i,j,k), &
+ uxx,uxy,uxz,uyy,uyz,uzz,det, &
+ Bvec(i,j,k,1), Bvec(i,j,k,2),Bvec(i,j,k,3),b2,&
+ epsnegative,GRHydro_C2P_failed(i,j,k))
end do
end do