aboutsummaryrefslogtreecommitdiff
path: root/src/SetTime.c
diff options
context:
space:
mode:
authorhawke <hawke@578cdeb0-5ea1-4b81-8215-5a3b8777ee0b>2003-06-19 09:54:07 +0000
committerhawke <hawke@578cdeb0-5ea1-4b81-8215-5a3b8777ee0b>2003-06-19 09:54:07 +0000
commitdaa55127b4b58ced31c6bd7248349e059177ded2 (patch)
tree422a89e149e2aac9c784335b80ff8bcf11117d21 /src/SetTime.c
parent49731c4b006be632a84b9a54284682fcb03bd896 (diff)
Patch for the new cctk_time / lev / levoffset etc. variables from Erik Schnetter.
git-svn-id: http://svn.cactuscode.org/arrangements/CactusNumerical/MoL/trunk@16 578cdeb0-5ea1-4b81-8215-5a3b8777ee0b
Diffstat (limited to 'src/SetTime.c')
-rw-r--r--src/SetTime.c28
1 files changed, 14 insertions, 14 deletions
diff --git a/src/SetTime.c b/src/SetTime.c
index 0516e27..7774102 100644
--- a/src/SetTime.c
+++ b/src/SetTime.c
@@ -71,18 +71,18 @@ int MoL_SetTime(CCTK_ARGUMENTS)
CCTK_REAL beta;
*Original_Time = cctkGH->cctk_time;
- *Original_Delta_Time = cctkGH->CCTK_DELTA_TIME;
- cctkGH->cctk_time -= cctkGH->CCTK_DELTA_TIME;
+ *Original_Delta_Time = cctkGH->cctk_delta_time;
+ cctkGH->cctk_time -= cctkGH->cctk_delta_time / cctkGH->cctk_timefac;
if (CCTK_EQUALS(ODE_Method,"ICN"))
{
- cctkGH->CCTK_DELTA_TIME = 0.5*(*Original_Delta_Time);
+ cctkGH->cctk_delta_time = 0.5*(*Original_Delta_Time);
}
else if (CCTK_EQUALS(ODE_Method,"Generic"))
{
beta = RKBetaCoefficients[0];
- cctkGH->CCTK_DELTA_TIME = beta*(*Original_Delta_Time);
+ cctkGH->cctk_delta_time = beta*(*Original_Delta_Time);
}
return 0;
@@ -123,26 +123,26 @@ int MoL_ResetTime(CCTK_ARGUMENTS)
if (*MoL_Intermediate_Step == 0)
{
cctkGH->cctk_time = (*Original_Time);
- cctkGH->CCTK_DELTA_TIME = (*Original_Delta_Time);
+ cctkGH->cctk_delta_time = (*Original_Delta_Time);
}
else if (CCTK_EQUALS(ODE_Method,"ICN"))
{
if (*MoL_Intermediate_Step == 1)
{
- cctkGH->CCTK_DELTA_TIME = (*Original_Delta_Time);
+ cctkGH->cctk_delta_time = (*Original_Delta_Time);
}
else
{
- cctkGH->CCTK_DELTA_TIME = 0.5*(*Original_Delta_Time);
+ cctkGH->cctk_delta_time = 0.5*(*Original_Delta_Time);
}
- cctkGH->cctk_time = (*Original_Time)-0.5*(*Original_Delta_Time);
+ cctkGH->cctk_time = (*Original_Time)-0.5*(*Original_Delta_Time)/cctkGH->cctk_timefac;
}
else if (CCTK_EQUALS(ODE_Method,"Generic"))
{
- cctkGH->CCTK_DELTA_TIME = RKBetaCoefficients[MoL_Intermediate_Steps -
+ cctkGH->cctk_delta_time = RKBetaCoefficients[MoL_Intermediate_Steps -
(*MoL_Intermediate_Step)] *
- (*Original_Delta_Time);
- previous_times[0] = (*Original_Time) - (*Original_Delta_Time);
+ (*Original_Delta_Time)/cctkGH->cctk_timefac;
+ previous_times[0] = (*Original_Time) - (*Original_Delta_Time)/cctkGH->cctk_timefac;
for (i = MoL_Intermediate_Steps - 1; i > *MoL_Intermediate_Step; i--)
{
for (j = MoL_Intermediate_Steps; j > i; j--)
@@ -151,7 +151,7 @@ int MoL_ResetTime(CCTK_ARGUMENTS)
MoL_Intermediate_Steps + j;
previous_times[MoL_Intermediate_Steps - i] =
RKBetaCoefficients[MoL_Intermediate_Steps - j] *
- (*Original_Delta_Time) +
+ (*Original_Delta_Time)/cctkGH->cctk_timefac +
RKAlphaCoefficients[alphaindex];
}
}
@@ -162,13 +162,13 @@ int MoL_ResetTime(CCTK_ARGUMENTS)
{
if (*MoL_Intermediate_Step == 1)
{
- cctkGH->CCTK_DELTA_TIME = 0.5*(*Original_Delta_Time);
+ cctkGH->cctk_delta_time = 0.5*(*Original_Delta_Time);
cctkGH->cctk_time = (*Original_Time);
}
}
#ifdef MOLDEBUG
printf("MoL has once more reset t: %f: and dt: %f.\n", cctkGH->cctk_time,
- cctkGH->CCTK_DELTA_TIME);
+ cctkGH->cctk_delta_time/cctkGH->cctk_timefac);
#endif
free(previous_times);