v0.13.2
Loading...
Searching...
No Matches
Public Member Functions | Private Member Functions | Private Attributes | List of all members
Poisson3DHomogeneous Struct Reference
Collaboration diagram for Poisson3DHomogeneous:
[legend]

Public Member Functions

 Poisson3DHomogeneous (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProgram ()
 

Private Member Functions

MoFEMErrorCode readMesh ()
 
MoFEMErrorCode setupProblem ()
 
MoFEMErrorCode boundaryCondition ()
 
MoFEMErrorCode assembleSystem ()
 
MoFEMErrorCode setIntegrationRules ()
 
MoFEMErrorCode solveSystem ()
 
MoFEMErrorCode outputResults ()
 

Private Attributes

MoFEM::InterfacemField
 
SimplesimpleInterface
 
std::string domainField
 
int oRder
 

Detailed Description

Examples
poisson_3d_homogeneous.cpp.

Definition at line 21 of file poisson_3d_homogeneous.cpp.

Constructor & Destructor Documentation

◆ Poisson3DHomogeneous()

Poisson3DHomogeneous::Poisson3DHomogeneous ( MoFEM::Interface m_field)
Examples
poisson_3d_homogeneous.cpp.

Definition at line 47 of file poisson_3d_homogeneous.cpp.

48 : domainField("U"), mField(m_field) {}

Member Function Documentation

◆ assembleSystem()

MoFEMErrorCode Poisson3DHomogeneous::assembleSystem ( )
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 94 of file poisson_3d_homogeneous.cpp.

94 {
96
97 auto pipeline_mng = mField.getInterface<PipelineManager>();
98
99 { // Push operators to the Pipeline that is responsible for calculating LHS
100 CHKERR AddHOOps<3, 3, 3>::add(pipeline_mng->getOpDomainLhsPipeline(), {H1});
101 pipeline_mng->getOpDomainLhsPipeline().push_back(
103 }
104
105 { // Push operators to the Pipeline that is responsible for calculating LHS
106
107 constexpr int space_dim = 3;
108
109 auto set_values_to_bc_dofs = [&](auto &fe) {
110 auto get_bc_hook = [&]() {
112 return hook;
113 };
114 fe->preProcessHook = get_bc_hook();
115 };
116
117 auto calculate_residual_from_set_values_on_bc = [&](auto &pipeline) {
118 using DomainEle =
123
124 auto grad_u_vals_ptr = boost::make_shared<MatrixDouble>();
125 pipeline_mng->getOpDomainRhsPipeline().push_back(
127 grad_u_vals_ptr));
128 pipeline_mng->getOpDomainRhsPipeline().push_back(
129 new OpInternal(domainField, grad_u_vals_ptr,
130 [](double, double, double) constexpr { return -1; }));
131 };
132
134 pipeline_mng->getOpDomainRhsPipeline(), {H1});
135 set_values_to_bc_dofs(pipeline_mng->getDomainRhsFE());
136 calculate_residual_from_set_values_on_bc(
137 pipeline_mng->getOpDomainRhsPipeline());
138 pipeline_mng->getOpDomainRhsPipeline().push_back(
140 }
141
143}
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
Add operators pushing bases from local to physical configuration.
Class (Function) to enforce essential constrains.
Definition: Essential.hpp:39
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
PipelineManager interface.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

◆ boundaryCondition()

MoFEMErrorCode Poisson3DHomogeneous::boundaryCondition ( )
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 75 of file poisson_3d_homogeneous.cpp.

75 {
77
78 auto bc_mng = mField.getInterface<BcManager>();
79
80 // Remove BCs from blockset name "BOUNDARY_CONDITION" or SETU, note that you
81 // can use regular expression to put list of blocksets;
82 CHKERR bc_mng->removeBlockDOFsOnEntities<BcScalarMeshsetType<BLOCKSET>>(
83 simpleInterface->getProblemName(), "(BOUNDARY_CONDITION|SETU)",
84 domainField, true);
85
86 // Remove BCs from cubit TEMPERATURESET, i.e. set by cubit, and meshsets named
87 // FIX_SCALAR (default name to name boundary conditions for scalar fields)
88 CHKERR bc_mng->removeBlockDOFsOnEntities<TemperatureCubitBcData>(
89 simpleInterface->getProblemName(), domainField, true, false, true);
90
92}
Simple interface for fast problem set-up.
Definition: BcManager.hpp:24
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:348
Definition of the temperature bc data structure.
Definition: BCData.hpp:301

◆ outputResults()

MoFEMErrorCode Poisson3DHomogeneous::outputResults ( )
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 183 of file poisson_3d_homogeneous.cpp.

183 {
185
186 auto pipeline_mng = mField.getInterface<PipelineManager>();
187 pipeline_mng->getDomainLhsFE().reset();
188
189 auto post_proc_fe = boost::make_shared<PostProcVolEle>(mField);
190
191 auto u_ptr = boost::make_shared<VectorDouble>();
192 post_proc_fe->getOpPtrVector().push_back(
194
196
197 post_proc_fe->getOpPtrVector().push_back(
198
199 new OpPPMap(
200
201 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
202
203 {{domainField, u_ptr}},
204
205 {},
206
207 {},
208
209 {})
210
211 );
212
213 pipeline_mng->getDomainRhsFE() = post_proc_fe;
214 CHKERR pipeline_mng->loopFiniteElements();
215 CHKERR post_proc_fe->writeFile("out_result.h5m");
216
218}
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
Get value at integration points for scalar field.
Post post-proc data at points from hash maps.
boost::shared_ptr< FEMethod > & getDomainLhsFE()

◆ readMesh()

MoFEMErrorCode Poisson3DHomogeneous::readMesh ( )
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 50 of file poisson_3d_homogeneous.cpp.

50 {
52
56
58}
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:180
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:194

◆ runProgram()

MoFEMErrorCode Poisson3DHomogeneous::runProgram ( )

◆ setIntegrationRules()

MoFEMErrorCode Poisson3DHomogeneous::setIntegrationRules ( )
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 145 of file poisson_3d_homogeneous.cpp.

145 {
147
148 auto rule_lhs = [](int, int, int p) -> int { return 2 * (p - 1); };
149 auto rule_rhs = [](int, int, int p) -> int { return p; };
150
151 auto pipeline_mng = mField.getInterface<PipelineManager>();
152 CHKERR pipeline_mng->setDomainLhsIntegrationRule(rule_lhs);
153 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule_rhs);
154
156}
static Index< 'p', 3 > p

◆ setupProblem()

MoFEMErrorCode Poisson3DHomogeneous::setupProblem ( )
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 60 of file poisson_3d_homogeneous.cpp.

60 {
62
65
66 int oRder = 3;
67 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &oRder, PETSC_NULL);
69
71
73}
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ H1
continuous field
Definition: definitions.h:85
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
MoFEMErrorCode addDomainField(const std::string &name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_ZERO, int verb=-1)
Add field on domain.
Definition: Simple.cpp:264
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:476
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:614

◆ solveSystem()

MoFEMErrorCode Poisson3DHomogeneous::solveSystem ( )
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 158 of file poisson_3d_homogeneous.cpp.

158 {
160
161 auto pipeline_mng = mField.getInterface<PipelineManager>();
162
163 auto ksp_solver = pipeline_mng->createKSP();
164 CHKERR KSPSetFromOptions(ksp_solver);
165 CHKERR KSPSetUp(ksp_solver);
166
167 // Create RHS and solution vectors
168 auto dm = simpleInterface->getDM();
169 auto F = smartCreateDMVector(dm);
170 auto D = smartVectorDuplicate(F);
171
172 // Solve the system
173 CHKERR KSPSolve(ksp_solver, F, D);
174
175 // Scatter result data on the mesh
176 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
177 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
178 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
179
181}
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMoFEM.cpp:511
auto smartCreateDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:995
SmartPetscObj< KSP > createKSP(SmartPetscObj< DM > dm=nullptr)
Create KSP (linear) solver.
double D
SmartPetscObj< Vec > smartVectorDuplicate(SmartPetscObj< Vec > &vec)
Create duplicate vector of smart vector.
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition: Simple.cpp:671

Member Data Documentation

◆ domainField

std::string Poisson3DHomogeneous::domainField
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 43 of file poisson_3d_homogeneous.cpp.

◆ mField

MoFEM::Interface& Poisson3DHomogeneous::mField
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 39 of file poisson_3d_homogeneous.cpp.

◆ oRder

int Poisson3DHomogeneous::oRder
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 44 of file poisson_3d_homogeneous.cpp.

◆ simpleInterface

Simple* Poisson3DHomogeneous::simpleInterface
private
Examples
poisson_3d_homogeneous.cpp.

Definition at line 40 of file poisson_3d_homogeneous.cpp.


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