v0.13.0
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:
auto dm_field = \
static_cast<DMCtx *>(DM->data)->getInterface<DMMGViaApproxOrdersCtx>(); \
NOT_USED(dm_field)
Structure for DM for multi-grid via approximation orders.

Definition at line 181 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 420 of file PCMGSetUpViaApproxOrders.cpp.

420  {
421  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
423  if (!dm->data) {
424  dm->data = new DMMGViaApproxOrdersCtx();
425  } else {
426  ((DMCtx *)(dm->data))->referenceNumber++;
427  }
428  // cerr << "Create " << ((DMCtx*)(dm->data))->referenceNumber << endl;
430  dm->ops->creatematrix = DMCreateMatrix_MGViaApproxOrders;
431  dm->ops->createglobalvector = DMCreateGlobalVector_MGViaApproxOrders;
432  dm->ops->coarsen = DMCoarsen_MGViaApproxOrders;
433  // dm->ops->destroy = DMDestroy_MGViaApproxOrders;
434  dm->ops->createinterpolation = DMCreateInterpolation_MGViaApproxOrders;
435  CHKERR DMKSPSetComputeOperators(dm, ksp_set_operators, NULL);
436  PetscInfo1(dm, "Create DMMGViaApproxOrders reference = %d\n",
437  ((DMCtx *)(dm->data))->referenceNumber);
439 }
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 ...
Definition: definitions.h:359
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
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: DMMMoFEM.cpp:65
PETSc Discrete Manager data structure.
Definition: DMMoFEM.hpp:894

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

587  {
588 
589  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
591  int leveldown = dm->leveldown;
592  GET_DM_FIELD(dm);
593  if (dm_field->kspOperators.empty()) {
595  } else {
596 #if PETSC_VERSION_GE(3, 5, 3)
597  CHKERR MatCreateVecs(
598  dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown],
599  g, NULL);
600 #else
601  CHKERR MatGetVecs(
602  dm_field->kspOperators[dm_field->kspOperators.size() - 1 - leveldown],
603  g, NULL);
604 #endif
605  CHKERR VecSetDM(*g, dm);
606  }
607  PetscInfo1(dm, "Create global vector DMMGViaApproxOrders leveldown = %d\n",
608  dm->leveldown);
610 }
#define GET_DM_FIELD(DM)
PetscErrorCode DMCreateGlobalVector_MoFEM(DM dm, Vec *g)
DMShellSetCreateGlobalVector.
Definition: DMMMoFEM.cpp:1105
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 488 of file PCMGSetUpViaApproxOrders.cpp.

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

209  {
210  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
212  GET_DM_FIELD(dm);
213  *size = dm_field->coarseningIS.size();
215 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453

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

379  {
380  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
382  GET_DM_FIELD(dm);
383  *ctx = dm_field;
385 }

◆ DMMGViaApproxOrdersGetCtx() [2/2]

MoFEMErrorCode DMMGViaApproxOrdersGetCtx ( DM  dm,
DMMGViaApproxOrdersCtx **  ctx 
)

Get DM Ctx

Definition at line 186 of file PCMGSetUpViaApproxOrders.cpp.

186  {
187  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
189  GET_DM_FIELD(dm);
190  *ctx = dm_field;
192 }

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

298  {
299  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
301  GET_DM_FIELD(dm);
302  int nb_no_changed = 0;
303  int nb_replaced = 0;
304  int nb_deleted = 0;
305  int nb_added = 0;
306  std::vector<IS>::iterator it;
307  it = dm_field->coarseningIS.begin();
308  int ii = 0;
309  for (; it != dm_field->coarseningIS.end(); it++, ii++) {
310  if (ii < nb_elems) {
311  PetscBool flg;
312  CHKERR ISEqual(*it, is_vec[ii], &flg);
313  if (!flg) {
314  CHKERR ISDestroy(&*it);
315  CHKERR MatDestroy(&dm_field->kspOperators[ii]);
316  *it = is_vec[ii];
317  CHKERR PetscObjectReference((PetscObject)is_vec[ii]);
318  if (ii < nb_elems - 1) {
319  IS is = is_vec[ii];
320  if (dm_field->aO) {
321  CHKERR ISDuplicate(is_vec[ii], &is);
322  CHKERR ISCopy(is_vec[ii], is);
323  CHKERR AOApplicationToPetscIS(dm_field->aO, is);
324  }
325  Mat subA;
326 #if PETSC_VERSION_GE(3, 8, 0)
327  CHKERR MatCreateSubMatrix(A, is, is, MAT_INITIAL_MATRIX, &subA);
328 #else
329  CHKERR MatGetSubMatrix(A, is, is, MAT_INITIAL_MATRIX, &subA);
330 #endif
331  CHKERR PetscObjectReference((PetscObject)subA);
332  dm_field->kspOperators[ii] = subA;
333  CHKERR MatDestroy(&subA);
334  if (dm_field->aO) {
335  CHKERR ISDestroy(&is);
336  }
337  } else {
338  CHKERR PetscObjectReference((PetscObject)A);
339  dm_field->kspOperators[ii] = A;
340  }
341  nb_replaced++;
342  }
343  } else {
344  nb_no_changed++;
345  continue;
346  }
347  }
348  if (static_cast<int>(dm_field->coarseningIS.size()) < nb_elems) {
349  for (; ii < nb_elems - 1; ii++) {
350  Mat subA;
351  CHKERR DMMGViaApproxOrdersPushBackCoarseningIS(dm, is_vec[ii], A, &subA,
352  true, false);
353  CHKERR MatDestroy(&subA);
354  nb_added++;
355  }
356  CHKERR DMMGViaApproxOrdersPushBackCoarseningIS(dm, is_vec[ii], A, &A, false,
357  false);
358  nb_added++;
359  } else {
360  for (; ii < static_cast<int>(dm_field->coarseningIS.size()); ii++) {
362  nb_deleted++;
363  }
364  }
365  MPI_Comm comm;
366  CHKERR PetscObjectGetComm((PetscObject)dm, &comm);
367  if (verb > 0) {
368  PetscPrintf(comm,
369  "DMMGViaApproxOrders nb_no_changed = %d, nb_replaced = %d, "
370  "nb_added = %d, nb_deleted = %d, size = %d\n",
371  nb_no_changed, nb_replaced, nb_added, nb_deleted,
372  dm_field->coarseningIS.size());
373  }
374  PetscInfo(dm, "Replace IS to DMMGViaApproxOrders\n");
376 }
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.
double 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 194 of file PCMGSetUpViaApproxOrders.cpp.

194  {
195  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
197  GET_DM_FIELD(dm);
198  if (dm_field->aO) {
199  // std::cerr << dm_field->aO << std::endl;
200  CHKERR AODestroy(&dm_field->aO);
201  // std::cerr << "destroy ao when adding\n";
202  }
203  dm_field->aO = ao;
204  CHKERR PetscObjectReference((PetscObject)ao);
205  // std::cerr << "add ao\n";
207 }

◆ ksp_set_operators()

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

Definition at line 415 of file PCMGSetUpViaApproxOrders.cpp.

415  {
418 }

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

740  {
741 
742  PetscValidHeaderSpecific(pc, PC_CLASSID, 1);
744 
745  MPI_Comm comm;
746  CHKERR PetscObjectGetComm((PetscObject)pc, &comm);
747  if (verb > 0) {
748  PetscPrintf(comm, "Start PCMGSetUpViaApproxOrders\n");
749  }
750 
751  CHKERR ctx->getOptions();
752  CHKERR ctx->buildProlongationOperator(true, verb);
753 
754 #if PETSC_VERSION_GE(3, 8, 0)
755  CHKERR PCMGSetGalerkin(pc, PC_MG_GALERKIN_NONE);
756 #else
757  CHKERR PCMGSetGalerkin(pc, PETSC_FALSE);
758 #endif
759 
760  CHKERR PCMGSetLevels(pc, ctx->nbLevels, NULL);
761 
762  if (verb > 0) {
763  PetscPrintf(comm, "End PCMGSetUpViaApproxOrders\n");
764  }
765 
767 }
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.

◆ 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 }
constexpr double a

◆ 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 }

◆ 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 }
constexpr double omega