From 91270e519a2413c00e39aee4688765bd0e4d9ab1 Mon Sep 17 00:00:00 2001 From: rhaas Date: Tue, 4 Sep 2012 12:28:40 +0000 Subject: Changes to make GRHydro MHD work with nuclear/hot equation of state: 1) Switch to EOSOmni pointwise C2P routine and modify where necessary. 2) Modify Con2PrimM.F90 to allow for the evolution of temperature and adjust the wrapper routine. 3) Create EigenProblemM_hot pointwise routine and call that from HLLEM.F90 when temperature is evolved. Additionally adjust HLLEM where necessary. 4) Adjust InterfacesM.h to incorporate the newly created functions. 5) Fix a loop problem (not taking into account constraint transport) in PPM reconstruction of Y_e 6) Introduce Prim2ConM_hot and call this pointwise routine from Prim2ConM.F90 when temperature is evolved. Additionally also make this routine available to initial data routine in GRHydro_InitData 7) Adjust loops in GRHydro_PoloidalMagFieldM.F90 to not set boundary points it cannot set but instead call boundary group afterwards! Pay attention as this will not work with boundary conditions set to "none" in MHD case anymore but is the correct thing to do. 8) Allow StarMapper to extend HydroBase::initial_hydro = "starmapper". 9) Smaller fixes. From: Philipp Moesta git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinInitialData/GRHydro_InitData/trunk@153 ac85fae7-cede-4708-beff-ae01c7fa1c26 --- src/GRHydro_PoloidalMagFieldM.F90 | 28 ++++++++++++++++++++++------ 1 file changed, 22 insertions(+), 6 deletions(-) (limited to 'src') diff --git a/src/GRHydro_PoloidalMagFieldM.F90 b/src/GRHydro_PoloidalMagFieldM.F90 index 47ce06c..16c6eef 100644 --- a/src/GRHydro_PoloidalMagFieldM.F90 +++ b/src/GRHydro_PoloidalMagFieldM.F90 @@ -75,9 +75,11 @@ subroutine GRHydro_PoloidalMagFieldM(CCTK_ARGUMENTS) dy = CCTK_DELTA_SPACE(2) dz = CCTK_DELTA_SPACE(3) - do i=1,nx - do j=1,ny - do k=1,nz + write(*,*)'GRHydro_InitData: Setting up initial poloidal magnetic field' + + do i=2,nx-1 + do j=2,ny-1 + do k=2,nz-1 rhofac = 1.0d0-rho(i,j,k)/poloidal_rho_max delPcut = press(i,j,k)-poloidal_P_cut @@ -121,6 +123,10 @@ subroutine GRHydro_PoloidalMagFieldM(CCTK_ARGUMENTS) Bvecy(i,j,k) = - Ax_dz/sdet Bvecz(i,j,k) = (Ax_dy-Ay_dx)/sdet + !Bvecx(i,j,k) = 0.0d0 + !Bvecy(i,j,k) = 0.0d0 + !Bvecz(i,j,k) = 0.00000001/sdet + if (CCTK_EQUALS(GRHydro_eos_type,"Polytype")) then call Prim2ConPolyM(GRHydro_eos_handle,gxx(i,j,k),gxy(i,j,k),& gxz(i,j,k),gyy(i,j,k),gyz(i,j,k),gzz(i,j,k),& @@ -130,13 +136,21 @@ subroutine GRHydro_PoloidalMagFieldM(CCTK_ARGUMENTS) eps(i,j,k),press(i,j,k),Bvecx(i,j,k),Bvecy(i,j,k),Bvecz(i,j,k),& w_lorentz(i,j,k)) else - call Prim2ConGenM(GRHydro_eos_handle,gxx(i,j,k),gxy(i,j,k),& + call Prim2ConGenM_hot(GRHydro_eos_handle,GRHydro_reflevel,& + i,j,k,x(i,j,k),y(i,j,k),z(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),& det, dens(i,j,k),sx(i,j,k),sy(i,j,k),sz(i,j,k),& tau(i,j,k),Bconsx(i,j,k),Bconsy(i,j,k),Bconsz(i,j,k),rho(i,j,k),& velx(i,j,k),vely(i,j,k),velz(i,j,k),& eps(i,j,k),press(i,j,k),Bvecx(i,j,k),Bvecy(i,j,k),Bvecz(i,j,k),& - w_lorentz(i,j,k)) + w_lorentz(i,j,k),temperature(i,j,k),y_e(i,j,k)) +! call Prim2ConGenM(GRHydro_eos_handle,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, dens(i,j,k),sx(i,j,k),sy(i,j,k),sz(i,j,k),& +! tau(i,j,k),Bconsx(i,j,k),Bconsy(i,j,k),Bconsz(i,j,k),rho(i,j,k),& +! velx(i,j,k),vely(i,j,k),velz(i,j,k),& +! eps(i,j,k),press(i,j,k),Bvecx(i,j,k),Bvecy(i,j,k),Bvecz(i,j,k),& +! w_lorentz(i,j,k)) end if enddo @@ -147,7 +161,9 @@ subroutine GRHydro_PoloidalMagFieldM(CCTK_ARGUMENTS) srhs = 0.d0 taurhs = 0.d0 Bconsrhs = 0.d0 - + if (clean_divergence .ne. 0) then + psidcrhs =0.0 + endif return end subroutine GRHydro_PoloidalMagFieldM -- cgit v1.2.3