v0.14.0
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. More...
 
 TSPrePostProc ()=default
 
virtual ~TSPrePostProc ()=default
 
MoFEMErrorCode tsSetUp (TS ts)
 Used to setup TS solver. More...
 
SmartPetscObj< VecScatter > getScatter (Vec x, Vec y, enum FR fr)
 
SmartPetscObj< Vec > getSubVector ()
 

Static Public Member Functions

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

Public Attributes

ExamplefsRawPtr
 
SmartPetscObj< DM > solverSubDM
 
SmartPetscObj< Vec > globSol
 
FreeSurfacefsRawPtr
 

Static Private Member Functions

static MoFEMErrorCode tsPreProc (TS ts)
 Pre process time step. More...
 
static MoFEMErrorCode tsPostProc (TS ts)
 Post process time step. More...
 
static MoFEMErrorCode tsPreStage (TS ts)
 
static MoFEMErrorCode tsSetIFunction (TS ts, PetscReal t, Vec u, Vec u_t, Vec f, void *ctx)
 
static MoFEMErrorCode tsSetIJacobian (TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
 Wrapper for SNES Lhs. More...
 
static MoFEMErrorCode tsMonitor (TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
 Wrapper for TS monitor. More...
 
static MoFEMErrorCode pcSetup (PC pc)
 
static MoFEMErrorCode pcApply (PC pc, Vec pc_f, Vec pc_x)
 

Private Attributes

SmartPetscObj< Vec > globRes
 
SmartPetscObj< Mat > subB
 
SmartPetscObj< KSP > subKSP
 
boost::shared_ptr< SnesCtxsnesCtxPtr
 
boost::shared_ptr< TsCtxtsCtxPtr
 

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
dynamic_first_order_con_law.cpp, and free_surface.cpp.

Definition at line 384 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 
)
Examples
free_surface.cpp.

Definition at line 3118 of file free_surface.cpp.

3118  {
3119  if (auto ptr = tsPrePostProc.lock()) {
3120  auto prb_ptr = ptr->fsRawPtr->mField.get_problem("SUB_SOLVER");
3121  if (auto sub_data = prb_ptr->getSubData()) {
3122  auto is = sub_data->getSmartColIs();
3123  VecScatter s;
3124  if (fr == R) {
3125  CHK_THROW_MESSAGE(VecScatterCreate(x, PETSC_NULL, y, is, &s),
3126  "crate scatter");
3127  } else {
3128  CHK_THROW_MESSAGE(VecScatterCreate(x, is, y, PETSC_NULL, &s),
3129  "crate scatter");
3130  }
3131  return SmartPetscObj<VecScatter>(s);
3132  }
3133  }
3134  CHK_THROW_MESSAGE(MOFEM_DATA_INCONSISTENCY, "No prb pinter");
3135  return SmartPetscObj<VecScatter>();
3136 }

◆ getSubVector()

SmartPetscObj< Vec > TSPrePostProc::getSubVector ( )
Examples
free_surface.cpp.

Definition at line 3138 of file free_surface.cpp.

3138  {
3139  return createDMVector(solverSubDM);
3140 }

◆ pcApply()

MoFEMErrorCode TSPrePostProc::pcApply ( PC  pc,
Vec  pc_f,
Vec  pc_x 
)
staticprivate
Examples
free_surface.cpp.

Definition at line 3098 of file free_surface.cpp.

3098  {
3100  if (auto ptr = tsPrePostProc.lock()) {
3101  auto sub_x = ptr->getSubVector();
3102  auto sub_f = vectorDuplicate(sub_x);
3103  auto scatter = ptr->getScatter(sub_x, pc_x, R);
3104  CHKERR VecScatterBegin(scatter, pc_f, sub_f, INSERT_VALUES,
3105  SCATTER_REVERSE);
3106  CHKERR VecScatterEnd(scatter, pc_f, sub_f, INSERT_VALUES, SCATTER_REVERSE);
3107  CHKERR KSPSetOperators(ptr->subKSP, ptr->subB, ptr->subB);
3108  MOFEM_LOG("FS", Sev::verbose) << "PCShell solve";
3109  CHKERR KSPSolve(ptr->subKSP, sub_f, sub_x);
3110  MOFEM_LOG("FS", Sev::verbose) << "PCShell solve <- done";
3111  CHKERR VecScatterBegin(scatter, sub_x, pc_x, INSERT_VALUES,
3112  SCATTER_FORWARD);
3113  CHKERR VecScatterEnd(scatter, sub_x, pc_x, INSERT_VALUES, SCATTER_FORWARD);
3114  }
3116 };

◆ pcSetup()

MoFEMErrorCode TSPrePostProc::pcSetup ( PC  pc)
staticprivate
Examples
free_surface.cpp.

Definition at line 3088 of file free_surface.cpp.

3088  {
3090  if (auto ptr = tsPrePostProc.lock()) {
3091  MOFEM_LOG("FS", Sev::verbose) << "SetUP sub PC";
3092  ptr->subKSP = createKSP(ptr->fsRawPtr->mField.get_comm());
3093  CHKERR KSPSetFromOptions(ptr->subKSP);
3094  }
3096 };

◆ tsMonitor()

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

Wrapper for TS monitor.

Examples
free_surface.cpp.

Definition at line 3063 of file free_surface.cpp.

3064  {
3066  if (auto ptr = tsPrePostProc.lock()) {
3067  auto get_norm = [&](auto x) {
3068  double nrm;
3069  CHKERR VecNorm(x, NORM_2, &nrm);
3070  return nrm;
3071  };
3072 
3073  auto sub_u = ptr->getSubVector();
3074  auto scatter = ptr->getScatter(sub_u, u, R);
3075  CHKERR VecScatterBegin(scatter, u, sub_u, INSERT_VALUES, SCATTER_REVERSE);
3076  CHKERR VecScatterEnd(scatter, u, sub_u, INSERT_VALUES, SCATTER_REVERSE);
3077  CHKERR VecGhostUpdateBegin(sub_u, INSERT_VALUES, SCATTER_FORWARD);
3078  CHKERR VecGhostUpdateEnd(sub_u, INSERT_VALUES, SCATTER_FORWARD);
3079 
3080  MOFEM_LOG("FS", Sev::verbose)
3081  << "u norm " << get_norm(u) << " u sub nom " << get_norm(sub_u);
3082 
3083  CHKERR TsMonitorSet(ts, step, t, sub_u, ptr->tsCtxPtr.get());
3084  }
3086 }

◆ tsPostProc()

MoFEMErrorCode TSPrePostProc::tsPostProc ( TS  ts)
staticprivate

Post process time step.

Currently that function do not make anything major

Parameters
ts
Returns
MoFEMErrorCode
Examples
free_surface.cpp.

Definition at line 2998 of file free_surface.cpp.

2998  {
2999  if (auto ptr = tsPrePostProc.lock()) {
3000  auto &m_field = ptr->fsRawPtr->mField;
3001  MOFEM_LOG_CHANNEL("SYNC");
3002  MOFEM_TAG_AND_LOG("SYNC", Sev::verbose, "FS") << "PostProc done";
3003  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::verbose);
3004  }
3005  return 0;
3006 }

◆ tsPostStage()

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

[Boundary condition]

Examples
dynamic_first_order_con_law.cpp.

Definition at line 581 of file dynamic_first_order_con_law.cpp.

582  {
584  // cerr << "tsPostStage " <<"\n";
585  if (auto ptr = tsPrePostProc.lock()) {
586  auto &m_field = ptr->fsRawPtr->mField;
587  auto *simple = m_field.getInterface<Simple>();
588  auto *pipeline_mng = m_field.getInterface<PipelineManager>();
589 
590  auto fb = m_field.getInterface<FieldBlas>();
591  double dt;
592  CHKERR TSGetTimeStep(ts, &dt);
593  double time;
594  CHKERR TSGetTime(ts, &time);
595  PetscInt num_stages;
596  Vec *stage_solutions;
597 
598  CHKERR TSGetStages(ts, &num_stages, &stage_solutions);
599  PetscPrintf(PETSC_COMM_WORLD, "Check timestep %d time %e dt %e\n",
600  num_stages, time, dt);
601 
602  const double inv_num_step = (double)num_stages;
603  CHKERR fb->fieldCopy(1., "x_1", "x_2");
604  CHKERR fb->fieldAxpy(dt, "V", "x_2");
605  CHKERR fb->fieldCopy(1., "x_2", "x_1");
606 
607  CHKERR fb->fieldCopy(-inv_num_step / dt, "F_0", "F_dot");
608  CHKERR fb->fieldAxpy(inv_num_step / dt, "F", "F_dot");
609  CHKERR fb->fieldCopy(1., "F", "F_0");
610  }
612 }

◆ tsPostStep()

MoFEMErrorCode TSPrePostProc::tsPostStep ( TS  ts)
static
Examples
dynamic_first_order_con_law.cpp.

Definition at line 614 of file dynamic_first_order_con_law.cpp.

614  {
616 
617  if (auto ptr = tsPrePostProc.lock()) {
618  auto &m_field = ptr->fsRawPtr->mField;
619  auto fb = m_field.getInterface<FieldBlas>();
620  double dt;
621  CHKERR TSGetTimeStep(ts, &dt);
622  double time;
623  CHKERR TSGetTime(ts, &time);
624  }
626 }

◆ tsPreProc()

MoFEMErrorCode TSPrePostProc::tsPreProc ( TS  ts)
staticprivate

Pre process time step.

Refine mesh and update fields

Parameters
ts
Returns
MoFEMErrorCode

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

Examples
free_surface.cpp.

Definition at line 2877 of file free_surface.cpp.

2877  {
2879 
2880  if (auto ptr = tsPrePostProc.lock()) {
2881 
2882  /**
2883  * @brief cut-off values at nodes, i.e. abs("H") <= 1
2884  *
2885  */
2886  auto cut_off_dofs = [&]() {
2888 
2889  auto &m_field = ptr->fsRawPtr->mField;
2890 
2891  Range current_verts;
2892  auto bit_mng = m_field.getInterface<BitRefManager>();
2893  CHKERR bit_mng->getEntitiesByTypeAndRefLevel(
2894  bit(get_current_bit()), BitRefLevel().set(), MBVERTEX, current_verts);
2895 
2896  auto cut_off_verts = [&](boost::shared_ptr<FieldEntity> ent_ptr) {
2898  for (auto &h : ent_ptr->getEntFieldData()) {
2899  h = cut_off(h);
2900  }
2902  };
2903 
2904  auto field_blas = m_field.getInterface<FieldBlas>();
2905  CHKERR field_blas->fieldLambdaOnEntities(cut_off_verts, "H",
2906  &current_verts);
2908  };
2909 
2910  CHKERR cut_off_dofs();
2911  }
2912 
2913  if (auto ptr = tsPrePostProc.lock()) {
2914  MOFEM_LOG("FS", Sev::inform) << "Run step pre proc";
2915 
2916  auto &m_field = ptr->fsRawPtr->mField;
2917  auto simple = m_field.getInterface<Simple>();
2918  auto bit_mng = m_field.getInterface<BitRefManager>();
2919  auto bc_mng = m_field.getInterface<BcManager>();
2920  auto field_blas = m_field.getInterface<FieldBlas>();
2921  auto opt = m_field.getInterface<OperatorsTester>();
2922  auto pip_mng = m_field.getInterface<PipelineManager>();
2923  auto prb_mng = m_field.getInterface<ProblemsManager>();
2924 
2925  // get vector norm
2926  auto get_norm = [&](auto x) {
2927  double nrm;
2928  CHKERR VecNorm(x, NORM_2, &nrm);
2929  return nrm;
2930  };
2931 
2932  // refine problem and project data, including theta data
2933  auto refine_problem = [&]() {
2935  MOFEM_LOG("FS", Sev::inform) << "Refine problem";
2936  CHKERR ptr->fsRawPtr->refineMesh(refine_overlap);
2937  CHKERR ptr->fsRawPtr->projectData();
2939  };
2940 
2941  // set new jacobin operator, since problem and thus tangent matrix size has
2942  // changed
2943  auto set_jacobian_operators = [&]() {
2945  ptr->subB = createDMMatrix(ptr->solverSubDM);
2946  CHKERR KSPReset(ptr->subKSP);
2948  };
2949 
2950  // set new solution
2951  auto set_solution = [&]() {
2953  MOFEM_LOG("FS", Sev::inform) << "Set solution";
2954 
2955  PetscObjectState state;
2956 
2957  // Record the state, and set it again. This is to fool PETSc that solution
2958  // vector is not updated. Otherwise PETSc will treat every step as a first
2959  // step.
2960 
2961  // globSol is updated as result mesh refinement - this is not really set
2962  // a new solution.
2963 
2964  CHKERR PetscObjectStateGet(getPetscObject(ptr->globSol.get()), &state);
2965  CHKERR DMoFEMMeshToLocalVector(simple->getDM(), ptr->globSol,
2966  INSERT_VALUES, SCATTER_FORWARD);
2967  CHKERR PetscObjectStateSet(getPetscObject(ptr->globSol.get()), state);
2968  MOFEM_LOG("FS", Sev::verbose)
2969  << "Set solution, vector norm " << get_norm(ptr->globSol);
2971  };
2972 
2973  PetscBool is_theta;
2974  PetscObjectTypeCompare((PetscObject)ts, TSTHETA, &is_theta);
2975  if (is_theta) {
2976 
2977  CHKERR refine_problem(); // refine problem
2978  CHKERR set_jacobian_operators(); // set new jacobian
2979  CHKERR set_solution(); // set solution
2980 
2981  } else {
2982  SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_IMPLEMENTED,
2983  "Sorry, only TSTheta handling is implemented");
2984  }
2985 
2986  // Need barriers, somme functions in TS solver need are called collectively
2987  // and requite the same state of variables
2988  PetscBarrier((PetscObject)ts);
2989 
2990  MOFEM_LOG_CHANNEL("SYNC");
2991  MOFEM_TAG_AND_LOG("SYNC", Sev::verbose, "FS") << "PreProc done";
2992  MOFEM_LOG_SEVERITY_SYNC(m_field.get_comm(), Sev::verbose);
2993  }
2994 
2996 }

◆ tsPreStage()

static MoFEMErrorCode TSPrePostProc::tsPreStage ( TS  ts)
staticprivate

◆ tsPreStep()

MoFEMErrorCode TSPrePostProc::tsPreStep ( TS  ts)
static
Examples
dynamic_first_order_con_law.cpp.

Definition at line 628 of file dynamic_first_order_con_law.cpp.

628  {
630 
631  if (auto ptr = tsPrePostProc.lock()) {
632  auto &m_field = ptr->fsRawPtr->mField;
633  auto *simple = m_field.getInterface<Simple>();
634  auto *pipeline_mng = m_field.getInterface<PipelineManager>();
635 
636  double dt;
637  CHKERR TSGetTimeStep(ts, &dt);
638  double time;
639  CHKERR TSGetTime(ts, &time);
640  int step_num;
641  CHKERR TSGetStepNumber(ts, &step_num);
642  }
644 }

◆ tsSetIFunction()

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

Definition at line 3008 of file free_surface.cpp.

3009  {
3011  if (auto ptr = tsPrePostProc.lock()) {
3012  auto sub_u = ptr->getSubVector();
3013  auto sub_u_t = vectorDuplicate(sub_u);
3014  auto sub_f = vectorDuplicate(sub_u);
3015  auto scatter = ptr->getScatter(sub_u, u, R);
3016 
3017  auto apply_scatter_and_update = [&](auto x, auto sub_x) {
3019  CHKERR VecScatterBegin(scatter, x, sub_x, INSERT_VALUES, SCATTER_REVERSE);
3020  CHKERR VecScatterEnd(scatter, x, sub_x, INSERT_VALUES, SCATTER_REVERSE);
3021  CHKERR VecGhostUpdateBegin(sub_x, INSERT_VALUES, SCATTER_FORWARD);
3022  CHKERR VecGhostUpdateEnd(sub_x, INSERT_VALUES, SCATTER_FORWARD);
3024  };
3025 
3026  CHKERR apply_scatter_and_update(u, sub_u);
3027  CHKERR apply_scatter_and_update(u_t, sub_u_t);
3028 
3029  CHKERR TsSetIFunction(ts, t, sub_u, sub_u_t, sub_f, ptr->tsCtxPtr.get());
3030  CHKERR VecScatterBegin(scatter, sub_f, f, INSERT_VALUES, SCATTER_FORWARD);
3031  CHKERR VecScatterEnd(scatter, sub_f, f, INSERT_VALUES, SCATTER_FORWARD);
3032  }
3034 }

◆ tsSetIJacobian()

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

Wrapper for SNES Lhs.

Examples
free_surface.cpp.

Definition at line 3036 of file free_surface.cpp.

3038  {
3040  if (auto ptr = tsPrePostProc.lock()) {
3041  auto sub_u = ptr->getSubVector();
3042  auto sub_u_t = vectorDuplicate(sub_u);
3043  auto scatter = ptr->getScatter(sub_u, u, R);
3044 
3045  auto apply_scatter_and_update = [&](auto x, auto sub_x) {
3047  CHKERR VecScatterBegin(scatter, x, sub_x, INSERT_VALUES, SCATTER_REVERSE);
3048  CHKERR VecScatterEnd(scatter, x, sub_x, INSERT_VALUES, SCATTER_REVERSE);
3049  CHKERR VecGhostUpdateBegin(sub_x, INSERT_VALUES, SCATTER_FORWARD);
3050  CHKERR VecGhostUpdateEnd(sub_x, INSERT_VALUES, SCATTER_FORWARD);
3052  };
3053 
3054  CHKERR apply_scatter_and_update(u, sub_u);
3055  CHKERR apply_scatter_and_update(u_t, sub_u_t);
3056 
3057  CHKERR TsSetIJacobian(ts, t, sub_u, sub_u_t, a, ptr->subB, ptr->subB,
3058  ptr->tsCtxPtr.get());
3059  }
3061 }

◆ tsSetUp() [1/2]

MoFEMErrorCode TSPrePostProc::tsSetUp ( TS  ts)

Used to setup TS solver.

Parameters
ts
Returns
MoFEMErrorCode
Examples
free_surface.cpp.

Definition at line 3142 of file free_surface.cpp.

3142  {
3143 
3144  auto &m_field = fsRawPtr->mField;
3145  auto simple = m_field.getInterface<Simple>();
3146 
3148 
3149  auto dm = simple->getDM();
3150 
3152  CHKERR TSSetIFunction(ts, globRes, tsSetIFunction, nullptr);
3153  CHKERR TSSetIJacobian(ts, PETSC_NULL, PETSC_NULL, tsSetIJacobian, nullptr);
3154  CHKERR TSMonitorSet(ts, tsMonitor, fsRawPtr, PETSC_NULL);
3155 
3156  SNES snes;
3157  CHKERR TSGetSNES(ts, &snes);
3158  auto snes_ctx_ptr = getDMSnesCtx(dm);
3159 
3160  auto set_section_monitor = [&](auto snes) {
3162  CHKERR SNESMonitorSet(snes,
3163  (MoFEMErrorCode(*)(SNES, PetscInt, PetscReal,
3164  void *))MoFEMSNESMonitorFields,
3165  (void *)(snes_ctx_ptr.get()), nullptr);
3167  };
3168 
3169  CHKERR set_section_monitor(snes);
3170 
3171  auto ksp = createKSP(m_field.get_comm());
3172  CHKERR KSPSetType(ksp, KSPPREONLY); // Run KSP internally in ShellPC
3173  auto sub_pc = createPC(fsRawPtr->mField.get_comm());
3174  CHKERR PCSetType(sub_pc, PCSHELL);
3175  CHKERR PCShellSetSetUp(sub_pc, pcSetup);
3176  CHKERR PCShellSetApply(sub_pc, pcApply);
3177  CHKERR KSPSetPC(ksp, sub_pc);
3178  CHKERR SNESSetKSP(snes, ksp);
3179 
3182 
3183  CHKERR TSSetPreStep(ts, tsPreProc);
3184  CHKERR TSSetPostStep(ts, tsPostProc);
3185 
3187 }

◆ tsSetUp() [2/2]

MoFEMErrorCode TSPrePostProc::tsSetUp ( TS  ts)

Used to setup TS solver.

Parameters
ts
Returns
MoFEMErrorCode

Member Data Documentation

◆ fsRawPtr [1/2]

Example* TSPrePostProc::fsRawPtr
Examples
free_surface.cpp.

Definition at line 398 of file dynamic_first_order_con_law.cpp.

◆ fsRawPtr [2/2]

FreeSurface* TSPrePostProc::fsRawPtr

Definition at line 421 of file free_surface.cpp.

◆ globRes

SmartPetscObj<Vec> TSPrePostProc::globRes
private
Examples
free_surface.cpp.

Definition at line 457 of file free_surface.cpp.

◆ globSol

SmartPetscObj<Vec> TSPrePostProc::globSol
Examples
free_surface.cpp.

Definition at line 420 of file free_surface.cpp.

◆ snesCtxPtr

boost::shared_ptr<SnesCtx> TSPrePostProc::snesCtxPtr
private
Examples
free_surface.cpp.

Definition at line 462 of file free_surface.cpp.

◆ solverSubDM

SmartPetscObj<DM> TSPrePostProc::solverSubDM
Examples
free_surface.cpp.

Definition at line 419 of file free_surface.cpp.

◆ subB

SmartPetscObj<Mat> TSPrePostProc::subB
private

Definition at line 458 of file free_surface.cpp.

◆ subKSP

SmartPetscObj<KSP> TSPrePostProc::subKSP
private

Definition at line 459 of file free_surface.cpp.

◆ tsCtxPtr

boost::shared_ptr<TsCtx> TSPrePostProc::tsCtxPtr
private
Examples
free_surface.cpp.

Definition at line 464 of file free_surface.cpp.


The documentation for this struct was generated from the following files:
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
TSPrePostProc::solverSubDM
SmartPetscObj< DM > solverSubDM
Definition: free_surface.cpp:419
MOFEM_LOG_SEVERITY_SYNC
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
Definition: LogManager.hpp:352
MOFEM_LOG_CHANNEL
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:284
MoFEM::ProblemsManager
Problem manager is used to build and partition problems.
Definition: ProblemsManager.hpp:21
cut_off
auto cut_off
Definition: free_surface.cpp:203
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
MoFEM::TsSetIFunction
PetscErrorCode TsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec F, void *ctx)
Set IFunction for TS solver.
Definition: TsCtx.cpp:56
CHK_THROW_MESSAGE
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:609
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::getPetscObject
PetscObject getPetscObject(T obj)
Definition: PetscSmartObj.hpp:9
MoFEM::PipelineManager
PipelineManager interface.
Definition: PipelineManager.hpp:24
MoFEM::DMoFEMMeshToLocalVector
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMoFEM.cpp:523
MoFEM::createKSP
auto createKSP(MPI_Comm comm)
Definition: PetscSmartObj.hpp:261
tsPrePostProc
static boost::weak_ptr< TSPrePostProc > tsPrePostProc
Definition: dynamic_first_order_con_law.cpp:405
MoFEM::Simple
Simple interface for fast problem set-up.
Definition: Simple.hpp:27
MoFEM::createDMMatrix
auto createDMMatrix(DM dm)
Get smart matrix from DM.
Definition: DMMoFEM.hpp:1056
h
double h
Definition: free_surface.cpp:169
MoFEM::TsMonitorSet
PetscErrorCode TsMonitorSet(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Set monitor for TS solver.
Definition: TsCtx.cpp:259
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::createDMVector
auto createDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:1099
TSPrePostProc::tsMonitor
static MoFEMErrorCode tsMonitor(TS ts, PetscInt step, PetscReal t, Vec u, void *ctx)
Wrapper for TS monitor.
Definition: free_surface.cpp:3063
a
constexpr double a
Definition: approx_sphere.cpp:30
MoFEM::BcManager
Simple interface for fast problem set-up.
Definition: BcManager.hpp:25
R
@ R
Definition: free_surface.cpp:394
get_current_bit
auto get_current_bit
dofs bit used to do calculations
Definition: free_surface.cpp:149
simple
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
double
TSPrePostProc::tsSetIFunction
static MoFEMErrorCode tsSetIFunction(TS ts, PetscReal t, Vec u, Vec u_t, Vec f, void *ctx)
Definition: free_surface.cpp:3008
TSPrePostProc::tsPostProc
static MoFEMErrorCode tsPostProc(TS ts)
Post process time step.
Definition: free_surface.cpp:2998
MoFEM::createPC
auto createPC(MPI_Comm comm)
Definition: PetscSmartObj.hpp:267
bit
auto bit
Definition: free_surface.cpp:312
MoFEM::TsSetIJacobian
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:165
t
constexpr double t
plate stiffness
Definition: plate.cpp:58
TSPrePostProc::snesCtxPtr
boost::shared_ptr< SnesCtx > snesCtxPtr
Definition: free_surface.cpp:462
MoFEM::OperatorsTester
Calculate directional derivative of the right hand side and compare it with tangent matrix derivative...
Definition: OperatorsTester.hpp:21
tsPrePostProc
static boost::weak_ptr< TSPrePostProc > tsPrePostProc
Definition: free_surface.cpp:467
TSPrePostProc::tsPreProc
static MoFEMErrorCode tsPreProc(TS ts)
Pre process time step.
Definition: free_surface.cpp:2877
Range
MOFEM_TAG_AND_LOG
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:362
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
TSPrePostProc::tsCtxPtr
boost::shared_ptr< TsCtx > tsCtxPtr
Definition: free_surface.cpp:464
MoFEM::vectorDuplicate
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
Definition: PetscSmartObj.hpp:221
TSPrePostProc::tsSetIJacobian
static MoFEMErrorCode tsSetIJacobian(TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a, Mat A, Mat B, void *ctx)
Wrapper for SNES Lhs.
Definition: free_surface.cpp:3036
HenckyOps::f
auto f
Definition: HenckyOps.hpp:15
refine_overlap
int refine_overlap
Definition: free_surface.cpp:142
EigenMatrix::Vec
const FTensor::Tensor2< T, Dim, Dim > Vec
Definition: MatrixFunction.hpp:66
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::MoFEMSNESMonitorFields
MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *snes_ctx)
Sens monitor printing residual field by field.
Definition: SnesCtx.cpp:232
TSPrePostProc::globSol
SmartPetscObj< Vec > globSol
Definition: free_surface.cpp:420
Example::mField
MoFEM::Interface & mField
Definition: plastic.cpp:184
TSPrePostProc::pcSetup
static MoFEMErrorCode pcSetup(PC pc)
Definition: free_surface.cpp:3088
MoFEM::BitRefManager
Managing BitRefLevels.
Definition: BitRefManager.hpp:21
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
MoFEM::getDMSnesCtx
auto getDMSnesCtx(DM dm)
Get SNES context data structure used by DM.
Definition: DMMoFEM.hpp:1127
dt
double dt
Definition: heat_method.cpp:26
MoFEM::SmartPetscObj< VecScatter >
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
TSPrePostProc::globRes
SmartPetscObj< Vec > globRes
Definition: free_surface.cpp:457
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEM::FieldBlas
Basic algebra on fields.
Definition: FieldBlas.hpp:21
MoFEM::getDMTsCtx
auto getDMTsCtx(DM dm)
Get TS context data structure used by DM.
Definition: DMMoFEM.hpp:1141
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
TSPrePostProc::fsRawPtr
Example * fsRawPtr
Definition: dynamic_first_order_con_law.cpp:398
TSPrePostProc::pcApply
static MoFEMErrorCode pcApply(PC pc, Vec pc_f, Vec pc_x)
Definition: free_surface.cpp:3098