v0.15.0
Loading...
Searching...
No Matches
MoFEM::KspCtx::KspCtx::KspCtxImpl Struct Reference
Collaboration diagram for MoFEM::KspCtx::KspCtx::KspCtxImpl:
[legend]

Public Member Functions

 KspCtxImpl (MoFEM::Interface &m_field, std::string problem_name)
 
virtual ~KspCtxImpl ()=default
 
FEMethodsSequencegetSetOperators ()
 
FEMethodsSequencegetComputeRhs ()
 
BasicMethodsSequencegetPreProcComputeRhs ()
 
BasicMethodsSequencegetPostProcComputeRhs ()
 
BasicMethodsSequencegetPreProcSetOperators ()
 
BasicMethodsSequencegetPostProcSetOperators ()
 
MoFEMErrorCode clearLoops ()
 Clear loops.
 

Protected Attributes

PetscLogEvent MOFEM_EVENT_KspRhs
 
PetscLogEvent MOFEM_EVENT_KspMat
 
boost::movelib::unique_ptr< boolvecAssembleSwitch
 
boost::movelib::unique_ptr< boolmatAssembleSwitch
 
MoFEM::InterfacemField
 
moab::Interface & moab
 
std::string problemName
 Problem name.
 
MoFEMTypes bH
 If set to MF_EXIST check if element exist.
 
FEMethodsSequence loopsOperator
 
FEMethodsSequence loopsRhs
 
BasicMethodsSequence preProcessOperator
 
BasicMethodsSequence postProcessOperator
 
BasicMethodsSequence preProcessRhs
 
BasicMethodsSequence postProcessRhs
 

Friends

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

Detailed Description

Definition at line 9 of file KspCtx.cpp.

Constructor & Destructor Documentation

◆ KspCtxImpl()

MoFEM::KspCtx::KspCtxImpl::KspCtxImpl ( MoFEM::Interface & m_field,
std::string problem_name )

Definition at line 128 of file KspCtx.cpp.

130 : mField(m_field), moab(m_field.get_moab()), problemName(problem_name),
131 bH(MF_EXIST) {
132 PetscLogEventRegister("LoopKSPRhs", 0, &MOFEM_EVENT_KspRhs);
133 PetscLogEventRegister("LoopKSPMat", 0, &MOFEM_EVENT_KspMat);
134}
@ MF_EXIST
virtual moab::Interface & get_moab()=0
MoFEM::Interface & mField
Definition KspCtx.cpp:76
std::string problemName
Problem name.
Definition KspCtx.cpp:79
MoFEMTypes bH
If set to MF_EXIST check if element exist.
Definition KspCtx.cpp:80

◆ ~KspCtxImpl()

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

Member Function Documentation

◆ clearLoops()

MoFEMErrorCode MoFEM::KspCtx::KspCtxImpl::clearLoops ( )

Clear loops.

Returns
MoFEMErrorCode

Definition at line 136 of file KspCtx.cpp.

136 {
138 loopsOperator.clear();
139 loopsRhs.clear();
140 preProcessOperator.clear();
141 postProcessOperator.clear();
142 preProcessRhs.clear();
144}
#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()
BasicMethodsSequence preProcessOperator
Definition KspCtx.cpp:86
BasicMethodsSequence postProcessOperator
Definition KspCtx.cpp:88
FEMethodsSequence loopsOperator
Definition KspCtx.cpp:82
BasicMethodsSequence preProcessRhs
Definition KspCtx.cpp:90

◆ getComputeRhs()

FEMethodsSequence & MoFEM::KspCtx::KspCtx::KspCtxImpl::getComputeRhs ( )
inline
Returns
return vector to vector with FEMethod to vector

Definition at line 23 of file KspCtx.cpp.

23{ return loopsRhs; }

◆ getPostProcComputeRhs()

BasicMethodsSequence & MoFEM::KspCtx::KspCtx::KspCtxImpl::getPostProcComputeRhs ( )
inline

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 41 of file KspCtx.cpp.

41{ return postProcessRhs; }
BasicMethodsSequence postProcessRhs
Definition KspCtx.cpp:92

◆ getPostProcSetOperators()

BasicMethodsSequence & MoFEM::KspCtx::KspCtx::KspCtxImpl::getPostProcSetOperators ( )
inline

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 55 of file KspCtx.cpp.

55 {
57 }

◆ getPreProcComputeRhs()

BasicMethodsSequence & MoFEM::KspCtx::KspCtx::KspCtxImpl::getPreProcComputeRhs ( )
inline

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 32 of file KspCtx.cpp.

32{ return preProcessRhs; }

◆ getPreProcSetOperators()

BasicMethodsSequence & MoFEM::KspCtx::KspCtx::KspCtxImpl::getPreProcSetOperators ( )
inline
Returns
reference to BasicMethod for preprocessing

Definition at line 46 of file KspCtx.cpp.

46{ return preProcessOperator; }

◆ getSetOperators()

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

Definition at line 18 of file KspCtx.cpp.

18{ return loopsOperator; }

Friends And Related Symbol Documentation

◆ KspMat

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

Run over elements in the list.

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

Definition at line 210 of file KspCtx.cpp.

210 {
211 // PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
213 KspCtx::KspCtxImpl *ksp_ctx = static_cast<KspCtx *>(ctx)->kspCtxImpl.get();
214 PetscLogEventBegin(ksp_ctx->MOFEM_EVENT_KspMat, 0, 0, 0, 0);
215
216 ksp_ctx->matAssembleSwitch = boost::movelib::make_unique<bool>(true);
217 auto cache_ptr = boost::make_shared<CacheTuple>();
218 CHKERR ksp_ctx->mField.cache_problem_entities(ksp_ctx->problemName,
219 cache_ptr);
220
221 auto set = [&](auto &fe) {
222 fe.ksp = ksp;
223 fe.ksp_A = A;
224 fe.ksp_B = B;
225 fe.ksp_ctx = KspMethod::CTX_OPERATORS;
226 fe.data_ctx = PetscData::CtxSetA | PetscData::CtxSetB;
227 fe.cacheWeakPtr = cache_ptr;
228 };
229
230 auto unset = [&](auto &fe) {
231 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
232 fe.data_ctx = PetscData::CtxSetNone;
233 };
234
235 auto ent_data_cache = boost::make_shared<std::vector<EntityCacheDofs>>();
236 auto ent_row_cache =
237 boost::make_shared<std::vector<EntityCacheNumeredDofs>>();
238 auto ent_col_cache =
239 boost::make_shared<std::vector<EntityCacheNumeredDofs>>();
240
241 // pre-procsess
242 for (auto &bit : ksp_ctx->preProcessOperator) {
243 bit->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
244 set(*bit);
245 CHKERR ksp_ctx->mField.problem_basic_method_preProcess(ksp_ctx->problemName,
246 *bit);
247 unset(*bit);
248 ksp_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
249 }
250
251 // operators
252 for (auto &lit : ksp_ctx->loopsOperator) {
253 lit.second->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
254 set(*lit.second);
255 CHKERR ksp_ctx->mField.loop_finite_elements(ksp_ctx->problemName, lit.first,
256 *(lit.second), nullptr,
257 ksp_ctx->bH, cache_ptr);
258 unset(*lit.second);
259 ksp_ctx->matAssembleSwitch = boost::move(lit.second->matAssembleSwitch);
260 }
261
262 // post-process
263 for (auto &bit : ksp_ctx->postProcessOperator) {
264 bit->matAssembleSwitch = boost::move(ksp_ctx->matAssembleSwitch);
265 set(*bit);
266 CHKERR ksp_ctx->mField.problem_basic_method_postProcess(
267 ksp_ctx->problemName, *bit);
268 unset(*bit);
269 ksp_ctx->matAssembleSwitch = boost::move(bit->matAssembleSwitch);
270 }
271
272 if (ksp_ctx->matAssembleSwitch) {
273 CHKERR MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
274 CHKERR MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
275 }
276 PetscLogEventEnd(ksp_ctx->MOFEM_EVENT_KspMat, 0, 0, 0, 0);
278}
#define CHKERR
Inline error check.
auto bit
set bit
constexpr AssemblyType A

◆ 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 146 of file KspCtx.cpp.

146 {
148 KspCtx::KspCtxImpl *ksp_ctx = static_cast<KspCtx *>(ctx)->kspCtxImpl.get();
149 PetscLogEventBegin(ksp_ctx->MOFEM_EVENT_KspRhs, 0, 0, 0, 0);
150
151 ksp_ctx->vecAssembleSwitch = boost::movelib::make_unique<bool>(true);
152
153 auto cache_ptr = boost::make_shared<CacheTuple>();
154 CHKERR ksp_ctx->mField.cache_problem_entities(ksp_ctx->problemName,
155 cache_ptr);
156
157 auto set = [&](auto &fe) {
158 fe.ksp = ksp;
159 fe.ksp_ctx = KspMethod::CTX_SETFUNCTION;
160 fe.data_ctx = PetscData::CtxSetF;
161 fe.ksp_f = f;
162 fe.cacheWeakPtr = cache_ptr;
163 };
164
165 auto unset = [&](auto &fe) {
166 fe.ksp_ctx = KspMethod::CTX_KSPNONE;
167 fe.data_ctx = PetscData::CtxSetNone;
168 };
169
170 // pre-process
171 for (auto &bit : ksp_ctx->preProcessRhs) {
172 bit->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
173 set(*bit);
174 CHKERR ksp_ctx->mField.problem_basic_method_preProcess(ksp_ctx->problemName,
175 *bit);
176 unset(*bit);
177 ksp_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
178 }
179
180 // operators
181 for (auto &lit : ksp_ctx->loopsRhs) {
182 lit.second->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
183 set(*lit.second);
184 CHKERR ksp_ctx->mField.loop_finite_elements(ksp_ctx->problemName, lit.first,
185 *(lit.second), nullptr,
186 ksp_ctx->bH, cache_ptr);
187 unset(*lit.second);
188 ksp_ctx->vecAssembleSwitch = boost::move(lit.second->vecAssembleSwitch);
189 }
190
191 // post-process
192 for (auto &bit : ksp_ctx->postProcessRhs) {
193 bit->vecAssembleSwitch = boost::move(ksp_ctx->vecAssembleSwitch);
194 set(*bit);
195 CHKERR ksp_ctx->mField.problem_basic_method_postProcess(
196 ksp_ctx->problemName, *bit);
197 unset(*bit);
198 ksp_ctx->vecAssembleSwitch = boost::move(bit->vecAssembleSwitch);
199 }
200
201 if (*ksp_ctx->vecAssembleSwitch) {
202 CHKERR VecGhostUpdateBegin(f, ADD_VALUES, SCATTER_REVERSE);
203 CHKERR VecGhostUpdateEnd(f, ADD_VALUES, SCATTER_REVERSE);
204 CHKERR VecAssemblyBegin(f);
205 CHKERR VecAssemblyEnd(f);
206 }
207 PetscLogEventEnd(ksp_ctx->MOFEM_EVENT_KspRhs, 0, 0, 0, 0);
209}

Member Data Documentation

◆ bH

MoFEMTypes MoFEM::KspCtx::KspCtx::KspCtxImpl::bH
protected

If set to MF_EXIST check if element exist.

Definition at line 80 of file KspCtx.cpp.

◆ loopsOperator

FEMethodsSequence MoFEM::KspCtx::KspCtx::KspCtxImpl::loopsOperator
protected

Sequence of finite elements instances assembling tangent matrix

Definition at line 82 of file KspCtx.cpp.

◆ loopsRhs

FEMethodsSequence MoFEM::KspCtx::KspCtx::KspCtxImpl::loopsRhs
protected

Sequence of finite elements instances assembling residual vector

Definition at line 84 of file KspCtx.cpp.

◆ matAssembleSwitch

boost::movelib::unique_ptr<bool> MoFEM::KspCtx::KspCtx::KspCtxImpl::matAssembleSwitch
protected

Definition at line 74 of file KspCtx.cpp.

◆ mField

MoFEM::Interface& MoFEM::KspCtx::KspCtx::KspCtxImpl::mField
protected

Definition at line 76 of file KspCtx.cpp.

◆ moab

moab::Interface& MoFEM::KspCtx::KspCtx::KspCtxImpl::moab
protected

Definition at line 77 of file KspCtx.cpp.

◆ MOFEM_EVENT_KspMat

PetscLogEvent MoFEM::KspCtx::KspCtx::KspCtxImpl::MOFEM_EVENT_KspMat
protected

Definition at line 71 of file KspCtx.cpp.

◆ MOFEM_EVENT_KspRhs

PetscLogEvent MoFEM::KspCtx::KspCtx::KspCtxImpl::MOFEM_EVENT_KspRhs
protected

Definition at line 70 of file KspCtx.cpp.

◆ postProcessOperator

BasicMethodsSequence MoFEM::KspCtx::KspCtx::KspCtxImpl::postProcessOperator
protected

Sequence of methods run after tangent matrix is assembled

Definition at line 88 of file KspCtx.cpp.

◆ postProcessRhs

BasicMethodsSequence MoFEM::KspCtx::KspCtx::KspCtxImpl::postProcessRhs
protected

Sequence of methods run after residual is assembled

Definition at line 92 of file KspCtx.cpp.

◆ preProcessOperator

BasicMethodsSequence MoFEM::KspCtx::KspCtx::KspCtxImpl::preProcessOperator
protected

Sequence of methods run before tangent matrix is assembled

Definition at line 86 of file KspCtx.cpp.

◆ preProcessRhs

BasicMethodsSequence MoFEM::KspCtx::KspCtx::KspCtxImpl::preProcessRhs
protected

Sequence of methods run before residual is assembled

Definition at line 90 of file KspCtx.cpp.

◆ problemName

std::string MoFEM::KspCtx::KspCtx::KspCtxImpl::problemName
protected

Problem name.

Definition at line 79 of file KspCtx.cpp.

◆ vecAssembleSwitch

boost::movelib::unique_ptr<bool> MoFEM::KspCtx::KspCtx::KspCtxImpl::vecAssembleSwitch
protected

Definition at line 73 of file KspCtx.cpp.


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