v0.13.2
Loading...
Searching...
No Matches
TSEPAdapt.hpp
Go to the documentation of this file.
1/**
2 * \file TSEPAdapt.hpp
3 * \brief Step adaptation
4 *
5 * \brief Problem implementation for mix element for large-strain elasticity
6 *
7 * \todo Implementation of plasticity
8 */
9
10#ifndef __TSEPADAPT_HPP__
11#define __TSEPADAPT_HPP__
12
13#include <petsc/private/tsimpl.h>
14
15#define TSADAPTEP "ep"
16
17struct TSAdapt_EP {};
18
19static PetscErrorCode TSAdaptChoose_EP(TSAdapt adapt, TS ts, PetscReal h,
20 PetscInt *next_sc, PetscReal *next_h,
21 PetscBool *accept, PetscReal *wlte,
22 PetscReal *wltea, PetscReal *wlter) {
23 PetscFunctionBegin;
24
25 TSAdapt_EP *basic = static_cast<TSAdapt_EP *>(adapt->data);
26
27 *next_sc = 0; /* Reuse the same order scheme */
28 *wlte = -1; /* Weighted local truncation error was not evaluated */
29 *wltea = -1; /* Weighted absolute local truncation error is not used */
30 *wlter = -1; /* Weighted relative local truncation error is not used */
31
32 *accept = PETSC_TRUE;
33 *next_h = h; /* Reuse the old step */
34
35 SNES snes;
36 ierr = TSGetSNES(ts,&snes);
38
39 SNESConvergedReason reason;
40 ierr = SNESGetConvergedReason(snes,&reason);
42
43 int it;
44 ierr = SNESGetIterationNumber(snes, &it);
46
47 if(reason < 0) {
48 h *= 0.75;
49 *next_h = h;
50 PetscPrintf(
51 PETSC_COMM_WORLD,
52 "\tDiverged set step length: it = %d, h = %3.4g set h = %3.4g \n", it,
53 h, *next_h);
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);
57 PetscPrintf(
58 PETSC_COMM_WORLD,
59 "\tConverged set step length: it = %d, h = %3.4g set h = %3.4g \n", it,
60 h, *next_h);
61 }
62
63 PetscFunctionReturn(0);
64}
65
66static PetscErrorCode TSAdaptReset_EP(TSAdapt adapt) {
67 TSAdapt_EP *basic = static_cast<TSAdapt_EP *>(adapt->data);
68 PetscFunctionBegin;
69 PetscFunctionReturn(0);
70}
71
72static PetscErrorCode TSAdaptDestroy_EP(TSAdapt adapt) {
73 PetscFunctionBegin;
74 ierr = TSAdaptReset_EP(adapt);
76 ierr = PetscFree(adapt->data);
78 PetscFunctionReturn(0);
79}
80
81PetscErrorCode TSAdaptCreate_EP(TSAdapt adapt) {
82 TSAdapt_EP *ep;
83 PetscFunctionBegin;
84 ierr = PetscNewLog(adapt, &ep);
86 adapt->data = (void *)ep;
87 adapt->ops->choose = TSAdaptChoose_EP;
88 adapt->ops->reset = TSAdaptReset_EP;
89 adapt->ops->destroy = TSAdaptDestroy_EP;
90 PetscFunctionReturn(0);
91}
92
93#endif // __TSEPADAPT_HPP__
static PetscErrorCode TSAdaptReset_EP(TSAdapt adapt)
Definition: TSEPAdapt.hpp:66
static PetscErrorCode TSAdaptChoose_EP(TSAdapt adapt, TS ts, PetscReal h, PetscInt *next_sc, PetscReal *next_h, PetscBool *accept, PetscReal *wlte, PetscReal *wltea, PetscReal *wlter)
Definition: TSEPAdapt.hpp:19
PetscErrorCode TSAdaptCreate_EP(TSAdapt adapt)
Definition: TSEPAdapt.hpp:81
static PetscErrorCode TSAdaptDestroy_EP(TSAdapt adapt)
Definition: TSEPAdapt.hpp:72
static PetscErrorCode ierr
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
double h