diff options
Diffstat (limited to 'src/GenericRK.c')
-rw-r--r-- | src/GenericRK.c | 87 |
1 files changed, 59 insertions, 28 deletions
diff --git a/src/GenericRK.c b/src/GenericRK.c index 6c052af..942866b 100644 --- a/src/GenericRK.c +++ b/src/GenericRK.c @@ -108,7 +108,7 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) #ifdef MOLDOESCOMPLEX - Complex_Delta_Time = CCTK_Cmplx((*Original_Delta_Time), 0); + Complex_Delta_Time = CCTK_Cmplx((*Original_Delta_Time) / cctkGH->cctk_timefac, 0); Complex_beta = CCTK_Cmplx(beta, 0); #endif @@ -122,22 +122,26 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) EvolvedVariableIndex[var]); RHSVar = (CCTK_REAL *)CCTK_VarDataPtrI(cctkGH, 0, RHSVariableIndex[var]); -/* #define MOLDEBUG */ +/* #define MOLDEBUG 1 */ #ifdef MOLDEBUG - printf("In generic RK. Variable %d (%s). RHS %d (%s).\n", + printf("In generic RK. Variable %d (%s). RHS %d (%s). beta %g.\n", EvolvedVariableIndex[var], CCTK_VarName(EvolvedVariableIndex[var]), RHSVariableIndex[var], - CCTK_VarName(RHSVariableIndex[var])); + CCTK_VarName(RHSVariableIndex[var]), + beta); #endif for (index = 0; index < totalsize; index++) { - UpdateVar[index] = (*Original_Delta_Time) * beta * RHSVar[index]; + UpdateVar[index] = (*Original_Delta_Time) / cctkGH->cctk_timefac * beta * RHSVar[index]; #ifdef MOLDEBUG - printf("Variable: %d. Index: %d. dt: %f. beta %f. RHS: %f. q: %f.\n", - var, index, (*Original_Delta_Time), beta, RHSVar[index], - UpdateVar[index]); + if (CCTK_EQUALS(verbose,"extreme")) + { + printf("Variable: %d. Index: %d. dt: %f. beta %f. RHS: %f. q: %f.\n", + var, index, (*Original_Delta_Time) / cctkGH->cctk_timefac, beta, RHSVar[index], + UpdateVar[index]); + } #endif } @@ -150,6 +154,15 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) scratchindex = scratchstep - 1; alpha = RKAlphaCoefficients[alphaindex]; +#ifdef MOLDEBUG + printf("In generic RK. Variable %d (%s). RHS %d (%s). step %d. alpha %g.\n", + EvolvedVariableIndex[var], + CCTK_VarName(EvolvedVariableIndex[var]), + RHSVariableIndex[var], + CCTK_VarName(RHSVariableIndex[var]), + scratchstep, + alpha); +#endif if (scratchstep) { @@ -166,9 +179,12 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) + var + MoL_Num_Evolved_Vars * scratchindex); #ifdef MOLDEBUG - printf("Reading from scratch space, initial address %ld index %d\n", - ScratchVar, (var * MoL_Num_Scratch_Levels + - scratchindex) * totalsize); + if (CCTK_EQUALS(verbose,"extreme")) + { + printf("Reading from scratch space, initial address %ld index %d\n", + ScratchVar, (var * MoL_Num_Scratch_Levels + + scratchindex) * totalsize); + } #endif } else @@ -183,10 +199,13 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) { UpdateVar[index] += alpha * ScratchVar[index]; #ifdef MOLDEBUG - printf("Variable: %d. Index: %d. step: %d. " - "alpha: %f. Scratch: %f. q: %f.\n", - var, index, (*MoL_Intermediate_Step), alpha, - ScratchVar[index], UpdateVar[index]); + if (CCTK_EQUALS(verbose,"extreme")) + { + printf("Variable: %d. Index: %d. step: %d. " + "alpha: %f. Scratch: %f. q: %f.\n", + var, index, (*MoL_Intermediate_Step), alpha, + ScratchVar[index], UpdateVar[index]); + } #endif } } @@ -225,9 +244,12 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) { ScratchVar[index] = UpdateVar[index]; #ifdef MOLDEBUG - printf("Variable: %d. Index: %d. step: %d. Scratch: %f.\n", - var, index, (*MoL_Intermediate_Step), ScratchVar[index]); - fflush(stdout); + if (CCTK_EQUALS(verbose,"extreme")) + { + printf("Variable: %d. Index: %d. step: %d. Scratch: %f.\n", + var, index, (*MoL_Intermediate_Step), ScratchVar[index]); + fflush(stdout); + } #endif } } @@ -344,7 +366,7 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) for (index = 0; index < arraytotalsize; index++) { - UpdateVar[index] = (*Original_Delta_Time) * beta * RHSVar[index]; + UpdateVar[index] = (*Original_Delta_Time) / cctkGH->cctk_timefac * beta * RHSVar[index]; } for (scratchstep = 0; @@ -363,9 +385,12 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) (MoL_Max_Evolved_Array_Size+1) + arrayscratchlocation]; #ifdef MOLDEBUG - printf("Reading from scratch space, initial address %ld index %d\n", - ScratchVar, scratchindex*(MoL_Max_Evolved_Array_Size+1) + - arrayscratchlocation); + if (CCTK_EQUALS(verbose,"extreme")) + { + printf("Reading from scratch space, initial address %ld index %d\n", + ScratchVar, scratchindex*(MoL_Max_Evolved_Array_Size+1) + + arrayscratchlocation); + } #endif } else @@ -380,10 +405,13 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) { UpdateVar[index] += alpha * ScratchVar[index]; #ifdef MOLDEBUG - printf("Variable: %d. Index: %d. step: %d. " - "alpha: %f. Scratch: %f. q: %f.\n", - var, index, (*MoL_Intermediate_Step), - alpha, ScratchVar[index], UpdateVar[index]); + if (CCTK_EQUALS(verbose,"extreme")) + { + printf("Variable: %d. Index: %d. step: %d. " + "alpha: %f. Scratch: %f. q: %f.\n", + var, index, (*MoL_Intermediate_Step), + alpha, ScratchVar[index], UpdateVar[index]); + } #endif } } @@ -433,8 +461,11 @@ void MoL_GenericRKAdd(CCTK_ARGUMENTS) { ScratchVar[index] = UpdateVar[index]; #ifdef MOLDEBUG - printf("Variable: %d. Index: %d. step: %d. Scratch: %f.\n", - var, index, (*MoL_Intermediate_Step), ScratchVar[index]); + if (CCTK_EQUALS(verbose,"extreme")) + { + printf("Variable: %d. Index: %d. step: %d. Scratch: %f.\n", + var, index, (*MoL_Intermediate_Step), ScratchVar[index]); + } #endif } arrayscratchlocation += arraytotalsize; |