User: bmundim Date: 2010/12/25 12:44 AM

Added: /trunk/src/ GRHydro_P2C2PM_polytype.F90

Removed: /trunk/ interface.ccl.omni

Modified: /trunk/ param.ccl, schedule.ccl /trunk/src/ CheckParam.c, GRHydro_C2P2C.F90, GRHydro_C2P2CM.F90, GRHydro_Con2Prim.F90, GRHydro_P2C2P.F90, GRHydro_P2C2PM.F90, GRHydro_ReadConformalData.F90, GRHydro_ShockTube.F90, GRHydro_ShockTubeM.F90, GRHydro_SimpleWave.F90, make.code.defn

Log: RIT MHD dev:

Update initial data routines to properly call out pressures and other hydro quantities using EOS_Omni rather than now-defunct interfaces.

Add C2P_polytype interface.

Un-hardwire fixed values of gamma for all MHD Con2Prim routines.

File Changes:

Directory: /trunk/src/ ======================

File [modified]: CheckParam.c Delta lines: +10 -1 =================================================================== --- trunk/src/CheckParam.c 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/CheckParam.c 2010-12-25 06:44:22 UTC (rev 117) @@ -23,6 +23,15 @@ { CCTK_PARAMWARN("That test not yet implemented in MHD!"); } - + + if(!CCTK_Equals(Bvec_evolution_method,"GRHydro") && + ((CCTK_Equals(shock_case,"Balsara1")) || + (CCTK_Equals(shock_case,"Balsara2")) || + (CCTK_Equals(shock_case,"Balsara3")) || + (CCTK_Equals(shock_case,"Balsara4")) || + (CCTK_Equals(shock_case,"Balsara5")))) + { + CCTK_PARAMWARN("That test requires MHD! Set Bvec_evolution_method=GRHYDRO!"); + } }

File [modified]: GRHydro_C2P2C.F90 Delta lines: +4 -16 =================================================================== --- trunk/src/GRHydro_C2P2C.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_C2P2C.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -34,13 +34,6 @@ DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS

-#if !USE_EOS_OMNI -#ifdef _EOS_BASE_INC_ -#undef _EOS_BASE_INC_ -#endif -#include "EOS_Base.inc" -#endif - integer didit,i,j,k,nx,ny,nz CCTK_REAL det CCTK_REAL uxx,uxy,uxz,uyy,uyz,uzz @@ -52,7 +45,6 @@ CCTK_INT C2P_failed logical epsnegative

-#if USE_EOS_OMNI ! begin EOS Omni vars integer :: n = 1 integer :: keytemp = 0 @@ -63,7 +55,6 @@ real*8 :: xtemp = 0.0d0 real*8 :: xye = 0.0d0 ! end EOS Omni vars -#endif

call CCTK_WARN(1,"This test works only with Ideal_Fluid EoS")

@@ -105,16 +96,13 @@ epsnegative = .false.

GRHydro_rho_min = 1.e-10 -#if USE_EOS_OMNI - call EOS_Omni_press(GRHydro_eos_handle,keytemp,n,& + + call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,1.0d0,xtemp,xye,pmin,keyerr,anyerr)

- call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,n,& + call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,xeps,xtemp,xye,pmin,epsmin,keyerr,anyerr) -#else - pmin = EOS_Pressure(GRHydro_eos_handle, GRHydro_rho_min, 1.0d0) - epsmin = EOS_SpecificIntEnergy(GRHydro_eos_handle, GRHydro_rho_min, pmin) -#endif + C2P_failed = 0

write(*,*) 'C2P2C test: initial values.'

File [modified]: GRHydro_C2P2CM.F90 Delta lines: +11 -18 =================================================================== --- trunk/src/GRHydro_C2P2CM.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_C2P2CM.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -34,13 +34,6 @@ DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS

-#if !USE_EOS_OMNI -#ifdef _EOS_BASE_INC_ -#undef _EOS_BASE_INC_ -#endif -#include "EOS_Base.inc" -#endif - integer didit,i,j,k,nx,ny,nz CCTK_REAL det CCTK_REAL uxx,uxy,uxz,uyy,uyz,uzz @@ -49,11 +42,10 @@ CCTK_REAL rho_send,velx_send,vely_send,velz_send,eps_send CCTK_REAL press_send,w_lorentz_send,x_send,y_send,z_send,r_send CCTK_REAL bvcx_send, bvcy_send, bvcz_send, b2_send - CCTK_REAL pmin, epsmin + CCTK_REAL pmin, epsmin, local_gam CCTK_INT C2P_failed logical epsnegative

-#if USE_EOS_OMNI ! begin EOS Omni vars integer :: n = 1 integer :: keytemp = 0 @@ -64,7 +56,6 @@ real*8 :: xtemp = 0.0d0 real*8 :: xye = 0.0d0 ! end EOS Omni vars -#endif

call CCTK_WARN(1,"This test works only with Ideal_Fluid EoS")

@@ -110,16 +101,18 @@ epsnegative = .false.

GRHydro_rho_min = 1.e-10 -#if USE_EOS_OMNI - call EOS_Omni_press(GRHydro_eos_handle,keytemp,n,& + + call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,1.0d0,xtemp,xye,pmin,keyerr,anyerr)

- call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,n,& + call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,xeps,xtemp,xye,pmin,epsmin,keyerr,anyerr) -#else - pmin = EOS_Pressure(GRHydro_eos_handle, GRHydro_rho_min, 1.0d0) - epsmin = EOS_SpecificIntEnergy(GRHydro_eos_handle, GRHydro_rho_min, pmin) -#endif + + local_gam=0.0 + 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 + C2P_failed = 0

write(*,*) 'C2P2CM test: initial values.' @@ -136,7 +129,7 @@ write(*,*) ' Bvecz : ',bvcz_send

write(*,*) 'C2P2CM test: getting the associated primitive variables.' - call GRHydro_Con2PrimM_pt(GRHydro_eos_handle,dens_send,sx_send,sy_send,sz_send, & + call GRHydro_Con2PrimM_pt(GRHydro_eos_handle,local_gam,dens_send,sx_send,sy_send,sz_send, & tau_send,rho_send,velx_send,vely_send,velz_send, & eps_send,press_send,w_lorentz_send, & gxx_send,gxy_send,gxz_send,gyy_send,gyz_send,gzz_send,&

File [modified]: GRHydro_Con2Prim.F90 Delta lines: +4 -17 =================================================================== --- trunk/src/GRHydro_Con2Prim.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_Con2Prim.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -52,14 +52,6 @@ CCTK_INT C2P_failed logical epsnegative

-#ifndef USE_EOS_OMNI -#ifdef _EOS_BASE_INC_ -#undef _EOS_BASE_INC_ -#endif -#include "EOS_Base.inc" -#endif - -#if USE_EOS_OMNI ! begin EOS Omni vars integer :: n = 1 integer :: keytemp = 0 @@ -70,7 +62,6 @@ real*8 :: xtemp = 0.0d0 real*8 :: xye = 0.0d0 ! end EOS Omni vars -#endif

call CCTK_WARN(1,"For this test, remember to use a polytropic EoS and to set eos_gamma = 2.0 and eos_k = 100.0")

@@ -107,16 +98,12 @@

epsnegative = .false.

-#if USE_EOS_OMNI - call EOS_Omni_press(GRHydro_eos_handle,keytemp,n,& - GRHydro_rho_min,1.0d0,xtemp,xye,pmin,keyerr,anyerr) + call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& + GRHydro_rho_min,xeps,xtemp,xye,pmin,keyerr,anyerr)

- call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,n,& + call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,xeps,xtemp,xye,pmin,epsmin,keyerr,anyerr) -#else - pmin = EOS_Pressure(GRHydro_eos_handle, GRHydro_rho_min, 1.0d0) - epsmin = EOS_SpecificIntEnergy(GRHydro_eos_handle, GRHydro_rho_min, pmin) -#endif + C2P_failed = 0

write(*,*) 'Con2Prim test: converting to primitive variables'

File [modified]: GRHydro_P2C2P.F90 Delta lines: +6 -19 =================================================================== --- trunk/src/GRHydro_P2C2P.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_P2C2P.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -34,13 +34,6 @@ DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS

-#if !USE_EOS_OMNI -#ifdef _EOS_BASE_INC_ -#undef _EOS_BASE_INC_ -#endif -#include "EOS_Base.inc" -#endif - integer didit,i,j,k,nx,ny,nz CCTK_REAL det CCTK_REAL uxx,uxy,uxz,uyy,uyz,uzz @@ -52,7 +45,6 @@ CCTK_INT C2P_failed logical epsnegative

-#if USE_EOS_OMNI ! begin EOS Omni vars integer :: n = 1 integer :: keytemp = 0 @@ -63,7 +55,6 @@ real*8 :: xtemp = 0.0d0 real*8 :: xye = 0.0d0 ! end EOS Omni vars -#endif

call CCTK_WARN(1,"This test works only with Ideal_Fluid EoS")

@@ -103,19 +94,15 @@ epsnegative = .false.

GRHydro_rho_min = 1.e-10 -#if USE_EOS_OMNI - call EOS_Omni_press(GRHydro_eos_handle,keytemp,n,& - rho_send,1.0d0,xtemp,eps_send,press_send,keyerr,anyerr) - call EOS_Omni_press(GRHydro_eos_handle,keytemp,n,& + + call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& + rho_send,eps_send,xtemp,xye,press_send,keyerr,anyerr) + call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,1.0d0,xtemp,xye,pmin,keyerr,anyerr)

- call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,n,& + call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,xeps,xtemp,xye,pmin,epsmin,keyerr,anyerr) -#else - press_send = EOS_Pressure(GRHydro_eos_handle, rho_send, eps_send) - pmin = EOS_Pressure(GRHydro_eos_handle, GRHydro_rho_min, 1.0d0) - epsmin = EOS_SpecificIntEnergy(GRHydro_eos_handle, GRHydro_rho_min, pmin) -#endif + C2P_failed = 0

write(*,*) 'P2C2P test: the primitive variables are'

File [modified]: GRHydro_P2C2PM.F90 Delta lines: +13 -21 =================================================================== --- trunk/src/GRHydro_P2C2PM.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_P2C2PM.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -34,13 +34,6 @@ DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS

-#if !USE_EOS_OMNI -#ifdef _EOS_BASE_INC_ -#undef _EOS_BASE_INC_ -#endif -#include "EOS_Base.inc" -#endif - integer didit,i,j,k,nx,ny,nz CCTK_REAL det CCTK_REAL uxx,uxy,uxz,uyy,uyz,uzz @@ -49,11 +42,10 @@ CCTK_REAL rho_send,velx_send,vely_send,velz_send,eps_send CCTK_REAL press_send,w_lorentz_send,x_send,y_send,z_send,r_send CCTK_REAL bvcx_send,bvcy_send,bvcz_send,b2_send - CCTK_REAL pmin, epsmin + CCTK_REAL pmin, epsmin, local_gam CCTK_INT C2P_failed logical epsnegative

-#if USE_EOS_OMNI ! begin EOS Omni vars integer :: n = 1 integer :: keytemp = 0 @@ -64,7 +56,6 @@ real*8 :: xtemp = 0.0d0 real*8 :: xye = 0.0d0 ! end EOS Omni vars -#endif

call CCTK_WARN(1,"This test works only with Ideal_Fluid EoS")

@@ -108,19 +99,19 @@ epsnegative = .false.

GRHydro_rho_min = 1.e-10 -#if USE_EOS_OMNI - call EOS_Omni_press(GRHydro_eos_handle,keytemp,n,& - rho_send,1.0d0,xtemp,eps_send,press_send,keyerr,anyerr) - call EOS_Omni_press(GRHydro_eos_handle,keytemp,n,& + + call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& + rho_send,eps_send,xtemp,xye,press_send,keyerr,anyerr) + call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,1.0d0,xtemp,xye,pmin,keyerr,anyerr)

- call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,n,& + call EOS_Omni_EpsFromPress(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& GRHydro_rho_min,xeps,xtemp,xye,pmin,epsmin,keyerr,anyerr) -#else - press_send = EOS_Pressure(GRHydro_eos_handle, rho_send, eps_send) - pmin = EOS_Pressure(GRHydro_eos_handle, GRHydro_rho_min, 1.0d0) - epsmin = EOS_SpecificIntEnergy(GRHydro_eos_handle, GRHydro_rho_min, pmin) -#endif + + call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& + 1.0d0,1.0d0,xtemp,xye,local_gam,keyerr,anyerr) + local_gam=local_gam+1.0 + C2P_failed = 0

write(*,*) 'P2C2PM test: the primitive variables are' @@ -155,7 +146,8 @@ write(*,*) ' Bvecz : ',bvcz_send

write(*,*) 'P2C2PM test: getting the associated primitive variables.' - call GRHydro_Con2PrimM_pt(GRHydro_eos_handle,dens_send,sx_send,sy_send,sz_send, & + + call GRHydro_Con2PrimM_pt(GRHydro_eos_handle,local_gam,dens_send,sx_send,sy_send,sz_send, & tau_send,rho_send,velx_send,vely_send,velz_send, & eps_send,press_send,w_lorentz_send, & gxx_send,gxy_send,gxz_send,gyy_send,gyz_send,gzz_send,&

File [added]: GRHydro_P2C2PM_polytype.F90 Delta lines: +181 -0 =================================================================== --- trunk/src/GRHydro_P2C2PM_polytype.F90 (rev 0) +++ trunk/src/GRHydro_P2C2PM_polytype.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -0,0 +1,181 @@ + /*@@ + @file GRHydro_P2C2PM_polytype.F90 + @date Sep 25, 2010 + @author Joshua Faber, Scott Noble, Bruno Mundim, Luca Baiotti + @desc + A test of the conservative <--> primitive variable exchange + @enddesc + @@*/ + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "cctk_Arguments.h" + + /*@@ + @routine p2c2pm + @date Sep 25, 2010 + @author Joshua Faber, Scott Noble, Bruno Mundim, Luca Baiotti + @desc + Testing the conservative <--> primitive variable transformations. + The values before and after should match. + @enddesc + @calls + @calledby + @history + + @endhistory + +@@*/ + +subroutine p2c2pm_polytype(CCTK_ARGUMENTS) + + implicit none + + DECLARE_CCTK_ARGUMENTS + DECLARE_CCTK_PARAMETERS + + integer didit,i,j,k,nx,ny,nz + CCTK_REAL det, local_K + CCTK_REAL uxx,uxy,uxz,uyy,uyz,uzz + CCTK_REAL gxx_send,gxy_send,gxz_send,gyy_send,gyz_send,gzz_send + CCTK_REAL dens_send,sx_send,sy_send,sz_send,tau_send + CCTK_REAL rho_send,velx_send,vely_send,velz_send,eps_send + CCTK_REAL press_send,w_lorentz_send,x_send,y_send,z_send,r_send + CCTK_REAL bvcx_send,bvcy_send,bvcz_send,b2_send + CCTK_REAL pmin, epsmin, local_gam, sc_send + CCTK_INT C2P_failed + logical epsnegative + +! begin EOS Omni vars + integer :: n = 1 + integer :: keytemp = 0 + integer :: anyerr = 0 + integer :: keyerr(1) = 0 + real*8 :: xpress = 0.0d0 + real*8 :: xeps = 0.0d0 + real*8 :: xtemp = 0.0d0 + real*8 :: xye = 0.0d0 +! end EOS Omni vars + + call CCTK_WARN(1,"This test works only with Ideal_Fluid EoS") + + nx = cctk_lsh(1) + ny = cctk_lsh(2) + nz = cctk_lsh(3) + + x_send = 0.0d0 + y_send = 0.0d0 + z_send = 0.0d0 + r_send = 0.0d0 + + gxx_send = 1.0d0 + gyy_send = 1.0d0 + gzz_send = 1.0d0 + gxy_send = 0.0d0 + gxz_send = 0.0d0 + gyz_send = 0.0d0 + + det = 1.0d0 + + uxx = 1.0d0 + uyy = 1.0d0 + uzz = 1.0d0 + uxy = 0.0d0 + uxz = 0.0d0 + uyz = 0.0d0 + + rho_send = 1.29047362d0 + velx_send = 0.166666658d0 + vely_send = 0.166666658d0 + velz_send = 0.166666658d0 +!!$ eps_send = 0.484123939d0 + + bvcx_send = Bx_init + bvcy_send = By_init + bvcz_send = Bz_init + + w_lorentz_send = 1.d0/sqrt(1.0d0-velx_send*velx_send-vely_send*vely_send-velz_send*velz_send) + + epsnegative = .false. + + GRHydro_rho_min = 1.e-10 + + call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& + rho_send,1.0d0,xtemp,xye,press_send,keyerr,anyerr) + call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& + GRHydro_rho_min,1.0d0,xtemp,xye,pmin,keyerr,anyerr) + + call EOS_Omni_EpsFromPress(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& + rho_send,xeps,xtemp,xye,press_send,eps_send,keyerr,anyerr) + 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_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& + 1.d0,1.d0,xtemp,xye,xpress,keyerr,anyerr) + local_K = xpress + + local_gam=log(press_send/local_K)/log(rho_send) + + C2P_failed = 0 + + write(*,*) 'P2C2PM test: the primitive variables are' + write(*,*) ' primitive variables: ' + write(*,*) ' rho : ',rho_send + write(*,*) ' velx : ',velx_send + write(*,*) ' vely : ',vely_send + write(*,*) ' velz : ',velz_send + write(*,*) ' press : ',press_send + write(*,*) ' eps : ',eps_send + write(*,*) ' W : ',w_lorentz_send + write(*,*) ' Bvecx : ',bvcx_send + write(*,*) ' Bvecy : ',bvcy_send + write(*,*) ' Bvecz : ',bvcz_send + + write(*,*) 'P2C2PM test: converting back to conserved variables.' + call prim2conpolytypeM(GRHydro_polytrope_handle,gxx_send, gxy_send, gxz_send, gyy_send, gyz_send, gzz_send, det, & + dens_send, sx_send, sy_send, sz_send, tau_send, bvcx_send, bvcy_send, bvcz_send, rho_send, & + velx_send, vely_send, velz_send, eps_send, press_send, w_lorentz_send) + + write(*,*) 'P2C2PM test: initial values.' + write(*,*) ' conservative variables: ' + write(*,*) ' dens: ',dens_send + write(*,*) ' sx : ',sx_send + write(*,*) ' sy : ',sy_send + write(*,*) ' sz : ',sz_send + write(*,*) ' tau : ',tau_send + write(*,*) ' eps : ',eps_send + write(*,*) ' W : ',w_lorentz_send + write(*,*) ' Bvecx : ',bvcx_send + write(*,*) ' Bvecy : ',bvcy_send + write(*,*) ' Bvecz : ',bvcz_send + + write(6,*)'local_K:',local_K + sc_send=local_K*dens_send + + write(*,*) 'P2C2PM test: getting the associated primitive variables.' + call GRHydro_Con2PrimM_Polytype_pt(GRHydro_polytrope_handle,local_gam,dens_send,sx_send,sy_send,sz_send, & + sc_send,rho_send,velx_send,vely_send,velz_send, & + eps_send,press_send,w_lorentz_send, & + gxx_send,gxy_send,gxz_send,gyy_send,gyz_send,gzz_send,& + uxx,uxy,uxz,uyy,uyz,uzz,det,& + bvcx_send,bvcy_send,bvcz_send,b2_send,& + epsnegative,C2P_failed) + + write(*,*) 'P2C2PM test: the primitive variables are' + write(*,*) ' primitive variables: ' + write(*,*) ' rho : ',rho_send + write(*,*) ' velx : ',velx_send + write(*,*) ' vely : ',vely_send + write(*,*) ' velz : ',velz_send + write(*,*) ' press : ',press_send + write(*,*) ' eps : ',eps_send + write(*,*) ' W : ',w_lorentz_send + write(*,*) ' Bvecx : ',bvcx_send + write(*,*) ' Bvecy : ',bvcy_send + write(*,*) ' Bvecz : ',bvcz_send + + STOP + + return + +end subroutine p2c2pm_polytype

File [modified]: GRHydro_ReadConformalData.F90 Delta lines: +1 -14 =================================================================== --- trunk/src/GRHydro_ReadConformalData.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_ReadConformalData.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -43,18 +43,9 @@ DECLARE_CCTK_PARAMETERS DECLARE_CCTK_FUNCTIONS

-!#include "EOS_Base.h" -#ifdef _EOS_BASE_INC_ -#undef _EOS_BASE_INC_ -#endif -#include "EOS_Base.inc" - - CCTK_INT :: i,j,k,handle,ierr -! CCTK_REAL :: eos_k, eos_gamma CCTK_REAL :: rho_min, det

-#if USE_EOS_OMNI ! begin EOS Omni vars integer :: n = 1 integer :: poly_eoskey = 0 @@ -68,7 +59,6 @@ real*8 :: xye(1) = 0.0d0 poly_eoskey = EOS_Omni_GetHandle("2D_Polytrope") ! end EOS Omni vars -#endif

! only gxx has been read in; copy it into gyy and gzz as well gyy = gxx @@ -155,14 +145,11 @@ do i=1,cctk_lsh(1) do j=1,cctk_lsh(2) do k=1,cctk_lsh(3) -#ifdef USE_EOS_OMNI + call EOS_Omni_EpsFromPress(poly_eoskey,keytemp,rf_precision,& n,rho(i,j,k),xeps,xtemp,xye,& press(i,j,k),eps(i,j,k),keyerr,anyerr)

-#else - eps(i,j,k) = EOS_SpecificIntEnergy(GRHydro_eos_handle,rho(i,j,k),press(i,j,k)) -#endif end do end do end do

File [modified]: GRHydro_ShockTube.F90 Delta lines: +2 -2 =================================================================== --- trunk/src/GRHydro_ShockTube.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_ShockTube.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -11,6 +11,7 @@ #include "cctk_Parameters.h" #include "cctk_Arguments.h" #include "cctk_Functions.h" +#include "GRHydro_Macros.h"

#define velx(i,j,k) vel(i,j,k,1) #define vely(i,j,k) vel(i,j,k,2) @@ -124,8 +125,7 @@ eps(i,j,k) = epsr end if

- call SpatialDet(gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),& - gyy(i,j,k),gyz(i,j,k),gzz(i,j,k),det) + det=SPATIAL_DETERMINANT(gxx(i,j,k),gxy(i,j,k),gxz(i,j,k),gyy(i,j,k),gyz(i,j,k),gzz(i,j,k))

if (CCTK_EQUALS(GRHydro_eos_type,"Polytype")) then call Prim2ConPoly(GRHydro_eos_handle,gxx(i,j,k),gxy(i,j,k),&

File [modified]: GRHydro_ShockTubeM.F90 Delta lines: +93 -8 =================================================================== --- trunk/src/GRHydro_ShockTubeM.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_ShockTubeM.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -54,6 +54,13 @@ velzl, velzr, epsl, epsr CCTK_REAL :: bvcxl,bvcyl,bvczl,bvcxr,bvcyr,bvczr

+ bvcxl = Bx_init + bvcyl = By_init + bvczl = Bz_init + bvcxr = Bx_init + bvcyr = By_init + bvczr = Bz_init + nx = cctk_lsh(1) ny = cctk_lsh(2) nz = cctk_lsh(3) @@ -110,17 +117,95 @@ velzr = 0.d0 epsl = 1500.d0 epsr = 1.5d-2 - else + else if (CCTK_EQUALS(shock_case,"Balsara1")) then + rhol = 1.d0 + rhor = 0.125d0 + velxl = 0.d0 + velxr = 0.d0 + velyl = 0.d0 + velyr = 0.d0 + velzl = 0.d0 + velzr = 0.d0 + bvcxl=0.5d0 + bvcxr=0.5d0 + bvcyl=1.d0 + bvcyr=-1.d0 + bvczl=0.d0 + bvczr=0.d0 + epsl = 1.0/rhol + epsr = 0.1/rhor + else if (CCTK_EQUALS(shock_case,"Balsara2")) then + rhol = 1.d0 + rhor = 1.d0 + velxl = 0.d0 + velxr = 0.d0 + velyl = 0.d0 + velyr = 0.d0 + velzl = 0.d0 + velzr = 0.d0 + bvcxl=5.0d0 + bvcxr=5.0d0 + bvcyl=6.d0 + bvcyr=0.7d0 + bvczl=6.d0 + bvczr=0.7d0 + epsl = 1.5d0*30.0d0/rhol + epsr = 1.5d0*1.0d0/rhor + else if (CCTK_EQUALS(shock_case,"Balsara3")) then + rhol = 1.d0 + rhor = 1.d0 + velxl = 0.d0 + velxr = 0.d0 + velyl = 0.d0 + velyr = 0.d0 + velzl = 0.d0 + velzr = 0.d0 + bvcxl=10.0d0 + bvcxr=10.0d0 + bvcyl=7.d0 + bvcyr=0.7d0 + bvczl=7.d0 + bvczr=0.7d0 + epsl = 1.5d0*1000.0d0/rhol + epsr = 1.5d0*0.1d0/rhor + else if (CCTK_EQUALS(shock_case,"Balsara4")) then + rhol = 1.d0 + rhor = 1.d0 + velxl = 0.999d0 + velxr = -0.999d0 + velyl = 0.d0 + velyr = 0.d0 + velzl = 0.d0 + velzr = 0.d0 + bvcxl=10.0d0 + bvcxr=10.0d0 + bvcyl=7.d0 + bvcyr=-7.d0 + bvczl=7.d0 + bvczr=-7.d0 + epsl = 1.5d0*0.1d0/rhol + epsr = 1.5d0*0.1d0/rhor + else if (CCTK_EQUALS(shock_case,"Balsara5")) then + rhol = 1.08d0 + rhor = 1.d0 + velxl = 0.4d0 + velxr = -0.45d0 + velyl = 0.3d0 + velyr = -0.2d0 + velzl = 0.2d0 + velzr = 0.2d0 + bvcxl=2.0d0 + bvcxr=2.0d0 + bvcyl=0.3d0 + bvcyr=-0.7d0 + bvczl=0.3d0 + bvczr=0.5d0 + epsl = 1.5d0*0.95d0/rhol + epsr = 1.5d0*1.0d0/rhor + else call CCTK_WARN(0,"Shock case not recognized") end if

- bvcxl = Bx_init - bvcyl = By_init - bvczl = Bz_init - bvcxr = Bx_init - bvcyr = By_init - bvczr = Bz_init - if ( ((change_shock_direction==0).and.(direction .lt. 0.0d0)).or.& ((change_shock_direction==1).and.(direction .gt. 0.0d0)) ) then rho(i,j,k) = rhol

File [modified]: GRHydro_SimpleWave.F90 Delta lines: +13 -24 =================================================================== --- trunk/src/GRHydro_SimpleWave.F90 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/GRHydro_SimpleWave.F90 2010-12-25 06:44:22 UTC (rev 117) @@ -46,17 +46,9 @@ DECLARE_CCTK_PARAMETERS DECLARE_CCTK_FUNCTIONS

-#ifndef USE_EOS_OMNI -#ifdef _EOS_BASE_INC_ -#undef _EOS_BASE_INC_ -#endif -#include "EOS_Base.inc" -#endif - CCTK_INT :: i, j, k, nx, ny, nz CCTK_REAL :: dr, k1, k2, k3, k4, in_data, old_data, source_data, new_data, c_0, det, pi

-#if USE_EOS_OMNI ! begin EOS Omni vars integer :: n = 1 integer :: keytemp = 0 @@ -68,7 +60,6 @@ real*8 :: xye(1) = 0.0d0 real*8 :: rf_precision = 1.0d-10 ! end EOS Omni vars -#endif

call CCTK_INFO("Setting initial data for a simple wave as Anile Miller Motta")

@@ -143,21 +134,19 @@

! atmosphere

- if ( (rho(i,1,1) < GRHydro_rho_min).OR.(velx(i,1,1) < 0) ) then - rho(i,1,1) = rho_abs_min -! rho(i,1,1) = 1.0 !the value of rho_min for the initial data - eps(i,1,1) = rho_abs_min**(1.d0/3.d0) - velx(i,1,1) = 0.d0 - w_lorentz(i,1,1) = 1.d0 -#if USE_EOS_OMNI - xeps = 1.0d0 - call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,rf_precision,n,& - rho(i,1,1),xeps,xtemp,xye,press(i,1,1),keyerr,anyerr) -#else - press(i,1,1) = EOS_Pressure(GRHydro_polytrope_handle, rho(i,1,1), 1.0d0) -#endif - ! polytrope only (initial data) - end if + if ( (rho(i,1,1) < GRHydro_rho_min).OR.(velx(i,1,1) < 0) ) then + rho(i,1,1) = rho_abs_min + ! rho(i,1,1) = 1.0 !the value of rho_min for the initial data + eps(i,1,1) = rho_abs_min**(1.d0/3.d0) + velx(i,1,1) = 0.d0 + w_lorentz(i,1,1) = 1.d0 + + xeps = 1.0d0 + call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,rf_precision,n,& + rho(i,1,1),xeps,xtemp,xye,press(i,1,1),keyerr,anyerr) + + ! polytrope only (initial data) + end if

! write(*,*) 'p',i, x(i,1,1), rho(i,1,1)**(4.d0/3.d0)/press(i,1,1)

File [modified]: make.code.defn Delta lines: +2 -1 =================================================================== --- trunk/src/make.code.defn 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/src/make.code.defn 2010-12-25 06:44:22 UTC (rev 117) @@ -13,7 +13,8 @@ GRHydro_SimpleWave.F90 \ CheckParam.c\ GRHydro_P2C2P.F90 \ - GRHydro_P2C2PM.F90 + GRHydro_P2C2PM.F90 \ + GRHydro_P2C2PM_polytype.F90

# Subdirectories containing source files SUBDIRS =

Directory: /trunk/ ==================

File [removed]: interface.ccl.omni Delta lines: +0 -84 =================================================================== --- trunk/interface.ccl.omni 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/interface.ccl.omni 2010-12-25 06:44:22 UTC (rev 117) @@ -1,84 +0,0 @@ -# Interface definition for thorn GRHydro_Init_Data -# $Header$ - -implements: GRHydro_init_data -inherits: GRHydro grid - -#USES INCLUDE: EOS_Base.inc -USES INCLUDE: SpaceMask.h -#USES INCLUDE: EOS_Base.h - - -SUBROUTINE SpatialDet(CCTK_REAL IN gxx, CCTK_REAL IN gxy, \ - CCTK_REAL IN gxz, CCTK_REAL IN gyy, \ - CCTK_REAL IN gyz, CCTK_REAL IN gzz, \ - CCTK_REAL OUT det) - - -SUBROUTINE Prim2ConPoly(CCTK_INT IN handle, \ - CCTK_REAL IN gxx, CCTK_REAL IN gxy, CCTK_REAL IN gxz, \ - CCTK_REAL IN gyy, CCTK_REAL IN gyz, CCTK_REAL IN gzz, \ - CCTK_REAL IN det, CCTK_REAL OUT dens, \ - CCTK_REAL OUT sx, CCTK_REAL OUT sy, \ - CCTK_REAL OUT sz, CCTK_REAL OUT tau, \ - CCTK_REAL IN rho, CCTK_REAL IN velx, \ - CCTK_REAL IN vely, \ - CCTK_REAL IN velz, CCTK_REAL OUT epsilon, \ - CCTK_REAL OUT press, CCTK_REAL OUT w_lorentz) - - -SUBROUTINE Prim2ConGen(CCTK_INT IN handle, \ - CCTK_REAL IN gxx, CCTK_REAL IN gxy, \ - CCTK_REAL IN gxz, CCTK_REAL IN gyy, \ - CCTK_REAL IN gyz, CCTK_REAL IN gzz, \ - CCTK_REAL IN det, CCTK_REAL OUT dens, \ - CCTK_REAL OUT sx, CCTK_REAL OUT sy, \ - CCTK_REAL OUT sz, CCTK_REAL OUT tau, \ - CCTK_REAL IN rho, CCTK_REAL IN velx, \ - CCTK_REAL IN vely, \ - CCTK_REAL IN velz, CCTK_REAL IN epsilon, \ - CCTK_REAL OUT press, CCTK_REAL OUT w_lorentz) - -SUBROUTINE Con2PrimPoly(CCTK_INT IN handle, CCTK_REAL OUT dens, \ - CCTK_REAL OUT sx, CCTK_REAL OUT sy, \ - CCTK_REAL OUT sz, CCTK_REAL OUT tau, \ - CCTK_REAL OUT rho, CCTK_REAL OUT velx, \ - CCTK_REAL OUT vely, CCTK_REAL OUT velz, \ - CCTK_REAL OUT epsilon, CCTK_REAL OUT press, \ - CCTK_REAL OUT w_lorentz, CCTK_REAL IN uxx, \ - CCTK_REAL IN uxy, CCTK_REAL IN uxz, CCTK_REAL IN uyy, \ - CCTK_REAL IN uyz, CCTK_REAL IN uzz, CCTK_REAL IN det, \ - CCTK_REAL IN x, CCTK_REAL IN y, CCTK_REAL IN z, \ - CCTK_REAL IN r, CCTK_REAL IN rho_min, \ - CCTK_INT IN GRHydro_reflevel, CCTK_REAL IN GRHydro_C2P_failed) - - -USES FUNCTION SpatialDet -USES FUNCTION Prim2ConPoly -USES FUNCTION Prim2ConGen -USES FUNCTION Con2PrimPoly - -protected: - -CCTK_REAL simple_wave_grid_functions TYPE=GF TIMELEVELS=1 tags='checkpoint="no"' -{ - simple_tmp - c_s -} "1D arrays for the simple-wave routine" - -CCTK_REAL simple_wave_scalars TYPE=scalar -{ - simple_rho_0 - simple_eps_0 -} "values at v=0" - -CCTK_REAL simple_wave_output TYPE=GF TIMELEVELS=1 tags='checkpoint="no"' -{ - simple_rho - simple_eps -# simple_entropy -} "output variables for the simple-wave routine" - -private: - -CCTK_INT GRHydro_init_data_reflevel type = SCALAR tags='checkpoint="no"' "Refinement level GRHydro is working on right now"

File [modified]: param.ccl Delta lines: +7 -0 =================================================================== --- trunk/param.ccl 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/param.ccl 2010-12-25 06:44:22 UTC (rev 117) @@ -22,6 +22,7 @@ "con2primtest" :: "Testing the con -> prim conversion" "con2prim2con_test" :: "Testing the con -> prim -> con conversion" "prim2con2prim_test" :: "Testing the prim -> con -> prim conversion" + "prim2con2prim_polytype_test" :: "Testing the prim -> con -> prim conversion - polytype version" "reconstruction_test" :: "Testing reconstruction" }

@@ -41,6 +42,11 @@ "Simple" :: "GRAstro_Hydro test case" "Sod" :: "Sod's problem" "Blast" :: "Strong blast wave" + "Balsara1" :: "Balsara Test #1" + "Balsara2" :: "Balsara Test #2" + "Balsara3" :: "Balsara Test #3" + "Balsara4" :: "Balsara Test #4" + "Balsara5" :: "Balsara Test #5" } "Sod"

REAL shock_xpos "Position of shock plane: x" @@ -106,6 +112,7 @@

shares:GRHydro

+USES real GRHydro_eos_rf_prec USES real GRHydro_rho_central USES real GRHydro_eps_min "" USES real GRHydro_perc_ptol ""

File [modified]: schedule.ccl Delta lines: +16 -0 =================================================================== --- trunk/schedule.ccl 2010-11-21 03:29:35 UTC (rev 116) +++ trunk/schedule.ccl 2010-12-25 06:44:22 UTC (rev 117) @@ -77,6 +77,22 @@ } }

+if (CCTK_Equals(initial_data,"prim2con2prim_polytype_test")) { + STORAGE:GRHydro_init_data_reflevel + schedule GRHydro_Init_Data_RefinementLevel IN HydroBase_Initial BEFORE p2c2p_call + { + LANG: Fortran + } "Calculate current refinement level" + + if(CCTK_Equals(Bvec_evolution_method,"GRHydro")) + { + schedule p2c2pM_polytype in HydroBase_Initial AS p2c2p_call + { + LANG: Fortran + } "Testing primitive to conservative to primitive - MHD polytype version" + } +} + if (CCTK_Equals(initial_data,"reconstruction_test")) { schedule GRHydro_reconstruction_test in HydroBase_Initial {