diff options
author | schnetter <schnetter@578cdeb0-5ea1-4b81-8215-5a3b8777ee0b> | 2005-12-11 18:49:18 +0000 |
---|---|---|
committer | schnetter <schnetter@578cdeb0-5ea1-4b81-8215-5a3b8777ee0b> | 2005-12-11 18:49:18 +0000 |
commit | d2581020e698490172cee2cdea2e2deb43bf9251 (patch) | |
tree | ab108f3656d349d0682026abbbd16ca4a6d3532b /src/GenericRK.c | |
parent | 7914bcb41d96ee8ab0213c91ebb763a6f3562fb0 (diff) |
Add const and restrict qualifiers to the pointers.
Simplify some complex arithmetic.
Initialise the error variable in the RK45 integrator only after the
last iteration.
git-svn-id: http://svn.cactuscode.org/arrangements/CactusNumerical/MoL/trunk@104 578cdeb0-5ea1-4b81-8215-5a3b8777ee0b
Diffstat (limited to 'src/GenericRK.c')
-rw-r--r-- | src/GenericRK.c | 25 |
1 files changed, 13 insertions, 12 deletions
diff --git a/src/GenericRK.c b/src/GenericRK.c index 7d278b5..196c83e 100644 --- a/src/GenericRK.c +++ b/src/GenericRK.c @@ -80,18 +80,19 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) #ifdef MOLDOESCOMPLEX CCTK_COMPLEX Complex_alpha, Complex_beta, Complex_Delta_Time; - CCTK_COMPLEX *UpdateComplexVar; - CCTK_COMPLEX *RHSComplexVar; - CCTK_COMPLEX *ScratchComplexVar; + CCTK_COMPLEX * restrict UpdateComplexVar; + CCTK_COMPLEX const * restrict RHSComplexVar; + CCTK_COMPLEX * restrict ScratchComplexVar; #endif + static CCTK_INT scratchspace_firstindex = -99; CCTK_INT index, var, scratchstep, alphaindex, scratchindex; CCTK_INT totalsize; CCTK_REAL alpha, beta; - CCTK_REAL *UpdateVar; - CCTK_REAL *RHSVar; - CCTK_REAL *ScratchVar; + CCTK_REAL * restrict UpdateVar; + CCTK_REAL const * restrict RHSVar; + CCTK_REAL * restrict ScratchVar; CCTK_INT arrayscratchlocation; @@ -120,8 +121,8 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) UpdateVar = (CCTK_REAL *)CCTK_VarDataPtrI(cctkGH, 0, EvolvedVariableIndex[var]); - RHSVar = (CCTK_REAL *)CCTK_VarDataPtrI(cctkGH, 0, - RHSVariableIndex[var]); + RHSVar = (CCTK_REAL const *)CCTK_VarDataPtrI(cctkGH, 0, + RHSVariableIndex[var]); /* #define MOLDEBUG 1 */ #ifdef MOLDEBUG printf("In generic RK. Variable %d (%s). RHS %d (%s). beta %g.\n", @@ -267,8 +268,8 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) UpdateComplexVar = (CCTK_COMPLEX *)CCTK_VarDataPtrI(cctkGH, 0, EvolvedComplexVariableIndex[var]); - RHSComplexVar = (CCTK_COMPLEX *)CCTK_VarDataPtrI(cctkGH, 0, - RHSVariableIndex[var]); + RHSComplexVar = (CCTK_COMPLEX const *)CCTK_VarDataPtrI(cctkGH, 0, + RHSVariableIndex[var]); for (index = 0; index < totalsize; index++) { @@ -345,8 +346,8 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) UpdateVar = (CCTK_REAL *)CCTK_VarDataPtrI(cctkGH, 0, EvolvedArrayVariableIndex[var]); - RHSVar = (CCTK_REAL *)CCTK_VarDataPtrI(cctkGH, 0, - RHSArrayVariableIndex[var]); + RHSVar = (CCTK_REAL const *)CCTK_VarDataPtrI(cctkGH, 0, + RHSArrayVariableIndex[var]); groupindex = CCTK_GroupIndexFromVarI(EvolvedArrayVariableIndex[var]); ierr = CCTK_GroupDynamicData(cctkGH, groupindex, |