From 8d3eff4200fba74f3b9260402117b0833f802c3a Mon Sep 17 00:00:00 2001 From: "Daniel R. Reynolds" Date: Wed, 11 Sep 2024 16:27:27 -0500 Subject: [PATCH] Bugfix: ARKODE Stepper SetDefaults memory leak fix (#560) Updated stepper `SetDefaults` routines to deallocate pre-existing structures before nullifying, and moved initial creation of the `SUNAdaptController` into stepper `SetDefaults` routines instead of `arkCreate`. This allows non-adaptive steppers to leave the SUNAdaptController object unallocated. --------- Co-authored-by: Steven Roberts Co-authored-by: David Gardner --- CHANGELOG.md | 3 ++ doc/shared/RecentChanges.rst | 3 ++ src/arkode/arkode.c | 46 +++++++----------- src/arkode/arkode_adapt.c | 12 ++++- src/arkode/arkode_arkstep_io.c | 62 +++++++++++++++++++++-- src/arkode/arkode_erkstep_io.c | 43 ++++++++++------ src/arkode/arkode_io.c | 89 +++++++++++++++++++++------------- src/arkode/arkode_mristep_io.c | 31 ++++++++---- 8 files changed, 198 insertions(+), 91 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 247f6f26ea..b73b8b7131 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -26,6 +26,9 @@ Fixed the loading of ARKStep's default first order explicit method. Fixed a CMake bug regarding usage of missing "print_warning" macro that was only triggered when the deprecated `CUDA_ARCH` option was used. +Fixed a memory leak that could occur if ``ARKodeSetDefaults`` is called +repeatedly. + Fixed compilation errors when building the Trilinos Teptra NVector with CUDA support. diff --git a/doc/shared/RecentChanges.rst b/doc/shared/RecentChanges.rst index ed8b0b13a6..5d4e6e9554 100644 --- a/doc/shared/RecentChanges.rst +++ b/doc/shared/RecentChanges.rst @@ -22,6 +22,9 @@ Fixed the loading of ARKStep's default first order explicit method. Fixed a CMake bug regarding usage of missing "print_warning" macro that was only triggered when the deprecated ``CUDA_ARCH`` option was used. +Fixed a memory leak that could occur if :c:func:`ARKodeSetDefaults` is called +repeatedly. + Fixed compilation errors when building the Trilinos Teptra NVector with CUDA support. diff --git a/src/arkode/arkode.c b/src/arkode/arkode.c index be15bd2b2f..082a2c1755 100644 --- a/src/arkode/arkode.c +++ b/src/arkode/arkode.c @@ -1364,7 +1364,6 @@ void ARKodePrintMem(void* arkode_mem, FILE* outfile) ARKodeMem arkCreate(SUNContext sunctx) { int iret; - long int lenrw, leniw; ARKodeMem ark_mem; if (!sunctx) @@ -1485,21 +1484,6 @@ ARKodeMem arkCreate(SUNContext sunctx) ark_mem->lrw += ARK_ADAPT_LRW; ark_mem->liw += ARK_ADAPT_LIW; - /* Allocate default step controller (PID) and note storage */ - ark_mem->hadapt_mem->hcontroller = SUNAdaptController_PID(sunctx); - if (ark_mem->hadapt_mem->hcontroller == NULL) - { - arkProcessError(NULL, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, - "Allocation of step controller object failed"); - ARKodeFree((void**)&ark_mem); - return (NULL); - } - ark_mem->hadapt_mem->owncontroller = SUNTRUE; - (void)SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, - &leniw); - ark_mem->lrw += lenrw; - ark_mem->liw += leniw; - /* Initialize the interpolation structure to NULL */ ark_mem->interp = NULL; ark_mem->interp_type = ARK_INTERP_HERMITE; @@ -1523,7 +1507,7 @@ ARKodeMem arkCreate(SUNContext sunctx) ark_mem->h = ZERO; ark_mem->h0u = ZERO; - /* Set default values for integrator optional inputs */ + /* Set default values for integrator and stepper optional inputs */ iret = ARKodeSetDefaults(ark_mem); if (iret != ARK_SUCCESS) { @@ -1706,12 +1690,15 @@ int arkInit(ARKodeMem ark_mem, sunrealtype t0, N_Vector y0, int init_type) ark_mem->tolsf = ONE; /* Reset error controller object */ - retval = SUNAdaptController_Reset(ark_mem->hadapt_mem->hcontroller); - if (retval != SUN_SUCCESS) + if (ark_mem->hadapt_mem->hcontroller) { - arkProcessError(ark_mem, ARK_CONTROLLER_ERR, __LINE__, __func__, __FILE__, - "Unable to reset error controller object"); - return (ARK_CONTROLLER_ERR); + retval = SUNAdaptController_Reset(ark_mem->hadapt_mem->hcontroller); + if (retval != SUN_SUCCESS) + { + arkProcessError(ark_mem, ARK_CONTROLLER_ERR, __LINE__, __func__, + __FILE__, "Unable to reset error controller object"); + return (ARK_CONTROLLER_ERR); + } } /* Adaptivity counters */ @@ -2543,13 +2530,16 @@ int arkCompleteStep(ARKodeMem ark_mem, sunrealtype dsm) ark_mem->fn_is_current = SUNFALSE; /* Notify time step controller object of successful step */ - retval = SUNAdaptController_UpdateH(ark_mem->hadapt_mem->hcontroller, - ark_mem->h, dsm); - if (retval != SUN_SUCCESS) + if (ark_mem->hadapt_mem->hcontroller) { - arkProcessError(ark_mem, ARK_CONTROLLER_ERR, __LINE__, __func__, __FILE__, - "Failure updating controller object"); - return (ARK_CONTROLLER_ERR); + retval = SUNAdaptController_UpdateH(ark_mem->hadapt_mem->hcontroller, + ark_mem->h, dsm); + if (retval != SUN_SUCCESS) + { + arkProcessError(ark_mem, ARK_CONTROLLER_ERR, __LINE__, __func__, __FILE__, + "Failure updating controller object"); + return (ARK_CONTROLLER_ERR); + } } /* update scalar quantities */ diff --git a/src/arkode/arkode_adapt.c b/src/arkode/arkode_adapt.c index 7fef0ec2ae..ada5cda343 100644 --- a/src/arkode/arkode_adapt.c +++ b/src/arkode/arkode_adapt.c @@ -82,7 +82,10 @@ void arkPrintAdaptMem(ARKodeHAdaptMem hadapt_mem, FILE* outfile) fprintf(outfile, " ark_hadapt: stability function data pointer = %p\n", hadapt_mem->estab_data); } - (void)SUNAdaptController_Write(hadapt_mem->hcontroller, outfile); + if (hadapt_mem->hcontroller != NULL) + { + (void)SUNAdaptController_Write(hadapt_mem->hcontroller, outfile); + } } } @@ -98,6 +101,13 @@ int arkAdapt(ARKodeMem ark_mem, ARKodeHAdaptMem hadapt_mem, N_Vector ycur, sunrealtype h_acc, h_cfl, int_dir; int controller_order; + /* Return with no stepsize adjustment if the controller is NULL */ + if (hadapt_mem->hcontroller == NULL) + { + ark_mem->eta = ONE; + return (ARK_SUCCESS); + } + /* Request error-based step size from adaptivity controller */ if (hadapt_mem->pq == 0) { diff --git a/src/arkode/arkode_arkstep_io.c b/src/arkode/arkode_arkstep_io.c index 7df74cb590..ae90b65387 100644 --- a/src/arkode/arkode_arkstep_io.c +++ b/src/arkode/arkode_arkstep_io.c @@ -653,6 +653,8 @@ int arkStep_SetUserData(ARKodeMem ark_mem, void* user_data) int arkStep_SetDefaults(ARKodeMem ark_mem) { ARKodeARKStepMem step_mem; + sunindextype Blrw, Bliw; + long int lenrw, leniw; int retval; /* access ARKodeARKStepMem structure */ @@ -677,12 +679,66 @@ int arkStep_SetDefaults(ARKodeMem ark_mem) step_mem->msbp = MSBP; /* max steps between updates to J or P */ step_mem->stages = 0; /* no stages */ step_mem->istage = 0; /* current stage */ - step_mem->Be = NULL; /* no Butcher tables */ - step_mem->Bi = NULL; - step_mem->NLS = NULL; /* no nonlinear solver object */ step_mem->jcur = SUNFALSE; step_mem->convfail = ARK_NO_FAILURES; step_mem->stage_predict = NULL; /* no user-supplied stage predictor */ + + /* Remove pre-existing Butcher tables */ + if (step_mem->Be) + { + ARKodeButcherTable_Space(step_mem->Be, &Bliw, &Blrw); + ark_mem->liw -= Bliw; + ark_mem->lrw -= Blrw; + ARKodeButcherTable_Free(step_mem->Be); + } + step_mem->Be = NULL; + if (step_mem->Bi) + { + ARKodeButcherTable_Space(step_mem->Bi, &Bliw, &Blrw); + ark_mem->liw -= Bliw; + ark_mem->lrw -= Blrw; + ARKodeButcherTable_Free(step_mem->Bi); + } + step_mem->Bi = NULL; + + /* Remove pre-existing nonlinear solver object */ + if (step_mem->NLS && step_mem->ownNLS) { SUNNonlinSolFree(step_mem->NLS); } + step_mem->NLS = NULL; + + /* Remove pre-existing SUNAdaptController object, and replace with "PID" */ + if (ark_mem->hadapt_mem->owncontroller) + { + retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, + &leniw); + if (retval == SUN_SUCCESS) + { + ark_mem->liw -= leniw; + ark_mem->lrw -= lenrw; + } + retval = SUNAdaptController_Destroy(ark_mem->hadapt_mem->hcontroller); + ark_mem->hadapt_mem->owncontroller = SUNFALSE; + if (retval != SUN_SUCCESS) + { + arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, + "SUNAdaptController_Destroy failure"); + return (ARK_MEM_FAIL); + } + } + ark_mem->hadapt_mem->hcontroller = SUNAdaptController_PID(ark_mem->sunctx); + if (ark_mem->hadapt_mem->hcontroller == NULL) + { + arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, + "SUNAdaptControllerPID allocation failure"); + return (ARK_MEM_FAIL); + } + ark_mem->hadapt_mem->owncontroller = SUNTRUE; + retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, + &leniw); + if (retval == SUN_SUCCESS) + { + ark_mem->liw += leniw; + ark_mem->lrw += lenrw; + } return (ARK_SUCCESS); } diff --git a/src/arkode/arkode_erkstep_io.c b/src/arkode/arkode_erkstep_io.c index ffbf9a18fb..77987b4245 100644 --- a/src/arkode/arkode_erkstep_io.c +++ b/src/arkode/arkode_erkstep_io.c @@ -259,23 +259,42 @@ int erkStep_SetRelaxFn(ARKodeMem ark_mem, ARKRelaxFn rfn, ARKRelaxJacFn rjac) int erkStep_SetDefaults(ARKodeMem ark_mem) { ARKodeERKStepMem step_mem; - int retval; + sunindextype Blrw, Bliw; long int lenrw, leniw; + int retval; /* access ARKodeERKStepMem structure */ retval = erkStep_AccessStepMem(ark_mem, __func__, &step_mem); if (retval != ARK_SUCCESS) { return (retval); } - /* Remove current SUNAdaptController object, and replace with "PI" */ - retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, - &leniw); - if (retval == SUN_SUCCESS) + /* Set default values for integrator optional inputs */ + step_mem->q = Q_DEFAULT; /* method order */ + step_mem->p = 0; /* embedding order */ + step_mem->stages = 0; /* no stages */ + ark_mem->hadapt_mem->etamxf = SUN_RCONST(0.3); /* max change on error-failed step */ + ark_mem->hadapt_mem->safety = SUN_RCONST(0.99); /* step adaptivity safety factor */ + ark_mem->hadapt_mem->growth = SUN_RCONST(25.0); /* step adaptivity growth factor */ + + /* Remove pre-existing Butcher table */ + if (step_mem->B) { - ark_mem->liw -= leniw; - ark_mem->lrw -= lenrw; + ARKodeButcherTable_Space(step_mem->B, &Bliw, &Blrw); + ark_mem->liw -= Bliw; + ark_mem->lrw -= Blrw; + ARKodeButcherTable_Free(step_mem->B); } + step_mem->B = NULL; + + /* Remove current SUNAdaptController object, and replace with "PI" */ if (ark_mem->hadapt_mem->owncontroller) { + retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, + &leniw); + if (retval == SUN_SUCCESS) + { + ark_mem->liw -= leniw; + ark_mem->lrw -= lenrw; + } retval = SUNAdaptController_Destroy(ark_mem->hadapt_mem->hcontroller); ark_mem->hadapt_mem->owncontroller = SUNFALSE; if (retval != SUN_SUCCESS) @@ -302,15 +321,7 @@ int erkStep_SetDefaults(ARKodeMem ark_mem) ark_mem->lrw += lenrw; } - /* Set default values for integrator optional inputs - (overwrite some adaptivity params for ERKStep use) */ - step_mem->q = Q_DEFAULT; /* method order */ - step_mem->p = 0; /* embedding order */ - step_mem->stages = 0; /* no stages */ - step_mem->B = NULL; /* no Butcher table */ - ark_mem->hadapt_mem->etamxf = SUN_RCONST(0.3); /* max change on error-failed step */ - ark_mem->hadapt_mem->safety = SUN_RCONST(0.99); /* step adaptivity safety factor */ - ark_mem->hadapt_mem->growth = SUN_RCONST(25.0); /* step adaptivity growth factor */ + /* Update some controller parameters */ (void)SUNAdaptController_SetErrorBias(ark_mem->hadapt_mem->hcontroller, SUN_RCONST(1.2)); (void)SUNAdaptController_SetParams_PI(ark_mem->hadapt_mem->hcontroller, diff --git a/src/arkode/arkode_io.c b/src/arkode/arkode_io.c index 357912409d..dd1325f794 100644 --- a/src/arkode/arkode_io.c +++ b/src/arkode/arkode_io.c @@ -55,13 +55,6 @@ int ARKodeSetDefaults(void* arkode_mem) } ark_mem = (ARKodeMem)arkode_mem; - /* Set stepper defaults (if provided) */ - if (ark_mem->step_setdefaults) - { - retval = ark_mem->step_setdefaults(arkode_mem); - if (retval != ARK_SUCCESS) { return retval; } - } - /* Set default values for integrator optional inputs */ ark_mem->use_compensated_sums = SUNFALSE; ark_mem->fixedstep = SUNFALSE; /* default to use adaptive steps */ @@ -106,6 +99,14 @@ int ARKodeSetDefaults(void* arkode_mem) ark_mem->hadapt_mem->p = 0; /* no default embedding order */ ark_mem->hadapt_mem->q = 0; /* no default method order */ ark_mem->hadapt_mem->adjust = ADJUST; /* controller order adjustment */ + + /* Set stepper defaults (if provided) */ + if (ark_mem->step_setdefaults) + { + retval = ark_mem->step_setdefaults(arkode_mem); + if (retval != ARK_SUCCESS) { return retval; } + } + return (ARK_SUCCESS); } @@ -910,15 +911,17 @@ int ARKodeSetAdaptController(void* arkode_mem, SUNAdaptController C) /* Remove current SUNAdaptController object (delete if owned, and then nullify pointer) */ - retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, - &leniw); - if (retval == SUN_SUCCESS) - { - ark_mem->liw -= leniw; - ark_mem->lrw -= lenrw; - } - if (ark_mem->hadapt_mem->owncontroller) + if (ark_mem->hadapt_mem->owncontroller && + (ark_mem->hadapt_mem->hcontroller != NULL)) { + retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, + &leniw); + if (retval == SUN_SUCCESS) + { + ark_mem->liw -= leniw; + ark_mem->lrw -= lenrw; + } + retval = SUNAdaptController_Destroy(ark_mem->hadapt_mem->hcontroller); ark_mem->hadapt_mem->owncontroller = SUNFALSE; if (retval != SUN_SUCCESS) @@ -1043,8 +1046,11 @@ int ARKodeSetInitStep(void* arkode_mem, sunrealtype hin) ark_mem->h0u = ZERO; /* Reset error controller (e.g., error and step size history) */ - retval = SUNAdaptController_Reset(ark_mem->hadapt_mem->hcontroller); - if (retval != SUN_SUCCESS) { return (ARK_CONTROLLER_ERR); } + if (ark_mem->hadapt_mem->hcontroller != NULL) + { + retval = SUNAdaptController_Reset(ark_mem->hadapt_mem->hcontroller); + if (retval != SUN_SUCCESS) { return (ARK_CONTROLLER_ERR); } + } return (ARK_SUCCESS); } @@ -1641,6 +1647,14 @@ int ARKodeSetErrorBias(void* arkode_mem, sunrealtype bias) return (ARK_STEPPER_UNSUPPORTED); } + /* Return an error if there is not a current SUNAdaptController */ + if (ark_mem->hadapt_mem->hcontroller == NULL) + { + arkProcessError(ark_mem, ARK_MEM_NULL, __LINE__, __func__, + __FILE__, "SUNAdaptController NULL -- must be set before setting the error bias"); + return (ARK_MEM_NULL); + } + /* set allowed value, otherwise set default */ if (bias < ONE) { @@ -2947,7 +2961,10 @@ int ARKodeWriteParameters(void* arkode_mem, FILE* fp) fprintf(fp, " Default explicit stability function\n"); } else { fprintf(fp, " User provided explicit stability function\n"); } - (void)SUNAdaptController_Write(ark_mem->hadapt_mem->hcontroller, fp); + if (ark_mem->hadapt_mem->hcontroller != NULL) + { + (void)SUNAdaptController_Write(ark_mem->hadapt_mem->hcontroller, fp); + } fprintf(fp, " Maximum number of error test failures = %i\n", ark_mem->maxnef); fprintf(fp, " Maximum number of convergence test failures = %i\n", @@ -3050,15 +3067,17 @@ int arkSetAdaptivityMethod(void* arkode_mem, int imethod, int idefault, int pq, /* Remove current SUNAdaptController object (delete if owned, and then nullify pointer) */ - retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, - &leniw); - if (retval == SUN_SUCCESS) - { - ark_mem->liw -= leniw; - ark_mem->lrw -= lenrw; - } - if (ark_mem->hadapt_mem->owncontroller) + if (ark_mem->hadapt_mem->owncontroller && + (ark_mem->hadapt_mem->hcontroller != NULL)) { + retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, + &leniw); + if (retval == SUN_SUCCESS) + { + ark_mem->liw -= leniw; + ark_mem->lrw -= lenrw; + } + retval = SUNAdaptController_Destroy(ark_mem->hadapt_mem->hcontroller); ark_mem->hadapt_mem->owncontroller = SUNFALSE; if (retval != SUN_SUCCESS) @@ -3251,15 +3270,17 @@ int arkSetAdaptivityFn(void* arkode_mem, ARKAdaptFn hfun, void* h_data) /* Remove current SUNAdaptController object (delete if owned, and then nullify pointer) */ - retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, - &leniw); - if (retval == SUN_SUCCESS) - { - ark_mem->liw -= leniw; - ark_mem->lrw -= lenrw; - } - if (ark_mem->hadapt_mem->owncontroller) + if (ark_mem->hadapt_mem->owncontroller && + (ark_mem->hadapt_mem->hcontroller != NULL)) { + retval = SUNAdaptController_Space(ark_mem->hadapt_mem->hcontroller, &lenrw, + &leniw); + if (retval == SUN_SUCCESS) + { + ark_mem->liw -= leniw; + ark_mem->lrw -= lenrw; + } + retval = SUNAdaptController_Destroy(ark_mem->hadapt_mem->hcontroller); ark_mem->hadapt_mem->owncontroller = SUNFALSE; if (retval != SUN_SUCCESS) diff --git a/src/arkode/arkode_mristep_io.c b/src/arkode/arkode_mristep_io.c index d803141044..1ad11ab59e 100644 --- a/src/arkode/arkode_mristep_io.c +++ b/src/arkode/arkode_mristep_io.c @@ -231,6 +231,7 @@ int mriStep_SetUserData(ARKodeMem ark_mem, void* user_data) int mriStep_SetDefaults(ARKodeMem ark_mem) { ARKodeMRIStepMem step_mem; + sunindextype lenrw, leniw; int retval; /* access ARKodeMRIStepMem structure */ @@ -248,15 +249,27 @@ int mriStep_SetDefaults(ARKodeMem ark_mem) step_mem->nlscoef = NLSCOEF; /* nonlinear tolerance coefficient */ step_mem->crdown = CRDOWN; /* nonlinear convergence estimate coeff. */ step_mem->rdiv = RDIV; /* nonlinear divergence tolerance */ - step_mem->dgmax = DGMAX; /* max gamma change before recomputing J or P */ - step_mem->msbp = MSBP; /* max steps between updates to J or P */ - step_mem->stages = 0; /* no stages */ - step_mem->istage = 0; /* current stage index */ - step_mem->MRIC = NULL; /* no slow->fast coupling */ - step_mem->NLS = NULL; /* no nonlinear solver object */ - step_mem->jcur = SUNFALSE; - step_mem->convfail = ARK_NO_FAILURES; - step_mem->stage_predict = NULL; /* no user-supplied stage predictor */ + step_mem->dgmax = DGMAX; /* max gamma change to recompute J or P */ + step_mem->msbp = MSBP; /* max steps between updating J or P */ + step_mem->stages = 0; /* no stages */ + step_mem->istage = 0; /* current stage index */ + step_mem->jcur = SUNFALSE; + step_mem->convfail = ARK_NO_FAILURES; + step_mem->stage_predict = NULL; /* no user-supplied stage predictor */ + + /* Remove pre-existing nonlinear solver object */ + if (step_mem->NLS && step_mem->ownNLS) { SUNNonlinSolFree(step_mem->NLS); } + step_mem->NLS = NULL; + + /* Remove pre-existing coupling table */ + if (step_mem->MRIC) + { + MRIStepCoupling_Space(step_mem->MRIC, &leniw, &lenrw); + ark_mem->lrw -= lenrw; + ark_mem->liw -= leniw; + MRIStepCoupling_Free(step_mem->MRIC); + } + step_mem->MRIC = NULL; return (ARK_SUCCESS); }