v0.9.1
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>

Inheritance diagram for PCMGSetUpViaApproxOrdersCtx:
[legend]
Collaboration diagram for PCMGSetUpViaApproxOrdersCtx:
[legend]

Public Member Functions

 PCMGSetUpViaApproxOrdersCtx (DM dm, Mat a, bool shell_sub_a)
 
virtual ~PCMGSetUpViaApproxOrdersCtx ()
 
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 198 of file PCMGSetUpViaApproxOrders.hpp.

Constructor & Destructor Documentation

◆ PCMGSetUpViaApproxOrdersCtx()

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

Definition at line 206 of file PCMGSetUpViaApproxOrders.hpp.

208  :
209  // mFieldPtr(mfield_ptr),
210  // problemName(problem_name),
211  dM(dm),
212  A(a),
213  nbLevels(2),
214  coarseOrder(2),
215  orderAtLastLevel(1000),
216  shellSubA(shell_sub_a),
217  verboseLevel(0) {
218  }
DM dM
Distributed mesh manager.
int nbLevels
number of multi-grid levels
int orderAtLastLevel
set maximal evaluated order
int coarseOrder
approximation order of coarse level

◆ ~PCMGSetUpViaApproxOrdersCtx()

virtual PCMGSetUpViaApproxOrdersCtx::~PCMGSetUpViaApproxOrdersCtx ( )
virtual

Definition at line 220 of file PCMGSetUpViaApproxOrders.hpp.

220  {
221  }

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

683  {
685  verb = verb > verboseLevel ? verb : verboseLevel;
686 
687  MPI_Comm comm;
688  CHKERR PetscObjectGetComm((PetscObject)dM, &comm);
689 
690  if (verb > QUIET) {
691  PetscPrintf(comm, "set MG levels %u\n", nbLevels);
692  }
693 
694  std::vector<IS> is_vec(nbLevels + 1);
695  std::vector<int> is_glob_size(nbLevels + 1), is_loc_size(nbLevels + 1);
696 
697  for (int kk = 0; kk < nbLevels; kk++) {
698 
699  // get indices up to up to give approximation order
700  CHKERR createIsAtLevel(kk, &is_vec[kk]);
701  CHKERR ISGetSize(is_vec[kk], &is_glob_size[kk]);
702  CHKERR ISGetLocalSize(is_vec[kk], &is_loc_size[kk]);
703 
704  if (verb > QUIET) {
705  PetscSynchronizedPrintf(comm,
706  "Nb. dofs at level [ %d ] global %u local %d\n",
707  kk, is_glob_size[kk], is_loc_size[kk]);
708  }
709 
710  // if no dofs on level kk finish here
711  if (is_glob_size[kk] == 0) {
712  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "no dofs at level");
713  }
714  }
715 
716  for (int kk = 0; kk != nbLevels; kk++) {
717  Mat subA;
718  if (kk == nbLevels - 1 && use_mat_a) {
719  subA = A;
721  false, false);
722  } else {
723  if (kk > 0) {
724  // Not coarse level
726  true, shellSubA);
727  } else {
728  // Coarse lave is compressed matrix allowing for factorization when
729  // needed
731  true, false);
732  }
733  if (subA) {
734  CHKERR MatDestroy(&subA);
735  }
736  }
737  }
738 
739  for (unsigned int kk = 0; kk < is_vec.size(); kk++) {
740  CHKERR destroyIsAtLevel(kk, &is_vec[kk]);
741  }
742 
743  if (verb > QUIET) {
744  PetscSynchronizedFlush(comm, PETSC_STDOUT);
745  }
746 
748 }
DM dM
Distributed mesh manager.
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:483
int nbLevels
number of multi-grid levels
#define CHKERR
Inline error check.
Definition: definitions.h:602
virtual MoFEMErrorCode createIsAtLevel(int kk, IS *is)
Set IS for levels.
virtual MoFEMErrorCode destroyIsAtLevel(int kk, IS *is)
Destroy IS if internally created.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ createIsAtLevel()

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

Set IS for levels.

Parameters
kklevel
ispointer to IS
Returns
error code

Reimplemented in SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx.

Definition at line 648 of file PCMGSetUpViaApproxOrders.cpp.

648  {
649  MoFEM::Interface *m_field_ptr;
650  MoFEM::ISManager *is_manager_ptr;
652  // if is last level, take all remaining orders dofs, if any left
653  CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
654  CHKERR m_field_ptr->getInterface(is_manager_ptr);
655  const Problem *problem_ptr;
656  CHKERR DMMoFEMGetProblemPtr(dM, &problem_ptr);
657  int order_at_next_level = kk + coarseOrder;
658  if (kk == nbLevels - 1) {
659  int first = problem_ptr->getNumeredDofsRows()
660  ->get<PetscLocalIdx_mi_tag>()
661  .find(0)
662  ->get()
663  ->getPetscGlobalDofIdx();
664  CHKERR ISCreateStride(PETSC_COMM_WORLD, problem_ptr->getNbLocalDofsRow(),
665  first, 1, is);
667  // order_at_next_level = orderAtLastLevel;
668  }
669  string problem_name = problem_ptr->getName();
670  CHKERR is_manager_ptr->isCreateProblemOrder(problem_name, ROW, 0,
671  order_at_next_level, is);
673 }
Deprecated interface functions.
DM dM
Distributed mesh manager.
Section manager is used to create indexes and sectionsFIXME: ISManager is not properly testsed by ato...
Definition: ISManager.hpp:36
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
DofIdx getNbLocalDofsRow() const
int nbLevels
number of multi-grid levels
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
keeps basic data about problemThis is low level structure with information about problem,...
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
#define CHKERR
Inline error check.
Definition: definitions.h:602
int coarseOrder
approximation order of coarse level
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition: DMMMoFEM.cpp:336
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMMoFEM.cpp:348
MoFEMErrorCode isCreateProblemOrder(const std::string &problem, RowColData rc, int min_order, int max_order, IS *is) const
create IS for given order range (collective)
Definition: ISManager.cpp:174
std::string getName() const
const boost::shared_ptr< NumeredDofEntity_multiIndex > & getNumeredDofsRows() const
get access to numeredDofsRows storing DOFs on rows

◆ destroyIsAtLevel()

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

Destroy IS if internally created.

Parameters
kklevel
ispointer to is
Returns
error code

Reimplemented in SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx.

Definition at line 675 of file PCMGSetUpViaApproxOrders.cpp.

675  {
677  CHKERR ISDestroy(is);
679 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ getOptions()

MoFEMErrorCode PCMGSetUpViaApproxOrdersCtx::getOptions ( )
virtual

get options from line command

Returns
error code

Reimplemented in SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx.

Definition at line 623 of file PCMGSetUpViaApproxOrders.cpp.

623  {
625  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
626  "MOFEM Multi-Grid (Orders) pre-conditioner", "none");
627  CHKERRG(ierr);
628 
629  CHKERR PetscOptionsInt("-mofem_mg_levels", "nb levels of multi-grid solver",
630  "", 2, &nbLevels, PETSC_NULL);
631  CHKERR PetscOptionsInt("-mofem_mg_coarse_order",
632  "approximation order of coarse level", "", 2,
633  &coarseOrder, PETSC_NULL);
634  CHKERR PetscOptionsInt("-mofem_mg_order_at_last_level", "order at last level",
635  "", 100, &orderAtLastLevel, PETSC_NULL);
636  CHKERR PetscOptionsInt("-mofem_mg_verbose", "nb levels of multi-grid solver",
637  "", 0, &verboseLevel, PETSC_NULL);
638  PetscBool shell_sub_a = shellSubA ? PETSC_TRUE : PETSC_FALSE;
639  CHKERR PetscOptionsBool("-mofem_mg_shell_a", "use shell matrix as sub matrix",
640  "", shell_sub_a, &shell_sub_a, NULL);
641  shellSubA = (shellSubA == PETSC_TRUE);
642 
643  ierr = PetscOptionsEnd();
644  CHKERRG(ierr);
646 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:507
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
int nbLevels
number of multi-grid levels
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
int orderAtLastLevel
set maximal evaluated order
#define CHKERR
Inline error check.
Definition: definitions.h:602
int coarseOrder
approximation order of coarse level

Member Data Documentation

◆ A

Mat PCMGSetUpViaApproxOrdersCtx::A

Matrix at fine level.

Definition at line 204 of file PCMGSetUpViaApproxOrders.hpp.

◆ coarseOrder

int PCMGSetUpViaApproxOrdersCtx::coarseOrder

approximation order of coarse level

Definition at line 224 of file PCMGSetUpViaApproxOrders.hpp.

◆ dM

DM PCMGSetUpViaApproxOrdersCtx::dM

Distributed mesh manager.

Definition at line 203 of file PCMGSetUpViaApproxOrders.hpp.

◆ nbLevels

int PCMGSetUpViaApproxOrdersCtx::nbLevels

number of multi-grid levels

Definition at line 223 of file PCMGSetUpViaApproxOrders.hpp.

◆ orderAtLastLevel

int PCMGSetUpViaApproxOrdersCtx::orderAtLastLevel

set maximal evaluated order

Definition at line 225 of file PCMGSetUpViaApproxOrders.hpp.

◆ shellSubA

bool PCMGSetUpViaApproxOrdersCtx::shellSubA

Definition at line 227 of file PCMGSetUpViaApproxOrders.hpp.

◆ verboseLevel

int PCMGSetUpViaApproxOrdersCtx::verboseLevel

Definition at line 228 of file PCMGSetUpViaApproxOrders.hpp.


The documentation for this struct was generated from the following files: