v0.8.16
Classes | Macros | Functions
PCMGSetUpViaApproxOrders.cpp File Reference

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

#include <MoFEM.hpp>
#include <UnknownInterface.hpp>
#include <PCMGSetUpViaApproxOrders.hpp>
#include <petsc-private/petscimpl.h>
#include <petsc-private/dmimpl.h>
#include <petsc-private/vecimpl.h>

Go to the source code of this file.

Classes

struct  PCMGSubMatrixCtx_private
 

Macros

#define PETSC_VERSION_RELEASE   1
 
#define GET_DM_FIELD(DM)
 

Functions

template<InsertMode MODE>
MoFEMErrorCode sub_mat_mult_generic (Mat a, Vec x, Vec f)
 
MoFEMErrorCode sub_mat_mult (Mat a, Vec x, Vec f)
 
MoFEMErrorCode sub_mat_mult_add (Mat a, Vec x, Vec f)
 
MoFEMErrorCode sub_mat_sor (Mat mat, Vec b, PetscReal omega, MatSORType flag, PetscReal shift, PetscInt its, PetscInt lits, Vec x)
 
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 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 dm)
 Pop is form MG via approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersClearCoarseningIS (DM dm)
 Clear approximation orders. More...
 
MoFEMErrorCode DMMGViaApproxOrdersReplaceCoarseningIS (DM dm, IS *is_vec, int nb_elems, Mat A, int verb)
 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...
 
static MoFEMErrorCode ksp_set_operators (KSP ksp, Mat A, Mat B, void *ctx)
 
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)
 Function build MG structure. More...
 

Detailed Description

implementation of multi-grid solver for p- adaptivity

MoFEM is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

MoFEM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with MoFEM. If not, see http://www.gnu.org/licenses/

Definition in file PCMGSetUpViaApproxOrders.cpp.

Macro Definition Documentation

◆ GET_DM_FIELD

#define GET_DM_FIELD (   DM)
Value:
CHKERR((DMCtx *)DM->data) \
DMMGViaApproxOrdersCtx *dm_field = \
static_cast<DMMGViaApproxOrdersCtx *>(iface); \
NOT_USED(dm_field)
base class for all interface classes
static const MOFEMuuid IDD_DMMGVIAAPPROXORDERSCTX
Structure for DM for multi-grid via approximation orders.
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:737
virtual MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const =0

Definition at line 191 of file PCMGSetUpViaApproxOrders.cpp.

◆ PETSC_VERSION_RELEASE

#define PETSC_VERSION_RELEASE   1

Definition at line 27 of file PCMGSetUpViaApproxOrders.cpp.

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 425 of file PCMGSetUpViaApproxOrders.cpp.

425  {
426  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
428  if (!dm->data) {
429  dm->data = new DMMGViaApproxOrdersCtx();
430  } else {
431  ((DMCtx *)(dm->data))->referenceNumber++;
432  }
433  // cerr << "Create " << ((DMCtx*)(dm->data))->referenceNumber << endl;
435  dm->ops->creatematrix = DMCreateMatrix_MGViaApproxOrders;
436  dm->ops->createglobalvector = DMCreateGlobalVector_MGViaApproxOrders;
437  dm->ops->coarsen = DMCoarsen_MGViaApproxOrders;
438  // dm->ops->destroy = DMDestroy_MGViaApproxOrders;
439  dm->ops->createinterpolation = DMCreateInterpolation_MGViaApproxOrders;
440  CHKERR DMKSPSetComputeOperators(dm, ksp_set_operators, NULL);
441  PetscInfo1(dm, "Create DMMGViaApproxOrders reference = %d\n",
442  ((DMCtx *)(dm->data))->referenceNumber);
444 }
MoFEMErrorCode DMCoarsen_MGViaApproxOrders(DM dm, MPI_Comm comm, DM *dmc)
Coarsen DM.
static MoFEMErrorCode ksp_set_operators(KSP ksp, Mat A, Mat B, void *ctx)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
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.
Structure for DM for multi-grid via approximation orders.
MoFEMErrorCode DMCreateMatrix_MGViaApproxOrders(DM dm, Mat *M)
Create matrix for Multi-Grid via approximation orders.
PetscErrorCode DMSetOperators_MoFEM(DM dm)
Set operators for MoFEM dm.
Definition: DMMMoFEM.cpp:97
#define CHKERR
Inline error check.
Definition: definitions.h:578
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:737
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ 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 591 of file PCMGSetUpViaApproxOrders.cpp.

591  {
592 
593  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
595  int leveldown = dm->leveldown;
596  GET_DM_FIELD(dm);
597  if (dm_field->kspOperators.empty()) {
599  } else {
600 #if PETSC_VERSION_GE(3, 5, 3)
601  CHKERR MatCreateVecs(
602  dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown],
603  g, NULL);
604 #else
605  CHKERR MatGetVecs(
606  dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown],
607  g, NULL);
608 #endif
609  }
610  PetscInfo1(dm, "Create global vector DMMGViaApproxOrders leveldown = %d\n",
611  dm->leveldown);
613 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVectorsets the routine to create a global vector associated with the shell DM...
Definition: DMMMoFEM.cpp:889
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ 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 492 of file PCMGSetUpViaApproxOrders.cpp.

493  {
494  PetscValidHeaderSpecific(dm1, DM_CLASSID, 1);
495  PetscValidHeaderSpecific(dm2, DM_CLASSID, 1);
497 
498  MPI_Comm comm;
499  CHKERR PetscObjectGetComm((PetscObject)dm1, &comm);
500 
501  int m, n, M, N;
502 
503  DM dm_down = dm1;
504  DM dm_up = dm2;
505 
506  int dm_down_leveldown = dm_down->leveldown;
507  int dm_up_leveldown = dm_up->leveldown;
508 
509  PetscInfo2(dm1,
510  "Create interpolation DMMGViaApproxOrders dm1_leveldown = %d "
511  "dm2_leveldown = %d\n",
512  dm_down_leveldown, dm_up_leveldown);
513 
514  IS is_down, is_up;
515  {
516  // Coarser mesh
517  GET_DM_FIELD(dm_down);
518  if (static_cast<int>(dm_field->coarseningIS.size()) < dm_down_leveldown) {
519  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
520  }
521  is_down = dm_field->coarseningIS[dm_field->coarseningIS.size() - 1 -
522  dm_down_leveldown];
523  CHKERR ISGetSize(is_down, &M);
524  CHKERR ISGetLocalSize(is_down, &m);
525  }
526  {
527  // Finer mesh
528  GET_DM_FIELD(dm_up);
529  if (static_cast<int>(dm_field->coarseningIS.size()) < dm_up_leveldown) {
530  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
531  }
532  is_up =
533  dm_field
534  ->coarseningIS[dm_field->coarseningIS.size() - 1 - dm_up_leveldown];
535  CHKERR ISGetSize(is_up, &N);
536  CHKERR ISGetLocalSize(is_up, &n);
537  }
538 
539  // is_dow rows
540  // is_up columns
541 
542  CHKERR MatCreate(comm, mat);
543  CHKERR MatSetSizes(*mat, m, n, M, N);
544  CHKERR MatSetType(*mat, MATMPIAIJ);
545  CHKERR MatMPIAIJSetPreallocation(*mat, 1, PETSC_NULL, 0, PETSC_NULL);
546 
547  // get matrix layout
548  PetscLayout rmap, cmap;
549  CHKERR MatGetLayouts(*mat, &rmap, &cmap);
550  int rstart, rend, cstart, cend;
551  CHKERR PetscLayoutGetRange(rmap, &rstart, &rend);
552  CHKERR PetscLayoutGetRange(cmap, &cstart, &cend);
553 
554  // if(verb>0) {
555  // PetscSynchronizedPrintf(comm,"level %d row start %d row end
556  // %d\n",kk,rstart,rend); PetscSynchronizedPrintf(comm,"level %d col start
557  // %d col end %d\n",kk,cstart,cend);
558  // }
559 
560  const int *row_indices_ptr, *col_indices_ptr;
561  CHKERR ISGetIndices(is_down, &row_indices_ptr);
562  CHKERR ISGetIndices(is_up, &col_indices_ptr);
563 
564  map<int, int> idx_map;
565  for (int ii = 0; ii < m; ii++) {
566  idx_map[row_indices_ptr[ii]] = rstart + ii;
567  }
568 
569  CHKERR MatZeroEntries(*mat);
570  // FIXME: Use MatCreateMPIAIJWithArrays and set array directly
571  for (int jj = 0; jj < n; jj++) {
572  map<int, int>::iterator mit = idx_map.find(col_indices_ptr[jj]);
573  if (mit != idx_map.end()) {
574  CHKERR MatSetValue(*mat, mit->second, cstart + jj, 1, INSERT_VALUES);
575  }
576  }
577 
578  CHKERR ISRestoreIndices(is_down, &row_indices_ptr);
579  CHKERR ISRestoreIndices(is_up, &col_indices_ptr);
580 
581  CHKERR MatAssemblyBegin(*mat, MAT_FINAL_ASSEMBLY);
582  CHKERR MatAssemblyEnd(*mat, MAT_FINAL_ASSEMBLY);
583 
584  if (vec != NULL) {
585  *vec = PETSC_NULL;
586  }
587 
589 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
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 222 of file PCMGSetUpViaApproxOrders.cpp.

222  {
223  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
225  GET_DM_FIELD(dm);
226  *size = dm_field->coarseningIS.size();
228 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
#define GET_DM_FIELD(DM)

◆ DMMGViaApproxOrdersGetCtx() [1/2]

MoFEMErrorCode DMMGViaApproxOrdersGetCtx ( DM  dm,
DMMGViaApproxOrdersCtx **  ctx 
)

Get DM Ctx

Definition at line 199 of file PCMGSetUpViaApproxOrders.cpp.

199  {
200  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
202  GET_DM_FIELD(dm);
203  *ctx = dm_field;
205 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
#define GET_DM_FIELD(DM)

◆ DMMGViaApproxOrdersGetCtx() [2/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 383 of file PCMGSetUpViaApproxOrders.cpp.

384  {
385  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
387  GET_DM_FIELD(dm);
388  *ctx = dm_field;
390 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
#define GET_DM_FIELD(DM)

◆ 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 305 of file PCMGSetUpViaApproxOrders.cpp.

307  {
308  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
310  GET_DM_FIELD(dm);
311  int nb_no_changed = 0;
312  int nb_replaced = 0;
313  int nb_deleted = 0;
314  int nb_added = 0;
315  std::vector<IS>::iterator it;
316  it = dm_field->coarseningIS.begin();
317  int ii = 0;
318  for (; it != dm_field->coarseningIS.end(); it++, ii++) {
319  if (ii < nb_elems) {
320  PetscBool flg;
321  CHKERR ISEqual(*it, is_vec[ii], &flg);
322  if (!flg) {
323  CHKERR ISDestroy(&*it);
324  CHKERR MatDestroy(&dm_field->kspOperators[ii]);
325  *it = is_vec[ii];
326  CHKERR PetscObjectReference((PetscObject)is_vec[ii]);
327  if (ii < nb_elems - 1) {
328  IS is = is_vec[ii];
329  if (dm_field->aO) {
330  CHKERR ISDuplicate(is_vec[ii], &is);
331  CHKERR ISCopy(is_vec[ii], is);
332  CHKERR AOApplicationToPetscIS(dm_field->aO, is);
333  }
334  Mat subA;
335  CHKERR MatGetSubMatrix(A, is, is, MAT_INITIAL_MATRIX, &subA);
336  CHKERR PetscObjectReference((PetscObject)subA);
337  dm_field->kspOperators[ii] = subA;
338  CHKERR MatDestroy(&subA);
339  if (dm_field->aO) {
340  CHKERR ISDestroy(&is);
341  }
342  } else {
343  CHKERR PetscObjectReference((PetscObject)A);
344  dm_field->kspOperators[ii] = A;
345  }
346  nb_replaced++;
347  }
348  } else {
349  nb_no_changed++;
350  continue;
351  }
352  }
353  if (static_cast<int>(dm_field->coarseningIS.size()) < nb_elems) {
354  for (; ii < nb_elems - 1; ii++) {
355  Mat subA;
356  CHKERR DMMGViaApproxOrdersPushBackCoarseningIS(dm, is_vec[ii], A, &subA,
357  true, false);
358  CHKERR MatDestroy(&subA);
359  nb_added++;
360  }
361  CHKERR DMMGViaApproxOrdersPushBackCoarseningIS(dm, is_vec[ii], A, &A, false,
362  false);
363  nb_added++;
364  } else {
365  for (; ii < static_cast<int>(dm_field->coarseningIS.size()); ii++) {
367  nb_deleted++;
368  }
369  }
370  MPI_Comm comm;
371  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
372  if (verb > 0) {
373  PetscPrintf(comm,
374  "DMMGViaApproxOrders nb_no_changed = %d, nb_replaced = %d, "
375  "nb_added = %d, nb_deleted = %d, size = %d\n",
376  nb_no_changed, nb_replaced, nb_added, nb_deleted,
377  dm_field->coarseningIS.size());
378  }
379  PetscInfo(dm, "Replace IS to DMMGViaApproxOrders\n");
381 }
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.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define GET_DM_FIELD(DM)
MoFEMErrorCode DMMGViaApproxOrdersPopBackCoarseningIS(DM dm)
Pop is form MG via approximation orders.
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ 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 207 of file PCMGSetUpViaApproxOrders.cpp.

207  {
208  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
210  GET_DM_FIELD(dm);
211  if (dm_field->aO) {
212  // std::cerr << dm_field->aO << std::endl;
213  CHKERR AODestroy(&dm_field->aO);
214  // std::cerr << "destroy ao when adding\n";
215  }
216  dm_field->aO = ao;
217  CHKERR PetscObjectReference((PetscObject)ao);
218  // std::cerr << "add ao\n";
220 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define GET_DM_FIELD(DM)
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ ksp_set_operators()

static MoFEMErrorCode ksp_set_operators ( KSP  ksp,
Mat  A,
Mat  B,
void *  ctx 
)
static

Definition at line 420 of file PCMGSetUpViaApproxOrders.cpp.

420  {
423 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490

◆ 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
Examples:
elasticity.cpp.

Definition at line 742 of file PCMGSetUpViaApproxOrders.cpp.

743  {
744 
745  PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
747 
748  MPI_Comm comm;
749  CHKERR PetscObjectGetComm((PetscObject)pc, &comm);
750  if (verb > 0) {
751  PetscPrintf(comm, "Start PCMGSetUpViaApproxOrders\n");
752  }
753 
754  CHKERR ctx->getOptions();
755  CHKERR ctx->buildProlongationOperator(true, verb);
756 
757 #if PETSC_VERSION_GE(3, 8, 0)
758  CHKERR PCMGSetGalerkin(pc, PC_MG_GALERKIN_NONE);
759 #else
760  CHKERR PCMGSetGalerkin(pc, PETSC_FALSE);
761 #endif
762 
763  CHKERR PCMGSetLevels(pc, ctx->nbLevels, NULL);
764 
765  if (verb > 0) {
766  PetscPrintf(comm, "End PCMGSetUpViaApproxOrders\n");
767  }
768 
770 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
int nbLevels
number of multi-grid levels
virtual MoFEMErrorCode buildProlongationOperator(bool use_mat_a, int verb=0)
Set up data structures for MG.
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
virtual MoFEMErrorCode getOptions()
get options from line command

◆ sub_mat_mult()

MoFEMErrorCode sub_mat_mult ( Mat  a,
Vec  x,
Vec  f 
)

Definition at line 117 of file PCMGSetUpViaApproxOrders.cpp.

117  {
118  return sub_mat_mult_generic<INSERT_VALUES>(a, x, f);
119 }

◆ sub_mat_mult_add()

MoFEMErrorCode sub_mat_mult_add ( Mat  a,
Vec  x,
Vec  f 
)

Definition at line 121 of file PCMGSetUpViaApproxOrders.cpp.

121  {
122  return sub_mat_mult_generic<ADD_VALUES>(a, x, f);
123 }

◆ sub_mat_mult_generic()

template<InsertMode MODE>
MoFEMErrorCode sub_mat_mult_generic ( Mat  a,
Vec  x,
Vec  f 
)

Definition at line 99 of file PCMGSetUpViaApproxOrders.cpp.

99  {
100  void *void_ctx;
102  CHKERR MatShellGetContext(a, &void_ctx);
104  if (!ctx->isInitisalised) {
105  CHKERR ctx->initData(x);
106  }
107  PetscLogEventBegin(ctx->MOFEM_EVENT_mult, 0, 0, 0, 0);
108  CHKERR VecScatterBegin(ctx->sCat, x, ctx->X, INSERT_VALUES, SCATTER_REVERSE);
109  CHKERR VecScatterEnd(ctx->sCat, x, ctx->X, INSERT_VALUES, SCATTER_REVERSE);
110  CHKERR MatMult(ctx->A, ctx->X, ctx->F);
111  CHKERR VecScatterBegin(ctx->sCat, ctx->F, f, MODE, SCATTER_FORWARD);
112  CHKERR VecScatterEnd(ctx->sCat, ctx->F, f, MODE, SCATTER_FORWARD);
113  PetscLogEventEnd(ctx->MOFEM_EVENT_mult, 0, 0, 0, 0);
115 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ sub_mat_sor()

MoFEMErrorCode sub_mat_sor ( Mat  mat,
Vec  b,
PetscReal  omega,
MatSORType  flag,
PetscReal  shift,
PetscInt  its,
PetscInt  lits,
Vec  x 
)

Definition at line 125 of file PCMGSetUpViaApproxOrders.cpp.

127  {
128  void *void_ctx;
130  CHKERR MatShellGetContext(mat, &void_ctx);
132  if (!ctx->isInitisalised) {
133  CHKERR ctx->initData(x);
134  }
135  PetscLogEventBegin(ctx->MOFEM_EVENT_sor, 0, 0, 0, 0);
136  CHKERR VecScatterBegin(ctx->sCat, b, ctx->X, INSERT_VALUES, SCATTER_REVERSE);
137  CHKERR VecScatterEnd(ctx->sCat, b, ctx->X, INSERT_VALUES, SCATTER_REVERSE);
138  CHKERR MatSOR(ctx->A, ctx->X, omega, flag, shift, its, lits, ctx->F);
139  CHKERR VecScatterBegin(ctx->sCat, ctx->F, x, INSERT_VALUES, SCATTER_FORWARD);
140  CHKERR VecScatterEnd(ctx->sCat, ctx->F, x, INSERT_VALUES, SCATTER_FORWARD);
141  PetscLogEventEnd(ctx->MOFEM_EVENT_sor, 0, 0, 0, 0);
143 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403