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

SimplesimpleInterface
 
MoFEM::InterfacemField
 
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}
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 290 of file poisson_2d_homogeneous.cpp.

290 {
292
293 auto check_result_fe_ptr = boost::make_shared<DomainEle>(mField);
294 auto petscVec =
296 (mField.get_comm_rank() == 0) ? LAST_NORM : 0, LAST_NORM);
297
299 check_result_fe_ptr->getOpPtrVector(), {H1})),
300 "Apply transform");
301
302 check_result_fe_ptr->getRuleHook = [](int, int, int p) { return p; };
303 auto analyticalFunction = [&](double x, double y, double z) {
304 return sin(M_PI * x) * sin(M_PI * y);
305 };
306
307 auto u_ptr = boost::make_shared<VectorDouble>();
308
309 check_result_fe_ptr->getOpPtrVector().push_back(
311 auto mValFuncPtr = boost::make_shared<VectorDouble>();
312 check_result_fe_ptr->getOpPtrVector().push_back(
313 new OpGetTensor0fromFunc(mValFuncPtr, analyticalFunction));
314 check_result_fe_ptr->getOpPtrVector().push_back(
315 new OpCalcNormL2Tensor0(u_ptr, petscVec, NORM, mValFuncPtr));
316 CHKERR VecZeroEntries(petscVec);
319 check_result_fe_ptr);
320 CHKERR VecAssemblyBegin(petscVec);
321 CHKERR VecAssemblyEnd(petscVec);
322 MOFEM_LOG_CHANNEL("SELF"); // Clear channel from old tags
323 // print norm in general log
324 if (mField.get_comm_rank() == 0) {
325 const double *norms;
326 CHKERR VecGetArrayRead(petscVec, &norms);
327 MOFEM_TAG_AND_LOG("SELF", Sev::inform, "Errors")
328 << "NORM: " << std::sqrt(norms[NORM]);
329 CHKERR VecRestoreArrayRead(petscVec, &norms);
330 }
331 // compare norm for ctest
332 if (atom_test && !mField.get_comm_rank()) {
333 const double *t_ptr;
334 CHKERR VecGetArrayRead(petscVec, &t_ptr);
335 double ref_norm = 2.2e-04;
336 double cal_norm;
337 switch (atom_test) {
338 case 1: // 2D
339 cal_norm = sqrt(t_ptr[0]);
340 break;
341 default:
342 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
343 "atom test %d does not exist", atom_test);
344 }
345 if (cal_norm > ref_norm) {
346 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID,
347 "atom test %d failed! Calculated Norm %3.16e is greater than "
348 "reference Norm %3.16e",
349 atom_test, cal_norm, ref_norm);
350 }
351 CHKERR VecRestoreArrayRead(petscVec, &t_ptr);
352 }
354}
#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:576
#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 244 of file poisson_2d_homogeneous.cpp.

244 {
246
247 auto pipeline_mng = mField.getInterface<PipelineManager>();
248 pipeline_mng->getDomainLhsFE().reset();
249
250 auto post_proc_fe = boost::make_shared<PostProcFaceEle>(mField);
252 post_proc_fe->getOpPtrVector(), {H1});
253
254 auto u_ptr = boost::make_shared<VectorDouble>();
255 auto grad_u_ptr = boost::make_shared<MatrixDouble>();
256 post_proc_fe->getOpPtrVector().push_back(
258
259 post_proc_fe->getOpPtrVector().push_back(
261
263
264 post_proc_fe->getOpPtrVector().push_back(
265
266 new OpPPMap(post_proc_fe->getPostProcMesh(),
267 post_proc_fe->getMapGaussPts(),
268
269 OpPPMap::DataMapVec{{"U", u_ptr}},
270
271 OpPPMap::DataMapMat{{"GRAD_U", grad_u_ptr}},
272
274
276
277 )
278
279 );
280
281 pipeline_mng->getDomainRhsFE() = post_proc_fe;
282 CHKERR pipeline_mng->loopFiniteElements();
283 CHKERR post_proc_fe->writeFile("out_result.h5m");
284
286}
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 > & getDomainLhsFE()
Get domain left-hand side 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 358 of file poisson_2d_homogeneous.cpp.

358 {
360
369
371}
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);
233
234 // Scatter result data on the mesh
235 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
236 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
237 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
238
240}
@ 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:514
auto createDMVector(DM dm, RowColData rc=RowColData::COL)
Get smart vector from DM.
Definition DMMoFEM.hpp:1237
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: