v0.16.3
Loading...
Searching...
No Matches
Public Member Functions | Static Public Member Functions | Public Attributes | Static Private Member Functions | Private Attributes | List of all members
TSPrePostProc Struct Reference

Set of functions called by PETSc solver used to refine and update mesh. More...

Collaboration diagram for TSPrePostProc:
[legend]

Public Member Functions

 TSPrePostProc ()=default
 
virtual ~TSPrePostProc ()=default
 
MoFEMErrorCode tsSetUp (TS ts)
 Used to setup TS solver.
 
 TSPrePostProc ()=default
 
virtual ~TSPrePostProc ()=default
 
MoFEMErrorCode tsSetUp (TS ts)
 Used to setup TS solver.
 
SmartPetscObj< VecScatter > getScatter (Vec x, Vec y, enum FR fr)
 Get scatter context for vector operations.
 
SmartPetscObj< Vec > getSubVector ()
 Create sub-problem vector.
 

Static Public Member Functions

static MoFEMErrorCode tsPostStage (TS ts, PetscReal stagetime, PetscInt stageindex, Vec *Y)
 [Boundary condition]
 
static MoFEMErrorCode tsPostStep (TS ts)
 
static MoFEMErrorCode tsPreStep (TS ts)
 

Public Attributes

Example * fsRawPtr
 
SmartPetscObj< DM > solverSubDM
 
SmartPetscObj< Vec > globSol
 
FreeSurface * fsRawPtr
 

Static Private Member Functions

static MoFEMErrorCode tsPreProc (TS ts)
 Pre process time step.
 
static MoFEMErrorCode tsPostProc (TS ts)
 Post process time step.
 
static MoFEMErrorCode tsPreStage (TS ts)
 Pre-stage processing for time stepping.
 
static MoFEMErrorCode tsSetIFunction (TS ts, PetscReal t, Vec u, Vec u_t, Vec f, void *ctx)
 Set implicit function for time stepping.
 
static MoFEMErrorCode tsSetIJacobian (TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
 Set implicit Jacobian for time stepping.
 
static MoFEMErrorCode tsMonitor (TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
 Monitor solution during time stepping.
 
static MoFEMErrorCode pcSetup (PC pc)
 Setup preconditioner.
 
static MoFEMErrorCode pcApply (PC pc, Vec pc_f, Vec pc_x)
 Apply preconditioner.
 

Private Attributes

SmartPetscObj< Vec > globRes
 
SmartPetscObj< Mat > subB
 
SmartPetscObj< KSP > subKSP
 
boost::shared_ptr< SnesCtx > snesCtxPtr
 
boost::shared_ptr< TsCtx > tsCtxPtr
 

Detailed Description

Set of functions called by PETSc solver used to refine and update mesh.

Note
Currently theta method is only handled by this code.
Examples
mofem/tutorials/adv-4_dynamic_first_order_con_law/dynamic_first_order_con_law.cpp, and mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 380 of file dynamic_first_order_con_law.cpp.

Constructor & Destructor Documentation

◆ TSPrePostProc() [1/2]

TSPrePostProc::TSPrePostProc ( )
default

◆ ~TSPrePostProc() [1/2]

virtual TSPrePostProc::~TSPrePostProc ( )
virtualdefault

◆ TSPrePostProc() [2/2]

TSPrePostProc::TSPrePostProc ( )
default

◆ ~TSPrePostProc() [2/2]

virtual TSPrePostProc::~TSPrePostProc ( )
virtualdefault

Member Function Documentation

◆ getScatter()

SmartPetscObj< VecScatter > TSPrePostProc::getScatter ( Vec  x,
Vec  y,
enum FR  fr 
)

Get scatter context for vector operations.

Parameters
xLocal sub-vector
yGlobal vector
frDirection flag (F=forward, R=reverse)
Returns
SmartPetscObj<VecScatter> Scatter context for vector operations

Creates scatter context to transfer data between global and sub-problem vectors

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3596 of file free_surface.cpp.

3596 {
3597 if (auto ptr = tsPrePostProc.lock()) {
3598 auto prb_ptr = ptr->fsRawPtr->mField.get_problem("SUB_SOLVER");
3599 if (auto sub_data = prb_ptr->getSubData()) {
3600 auto is = sub_data->getSmartColIs();
3601 VecScatter s;
3602 if (fr == R) {
3603 CHK_THROW_MESSAGE(VecScatterCreate(x, PETSC_NULLPTR, y, is, &s),
3604 "crate scatter");
3605 } else {
3606 CHK_THROW_MESSAGE(VecScatterCreate(x, is, y, PETSC_NULLPTR, &s),
3607 "crate scatter");
3608 }
3609 return SmartPetscObj<VecScatter>(s);
3610 }
3611 }
3614}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
static boost::weak_ptr< TSPrePostProc > tsPrePostProc
@ R
intrusive_ptr for managing petsc objects

◆ getSubVector()

SmartPetscObj< Vec > TSPrePostProc::getSubVector ( )

Create sub-problem vector.

Returns
SmartPetscObj<Vec> Vector compatible with sub-problem DM

Creates a vector with proper size and ghost structure for sub-problem

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3616 of file free_surface.cpp.

3616 {
3618}
auto createDMVector(DM dm, RowColData rc=RowColData::COL)
Get smart vector from DM.
Definition DMMoFEM.hpp:1250
SmartPetscObj< DM > solverSubDM

◆ pcApply()

MoFEMErrorCode TSPrePostProc::pcApply ( PC  pc,
Vec  pc_f,
Vec  pc_x 
)
staticprivate

Apply preconditioner.

Parameters
pcPETSc preconditioner object
pc_fInput vector (right-hand side)
pc_xOutput vector (preconditioned solution)
Returns
MoFEMErrorCode Success/failure code

Applies preconditioner by solving sub-problem with KSP

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3576 of file free_surface.cpp.

3576 {
3578 if (auto ptr = tsPrePostProc.lock()) {
3579 auto sub_x = ptr->getSubVector();
3580 auto sub_f = vectorDuplicate(sub_x);
3581 auto scatter = ptr->getScatter(sub_x, pc_x, R);
3582 CHKERR VecScatterBegin(scatter, pc_f, sub_f, INSERT_VALUES,
3583 SCATTER_REVERSE);
3584 CHKERR VecScatterEnd(scatter, pc_f, sub_f, INSERT_VALUES, SCATTER_REVERSE);
3585 CHKERR KSPSetOperators(ptr->subKSP, ptr->subB, ptr->subB);
3586 MOFEM_LOG("FS", Sev::verbose) << "PCShell solve";
3587 CHKERR KSPSolve(ptr->subKSP, sub_f, sub_x);
3588 MOFEM_LOG("FS", Sev::verbose) << "PCShell solve <- done";
3589 CHKERR VecScatterBegin(scatter, sub_x, pc_x, INSERT_VALUES,
3590 SCATTER_FORWARD);
3591 CHKERR VecScatterEnd(scatter, sub_x, pc_x, INSERT_VALUES, SCATTER_FORWARD);
3592 }
3594};
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MOFEM_LOG(channel, severity)
Log.
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.

◆ pcSetup()

MoFEMErrorCode TSPrePostProc::pcSetup ( PC  pc)
staticprivate

Setup preconditioner.

Parameters
pcPETSc preconditioner object
Returns
MoFEMErrorCode Success/failure code

Initializes KSP solver for shell preconditioner

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3566 of file free_surface.cpp.

3566 {
3568 if (auto ptr = tsPrePostProc.lock()) {
3569 MOFEM_LOG("FS", Sev::verbose) << "SetUP sub PC";
3570 ptr->subKSP = createKSP(ptr->fsRawPtr->mField.get_comm());
3571 CHKERR KSPSetFromOptions(ptr->subKSP);
3572 }
3574};
auto createKSP(MPI_Comm comm)

◆ tsMonitor()

MoFEMErrorCode TSPrePostProc::tsMonitor ( TS  ts,
PetscInt  step,
PetscReal  t,
Vec  u,
void *  ctx 
)
staticprivate

Monitor solution during time stepping.

Parameters
tsPETSc time stepping object
stepCurrent time step number
tCurrent time value
uCurrent solution vector
ctxUser context (pointer to FreeSurface)
Returns
MoFEMErrorCode Success/failure code

Called after each time step to monitor solution and save output Wrapper for TS monitor

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3541 of file free_surface.cpp.

3542 {
3544 if (auto ptr = tsPrePostProc.lock()) {
3545 auto get_norm = [&](auto x) {
3546 double nrm;
3547 CHKERR VecNorm(x, NORM_2, &nrm);
3548 return nrm;
3549 };
3550
3551 auto sub_u = ptr->getSubVector();
3552 auto scatter = ptr->getScatter(sub_u, u, R);
3553 CHKERR VecScatterBegin(scatter, u, sub_u, INSERT_VALUES, SCATTER_REVERSE);
3554 CHKERR VecScatterEnd(scatter, u, sub_u, INSERT_VALUES, SCATTER_REVERSE);
3555 CHKERR VecGhostUpdateBegin(sub_u, INSERT_VALUES, SCATTER_FORWARD);
3556 CHKERR VecGhostUpdateEnd(sub_u, INSERT_VALUES, SCATTER_FORWARD);
3557
3558 MOFEM_LOG("FS", Sev::verbose)
3559 << "u norm " << get_norm(u) << " u sub nom " << get_norm(sub_u);
3560
3561 CHKERR TsMonitorSet(ts, step, t, sub_u, ptr->tsCtxPtr.get());
3562 }
3564}
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition TsCtx.cpp:263
constexpr double t
plate stiffness
Definition plate.cpp:58

◆ tsPostProc()

MoFEMErrorCode TSPrePostProc::tsPostProc ( TS  ts)
staticprivate

Post process time step.

Currently this function does not make anything major

Parameters
tsPETSc time stepping object
Returns
MoFEMErrorCode Success/failure code

Called after each time step completion for cleanup operations

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3476 of file free_surface.cpp.

3476 {
3477 if (auto ptr = tsPrePostProc.lock()) {
3478 auto &m_field = ptr->fsRawPtr->mField;
3479 MOFEM_LOG_CHANNEL("SYNC");
3480 MOFEM_TAG_AND_LOG("SYNC", Sev::verbose, "FS") << "PostProc done";
3481 MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::verbose);
3482 }
3483 return 0;
3484}
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.

◆ tsPostStage()

MoFEMErrorCode TSPrePostProc::tsPostStage ( TS  ts,
PetscReal  stagetime,
PetscInt  stageindex,
Vec *  Y 
)
static

[Boundary condition]

Examples
mofem/tutorials/adv-4_dynamic_first_order_con_law/dynamic_first_order_con_law.cpp.

Definition at line 577 of file dynamic_first_order_con_law.cpp.

578 {
580 // cerr << "tsPostStage " <<"\n";
581 if (auto ptr = tsPrePostProc.lock()) {
582 auto &m_field = ptr->fsRawPtr->mField;
583
584 auto fb = m_field.getInterface<FieldBlas>();
585 double dt;
586 CHKERR TSGetTimeStep(ts, &dt);
587 double time;
588 CHKERR TSGetTime(ts, &time);
589 PetscInt num_stages;
590 Vec *stage_solutions;
591
592 CHKERR TSGetStages(ts, &num_stages, &stage_solutions);
593 PetscPrintf(PETSC_COMM_WORLD, "Check timestep %d time %e dt %e\n",
594 num_stages, time, dt);
595
596 const double inv_num_step = (double)num_stages;
597 CHKERR fb->fieldCopy(1., "x_1", "x_2");
598 CHKERR fb->fieldAxpy(dt, "V", "x_2");
599 CHKERR fb->fieldCopy(1., "x_2", "x_1");
600
601 CHKERR fb->fieldCopy(-inv_num_step / dt, "F_0", "F_dot");
602 CHKERR fb->fieldAxpy(inv_num_step / dt, "F", "F_dot");
603 CHKERR fb->fieldCopy(1., "F", "F_0");
604 }
606}
double dt
const FTensor::Tensor2< T, Dim, Dim > Vec
Basic algebra on fields.
Definition FieldBlas.hpp:21

◆ tsPostStep()

MoFEMErrorCode TSPrePostProc::tsPostStep ( TS  ts)
static
Examples
mofem/tutorials/adv-4_dynamic_first_order_con_law/dynamic_first_order_con_law.cpp.

Definition at line 608 of file dynamic_first_order_con_law.cpp.

608 {
610
611 if (auto ptr = tsPrePostProc.lock()) {
612 double dt;
613 CHKERR TSGetTimeStep(ts, &dt);
614 double time;
615 CHKERR TSGetTime(ts, &time);
616 }
618}

◆ tsPreProc()

MoFEMErrorCode TSPrePostProc::tsPreProc ( TS  ts)
staticprivate

Pre process time step.

Refine mesh and update fields before each time step

Parameters
tsPETSc time stepping object
Returns
MoFEMErrorCode Success/failure code

This function is called before each time step to:

  • Apply phase field cutoff constraints
  • Refine mesh based on interface location
  • Project solution data to new mesh
  • Update solver operators

cut-off values at nodes, i.e. abs("H") <= 1

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3361 of file free_surface.cpp.

3361 {
3363
3364 if (auto ptr = tsPrePostProc.lock()) {
3365
3366 /**
3367 * @brief cut-off values at nodes, i.e. abs("H") <= 1
3368 *
3369 */
3370 auto cut_off_dofs = [&]() {
3372
3373 auto &m_field = ptr->fsRawPtr->mField;
3374
3375 Range current_verts;
3376 auto bit_mng = m_field.getInterface<BitRefManager>();
3378 bit(get_current_bit()), BitRefLevel().set(), MBVERTEX, current_verts);
3379
3380 auto cut_off_verts = [&](boost::shared_ptr<FieldEntity> ent_ptr) {
3382 for (auto &h : ent_ptr->getEntFieldData()) {
3383 h = cut_off(h);
3384 }
3386 };
3387
3388 auto field_blas = m_field.getInterface<FieldBlas>();
3389 CHKERR field_blas->fieldLambdaOnEntities(cut_off_verts, "H",
3390 &current_verts);
3392 };
3393
3394 CHKERR cut_off_dofs();
3395 }
3396
3397 if (auto ptr = tsPrePostProc.lock()) {
3398 MOFEM_LOG("FS", Sev::inform) << "Run step pre proc";
3399
3400 auto &m_field = ptr->fsRawPtr->mField;
3401 auto simple = m_field.getInterface<Simple>();
3402
3403 // get vector norm
3404 auto get_norm = [&](auto x) {
3405 double nrm;
3406 CHKERR VecNorm(x, NORM_2, &nrm);
3407 return nrm;
3408 };
3409
3410 // refine problem and project data, including theta data
3411 auto refine_problem = [&]() {
3413 MOFEM_LOG("FS", Sev::inform) << "Refine problem";
3414 CHKERR ptr->fsRawPtr->refineMesh(refine_overlap);
3415 CHKERR ptr->fsRawPtr->projectData();
3417 };
3418
3419 // set new jacobin operator, since problem and thus tangent matrix size has
3420 // changed
3421 auto set_jacobian_operators = [&]() {
3423 ptr->subB = createDMMatrix(ptr->solverSubDM);
3424 CHKERR KSPReset(ptr->subKSP);
3426 };
3427
3428 // set new solution
3429 auto set_solution = [&]() {
3431 MOFEM_LOG("FS", Sev::inform) << "Set solution";
3432
3433 PetscObjectState state;
3434
3435 // Record the state, and set it again. This is to fool PETSc that solution
3436 // vector is not updated. Otherwise PETSc will treat every step as a first
3437 // step.
3438
3439 // globSol is updated as result mesh refinement - this is not really set
3440 // a new solution.
3441
3442 CHKERR PetscObjectStateGet(getPetscObject(ptr->globSol.get()), &state);
3443 CHKERR DMoFEMMeshToLocalVector(simple->getDM(), ptr->globSol,
3444 INSERT_VALUES, SCATTER_FORWARD);
3445 CHKERR PetscObjectStateSet(getPetscObject(ptr->globSol.get()), state);
3446 MOFEM_LOG("FS", Sev::verbose)
3447 << "Set solution, vector norm " << get_norm(ptr->globSol);
3449 };
3450
3451 PetscBool is_theta;
3452 PetscObjectTypeCompare((PetscObject)ts, TSTHETA, &is_theta);
3453 if (is_theta) {
3454
3455 CHKERR refine_problem(); // refine problem
3456 CHKERR set_jacobian_operators(); // set new jacobian
3457 CHKERR set_solution(); // set solution
3458
3459 } else {
3460 SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_IMPLEMENTED,
3461 "Sorry, only TSTheta handling is implemented");
3462 }
3463
3464 // Need barriers, some functions in TS solver need are called collectively
3465 // and require the same state of variables
3466 PetscBarrier((PetscObject)ts);
3467
3468 MOFEM_LOG_CHANNEL("SYNC");
3469 MOFEM_TAG_AND_LOG("SYNC", Sev::verbose, "FS") << "PreProc done";
3470 MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::verbose);
3471 }
3472
3474}
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition acoustic.cpp:69
@ MOFEM_NOT_IMPLEMENTED
Definition definitions.h:32
double h
int refine_overlap
auto cut_off
Phase field cutoff function.
auto get_current_bit
dofs bit used to do calculations
auto bit
Create bit reference level.
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode, RowColData rc=RowColData::COL)
set local (or ghosted) vector values on mesh for partition only
Definition DMMoFEM.cpp:543
auto createDMMatrix(DM dm)
Get smart matrix from DM.
Definition DMMoFEM.hpp:1207
MoFEMErrorCode getEntitiesByTypeAndRefLevel(const BitRefLevel bit, const BitRefLevel mask, const EntityType type, const EntityHandle meshset, int verb=0) const
add all ents from ref level given by bit to meshset
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition Types.hpp:40
PetscObject getPetscObject(T obj)
Managing BitRefLevels.
MoFEMErrorCode fieldLambdaOnEntities(OneFieldFunctionOnEntities lambda, const std::string field_name, Range *ents_ptr=nullptr)
field lambda
Definition FieldBlas.cpp:50
Simple interface for fast problem set-up.
Definition Simple.hpp:27

◆ tsPreStage()

static MoFEMErrorCode TSPrePostProc::tsPreStage ( TS  ts)
staticprivate

Pre-stage processing for time stepping.

Parameters
tsPETSc time stepping object
Returns
MoFEMErrorCode Success/failure code
Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

◆ tsPreStep()

MoFEMErrorCode TSPrePostProc::tsPreStep ( TS  ts)
static
Examples
mofem/tutorials/adv-4_dynamic_first_order_con_law/dynamic_first_order_con_law.cpp.

Definition at line 620 of file dynamic_first_order_con_law.cpp.

620 {
622
623 if (auto ptr = tsPrePostProc.lock()) {
624 double dt;
625 CHKERR TSGetTimeStep(ts, &dt);
626 double time;
627 CHKERR TSGetTime(ts, &time);
628 int step_num;
629 CHKERR TSGetStepNumber(ts, &step_num);
630 }
632}

◆ tsSetIFunction()

MoFEMErrorCode TSPrePostProc::tsSetIFunction ( TS  ts,
PetscReal  t,
Vec  u,
Vec  u_t,
Vec  f,
void *  ctx 
)
staticprivate

Set implicit function for time stepping.

Parameters
tsPETSc time stepping object
tCurrent time
uSolution vector at current time
u_tTime derivative of solution
fOutput residual vector F(t,u,u_t)
ctxUser context (unused)
Returns
MoFEMErrorCode Success/failure code

Wrapper that scatters global vectors to sub-problem and evaluates residual

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3486 of file free_surface.cpp.

3487 {
3489 if (auto ptr = tsPrePostProc.lock()) {
3490 auto sub_u = ptr->getSubVector();
3491 auto sub_u_t = vectorDuplicate(sub_u);
3492 auto sub_f = vectorDuplicate(sub_u);
3493 auto scatter = ptr->getScatter(sub_u, u, R);
3494
3495 auto apply_scatter_and_update = [&](auto x, auto sub_x) {
3497 CHKERR VecScatterBegin(scatter, x, sub_x, INSERT_VALUES, SCATTER_REVERSE);
3498 CHKERR VecScatterEnd(scatter, x, sub_x, INSERT_VALUES, SCATTER_REVERSE);
3499 CHKERR VecGhostUpdateBegin(sub_x, INSERT_VALUES, SCATTER_FORWARD);
3500 CHKERR VecGhostUpdateEnd(sub_x, INSERT_VALUES, SCATTER_FORWARD);
3502 };
3503
3504 CHKERR apply_scatter_and_update(u, sub_u);
3505 CHKERR apply_scatter_and_update(u_t, sub_u_t);
3506
3507 CHKERR TsSetIFunction(ts, t, sub_u, sub_u_t, sub_f, ptr->tsCtxPtr.get());
3508 CHKERR VecScatterBegin(scatter, sub_f, f, INSERT_VALUES, SCATTER_FORWARD);
3509 CHKERR VecScatterEnd(scatter, sub_f, f, INSERT_VALUES, SCATTER_FORWARD);
3510 }
3512}
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition TsCtx.cpp:56

◆ tsSetIJacobian()

MoFEMErrorCode TSPrePostProc::tsSetIJacobian ( TS  ts,
PetscReal  t,
Vec  u,
Vec  u_t,
PetscReal  a,
Mat  A,
Mat  B,
void *  ctx 
)
staticprivate

Set implicit Jacobian for time stepping.

Parameters
tsPETSc time stepping object
tCurrent time
uSolution vector at current time
u_tTime derivative of solution
aShift parameter for implicit methods
AJacobian matrix (input)
BJacobian matrix (output, often same as A)
ctxUser context (unused)
Returns
MoFEMErrorCode Success/failure code

Wrapper that assembles Jacobian matrix for sub-problem Wrapper for SNES Lhs

Examples
mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3514 of file free_surface.cpp.

3516 {
3518 if (auto ptr = tsPrePostProc.lock()) {
3519 auto sub_u = ptr->getSubVector();
3520 auto sub_u_t = vectorDuplicate(sub_u);
3521 auto scatter = ptr->getScatter(sub_u, u, R);
3522
3523 auto apply_scatter_and_update = [&](auto x, auto sub_x) {
3525 CHKERR VecScatterBegin(scatter, x, sub_x, INSERT_VALUES, SCATTER_REVERSE);
3526 CHKERR VecScatterEnd(scatter, x, sub_x, INSERT_VALUES, SCATTER_REVERSE);
3527 CHKERR VecGhostUpdateBegin(sub_x, INSERT_VALUES, SCATTER_FORWARD);
3528 CHKERR VecGhostUpdateEnd(sub_x, INSERT_VALUES, SCATTER_FORWARD);
3530 };
3531
3532 CHKERR apply_scatter_and_update(u, sub_u);
3533 CHKERR apply_scatter_and_update(u_t, sub_u_t);
3534
3535 CHKERR TsSetIJacobian(ts, t, sub_u, sub_u_t, a, ptr->subB, ptr->subB,
3536 ptr->tsCtxPtr.get());
3537 }
3539}
constexpr double a
PetscErrorCode TsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Set function evaluating jacobian in TS solver.
Definition TsCtx.cpp:169

◆ tsSetUp() [1/2]

MoFEMErrorCode TSPrePostProc::tsSetUp ( TS  ts)

Used to setup TS solver.

Parameters
ts
Returns
MoFEMErrorCode
Examples
mofem/tutorials/adv-4_dynamic_first_order_con_law/dynamic_first_order_con_law.cpp, and mofem/tutorials/vec-5_free_surface/free_surface.cpp.

Definition at line 3620 of file free_surface.cpp.

3620 {
3621
3622 auto &m_field = fsRawPtr->mField;
3623 auto simple = m_field.getInterface<Simple>();
3624
3626
3627 auto dm = simple->getDM();
3628
3630 CHKERR TSSetIFunction(ts, globRes, tsSetIFunction, nullptr);
3631 CHKERR TSSetIJacobian(ts, PETSC_NULLPTR, PETSC_NULLPTR, tsSetIJacobian, nullptr);
3632 CHKERR TSMonitorSet(ts, tsMonitor, fsRawPtr, PETSC_NULLPTR);
3633
3634 SNES snes;
3635 CHKERR TSGetSNES(ts, &snes);
3636 auto snes_ctx_ptr = getDMSnesCtx(dm);
3637
3638 auto set_section_monitor = [&](auto snes) {
3640 CHKERR SNESMonitorSet(snes,
3641 (MoFEMErrorCode(*)(SNES, PetscInt, PetscReal,
3642 void *))MoFEMSNESMonitorFields,
3643 (void *)(snes_ctx_ptr.get()), nullptr);
3645 };
3646
3647 CHKERR set_section_monitor(snes);
3648
3649 auto ksp = createKSP(m_field.get_comm());
3650 CHKERR KSPSetType(ksp, KSPPREONLY); // Run KSP internally in ShellPC
3651 auto sub_pc = createPC(fsRawPtr->mField.get_comm());
3652 CHKERR PCSetType(sub_pc, PCSHELL);
3653 CHKERR PCShellSetSetUp(sub_pc, pcSetup);
3654 CHKERR PCShellSetApply(sub_pc, pcApply);
3655 CHKERR KSPSetPC(ksp, sub_pc);
3656 CHKERR SNESSetKSP(snes, ksp);
3657
3660
3661 CHKERR TSSetPreStep(ts, tsPreProc);
3662 CHKERR TSSetPostStep(ts, tsPostProc);
3663
3665}
SmartPetscObj< KSP > ksp
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
auto getDMTsCtx(DM dm)
Get TS context data structure used by DM.
Definition DMMoFEM.hpp:1292
MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *ctx)
Sens monitor printing residual field by field.
Definition SnesCtx.cpp:600
auto createPC(MPI_Comm comm)
auto getDMSnesCtx(DM dm)
Get SNES context data structure used by DM.
Definition DMMoFEM.hpp:1278
MoFEM::Interface & mField
Reference to MoFEM interface.
Definition plastic.cpp:227
virtual MPI_Comm & get_comm() const =0
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition Simple.cpp:799
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
boost::shared_ptr< TsCtx > tsCtxPtr
SmartPetscObj< Vec > globRes
SmartPetscObj< Vec > globSol
static MoFEMErrorCode tsMonitor(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Monitor solution during time stepping.
static MoFEMErrorCode pcApply(PC pc, Vec pc_f, Vec pc_x)
Apply preconditioner.
static MoFEMErrorCode tsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Set implicit Jacobian for time stepping.
boost::shared_ptr< SnesCtx > snesCtxPtr
static MoFEMErrorCode tsPreProc(TS ts)
Pre process time step.
static MoFEMErrorCode tsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec f, void *ctx)
Set implicit function for time stepping.
static MoFEMErrorCode pcSetup(PC pc)
Setup preconditioner.
static MoFEMErrorCode tsPostProc(TS ts)
Post process time step.

◆ tsSetUp() [2/2]

MoFEMErrorCode TSPrePostProc::tsSetUp ( TS  ts)

Used to setup TS solver.

Parameters
tsPETSc time stepping object to configure
Returns
MoFEMErrorCode Success/failure code

Sets up time stepping solver with custom preconditioner, monitors, and callback functions for mesh refinement and solution projection

Member Data Documentation

◆ fsRawPtr [1/2]

Example* TSPrePostProc::fsRawPtr

◆ fsRawPtr [2/2]

FreeSurface* TSPrePostProc::fsRawPtr

Definition at line 678 of file free_surface.cpp.

◆ globRes

SmartPetscObj<Vec> TSPrePostProc::globRes
private

◆ globSol

SmartPetscObj<Vec> TSPrePostProc::globSol

◆ snesCtxPtr

boost::shared_ptr<SnesCtx> TSPrePostProc::snesCtxPtr
private

◆ solverSubDM

SmartPetscObj<DM> TSPrePostProc::solverSubDM

◆ subB

SmartPetscObj<Mat> TSPrePostProc::subB
private

◆ subKSP

SmartPetscObj<KSP> TSPrePostProc::subKSP
private

◆ tsCtxPtr

boost::shared_ptr<TsCtx> TSPrePostProc::tsCtxPtr
private

The documentation for this struct was generated from the following files: