From 29091dc1b2dfc37f7e681b91a7dbfcab1c01f9f1 Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Sat, 7 Apr 2018 10:06:12 +0200 Subject: Make recovery from checkpoints work seamlessly. --- interface.ccl | 3 ++- schedule.ccl | 6 ++---- src/qms.c | 56 ++++++++++++++++++++++++-------------------------------- 3 files changed, 28 insertions(+), 37 deletions(-) diff --git a/interface.ccl b/interface.ccl index 324fc5f..7216b74 100644 --- a/interface.ccl +++ b/interface.ccl @@ -13,4 +13,5 @@ REQUIRES FUNCTION MoLRegisterSaveAndRestoreGroup public: CCTK_REAL W TYPE=GF -CCTK_REAL w_coeffs TYPE=array DIM=2 SIZE=basis_order_z,basis_order_r DISTRIB=constant +CCTK_REAL W_coeffs TYPE=array DIM=3 SIZE=2,basis_order_z,basis_order_r DISTRIB=constant +CCTK_REAL W_pred TYPE=array DIM=3 SIZE=2,basis_order_z,basis_order_r DISTRIB=constant diff --git a/schedule.ccl b/schedule.ccl index b4fcda2..d8eb719 100644 --- a/schedule.ccl +++ b/schedule.ccl @@ -26,7 +26,5 @@ SCHEDULE qms_init IN ADMBase_InitialData { } "" STORAGE: W - -if (export_coeffs) { - STORAGE: w_coeffs -} +STORAGE: W_pred +STORAGE: W_coeffs diff --git a/src/qms.c b/src/qms.c index 3d3fb64..4be70c6 100644 --- a/src/qms.c +++ b/src/qms.c @@ -256,24 +256,16 @@ static int context_init(cGH *cctkGH) return -ENOMEM; for (int i = 0; i < ARRAY_ELEMS(qms->solutions); i++) { - ret = posix_memalign((void**)&qms->solutions[i].coeffs, 32, - basis_order_r * basis_order_z * sizeof(*qms->solutions[i].coeffs)); - if (ret) - return -ENOMEM; - - memset(qms->solutions[i].coeffs, 0, basis_order_r * basis_order_z * sizeof(*qms->solutions[i].coeffs)); - qms->solutions[i].time = -i * cctkGH->cctk_delta_time; + qms->solutions[i].coeffs = W_coeffs + i * basis_order_r * basis_order_z; + //memset(qms->solutions[i].coeffs, 0, basis_order_r * basis_order_z * sizeof(*qms->solutions[i].coeffs)); + qms->solutions[i].time = cctkGH->cctk_time - (2 - i) * cctkGH->cctk_delta_time / (1.0 * (1 << solve_level)); fprintf(stderr, "init time %d %f\n", i, qms->solutions[i].time); } for (int i = 0; i < ARRAY_ELEMS(qms->pred); i++) { - ret = posix_memalign((void**)&qms->pred[i].coeffs, 32, - basis_order_r * basis_order_z * sizeof(*qms->pred[i].coeffs)); - if (ret) - return -ENOMEM; - - memset(qms->pred[i].coeffs, 0, basis_order_r * basis_order_z * sizeof(*qms->pred[i].coeffs)); - qms->pred[i].time = -i * cctkGH->cctk_delta_time; + qms->pred[i].coeffs = W_pred + i * basis_order_r * basis_order_z; + //memset(qms->pred[i].coeffs, 0, basis_order_r * basis_order_z * sizeof(*qms->pred[i].coeffs)); + qms->pred[i].time = cctkGH->cctk_time - (1 - i) * cctkGH->cctk_delta_time / (1.0 * (1 << solve_level)); fprintf(stderr, "init time %d %f\n", i, qms->pred[i].time); } @@ -316,30 +308,33 @@ void quasimaximal_slicing_axi_solve(CCTK_ARGUMENTS) /* add the solution to the list of past solutions */ { - double *tmp = ms->solutions[0].coeffs; + int N = ms->solver->nb_coeffs[0] * ms->solver->nb_coeffs[1]; - memmove(ms->solutions, ms->solutions + 1, sizeof(ms->solutions[0]) * (ARRAY_ELEMS(ms->solutions) - 1)); - ms->solutions[ARRAY_ELEMS(ms->solutions) - 1].coeffs = ms->solver->coeffs; - ms->solutions[ARRAY_ELEMS(ms->solutions) - 1].time = ms->gh->cctk_time; + memcpy(ms->solutions[0].coeffs, ms->solutions[1].coeffs, sizeof(*ms->solutions[0].coeffs) * N); + ms->solutions[0].time = ms->solutions[1].time; - ms->solver->coeffs = tmp; + memcpy(ms->solutions[1].coeffs, ms->solver->coeffs, sizeof(*ms->solutions[1].coeffs) * N); + ms->solutions[1].time = ms->gh->cctk_time; } /* linearly extrapolate the past two solution to predict the next one */ { - double time0 = ms->solutions[ARRAY_ELEMS(ms->solutions) - 2].time; - double time1 = ms->solutions[ARRAY_ELEMS(ms->solutions) - 1].time; - double time = time1 + ms->gh->cctk_delta_time; + int N = ms->solver->nb_coeffs[0] * ms->solver->nb_coeffs[1]; - double *coeffs0 = ms->solutions[ARRAY_ELEMS(ms->solutions) - 2].coeffs;; - double *coeffs1 = ms->solutions[ARRAY_ELEMS(ms->solutions) - 1].coeffs;; - double *pred = ms->pred[0].coeffs; + double time0 = ms->solutions[0].time; + double time1 = ms->solutions[1].time; + double time = time1 + ms->gh->cctk_delta_time / (1.0 * (1 << solve_level)); - for (int i = 0; i < ms->solver->nb_coeffs[0] * ms->solver->nb_coeffs[1]; i++) - pred[i] = (coeffs1[i] * (time - time0) / (time1 - time0) + coeffs0[i] * (time - time1) / (time0 - time1)); + double *coeffs0 = ms->solutions[0].coeffs; + double *coeffs1 = ms->solutions[1].coeffs; + double *pred = ms->pred[1].coeffs; + + memcpy(ms->pred[0].coeffs, ms->pred[1].coeffs, sizeof(*ms->pred[0].coeffs) * N); + ms->pred[0].time = ms->pred[1].time; + + for (int i = 0; i < N; i++) + ms->pred[1].coeffs[i] = (coeffs1[i] * (time - time0) / (time1 - time0) + coeffs0[i] * (time - time1) / (time0 - time1)); - memmove(ms->pred, ms->pred + 1, sizeof(ms->pred[0]) * (ARRAY_ELEMS(ms->pred) - 1)); - ms->pred[1].coeffs = pred; ms->pred[1].time = time; } } @@ -397,9 +392,6 @@ void quasimaximal_slicing_axi_eval(CCTK_ARGUMENTS) } #endif - if (export_coeffs) - memcpy(w_coeffs, coeffs, sizeof(*w_coeffs) * ms->solver->nb_coeffs[0] * ms->solver->nb_coeffs[1]); - CCTK_TimerStart("QuasiMaximalSlicing_Expand"); expand_start = gettime(); #if 0 -- cgit v1.2.3