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