v0.13.0
Public Types | Public Member Functions | Public Attributes | Private Attributes | Friends | List of all members
MoFEM::KspCtx Struct Reference

Interface for linear (KSP) solver. More...

#include <src/petsc/KspCtx.hpp>

Collaboration diagram for MoFEM::KspCtx:
[legend]

Public Types

typedef MoFEM::PairNameFEMethodPtr PairNameFEMethodPtr
 
typedef MoFEM::FEMethodsSequence FEMethodsSequence
 
typedef MoFEM::BasicMethodsSequence BasicMethodsSequence
 

Public Member Functions

 KspCtx (MoFEM::Interface &m_field, const std::string &_problem_name)
 
virtual ~KspCtx ()=default
 
FEMethodsSequenceget_loops_to_do_Mat ()
 
FEMethodsSequenceget_loops_to_do_Rhs ()
 
BasicMethodsSequenceget_preProcess_to_do_Rhs ()
 
BasicMethodsSequenceget_postProcess_to_do_Rhs ()
 
BasicMethodsSequenceget_preProcess_to_do_Mat ()
 
BasicMethodsSequenceget_postProcess_to_do_Mat ()
 
MoFEMErrorCode clearLoops ()
 Clear loops. More...
 

Public Attributes

MoFEM::InterfacemField
 
moab::Interface & moab
 
std::string problemName
 Problem name. More...
 
MoFEMTypes bH
 If set to MF_EXIST check if element exist. More...
 
FEMethodsSequence loops_to_do_Mat
 
FEMethodsSequence loops_to_do_Rhs
 
BasicMethodsSequence preProcess_Mat
 
BasicMethodsSequence postProcess_Mat
 
BasicMethodsSequence preProcess_Rhs
 
BasicMethodsSequence postProcess_Rhs
 

Private Attributes

PetscLogEvent MOFEM_EVENT_KspRhs
 
PetscLogEvent MOFEM_EVENT_KspMat
 
boost::movelib::unique_ptr< boolvecAssembleSwitch
 
boost::movelib::unique_ptr< boolmatAssembleSwitch
 

Friends

PetscErrorCode KspRhs (KSP ksp, Vec f, void *ctx)
 Run over elements in the lists. More...
 
PetscErrorCode KspMat (KSP ksp, Mat A, Mat B, void *ctx)
 Run over elenents in the list. More...
 

Detailed Description

Interface for linear (KSP) solver.

Definition at line 27 of file KspCtx.hpp.

Member Typedef Documentation

◆ BasicMethodsSequence

Definition at line 37 of file KspCtx.hpp.

◆ FEMethodsSequence

Definition at line 36 of file KspCtx.hpp.

◆ PairNameFEMethodPtr

Definition at line 35 of file KspCtx.hpp.

Constructor & Destructor Documentation

◆ KspCtx()

MoFEM::KspCtx::KspCtx ( MoFEM::Interface m_field,
const std::string &  _problem_name 
)

Definition at line 52 of file KspCtx.hpp.

53  : mField(m_field), moab(m_field.get_moab()), problemName(_problem_name),
54  bH(MF_EXIST) {
55  PetscLogEventRegister("LoopKSPRhs", 0, &MOFEM_EVENT_KspRhs);
56  PetscLogEventRegister("LoopKSPMat", 0, &MOFEM_EVENT_KspMat);
57  }
@ MF_EXIST
Definition: definitions.h:113
virtual moab::Interface & get_moab()=0
moab::Interface & moab
Definition: KspCtx.hpp:30
PetscLogEvent MOFEM_EVENT_KspRhs
Definition: KspCtx.hpp:113
PetscLogEvent MOFEM_EVENT_KspMat
Definition: KspCtx.hpp:114
MoFEM::Interface & mField
Definition: KspCtx.hpp:29
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition: KspCtx.hpp:33
std::string problemName
Problem name.
Definition: KspCtx.hpp:32

◆ ~KspCtx()

virtual MoFEM::KspCtx::~KspCtx ( )
virtualdefault

Member Function Documentation

◆ clearLoops()

MoFEMErrorCode MoFEM::KspCtx::clearLoops ( )

Clear loops.

Returns
MoFEMErrorCode

Definition at line 23 of file KspCtx.cpp.

23  {
25  loops_to_do_Mat.clear();
26  loops_to_do_Rhs.clear();
27  preProcess_Mat.clear();
28  postProcess_Mat.clear();
29  preProcess_Rhs.clear();
31 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
FEMethodsSequence loops_to_do_Rhs
Definition: KspCtx.hpp:41
BasicMethodsSequence preProcess_Rhs
Definition: KspCtx.hpp:47
BasicMethodsSequence postProcess_Mat
Definition: KspCtx.hpp:45
BasicMethodsSequence preProcess_Mat
Definition: KspCtx.hpp:43
FEMethodsSequence loops_to_do_Mat
Definition: KspCtx.hpp:39

◆ get_loops_to_do_Mat()

FEMethodsSequence& MoFEM::KspCtx::get_loops_to_do_Mat ( )
Returns
return reference to vector with FEMethod to calculate matrix

Definition at line 63 of file KspCtx.hpp.

63 { return loops_to_do_Mat; }

◆ get_loops_to_do_Rhs()

FEMethodsSequence& MoFEM::KspCtx::get_loops_to_do_Rhs ( )
Returns
return vector to vector with FEMethod to vector

Definition at line 68 of file KspCtx.hpp.

68 { return loops_to_do_Rhs; }

◆ get_postProcess_to_do_Mat()

BasicMethodsSequence& MoFEM::KspCtx::get_postProcess_to_do_Mat ( )

The sequence of BasicMethod is executed after tangent matrix is calculated. It can be used to setup data structures, e.g. aggregate data from processors or to apply essential boundary conditions.

Returns
reference to BasicMethod for postprocessing

Definition at line 100 of file KspCtx.hpp.

100 { return postProcess_Mat; }

◆ get_postProcess_to_do_Rhs()

BasicMethodsSequence& MoFEM::KspCtx::get_postProcess_to_do_Rhs ( )

The sequence of BasicMethod is executed after residual is calculated. It can be used to setup data structures, e.g. aggregate data from processors or to apply essential boundary conditions.

Returns
reference to BasicMethod for postprocessing

Definition at line 86 of file KspCtx.hpp.

86 { return postProcess_Rhs; }
BasicMethodsSequence postProcess_Rhs
Definition: KspCtx.hpp:49

◆ get_preProcess_to_do_Mat()

BasicMethodsSequence& MoFEM::KspCtx::get_preProcess_to_do_Mat ( )
Returns
reference to BasicMethod for preprocessing

Definition at line 91 of file KspCtx.hpp.

91 { return preProcess_Mat; }

◆ get_preProcess_to_do_Rhs()

BasicMethodsSequence& MoFEM::KspCtx::get_preProcess_to_do_Rhs ( )

The sequence of BasicMethod is executed before residual is calculated. It can be used to setup data structures, e.g. zero global variable which is integrated in domain, e.g. for calculation of strain energy.

Returns
reference to BasicMethod for preprocessing

Definition at line 77 of file KspCtx.hpp.

77 { return preProcess_Rhs; }

Friends And Related Function Documentation

◆ KspMat

PetscErrorCode KspMat ( KSP  ksp,
Mat  A,
Mat  B,
void *  ctx 
)
friend

Run over elenents in the list.

Parameters
kspKSP solver
Amatrix
BPreconditioner matrix
ctxdata context, i.e. KspCtx
Returns
error code

Definition at line 97 of file KspCtx.cpp.

97  {
98  // PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
100  KspCtx *ksp_ctx = static_cast<KspCtx *>(ctx);
101  PetscLogEventBegin(ksp_ctx->MOFEM_EVENT_KspMat, 0, 0, 0, 0);
102 
103  ksp_ctx->matAssembleSwitch = boost::movelib::make_unique<bool>(true);
104 
105  auto set = [&](auto &fe) {
106  fe.ksp = ksp;
107  fe.ksp_A = A;
108  fe.ksp_B = B;
109  fe.ksp_ctx = KspMethod::CTX_OPERATORS;
110  fe.data_ctx = PetscData::CtxSetA | PetscData::CtxSetB;
111  };
112 
113  auto unset = [&](auto &fe) {
114  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
115  fe.data_ctx = PetscData::CtxSetNone;
116  };
117 
118  auto ent_data_cache = boost::make_shared<std::vector<EntityCacheDofs>>();
119  auto ent_row_cache =
120  boost::make_shared<std::vector<EntityCacheNumeredDofs>>();
121  auto ent_col_cache =
122  boost::make_shared<std::vector<EntityCacheNumeredDofs>>();
123 
124  // pre-procsess
125  for (auto &bit : ksp_ctx->preProcess_Mat) {
126  bit->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
127  set(*bit);
128  CHKERR ksp_ctx->mField.problem_basic_method_preProcess(ksp_ctx->problemName,
129  *bit);
130  unset(*bit);
131  ksp_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
132  }
133 
134  auto cache_ptr = boost::make_shared<CacheTuple>();
135  CHKERR ksp_ctx->mField.cache_problem_entities(ksp_ctx->problemName,
136  cache_ptr);
137 
138  // operators
139  for (auto &lit : ksp_ctx->loops_to_do_Mat) {
140  lit.second->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
141  set(*lit.second);
142  CHKERR ksp_ctx->mField.loop_finite_elements(ksp_ctx->problemName, lit.first,
143  *(lit.second), nullptr,
144  ksp_ctx->bH, cache_ptr);
145  unset(*lit.second);
146  ksp_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
147  }
148 
149  // post-process
150  for (auto &bit : ksp_ctx->postProcess_Mat) {
151  bit->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
152  set(*bit);
153  CHKERR ksp_ctx->mField.problem_basic_method_postProcess(
154  ksp_ctx->problemName, *bit);
155  unset(*bit);
156  ksp_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
157  }
158 
159  if (ksp_ctx->matAssembleSwitch) {
160  CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
161  CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
162  }
163  PetscLogEventEnd(ksp_ctx->MOFEM_EVENT_KspMat, 0, 0, 0, 0);
165 }
#define CHKERR
Inline error check.
Definition: definitions.h:548
auto bit
set bit
double A
KspCtx(MoFEM::Interface &m_field, const std::string &_problem_name)
Definition: KspCtx.hpp:52
static constexpr Switches CtxSetA
Definition: LoopMethods.hpp:47
static constexpr Switches CtxSetNone
Definition: LoopMethods.hpp:45
static constexpr Switches CtxSetB
Definition: LoopMethods.hpp:48

◆ KspRhs

PetscErrorCode KspRhs ( KSP  ksp,
Vec  f,
void *  ctx 
)
friend

Run over elements in the lists.

Parameters
kspKSP solver
fthe right hand side vector
ctxdata context, i.e. KspCtx
Returns
error code

Definition at line 33 of file KspCtx.cpp.

33  {
34  // PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
36  KspCtx *ksp_ctx = static_cast<KspCtx *>(ctx);
37  PetscLogEventBegin(ksp_ctx->MOFEM_EVENT_KspRhs, 0, 0, 0, 0);
38 
39  ksp_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
40 
41  auto set = [&](auto &fe) {
42  fe.ksp = ksp;
43  fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
44  fe.data_ctx = PetscData::CtxSetF;
45  fe.ksp_f = f;
46  };
47 
48  auto unset = [&](auto &fe) {
49  fe.ksp_ctx = KspMethod::CTX_KSPNONE;
50  fe.data_ctx = PetscData::CtxSetNone;
51  };
52 
53  // pre-process
54  for (auto &bit : ksp_ctx->preProcess_Rhs) {
55  bit->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
56  set(*bit);
57  CHKERR ksp_ctx->mField.problem_basic_method_preProcess(ksp_ctx->problemName,
58  *bit);
59  unset(*bit);
60  ksp_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
61  }
62 
63  auto cache_ptr = boost::make_shared<CacheTuple>();
64  CHKERR ksp_ctx->mField.cache_problem_entities(ksp_ctx->problemName,
65  cache_ptr);
66 
67  // operators
68  for (auto &lit : ksp_ctx->loops_to_do_Rhs) {
69  lit.second->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
70  set(*lit.second);
71  CHKERR ksp_ctx->mField.loop_finite_elements(ksp_ctx->problemName, lit.first,
72  *(lit.second), nullptr,
73  ksp_ctx->bH, cache_ptr);
74  unset(*lit.second);
75  ksp_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
76  }
77 
78  // post-process
79  for (auto &bit : ksp_ctx->postProcess_Rhs) {
80  bit->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
81  set(*bit);
82  CHKERR ksp_ctx->mField.problem_basic_method_postProcess(
83  ksp_ctx->problemName, *bit);
84  unset(*bit);
85  ksp_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
86  }
87 
88  if (*ksp_ctx->vecAssembleSwitch) {
89  CHKERR VecGhostUpdateBegin(f, ADD_VALUES, SCATTER_REVERSE);
90  CHKERR VecGhostUpdateEnd(f, ADD_VALUES, SCATTER_REVERSE);
91  CHKERR VecAssemblyBegin(f);
92  CHKERR VecAssemblyEnd(f);
93  }
94  PetscLogEventEnd(ksp_ctx->MOFEM_EVENT_KspRhs, 0, 0, 0, 0);
96 }
static constexpr Switches CtxSetF
Definition: LoopMethods.hpp:46

Member Data Documentation

◆ bH

MoFEMTypes MoFEM::KspCtx::bH

If set to MF_EXIST check if element exist.

Definition at line 33 of file KspCtx.hpp.

◆ loops_to_do_Mat

FEMethodsSequence MoFEM::KspCtx::loops_to_do_Mat

Sequence of finite elements instances assembling tangent matrix

Definition at line 39 of file KspCtx.hpp.

◆ loops_to_do_Rhs

FEMethodsSequence MoFEM::KspCtx::loops_to_do_Rhs

Sequence of finite elements instances assembling residual vector

Definition at line 41 of file KspCtx.hpp.

◆ matAssembleSwitch

boost::movelib::unique_ptr<bool> MoFEM::KspCtx::matAssembleSwitch
private

Definition at line 117 of file KspCtx.hpp.

◆ mField

MoFEM::Interface& MoFEM::KspCtx::mField

Definition at line 29 of file KspCtx.hpp.

◆ moab

moab::Interface& MoFEM::KspCtx::moab

Definition at line 30 of file KspCtx.hpp.

◆ MOFEM_EVENT_KspMat

PetscLogEvent MoFEM::KspCtx::MOFEM_EVENT_KspMat
private

Definition at line 114 of file KspCtx.hpp.

◆ MOFEM_EVENT_KspRhs

PetscLogEvent MoFEM::KspCtx::MOFEM_EVENT_KspRhs
private

Definition at line 113 of file KspCtx.hpp.

◆ postProcess_Mat

BasicMethodsSequence MoFEM::KspCtx::postProcess_Mat

Sequence of methods run after tangent matrix is assembled

Definition at line 45 of file KspCtx.hpp.

◆ postProcess_Rhs

BasicMethodsSequence MoFEM::KspCtx::postProcess_Rhs

Sequence of methods run after residual is assembled

Definition at line 49 of file KspCtx.hpp.

◆ preProcess_Mat

BasicMethodsSequence MoFEM::KspCtx::preProcess_Mat

Sequence of methods run before tangent matrix is assembled

Definition at line 43 of file KspCtx.hpp.

◆ preProcess_Rhs

BasicMethodsSequence MoFEM::KspCtx::preProcess_Rhs

Sequence of methods run before residual is assembled

Definition at line 47 of file KspCtx.hpp.

◆ problemName

std::string MoFEM::KspCtx::problemName

Problem name.

Definition at line 32 of file KspCtx.hpp.

◆ vecAssembleSwitch

boost::movelib::unique_ptr<bool> MoFEM::KspCtx::vecAssembleSwitch
private

Definition at line 116 of file KspCtx.hpp.


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