User: rhaas Date: 2012/05/13 02:01 PM
Modified: /trunk/src/ GRHydro_Con2PrimM.F90
Log: Error checking for C2P_failed > 1, not only equals 1 and check for atmosphere using density after Con2PrimM_pt was performed.
File Changes:
Directory: /trunk/src/ ======================
File [modified]: GRHydro_Con2PrimM.F90 Delta lines: +28 -1 =================================================================== --- trunk/src/GRHydro_Con2PrimM.F90 2012-05-11 22:46:36 UTC (rev 319) +++ trunk/src/GRHydro_Con2PrimM.F90 2012-05-13 19:01:37 UTC (rev 320) @@ -243,6 +243,34 @@ #endif
end if + + if ( rho(i,j,k) .le. GRHydro_rho_min*(1.d0+GRHydro_atmo_tolerance) .or. GRHydro_C2P_failed(i,j,k) .ge. 1) then + + b2=gxx(i,j,k)*Bvec(i,j,k,1)**2+gyy(i,j,k)*Bvec(i,j,k,2)**2+gzz(i,j,k)*Bvec(i,j,k,3)**2+ & + 2.0*(gxy(i,j,k)*Bvec(i,j,k,1)*Bvec(i,j,k,2)+gxz(i,j,k)*Bvec(i,j,k,1)*Bvec(i,j,k,3)+ & + gyz(i,j,k)*Bvec(i,j,k,2)*Bvec(i,j,k,3)) + + dens(i,j,k) = sqrt(det)*GRHydro_rho_min !/(1.d0+GRHydro_atmo_tolerance) + rho(i,j,k) = GRHydro_rho_min + scon(i,j,k,:) = 0.d0 + vel(i,j,k,:) = 0.d0 + w_lorentz(i,j,k) = 1.d0 + xtemp = 0.0d0 + + 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) + + ! w_lorentz=1, so the expression for tau reduces to: + + !!$ 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) = sqrt(det) * (rho(i,j,k)*(1.0+eps(i,j,k)+b2/2.0)) - dens(i,j,k) + + end if + end do end do end do @@ -487,7 +515,6 @@ !$OMP END CRITICAL endif endif - end do end do end do
commits@lists.einsteintoolkit.org