diff options
Diffstat (limited to 'src')
| -rw-r--r-- | src/GRHydro_PoloidalMagFieldM.F90 | 28 |
1 files changed, 22 insertions, 6 deletions
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 |