[Commits] [svn:einsteintoolkit] GRHydro/trunk/ (Rev. 196)
bcmsma at astro.rit.edu
bcmsma at astro.rit.edu
Tue Dec 21 23:35:44 CST 2010
User: bmundim
Date: 2010/12/21 11:35 PM
Added:
/trunk/src/
GRHydro_DivergenceClean.F90
Removed:
/trunk/src/
GRHydro_ParamCheckM.F90
Modified:
/trunk/
interface.ccl, param.ccl, schedule.ccl
/trunk/src/
GRHydro_CalcUpdateM.F90, GRHydro_HLLEM.F90, GRHydro_HLLE_EOSGM.F90, GRHydro_ParamCheck.F90, GRHydro_RiemannSolveM.F90, GRHydro_SourceM.F90, make.code.defn
Log:
RIT MHD development: an update.
Several bugfixes, including properly isolating divergence cleaning
code and initialization of the "psidc" gridfunction in
GRHydro_DivergenceClean.F90, as well as properly initializing rhoenth
in GRHydro_HLLEM (thanks to Roland for catching that). Minor
changes to scheduler to properly call MHD and non-MHD routines in
parallel for Riemann solving. Remove GRHydro_ParamCheckM.F90 and
add GRHydro_DivergenceClean.F90.
File Changes:
Directory: /trunk/src/
======================
File [modified]: GRHydro_CalcUpdateM.F90
Delta lines: +10 -7
===================================================================
--- trunk/src/GRHydro_CalcUpdateM.F90 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/src/GRHydro_CalcUpdateM.F90 2010-12-22 05:35:43 UTC (rev 196)
@@ -89,10 +89,11 @@
Bvecrhs(i,j,k,3) = Bvecrhs(i,j,k,3) + &
(alp_l * Bveczflux(i-xoffset,j-yoffset,k-zoffset) - &
alp_r * Bveczflux(i,j,k)) * idx
- psidcrhs(i,j,k) = psidcrhs(i,j,k) + &
- (alp_l * psidcflux(i-xoffset,j-yoffset,k-zoffset) - &
- alp_r * psidcflux(i,j,k)) * idx
-
+ if(clean_divergence.ne.0) then
+ psidcrhs(i,j,k) = psidcrhs(i,j,k) + &
+ (alp_l * psidcflux(i-xoffset,j-yoffset,k-zoffset) - &
+ alp_r * psidcflux(i,j,k)) * idx
+ endif
if (evolve_tracer .ne. 0) then
do itracer=1,number_of_tracers
@@ -236,9 +237,11 @@
Bvecrhs(i,j,k,3) = Bvecrhs(i,j,k,3) + &
(Bveczflux(i-xoffset,j-yoffset,k-zoffset) - &
Bveczflux(i,j,k)) * idx
- psidcrhs(i,j,k) = psidcrhs(i,j,k) + &
- (psidcflux(i-xoffset,j-yoffset,k-zoffset) - &
- psidcflux(i,j,k)) * idx
+ if(clean_divergence.ne.0) then
+ psidcrhs(i,j,k) = psidcrhs(i,j,k) + &
+ (psidcflux(i-xoffset,j-yoffset,k-zoffset) - &
+ psidcflux(i,j,k)) * idx
+ endif
enddo
enddo
File [added]: GRHydro_DivergenceClean.F90
Delta lines: +40 -0
===================================================================
--- trunk/src/GRHydro_DivergenceClean.F90 (rev 0)
+++ trunk/src/GRHydro_DivergenceClean.F90 2010-12-22 05:35:43 UTC (rev 196)
@@ -0,0 +1,40 @@
+ /*@@
+ @file GRHydro_DivergenceClean.F90
+ @date Nov 24, 2010
+ @author
+ @desc
+ Routines controlling divergence cleaning
+ @enddesc
+ @@*/
+
+#include "cctk.h"
+#include "cctk_Parameters.h"
+#include "cctk_Arguments.h"
+
+ /*@@
+ @routine GRHydro_InitDivergenceClean
+ @date Nov 24, 2010
+ @author Joshua Faber
+ @desc
+ Set Psi=0 initially for use with divergence cleaning
+ @enddesc
+ @calls
+ @calledby
+ @history
+
+ @endhistory
+
+@@*/
+
+subroutine GRHydro_InitDivergenceClean(CCTK_ARGUMENTS)
+
+ implicit none
+
+ DECLARE_CCTK_ARGUMENTS
+
+ psidc=0.0
+
+end subroutine GRHydro_InitDivergenceClean
+
+
+
File [modified]: GRHydro_HLLEM.F90
Delta lines: +6 -1
===================================================================
--- trunk/src/GRHydro_HLLEM.F90 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/src/GRHydro_HLLEM.F90 2010-12-22 05:35:43 UTC (rev 196)
@@ -187,6 +187,9 @@
prim_m(2),prim_m(3),prim_m(4),cons_m(6),cons_m(7),cons_m(8), &
velxlowm,velylowm,velzlowm,Bvecxlowm,Bvecylowm,Bveczlowm, &
bdotvm,b2m,v2m,wm,bxlowm,bylowm,bzlowm)
+
+ rhoenth_p = prim_p(1)*(1.d0+prim_p(5))+prim_p(6)
+ rhoenth_m = prim_m(1)*(1.d0+prim_m(5))+prim_m(6)
ab0p = wp*bdotvp
ab0m = wm*bdotvm
@@ -400,7 +403,9 @@
Bvecxflux(i, j, k) = f1(6)
Bvecyflux(i, j, k) = f1(7)
Bveczflux(i, j, k) = f1(8)
- psidcflux(i,j,k) = psidcf
+ if(clean_divergence.ne.0) then
+ psidcflux(i,j,k) = psidcf
+ endif
if(evolve_Y_e.ne.0) then
if (densflux(i, j, k) > 0.d0) then
File [modified]: GRHydro_HLLE_EOSGM.F90
Delta lines: +4 -1
===================================================================
--- trunk/src/GRHydro_HLLE_EOSGM.F90 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/src/GRHydro_HLLE_EOSGM.F90 2010-12-22 05:35:43 UTC (rev 196)
@@ -438,7 +438,10 @@
Bvecxflux(i, j, k) = f1(6)
Bvecyflux(i, j, k) = f1(7)
Bveczflux(i, j, k) = f1(8)
- psidcflux(i,j,k) = psidcf
+ if(clean_divergence.ne.0) then
+ psidcflux(i,j,k) = psidcf
+ endif
+
if(evolve_Y_e.ne.0) then
if (densflux(i, j, k) > 0.d0) then
File [modified]: GRHydro_ParamCheck.F90
Delta lines: +9 -2
===================================================================
--- trunk/src/GRHydro_ParamCheck.F90 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/src/GRHydro_ParamCheck.F90 2010-12-22 05:35:43 UTC (rev 196)
@@ -87,9 +87,9 @@
end if
if (CCTK_EQUALS(Bvec_evolution_method,"GRHydro")) then
- MHD = 1
+ evolve_MHD = 1
else
- MHD = 0
+ evolve_MHD = 0
endif
if (CCTK_EQUALS(Y_e_evolution_method,"GRHydro")) then
@@ -104,6 +104,13 @@
evolve_temper = 0
endif
+ if (CCTK_EQUALS(riemann_solver,"Roe").and.CCTK_EQUALS(Bvec_evolution_method,"GRHydro")) then
+ call CCTK_PARAMWARN("Roe solver is not implemented yet for MHD")
+ end if
+ if (CCTK_EQUALS(riemann_solver,"Marquina").and.CCTK_EQUALS(Bvec_evolution_method,"GRHydro")) then
+ call CCTK_PARAMWARN("Marquina solver is not implemented yet for MHD")
+ end if
+
end subroutine GRHydro_ParamCheck
File [removed]: GRHydro_ParamCheckM.F90
Delta lines: +0 -50
===================================================================
--- trunk/src/GRHydro_ParamCheckM.F90 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/src/GRHydro_ParamCheckM.F90 2010-12-22 05:35:43 UTC (rev 196)
@@ -1,50 +0,0 @@
- /*@@
- @file GRHydro_ParamCheckM.F90
- @date Sep 23, 2010
- @author
- @desc
- Parameter checking routine.
- @enddesc
- @@*/
-
-#include "cctk.h"
-#include "cctk_Parameters.h"
-#include "cctk_Arguments.h"
-#include "cctk_Functions.h"
-
- /*@@
- @routine GRHydro_ParamCheckM
- @date Sep 23, 2010
- @author Joshua Faber, Scott Noble, Bruno Mundim
- @desc
- Checks the parameters.
- @enddesc
- @calls
- @calledby
- @history
-
- @endhistory
-
-@@*/
-
-subroutine GRHydro_ParamCheckM(CCTK_ARGUMENTS)
-
- implicit none
-
- DECLARE_CCTK_ARGUMENTS
- DECLARE_CCTK_PARAMETERS
- DECLARE_CCTK_FUNCTIONS
-
-
-
- if (CCTK_EQUALS(riemann_solver,"Roe").and.CCTK_EQUALS(Bvec_evolution_method,"GRHydro")) then
- call CCTK_PARAMWARN("Roe solver is not implemented yet for MHD")
- end if
-
- if (CCTK_EQUALS(riemann_solver,"Marquina").and.CCTK_EQUALS(Bvec_evolution_method,"GRHydro")) then
- call CCTK_PARAMWARN("Marquina solver is not implemented yet for MHD")
- end if
-
-
-end subroutine GRHydro_ParamCheckM
-
File [modified]: GRHydro_RiemannSolveM.F90
Delta lines: +0 -1
===================================================================
--- trunk/src/GRHydro_RiemannSolveM.F90 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/src/GRHydro_RiemannSolveM.F90 2010-12-22 05:35:43 UTC (rev 196)
@@ -362,7 +362,6 @@
Bvecxflux(i,j,k)= Bvecxflux(i,j,k)+ 0.5d0 * tmp_flux(6)
Bvecyflux(i,j,k)= Bvecyflux(i,j,k)+ 0.5d0 * tmp_flux(7)
Bveczflux(i,j,k)= Bveczflux(i,j,k)+ 0.5d0 * tmp_flux(8)
- psidcflux(i,j,k)= psidcflux(i,j,k)+ 0.5d0 * psidcf
if(clean_divergence.ne.0) then
psidcf = cons_m(6)
File [modified]: GRHydro_SourceM.F90
Delta lines: +32 -24
===================================================================
--- trunk/src/GRHydro_SourceM.F90 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/src/GRHydro_SourceM.F90 2010-12-22 05:35:43 UTC (rev 196)
@@ -9,18 +9,18 @@
! Second order f.d.
-#define DIFF_X_2(q) 0.5d0 * (q(i+1,j,k) - q(i-1,j,k)) * idx
-#define DIFF_Y_2(q) 0.5d0 * (q(i,j+1,k) - q(i,j-1,k)) * idy
-#define DIFF_Z_2(q) 0.5d0 * (q(i,j,k+1) - q(i,j,k-1)) * idz
+#define DIFF_X_2(q) (0.5d0 * (q(i+1,j,k) - q(i-1,j,k)) * idx)
+#define DIFF_Y_2(q) (0.5d0 * (q(i,j+1,k) - q(i,j-1,k)) * idy)
+#define DIFF_Z_2(q) (0.5d0 * (q(i,j,k+1) - q(i,j,k-1)) * idz)
! Fourth order f.d.
-#define DIFF_X_4(q) (-q(i+2,j,k) + 8.d0 * q(i+1,j,k) - 8.d0 * q(i-1,j,k) + \
- q(i-2,j,k)) / 12.d0 * idx
-#define DIFF_Y_4(q) (-q(i,j+2,k) + 8.d0 * q(i,j+1,k) - 8.d0 * q(i,j-1,k) + \
- q(i,j-2,k)) / 12.d0 * idy
-#define DIFF_Z_4(q) (-q(i,j,k+2) + 8.d0 * q(i,j,k+1) - 8.d0 * q(i,j,k-1) + \
- q(i,j,k-2)) / 12.d0 * idz
+#define DIFF_X_4(q) ((-q(i+2,j,k) + 8.d0 * q(i+1,j,k) - 8.d0 * q(i-1,j,k) + \
+ q(i-2,j,k)) / 12.d0 * idx)
+#define DIFF_Y_4(q) ((-q(i,j+2,k) + 8.d0 * q(i,j+1,k) - 8.d0 * q(i,j-1,k) + \
+ q(i,j-2,k)) / 12.d0 * idy)
+#define DIFF_Z_4(q) ((-q(i,j,k+2) + 8.d0 * q(i,j,k+1) - 8.d0 * q(i,j,k-1) + \
+ q(i,j,k-2)) / 12.d0 * idz)
#include "cctk.h"
#include "cctk_Parameters.h"
@@ -103,6 +103,7 @@
densrhs = 0.d0
srhs = 0.d0
taurhs = 0.d0
+ Bvecrhs=0.d0
if (evolve_tracer .ne. 0) then
cons_tracerrhs = 0.d0
@@ -113,7 +114,6 @@
endif
if (clean_divergence .ne. 0) then
- Bvecrhs=0.d0
psidcrhs=0.d0
endif
@@ -294,19 +294,23 @@
dx_alp = DIFF_X_2(alp)
dy_alp = DIFF_Y_2(alp)
dz_alp = DIFF_Z_2(alp)
- dx_psidc = DIFF_X_2(psidc)
- dy_psidc = DIFF_Y_2(psidc)
- dz_psidc = DIFF_Z_2(psidc)
+ if(clean_divergence.ne.0) then
+ dx_psidc = DIFF_X_2(psidc)
+ dy_psidc = DIFF_Y_2(psidc)
+ dz_psidc = DIFF_Z_2(psidc)
+ endif
+
else
dx_alp = DIFF_X_4(alp)
dy_alp = DIFF_Y_4(alp)
dz_alp = DIFF_Z_4(alp)
- dx_psidc = DIFF_X_4(psidc)
- dy_psidc = DIFF_Y_4(psidc)
- dz_psidc = DIFF_Z_4(psidc)
-
+ if(clean_divergence.ne.0) then
+ dx_psidc = DIFF_X_4(psidc)
+ dy_psidc = DIFF_Y_4(psidc)
+ dz_psidc = DIFF_Z_4(psidc)
+ endif
end if
if (local_spatial_order .eq. 2) then
@@ -432,19 +436,23 @@
rhohstarW2*invalp*(vlowx*dz_betax + vlowy*dz_betay + vlowz*dz_betaz) -&
bt*(bxlow*dz_betax + bylow*dz_betay + bzlow*dz_betaz)
- bvcx_source = uxx*dx_psidc+uxy*dy_psidc+uxz*dz_psidc
- bvcy_source = uxy*dx_psidc+uyy*dy_psidc+uyz*dz_psidc
- bvcz_source = uxz*dx_psidc+uyz*dy_psidc+uzz*dz_psidc
+ if(clean_divergence.ne.0) then
+ bvcx_source = uxx*dx_psidc+uxy*dy_psidc+uxz*dz_psidc
+ bvcy_source = uxy*dx_psidc+uyy*dy_psidc+uyz*dz_psidc
+ bvcz_source = uxz*dx_psidc+uyz*dy_psidc+uzz*dz_psidc
+ endif
densrhs(i,j,k) = 0.d0
srhs(i,j,k,1) = alp(i,j,k)*sqrtdet*sx_source
srhs(i,j,k,2) = alp(i,j,k)*sqrtdet*sy_source
srhs(i,j,k,3) = alp(i,j,k)*sqrtdet*sz_source
taurhs(i,j,k) = alp(i,j,k)*sqrtdet*tau_source
- Bvecrhs(i,j,k,1) = alp(i,j,k)*sqrtdet*bvcx_source
- Bvecrhs(i,j,k,2) = alp(i,j,k)*sqrtdet*bvcy_source
- Bvecrhs(i,j,k,3) = alp(i,j,k)*sqrtdet*bvcz_source
- psidcrhs(i,j,k) = -1.d0*ch_dc*ch_dc/cp_dc/cp_dc*psidc(i,j,k)
+ if(clean_divergence.ne.0) then
+ Bvecrhs(i,j,k,1) = alp(i,j,k)*sqrtdet*bvcx_source
+ Bvecrhs(i,j,k,2) = alp(i,j,k)*sqrtdet*bvcy_source
+ Bvecrhs(i,j,k,3) = alp(i,j,k)*sqrtdet*bvcz_source
+ psidcrhs(i,j,k) = -1.d0*ch_dc*ch_dc/cp_dc/cp_dc*psidc(i,j,k)
+ endif
enddo
enddo
File [modified]: make.code.defn
Delta lines: +1 -1
===================================================================
--- trunk/src/make.code.defn 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/src/make.code.defn 2010-12-22 05:35:43 UTC (rev 196)
@@ -7,6 +7,7 @@
GRHydro_Boundaries.F90 \
GRHydro_CalcUpdate.F90 \
GRHydro_Con2Prim.F90 \
+ GRHydro_DivergenceClean.F90 \
GRHydro_Eigenproblem.F90 \
GRHydro_Eigenproblem_Marquina.F90 \
GRHydro_ENOReconstruct.F90 \
@@ -55,7 +56,6 @@
GRHydro_HLLE_EOSGM.F90 \
GRHydro_HLLEM.F90 \
GRHydro_PPMM.F90 \
- GRHydro_ParamCheckM.F90 \
GRHydro_Prim2ConM.F90 \
GRHydro_RiemannSolveM.F90 \
GRHydro_ReconstructM.F90 \
Directory: /trunk/
==================
File [modified]: interface.ccl
Delta lines: +1 -1
===================================================================
--- trunk/interface.ccl 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/interface.ccl 2010-12-22 05:35:43 UTC (rev 196)
@@ -392,7 +392,7 @@
int evolve_temper type = SCALAR tags='checkpoint="no"' "Are we evolving temperature? Set in Paramcheck"
-int MHD type = SCALAR tags='checkpoint="no"' "Are we doing MHD? Set in ParamCheck"
+int evolve_MHD type = SCALAR tags='checkpoint="no"' "Are we doing MHD? Set in ParamCheck"
int GRHydro_reflevel type = SCALAR tags='checkpoint="no"' "Refinement level GRHydro is working on right now"
File [modified]: param.ccl
Delta lines: +1 -1
===================================================================
--- trunk/param.ccl 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/param.ccl 2010-12-22 05:35:43 UTC (rev 196)
@@ -76,7 +76,7 @@
CCTK_INT GRHydro_MaxNumEvolvedVars "The maximum number of evolved variables used by GRHydro" ACCUMULATOR-BASE=MethodofLines::MoL_Num_Evolved_Vars
{
- 5:9 :: "dens scon[3] tau Bvec[3] ye"
+ 5:10 :: "dens scon[3] tau Bvec[3] psidc ye"
} 5
CCTK_INT GRHydro_MaxNumConstrainedVars "The maximum number of constrained variables used by GRHydro" ACCUMULATOR-BASE=MethodofLines::MoL_Num_Constrained_Vars
File [modified]: schedule.ccl
Delta lines: +18 -18
===================================================================
--- trunk/schedule.ccl 2010-12-21 10:58:39 UTC (rev 195)
+++ trunk/schedule.ccl 2010-12-22 05:35:43 UTC (rev 196)
@@ -31,11 +31,11 @@
if(CCTK_Equals(Bvec_evolution_method,"GRHydro"))
{
STORAGE: HydroBase::Bvec[3]
+ if (clean_divergence)
+ {
+ STORAGE:psidc[3]
+ }
}
- if (clean_divergence)
- {
- STORAGE:psidc[3]
- }
}
else
{
@@ -57,13 +57,13 @@
if(CCTK_Equals(Bvec_evolution_method,"GRHydro"))
{
STORAGE: HydroBase::Bvec[2]
+ if (clean_divergence)
+ {
+ STORAGE:psidc[2]
+ }
}
- if (clean_divergence)
- {
- STORAGE:psidc[2]
- }
}
-STORAGE:MHD
+STORAGE:evolve_MHD
STORAGE:evolve_Y_e
STORAGE:evolve_temper
STORAGE:GRHydro_reflevel
@@ -99,6 +99,14 @@
{
} "GRHydro initial data group"
+if(CCTK_Equals(Bvec_evolution_method,"GRHydro") && clean_divergence)
+{
+ schedule GRHydro_InitDivergenceClean IN HydroBase_Initial
+ {
+ LANG: Fortran
+ } "Set psi for divergence cleaning initially to zero"
+}
+
#################################################
### Storage for the extra timelevels for the ###
### use of MoL with Einstein ###
@@ -188,14 +196,6 @@
LANG: Fortran
} "Check parameters"
-if(CCTK_Equals(Bvec_evolution_method,"GRHydro"))
-{
- schedule GRHydro_ParamCheckM AT PARAMCHECK AFTER GRHydro_ParamCheck
- {
- LANG: Fortran
- } "Check parameters - MHD version"
-}
-
######################################
### Standard symmetry registration ###
######################################
@@ -733,7 +733,7 @@
else if (CCTK_Equals(GRHydro_eos_type,"General")) {
if(CCTK_Equals(Bvec_evolution_method,"GRHydro"))
{
- schedule RiemannSolve IN FluxTerms AFTER Reconstruct AS Riemann
+ schedule RiemannSolveM IN FluxTerms AFTER Reconstruct AS Riemann
{
LANG: Fortran
STORAGE: EOS_temps
More information about the Commits
mailing list