diff options
author | eschnett <eschnett@c83d129a-5a75-4d5a-9c4d-ed3a5855bf45> | 2012-05-11 22:46:36 +0000 |
---|---|---|
committer | eschnett <eschnett@c83d129a-5a75-4d5a-9c4d-ed3a5855bf45> | 2012-05-11 22:46:36 +0000 |
commit | 358d5165056ced200988620b63d158a14c1f0d3f (patch) | |
tree | 084e704a6708b5b548f6c0562bf257adb146ce4a /src/GRHydro_FluxSplit.F90 | |
parent | af6dfaa5b36603dc81d04e3c5a8e054f64c37293 (diff) |
Support real*16 (and real*4) in GRHydro
Support using precisions other than real*8 in GRHydro by removing all
explicit references to double precision functions and constants, and
using type-generic functions and constants instead.
In particular, use "one" and "half" as constants in some places, and
use "abs", "max" etc. instead of "dabs", "dmax" etc.
git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinEvolve/GRHydro/trunk@319 c83d129a-5a75-4d5a-9c4d-ed3a5855bf45
Diffstat (limited to 'src/GRHydro_FluxSplit.F90')
-rw-r--r-- | src/GRHydro_FluxSplit.F90 | 32 |
1 files changed, 17 insertions, 15 deletions
diff --git a/src/GRHydro_FluxSplit.F90 b/src/GRHydro_FluxSplit.F90 index 1c701e8..16777c3 100644 --- a/src/GRHydro_FluxSplit.F90 +++ b/src/GRHydro_FluxSplit.F90 @@ -437,6 +437,8 @@ subroutine GRHydro_SplitFlux_1D(handle, nx, & DECLARE_CCTK_PARAMETERS DECLARE_CCTK_FUNCTIONS + CCTK_REAL, parameter :: half = 0.5d0 + CCTK_INT :: i, nx, handle, ll CCTK_REAL, dimension(nx) :: gxx, gxy, gxz, gyy, gyz, gzz, & u, det, alp, beta, & @@ -493,21 +495,21 @@ subroutine GRHydro_SplitFlux_1D(handle, nx, & !!$ Calculate the maximum eigenvalue and put it here call eigenproblem_leftright(handle, & - 0.5d0 * (rho (i) + rho (i+1)), & - 0.5d0 * (velx1 (i) + velx1 (i+1)), & - 0.5d0 * (vely1 (i) + vely1 (i+1)), & - 0.5d0 * (velz1 (i) + velz1 (i+1)), & - 0.5d0 * (eps (i) + eps (i+1)), & - 0.5d0 * (w_lorentz(i) + w_lorentz(i+1)), & - 0.5d0 * (gxx (i) + gxx (i+1)), & - 0.5d0 * (gxy (i) + gxy (i+1)), & - 0.5d0 * (gxz (i) + gxz (i+1)), & - 0.5d0 * (gyy (i) + gyy (i+1)), & - 0.5d0 * (gyz (i) + gyz (i+1)), & - 0.5d0 * (gzz (i) + gzz (i+1)), & - 0.5d0 * (u (i) + u (i+1)), & - 0.5d0 * (alp (i) + alp (i+1)), & - 0.5d0 * (beta (i) + beta (i+1)), & + half * (rho (i) + rho (i+1)), & + half * (velx1 (i) + velx1 (i+1)), & + half * (vely1 (i) + vely1 (i+1)), & + half * (velz1 (i) + velz1 (i+1)), & + half * (eps (i) + eps (i+1)), & + half * (w_lorentz(i) + w_lorentz(i+1)), & + half * (gxx (i) + gxx (i+1)), & + half * (gxy (i) + gxy (i+1)), & + half * (gxz (i) + gxz (i+1)), & + half * (gyy (i) + gyy (i+1)), & + half * (gyz (i) + gyz (i+1)), & + half * (gzz (i) + gzz (i+1)), & + half * (u (i) + u (i+1)), & + half * (alp (i) + alp (i+1)), & + half * (beta (i) + beta (i+1)), & lambda,& levec,& revec) |