v0.14.0
Loading...
Searching...
No Matches
Classes | Functions
PCMGSetUpViaApproxOrders.hpp File Reference

header of multi-grid solver for p- adaptivity More...

Go to the source code of this file.

Classes

struct  PCMGSubMatrixCtx
 
struct  DMMGViaApproxOrdersCtx
 Structure for DM for multi-grid via approximation orders. More...
 
struct  PCMGSetUpViaApproxOrdersCtx
 Set data structures of MG pre-conditioner via approximation orders. More...
 

Functions

MoFEMErrorCode DMMGViaApproxOrdersGetCtx (DM dm, DMMGViaApproxOrdersCtx **ctx)
 
MoFEMErrorCode DMMGViaApproxOrdersSetAO (DM dm, AO ao)
 Set DM ordering.
 
MoFEMErrorCode DMMGViaApproxOrdersGetCoarseningISSize (DM dm, int *size)
 Gets size of coarseningIS in internal data struture DMMGViaApproxOrdersCtx.
 
MoFEMErrorCode DMMGViaApproxOrdersPushBackCoarseningIS (DM, IS is, Mat A, Mat *subA, bool create_sub_matrix, bool shell_sub_a)
 Push back coarsening level to MG via approximation orders.
 
MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS (DM)
 Pop is form MG via approximation orders.
 
MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS (DM)
 Clear approximation orders.
 
MoFEMErrorCode DMMGViaApproxOrdersReplaceCoarseningIS (DM dm, IS *is_vec, int nb_elems, Mat A, int verb=0)
 Replace coarsening IS in DM via approximation orders.
 
MoFEMErrorCode DMMGViaApproxOrdersGetCtx (DM dm, const DMMGViaApproxOrdersCtx **ctx)
 Get context for DM via approximation orders.
 
MoFEMErrorCode DMRegister_MGViaApproxOrders (const char sname[])
 Register DM for Multi-Grid via approximation orders.
 
MoFEMErrorCode DMCreate_MGViaApproxOrders (DM dm)
 Create DM data structure for Multi-Grid via approximation orders.
 
MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders (DM dm, Mat *M)
 Create matrix for Multi-Grid via approximation orders.
 
MoFEMErrorCode DMCoarsen_MGViaApproxOrders (DM dm, MPI_Comm comm, DM *dmc)
 Coarsen DM.
 
MoFEMErrorCode DMCreateInterpolation_MGViaApproxOrders (DM dm1, DM dm2, Mat *mat, Vec *vec)
 Create interpolation matrix between data managers dm1 and dm2.
 
MoFEMErrorCode DMCreateGlobalVector_MGViaApproxOrders (DM dm, Vec *g)
 Create global vector for DMGViaApproxOrders.
 
MoFEMErrorCode PCMGSetUpViaApproxOrders (PC pc, PCMGSetUpViaApproxOrdersCtx *ctx, int verb=0)
 Function build MG structure.
 

Detailed Description

header of multi-grid solver for p- adaptivity

Definition in file PCMGSetUpViaApproxOrders.hpp.

Function Documentation

◆ DMCreate_MGViaApproxOrders()

MoFEMErrorCode DMCreate_MGViaApproxOrders ( DM dm)

Create DM data structure for Multi-Grid via approximation orders.

It set data structure and operators needed

Parameters
dmDiscrete manager
Returns
Error code

Definition at line 409 of file PCMGSetUpViaApproxOrders.cpp.

409 {
410 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
412 if (!dm->data) {
413 dm->data = new DMMGViaApproxOrdersCtx();
414 } else {
415 ((DMCtx *)(dm->data))->referenceNumber++;
416 }
417 // cerr << "Create " << ((DMCtx*)(dm->data))->referenceNumber << endl;
419 dm->ops->creatematrix = DMCreateMatrix_MGViaApproxOrders;
420 dm->ops->createglobalvector = DMCreateGlobalVector_MGViaApproxOrders;
421 dm->ops->coarsen = DMCoarsen_MGViaApproxOrders;
422 // dm->ops->destroy = DMDestroy_MGViaApproxOrders;
423 dm->ops->createinterpolation = DMCreateInterpolation_MGViaApproxOrders;
424 CHKERR DMKSPSetComputeOperators(dm, ksp_set_operators, NULL);
425 PetscInfo1(dm, "Create DMMGViaApproxOrders reference = %d\n",
426 ((DMCtx *)(dm->data))->referenceNumber);
428}
MoFEMErrorCode DMCreateGlobalVector_MGViaApproxOrders(DM dm, Vec *g)
Create global vector for DMGViaApproxOrders.
static MoFEMErrorCode ksp_set_operators(KSP ksp, Mat A, Mat B, void *ctx)
MoFEMErrorCode DMCreateInterpolation_MGViaApproxOrders(DM dm1, DM dm2, Mat *mat, Vec *vec)
Create interpolation matrix between data managers dm1 and dm2.
#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.
MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders(DM dm, Mat *M)
Create matrix for Multi-Grid via approximation orders.
MoFEMErrorCode DMCoarsen_MGViaApproxOrders(DM dm, MPI_Comm comm, DM *dmc)
Coarsen DM.
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition DMMoFEM.cpp:53
Structure for DM for multi-grid via approximation orders.
PETSc Discrete Manager data structure.
Definition DMMoFEM.hpp:929

◆ DMCreateGlobalVector_MGViaApproxOrders()

MoFEMErrorCode DMCreateGlobalVector_MGViaApproxOrders ( DM dm,
Vec * g )

Create global vector for DMGViaApproxOrders.

Parameters
dmDistributed mesh data structure
greturned pointer to vector
Returns
Error code

Definition at line 576 of file PCMGSetUpViaApproxOrders.cpp.

576 {
577
578 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
580 int leveldown = dm->leveldown;
581 GET_DM_FIELD(dm);
582 if (dm_field->kspOperators.empty()) {
584 } else {
585#if PETSC_VERSION_GE(3, 5, 3)
586 CHKERR MatCreateVecs(
587 dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown],
588 g, NULL);
589#else
590 CHKERR MatGetVecs(
591 dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown],
592 g, NULL);
593#endif
594 CHKERR VecSetDM(*g, dm);
595 }
596 PetscInfo1(dm, "Create global vector DMMGViaApproxOrders leveldown = %d\n",
597 dm->leveldown);
599}
#define GET_DM_FIELD(DM)
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition DMMoFEM.cpp:1153
constexpr double g

◆ DMCreateInterpolation_MGViaApproxOrders()

MoFEMErrorCode DMCreateInterpolation_MGViaApproxOrders ( DM dm1,
DM dm2,
Mat * mat,
Vec * vec )

Create interpolation matrix between data managers dm1 and dm2.

Parameters
dm1Distributed mesh data structure
dm2Distributed mesh data structure
matPointer to returned interpolation matrix
vecPointer to scaling vector here returned NULL
Returns
Error code

Definition at line 477 of file PCMGSetUpViaApproxOrders.cpp.

478 {
479 PetscValidHeaderSpecific(dm1, DM_CLASSID, 1);
480 PetscValidHeaderSpecific(dm2, DM_CLASSID, 1);
482
483 MPI_Comm comm;
484 CHKERR PetscObjectGetComm((PetscObject)dm1, &comm);
485
486 int m, n, M, N;
487
488 DM dm_down = dm1;
489 DM dm_up = dm2;
490
491 int dm_down_leveldown = dm_down->leveldown;
492 int dm_up_leveldown = dm_up->leveldown;
493
494 PetscInfo2(dm1,
495 "Create interpolation DMMGViaApproxOrders dm1_leveldown = %d "
496 "dm2_leveldown = %d\n",
497 dm_down_leveldown, dm_up_leveldown);
498
499 IS is_down, is_up;
500 {
501 // Coarser mesh
502 GET_DM_FIELD(dm_down);
503 if (static_cast<int>(dm_field->coarseningIS.size()) < dm_down_leveldown) {
504 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
505 }
506 is_down = dm_field->coarseningIS[dm_field->coarseningIS.size() - 1 -
507 dm_down_leveldown];
508 CHKERR ISGetSize(is_down, &M);
509 CHKERR ISGetLocalSize(is_down, &m);
510 }
511 {
512 // Finer mesh
513 GET_DM_FIELD(dm_up);
514 if (static_cast<int>(dm_field->coarseningIS.size()) < dm_up_leveldown) {
515 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
516 }
517 is_up =
518 dm_field
519 ->coarseningIS[dm_field->coarseningIS.size() - 1 - dm_up_leveldown];
520 CHKERR ISGetSize(is_up, &N);
521 CHKERR ISGetLocalSize(is_up, &n);
522 }
523
524 // is_dow rows
525 // is_up columns
526
527 CHKERR MatCreate(comm, mat);
528 CHKERR MatSetSizes(*mat, m, n, M, N);
529 CHKERR MatSetType(*mat, MATMPIAIJ);
530 CHKERR MatMPIAIJSetPreallocation(*mat, 1, PETSC_NULL, 0, PETSC_NULL);
531
532 // get matrix layout
533 PetscLayout rmap, cmap;
534 CHKERR MatGetLayouts(*mat, &rmap, &cmap);
535 int rstart, rend, cstart, cend;
536 CHKERR PetscLayoutGetRange(rmap, &rstart, &rend);
537 CHKERR PetscLayoutGetRange(cmap, &cstart, &cend);
538
539 // if(verb>0) {
540 // PetscSynchronizedPrintf(comm,"level %d row start %d row end
541 // %d\n",kk,rstart,rend); PetscSynchronizedPrintf(comm,"level %d col start
542 // %d col end %d\n",kk,cstart,cend);
543 // }
544
545 const int *row_indices_ptr, *col_indices_ptr;
546 CHKERR ISGetIndices(is_down, &row_indices_ptr);
547 CHKERR ISGetIndices(is_up, &col_indices_ptr);
548
549 map<int, int> idx_map;
550 for (int ii = 0; ii < m; ii++) {
551 idx_map[row_indices_ptr[ii]] = rstart + ii;
552 }
553
554 CHKERR MatZeroEntries(*mat);
555 // FIXME: Use MatCreateMPIAIJWithArrays and set array directly
556 for (int jj = 0; jj < n; jj++) {
557 map<int, int>::iterator mit = idx_map.find(col_indices_ptr[jj]);
558 if (mit != idx_map.end()) {
559 CHKERR MatSetValue(*mat, mit->second, cstart + jj, 1, INSERT_VALUES);
560 }
561 }
562
563 CHKERR ISRestoreIndices(is_down, &row_indices_ptr);
564 CHKERR ISRestoreIndices(is_up, &col_indices_ptr);
565
566 CHKERR MatAssemblyBegin(*mat, MAT_FINAL_ASSEMBLY);
567 CHKERR MatAssemblyEnd(*mat, MAT_FINAL_ASSEMBLY);
568
569 if (vec != NULL) {
570 *vec = PETSC_NULL;
571 }
572
574}
static Index< 'M', 3 > M
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
const int N
Definition speed_test.cpp:3

◆ DMMGViaApproxOrdersGetCoarseningISSize()

MoFEMErrorCode DMMGViaApproxOrdersGetCoarseningISSize ( DM dm,
int * size )

Gets size of coarseningIS in internal data struture DMMGViaApproxOrdersCtx.

Parameters
dmDM
sizesize of coarseningIS
Returns
Error code

Definition at line 198 of file PCMGSetUpViaApproxOrders.cpp.

198 {
199 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
201 GET_DM_FIELD(dm);
202 *size = dm_field->coarseningIS.size();
204}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...

◆ DMMGViaApproxOrdersGetCtx() [1/2]

MoFEMErrorCode DMMGViaApproxOrdersGetCtx ( DM dm,
const DMMGViaApproxOrdersCtx ** ctx )

Get context for DM via approximation orders.

Parameters
dmthe DM object
ctxdata context
Returns
error code

Definition at line 367 of file PCMGSetUpViaApproxOrders.cpp.

368 {
369 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
371 GET_DM_FIELD(dm);
372 *ctx = dm_field;
374}

◆ DMMGViaApproxOrdersGetCtx() [2/2]

MoFEMErrorCode DMMGViaApproxOrdersGetCtx ( DM dm,
DMMGViaApproxOrdersCtx ** ctx )

Get DM Ctx

Definition at line 175 of file PCMGSetUpViaApproxOrders.cpp.

175 {
176 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
178 GET_DM_FIELD(dm);
179 *ctx = dm_field;
181}

◆ DMMGViaApproxOrdersReplaceCoarseningIS()

MoFEMErrorCode DMMGViaApproxOrdersReplaceCoarseningIS ( DM dm,
IS * is_vec,
int nb_elems,
Mat A,
int verb = 0 )

Replace coarsening IS in DM via approximation orders.

Parameters
dmdm
is_vecPointer to vector of is
nb_elemsNumber of elements
AFine matrix
Returns
Error code

Definition at line 285 of file PCMGSetUpViaApproxOrders.cpp.

287 {
288 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
290 GET_DM_FIELD(dm);
291 int nb_no_changed = 0;
292 int nb_replaced = 0;
293 int nb_deleted = 0;
294 int nb_added = 0;
295 std::vector<IS>::iterator it;
296 it = dm_field->coarseningIS.begin();
297 int ii = 0;
298 for (; it != dm_field->coarseningIS.end(); it++, ii++) {
299 if (ii < nb_elems) {
300 PetscBool flg;
301 CHKERR ISEqual(*it, is_vec[ii], &flg);
302 if (!flg) {
303 CHKERR ISDestroy(&*it);
304 CHKERR MatDestroy(&dm_field->kspOperators[ii]);
305 *it = is_vec[ii];
306 CHKERR PetscObjectReference((PetscObject)is_vec[ii]);
307 if (ii < nb_elems - 1) {
308 IS is = is_vec[ii];
309 if (dm_field->aO) {
310 CHKERR ISDuplicate(is_vec[ii], &is);
311 CHKERR ISCopy(is_vec[ii], is);
312 CHKERR AOApplicationToPetscIS(dm_field->aO, is);
313 }
314 Mat subA;
315#if PETSC_VERSION_GE(3, 8, 0)
316 CHKERR MatCreateSubMatrix(A, is, is, MAT_INITIAL_MATRIX, &subA);
317#else
318 CHKERR MatGetSubMatrix(A, is, is, MAT_INITIAL_MATRIX, &subA);
319#endif
320 CHKERR PetscObjectReference((PetscObject)subA);
321 dm_field->kspOperators[ii] = subA;
322 CHKERR MatDestroy(&subA);
323 if (dm_field->aO) {
324 CHKERR ISDestroy(&is);
325 }
326 } else {
327 CHKERR PetscObjectReference((PetscObject)A);
328 dm_field->kspOperators[ii] = A;
329 }
330 nb_replaced++;
331 }
332 } else {
333 nb_no_changed++;
334 continue;
335 }
336 }
337 if (static_cast<int>(dm_field->coarseningIS.size()) < nb_elems) {
338 for (; ii < nb_elems - 1; ii++) {
339 Mat subA;
340 CHKERR DMMGViaApproxOrdersPushBackCoarseningIS(dm, is_vec[ii], A, &subA,
341 true, false);
342 CHKERR MatDestroy(&subA);
343 nb_added++;
344 }
345 CHKERR DMMGViaApproxOrdersPushBackCoarseningIS(dm, is_vec[ii], A, &A, false,
346 false);
347 nb_added++;
348 } else {
349 for (; ii < static_cast<int>(dm_field->coarseningIS.size()); ii++) {
351 nb_deleted++;
352 }
353 }
354 MPI_Comm comm;
355 CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
356 if (verb > 0) {
357 PetscPrintf(comm,
358 "DMMGViaApproxOrders nb_no_changed = %d, nb_replaced = %d, "
359 "nb_added = %d, nb_deleted = %d, size = %d\n",
360 nb_no_changed, nb_replaced, nb_added, nb_deleted,
361 dm_field->coarseningIS.size());
362 }
363 PetscInfo(dm, "Replace IS to DMMGViaApproxOrders\n");
365}
MoFEMErrorCode DMMGViaApproxOrdersPushBackCoarseningIS(DM dm, IS is, Mat A, Mat *subA, bool create_sub_matrix, bool shell_sub_a)
Push back coarsening level to MG via approximation orders.
MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS(DM dm)
Pop is form MG via approximation orders.
constexpr AssemblyType A

◆ DMMGViaApproxOrdersSetAO()

MoFEMErrorCode DMMGViaApproxOrdersSetAO ( DM dm,
AO ao )

Set DM ordering.

IS can be given is some other ordering, AO will transform indices from coarseningIS ordering to ordering used to construct fine matrix.

Parameters
dm[description]
ao[description]
Returns
[description]

Definition at line 183 of file PCMGSetUpViaApproxOrders.cpp.

183 {
184 PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
186 GET_DM_FIELD(dm);
187 if (dm_field->aO) {
188 // std::cerr << dm_field->aO << std::endl;
189 CHKERR AODestroy(&dm_field->aO);
190 // std::cerr << "destroy ao when adding\n";
191 }
192 dm_field->aO = ao;
193 CHKERR PetscObjectReference((PetscObject)ao);
194 // std::cerr << "add ao\n";
196}

◆ PCMGSetUpViaApproxOrders()

MoFEMErrorCode PCMGSetUpViaApproxOrders ( PC pc,
PCMGSetUpViaApproxOrdersCtx * ctx,
int verb = 0 )

Function build MG structure.

Parameters
pcMG pre-conditioner http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/PC/PCMG.html
ctxMoFEM data structure for MG
verbverbosity level
Returns
error code

Definition at line 728 of file PCMGSetUpViaApproxOrders.cpp.

729 {
730
731 PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
733
734 MPI_Comm comm;
735 CHKERR PetscObjectGetComm((PetscObject)pc, &comm);
736 if (verb > 0) {
737 PetscPrintf(comm, "Start PCMGSetUpViaApproxOrders\n");
738 }
739
740 CHKERR ctx->getOptions();
741 CHKERR ctx->buildProlongationOperator(true, verb);
742
743#if PETSC_VERSION_GE(3, 8, 0)
744 CHKERR PCMGSetGalerkin(pc, PC_MG_GALERKIN_NONE);
745#else
746 CHKERR PCMGSetGalerkin(pc, PETSC_FALSE);
747#endif
748
749 CHKERR PCMGSetLevels(pc, ctx->nbLevels, NULL);
750
751 if (verb > 0) {
752 PetscPrintf(comm, "End PCMGSetUpViaApproxOrders\n");
753 }
754
756}
virtual MoFEMErrorCode getOptions()
get options from line command
int nbLevels
number of multi-grid levels
virtual MoFEMErrorCode buildProlongationOperator(bool use_mat_a, int verb=0)
Set up data structures for MG.