v0.14.0
Public Member Functions | Public Attributes | List of all members
PCMGSetUpViaApproxOrdersCtx Struct Reference

Set data structures of MG pre-conditioner via approximation orders. More...

#include <users_modules/basic_finite_elements/src/PCMGSetUpViaApproxOrders.hpp>

Collaboration diagram for PCMGSetUpViaApproxOrdersCtx:
[legend]

Public Member Functions

 PCMGSetUpViaApproxOrdersCtx (DM dm, Mat a, bool shell_sub_a)
 
virtual ~PCMGSetUpViaApproxOrdersCtx ()=default
 
virtual MoFEMErrorCode getOptions ()
 get options from line command More...
 
virtual MoFEMErrorCode createIsAtLevel (int kk, IS *is)
 Set IS for levels. More...
 
virtual MoFEMErrorCode destroyIsAtLevel (int kk, IS *is)
 Destroy IS if internally created. More...
 
virtual MoFEMErrorCode buildProlongationOperator (bool use_mat_a, int verb=0)
 Set up data structures for MG. More...
 

Public Attributes

DM dM
 Distributed mesh manager. More...
 
Mat A
 Matrix at fine level. More...
 
int nbLevels
 number of multi-grid levels More...
 
int coarseOrder
 approximation order of coarse level More...
 
int orderAtLastLevel
 set maximal evaluated order More...
 
bool shellSubA
 
int verboseLevel
 

Detailed Description

Set data structures of MG pre-conditioner via approximation orders.

Examples
elasticity.cpp.

Definition at line 190 of file PCMGSetUpViaApproxOrders.hpp.

Constructor & Destructor Documentation

◆ PCMGSetUpViaApproxOrdersCtx()

PCMGSetUpViaApproxOrdersCtx::PCMGSetUpViaApproxOrdersCtx ( DM  dm,
Mat  a,
bool  shell_sub_a 
)
inline

Definition at line 198 of file PCMGSetUpViaApproxOrders.hpp.

199  : dM(dm), A(a), nbLevels(2), coarseOrder(2), orderAtLastLevel(1000),
200  shellSubA(shell_sub_a), verboseLevel(0) {}

◆ ~PCMGSetUpViaApproxOrdersCtx()

virtual PCMGSetUpViaApproxOrdersCtx::~PCMGSetUpViaApproxOrdersCtx ( )
virtualdefault

Member Function Documentation

◆ buildProlongationOperator()

MoFEMErrorCode PCMGSetUpViaApproxOrdersCtx::buildProlongationOperator ( bool  use_mat_a,
int  verb = 0 
)
virtual

Set up data structures for MG.

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

Definition at line 659 of file PCMGSetUpViaApproxOrders.cpp.

660  {
662  verb = verb > verboseLevel ? verb : verboseLevel;
663 
664  MPI_Comm comm;
665  CHKERR PetscObjectGetComm((PetscObject)dM, &comm);
666 
667  if (verb > QUIET) {
668  PetscPrintf(comm, "set MG levels %u\n", nbLevels);
669  }
670 
671  std::vector<IS> is_vec(nbLevels + 1);
672  std::vector<int> is_glob_size(nbLevels + 1), is_loc_size(nbLevels + 1);
673 
674  for (int kk = 0; kk < nbLevels; kk++) {
675 
676  // get indices up to up to give approximation order
677  CHKERR createIsAtLevel(kk, &is_vec[kk]);
678  CHKERR ISGetSize(is_vec[kk], &is_glob_size[kk]);
679  CHKERR ISGetLocalSize(is_vec[kk], &is_loc_size[kk]);
680 
681  if (verb > QUIET) {
682  PetscSynchronizedPrintf(comm,
683  "Nb. dofs at level [ %d ] global %u local %d\n",
684  kk, is_glob_size[kk], is_loc_size[kk]);
685  }
686 
687  // if no dofs on level kk finish here
688  if (is_glob_size[kk] == 0) {
689  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "no dofs at level");
690  }
691  }
692 
693  for (int kk = 0; kk != nbLevels; kk++) {
694  Mat subA;
695  if (kk == nbLevels - 1 && use_mat_a) {
696  subA = A;
698  false, false);
699  } else {
700  if (kk > 0) {
701  // Not coarse level
703  true, shellSubA);
704  } else {
705  // Coarse lave is compressed matrix allowing for factorization when
706  // needed
708  true, false);
709  }
710  if (subA) {
711  CHKERR MatDestroy(&subA);
712  }
713  }
714  }
715 
716  for (unsigned int kk = 0; kk < is_vec.size(); kk++) {
717  CHKERR destroyIsAtLevel(kk, &is_vec[kk]);
718  }
719 
720  if (verb > QUIET) {
721  PetscSynchronizedFlush(comm, PETSC_STDOUT);
722  }
723 
725 }

◆ createIsAtLevel()

MoFEMErrorCode PCMGSetUpViaApproxOrdersCtx::createIsAtLevel ( int  kk,
IS *  is 
)
virtual

Set IS for levels.

Parameters
kklevel
ispointer to IS
Returns
error code

Definition at line 625 of file PCMGSetUpViaApproxOrders.cpp.

625  {
626  MoFEM::Interface *m_field_ptr;
627  MoFEM::ISManager *is_manager_ptr;
629  // if is last level, take all remaining orders dofs, if any left
630  CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
631  CHKERR m_field_ptr->getInterface(is_manager_ptr);
632  const Problem *problem_ptr;
633  CHKERR DMMoFEMGetProblemPtr(dM, &problem_ptr);
634  int order_at_next_level = kk + coarseOrder;
635  if (kk == nbLevels - 1) {
636  int first = problem_ptr->getNumeredRowDofsPtr()
637  ->get<PetscLocalIdx_mi_tag>()
638  .find(0)
639  ->get()
640  ->getPetscGlobalDofIdx();
641  CHKERR ISCreateStride(PETSC_COMM_WORLD, problem_ptr->getNbLocalDofsRow(),
642  first, 1, is);
644  // order_at_next_level = orderAtLastLevel;
645  }
646  string problem_name = problem_ptr->getName();
647  CHKERR is_manager_ptr->isCreateProblemOrder(problem_name, ROW, 0,
648  order_at_next_level, is);
650 }

◆ destroyIsAtLevel()

MoFEMErrorCode PCMGSetUpViaApproxOrdersCtx::destroyIsAtLevel ( int  kk,
IS *  is 
)
virtual

Destroy IS if internally created.

Parameters
kklevel
ispointer to is
Returns
error code

Definition at line 652 of file PCMGSetUpViaApproxOrders.cpp.

652  {
654  CHKERR ISDestroy(is);
656 }

◆ getOptions()

MoFEMErrorCode PCMGSetUpViaApproxOrdersCtx::getOptions ( )
virtual

get options from line command

Returns
error code

Definition at line 600 of file PCMGSetUpViaApproxOrders.cpp.

600  {
602  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
603  "MOFEM Multi-Grid (Orders) pre-conditioner", "none");
604  CHKERRG(ierr);
605 
606  CHKERR PetscOptionsInt("-mofem_mg_levels", "nb levels of multi-grid solver",
607  "", 2, &nbLevels, PETSC_NULL);
608  CHKERR PetscOptionsInt("-mofem_mg_coarse_order",
609  "approximation order of coarse level", "", 2,
610  &coarseOrder, PETSC_NULL);
611  CHKERR PetscOptionsInt("-mofem_mg_order_at_last_level", "order at last level",
612  "", 100, &orderAtLastLevel, PETSC_NULL);
613  CHKERR PetscOptionsInt("-mofem_mg_verbose", "nb levels of multi-grid solver",
614  "", 0, &verboseLevel, PETSC_NULL);
615  PetscBool shell_sub_a = shellSubA ? PETSC_TRUE : PETSC_FALSE;
616  CHKERR PetscOptionsBool("-mofem_mg_shell_a", "use shell matrix as sub matrix",
617  "", shell_sub_a, &shell_sub_a, NULL);
618  shellSubA = (shellSubA == PETSC_TRUE);
619 
620  ierr = PetscOptionsEnd();
621  CHKERRG(ierr);
623 }

Member Data Documentation

◆ A

Mat PCMGSetUpViaApproxOrdersCtx::A

Matrix at fine level.

Definition at line 196 of file PCMGSetUpViaApproxOrders.hpp.

◆ coarseOrder

int PCMGSetUpViaApproxOrdersCtx::coarseOrder

approximation order of coarse level

Definition at line 205 of file PCMGSetUpViaApproxOrders.hpp.

◆ dM

DM PCMGSetUpViaApproxOrdersCtx::dM

Distributed mesh manager.

Definition at line 195 of file PCMGSetUpViaApproxOrders.hpp.

◆ nbLevels

int PCMGSetUpViaApproxOrdersCtx::nbLevels

number of multi-grid levels

Definition at line 204 of file PCMGSetUpViaApproxOrders.hpp.

◆ orderAtLastLevel

int PCMGSetUpViaApproxOrdersCtx::orderAtLastLevel

set maximal evaluated order

Definition at line 206 of file PCMGSetUpViaApproxOrders.hpp.

◆ shellSubA

bool PCMGSetUpViaApproxOrdersCtx::shellSubA

Definition at line 208 of file PCMGSetUpViaApproxOrders.hpp.

◆ verboseLevel

int PCMGSetUpViaApproxOrdersCtx::verboseLevel

Definition at line 209 of file PCMGSetUpViaApproxOrders.hpp.


The documentation for this struct was generated from the following files:
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
PCMGSetUpViaApproxOrdersCtx::shellSubA
bool shellSubA
Definition: PCMGSetUpViaApproxOrders.hpp:208
MoFEM::PetscLocalIdx_mi_tag
Definition: TagMultiIndices.hpp:45
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:205
PCMGSetUpViaApproxOrdersCtx::nbLevels
int nbLevels
number of multi-grid levels
Definition: PCMGSetUpViaApproxOrders.hpp:204
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
ROW
@ ROW
Definition: definitions.h:123
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
MoFEM::ISManager
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:23
a
constexpr double a
Definition: approx_sphere.cpp:30
PCMGSetUpViaApproxOrdersCtx::dM
DM dM
Distributed mesh manager.
Definition: PCMGSetUpViaApproxOrders.hpp:195
MoFEM::Problem::getNumeredRowDofsPtr
auto & getNumeredRowDofsPtr() const
get access to numeredRowDofsPtr storing DOFs on rows
Definition: ProblemsMultiIndices.hpp:82
PCMGSetUpViaApproxOrdersCtx::createIsAtLevel
virtual MoFEMErrorCode createIsAtLevel(int kk, IS *is)
Set IS for levels.
Definition: PCMGSetUpViaApproxOrders.cpp:625
MoFEM::Problem::getNbLocalDofsRow
DofIdx getNbLocalDofsRow() const
Definition: ProblemsMultiIndices.hpp:378
PCMGSetUpViaApproxOrdersCtx::orderAtLastLevel
int orderAtLastLevel
set maximal evaluated order
Definition: PCMGSetUpViaApproxOrders.hpp:206
PCMGSetUpViaApproxOrdersCtx::destroyIsAtLevel
virtual MoFEMErrorCode destroyIsAtLevel(int kk, IS *is)
Destroy IS if internally created.
Definition: PCMGSetUpViaApproxOrders.cpp:652
PCMGSetUpViaApproxOrdersCtx::verboseLevel
int verboseLevel
Definition: PCMGSetUpViaApproxOrders.hpp:209
MoFEM::DMoFEMGetInterfacePtr
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMoFEM.cpp:414
PCMGSetUpViaApproxOrdersCtx::coarseOrder
int coarseOrder
approximation order of coarse level
Definition: PCMGSetUpViaApproxOrders.hpp:205
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::DMMoFEMGetProblemPtr
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMoFEM.cpp:426
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
QUIET
@ QUIET
Definition: definitions.h:208
MoFEM::Problem
keeps basic data about problem
Definition: ProblemsMultiIndices.hpp:54
MoFEM::ISManager::isCreateProblemOrder
MoFEMErrorCode isCreateProblemOrder(const std::string problem_name, RowColData rc, int min_order, int max_order, IS *is) const
create IS for given order range (collective)
Definition: ISManager.cpp:204
MoFEM::Problem::getName
auto getName() const
Definition: ProblemsMultiIndices.hpp:372
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
CHKERRG
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
PCMGSetUpViaApproxOrdersCtx::A
Mat A
Matrix at fine level.
Definition: PCMGSetUpViaApproxOrders.hpp:196
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346