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

Public Member Functions

 Poisson2DHomogeneous (MoFEM::Interface &m_field)
 
MoFEMErrorCode runProgram ()
 [Check]
 

Private Types

enum  NORMS { NORM = 0 , LAST_NORM }
 

Private Member Functions

MoFEMErrorCode readMesh ()
 [Read mesh]
 
MoFEMErrorCode setupProblem ()
 [Read mesh]
 
MoFEMErrorCode boundaryCondition ()
 [Setup problem]
 
MoFEMErrorCode assembleSystem ()
 [Boundary condition]
 
MoFEMErrorCode setIntegrationRules ()
 [Assemble system]
 
MoFEMErrorCode solveSystem ()
 [Set integration rules]
 
MoFEMErrorCode outputResults ()
 [Solve system]
 
MoFEMErrorCode checkResults ()
 [Output results]
 

Static Private Member Functions

static double sourceTermFunction (const double x, const double y, const double z)
 

Private Attributes

Simple * simpleInterface
 
MoFEM::Interface & mField
 
int oRder
 
SmartPetscObj< Vec > petscVec
 
int atom_test = 0
 

Detailed Description

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 27 of file poisson_2d_homogeneous.cpp.

Member Enumeration Documentation

◆ NORMS

Constructor & Destructor Documentation

◆ Poisson2DHomogeneous()

Poisson2DHomogeneous::Poisson2DHomogeneous ( MoFEM::Interface &  m_field)

Member Function Documentation

◆ assembleSystem()

MoFEMErrorCode Poisson2DHomogeneous::assembleSystem ( )
private

[Boundary condition]

[Assemble system]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 149 of file poisson_2d_homogeneous.cpp.

149 {
151
152 auto pipeline_mng = mField.getInterface<PipelineManager>();
153
154 { // Push operators to the Pipeline that is responsible for calculating LHS
156 pipeline_mng->getOpDomainLhsPipeline(), {H1});
157 pipeline_mng->getOpDomainLhsPipeline().push_back(
159 }
160
161 { // Push operators to the Pipeline that is responsible for calculating RHS
162
163 auto set_values_to_bc_dofs = [&](auto &fe) {
164 auto get_bc_hook = [&]() {
166 return hook;
167 };
168 fe->preProcessHook = get_bc_hook();
169 };
170
171 // you can skip that if boundary condition is prescribing zero
172 auto calculate_residual_from_set_values_on_bc = [&](auto &pipeline) {
173 using DomainEle =
175 using DomainEleOp = DomainEle::UserDataOperator;
178
179 auto grad_u_vals_ptr = boost::make_shared<MatrixDouble>();
180 pipeline_mng->getOpDomainRhsPipeline().push_back(
182 grad_u_vals_ptr));
183 pipeline_mng->getOpDomainRhsPipeline().push_back(
184 new OpInternal(field_name, grad_u_vals_ptr,
185 [](double, double, double) constexpr { return -1; }));
186 };
187
189 pipeline_mng->getOpDomainRhsPipeline(), {H1});
190 set_values_to_bc_dofs(pipeline_mng->getDomainRhsFE());
191 calculate_residual_from_set_values_on_bc(
192 pipeline_mng->getOpDomainRhsPipeline());
193 pipeline_mng->getOpDomainRhsPipeline().push_back(
195 }
196
198}
TestHook hook
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
@ PETSC
Standard PETSc assembly.
constexpr auto field_name
constexpr int SPACE_DIM
Add operators pushing bases from local to physical configuration.
Class (Function) to enforce essential constrains.
Definition Essential.hpp:25
Get field gradients at integration pts for scalar field rank 0, i.e. vector field.
Template struct for dimension-specific finite element types.
PipelineManager interface.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
static double sourceTermFunction(const double x, const double y, const double z)

◆ boundaryCondition()

MoFEMErrorCode Poisson2DHomogeneous::boundaryCondition ( )
private

[Setup problem]

[Boundary condition]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 133 of file poisson_2d_homogeneous.cpp.

133 {
135
136 auto bc_mng = mField.getInterface<BcManager>();
137
138 // Remove BCs from blockset name "BOUNDARY_CONDITION" or SETU, note that you
139 // can use regular expression to put list of blocksets;
141 simpleInterface->getProblemName(), "(BOUNDARY_CONDITION|SETU)",
142 std::string(field_name), true);
143
145}
MoFEMErrorCode removeBlockDOFsOnEntities(const std::string problem_name, const std::string block_name, const std::string field_name, int lo, int hi, bool get_low_dim_ents=true, bool is_distributed_mesh=true)
Remove DOFs from problem based on block entities.
Definition BcManager.cpp:72
Boundary condition manager for finite element problem setup.
Template specialization for scalar field boundary conditions.
const std::string getProblemName() const
Get the Problem Name.
Definition Simple.hpp:450

◆ checkResults()

MoFEMErrorCode Poisson2DHomogeneous::checkResults ( )
private

[Output results]

[Check]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 291 of file poisson_2d_homogeneous.cpp.

291 {
293
294 auto check_result_fe_ptr = boost::make_shared<DomainEle>(mField);
295 auto petscVec =
297 (mField.get_comm_rank() == 0) ? LAST_NORM : 0, LAST_NORM);
298
300 check_result_fe_ptr->getOpPtrVector(), {H1})),
301 "Apply transform");
302
303 check_result_fe_ptr->getRuleHook = [](int, int, int p) { return p; };
304 auto analyticalFunction = [&](double x, double y, double z) {
305 return sin(M_PI * x) * sin(M_PI * y);
306 };
307
308 auto u_ptr = boost::make_shared<VectorDouble>();
309
310 check_result_fe_ptr->getOpPtrVector().push_back(
312 auto mValFuncPtr = boost::make_shared<VectorDouble>();
313 check_result_fe_ptr->getOpPtrVector().push_back(
314 new OpGetTensor0fromFunc(mValFuncPtr, analyticalFunction));
315 check_result_fe_ptr->getOpPtrVector().push_back(
316 new OpCalcNormL2Tensor0(u_ptr, petscVec, NORM, mValFuncPtr));
317 CHKERR VecZeroEntries(petscVec);
320 check_result_fe_ptr);
321 CHKERR VecAssemblyBegin(petscVec);
322 CHKERR VecAssemblyEnd(petscVec);
323 MOFEM_LOG_CHANNEL("SELF"); // Clear channel from old tags
324 // print norm in general log
325 if (mField.get_comm_rank() == 0) {
326 const double *norms;
327 CHKERR VecGetArrayRead(petscVec, &norms);
328 MOFEM_TAG_AND_LOG("SELF", Sev::inform, "Errors")
329 << "NORM: " << std::sqrt(norms[NORM]);
330 CHKERR VecRestoreArrayRead(petscVec, &norms);
331 }
332 // compare norm for ctest
333 if (atom_test && !mField.get_comm_rank()) {
334 const double *t_ptr;
335 CHKERR VecGetArrayRead(petscVec, &t_ptr);
336 double ref_norm = 2.2e-04;
337 double cal_norm;
338 switch (atom_test) {
339 case 1: // 2D
340 cal_norm = sqrt(t_ptr[0]);
341 break;
342 default:
343 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
344 "atom test %d does not exist", atom_test);
345 }
346 if (cal_norm > ref_norm) {
347 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
348 "atom test %d failed! Calculated Norm %3.16e is greater than "
349 "reference Norm %3.16e",
350 atom_test, cal_norm, ref_norm);
351 }
352 CHKERR VecRestoreArrayRead(petscVec, &t_ptr);
353 }
355}
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
@ MOFEM_ATOM_TEST_INVALID
Definition definitions.h:40
@ MOFEM_INVALID_DATA
Definition definitions.h:36
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
Definition DMMoFEM.cpp:605
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
auto createVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
Get norm of input VectorDouble for Tensor0.
Specialization for double precision scalar field values calculation.
Get values from scalar function at integration points and save them to VectorDouble for Tensor0.
MoFEMErrorCode getDM(DM *dm)
Get DM.
Definition Simple.cpp:799
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition Simple.hpp:429
SmartPetscObj< Vec > petscVec

◆ outputResults()

MoFEMErrorCode Poisson2DHomogeneous::outputResults ( )
private

[Solve system]

[Output results]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 246 of file poisson_2d_homogeneous.cpp.

246 {
248
249 auto pipeline_mng = mField.getInterface<PipelineManager>();
250
251 auto post_proc_fe = boost::make_shared<PostProcFaceEle>(mField);
253 post_proc_fe->getOpPtrVector(), {H1});
254
255 auto u_ptr = boost::make_shared<VectorDouble>();
256 auto grad_u_ptr = boost::make_shared<MatrixDouble>();
257 post_proc_fe->getOpPtrVector().push_back(
259
260 post_proc_fe->getOpPtrVector().push_back(
262
264
265 post_proc_fe->getOpPtrVector().push_back(
266
267 new OpPPMap(post_proc_fe->getPostProcMesh(),
268 post_proc_fe->getMapGaussPts(),
269
270 OpPPMap::DataMapVec{{"U", u_ptr}},
271
272 OpPPMap::DataMapMat{{"GRAD_U", grad_u_ptr}},
273
275
277
278 )
279
280 );
281
282 pipeline_mng->getDomainPostProcFE() = post_proc_fe;
283 CHKERR pipeline_mng->loopFiniteElementsPostProc();
284 CHKERR post_proc_fe->writeFile("out_result.h5m");
285
287}
MoFEMErrorCode loopFiniteElementsPostProc(SmartPetscObj< DM > dm=nullptr)
Iterate postprocessing finite elements.
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
Post post-proc data at points from hash maps.
std::map< std::string, ScalarDataPtr > DataMapVec
std::map< std::string, boost::shared_ptr< MatrixDouble > > DataMapMat
boost::shared_ptr< FEMethod > & getDomainPostProcFE()
Get domain postprocessing finite element.

◆ readMesh()

MoFEMErrorCode Poisson2DHomogeneous::readMesh ( )
private

[Read mesh]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 66 of file poisson_2d_homogeneous.cpp.

66 {
68
72 CHKERR mField.getInterface<MeshsetsManager>()->setMeshsetFromFile();
73
75}
Interface for managing meshsets containing materials and boundary conditions.
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
Definition Simple.cpp:191
MoFEMErrorCode getOptions()
get options
Definition Simple.cpp:180

◆ runProgram()

MoFEMErrorCode Poisson2DHomogeneous::runProgram ( )

[Check]

[Run program]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 359 of file poisson_2d_homogeneous.cpp.

359 {
361
370
372}
MoFEMErrorCode setupProblem()
[Read mesh]
MoFEMErrorCode boundaryCondition()
[Setup problem]
MoFEMErrorCode readMesh()
[Read mesh]
MoFEMErrorCode checkResults()
[Output results]
MoFEMErrorCode outputResults()
[Solve system]
MoFEMErrorCode assembleSystem()
[Boundary condition]
MoFEMErrorCode solveSystem()
[Set integration rules]
MoFEMErrorCode setIntegrationRules()
[Assemble system]

◆ setIntegrationRules()

MoFEMErrorCode Poisson2DHomogeneous::setIntegrationRules ( )
private

[Assemble system]

[Set integration rules]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 202 of file poisson_2d_homogeneous.cpp.

202 {
204
205 auto rule_lhs = [](int, int, int p) -> int { return 2 * (p - 1); };
206 auto rule_rhs = [](int, int, int p) -> int { return p; };
207
208 auto pipeline_mng = mField.getInterface<PipelineManager>();
209 CHKERR pipeline_mng->setDomainLhsIntegrationRule(rule_lhs);
210 CHKERR pipeline_mng->setDomainRhsIntegrationRule(rule_rhs);
211
213}
MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule)
Set integration rule for domain left-hand side finite element.

◆ setupProblem()

MoFEMErrorCode Poisson2DHomogeneous::setupProblem ( )
private

[Read mesh]

[Setup problem]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 79 of file poisson_2d_homogeneous.cpp.

79 {
81 Range domain_ents;
82 CHKERR mField.get_moab().get_entities_by_dimension(0, SPACE_DIM, domain_ents,
83 true);
84 auto get_ents_by_dim = [&](const auto dim) {
85 if (dim == SPACE_DIM) {
86 return domain_ents;
87 } else {
88 Range ents;
89 if (dim == 0)
90 CHKERR mField.get_moab().get_connectivity(domain_ents, ents, true);
91 else
92 CHKERR mField.get_moab().get_entities_by_dimension(0, dim, ents, true);
93 return ents;
94 }
95 };
96
97 // Select base
98 auto get_base = [&]() {
99 auto domain_ents = get_ents_by_dim(SPACE_DIM);
100 if (domain_ents.empty())
101 CHK_THROW_MESSAGE(MOFEM_NOT_FOUND, "Empty mesh");
102 const auto type = type_from_handle(domain_ents[0]);
103 switch (type) {
104 case MBQUAD:
106 case MBHEX:
108 case MBTRI:
110 case MBTET:
112 default:
113 CHK_THROW_MESSAGE(MOFEM_NOT_FOUND, "Element type not handled");
114 }
115 return NOBASE;
116 };
117 auto base = get_base();
119
120 int oRder = 3;
121 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, "", "-order", &oRder, PETSC_NULLPTR);
122 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, "", "-atom_test", &atom_test,
123 PETSC_NULLPTR);
125
127
129}
std::string type
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition definitions.h:60
@ NOBASE
Definition definitions.h:59
@ DEMKOWICZ_JACOBI_BASE
Definition definitions.h:66
@ H1
continuous field
Definition definitions.h:85
@ MOFEM_NOT_FOUND
Definition definitions.h:33
auto type_from_handle(const EntityHandle h)
get type from entity handle
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
virtual moab::Interface & get_moab()=0
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:261
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition Simple.cpp:575
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition Simple.cpp:735

◆ solveSystem()

MoFEMErrorCode Poisson2DHomogeneous::solveSystem ( )
private

[Set integration rules]

[Solve system]

Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 217 of file poisson_2d_homogeneous.cpp.

217 {
219
220 auto pipeline_mng = mField.getInterface<PipelineManager>();
221
222 auto ksp_solver = pipeline_mng->createKSP();
223 CHKERR KSPSetFromOptions(ksp_solver);
224 CHKERR KSPSetUp(ksp_solver);
225
226 // Create RHS and solution vectors
227 auto dm = simpleInterface->getDM();
228 auto F = createDMVector(dm);
229 auto D = vectorDuplicate(F);
230
231 // Solve the system
232 CHKERR KSPSolve(ksp_solver, F, D);
234 ksp_solver, "poisson 2d homogeneous");
235
236 // Scatter result data on the mesh
237 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
238 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
239 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
240
242}
@ F
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode, RowColData rc=RowColData::COL)
set local (or ghosted) vector values on mesh for partition only
Definition DMMoFEM.cpp:543
auto createDMVector(DM dm, RowColData rc=RowColData::COL)
Get smart vector from DM.
Definition DMMoFEM.hpp:1250
SmartPetscObj< KSP > createKSP(SmartPetscObj< DM > dm=nullptr)
Create KSP (linear) solver.
double D
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.

◆ sourceTermFunction()

static double Poisson2DHomogeneous::sourceTermFunction ( const double  x,
const double  y,
const double  z 
)
inlinestaticprivate
Examples
mofem/tutorials/scl-1_poisson_2d_homogeneous/poisson_2d_homogeneous.cpp.

Definition at line 52 of file poisson_2d_homogeneous.cpp.

53 {
54 return 2. * M_PI * M_PI * sin(M_PI * x) * sin(M_PI * y);
55 }

Member Data Documentation

◆ atom_test

int Poisson2DHomogeneous::atom_test = 0
private

◆ mField

MoFEM::Interface& Poisson2DHomogeneous::mField
private

◆ oRder

int Poisson2DHomogeneous::oRder
private

◆ petscVec

SmartPetscObj<Vec> Poisson2DHomogeneous::petscVec
private

◆ simpleInterface

Simple* Poisson2DHomogeneous::simpleInterface
private

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