diff options
| -rw-r--r-- | interface.ccl | 17 | ||||
| -rw-r--r-- | schedule.ccl | 6 | ||||
| -rw-r--r-- | src/GRHydro_PoloidalMagFieldM.F90 | 28 |
3 files changed, 44 insertions, 7 deletions
diff --git a/interface.ccl b/interface.ccl index 98fcd8d..ec711c4 100644 --- a/interface.ccl +++ b/interface.ccl @@ -66,6 +66,22 @@ SUBROUTINE Prim2ConGenM(CCTK_INT IN handle, \ CCTK_REAL OUT press, CCTK_REAL IN Bvecx, CCTK_REAL IN Bvecy, \ CCTK_REAL IN Bvecz, CCTK_REAL OUT w_lorentz) +SUBROUTINE Prim2ConGenM_hot(CCTK_INT IN handle, CCTK_INT IN GRHydro_reflevel, CCTK_INT IN i, CCTK_INT IN j, CCTK_INT IN k, \ + CCTK_REAL IN x, CCTK_REAL IN y, CCTK_REAL IN z, \ + CCTK_REAL IN gxx, CCTK_REAL IN gxy, \ + CCTK_REAL IN gxz, CCTK_REAL IN gyy, \ + CCTK_REAL IN gyz, CCTK_REAL IN gzz, \ + CCTK_REAL IN det, CCTK_REAL OUT dens, \ + CCTK_REAL OUT sx, CCTK_REAL OUT sy, \ + CCTK_REAL OUT sz, CCTK_REAL OUT tau, \ + CCTK_REAL OUT Bconsx, CCTK_REAL OUT Bconsy, \ + CCTK_REAL OUT Bconsz, CCTK_REAL IN rho, CCTK_REAL IN velx, \ + CCTK_REAL IN vely, \ + CCTK_REAL IN velz, CCTK_REAL IN epsilon, \ + CCTK_REAL OUT press, CCTK_REAL IN Bvecx, CCTK_REAL IN Bvecy, \ + CCTK_REAL IN Bvecz, CCTK_REAL OUT w_lorentz, \ + CCTK_REAL INOUT temperature, CCTK_REAL IN Y_e) + SUBROUTINE Con2PrimPoly(CCTK_INT IN handle, CCTK_REAL OUT dens, \ CCTK_REAL OUT sx, CCTK_REAL OUT sy, \ @@ -154,6 +170,7 @@ USES FUNCTION Prim2ConPoly USES FUNCTION Prim2ConGen USES FUNCTION Prim2ConPolyM USES FUNCTION Prim2ConGenM +USES FUNCTION Prim2ConGenM_hot USES FUNCTION Con2PrimPoly USES FUNCTION Con2PrimGen USES FUNCTION Con2PrimGenM diff --git a/schedule.ccl b/schedule.ccl index 2b72ab4..4bca6fc 100644 --- a/schedule.ccl +++ b/schedule.ccl @@ -245,9 +245,13 @@ if (CCTK_EQUALS(initial_hydro, "magnetized_bondi_solution")) if (CCTK_EQUALS(initial_Bvec, "poloidalmagfield")) { - SCHEDULE GRHydro_PoloidalMagFieldM IN HydroBase_Initial AFTER TOV_Initial_Data +# SCHEDULE GRHydro_PoloidalMagFieldM AT CCTK_INITIAL AFTER IN HydroBase_Initial AFTER rnsid_init AFTER TOV_Initial_Data after CCCC_StarMapper_InitialData + SCHEDULE GRHydro_PoloidalMagFieldM AT CCTK_POSTINITIAL { LANGUAGE: Fortran } "Set up a poloidal magnetic field. It expects the other fluid variables already to be set, as for example in the TOV solution" + SCHEDULE group HydroBase_Boundaries IN HydroBase_Initial AFTER GRHydro_PoloidalMagFieldM + { + } "Call boundary conditions after magnetic field initial data setup" } 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 |