v0.14.0
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. More...
 
MoFEMErrorCode DMMGViaApproxOrdersGetCoarseningISSize (DM dm, int *size)
 Gets size of coarseningIS in internal data struture DMMGViaApproxOrdersCtx. More...
 
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. More...
 
MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS (DM)
 Pop is form MG via approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS (DM)
 Clear approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersReplaceCoarseningIS (DM dm, IS *is_vec, int nb_elems, Mat A, int verb=0)
 Replace coarsening IS in DM via approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersGetCtx (DM dm, const DMMGViaApproxOrdersCtx **ctx)
 Get context for DM via approximation orders. More...
 
MoFEMErrorCode DMRegister_MGViaApproxOrders (const char sname[])
 Register DM for Multi-Grid via approximation orders. More...
 
MoFEMErrorCode DMCreate_MGViaApproxOrders (DM dm)
 Create DM data structure for Multi-Grid via approximation orders. More...
 
MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders (DM dm, Mat *M)
 Create matrix for Multi-Grid via approximation orders. More...
 
MoFEMErrorCode DMCoarsen_MGViaApproxOrders (DM dm, MPI_Comm comm, DM *dmc)
 Coarsen DM. More...
 
MoFEMErrorCode DMCreateInterpolation_MGViaApproxOrders (DM dm1, DM dm2, Mat *mat, Vec *vec)
 Create interpolation matrix between data managers dm1 and dm2. More...
 
MoFEMErrorCode DMCreateGlobalVector_MGViaApproxOrders (DM dm, Vec *g)
 Create global vector for DMGViaApproxOrders. More...
 
MoFEMErrorCode PCMGSetUpViaApproxOrders (PC pc, PCMGSetUpViaApproxOrdersCtx *ctx, int verb=0)
 Function build MG structure. More...
 

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);
411  MoFEMFunctionBegin;
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;
418  CHKERR DMSetOperators_MoFEM(dm);
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);
427  MoFEMFunctionReturn(0);
428 }

◆ 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);
579  MoFEMFunctionBegin;
580  int leveldown = dm->leveldown;
581  GET_DM_FIELD(dm);
582  if (dm_field->kspOperators.empty()) {
583  CHKERR DMCreateGlobalVector_MoFEM(dm, g);
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);
598  MoFEMFunctionReturn(0);
599 }

◆ 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);
481  MoFEMFunctionBegin;
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 
573  MoFEMFunctionReturn(0);
574 }

◆ 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);
200  MoFEMFunctionBeginHot;
201  GET_DM_FIELD(dm);
202  *size = dm_field->coarseningIS.size();
203  MoFEMFunctionReturnHot(0);
204 }

◆ 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);
370  MoFEMFunctionBeginHot;
371  GET_DM_FIELD(dm);
372  *ctx = dm_field;
373  MoFEMFunctionReturnHot(0);
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);
177  MoFEMFunctionBeginHot;
178  GET_DM_FIELD(dm);
179  *ctx = dm_field;
180  MoFEMFunctionReturnHot(0);
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);
289  MoFEMFunctionBegin;
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++) {
350  CHKERR DMMGViaApproxOrdersPopBackCoarseningIS(dm);
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");
364  MoFEMFunctionReturn(0);
365 }

◆ 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);
185  MoFEMFunctionBegin;
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";
195  MoFEMFunctionReturn(0);
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);
732  MoFEMFunctionBegin;
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 
755  MoFEMFunctionReturn(0);
756 }
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
g
constexpr double g
Definition: shallow_wave.cpp:63
PCMGSetUpViaApproxOrdersCtx::getOptions
virtual MoFEMErrorCode getOptions()
get options from line command
Definition: PCMGSetUpViaApproxOrders.cpp:601
DMCoarsen_MGViaApproxOrders
MoFEMErrorCode DMCoarsen_MGViaApproxOrders(DM dm, MPI_Comm comm, DM *dmc)
Coarsen DM.
Definition: PCMGSetUpViaApproxOrders.cpp:462
DMMGViaApproxOrdersCtx
Structure for DM for multi-grid via approximation orders.
Definition: PCMGSetUpViaApproxOrders.hpp:22
A
constexpr AssemblyType A
Definition: operators_tests.cpp:30
DMMGViaApproxOrdersPushBackCoarseningIS
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.
Definition: PCMGSetUpViaApproxOrders.cpp:206
PCMGSetUpViaApproxOrdersCtx::nbLevels
int nbLevels
number of multi-grid levels
Definition: PCMGSetUpViaApproxOrders.hpp:204
GET_DM_FIELD
#define GET_DM_FIELD(DM)
Definition: PCMGSetUpViaApproxOrders.cpp:170
MoFEM::DMCreateGlobalVector_MoFEM
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMoFEM.cpp:1171
PlasticOps::M
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:117
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
MoFEM::DMCtx
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:944
DMCreateMatrix_MGViaApproxOrders
MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders(DM dm, Mat *M)
Create matrix for Multi-Grid via approximation orders.
Definition: PCMGSetUpViaApproxOrders.cpp:430
DMMGViaApproxOrdersPopBackCoarseningIS
MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS(DM dm)
Pop is form MG via approximation orders.
Definition: PCMGSetUpViaApproxOrders.cpp:260
MoFEM::DMSetOperators_MoFEM
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMoFEM.cpp:53
DMCreateGlobalVector_MGViaApproxOrders
MoFEMErrorCode DMCreateGlobalVector_MGViaApproxOrders(DM dm, Vec *g)
Create global vector for DMGViaApproxOrders.
Definition: PCMGSetUpViaApproxOrders.cpp:576
convert.n
n
Definition: convert.py:82
N
const int N
Definition: speed_test.cpp:3
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
ksp_set_operators
static MoFEMErrorCode ksp_set_operators(KSP ksp, Mat A, Mat B, void *ctx)
Definition: PCMGSetUpViaApproxOrders.cpp:404
PCMGSetUpViaApproxOrdersCtx::buildProlongationOperator
virtual MoFEMErrorCode buildProlongationOperator(bool use_mat_a, int verb=0)
Set up data structures for MG.
Definition: PCMGSetUpViaApproxOrders.cpp:660
DMCreateInterpolation_MGViaApproxOrders
MoFEMErrorCode DMCreateInterpolation_MGViaApproxOrders(DM dm1, DM dm2, Mat *mat, Vec *vec)
Create interpolation matrix between data managers dm1 and dm2.
Definition: PCMGSetUpViaApproxOrders.cpp:477
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
Generated by Doxygen 1.8.17 and hosted at University of Glasgow