v0.15.0
Loading...
Searching...
No Matches
SetUpSchurImpl.cpp
Go to the documentation of this file.
1/** @file SetUpSchurImpl
2 * @brief
3 * @date 2023-05-13
4 *
5 * @license{This project is released under the MIT License.}
6 *
7 */
8
10
12 : SetUpSchur(), mField(m_field), epCorePtr(ep_core_ptr) {}
13 virtual ~SetUpSchurImpl() {}
14
15 MoFEMErrorCode setUp(TS ts);
16 MoFEMErrorCode preProc();
17 MoFEMErrorCode postProc();
18
19private:
22
23 SmartPetscObj<Mat> S;
24 SmartPetscObj<AO> aoS;
25
26 boost::shared_ptr<std::vector<boost::weak_ptr<NumeredDofEntity>>>
27 piolaZeroDofsVec; //< Dofs on crack surface
28 boost::shared_ptr<std::vector<unsigned char>>
29 piolaZeroDofsMarker; //< marker for crack dofs on surface
30
32
33 P_MultiGridData(SmartPetscObj<DM> dm, SmartPetscObj<PC> pc,
34 SmartPetscObj<Mat> S)
35 : dm(dm), pc(pc), S(S) {}
36
37 MoFEMErrorCode setUP() {
39 CHKERR PCMGSetUpViaApproxOrders(
40 pc, createPCMGSetUpViaApproxOrdersCtx(dm, S, true));
41 CHKERR PCSetFromOptions(pc);
43 };
44
45 private:
46 SmartPetscObj<DM> dm;
47 SmartPetscObj<PC> pc;
48 SmartPetscObj<Mat> S;
49 };
50
51 boost::shared_ptr<P_MultiGridData> pMGPtr;
52
54 std::vector<std::string> schur_field_list{epCorePtr->hybridSpatialDisp,
56 std::vector<boost::shared_ptr<Range>> dm_range_list{nullptr, nullptr};
57 return std::make_pair(schur_field_list, dm_range_list);
58 };
59
60 auto getA00Fields() {
61 std::vector<std::string> a00_field_list{
62
64
66
68
70
72
73 };
74 std::vector<boost::shared_ptr<Range>> range_list_ptr(a00_field_list.size(),
75 nullptr);
76 return std::make_pair(a00_field_list, range_list_ptr);
77 }
78};
79
82
83 auto create_schur_dm = [&](SmartPetscObj<DM> &dm_sub) {
85
86 dm_sub = createDM(mField.get_comm(), "DMMOFEM_MG");
87 CHKERR DMMoFEMCreateSubDM(dm_sub, epCorePtr->dmElastic, "SUB_SCHUR");
88 CHKERR DMMoFEMSetSquareProblem(dm_sub, PETSC_TRUE);
89 CHKERR DMMoFEMSetIsPartitioned(dm_sub, PETSC_TRUE);
92
93 int r_idx = 0;
94 auto [schur_field_list, schur_range_list] = getSchurFields();
95 for (auto f : schur_field_list) {
96 MOFEM_LOG("EP", Sev::inform) << "Add schur field: " << f;
97 CHKERR DMMoFEMAddSubFieldRow(dm_sub, f, schur_range_list[r_idx]);
98 CHKERR DMMoFEMAddSubFieldCol(dm_sub, f, schur_range_list[r_idx]);
99 ++r_idx;
100 }
101 CHKERR DMSetUp(dm_sub);
103 };
104
105 auto create_a00_dm = [&](SmartPetscObj<DM> &dm_sub) {
107 dm_sub = createDM(mField.get_comm(), "DMMOFEM");
108 CHKERR DMMoFEMCreateSubDM(dm_sub, epCorePtr->dmElastic, "SUB_A00");
109 CHKERR DMMoFEMSetSquareProblem(dm_sub, PETSC_TRUE);
110 CHKERR DMMoFEMSetIsPartitioned(dm_sub, PETSC_TRUE);
114
115 int r_idx = 0;
116 auto [a00_field_list, a00_range_list] = getA00Fields();
117 for (auto f : a00_field_list) {
118 MOFEM_LOG("EP", Sev::inform) << "Add a00 field: " << f;
119 CHKERR DMMoFEMAddSubFieldRow(dm_sub, f, a00_range_list[r_idx]);
120 CHKERR DMMoFEMAddSubFieldCol(dm_sub, f, a00_range_list[r_idx]);
121 ++r_idx;
122 }
123 CHKERR DMSetUp(dm_sub);
125 };
126
127 auto get_snes = [&](TS ts) {
128 SNES snes;
129 CHKERR TSGetSNES(ts, &snes);
130 return snes;
131 };
132
133 auto get_ksp = [&](SNES snes) {
134 KSP ksp;
135 CHKERR SNESGetKSP(snes, &ksp);
136 CHKERR KSPSetFromOptions(ksp);
137 return ksp;
138 };
139
140 auto get_pc = [&](KSP ksp) {
141 PC pc;
142 CHKERR KSPGetPC(ksp, &pc);
143 return pc;
144 };
145
146 auto ksp = get_ksp(get_snes(ts));
147 auto pc = get_pc(ksp);
148
149 PetscBool is_pcfs = PETSC_FALSE;
150 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
151 if (is_pcfs) {
152
153 MOFEM_LOG("EP", Sev::inform) << "SetUpSchurImpl::setUp: PCFIELDSPLIT";
154
155 SmartPetscObj<DM> schur_dm, a00_dm;
156 CHKERR create_schur_dm(schur_dm);
157 CHKERR create_a00_dm(a00_dm);
158
159 auto dm_elastic = epCorePtr->dmElastic;
160 auto vol_elem_name = epCorePtr->elementVolumeName;
161 auto skel_elem_name = epCorePtr->skeletonElement;
162 auto contact_elem_name = epCorePtr->contactElement;
163 auto natural_bc_element_name = epCorePtr->naturalBcElement;
164
165 std::vector<std::pair<std::string, std::string>> mat_block_list = {
166
173
181
182 };
183
184 if (epCorePtr->noStretch) {
185 mat_block_list.push_back(
187 mat_block_list.push_back(
189 mat_block_list.push_back(
191 }
192
193 // if (epCorePtr->gradApproximator > EshelbianPlasticity::SMALL_ROT ||
194 // epCorePtr->rotSelector > EshelbianPlasticity::SMALL_ROT) {
195 mat_block_list.push_back({epCorePtr->rotAxis, epCorePtr->rotAxis});
196 mat_block_list.push_back({epCorePtr->stretchTensor, epCorePtr->rotAxis});
197 mat_block_list.push_back({epCorePtr->rotAxis, epCorePtr->stretchTensor});
198 // }
199
200 auto get_nested_mat_data = [&](auto schur_dm, auto block_dm) {
201 auto block_mat_data = createBlockMatStructure(
202 dm_elastic,
203
204 {
205
206 {vol_elem_name, mat_block_list},
207
208 {skel_elem_name,
209
210 {
211
214
215 }},
216
217 {contact_elem_name,
218
219 {
220
224
225 }},
226
227 {natural_bc_element_name,
228
229 {
230
234
235 }}
236
237 }
238
239 );
240
241 auto [a00_field_list, a00_range_list] = getA00Fields();
242
244
245 {schur_dm, a00_dm}, block_mat_data,
246
247 a00_field_list,
248
249 a00_range_list,
250
251 false
252
253 );
254 };
255
256 auto nested_mat_data = get_nested_mat_data(schur_dm, a00_dm);
258 CHKERR DMSetMatType(epCorePtr->dmElastic, MATSHELL);
259
262
263 if (std::abs(epCorePtr->alphaRho) >
264 std::numeric_limits<double>::epsilon()) {
265 auto swap_assemble = [](TS ts, PetscReal t, Vec u, Vec u_t, Vec utt,
266 PetscReal a, PetscReal aa, Mat A, Mat B,
267 void *ctx) {
268 return TsSetI2Jacobian(ts, t, u, u_t, utt, a, aa, B, A, ctx);
269 };
270 auto ts_ctx_ptr = getDMTsCtx(epCorePtr->dmElastic);
271 CHKERR TSSetI2Jacobian(ts, m, p, swap_assemble, ts_ctx_ptr.get());
272 } else {
273 auto swap_assemble = [](TS ts, PetscReal t, Vec u, Vec u_t, PetscReal a,
274 Mat A, Mat B, void *ctx) {
275 return TsSetIJacobian(ts, t, u, u_t, a, B, A, ctx);
276 };
277 auto ts_ctx_ptr = getDMTsCtx(epCorePtr->dmElastic);
278 CHKERR TSSetIJacobian(ts, m, p, swap_assemble, ts_ctx_ptr.get());
279 }
280 CHKERR KSPSetOperators(ksp, m, p);
281
282 auto set_assembly = [&]() {
284
285 aoS = getDMSubData(schur_dm)->getSmartRowMap();
286 S = createDMHybridisedL2Matrix(schur_dm);
287 CHKERR MatSetBlockSize(S, SPACE_DIM *
289 epCorePtr->S = S;
290 epCorePtr->aoS = aoS;
291
292 auto set_assemble = [&]() {
294 auto schur_asmb_pre_proc_lhs = boost::make_shared<FEMethod>();
295 auto schur_asmb_pre_proc_rhs = boost::make_shared<FEMethod>();
296
297 schur_asmb_pre_proc_lhs->preProcessHook = [this]() {
299 CHKERR MatZeroEntries(S);
301 };
302
303 schur_asmb_pre_proc_rhs->preProcessHook = [this,
304 schur_asmb_pre_proc_rhs]() {
306 auto prb_ptr = schur_asmb_pre_proc_rhs->problemPtr;
307 auto dofs_prb = prb_ptr->getNumeredRowDofsPtr();
308
309 auto crack_faces = epCorePtr->crackFaces;
310 piolaZeroDofsVec->clear();
311 CHKERR mField.getInterface<ProblemsManager>()
312 ->getSideDofsOnBrokenSpaceEntities(
313 *piolaZeroDofsVec, prb_ptr->getName(), ROW,
314 epCorePtr->piolaStress, *crack_faces, SPACE_DIM, 0,
315 SPACE_DIM);
316
317 piolaZeroDofsMarker->resize(dofs_prb->size(), 0);
318 piolaZeroDofsMarker->clear();
319 for (auto &dof : *piolaZeroDofsVec) {
320 if (auto dof_ptr = dof.lock()) {
321 auto idx = dof_ptr->getPetscLocalDofIdx();
322 (*piolaZeroDofsMarker)[idx] = 1;
323 }
324 }
325
326 constexpr bool hard_coded_set_bc_debug = false;
327 if constexpr (hard_coded_set_bc_debug) {
328
329 auto problem_name = schur_asmb_pre_proc_rhs->problemPtr->getName();
330 auto crack_faces = epCorePtr->crackFaces;
331
332 SmartPetscObj<IS> crack_hybrid_is;
333 CHKERR epCorePtr->mField.getInterface<ISManager>()
334 ->isCreateProblemFieldAndRankLocal(
335 problem_name, ROW, epCorePtr->hybridSpatialDisp, 0,
336 SPACE_DIM, crack_hybrid_is, &*crack_faces);
337
338 SmartPetscObj<IS> crack_piola_is;
339 CHKERR epCorePtr->mField.getInterface<ISManager>()
340 ->isCreateProblemBrokenFieldAndRankLocal(*piolaZeroDofsVec,
341 crack_piola_is);
342
343 const double *a_x;
344 CHKERR VecGetArrayRead(schur_asmb_pre_proc_rhs->x, &a_x);
345 auto zero_by_is = [&](auto is) {
347 const PetscInt *is_array;
348 PetscInt is_size;
349 CHKERR ISGetLocalSize(is, &is_size);
350 CHKERR ISGetIndices(is, &is_array);
351 for (int i = 0; i != is_size; ++i) {
352 // FIXME: That is irregular, need investigate if needed
353 const_cast<double *>(a_x)[is_array[i]] = 0;
354 }
355 CHKERR ISRestoreIndices(is, &is_array);
357 };
358
359 CHKERR zero_by_is(crack_hybrid_is);
360 CHKERR zero_by_is(crack_piola_is);
361
362 CHKERR VecRestoreArrayRead(schur_asmb_pre_proc_rhs->x, &a_x);
363
364 CHKERR epCorePtr->mField.getInterface<VecManager>()
365 ->setLocalGhostVector(problem_name, COL,
366 schur_asmb_pre_proc_rhs->x, INSERT_VALUES,
367 SCATTER_REVERSE);
368 }
369
371 };
372
373 auto schur_asmb_post_proc_lhs = boost::make_shared<FEMethod>();
374 auto schur_asmb_post_proc_rhs = boost::make_shared<FEMethod>();
375
376 schur_asmb_post_proc_rhs->postProcessHook =
377 [this, schur_asmb_post_proc_rhs]() {
379
380 CHKERR VecGhostUpdateBegin(schur_asmb_post_proc_rhs->f,
381 ADD_VALUES, SCATTER_REVERSE);
382 CHKERR VecGhostUpdateEnd(schur_asmb_post_proc_rhs->f, ADD_VALUES,
383 SCATTER_REVERSE);
384 CHKERR VecAssemblyBegin(schur_asmb_post_proc_rhs->f);
385 CHKERR VecAssemblyEnd(schur_asmb_post_proc_rhs->f);
386 *(schur_asmb_post_proc_rhs->vecAssembleSwitch) = false;
387
388 {
389
390 auto problem_name =
391 schur_asmb_post_proc_rhs->problemPtr->getName();
392
393 auto crack_faces = epCorePtr->crackFaces;
394
395 SmartPetscObj<IS> crack_hybrid_is;
396 CHKERR epCorePtr->mField.getInterface<ISManager>()
397 ->isCreateProblemFieldAndRankLocal(
398 problem_name, ROW, epCorePtr->hybridSpatialDisp, 0,
399 SPACE_DIM, crack_hybrid_is, &*crack_faces);
400
401 SmartPetscObj<IS> crack_piola_is;
402 CHKERR epCorePtr->mField.getInterface<ISManager>()
403 ->isCreateProblemBrokenFieldAndRankLocal(*piolaZeroDofsVec,
404 crack_piola_is);
405
406 double *a_f;
407 CHKERR VecGetArray(schur_asmb_post_proc_rhs->f, &a_f);
408 const double *a_x;
409 CHKERR VecGetArrayRead(schur_asmb_post_proc_rhs->x, &a_x);
410 auto zero_by_is = [&](auto is) {
412 const PetscInt *is_array;
413 PetscInt is_size;
414 CHKERR ISGetLocalSize(is, &is_size);
415 CHKERR ISGetIndices(is, &is_array);
416 for (int i = 0; i != is_size; ++i) {
417 a_f[is_array[i]] = -a_x[is_array[i]];
418 }
419 CHKERR ISRestoreIndices(is, &is_array);
421 };
422
423 CHKERR zero_by_is(crack_hybrid_is);
424 CHKERR zero_by_is(crack_piola_is);
425
426 CHKERR VecRestoreArray(schur_asmb_post_proc_rhs->f, &a_f);
427 CHKERR VecRestoreArrayRead(schur_asmb_post_proc_rhs->x, &a_x);
428 }
429
431 };
432
433 schur_asmb_post_proc_lhs->postProcessHook =
434 [this, schur_asmb_post_proc_lhs]() {
436
437 if (pMGPtr) {
438 CHKERR pMGPtr->setUP();
439 pMGPtr.reset();
440 }
441
442 auto crack_faces = epCorePtr->crackFaces;
443
444 // Assemble matrix
445 CHKERR MatAssemblyBegin(schur_asmb_post_proc_lhs->B,
446 MAT_FINAL_ASSEMBLY);
447 CHKERR MatAssemblyEnd(schur_asmb_post_proc_lhs->B,
448 MAT_FINAL_ASSEMBLY);
449 *(schur_asmb_post_proc_lhs->matAssembleSwitch) = false;
450 {
451 SmartPetscObj<IS> crack_hybrid_is;
452 CHKERR epCorePtr->mField.getInterface<ISManager>()
453 ->isCreateProblemFieldAndRank(
454 "ELASTIC_PROBLEM", ROW, epCorePtr->hybridSpatialDisp, 0,
455 SPACE_DIM, crack_hybrid_is, &*crack_faces);
456 CHKERR MatZeroRowsColumnsIS(schur_asmb_post_proc_lhs->B,
457 crack_hybrid_is, 1, PETSC_NULLPTR,
458 PETSC_NULLPTR);
459 }
460 {
461 SmartPetscObj<IS> crack_piola_is;
462 CHKERR epCorePtr->mField.getInterface<ISManager>()
463 ->isCreateProblemBrokenFieldAndRank(*piolaZeroDofsVec,
464 crack_piola_is);
465 CHKERR MatZeroRowsColumnsIS(schur_asmb_post_proc_lhs->B,
466 crack_piola_is, 1, PETSC_NULLPTR,
467 PETSC_NULLPTR);
468 }
469
470 auto [a00_field_list, a00_range_list] = getA00Fields();
472 schur_asmb_post_proc_lhs->B, S,
473 a00_field_list, a00_range_list, aoS);
474 epCorePtr->a00FieldList = a00_field_list;
475 epCorePtr->a00RangeList = a00_range_list;
476
477 // Apply essential constrains to Schur complement
478 CHKERR MatAssemblyBegin(S, MAT_FINAL_ASSEMBLY);
479 CHKERR MatAssemblyEnd(S, MAT_FINAL_ASSEMBLY);
480
481 SmartPetscObj<IS> crack_hybrid_is;
482 CHKERR epCorePtr->mField.getInterface<ISManager>()
483 ->isCreateProblemFieldAndRank(
484 "SUB_SCHUR", ROW, epCorePtr->hybridSpatialDisp, 0,
485 SPACE_DIM, crack_hybrid_is, &*crack_faces);
486 epCorePtr->crackHybridIs = crack_hybrid_is;
487 CHKERR MatZeroRowsColumnsIS(S, crack_hybrid_is, 1, PETSC_NULLPTR,
488 PETSC_NULLPTR);
489
491 };
492
493 auto ts_ctx_ptr = getDMTsCtx(epCorePtr->dmElastic);
494 ts_ctx_ptr->getPreProcessIFunction().push_front(
495 schur_asmb_pre_proc_rhs);
496 ts_ctx_ptr->getPostProcessIFunction().push_back(
497 schur_asmb_post_proc_rhs);
498 ts_ctx_ptr->getPreProcessIJacobian().push_front(
499 schur_asmb_pre_proc_lhs);
500 ts_ctx_ptr->getPostProcessIJacobian().push_back(
501 schur_asmb_post_proc_lhs);
503 };
504
506 boost::make_shared<std::vector<boost::weak_ptr<NumeredDofEntity>>>();
507 piolaZeroDofsMarker = boost::make_shared<std::vector<unsigned char>>();
508 CHKERR set_assemble();
509
511 };
512
513 auto set_pc = [&]() {
515 auto a00_is = getDMSubData(a00_dm)->getSmartRowIs();
516 auto schur_is = getDMSubData(schur_dm)->getSmartRowIs();
517 CHKERR PCFieldSplitSetIS(pc, NULL, a00_is);
518 CHKERR PCFieldSplitSetIS(pc, NULL, schur_is);
519 CHKERR PCFieldSplitSetSchurPre(pc, PC_FIELDSPLIT_SCHUR_PRE_USER, S);
521 };
522
523 auto set_diagonal_pc = [&]() {
525 KSP *subksp;
526 CHKERR PCFieldSplitSchurGetSubKSP(pc, PETSC_NULLPTR, &subksp);
527 auto get_pc = [](auto ksp) {
528 PC pc_raw;
529 CHKERR KSPGetPC(ksp, &pc_raw);
530 return SmartPetscObj<PC>(pc_raw, true); // bump reference
531 };
532 CHKERR setSchurA00MatSolvePC(get_pc(subksp[0]));
533
534 auto set_pc_p_mg = [&](auto dm, auto pc, auto S) {
536 CHKERR PCSetDM(pc, dm);
537 PetscBool same = PETSC_FALSE;
538 PetscObjectTypeCompare((PetscObject)pc, PCMG, &same);
539 if (same) {
540 auto smart_pc = SmartPetscObj<PC>(pc, true);
541 pMGPtr = boost::make_shared<P_MultiGridData>(dm, smart_pc, S);
542 }
543 PetscObjectTypeCompare((PetscObject)pc, PCKSP, &same);
544 if (same) {
545 MOFEM_LOG("EP", Sev::inform)
546 << "SetUpSchurImpl::setUp: fieldsplit 1 PCKSP";
547 CHKERR PCSetFromOptions(pc);
548 KSP ksp;
549 CHKERR PCKSPGetKSP(pc, &ksp);
550 CHKERR KSPSetFromOptions(ksp);
551 PC ksp_pc;
552 CHKERR KSPGetPC(ksp, &ksp_pc);
553 CHKERR PCSetFromOptions(ksp_pc);
554 PetscObjectTypeCompare((PetscObject)ksp_pc, PCMG, &same);
555 if (same) {
556 auto smart_pc = SmartPetscObj<PC>(ksp_pc, true);
557 pMGPtr = boost::make_shared<P_MultiGridData>(dm, smart_pc, S);
558 }
559 }
561 };
562
563 CHKERR set_pc_p_mg(schur_dm, get_pc(subksp[1]), S);
564
565 CHKERR PetscFree(subksp);
567 };
568
569 CHKERR set_assembly();
570 CHKERR set_pc();
571 CHKERR TSSetUp(ts);
572 CHKERR KSPSetUp(ksp);
573 CHKERR set_diagonal_pc();
574
575 } else {
576 MOFEM_LOG("EP", Sev::inform) << "SetUpSchurImpl::setUp: PCLU or other";
577
578 epCorePtr->elasticFeLhs->getOpPtrVector().push_front(
580 epCorePtr->elasticFeLhs->getOpPtrVector().push_back(
582 epCorePtr->elasticBcLhs->getOpPtrVector().push_front(
584 epCorePtr->elasticBcLhs->getOpPtrVector().push_back(
586 }
587
589}
590
591boost::shared_ptr<EshelbianCore::SetUpSchur>
593 EshelbianCore *ep_core_ptr) {
594 return boost::shared_ptr<SetUpSchur>(
595 new SetUpSchurImpl(m_field, ep_core_ptr));
596}
constexpr double a
constexpr int SPACE_DIM
@ COL
@ ROW
#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.
PetscErrorCode DMMoFEMSetIsPartitioned(DM dm, PetscBool is_partitioned)
Definition DMMoFEM.cpp:1113
PetscErrorCode DMMoFEMCreateSubDM(DM subdm, DM dm, const char problem_name[])
Must be called by user to set Sub DM MoFEM data structures.
Definition DMMoFEM.cpp:215
PetscErrorCode DMMoFEMAddElement(DM dm, std::string fe_name)
add element to dm
Definition DMMoFEM.cpp:488
PetscErrorCode DMMoFEMSetSquareProblem(DM dm, PetscBool square_problem)
set squared problem
Definition DMMoFEM.cpp:450
PetscErrorCode DMMoFEMAddSubFieldRow(DM dm, const char field_name[])
Definition DMMoFEM.cpp:238
PetscErrorCode DMMoFEMAddSubFieldCol(DM dm, const char field_name[])
Definition DMMoFEM.cpp:280
#define MOFEM_LOG(channel, severity)
Log.
#define NBFACETRI_L2(P)
Number of base functions on triangle for L2 space.
FTensor::Index< 'i', SPACE_DIM > i
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
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
auto getDMTsCtx(DM dm)
Get TS context data structure used by DM.
Definition DMMoFEM.hpp:1144
OpSchurAssembleBase * createOpSchurAssembleEnd(std::vector< std::string > fields_name, std::vector< boost::shared_ptr< Range > > field_ents, SmartPetscObj< AO > ao, SmartPetscObj< Mat > schur, bool sym_schur, bool symm_op)
Construct a new Op Schur Assemble End object.
Definition Schur.cpp:2585
MoFEMErrorCode setSchurA00MatSolvePC(SmartPetscObj< PC > pc)
Set PC for A00 block.
Definition Schur.cpp:2627
auto createDMHybridisedL2Matrix(DM dm)
Get smart hybridised L2 matrix from DM.
Definition DMMoFEM.hpp:1072
auto getDMSubData(DM dm)
Get sub problem data structure.
Definition DMMoFEM.hpp:1160
PetscErrorCode TsSetI2Jacobian(TS ts, PetscReal t, Vec u, Vec u_t, Vec u_tt, PetscReal a, PetscReal aa, Mat A, Mat B, void *ctx)
Calculation Jacobian for second order PDE in time.
Definition TsCtx.cpp:519
boost::shared_ptr< BlockStructure > createBlockMatStructure(DM dm, SchurFEOpsFEandFields schur_fe_op_vec)
Create a Mat Diag Blocks object.
Definition Schur.cpp:1082
boost::shared_ptr< NestSchurData > createSchurNestedMatrixStruture(std::pair< SmartPetscObj< DM >, SmartPetscObj< DM > > dms, boost::shared_ptr< BlockStructure > block_mat_data_ptr, std::vector< std::string > fields_names, std::vector< boost::shared_ptr< Range > > field_ents, bool add_preconditioner_block)
Get the Schur Nest Mat Array object.
Definition Schur.cpp:2343
MoFEMErrorCode DMMoFEMSetNestSchurData(DM dm, boost::shared_ptr< NestSchurData >)
Definition DMMoFEM.cpp:1554
auto createDMNestSchurMat(DM dm)
Definition DMMoFEM.hpp:1086
auto createDM(MPI_Comm comm, const std::string dm_type_name)
Creates smart DM object.
MoFEMErrorCode assembleBlockMatSchur(MoFEM::Interface &m_field, Mat B, Mat S, std::vector< std::string > fields_name, std::vector< boost::shared_ptr< Range > > field_ents, SmartPetscObj< AO > ao)
Assemble Schur matrix.
Definition Schur.cpp:1817
OpSchurAssembleBase * createOpSchurAssembleBegin()
Definition Schur.cpp:2580
auto createDMBlockMat(DM dm)
Definition DMMoFEM.hpp:1079
constexpr AssemblyType A
constexpr double t
plate stiffness
Definition plate.cpp:58
FTensor::Index< 'm', 3 > m
static boost::shared_ptr< SetUpSchur > createSetUpSchur(MoFEM::Interface &m_field, EshelbianCore *ep_core_ptr)
std::vector< boost::shared_ptr< Range > > a00RangeList
const std::string skeletonElement
MoFEM::Interface & mField
const std::string spatialL2Disp
SmartPetscObj< IS > crackHybridIs
boost::shared_ptr< FaceElementForcesAndSourcesCore > elasticBcLhs
const std::string elementVolumeName
const std::string piolaStress
std::vector< std::string > a00FieldList
const std::string bubbleField
boost::shared_ptr< VolumeElementForcesAndSourcesCore > elasticFeLhs
static PetscBool noStretch
const std::string rotAxis
const std::string contactDisp
const std::string naturalBcElement
boost::shared_ptr< Range > crackFaces
const std::string hybridSpatialDisp
SmartPetscObj< DM > dmElastic
Elastic problem.
const std::string stretchTensor
const std::string contactElement
virtual MPI_Comm & get_comm() const =0
Deprecated interface functions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
P_MultiGridData(SmartPetscObj< DM > dm, SmartPetscObj< PC > pc, SmartPetscObj< Mat > S)
EshelbianCore * epCorePtr
boost::shared_ptr< std::vector< unsigned char > > piolaZeroDofsMarker
SmartPetscObj< Mat > S
boost::shared_ptr< std::vector< boost::weak_ptr< NumeredDofEntity > > > piolaZeroDofsVec
boost::shared_ptr< P_MultiGridData > pMGPtr
MoFEMErrorCode setUp(TS ts)
MoFEMErrorCode setUp(SmartPetscObj< KSP >)
virtual ~SetUpSchurImpl()
MoFEMErrorCode postProc()
SmartPetscObj< AO > aoS
MoFEMErrorCode preProc()
MoFEM::Interface & mField
SetUpSchurImpl(MoFEM::Interface &m_field, EshelbianCore *ep_core_ptr)
[Push operators to pipeline]