diff options
-rw-r--r-- | src/GRHydro_Con2PrimM.F90 | 67 |
1 files changed, 35 insertions, 32 deletions
diff --git a/src/GRHydro_Con2PrimM.F90 b/src/GRHydro_Con2PrimM.F90 index 43c8374..67774ca 100644 --- a/src/GRHydro_Con2PrimM.F90 +++ b/src/GRHydro_Con2PrimM.F90 @@ -48,20 +48,22 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS + DECLARE_CCTK_FUNCTIONS integer :: i, j, k, itracer, nx, ny, nz - CCTK_REAL :: uxx, uxy, uxz, uyy, uyz, uzz, det, pmin, epsmin + CCTK_REAL :: uxx, uxy, uxz, uyy, uyz, uzz, det, pmin(1), epsmin(1) CCTK_REAL :: b2 CCTK_INT :: epsnegative character(len=100) warnline - CCTK_REAL :: local_min_tracer, local_gam, local_pgam,local_K,sc + CCTK_REAL :: local_min_tracer, local_gam(1), local_pgam,local_K, sc ! begin EOS Omni vars integer :: n,keytemp,anyerr,keyerr(1) - real*8 :: xpress,xtemp,xye,xeps,xrho + real*8 :: xpress(1),xtemp(1),xye(1),xeps(1),xrho(1),one(1)=1.0d0 + n=1;keytemp=0;anyerr=0;keyerr(1)=0 - xpress=0.0d0;xtemp=0.0d0;xye=0.0d0;xeps=0.0d0 + xpress(1)=0.0d0;xtemp(1)=0.0d0;xye(1)=0.0d0;xeps(1)=0.0d0 ! end EOS Omni vars nx = cctk_lsh(1) @@ -75,19 +77,19 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) end if ! this is a poly call + xrho(1)=GRHydro_rho_min call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& - GRHydro_rho_min,xeps,xtemp,xye,pmin,keyerr,anyerr) + xrho,xeps,xtemp,xye,pmin,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) + xrho,xeps,xtemp,xye,xpress,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) + one,one,xtemp,xye,local_gam,keyerr,anyerr) local_gam = local_gam+1.d0 - !$OMP PARALLEL DO PRIVATE(i,j,k,itracer,& !$OMP uxx, uxy, uxz, uyy, uyz, uzz, det, epsnegative, & - !$OMP b2,xrho,xeps,xpress,xtemp,local_K,local_pgam,sc,keyerr,anyerr) + !$OMP b2,xrho,xeps,xpress,xtemp,local_K,local_pgam,sc,keyerr,anyerr ) do k = 1, nz do j = 1, ny do i = 1, nx @@ -151,11 +153,11 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) cycle end if - + if(evolve_temper.eq.0) then - call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam, dens(i,j,k), & + call GRHydro_Con2PrimM_pt(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), & - Bcons(i,j,k,1), Bcons(i,j,k,2), Bcons(i,j,k,3), rho(i,j,k),& + Bcons(i,j,k,1),Bcons(i,j,k,2),Bcons(i,j,k,3), 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),& Bvec(i,j,k,1), Bvec(i,j,k,2), Bvec(i,j,k,3),b2, 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), & @@ -221,12 +223,12 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) xrho=1.0d0; xtemp=0.0d0; xeps=1.0d0 call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& xrho,xeps,xtemp,xye,xpress,keyerr,anyerr) - local_K = xpress; + local_K = xpress(1); xrho=10.0d0; xeps=1.0d0 call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& xrho,xeps,xtemp,xye,xpress,keyerr,anyerr) - local_pgam=log(xpress/local_K)/log(xrho) + 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, local_pgam, dens(i,j,k), & @@ -280,13 +282,14 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS) DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS + DECLARE_CCTK_FUNCTIONS integer :: i, j, k, itracer, nx, ny, nz CCTK_REAL :: uxxl, uxyl, uxzl, uyyl, uyzl, uzzl,& - uxxr, uxyr, uxzr, uyyr, uyzr, uzzr, pmin, epsmin + uxxr, uxyr, uxzr, uyyr, uyzr, uzzr, pmin(1), epsmin(1) CCTK_REAL :: gxxl,gxyl,gxzl,gyyl,gyzl,gzzl,avg_detl,& gxxr,gxyr,gxzr,gyyr,gyzr,gzzr,avg_detr - CCTK_REAL :: b2minus, b2plus, local_gam, local_pgam,local_K,scminus,scplus + CCTK_REAL :: b2minus, b2plus, local_gam(1), local_pgam,local_K,scminus,scplus CCTK_INT :: epsnegative character(len=100) warnline @@ -294,20 +297,21 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS) ! begin EOS Omni vars integer :: n,keytemp,anyerr,keyerr(1) - real*8 :: xpress,xtemp,xye,xeps,xrho + real*8 :: xpress(1),xtemp(1),xye(1),xeps(1),xrho(1),one(1)=1.0d0 n=1;keytemp=0;anyerr=0;keyerr(1)=0 - xpress=0.0d0;xeps=0.0d0;xtemp=0.0d0;xye=0.0d0 + xpress(1)=0.0d0;xeps(1)=0.0d0;xtemp(1)=0.0d0;xye(1)=0.0d0 ! end EOS Omni vars ! this is a poly call + xrho(1)=GRHydro_rho_min call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& - GRHydro_rho_min,1.0d0,xtemp,xye,pmin,keyerr,anyerr) + xrho,one,xtemp,xye,pmin,keyerr,anyerr) call EOS_Omni_EpsFromPress(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& - GRHydro_rho_min,epsmin,xtemp,xye,pmin,epsmin,keyerr,anyerr) + xrho,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) + one,one,xtemp,xye,local_gam,keyerr,anyerr) local_gam=local_gam+1.0 nx = cctk_lsh(1) @@ -370,7 +374,7 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS) endif - call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam,densminus(i,j,k), & + call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam(1),densminus(i,j,k), & sxminus(i,j,k),syminus(i,j,k),szminus(i,j,k), tauminus(i,j,k), & Bconsxminus(i,j,k), Bconsyminus(i,j,k), Bconszminus(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),& @@ -386,12 +390,12 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS) 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 + one,one,xtemp,xye,xpress,keyerr,anyerr) + local_K = xpress(1) 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) + xrho,one,xtemp,xye,xpress,keyerr,anyerr) + local_pgam=log(xpress(1)/local_K)/log(xrho(1)) scminus = local_K*densminus(i,j,k) call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, local_pgam, densminus(i,j,k), & @@ -427,7 +431,7 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS) endif epsnegative = 0 - call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam, densplus(i,j,k), & + call GRHydro_Con2PrimM_pt(GRHydro_eos_handle, local_gam(1), densplus(i,j,k), & sxplus(i,j,k),syplus(i,j,k),szplus(i,j,k), tauplus(i,j,k),& Bconsxplus(i,j,k), Bconsyplus(i,j,k), Bconszplus(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),& @@ -443,12 +447,12 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS) 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 + one,one,xtemp,xye,xpress,keyerr,anyerr) + local_K = xpress(1) 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) + xrho,one,xtemp,xye,xpress,keyerr,anyerr) + local_pgam=log(xpress(1)/local_K)/log(xrho(1)) scplus = local_K*densplus(i,j,k) call GRHydro_Con2PrimM_Polytype_pt(GRHydro_eos_handle, local_pgam, densplus(i,j,k), & @@ -490,7 +494,6 @@ subroutine Conservative2PrimitiveBoundsM(CCTK_ARGUMENTS) end subroutine Conservative2PrimitiveBoundsM - /*@@ @routine Con2PrimPolytypeM @date Sep 16, 2010 |