v0.15.0
Loading...
Searching...
No Matches
MoFEM::PCMGSetUpViaApproxOrdersCtx Struct Reference
Collaboration diagram for MoFEM::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
 
virtual MoFEMErrorCode createIsAtLevel (int kk, IS *is)
 Set IS for levels.
 
virtual MoFEMErrorCode destroyIsAtLevel (int kk, IS *is)
 Destroy IS if internally created.
 
virtual MoFEMErrorCode buildProlongationOperator (bool use_mat_a, int verb=0)
 Set up data structures for MG.
 

Public Attributes

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

Detailed Description

Definition at line 92 of file PCMGSetUpViaApproxOrders.cpp.

Constructor & Destructor Documentation

◆ PCMGSetUpViaApproxOrdersCtx()

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

◆ ~PCMGSetUpViaApproxOrdersCtx()

virtual MoFEM::PCMGSetUpViaApproxOrdersCtx::~PCMGSetUpViaApproxOrdersCtx ( )
virtualdefault

Member Function Documentation

◆ buildProlongationOperator()

MoFEMErrorCode MoFEM::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 590 of file PCMGSetUpViaApproxOrders.cpp.

591 {
593 verb = verb > verboseLevel ? verb : verboseLevel;
594
595 MPI_Comm comm;
596 CHKERR PetscObjectGetComm((PetscObject)dM, &comm);
597
598 MOFEM_LOG_C("PCMGViaApproximationOrders", Sev::inform, "set MG levels %u\n",
599 nbLevels);
600
601 MOFEM_LOG_CHANNEL("SYNC");
602 MOFEM_LOG_TAG("SYNC", "PCMGViaApproximationOrders")
603
604 std::vector<IS> is_vec(nbLevels + 1);
605 std::vector<int> is_glob_size(nbLevels + 1), is_loc_size(nbLevels + 1);
606
607 for (int kk = 0; kk < nbLevels; kk++) {
608
609 // get indices up to up to give approximation order
610 CHKERR createIsAtLevel(kk, &is_vec[kk]);
611 CHKERR ISGetSize(is_vec[kk], &is_glob_size[kk]);
612 CHKERR ISGetLocalSize(is_vec[kk], &is_loc_size[kk]);
613
614 MOFEM_LOG_C("SYNC", Sev::inform,
615 "Nb. dofs at level [ %d ] global %u local %d\n", kk,
616 is_glob_size[kk], is_loc_size[kk]);
617
618 // if no dofs on level kk finish here
619 if (is_glob_size[kk] == 0) {
620 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "no dofs at level");
621 }
622 }
623
624 for (int kk = 0; kk != nbLevels; kk++) {
625 if (kk == nbLevels - 1 && use_mat_a) {
627 false);
628 } else {
629 if (kk > 0) {
630 // Not coarse level
632 shellSubA);
633 } else {
634 // Coarse lave is compressed matrix allowing for factorization when
635 // needed
637 false);
638 }
639 }
640 }
641
642 for (unsigned int kk = 0; kk < is_vec.size(); kk++) {
643 CHKERR destroyIsAtLevel(kk, &is_vec[kk]);
644 }
645
646 MOFEM_LOG_SEVERITY_SYNC(comm, Sev::inform);
647 MOFEM_LOG_CHANNEL("SYNC");
648
650}
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
#define MOFEM_LOG_C(channel, severity, format,...)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
MoFEMErrorCode DMMGViaApproxOrdersPushBackCoarseningIS(DM dm, IS is, Mat A, bool create_sub_matrix, bool shell_sub_a)
Push back coarsening level to MG via approximation orders.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
virtual MoFEMErrorCode destroyIsAtLevel(int kk, IS *is)
Destroy IS if internally created.
virtual MoFEMErrorCode createIsAtLevel(int kk, IS *is)
Set IS for levels.

◆ createIsAtLevel()

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

Set IS for levels.

Parameters
kklevel
ispointer to IS
Returns
error code

Definition at line 556 of file PCMGSetUpViaApproxOrders.cpp.

556 {
557 MoFEM::Interface *m_field_ptr;
558 MoFEM::ISManager *is_manager_ptr;
560 // if is last level, take all remaining orders dofs, if any left
561 CHKERR DMoFEMGetInterfacePtr(dM, &m_field_ptr);
562 CHKERR m_field_ptr->getInterface(is_manager_ptr);
563 const Problem *problem_ptr;
564 CHKERR DMMoFEMGetProblemPtr(dM, &problem_ptr);
565 int order_at_next_level = kk + coarseOrder;
566 if (kk == nbLevels - 1) {
567 int first = problem_ptr->getNumeredRowDofsPtr()
568 ->get<PetscGlobalIdx_mi_tag>()
569 .lower_bound(0)
570 ->get()
571 ->getPetscGlobalDofIdx();
572 CHKERR ISCreateStride(PETSC_COMM_WORLD, problem_ptr->getNbLocalDofsRow(),
573 first, 1, is);
575 // order_at_next_level = orderAtLastLevel;
576 }
577 string problem_name = problem_ptr->getName();
578 CHKERR is_manager_ptr->isCreateProblemOrder(problem_name, ROW, 0,
579 order_at_next_level, is);
581}
@ ROW
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition DMMoFEM.cpp:422
PetscErrorCode DMoFEMGetInterfacePtr(DM dm, MoFEM::Interface **m_field_ptr)
Get pointer to MoFEM::Interface.
Definition DMMoFEM.cpp:410
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)
Deprecated interface functions.
Section manager is used to create indexes and sections.
Definition ISManager.hpp:23
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.

◆ destroyIsAtLevel()

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

Destroy IS if internally created.

Parameters
kklevel
ispointer to is
Returns
error code

Definition at line 583 of file PCMGSetUpViaApproxOrders.cpp.

583 {
585 CHKERR ISDestroy(is);
587}

◆ getOptions()

MoFEMErrorCode MoFEM::PCMGSetUpViaApproxOrdersCtx::getOptions ( )
virtual

get options from line command

Returns
error code

Definition at line 533 of file PCMGSetUpViaApproxOrders.cpp.

533 {
535 PetscOptionsBegin(PETSC_COMM_WORLD, "",
536 "MOFEM Multi-Grid (Orders) pre-conditioner", "none");
537
538 CHKERR PetscOptionsInt("-mofem_mg_levels", "nb levels of multi-grid solver",
539 "", 2, &nbLevels, PETSC_NULLPTR);
540 CHKERR PetscOptionsInt("-mofem_mg_coarse_order",
541 "approximation order of coarse level", "", 2,
542 &coarseOrder, PETSC_NULLPTR);
543 CHKERR PetscOptionsInt("-mofem_mg_order_at_last_level", "order at last level",
544 "", 100, &orderAtLastLevel, PETSC_NULLPTR);
545 CHKERR PetscOptionsInt("-mofem_mg_verbose", "nb levels of multi-grid solver",
546 "", 0, &verboseLevel, PETSC_NULLPTR);
547 PetscBool shell_sub_a = shellSubA ? PETSC_TRUE : PETSC_FALSE;
548 CHKERR PetscOptionsBool("-mofem_mg_shell_a", "use shell matrix as sub matrix",
549 "", shell_sub_a, &shell_sub_a, NULL);
550 shellSubA = (shellSubA == PETSC_TRUE);
551
552 PetscOptionsEnd();
554}
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...

Member Data Documentation

◆ A

Mat MoFEM::PCMGSetUpViaApproxOrdersCtx::A

Matrix at fine level.

Definition at line 133 of file PCMGSetUpViaApproxOrders.cpp.

◆ coarseOrder

int MoFEM::PCMGSetUpViaApproxOrdersCtx::coarseOrder

approximation order of coarse level

Definition at line 136 of file PCMGSetUpViaApproxOrders.cpp.

◆ dM

DM MoFEM::PCMGSetUpViaApproxOrdersCtx::dM

Distributed mesh manager.

Definition at line 132 of file PCMGSetUpViaApproxOrders.cpp.

◆ nbLevels

int MoFEM::PCMGSetUpViaApproxOrdersCtx::nbLevels

number of multi-grid levels

Definition at line 135 of file PCMGSetUpViaApproxOrders.cpp.

◆ orderAtLastLevel

int MoFEM::PCMGSetUpViaApproxOrdersCtx::orderAtLastLevel

set maximal evaluated order

Definition at line 137 of file PCMGSetUpViaApproxOrders.cpp.

◆ shellSubA

bool MoFEM::PCMGSetUpViaApproxOrdersCtx::shellSubA

Definition at line 139 of file PCMGSetUpViaApproxOrders.cpp.

◆ verboseLevel

int MoFEM::PCMGSetUpViaApproxOrdersCtx::verboseLevel

Definition at line 140 of file PCMGSetUpViaApproxOrders.cpp.


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