10#ifndef __TSEPADAPT_HPP__
11#define __TSEPADAPT_HPP__
13#include <petsc/private/tsimpl.h>
20 PetscInt *next_sc, PetscReal *next_h,
21 PetscBool *accept, PetscReal *wlte,
22 PetscReal *wltea, PetscReal *wlter) {
36 ierr = TSGetSNES(ts,&snes);
39 SNESConvergedReason reason;
40 ierr = SNESGetConvergedReason(snes,&reason);
44 ierr = SNESGetIterationNumber(snes, &it);
52 "\tDiverged set step length: it = %d, h = %3.4g set h = %3.4g \n", it,
54 }
else if (reason > 0) {
55 h *= sqrt(
static_cast<double>(6) /
static_cast<double>(it + 1));
56 *next_h = PetscClipInterval(
h, adapt->dt_min, adapt->dt_max);
59 "\tConverged set step length: it = %d, h = %3.4g set h = %3.4g \n", it,
63 PetscFunctionReturn(0);
69 PetscFunctionReturn(0);
76 ierr = PetscFree(adapt->data);
78 PetscFunctionReturn(0);
84 ierr = PetscNewLog(adapt, &ep);
86 adapt->data = (
void *)ep;
90 PetscFunctionReturn(0);
static PetscErrorCode TSAdaptReset_EP(TSAdapt adapt)
static PetscErrorCode TSAdaptChoose_EP(TSAdapt adapt, TS ts, PetscReal h, PetscInt *next_sc, PetscReal *next_h, PetscBool *accept, PetscReal *wlte, PetscReal *wltea, PetscReal *wlter)
PetscErrorCode TSAdaptCreate_EP(TSAdapt adapt)
static PetscErrorCode TSAdaptDestroy_EP(TSAdapt adapt)
static PetscErrorCode ierr
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.