diff options
Diffstat (limited to 'src/GRHydro_Con2Prim.F90')
| -rw-r--r-- | src/GRHydro_Con2Prim.F90 | 117 |
1 files changed, 57 insertions, 60 deletions
diff --git a/src/GRHydro_Con2Prim.F90 b/src/GRHydro_Con2Prim.F90 index 9d33962..176ba77 100644 --- a/src/GRHydro_Con2Prim.F90 +++ b/src/GRHydro_Con2Prim.F90 @@ -54,7 +54,7 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) DECLARE_CCTK_PARAMETERS integer :: i, j, k, itracer, nx, ny, nz - CCTK_REAL :: uxx, uxy, uxz, uyy, uyz, uzz, det, pmin, epsmin, dummy1, dummy2 + CCTK_REAL :: uxx, uxy, uxz, uyy, uyz, uzz, pmin, epsmin, dummy1, dummy2 logical :: epsnegative character*256 :: warnline @@ -115,7 +115,7 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) GRHydro_rho_min,xeps,xtemp,xye,pmin,epsmin,keyerr,anyerr) !$OMP PARALLEL DO PRIVATE(i,j,k,itracer,& - !$OMP uxx, uxy, uxz, uyy, uyz, uzz, det, epsnegative, anyerr, keyerr, keytemp,& + !$OMP uxx, uxy, uxz, uyy, uyz, uzz, epsnegative, anyerr, keyerr, keytemp,& !$OMP warnline, dummy1, dummy2) do k = 1, nz do j = 1, ny @@ -160,15 +160,13 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) epsnegative = .false. - det = SPATIAL_DETERMINANT(g11(i,j,k),g12(i,j,k),g13(i,j,k), \ - g22(i,j,k),g23(i,j,k),g33(i,j,k)) - call UpperMetric(uxx,uxy,uxz,uyy,uyz,uzz,det,& + call UpperMetric(uxx,uxy,uxz,uyy,uyz,uzz,sdetg(i,j,k)*sdetg(i,j,k),& g11(i,j,k),g12(i,j,k),g13(i,j,k),g22(i,j,k),& g23(i,j,k),g33(i,j,k)) ! In excised region, set to atmosphere! if (GRHydro_enable_internal_excision /= 0 .and. (hydro_excision_mask(i,j,k) .gt. 0)) then - SET_ATMO_MIN(dens(i,j,k), sqrt(det)*GRHydro_rho_min, r(i,j,k)) !sqrt(det)*GRHydro_rho_min !/(1.d0+GRHydro_atmo_tolerance) + SET_ATMO_MIN(dens(i,j,k), sdetg(i,j,k)*GRHydro_rho_min, r(i,j,k)) !sqrt(det)*GRHydro_rho_min !/(1.d0+GRHydro_atmo_tolerance) SET_ATMO_MIN(rho(i,j,k), GRHydro_rho_min, r(i,j,k)) !GRHydro_rho_min scon(i,j,k,:) = 0.d0 vup(i,j,k,:) = 0.d0 @@ -193,7 +191,7 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) endif ! w_lorentz=1, so the expression for tau reduces to: - tau(i,j,k) = sqrt(det) * (rho(i,j,k)+rho(i,j,k)*eps(i,j,k)) - dens(i,j,k) + tau(i,j,k) = sdetg(i,j,k) * (rho(i,j,k)+rho(i,j,k)*eps(i,j,k)) - dens(i,j,k) cycle endif @@ -223,8 +221,8 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) endif !if ( dens(i,j,k) .le. sqrt(det)*GRHydro_rho_min*(1.d0+GRHydro_atmo_tolerance) ) then - IF_BELOW_ATMO(dens(i,j,k), sqrt(det)*GRHydro_rho_min, GRHydro_atmo_tolerance, r(i,j,k)) then - SET_ATMO_MIN(dens(i,j,k), sqrt(det)*GRHydro_rho_min, r(i,j,k)) !sqrt(det)*GRHydro_rho_min !/(1.d0+GRHydro_atmo_tolerance) + IF_BELOW_ATMO(dens(i,j,k), sdetg(i,j,k)*GRHydro_rho_min, GRHydro_atmo_tolerance, r(i,j,k)) then + SET_ATMO_MIN(dens(i,j,k), sdetg(i,j,k)*GRHydro_rho_min, r(i,j,k)) !sqrt(det)*GRHydro_rho_min !/(1.d0+GRHydro_atmo_tolerance) SET_ATMO_MIN(rho(i,j,k), GRHydro_rho_min, r(i,j,k)) !GRHydro_rho_min scon(i,j,k,:) = 0.d0 vup(i,j,k,:) = 0.d0 @@ -249,7 +247,7 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) endif ! w_lorentz=1, so the expression for tau reduces to: - tau(i,j,k) = sqrt(det) * (rho(i,j,k)+rho(i,j,k)*eps(i,j,k)) - dens(i,j,k) + tau(i,j,k) = sdetg(i,j,k) * (rho(i,j,k)+rho(i,j,k)*eps(i,j,k)) - dens(i,j,k) cycle @@ -260,7 +258,7 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) 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),vup(i,j,k,1),vup(i,j,k,2), & vup(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k), & - uxx,uxy,uxz,uyy,uyz,uzz,det,x(i,j,k),y(i,j,k), & + uxx,uxy,uxz,uyy,uyz,uzz,sdetg(i,j,k),x(i,j,k),y(i,j,k), & z(i,j,k),r(i,j,k),epsnegative,GRHydro_rho_min,pmin, epsmin, & GRHydro_reflevel, GRHydro_C2P_failed(i,j,k)) else @@ -292,7 +290,7 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) !$OMP END CRITICAL ! for safety, let's set the point to atmosphere - dens(i,j,k) = ATMOMIN(sqrt(det)*GRHydro_rho_min, r(i,j,k)) !sqrt(det)*GRHydro_rho_min !/(1.d0+GRHydro_atmo_tolerance) + dens(i,j,k) = ATMOMIN(sdetg(i,j,k)*GRHydro_rho_min, r(i,j,k)) !sqrt(det)*GRHydro_rho_min !/(1.d0+GRHydro_atmo_tolerance) rho(i,j,k) = ATMOMIN(GRHydro_rho_min, r(i,j,k)) !GRHydro_rho_min scon(i,j,k,:) = 0.d0 vup(i,j,k,:) = 0.d0 @@ -305,7 +303,7 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) rho(i,j,k),eps(i,j,k),xtemp,xye,press(i,j,k),eps(i,j,k),keyerr,anyerr) ! w_lorentz=1, so the expression for tau reduces to: - tau(i,j,k) = sqrt(det) * (rho(i,j,k)+rho(i,j,k)*eps(i,j,k)) - dens(i,j,k) + tau(i,j,k) = sdetg(i,j,k) * (rho(i,j,k)+rho(i,j,k)*eps(i,j,k)) - dens(i,j,k) #else ! cott 2010/03/27: @@ -319,7 +317,7 @@ subroutine Conservative2Primitive(CCTK_ARGUMENTS) 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), vup(i,j,k,1), vup(i,j,k,2), & vup(i,j,k,3), eps(i,j,k), press(i,j,k), w_lorentz(i,j,k), & - uxx, uxy, uxz, uyy, uyz, uzz, det, x(i,j,k), y(i,j,k), & + uxx, uxy, uxz, uyy, uyz, uzz, sdetg(i,j,k), x(i,j,k), y(i,j,k), & z(i,j,k), r(i,j,k),GRHydro_rho_min, GRHydro_reflevel, GRHydro_C2P_failed(i,j,k)) #endif @@ -393,7 +391,7 @@ a single point. subroutine Con2Prim_pt(cctk_iteration,ii,jj,kk,& handle, dens, sx, sy, sz, tau, rho, velx, vely, & velz, epsilon, press, w_lorentz, uxx, uxy, uxz, uyy, & - uyz, uzz, det, x, y, z, r, epsnegative, GRHydro_rho_min, pmin, epsmin, & + uyz, uzz, sdetg, x, y, z, r, epsnegative, GRHydro_rho_min, pmin, epsmin, & GRHydro_reflevel, GRHydro_C2P_failed) implicit none @@ -401,7 +399,7 @@ subroutine Con2Prim_pt(cctk_iteration,ii,jj,kk,& DECLARE_CCTK_PARAMETERS CCTK_REAL dens, sx, sy, sz, tau, rho, velx, vely, velz, epsilon, & - press, uxx, uxy, uxz, uyy, uyz, uzz, det, w_lorentz, x, & + press, uxx, uxy, uxz, uyy, uyz, uzz, sdetg, isdetg, w_lorentz, x, & y, z, r, GRHydro_rho_min CCTK_REAL s2, f, df, vlowx, vlowy, vlowz @@ -424,11 +422,12 @@ subroutine Con2Prim_pt(cctk_iteration,ii,jj,kk,& !!$ Undensitize the variables - udens = dens /sqrt(det) - usx = sx /sqrt(det) - usy = sy /sqrt(det) - usz = sz /sqrt(det) - utau = tau /sqrt(det) + isdetg = 1.0d0/sdetg + udens = dens * isdetg + usx = sx * isdetg + usy = sy * isdetg + usz = sz * isdetg + utau = tau * isdetg s2 = usx*usx*uxx + usy*usy*uyy + usz*usz*uzz + 2.*usx*usy*uxy + & 2.*usx*usz*uxz + 2.*usy*usz*uyz @@ -512,7 +511,7 @@ subroutine Con2Prim_pt(cctk_iteration,ii,jj,kk,& ! for safety, let's set the point to atmosphere SET_ATMO_MIN(rho, GRHydro_rho_min, r) udens = rho - dens = sqrt(det) * rho + dens = sdetg * rho pnew = pmin epsilon = epsmin ! w_lorentz=1, so the expression for utau reduces to: @@ -685,7 +684,7 @@ subroutine Con2Prim_pt(cctk_iteration,ii,jj,kk,& IF_BELOW_ATMO(rho, GRHydro_rho_min, GRHydro_atmo_tolerance, r) then SET_ATMO_MIN(rho, GRHydro_rho_min, r) !GRHydro_rho_min udens = rho - dens = sqrt(det) * rho + dens = isdetg * rho ! epsilon = (sqrt( (utau + pnew + udens)**2) - pnew - udens)/udens epsilon = epsmin ! w_lorentz=1, so the expression for utau reduces to: @@ -874,11 +873,11 @@ subroutine Conservative2PrimitiveBounds(CCTK_ARGUMENTS) epsnegative = .false. - avg_detl = SPATIAL_DETERMINANT(gxxl,gxyl,gxzl,gyyl, gyzl,gzzl) - avg_detr = SPATIAL_DETERMINANT(gxxr,gxyr,gxzr,gyyr, gyzr,gzzr) - call UpperMetric(uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl,& + avg_detl = sqrt(SPATIAL_DETERMINANT(gxxl,gxyl,gxzl,gyyl, gyzl,gzzl)) + avg_detr = sqrt(SPATIAL_DETERMINANT(gxxr,gxyr,gxzr,gyyr, gyzr,gzzr)) + call UpperMetric(uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl*avg_detl,& gxxl,gxyl,gxzl,gyyl,gyzl,gzzl) - call UpperMetric(uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr,& + call UpperMetric(uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr*avg_detr,& gxxr,gxyr,gxzr,gyyr,gyzr,gzzr) if (evolve_tracer .ne. 0) then @@ -1035,7 +1034,7 @@ subroutine Conservative2PrimitivePolytype(CCTK_ARGUMENTS) DECLARE_CCTK_PARAMETERS integer :: i, j, k, itracer, nx, ny, nz - CCTK_REAL :: uxx, uxy, uxz, uyy, uyz, uzz, det + CCTK_REAL :: uxx, uxy, uxz, uyy, uyz, uzz CCTK_REAL :: local_min_tracer ! character(len=400) :: warnline @@ -1080,7 +1079,7 @@ subroutine Conservative2PrimitivePolytype(CCTK_ARGUMENTS) !!$ do j = GRHydro_stencil + 1, ny - GRHydro_stencil !!$ do i = GRHydro_stencil + 1, nx - GRHydro_stencil !$OMP PARALLEL DO PRIVATE(i,j,k,itracer,& - !$OMP uxx, uxy, uxz, uyy, uyz, uzz, det) + !$OMP uxx, uxy, uxz, uyy, uyz, uzz) do k = 1, nz do j = 1, ny do i = 1, nx @@ -1089,9 +1088,7 @@ subroutine Conservative2PrimitivePolytype(CCTK_ARGUMENTS) if ((atmosphere_mask(i,j,k) .ne. 0) .or. & GRHydro_enable_internal_excision /= 0 .and. (hydro_excision_mask(i,j,k) .ne. 0)) cycle - det = SPATIAL_DETERMINANT(g11(i,j,k),g12(i,j,k),g13(i,j,k),\ - g22(i,j,k),g23(i,j,k),g33(i,j,k)) - call UpperMetric(uxx,uxy,uxz,uyy,uyz,uzz,det,& + call UpperMetric(uxx,uxy,uxz,uyy,uyz,uzz,sdetg(i,j,k)*sdetg(i,j,k),& g11(i,j,k),g12(i,j,k),g13(i,j,k),g22(i,j,k),& g23(i,j,k),g33(i,j,k)) @@ -1117,7 +1114,7 @@ subroutine Conservative2PrimitivePolytype(CCTK_ARGUMENTS) scon(i,j,k,1),scon(i,j,k,2), & scon(i,j,k,3),tau(i,j,k),rho(i,j,k),vup(i,j,k,1),vup(i,j,k,2), & vup(i,j,k,3),eps(i,j,k),press(i,j,k),w_lorentz(i,j,k), & - uxx,uxy,uxz,uyy,uyz,uzz,det,x(i,j,k),y(i,j,k), & + uxx,uxy,uxz,uyy,uyz,uzz,sdetg(i,j,k),x(i,j,k),y(i,j,k), & z(i,j,k),r(i,j,k),GRHydro_rho_min, GRHydro_reflevel, GRHydro_C2P_failed(i,j,k)) end do @@ -1149,20 +1146,20 @@ a single point. subroutine Con2Prim_ptPolytype(handle, dens, sx, sy, sz, tau, rho, & velx, vely, velz, epsilon, press, w_lorentz, uxx, uxy, uxz, uyy, & - uyz, uzz, det, x, y, z, r, GRHydro_rho_min, GRHydro_reflevel, GRHydro_C2P_failed) + uyz, uzz, sdetg, x, y, z, r, GRHydro_rho_min, GRHydro_reflevel, GRHydro_C2P_failed) implicit none DECLARE_CCTK_PARAMETERS CCTK_REAL dens, sx, sy, sz, tau, rho, velx, vely, velz, epsilon, & - press, uxx, uxy, uxz, uyy, uyz, uzz, det, w_lorentz, x, & + press, uxx, uxy, uxz, uyy, uyz, uzz, sdetg, isdetg, w_lorentz, x, & y, z, r, GRHydro_rho_min CCTK_REAL s2, f, df, vlowx, vlowy, vlowz CCTK_INT count, handle, GRHydro_reflevel CCTK_REAL udens, usx, usy, usz, rhoold, rhonew, & - enthalpy, denthalpy, sqrtdet, invsqrtdet, invfac, GRHydro_C2P_failed, dummy1, dummy2 + enthalpy, denthalpy, invfac, GRHydro_C2P_failed, dummy1, dummy2 character(len=200) warnline ! begin EOS Omni vars @@ -1175,12 +1172,11 @@ subroutine Con2Prim_ptPolytype(handle, dens, sx, sy, sz, tau, rho, & !!$ Undensitize the variables - sqrtdet = sqrt(det) - invsqrtdet = 1.d0/sqrtdet - udens = dens*invsqrtdet - usx = sx*invsqrtdet - usy = sy*invsqrtdet - usz = sz*invsqrtdet + isdetg = 1.0d0/sdetg + udens = dens*isdetg + usx = sx*isdetg + usy = sy*isdetg + usz = sz*isdetg s2 = usx*usx*uxx + usy*usy*uyy + usz*usz*uzz + 2.d0*usx*usy*uxy + & 2.d0*usx*usz*uxz + 2.d0*usy*usz*uyz @@ -1231,7 +1227,7 @@ subroutine Con2Prim_ptPolytype(handle, dens, sx, sy, sz, tau, rho, & ! for safety, let's set the point to atmosphere SET_ATMO_MIN(rhonew, GRHydro_rho_min, r) udens = rhonew - dens = sqrtdet * rhonew + dens = sdetg * rhonew call EOS_Omni_press(handle,keytemp,GRHydro_eos_rf_prec,n,& rhonew,1.0d0,xtemp,xye,press,keyerr,anyerr) @@ -1286,7 +1282,7 @@ subroutine Con2Prim_ptPolytype(handle, dens, sx, sy, sz, tau, rho, & enthalpy = 1.0d0 + epsilon + press / rhonew w_lorentz = sqrt(1.d0 + s2 / ( (udens*enthalpy)**2 )) - tau = sqrtdet * ((rho * enthalpy) * w_lorentz**2 - press) - dens + tau = sdetg * ((rho * enthalpy) * w_lorentz**2 - press) - dens f = rhonew * w_lorentz - udens @@ -1302,7 +1298,7 @@ subroutine Con2Prim_ptPolytype(handle, dens, sx, sy, sz, tau, rho, & ! before 2009/10/07 dens was not reset and was used with its negative ! value below; this has since been changed, but the change produces ! significant changes in the results - dens = sqrtdet * rhonew + dens = sdetg * rhonew call EOS_Omni_press(handle,keytemp,GRHydro_eos_rf_prec,n,& rhonew,1.0d0,xtemp,xye,press,keyerr,anyerr) @@ -1338,7 +1334,7 @@ subroutine Con2Prim_ptPolytype(handle, dens, sx, sy, sz, tau, rho, & call EOS_Omni_EpsFromPress(handle,keytemp,GRHydro_eos_rf_prec,n,& rhonew,xeps,xtemp,xye,press,epsilon,keyerr,anyerr) - tau = sqrtdet * ((rho * enthalpy) * w_lorentz**2 - press) - dens + tau = sdetg * ((rho * enthalpy) * w_lorentz**2 - press) - dens if (epsilon .lt. 0.0d0) then GRHydro_C2P_failed = 1 @@ -1355,7 +1351,7 @@ subroutine Con2Prim_ptPolytype(handle, dens, sx, sy, sz, tau, rho, & !$OMP END CRITICAL rho = GRHydro_rho_min - dens = sqrtdet * rho + dens = sdetg * rho call EOS_Omni_press(handle,keytemp,GRHydro_eos_rf_prec,n,& rhonew,1.d0,xtemp,xye,press,keyerr,anyerr) @@ -1363,7 +1359,7 @@ subroutine Con2Prim_ptPolytype(handle, dens, sx, sy, sz, tau, rho, & rhonew,xeps,xtemp,xye,press,epsilon,keyerr,anyerr) ! w_lorentz=1, so the expression for tau reduces to: - tau = sqrtdet * (rho+rho*epsilon) - dens + tau = sdetg * (rho+rho*epsilon) - dens usx = 0.d0 usy = 0.d0 usz = 0.d0 @@ -1469,11 +1465,11 @@ subroutine Con2PrimBoundsPolytype(CCTK_ARGUMENTS) gyzr = 0.5d0 * (g23(i,j,k) + g23(i+xoffset,j+yoffset,k+zoffset)) gzzr = 0.5d0 * (g33(i,j,k) + g33(i+xoffset,j+yoffset,k+zoffset)) - avg_detl = SPATIAL_DETERMINANT(gxxl,gxyl,gxzl,gyyl, gyzl,gzzl) - avg_detr = SPATIAL_DETERMINANT(gxxr,gxyr,gxzr,gyyr, gyzr,gzzr) - call UpperMetric(uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl,& + avg_detl = sqrt(SPATIAL_DETERMINANT(gxxl,gxyl,gxzl,gyyl, gyzl,gzzl)) + avg_detr = sqrt(SPATIAL_DETERMINANT(gxxr,gxyr,gxzr,gyyr, gyzr,gzzr)) + call UpperMetric(uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl*avg_detl,& gxxl,gxyl,gxzl,gyyl,gyzl,gzzl) - call UpperMetric(uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr,& + call UpperMetric(uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr*avg_detr,& gxxr,gxyr,gxzr,gyyr,gyzr,gzzr) call Con2Prim_ptPolytype(GRHydro_eos_handle, densminus(i,j,k),& @@ -1484,7 +1480,8 @@ subroutine Con2PrimBoundsPolytype(CCTK_ARGUMENTS) uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl,& x(i,j,k)-0.5d0*CCTK_DELTA_SPACE(1),& y(i,j,k)-0.5d0*CCTK_DELTA_SPACE(2), & - z(i,j,k)-0.5d0*CCTK_DELTA_SPACE(3),r(i,j,k),GRHydro_rho_min, GRHydro_reflevel, GRHydro_C2P_failed(i,j,k)) + z(i,j,k)-0.5d0*CCTK_DELTA_SPACE(3),& + r(i,j,k),GRHydro_rho_min, GRHydro_reflevel, GRHydro_C2P_failed(i,j,k)) call Con2Prim_ptPolytype(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),& @@ -1582,11 +1579,11 @@ subroutine Con2PrimBoundsTracer(CCTK_ARGUMENTS) gyzr = 0.5d0 * (g23(i,j,k) + g23(i+xoffset,j+yoffset,k+zoffset)) gzzr = 0.5d0 * (g33(i,j,k) + g33(i+xoffset,j+yoffset,k+zoffset)) - avg_detl = SPATIAL_DETERMINANT(gxxl,gxyl,gxzl,gyyl, gyzl,gzzl) - avg_detr = SPATIAL_DETERMINANT(gxxr,gxyr,gxzr,gyyr, gyzr,gzzr) - call UpperMetric(uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl,& + avg_detl = sqrt(SPATIAL_DETERMINANT(gxxl,gxyl,gxzl,gyyl, gyzl,gzzl)) + avg_detr = sqrt(SPATIAL_DETERMINANT(gxxr,gxyr,gxzr,gyyr, gyzr,gzzr)) + call UpperMetric(uxxl,uxyl,uxzl,uyyl,uyzl,uzzl,avg_detl*avg_detl,& gxxl,gxyl,gxzl,gyyl,gyzl,gzzl) - call UpperMetric(uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr,& + call UpperMetric(uxxr,uxyr,uxzr,uyyr,uyzr,uzzr,avg_detr*avg_detr,& gxxr,gxyr,gxzr,gyyr,gyzr,gzzr) do itracer=1,number_of_tracers @@ -1681,15 +1678,15 @@ end subroutine Con2Prim_ptTracer @@*/ -subroutine Con2Prim_ptBoundsTracer(cons_tracer, tracer, rho, one_over_w_lorentz, det) +subroutine Con2Prim_ptBoundsTracer(cons_tracer, tracer, rho, one_over_w_lorentz, sdet) implicit none DECLARE_CCTK_PARAMETERS - CCTK_REAL cons_tracer, tracer, rho, one_over_w_lorentz, det + CCTK_REAL cons_tracer, tracer, rho, one_over_w_lorentz, sdet - tracer = cons_tracer / sqrt(det) / rho * one_over_w_lorentz + tracer = cons_tracer / sdet / rho * one_over_w_lorentz end subroutine Con2Prim_ptBoundsTracer |