v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | List of all members
SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx Struct Reference

#include <users_modules/solid_shell_prism_element/src/RunAdaptivity.hpp>

Inheritance diagram for SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx:
[legend]
Collaboration diagram for SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx:
[legend]

Public Member Functions

 ShellPCMGSetUpViaApproxOrdersCtx (DM dm, Mat a, vector< IS > &mg_levels)
 
PetscErrorCode createIsAtLevel (int kk, IS *is)
 Set IS for levels. More...
 
PetscErrorCode destroyIsAtLevel (int kk, IS *is)
 Destroy IS if internally created. More...
 
PetscErrorCode getOptions ()
 get options from line command More...
 
PetscErrorCode buildProlongationOperator (PC pc, int verb=0)
 
- Public Member Functions inherited from PCMGSetUpViaApproxOrdersCtx
 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

vector< IS > & mgLevels
 
AO aO
 
- Public Attributes inherited from PCMGSetUpViaApproxOrdersCtx
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

Multi grid pre-conditioner

Definition at line 60 of file RunAdaptivity.hpp.

Constructor & Destructor Documentation

◆ ShellPCMGSetUpViaApproxOrdersCtx()

SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx::ShellPCMGSetUpViaApproxOrdersCtx ( DM  dm,
Mat  a,
vector< IS > &  mg_levels 
)
inline

Definition at line 63 of file RunAdaptivity.hpp.

64 : PCMGSetUpViaApproxOrdersCtx(dm, a, false), mgLevels(mg_levels) {}
Set data structures of MG pre-conditioner via approximation orders.

Member Function Documentation

◆ buildProlongationOperator()

PetscErrorCode SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx::buildProlongationOperator ( PC  pc,
int  verb = 0 
)

Definition at line 89 of file RunAdaptivity.cpp.

90 {
91 PetscFunctionBegin;
92 verb = verb > verboseLevel ? verb : verboseLevel;
93
94 MPI_Comm comm;
95 CHKERR PetscObjectGetComm((PetscObject)dM, &comm);
96
97 if (verb > 0) {
98 PetscPrintf(comm, "set MG levels %u\n", nbLevels);
99 }
100
101 vector<IS> is_vec(nbLevels + 1);
102 vector<int> is_glob_size(nbLevels + 1), is_loc_size(nbLevels + 1);
103
104 for (int kk = 0; kk < nbLevels; kk++) {
105
106 // get indices up to up to give approximation order
107 CHKERR createIsAtLevel(kk, &is_vec[kk]);
108 CHKERR ISGetSize(is_vec[kk], &is_glob_size[kk]);
109 CHKERR ISGetLocalSize(is_vec[kk], &is_loc_size[kk]);
110
111 if (verb > 0) {
112 PetscSynchronizedPrintf(comm,
113 "Nb. dofs at level [ %d ] global %u local %d\n",
114 kk, is_glob_size[kk], is_loc_size[kk]);
115 }
116
117 // if no dofs on level kk finish here
118 if (is_glob_size[kk] == 0) {
119 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "no dofs at level");
120 }
121 }
122
123 // #if PETSC_VERSION_GE(3, 8, 0)
124 // #error "This is not working downgrade petsc-3.7.* or optimally petsc-3.6.*"
125 // #else
126 // CHKERR PCMGSetGalerkin(pc, PETSC_FALSE);
127 // CHKERRQ(ierr);
128 // #endif
129
130#if PETSC_VERSION_GE(3, 8, 0)
131 CHKERR PCMGSetGalerkin(pc, PC_MG_GALERKIN_PMAT);
132#warning "This is not working downgrade petsc-3.7.* or optimally petsc-3.6.*"
133#else
134 CHKERR PCMGSetGalerkin(pc, PETSC_TRUE);
135#endif
136
137 // CHKERR PCMGSetLevels(pc,nbLevels,NULL); CHKERRQ(ierr);
139
141
142 for (unsigned int kk = 0; kk < is_vec.size(); kk++) {
143 CHKERR destroyIsAtLevel(kk, &is_vec[kk]);
144 }
145
146 if (verb > 0) {
147 PetscSynchronizedFlush(comm, PETSC_STDOUT);
148 }
149
150 PetscFunctionReturn(0);
151}
MoFEMErrorCode DMMGViaApproxOrdersReplaceCoarseningIS(DM dm, IS *is_vec, int nb_elems, Mat A, int verb)
Replace coarsening IS in DM via approximation orders.
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#define CHKERR
Inline error check.
Definition: definitions.h:535
PetscErrorCode pc_mg_set_last_level(PC pc, PetscInt levels, MPI_Comm *comms)
Definition: HackMG.cpp:73
DM dM
Distributed mesh manager.
int nbLevels
number of multi-grid levels
PetscErrorCode createIsAtLevel(int kk, IS *is)
Set IS for levels.
PetscErrorCode destroyIsAtLevel(int kk, IS *is)
Destroy IS if internally created.

◆ createIsAtLevel()

PetscErrorCode SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx::createIsAtLevel ( int  kk,
IS *  is 
)
virtual

Set IS for levels.

Parameters
kklevel
ispointer to IS
Returns
error code

Reimplemented from PCMGSetUpViaApproxOrdersCtx.

Definition at line 59 of file RunAdaptivity.cpp.

60 {
61 PetscFunctionBegin;
62 *is = mgLevels[kk];
63 PetscFunctionReturn(0);
64}

◆ destroyIsAtLevel()

PetscErrorCode SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx::destroyIsAtLevel ( int  kk,
IS *  is 
)
virtual

Destroy IS if internally created.

Parameters
kklevel
ispointer to is
Returns
error code

Reimplemented from PCMGSetUpViaApproxOrdersCtx.

Definition at line 67 of file RunAdaptivity.cpp.

68 {
69 PetscFunctionBegin;
70 PetscFunctionReturn(0);
71}

◆ getOptions()

PetscErrorCode SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx::getOptions ( )
virtual

get options from line command

Returns
error code

Reimplemented from PCMGSetUpViaApproxOrdersCtx.

Definition at line 73 of file RunAdaptivity.cpp.

73 {
74 PetscFunctionBegin;
75 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "",
76 "MOFEM Shell Multi-Grid (Orders) pre-conditioner",
77 "none");
78
79 CHKERR PetscOptionsInt("-mofem_mg_verbose", "nb levels of multi-grid solver",
80 "", 0, &verboseLevel, PETSC_NULL);
81
82 ierr = PetscOptionsEnd();
83 CHKERRQ(ierr);
84
85 PetscFunctionReturn(0);
86}
static PetscErrorCode ierr
Definition: HackMG.cpp:17

Member Data Documentation

◆ aO

AO SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx::aO

Definition at line 62 of file RunAdaptivity.hpp.

◆ mgLevels

vector<IS>& SolidShellModule::RunAdaptivity::ShellPCMGSetUpViaApproxOrdersCtx::mgLevels

Definition at line 61 of file RunAdaptivity.hpp.


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