v0.15.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 speed up problem setup and analysis.
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 *
38 * Retrieves command-line options and configuration settings for the Simple interface
39 *
40 * @return MoFEMErrorCode Error code indicating success or failure
41 */
43
44 typedef boost::function<MoFEMErrorCode(Interface &, const char *,
45 const char *)>
47 /**
48 * \brief Set function for loading mesh file
49 *
50 * Default function for loading mesh files using MOAB interface
51 *
52 * @param m_field MoFEM interface object
53 * @param file_name Name of the mesh file to load
54 * @param options Loader options string for MOAB
55 * @return MoFEMErrorCode Error code indicating success or failure
56 */
58 const char *file_name,
59 const char *options) {
60 // Default behavior using the member moab.load_file
61 return m_field.get_moab().load_file(file_name, 0, options);
62 }
63
64 /**
65 * \brief Load mesh file
66 *
67 * Loads mesh file with specified options and custom loading function
68 *
69 * @param options File load options string passed to the loader
70 * @param mesh_file_name Name of mesh file to load (if empty, uses default or command line)
71 * @param loadFunc Function for loading mesh file (default: defaultLoadFileFunc)
72 *
73 * \note If bitRefLevel is set to any, bit ref level of loaded entities is not
74 * changed. After mesh is load, bit ref level should be set to create problem.
75 * Default setting of bit ref level is on first bit, and if is set all mesh
76 * entities on load are set to set level.
77 *
78 * @return MoFEMErrorCode Error code indicating success or failure
79 */
80 MoFEMErrorCode loadFile(const std::string options,
81 const std::string mesh_file_name,
83 /**
84 * \brief Load mesh file with parallel options if number of cores > 1
85 *
86 * Automatically selects appropriate loading options for parallel execution
87 *
88 * @param mesh_file_name Name of mesh file to load (if empty, uses default or command line value)
89 * @return MoFEMErrorCode Error code indicating success or failure
90 */
91 MoFEMErrorCode loadFile(const std::string mesh_file_name = "");
92
93 /**
94 * \brief Add field on domain
95 *
96 * Adds a field to the domain entities (highest dimension entities in the mesh)
97 *
98 * @param name Name of the field
99 * @param space Function space (L2, H1, Hdiv, Hcurl)
100 * @param base Approximation base (see FieldApproximationBase)
101 * @param nb_of_coefficients Number of field coefficients per DOF
102 * @param tag_type Tag storage type (MB_TAG_DENSE is faster, MB_TAG_SPARSE more flexible)
103 * @param bh Behavior if field exists (MF_EXCL throws error, MF_ZERO ignores)
104 * @param verb Verbosity level (-1 for default)
105 * @return MoFEMErrorCode Error code indicating success or failure
106 */
108 addDomainField(const std::string name, const FieldSpace space,
109 const FieldApproximationBase base,
110 const FieldCoefficientsNumber nb_of_coefficients,
111 const TagType tag_type = MB_TAG_SPARSE,
112 const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
113
114 /**
115 * \brief Add broken field on domain
116 *
117 * Adds a discontinuous field to the domain entities with element-wise DOF mapping
118 *
119 * @param name Name of the field
120 * @param space Function space (L2, H1, Hdiv, Hcurl)
121 * @param base Approximation base (see FieldApproximationBase)
122 * @param nb_of_coefficients Number of field coefficients per DOF
123 * @param tag_type Tag storage type (MB_TAG_DENSE is faster, MB_TAG_SPARSE more flexible)
124 * @param bh Behavior if field exists (MF_EXCL throws error, MF_ZERO ignores)
125 * @param verb Verbosity level (-1 for default)
126 * @return MoFEMErrorCode Error code indicating success or failure
127 */
129 addDomainBrokenField(const std::string name, const FieldSpace space,
130 const FieldApproximationBase base,
131 const FieldCoefficientsNumber nb_of_coefficients,
132 const TagType tag_type = MB_TAG_SPARSE,
133 const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
134
135 /**
136 * \brief Add field on boundary
137 *
138 * Adds a field to the boundary entities (entities on the mesh boundary)
139 *
140 * @param name Name of the field
141 * @param space Function space (L2, H1, Hdiv, Hcurl)
142 * @param base Approximation base (see FieldApproximationBase)
143 * @param nb_of_coefficients Number of field coefficients per DOF
144 * @param tag_type Tag storage type (MB_TAG_DENSE is faster, MB_TAG_SPARSE more flexible)
145 * @param bh Behavior if field exists (MF_EXCL throws error, MF_ZERO ignores)
146 * @param verb Verbosity level (-1 for default)
147 * @return MoFEMErrorCode Error code indicating success or failure
148 */
150 addBoundaryField(const std::string name, const FieldSpace space,
151 const FieldApproximationBase base,
152 const FieldCoefficientsNumber nb_of_coefficients,
153 const TagType tag_type = MB_TAG_SPARSE,
154 const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
155
156 /**
157 * \brief Add field on skeleton
158 *
159 * Adds a field to the skeleton entities (internal interfaces between domain elements)
160 *
161 * @param name Name of the field
162 * @param space Function space (L2, H1, Hdiv, Hcurl)
163 * @param base Approximation base (see FieldApproximationBase)
164 * @param nb_of_coefficients Number of field coefficients per DOF
165 * @param tag_type Tag storage type (MB_TAG_DENSE is faster, MB_TAG_SPARSE more flexible)
166 * @param bh Behavior if field exists (MF_EXCL throws error, MF_ZERO ignores)
167 * @param verb Verbosity level (-1 for default)
168 * @return MoFEMErrorCode Error code indicating success or failure
169 */
171 addSkeletonField(const std::string name, const FieldSpace space,
172 const FieldApproximationBase base,
173 const FieldCoefficientsNumber nb_of_coefficients,
174 const TagType tag_type = MB_TAG_SPARSE,
175 const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
176
177 /**
178 * \brief Add data field
179 *
180 * Adds a field for storing data values (typically for material properties or constants)
181 *
182 * @param name Name of the field
183 * @param space Function space (L2, H1, Hdiv, Hcurl)
184 * @param base Approximation base (see FieldApproximationBase)
185 * @param nb_of_coefficients Number of field coefficients per DOF
186 * @param tag_type Tag storage type (MB_TAG_DENSE is faster, MB_TAG_SPARSE more flexible)
187 * @param bh Behavior if field exists (MF_EXCL throws error, MF_ZERO ignores)
188 * @param verb Verbosity level (-1 for default)
189 * @return MoFEMErrorCode Error code indicating success or failure
190 */
191 MoFEMErrorCode addDataField(const std::string name, const FieldSpace space,
192 const FieldApproximationBase base,
193 const FieldCoefficientsNumber nb_of_coefficients,
194 const TagType tag_type = MB_TAG_SPARSE,
195 const enum MoFEMTypes bh = MF_ZERO,
196 int verb = -1);
197
198 /**
199 * @brief Add meshset field
200 *
201 * Add fields managed by meshset finite element for global DOFs
202 *
203 * @param name Name of the field
204 * @param space Function space (L2, H1, Hdiv, Hcurl)
205 * @param base Approximation base (see FieldApproximationBase)
206 * @param nb_of_coefficients Number of field coefficients per DOF
207 * @param tag_type Tag storage type (MB_TAG_DENSE is faster, MB_TAG_SPARSE more flexible)
208 * @param bh Behavior if field exists (MF_EXCL throws error, MF_ZERO ignores)
209 * @param verb Verbosity level (-1 for default)
210 * @return MoFEMErrorCode Error code indicating success or failure
211 */
213 addMeshsetField(const std::string name, const FieldSpace space,
214 const FieldApproximationBase base,
215 const FieldCoefficientsNumber nb_of_coefficients,
216 const TagType tag_type = MB_TAG_SPARSE,
217 const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
218
219 /**
220 * @brief Remove field from domain
221 *
222 * @param name Name of the field to remove from domain
223 * @return MoFEMErrorCode Error code indicating success or failure
224 */
225 MoFEMErrorCode removeDomainField(const std::string name);
226
227 /**
228 * @brief Remove field from boundary
229 *
230 * @param name Name of the field to remove from boundary
231 * @return MoFEMErrorCode Error code indicating success or failure
232 */
233 MoFEMErrorCode removeBoundaryField(const std::string name);
234
235 /**
236 * @brief Remove field from skeleton
237 *
238 * @param name Name of the field to remove from skeleton
239 * @return MoFEMErrorCode Error code indicating success or failure
240 */
241 MoFEMErrorCode removeSkeletonField(const std::string name);
242
243 /**
244 * \brief Define finite elements
245 *
246 * Creates and defines domain, boundary, and skeleton finite elements based on added fields
247 *
248 * @return MoFEMErrorCode Error code indicating success or failure
249 */
251
252 /**
253 * \brief define problem
254 *
255 * Creates the problem and associates finite elements with it
256 *
257 * @param is_partitioned Whether the mesh should be partitioned for parallel execution
258 * @return MoFEMErrorCode Error code indicating success or failure
259 */
260 MoFEMErrorCode defineProblem(const PetscBool is_partitioned = PETSC_TRUE);
261
262 /**
263 * \brief Set field order
264 *
265 * Sets the approximation order for a field on specified entities
266 *
267 * @param field_name Name of the field
268 * @param order Approximation order to set
269 * @param ents Range of entities to apply order to (NULL for all entities)
270 * @return MoFEMErrorCode Error code indicating success or failure
271 */
272 MoFEMErrorCode setFieldOrder(const std::string field_name, const int order,
273 const Range *ents = NULL);
274
275 /**
276 * \brief Build fields
277 *
278 * Builds the field data structures and creates DOFs
279 *
280 * @return MoFEMErrorCode Error code indicating success or failure
281 */
283
284 /**
285 * \brief Build finite elements
286 *
287 * Builds the finite element data structures and adjacencies
288 *
289 * @return MoFEMErrorCode Error code indicating success or failure
290 */
292
293 /**
294 * @brief Set the skeleton adjacency object
295 *
296 * Sets up adjacency relationships for skeleton finite elements
297 *
298 * @param dim Dimension for skeleton adjacency (-1 for automatic detection)
299 * @param fe_name Name of finite element to set adjacency for (empty for default)
300 * @return MoFEMErrorCode Error code indicating success or failure
301 */
302 MoFEMErrorCode setSkeletonAdjacency(int dim = -1, std::string fe_name = "");
303
304 /**
305 * \brief Build problem
306 *
307 * Builds the problem data structures and sets up DOF numbering
308 *
309 * @return MoFEMErrorCode Error code indicating success or failure
310 */
312
313 /**
314 * \brief Setup problem
315 *
316 * Complete setup of the problem including fields, finite elements, and problem definition
317 *
318 * @param is_partitioned Whether the mesh should be partitioned for parallel execution
319 * @return MoFEMErrorCode Error code indicating success or failure
320 */
321 MoFEMErrorCode setUp(const PetscBool is_partitioned = PETSC_TRUE);
322
323 /**
324 * @brief Rebuild internal MoFEM data structures
325 *
326 * Call this function after you add field or remove it.
327 *
328 * \note If you add field, or remove it, finite element and problem needs to
329 * be rebuild. However DM can remain the same.
330 *
331 * @param only_dm If true, only rebuild DM without rebuilding fields and finite elements
332 * @return MoFEMErrorCode Error code indicating success or failure
333 */
334 MoFEMErrorCode reSetUp(bool only_dm = false);
335
336 /**
337 * \brief Get DM
338 *
339 * Retrieves the PETSc Discrete Manager (DM) object
340 *
341 * @param dm Pointer to store the DM object
342 * @return MoFEMErrorCode Error code indicating success or failure
343 */
344 MoFEMErrorCode getDM(DM *dm);
345
346 /**
347 * @brief Return smart DM object
348 *
349 * Returns a smart pointer to the DM object that handles automatic cleanup
350 *
351 * \code
352 * {
353 * auto dm = simple_interface->getDM();
354 *
355 * // ...
356 *
357 * // dm is automatically destroyed when out of the scope
358 * }
359 * \endcode
360 *
361 * @return SmartPetscObj<DM> Smart pointer to DM object
362 */
363 inline SmartPetscObj<DM> getDM() { return dM; }
364
365 /**
366 * @brief Get the problem dimension
367 *
368 * Problem dimension is determined by highest dimension of entities on the
369 * mesh.
370 *
371 * @return int Problem dimension (1D, 2D, or 3D)
372 */
373 inline int getDim() const { return dIm; }
374
375 /**
376 * @brief Set the problem dimension
377 *
378 * @param dim Problem dimension to set (1, 2, or 3)
379 */
380 void setDim(int dim) { dIm = dim; };
381
382 /**
383 * @deprecated Use getMeshset
384 *
385 * @return EntityHandle& Reference to domain meshset handle
386 */
388
389 /**
390 * @brief Get the MeshSet object
391 *
392 * @return EntityHandle& Reference to domain meshset handle
393 */
394 inline EntityHandle &getMeshset() { return meshSet; }
395
396 /**
397 * @brief Get the BoundaryMeshSet object
398 *
399 * @return EntityHandle& Reference to boundary meshset handle
400 */
402
403 /**
404 * @brief Get the SkeletonMeshSet object
405 *
406 * @return EntityHandle& Reference to skeleton meshset handle
407 */
409
410 /**
411 * @brief Get the BitRefLevel
412 *
413 * @return BitRefLevel& Reference to bit refinement level
414 */
415 inline BitRefLevel &getBitRefLevel() { return bitLevel; }
416
417 /**
418 * @brief Get the BitRefLevelMask
419 *
420 * @return BitRefLevel& Reference to bit refinement level mask
421 */
423
424 /**
425 * @brief Get the Domain FE Name
426 *
427 * @return const std::string
428 */
429 inline const std::string getDomainFEName() const { return domainFE; }
430
431 /**
432 * @brief Get the Boundary FE Name
433 *
434 * @return const std::string
435 */
436 inline const std::string getBoundaryFEName() const { return boundaryFE; }
437
438 /**
439 * @brief Get the Skeleton FE Name
440 *
441 * @return const std::string
442 */
443 inline const std::string getSkeletonFEName() const { return skeletonFE; }
444
445 /**
446 * @brief Get the Problem Name
447 *
448 * @return const std::string
449 */
450 inline const std::string getProblemName() const { return nameOfProblem; }
451
452 /**
453 * @brief Get the Domain FE Name
454 *
455 * @return std::string&
456 */
457 inline std::string &getDomainFEName() { return domainFE; }
458
459 /**
460 * @brief Get the Boundary FE Name
461 *
462 * @return std::string&
463 */
464 inline std::string &getBoundaryFEName() { return boundaryFE; }
465
466 /**
467 * @brief Get the Skeleton FE Name
468 *
469 * @return std::string&
470 */
471 inline std::string &getSkeletonFEName() { return skeletonFE; }
472
473 /**
474 * @brief Get the Skeleton FE Name
475 *
476 * @return std::string&
477 */
478 inline std::string &getMeshsetFEName() { return meshsetFE; }
479
480 /**
481 * @brief Range of entities added to meshset finite element
482 *
483 * Meshset element have not defined integration rule, if you like to integrate
484 * you have to add mesh entity to meshset finite element. In principle
485 * meshset element become standard FE.
486 *
487 * @return std::vector<std::string>&
488 */
492
493 /**
494 * @brief Get the Problem Name
495 *
496 * @return std::string&
497 */
498 inline std::string &getProblemName() { return nameOfProblem; }
499
500 /**
501 * @brief Get the Other Finite Elements
502 *
503 * User can create finite elements using directly core interface and
504 * and them to simple problem by this function
505 *
506 * @return std::vector<std::string>&
507 */
508 inline std::vector<std::string> &getOtherFiniteElements() { return otherFEs; }
509
510 /**
511 * @brief Delete dm
512 *
513 * Destroys the PETSc DM object and cleans up resources
514 *
515 * @return MoFEMErrorCode Error code indicating success or failure
516 */
518
519 /**
520 * @brief Delete finite elements
521 *
522 * Removes all finite elements created by the Simple interface
523 *
524 * @return MoFEMErrorCode Error code indicating success or failure
525 */
527
528 /**
529 * @brief Get the addSkeletonFE flag
530 *
531 * If variable set to true, skeleton element is created regardless of whether
532 * a field on skeleton is added or not.
533 *
534 * @return bool& Reference to flag indicating if skeleton FE should be created
535 */
536 bool &getAddSkeletonFE() { return addSkeletonFE; }
537
538 /**
539 * @brief Get the addBoundaryFE flag
540 *
541 * If variable set to true, boundary element is created regardless of whether
542 * a field on boundary is added or not.
543 *
544 * @return bool& Reference to flag indicating if boundary FE should be created
545 */
546 bool &getAddBoundaryFE() { return addBoundaryFE; }
547
548 /**
549 * @brief Get the addParentAdjacencies flag
550 *
551 * If set true, parent adjacencies are added to finite elements
552 *
553 * @return bool& Reference to flag indicating if parent adjacencies should be added
554 */
556
557 /**
558 * @brief Get bit ref level for parent
559 *
560 * @return auto& Reference to bit refinement level for parent adjacencies
561 */
562 auto &getBitAdjParent() { return bitAdjParent; }
563
564 /**
565 * @brief Get bit ref level mask for parent
566 *
567 * @return auto& Reference to bit refinement level mask for parent adjacencies
568 */
570
571 /**
572 * @brief Get bit ref level for entities
573 *
574 * @return auto& Reference to bit refinement level for entity adjacencies
575 */
576 auto &getBitAdjEnt() { return bitAdjEnt; }
577
578 /**
579 * @brief Get bit ref level mask for entities
580 *
581 * @return auto& Reference to bit refinement level mask for entity adjacencies
582 */
583 auto &getBitAdjEntMask() { return bitAdjEntMask; }
584
585 /**
586 * @brief Add empty block to problem
587 *
588 * MatrixManager assumes that all blocks, i.e. all field combinations are non-zero.
589 * This is not always the case, to optimize code and reduce memory usage
590 * user can specify which blocks are empty.
591 *
592 * @param row_field Name of the row field
593 * @param col_field Name of the column field
594 * @return MoFEMErrorCode Error code indicating success or failure
595 */
596 MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field,
597 const std::string col_field) const;
598
599private:
601
602 BitRefLevel bitLevel; ///< BitRefLevel of the problem
603 BitRefLevel bitLevelMask; ///< BitRefLevel of the problem
604
611
612 EntityHandle meshSet; ///< domain meshset
613 EntityHandle boundaryMeshset; ///< meshset with boundary
614 EntityHandle skeletonMeshset; ///< skeleton meshset with boundary
615
616 bool addSkeletonFE; ///< Add skeleton FE
617 bool addBoundaryFE; ///< Add boundary FE
618 bool addParentAdjacencies; ///< If set to true parent adjacencies are build
619
620 BitRefLevel bitAdjParent; ///< bit ref level for parent
621 BitRefLevel bitAdjParentMask; ///< bit ref level for parent parent
622 BitRefLevel bitAdjEnt; ///< bit ref level for parent
623 BitRefLevel bitAdjEntMask; ///< bit ref level for parent parent
624
625 std::vector<std::string> domainFields; ///< domain fields
626 std::vector<std::string> boundaryFields; ///< boundary fields
627 std::vector<std::string> skeletonFields; ///< fields on the skeleton
628 std::vector<std::string> dataFields; ///< Data fields
629 std::vector<std::string> meshsetFields; ///< meshset fields
630 std::vector<std::string> noFieldFields; ///< NOFIELD field name
631 std::vector<std::string> noFieldDataFields; ///< NOFIELD field name
632
633 std::list<std::tuple<std::string, int, Range, bool>>
634 fieldsOrder; ///< fields order. 1: field name, order, range, set by range
635 ///< if true
636
637 std::string nameOfProblem = "SimpleProblem"; ///< problem name
638 std::string domainFE = "dFE"; ///< domain finite element
639 std::string boundaryFE = "bFE"; ///< boundary finite element
640 std::string skeletonFE = "sFE"; ///< skeleton finite element
641 std::string meshsetFE = "mFE"; ///< meshset finite element
642
643 std::vector<std::string> otherFEs; ///< Other finite elements
644
645 std::vector<Range> meshsetFiniteElementEntities; ///< Meshset element entities
646
647 char meshFileName[255]; ///< mesh file name
648 int dIm; ///< dimension of problem
649
651 dM; ///< Discrete manager (interface to PETSc/MoFEM functions)
652
656
657 template <int DIM = -1>
658 MoFEMErrorCode setSkeletonAdjacency(std::string fe_name);
659
660 template <int DIM = -1> MoFEMErrorCode setParentAdjacency();
661
662 boost::shared_ptr<ParentFiniteElementAdjacencyFunction<3>>
664 boost::shared_ptr<ParentFiniteElementAdjacencyFunction<2>>
666 boost::shared_ptr<ParentFiniteElementAdjacencyFunction<1>>
668
669 boost::shared_ptr<ParentFiniteElementAdjacencyFunctionSkeleton<2>>
671 boost::shared_ptr<ParentFiniteElementAdjacencyFunctionSkeleton<1>>
673};
674
675} // namespace MoFEM
676
677#endif // __SIMPLE_HPP__
678
679/**
680 * \defgroup mofem_simple_interface Simple interface
681 * \brief Implementation of simple interface for fast problem set-up.
682 *
683 * \ingroup mofem
684 **/
MoFEM interface.
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
@ MF_ZERO
FieldApproximationBase
approximation base
Definition definitions.h:58
FieldSpace
approximation spaces
Definition definitions.h:82
#define DEPRECATED
Definition definitions.h:17
constexpr int order
char mesh_file_name[255]
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
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:665
EntityHandle boundaryMeshset
meshset with boundary
Definition Simple.hpp:613
MoFEMErrorCode buildProblem()
Build problem.
Definition Simple.cpp:722
std::string & getDomainFEName()
Get the Domain FE Name.
Definition Simple.hpp:457
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
bool & getAddBoundaryFE()
Get the addBoundaryFE flag.
Definition Simple.hpp:546
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
Definition Simple.hpp:436
std::vector< std::string > dataFields
Data fields.
Definition Simple.hpp:628
MoFEMErrorCode setParentAdjacency()
Definition Simple.cpp:67
std::list< std::tuple< std::string, int, Range, bool > > fieldsOrder
Definition Simple.hpp:634
std::vector< std::string > domainFields
domain fields
Definition Simple.hpp:625
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition Simple.cpp:471
MoFEM::Core & cOre
Definition Simple.hpp:600
BitRefLevel bitAdjEnt
bit ref level for parent
Definition Simple.hpp:622
EntityHandle & getMeshset()
Get the MeshSet object.
Definition Simple.hpp:394
char meshFileName[255]
mesh file name
Definition Simple.hpp:647
EntityHandle meshSet
domain meshset
Definition Simple.hpp:612
auto & getBitAdjParent()
Get bit ref level for parent.
Definition Simple.hpp:562
int getDim() const
Get the problem dimension.
Definition Simple.hpp:373
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name, LoadFileFunc loadFunc=defaultLoadFileFunc)
Load mesh file.
Definition Simple.cpp:191
std::vector< std::string > otherFEs
Other finite elements.
Definition Simple.hpp:643
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition Simple.cpp:660
MoFEMErrorCode addMeshsetField(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 meshset field.
Definition Simple.cpp:411
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 2 > > parentAdjSkeletonFunctionDim2
Definition Simple.hpp:670
MoFEMErrorCode exchangeGhostCells()
Definition Simple.cpp:925
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition Simple.hpp:607
const std::string getSkeletonFEName() const
Get the Skeleton FE Name.
Definition Simple.hpp:443
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition Simple.hpp:606
static MoFEMErrorCode defaultLoadFileFunc(Interface &m_field, const char *file_name, const char *options)
Set function for loading mesh file.
Definition Simple.hpp:57
MoFEMErrorCode reSetUp(bool only_dm=false)
Rebuild internal MoFEM data structures.
Definition Simple.cpp:762
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:355
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition Simple.hpp:630
std::string meshsetFE
meshset finite element
Definition Simple.hpp:641
bool addSkeletonFE
Add skeleton FE.
Definition Simple.hpp:616
BitRefLevel bitAdjParent
bit ref level for parent
Definition Simple.hpp:620
bool addBoundaryFE
Add boundary FE.
Definition Simple.hpp:617
std::vector< std::string > & getOtherFiniteElements()
Get the Other Finite Elements.
Definition Simple.hpp:508
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
Add empty block to problem.
Definition Simple.cpp:835
std::string skeletonFE
skeleton finite element
Definition Simple.hpp:640
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition Simple.hpp:631
std::vector< std::string > meshsetFields
meshset fields
Definition Simple.hpp:629
EntityHandle & getSkeletonMeshSet()
Get the SkeletonMeshSet object.
Definition Simple.hpp:408
std::string boundaryFE
boundary finite element
Definition Simple.hpp:639
std::string nameOfProblem
problem name
Definition Simple.hpp:637
void setDim(int dim)
Set the problem dimension.
Definition Simple.hpp:380
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 3 > > parentAdjFunctionDim3
Definition Simple.hpp:663
MoFEMErrorCode getOptions()
get options
Definition Simple.cpp:180
virtual ~Simple()=default
std::string domainFE
domain finite element
Definition Simple.hpp:638
std::vector< std::string > boundaryFields
boundary fields
Definition Simple.hpp:626
boost::shared_ptr< ParentFiniteElementAdjacencyFunctionSkeleton< 1 > > parentAdjSkeletonFunctionDim1
Definition Simple.hpp:672
std::string & getBoundaryFEName()
Get the Boundary FE Name.
Definition Simple.hpp:464
MoFEMErrorCode removeSkeletonField(const std::string name)
Remove field from skeleton.
Definition Simple.cpp:457
BitRefLevel bitAdjEntMask
bit ref level for parent parent
Definition Simple.hpp:623
auto & getBitAdjEnt()
Get bit ref level for entities.
Definition Simple.hpp:576
DEPRECATED EntityHandle & getMeshSet()
Definition Simple.hpp:387
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition Simple.hpp:651
MoFEMErrorCode removeDomainField(const std::string name)
Remove field from domain.
Definition Simple.cpp:429
MoFEMErrorCode buildFields()
Build fields.
Definition Simple.cpp:584
SmartPetscObj< DM > getDM()
Return smart DM object.
Definition Simple.hpp:363
std::string & getMeshsetFEName()
Get the Skeleton FE Name.
Definition Simple.hpp:478
std::vector< Range > meshsetFiniteElementEntities
Meshset element entities.
Definition Simple.hpp:645
BitRefLevel bitLevel
BitRefLevel of the problem.
Definition Simple.hpp:602
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition Simple.cpp:575
auto & getMeshsetFiniteElementEntities()
Range of entities added to meshset finite element.
Definition Simple.hpp:489
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition Simple.cpp:551
BitRefLevel & getBitRefLevelMask()
Get the BitRefLevelMask.
Definition Simple.hpp:422
std::vector< std::string > skeletonFields
fields on the skeleton
Definition Simple.hpp:627
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 data field.
Definition Simple.cpp:393
EntityHandle & getBoundaryMeshSet()
Get the BoundaryMeshSet object.
Definition Simple.hpp:401
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition Simple.cpp:8
auto & getBitAdjEntMask()
Get bit ref level mask for entities.
Definition Simple.hpp:583
bool & getAddSkeletonFE()
Get the addSkeletonFE flag.
Definition Simple.hpp:536
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition Simple.cpp:736
boost::function< MoFEMErrorCode(Interface &, const char *, const char *)> LoadFileFunc
Definition Simple.hpp:46
MoFEMErrorCode removeBoundaryField(const std::string name)
Remove field from boundary.
Definition Simple.cpp:443
std::string & getSkeletonFEName()
Get the Skeleton FE Name.
Definition Simple.hpp:471
const std::string getProblemName() const
Get the Problem Name.
Definition Simple.hpp:450
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition Simple.hpp:429
int dIm
dimension of problem
Definition Simple.hpp:648
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:812
std::string & getProblemName()
Get the Problem Name.
Definition Simple.hpp:498
PetscLogEvent MOFEM_EVENT_SimpleSetUP
Definition Simple.hpp:605
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition Simple.hpp:610
MoFEMErrorCode createBoundaryMeshset()
Definition Simple.cpp:844
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition Simple.hpp:608
BitRefLevel bitLevelMask
BitRefLevel of the problem.
Definition Simple.hpp:603
auto & getBitAdjParentMask()
Get bit ref level mask for parent.
Definition Simple.hpp:569
BitRefLevel bitAdjParentMask
bit ref level for parent parent
Definition Simple.hpp:621
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition Simple.hpp:609
boost::shared_ptr< ParentFiniteElementAdjacencyFunction< 1 > > parentAdjFunctionDim1
Definition Simple.hpp:667
BitRefLevel & getBitRefLevel()
Get the BitRefLevel.
Definition Simple.hpp:415
bool & getParentAdjacencies()
Get the addParentAdjacencies flag.
Definition Simple.hpp:555
MoFEMErrorCode deleteFiniteElements()
Delete finite elements.
Definition Simple.cpp:823
EntityHandle skeletonMeshset
skeleton meshset with boundary
Definition Simple.hpp:614
MoFEMErrorCode createSkeletonMeshset()
Definition Simple.cpp:891
bool addParentAdjacencies
If set to true parent adjacencies are build.
Definition Simple.hpp:618
MoFEMErrorCode addDomainBrokenField(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 broken field on domain.
Definition Simple.cpp:282
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:373
intrusive_ptr for managing petsc objects
base class for all interface classes