diff options
Diffstat (limited to 'src/GRHydro_Con2PrimM.F90')
| -rw-r--r-- | src/GRHydro_Con2PrimM.F90 | 51 |
1 files changed, 8 insertions, 43 deletions
diff --git a/src/GRHydro_Con2PrimM.F90 b/src/GRHydro_Con2PrimM.F90 index 79000ef..8a5c93b 100644 --- a/src/GRHydro_Con2PrimM.F90 +++ b/src/GRHydro_Con2PrimM.F90 @@ -164,10 +164,6 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) !do not compute if in atmosphere region! if (atmosphere_mask(i,j,k) .gt. 0) cycle - !do not compute if in atmosphere or in excised region - !if ((atmosphere_mask(i,j,k) .ne. 0) .or. & - ! (hydro_excision_mask(i,j,k) .ne. 0)) cycle - epsnegative = 0 sdet = sdetg(i,j,k) @@ -192,43 +188,6 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) endif - ! 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), sdet*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 - w_lorentz(i,j,k) = 1.d0 - - if(evolve_temper.ne.0) then - ! set hot atmosphere values - temperature(i,j,k) = grhydro_hot_atmo_temp - y_e(i,j,k) = grhydro_hot_atmo_Y_e - y_e_con(i,j,k) = y_e(i,j,k) * dens(i,j,k) - keytemp = 1 - call EOS_Omni_press(GRHydro_eos_handle,keytemp,GRHydro_eos_rf_prec,n,& - rho(i,j,k),eps(i,j,k),temperature(i,j,k),y_e(i,j,k),& - press(i,j,k),keyerr,anyerr) - else - keytemp = 0 - call EOS_Omni_press(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& - rho(i,j,k),eps(i,j,k),xtemp,xye,press(i,j,k),keyerr,anyerr) - call EOS_Omni_EpsFromPress(GRHydro_polytrope_handle,keytemp,GRHydro_eos_rf_prec,n,& - rho(i,j,k),eps(i,j,k),xtemp,xye,press(i,j,k),eps(i,j,k),keyerr,anyerr) - endif - - !!$ 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 - - cycle - - endif - if(evolve_Y_e.ne.0) then Y_e(i,j,k) = max(min(Y_e_con(i,j,k) / dens(i,j,k),GRHydro_Y_e_max),& GRHydro_Y_e_min) @@ -239,8 +198,8 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) 2.0*(g12(i,j,k)*Bprim(i,j,k,1)*Bprim(i,j,k,2)+g13(i,j,k)*Bprim(i,j,k,1)*Bprim(i,j,k,3)+ & g23(i,j,k)*Bprim(i,j,k,2)*Bprim(i,j,k,3)) - - if ( dens(i,j,k) .le. sdet*GRHydro_rho_min*(1.d0+GRHydro_atmo_tolerance) ) then + + if ( (GRHydro_enable_internal_excision /= 0 .and. (hydro_excision_mask(i,j,k) .gt. 0)) .or. (dens(i,j,k) .le. sdet*GRHydro_rho_min*(1.d0+GRHydro_atmo_tolerance))) then !call CCTK_WARN(1,"Con2Prim: Resetting to atmosphere") !write(warnline,"(3i5,1P10E15.6)") i,j,k,x(i,j,k),y(i,j,k),z(i,j,k) @@ -256,6 +215,10 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) vup(i,j,k,:) = 0.d0 w_lorentz(i,j,k) = 1.d0 + Bprim(i,j,k,:) = 0.0d0 + Bcons(i,j,k,:) = 0.0d0 + b2 = 0 + if(evolve_temper.ne.0) then ! set hot atmosphere values temperature(i,j,k) = grhydro_hot_atmo_temp @@ -294,6 +257,8 @@ subroutine Conservative2PrimitiveM(CCTK_ARGUMENTS) end if + + if(evolve_temper.eq.0) then |