v0.13.2
Loading...
Searching...
No Matches
Simple.hpp
Go to the documentation of this file.
1/** \file Simple.hpp
2 * \brief Header file for simple interface
3 * \ingroup mofem_simple_interface
4 *
5 * Make simplified interface, to speedup problem setup and analysts.
6 * See discussion here
7 * <a
8 * href=https://groups.google.com/d/msg/mofem-group/Vkc00aia4dU/o9RF3ZmPAAAJ>link
9 * to google groups</a>
10 *
11 */
12
13#ifndef __SIMPLE_HPP__
14#define __SIMPLE_HPP__
15
16#include "UnknownInterface.hpp"
17
18namespace MoFEM {
19
20template <int DIM> struct ParentFiniteElementAdjacencyFunction;
21template <int DIM> struct ParentFiniteElementAdjacencyFunctionSkeleton;
22
23/**
24 * \brief Simple interface for fast problem set-up
25 * \ingroup mofem_simple_interface
26 */
27struct Simple : public UnknownInterface {
28
29 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
30 UnknownInterface **iface) const;
31
32 Simple(const MoFEM::Core &core);
33 virtual ~Simple() = default;
34
35 /**
36 * \brief get options
37 * @return error code
38 */
40
41 /**
42 * \brief Load mesh file
43 * @param options file load options
44 * @param mesh_file_name file name if not set default or set by command line
45 * is used.
46 *
47 * \note If bitRefLevel is set to any, bit ref level of loaded entities is not
48 * changed. After mesh is load, bit ref level should be set to create problem.
49 * Default setting of bit ref level is on first bit, and if is set all mesh
50 * entities on load are set to set level.
51 *
52 * @return error code
53 */
54 MoFEMErrorCode loadFile(const std::string options,
55 const std::string mesh_file_name);
56 /**
57 * \brief Load mesh file with parallel options if number of cores > 1
58 * @param mesh_file_name file name if not set default or set by command line
59 * is used.
60 * @return error code
61 */
62 MoFEMErrorCode loadFile(const std::string mesh_file_name = "");
63 /**
64 * \brief Add field on domain
65 * @param name name of the filed
66 * @param space space (L2,H1,Hdiv,Hcurl)
67 * @param base approximation base, see FieldApproximationBase
68 * @param nb_of_coefficients number of field coefficients
69 * @param tag_type type of the tag MB_TAG_DENSE or MB_TAG_SPARSE
70 * (DENSE is faster and uses less memory, SPARSE is more flexible if you
71 * define field on subdomains)
72 * @param bh if MF_EXCL throws error if field exits, MF_ZERO
73 * no error if field exist
74 * @param verb verbosity level
75 * @return error code
76 */
78 addDomainField(const std::string &name, const FieldSpace space,
79 const FieldApproximationBase base,
80 const FieldCoefficientsNumber nb_of_coefficients,
81 const TagType tag_type = MB_TAG_SPARSE,
82 const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
83
84 /**
85 * \brief Add field on boundary
86 * @param name name of the filed
87 * @param space space (L2,H1,Hdiv,Hcurl)
88 * @param base approximation base, see FieldApproximationBase
89 * @param nb_of_coefficients number of field coefficients
90 * @param tag_type type of the tag MB_TAG_DENSE or MB_TAG_SPARSE
91 * (DENSE is faster and uses less memory, SPARSE is more flexible if you
92 * define field on subdomains)
93 * @param bh if MF_EXCL throws error if field exits, MF_ZERO
94 * no error if field exist
95 * @param verb verbosity level
96 * @return error code
97 */
99 addBoundaryField(const std::string &name, const FieldSpace space,
100 const FieldApproximationBase base,
101 const FieldCoefficientsNumber nb_of_coefficients,
102 const TagType tag_type = MB_TAG_SPARSE,
103 const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
104
105 /**
106 * \brief Add field on skeleton
107 * @param name name of the filed
108 * @param space space (L2,H1,Hdiv,Hcurl)
109 * @param base approximation base, see FieldApproximationBase
110 * @param nb_of_coefficients number of field coefficients
111 * @param tag_type type of the tag MB_TAG_DENSE or MB_TAG_SPARSE
112 * (DENSE is faster and uses less memory, SPARSE is more flexible if you
113 * define field on subdomains)
114 * @param bh if MF_EXCL throws error if field exits, MF_ZERO
115 * no error if field exist
116 * @param verb verbosity level
117 * @return error code
118 */
120 addSkeletonField(const std::string &name, const FieldSpace space,
121 const FieldApproximationBase base,
122 const FieldCoefficientsNumber nb_of_coefficients,
123 const TagType tag_type = MB_TAG_SPARSE,
124 const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
125
126 /**
127 * \brief Add field on domain
128 * @param name name of the filed
129 * @param space space (L2,H1,Hdiv,Hcurl)
130 * @param base approximation base, see FieldApproximationBase
131 * @param nb_of_coefficients number of field coefficients
132 * @param tag_type type of the tag MB_TAG_DENSE or MB_TAG_SPARSE
133 * (DENSE is faster and uses less memory, SPARSE is more flexible if you
134 * define field on subdomains)
135 * @param bh if MF_EXCL throws error if field exits, MF_ZERO
136 * no error if field exist
137 * @param verb verbosity level
138 * @return error code
139 */
140 MoFEMErrorCode addDataField(const std::string &name, const FieldSpace space,
141 const FieldApproximationBase base,
142 const FieldCoefficientsNumber nb_of_coefficients,
143 const TagType tag_type = MB_TAG_SPARSE,
144 const enum MoFEMTypes bh = MF_ZERO,
145 int verb = -1);
146
147 /**
148 * @brief Remove field form domain
149 *
150 * @param name
151 * @return MoFEMErrorCode
152 */
153 MoFEMErrorCode removeDomainField(const std::string &name);
154
155 /**
156 * @brief Remove field form boundary
157 *
158 * @param name
159 * @return MoFEMErrorCode
160 */
161 MoFEMErrorCode removeBoundaryField(const std::string &name);
162
163 /**
164 * @brief Remove field form skeleton
165 *
166 * @param name
167 * @return MoFEMErrorCode
168 */
169 MoFEMErrorCode removeSkeletonField(const std::string &name);
170
171 /**
172 * \brief Define finite elements
173 * @return Error code
174 */
176
177 /**
178 * \brief define problem
179 * @return error code
180 */
181 MoFEMErrorCode defineProblem(const PetscBool is_partitioned = PETSC_TRUE);
182
183 /**
184 * \brief Set field order
185 * @param std::field_name field name
186 * @param order order
187 * @param range of entities to which order is set (If null it sat to all
188 * entities)
189 * @return error code
190 */
191 MoFEMErrorCode setFieldOrder(const std::string field_name, const int order,
192 const Range *ents = NULL);
193
194 /**
195 * \brief Build fields
196 * @return error code
197 */
199
200 /**
201 * \brief Build finite elements
202 * @return error code
203 */
205
206 /**
207 * @brief Set the skeleton adjacency object
208 *
209 * @param dim
210 * @param fe_name
211 * @return MoFEMErrorCode
212 */
213 MoFEMErrorCode setSkeletonAdjacency(int dim = -1, std::string fe_name = "");
214
215 /**
216 * \brief Build problem
217 * @return error code
218 */
220
221 /**
222 * \brief Setup problem
223 * @return error code
224 */
225 MoFEMErrorCode setUp(const PetscBool is_partitioned = PETSC_TRUE);
226
227 /**
228 * @brief Rebuild internal MoFEM data structures
229 *
230 * Call this function after you add field or remove it.
231 *
232 * \note If you add field, or remove it, finite element and problem needs to
233 * be rebuild. However DM can remain the same.
234 *
235 * @return MoFEMErrorCode
236 */
237 MoFEMErrorCode reSetUp(bool only_dm = false);
238
239 /**
240 * \brief Get DM
241 * @param dm discrete manager
242 * @return error code
243 */
244 MoFEMErrorCode getDM(DM *dm);
245
246 /**
247 * @brief Return smart DM object
248 *
249 * \code
250 * {
251 * auto dm = simple_interface->getDM();
252 *
253 * // ...
254 *
255 * // dm is automatically destroyed when out of the scope
256 * }
257 * \endcode
258 *
259 * @return SmartPetscObj<DM>
260 */
261 inline SmartPetscObj<DM> getDM() { return dM; }
262
263 /**
264 * @brief Get the problem dimension
265 *
266 * Problem dimension is determined by highest dimension of entities on the
267 * mesh.
268 *
269 * @return int
270 */
271 inline int getDim() const { return dIm; }
272
273 /**
274 * @brief Set the problem dimension
275 *
276 * @return int
277 */
278 void setDim(int dim) { dIm = dim; };
279
280 /**
281 * @deprecated Use getMeshset
282 *
283 * @return EntityHandle&
284 */
286
287 /**
288 * @brief Get the MeshSet object
289 *
290 * @return EntityHandle&
291 */
292 inline EntityHandle &getMeshset() { return meshSet; }
293
294 /**
295 * @brief Get the BoundaryMeshSet object
296 *
297 * @return EntityHandle&
298 */
300
301 /**
302 * @brief Get the SkeletonMeshSet object
303 *
304 * @return EntityHandle&
305 */
307
308 /**
309 * @brief Get the BitRefLevel
310 *
311 * @return BitRefLevel
312 */
313 inline BitRefLevel &getBitRefLevel() { return bitLevel; }
314
315 /**
316 * @brief Get the BitRefLevel
317 *
318 * @return BitRefLevel
319 */
321
322 /**
323 * @brief Get the Domain FE Name
324 *
325 * @return const std::string
326 */
327 inline const std::string getDomainFEName() const { return domainFE; }
328
329 /**
330 * @brief Get the Boundary FE Name
331 *
332 * @return const std::string
333 */
334 inline const std::string getBoundaryFEName() const { return boundaryFE; }
335
336 /**
337 * @brief Get the Skeleton FE Name
338 *
339 * @return const std::string
340 */
341 inline const std::string getSkeletonFEName() const { return skeletonFE; }
342
343 /**
344 * @brief Get the Problem Name
345 *
346 * @return const std::string
347 */
348 inline const std::string getProblemName() const { return nameOfProblem; }
349
350 /**
351 * @brief Get the Domain FE Name
352 *
353 * @return std::string&
354 */
355 inline std::string &getDomainFEName() { return domainFE; }
356
357 /**
358 * @brief Get the Boundary FE Name
359 *
360 * @return std::string&
361 */
362 inline std::string &getBoundaryFEName() { return boundaryFE; }
363
364 /**
365 * @brief Get the Skeleton FE Name
366 *
367 * @return std::string&
368 */
369 inline std::string &getSkeletonFEName() { return skeletonFE; }
370
371 /**
372 * @brief Get the Problem Name
373 *
374 * @return std::string&
375 */
376 inline std::string &getProblemName() { return nameOfProblem; }
377
378 /**
379 * @brief Get the Other Finite Elements
380 *
381 * User can create finite elements using directly core interface and
382 * and them to simple problem by this function
383 *
384 * @return std::vector<std::string>&
385 */
386 inline std::vector<std::string> &getOtherFiniteElements() { return otherFEs; }
387
388 /**
389 * @brief Delete dm
390 *
391 * @return MoFEMErrorCode
392 */
394
395 /**
396 * @brief Delete finite elements
397 *
398 * @return MoFEMErrorCode
399 */
401
402 /**
403 * @brief Get the addSkeletonFE
404 *
405 * If variable set to true, skeleton element is created regardless field on
406 * skelton is added or not.
407 *
408 * @return true
409 * @return false
410 */
411 bool &getAddSkeletonFE() { return addSkeletonFE; }
412
413 /**
414 * @brief Get the addSkeletonFE
415 *
416 * If variable set to true, boundary element is created regardless field on
417 * skelton is added or not.
418 *
419 * @return true
420 * @return false
421 */
422 bool &getAddBoundaryFE() { return addBoundaryFE; }
423
424 /**
425 * @brief Get the addParentAdjacencies
426 *
427 * If set true add parent adjacencies
428 *
429 * @return true
430 * @return false
431 */
433
434 /**
435 * @brief bit ref level for parent
436 *
437 * @return auto&
438 */
439 auto &getBitAdjParent() { return bitAdjParent; }
440
441 /**
442 * @brief bit ref level for parent parent
443 *
444 * @return auto&
445 */
447
448 /**
449 * @brief bit ref level for parent
450 *
451 * @return auto&
452 */
453 auto &getBitAdjEnt() { return bitAdjEnt; }
454
455 /**
456 * @brief bit ref level for parent parent
457 *
458 * @return auto&
459 */
460 auto &getBitAdjEntMask() { return bitAdjEntMask; }
461
462 /**
463 * @brief add empty block to problem
464 *
465 * MatrixManager assumes that all blocks, i.e. all fields combinations are non
466 * zero. This is not always the case, to optimise code and reduce memory usage
467 * user can specifi which blocks are empty.
468 *
469 * @param row_field row filed name
470 * @param col_field col field name
471 * @return MoFEMErrorCode
472 */
473 MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field,
474 const std::string col_field) const;
475
476
477private:
479
480 BitRefLevel bitLevel; ///< BitRefLevel of the problem
481 BitRefLevel bitLevelMask; ///< BitRefLevel of the problem
482
489
490 EntityHandle meshSet; ///< domain meshset
491 EntityHandle boundaryMeshset; ///< meshset with boundary
492 EntityHandle skeletonMeshset; ///< skeleton meshset with boundary
493
494 bool addSkeletonFE; ///< Add skeleton FE
495 bool addBoundaryFE; ///< Add boundary FE
496 bool addParentAdjacencies; ///< If set to true parent adjacencies are build
497
498 BitRefLevel bitAdjParent; ///< bit ref level for parent
499 BitRefLevel bitAdjParentMask; ///< bit ref level for parent parent
500 BitRefLevel bitAdjEnt; ///< bit ref level for parent
501 BitRefLevel bitAdjEntMask; ///< bit ref level for parent parent
502
503 std::vector<std::string> domainFields; ///< domain fields
504 std::vector<std::string> boundaryFields; ///< boundary fields
505 std::vector<std::string> skeletonFields; ///< fields on the skeleton
506 std::vector<std::string> dataFields; ///< Data fields
507 std::vector<std::string> noFieldFields; ///< NOFIELD field name
508 std::vector<std::string> noFieldDataFields; ///< NOFIELD field name
509
510 std::list<std::tuple<std::string, int, Range, bool>>
511 fieldsOrder; ///< fields order. 1: field name, order, range, set by range
512 ///< if true
513
514 std::string nameOfProblem; ///< problem name
515 std::string domainFE; ///< domain finite element
516 std::string boundaryFE; ///< boundary finite element
517 std::string skeletonFE; ///< skeleton finite element
518
519 std::vector<std::string> otherFEs; ///< Other finite elements
520
521 char meshFileName[255]; ///< mesh file name
522 int dIm; ///< dimension of problem
523
525 dM; ///< Discrete manager (interface to PETSc/MoFEM functions)
526
530
531 template <int DIM = -1>
532 MoFEMErrorCode setSkeletonAdjacency(std::string fe_name);
533
534 template <int DIM = -1> MoFEMErrorCode setParentAdjacency();
535
536 boost::shared_ptr<ParentFiniteElementAdjacencyFunction<3>>
538 boost::shared_ptr<ParentFiniteElementAdjacencyFunction<2>>
540 boost::shared_ptr<ParentFiniteElementAdjacencyFunction<1>>
542
543 boost::shared_ptr<ParentFiniteElementAdjacencyFunctionSkeleton<2>>
545 boost::shared_ptr<ParentFiniteElementAdjacencyFunctionSkeleton<1>>
547};
548
549} // namespace MoFEM
550
551#endif // __SIMPLE_HPP__
552
553/**
554 * \defgroup mofem_simple_interface Simple interface
555 * \brief Implementation of simple interface for fast problem set-up.
556 *
557 * \ingroup mofem
558 **/
MoFEM interface.
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:97
@ MF_ZERO
Definition: definitions.h:98
FieldApproximationBase
approximation base
Definition: definitions.h:58
FieldSpace
approximation spaces
Definition: definitions.h:82
#define DEPRECATED
Definition: definitions.h:17
const int dim
char mesh_file_name[255]
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
int FieldCoefficientsNumber
Number of field coefficients.
Definition: Types.hpp:27
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
constexpr auto field_name
Core (interface) class.
Definition: Core.hpp:82
Simple interface for fast problem set-up.
Definition: Simple.hpp:27
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 2 > > parentAdjFunctionDim2
Definition: Simple.hpp:539
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:491
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:597
std::string & getDomainFEName()
Get the Domain FE Name.
Definition: Simple.hpp:355
bool & getAddBoundaryFE()
Get the addSkeletonFE.
Definition: Simple.hpp:422
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
Definition: Simple.hpp:334
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:506
MoFEMErrorCode setParentAdjacency()
Definition: Simple.cpp:67
std::list< std::tuple< std::string, int, Range, bool > > fieldsOrder
Definition: Simple.hpp:511
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:503
MoFEMErrorCode removeDomainField(const std::string &name)
Remove field form domain.
Definition: Simple.cpp:337
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:380
MoFEM::Core & cOre
Definition: Simple.hpp:478
BitRefLevel bitAdjEnt
bit ref level for parent
Definition: Simple.hpp:500
EntityHandle & getMeshset()
Get the MeshSet object.
Definition: Simple.hpp:292
char meshFileName[255]
mesh file name
Definition: Simple.hpp:521
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:490
auto & getBitAdjParent()
bit ref level for parent
Definition: Simple.hpp:439
int getDim() const
Get the problem dimension.
Definition: Simple.hpp:271
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:519
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:571
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 2 > > parentAdjSkeletonFunctionDim2
Definition: Simple.hpp:544
MoFEMErrorCode addDataField(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:320
MoFEMErrorCode exchangeGhostCells()
Definition: Simple.cpp:792
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:485
const std::string getSkeletonFEName() const
Get the Skeleton FE Name.
Definition: Simple.hpp:341
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:484
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 reSetUp(bool only_dm=false)
Rebuild internal MoFEM data structures.
Definition: Simple.cpp:633
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:507
bool addSkeletonFE
Add skeleton FE.
Definition: Simple.hpp:494
BitRefLevel bitAdjParent
bit ref level for parent
Definition: Simple.hpp:498
bool addBoundaryFE
Add boundary FE.
Definition: Simple.hpp:495
std::vector< std::string > & getOtherFiniteElements()
Get the Other Finite Elements.
Definition: Simple.hpp:386
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
add empty block to problem
Definition: Simple.cpp:702
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:517
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:508
EntityHandle & getSkeletonMeshSet()
Get the SkeletonMeshSet object.
Definition: Simple.hpp:306
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:516
std::string nameOfProblem
problem name
Definition: Simple.hpp:514
void setDim(int dim)
Set the problem dimension.
Definition: Simple.hpp:278
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
Definition: Simple.hpp:537
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:180
virtual ~Simple()=default
std::string domainFE
domain finite element
Definition: Simple.hpp:515
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:504
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 1 > > parentAdjSkeletonFunctionDim1
Definition: Simple.hpp:546
std::string & getBoundaryFEName()
Get the Boundary FE Name.
Definition: Simple.hpp:362
BitRefLevel bitAdjEntMask
bit ref level for parent parent
Definition: Simple.hpp:501
auto & getBitAdjEnt()
bit ref level for parent
Definition: Simple.hpp:453
DEPRECATED EntityHandle & getMeshSet()
Definition: Simple.hpp:285
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:525
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:482
SmartPetscObj< DM > getDM()
Return smart DM object.
Definition: Simple.hpp:261
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:194
BitRefLevel bitLevel
BitRefLevel of the problem.
Definition: Simple.hpp:480
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:473
MoFEMErrorCode addBoundaryField(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 boundary.
Definition: Simple.cpp:282
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:452
BitRefLevel & getBitRefLevelMask()
Get the BitRefLevel.
Definition: Simple.hpp:320
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:505
MoFEMErrorCode addSkeletonField(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 skeleton.
Definition: Simple.cpp:300
EntityHandle & getBoundaryMeshSet()
Get the BoundaryMeshSet object.
Definition: Simple.hpp:299
MoFEMErrorCode removeSkeletonField(const std::string &name)
Remove field form skeleton.
Definition: Simple.cpp:366
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: Simple.cpp:8
auto & getBitAdjEntMask()
bit ref level for parent parent
Definition: Simple.hpp:460
bool & getAddSkeletonFE()
Get the addSkeletonFE.
Definition: Simple.hpp:411
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:611
MoFEMErrorCode removeBoundaryField(const std::string &name)
Remove field form boundary.
Definition: Simple.cpp:352
std::string & getSkeletonFEName()
Get the Skeleton FE Name.
Definition: Simple.hpp:369
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:348
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:327
int dIm
dimension of problem
Definition: Simple.hpp:522
MoFEMErrorCode setSkeletonAdjacency(int dim=-1, std::string fe_name="")
Set the skeleton adjacency object.
Definition: Simple.cpp:143
MoFEMErrorCode deleteDM()
Delete dm.
Definition: Simple.cpp:679
std::string & getProblemName()
Get the Problem Name.
Definition: Simple.hpp:376
PetscLogEvent MOFEM_EVENT_SimpleSetUP
Definition: Simple.hpp:483
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:488
MoFEMErrorCode createBoundaryMeshset()
Definition: Simple.cpp:711
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:486
BitRefLevel bitLevelMask
BitRefLevel of the problem.
Definition: Simple.hpp:481
auto & getBitAdjParentMask()
bit ref level for parent parent
Definition: Simple.hpp:446
BitRefLevel bitAdjParentMask
bit ref level for parent parent
Definition: Simple.hpp:499
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:487
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 1 > > parentAdjFunctionDim1
Definition: Simple.hpp:541
BitRefLevel & getBitRefLevel()
Get the BitRefLevel.
Definition: Simple.hpp:313
bool & getParentAdjacencies()
Get the addParentAdjacencies.
Definition: Simple.hpp:432
MoFEMErrorCode deleteFiniteElements()
Delete finite elements.
Definition: Simple.cpp:690
EntityHandle skeletonMeshset
skeleton meshset with boundary
Definition: Simple.hpp:492
MoFEMErrorCode createSkeletonMeshset()
Definition: Simple.cpp:758
bool addParentAdjacencies
If set to true parent adjacencies are build.
Definition: Simple.hpp:496
intrusive_ptr for managing petsc objects
base class for all interface classes