v0.13.0
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 /* MoFEM is distributed in the hope that it will be useful, but WITHOUT
14  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  * License for more details.
17  *
18  * You should have received a copy of the GNU Lesser General Public
19  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
20  */
21 
22 #ifndef __SIMPLE_HPP__
23 #define __SIMPLE_HPP__
24 
25 #include "UnknownInterface.hpp"
26 
27 namespace MoFEM {
28 
29 /**
30  * \brief Simple interface for fast problem set-up
31  * \ingroup mofem_simple_interface
32  */
33 struct Simple : public UnknownInterface {
34 
35  MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
36  UnknownInterface **iface) const;
37 
38  Simple(const MoFEM::Core &core);
39  virtual ~Simple() = default;
40 
41  /**
42  * \brief get options
43  * @return error code
44  */
46 
47  /**
48  * \brief Load mesh file
49  * @param options file load options
50  * @param mesh_file_name file name if not set default or set by command line
51  * is used.
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  * @return MoFEMErrorCode
211  */
213 
214  /**
215  * \brief Build problem
216  * @return error code
217  */
219 
220  /**
221  * \brief Setup problem
222  * @return error code
223  */
224  MoFEMErrorCode setUp(const PetscBool is_partitioned = PETSC_TRUE);
225 
226  /**
227  * @brief Rebuild internal MoFEM data structures
228  *
229  * Call this function after you add field or remove it.
230  *
231  * \note If you add field, or remove it, finite element and problem needs to
232  * be rebuild. However DM can remain the same.
233  *
234  * @return MoFEMErrorCode
235  */
237 
238  /**
239  * \brief Get DM
240  * @param dm discrete manager
241  * @return error code
242  */
243  MoFEMErrorCode getDM(DM *dm);
244 
245  /**
246  * @brief Return smart DM object
247  *
248  * \code
249  * {
250  * auto dm = simple_interface->getDM();
251  *
252  * // ...
253  *
254  * // dm is automatically destroyed when out of the scope
255  * }
256  * \endcode
257  *
258  * @return SmartPetscObj<DM>
259  */
260  inline SmartPetscObj<DM> getDM() { return dM; }
261 
262  /**
263  * @brief Get the problem dimension
264  *
265  * Problem dimension is determined by highest dimension of entities on the
266  * mesh.
267  *
268  * @return int
269  */
270  inline int getDim() const { return dIm; }
271 
272  /**
273  * @brief Set the problem dimension
274  *
275  * @return int
276  */
277  void setDim(int dim) { dIm = dim; };
278 
279  /**
280  * @brief Get the MeshSet object
281  *
282  * @return EntityHandle&
283  */
284  inline EntityHandle &getMeshSet() { return meshSet; }
285 
286  /**
287  * @brief Get the BoundaryMeshSet object
288  *
289  * @return EntityHandle&
290  */
292 
293  /**
294  * @brief Get the SkeletonMeshSet object
295  *
296  * @return EntityHandle&
297  */
299 
300  /**
301  * @brief Get the BitRefLevel
302  *
303  * @return BitRefLevel
304  */
305  inline BitRefLevel &getBitRefLevel() { return bitLevel; }
306 
307  /**
308  * @brief Get the BitRefLevel
309  *
310  * @return BitRefLevel
311  */
313 
314  /**
315  * @brief Get the Domain FE Name
316  *
317  * @return const std::string
318  */
319  inline const std::string getDomainFEName() const { return domainFE; }
320 
321  /**
322  * @brief Get the Boundary FE Name
323  *
324  * @return const std::string
325  */
326  inline const std::string getBoundaryFEName() const { return boundaryFE; }
327 
328  /**
329  * @brief Get the Skeleton FE Name
330  *
331  * @return const std::string
332  */
333  inline const std::string getSkeletonFEName() const { return skeletonFE; }
334 
335  /**
336  * @brief Get the Problem Name
337  *
338  * @return const std::string
339  */
340  inline const std::string getProblemName() const { return nameOfProblem; }
341 
342  /**
343  * @brief Get the Domain FE Name
344  *
345  * @return std::string&
346  */
347  inline std::string &getDomainFEName() { return domainFE; }
348 
349  /**
350  * @brief Get the Boundary FE Name
351  *
352  * @return std::string&
353  */
354  inline std::string &getBoundaryFEName() { return boundaryFE; }
355 
356  /**
357  * @brief Get the Skeleton FE Name
358  *
359  * @return std::string&
360  */
361  inline std::string &getSkeletonFEName() { return skeletonFE; }
362 
363  /**
364  * @brief Get the Problem Name
365  *
366  * @return std::string&
367  */
368  inline std::string &getProblemName() { return nameOfProblem; }
369 
370  /**
371  * @brief Get the Other Finite Elements
372  *
373  * User can create finite elements using directly core interface and
374  * and them to simple problem by this function
375  *
376  * @return std::vector<std::string>&
377  */
378  inline std::vector<std::string> &getOtherFiniteElements() { return otherFEs; }
379 
380  /**
381  * @brief Delete dm
382  *
383  * @return MoFEMErrorCode
384  */
386 
387  /**
388  * @brief Delete finite elements
389  *
390  * @return MoFEMErrorCode
391  */
393 
394  /**
395  * @brief Get the Add Skeleton F E object
396  *
397  * If variable set to true, skeleton element is created regardless field on
398  * skelton is added or not.
399  *
400  * @return true
401  * @return false
402  */
403  bool &getAddSkeletonFE() { return addSkeletonFE; }
404 
405  /**
406  * @brief Get the Add Boundary F E object
407  *
408  * If variable set to true, boundary element is created regardless field on
409  * skelton is added or not.
410  *
411  * @return true
412  * @return false
413  */
414  bool &getAddBoundaryFE() { return addBoundaryFE; }
415 
416  /**
417  * @brief add empty block to problem
418  *
419  * MatrixManager assumes that all blocks, i.e. all fields combinations are non
420  * zero. This is not always the case, to optimise code and reduce memory usage
421  * user can specifi which blocks are empty.
422  *
423  * @param row_field row filed name
424  * @param col_field col field name
425  * @return MoFEMErrorCode
426  */
427  MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field,
428  const std::string col_field) const;
429 
430 private:
432 
433  BitRefLevel bitLevel; ///< BitRefLevel of the probelm
434  BitRefLevel bitLevelMask; ///< BitRefLevel of the probelm
435 
436  PetscLogEvent MOFEM_EVENT_SimpleSetUP;
442 
443  EntityHandle meshSet; ///< domain meshset
444  EntityHandle boundaryMeshset; ///< meshset with boundary
445  EntityHandle skeletonMeshset; ///< skeleton meshset with boundary
446 
447  bool addSkeletonFE; ///< Add skeleton FE
448  bool addBoundaryFE; ///< Add boundary FE
449 
450  std::vector<std::string> domainFields; ///< domain fields
451  std::vector<std::string> boundaryFields; ///< boundary fields
452  std::vector<std::string> skeletonFields; ///< fields on the skeleton
453  std::vector<std::string> dataFields; ///< Data fields
454  std::vector<std::string> noFieldFields; ///< NOFIELD field name
455  std::vector<std::string> noFieldDataFields; ///< NOFIELD field name
456 
457  std::multimap<std::string, std::pair<int, Range>>
458  fieldsOrder; ///< fields order
459 
460  std::string nameOfProblem; ///< problem name
461  std::string domainFE; ///< domain finite element
462  std::string boundaryFE; ///< boundary finite element
463  std::string skeletonFE; ///< skeleton finite element
464 
465  std::vector<std::string> otherFEs; ///< Other finite elements
466 
467  char meshFileName[255]; ///< mesh file name
468  int dIm; ///< dimension of problem
469 
471  dM; ///< Discrete manager (interface to PETSc/MoFEM functions)
472 
476 
477  template <int DIM = -1> MoFEMErrorCode setSkeletonAdjacency();
478 };
479 
480 } // namespace MoFEM
481 
482 #endif // __SIMPLE_HPP__
483 
484 /**
485  * \defgroup mofem_simple_interface Simple interface
486  * \brief Implementation of simple interface for fast problem set-up.
487  *
488  * \ingroup mofem
489  **/
MoFEM interface.
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:110
@ MF_ZERO
Definition: definitions.h:111
FieldApproximationBase
approximation base
Definition: definitions.h:71
FieldSpace
approximation spaces
Definition: definitions.h:95
const int dim
char mesh_file_name[255]
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
int FieldCoefficientsNumber
Number of field coefficients.
Definition: Types.hpp:38
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
Core (interface) class.
Definition: Core.hpp:92
Simple interface for fast problem set-up.
Definition: Simple.hpp:33
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:444
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:654
EntityHandle & getMeshSet()
Get the MeshSet object.
Definition: Simple.hpp:284
std::string & getSkeletonFEName()
Get the Skeleton FE Name.
Definition: Simple.hpp:361
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
Definition: Simple.hpp:326
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:453
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:450
EntityHandle & getSkeletonMeshSet()
Get the SkeletonMeshSet object.
Definition: Simple.hpp:298
MoFEMErrorCode removeDomainField(const std::string &name)
Remove field form domain.
Definition: Simple.cpp:366
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:412
MoFEM::Core & cOre
Definition: Simple.hpp:431
char meshFileName[255]
mesh file name
Definition: Simple.hpp:467
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:443
bool & getAddBoundaryFE()
Get the Add Boundary F E object.
Definition: Simple.hpp:414
int getDim() const
Get the problem dimension.
Definition: Simple.hpp:270
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:465
MoFEMErrorCode setSkeletonAdjacency()
Definition: Simple.cpp:28
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:628
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:349
MoFEMErrorCode exchangeGhostCells()
Definition: Simple.cpp:840
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:458
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:438
const std::string getSkeletonFEName() const
Get the Skeleton FE Name.
Definition: Simple.hpp:333
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:437
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:293
BitRefLevel & getBitRefLevel()
Get the BitRefLevel.
Definition: Simple.hpp:305
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:454
bool addSkeletonFE
Add skeleton FE.
Definition: Simple.hpp:447
bool addBoundaryFE
Add boundary FE.
Definition: Simple.hpp:448
MoFEMErrorCode addFieldToEmptyFieldBlocks(const std::string row_field, const std::string col_field) const
add empty block to problem
Definition: Simple.cpp:750
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:463
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:455
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:462
std::string nameOfProblem
problem name
Definition: Simple.hpp:460
SmartPetscObj< DM > getDM()
Return smart DM object.
Definition: Simple.hpp:260
void setDim(int dim)
Set the problem dimension.
Definition: Simple.hpp:277
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:216
virtual ~Simple()=default
std::string domainFE
domain finite element
Definition: Simple.hpp:461
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:451
MoFEMErrorCode reSetUp()
Rebuild internal MoFEM data structures.
Definition: Simple.cpp:686
std::string & getDomainFEName()
Get the Domain FE Name.
Definition: Simple.hpp:347
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:471
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:516
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:230
BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:433
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:508
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:311
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:484
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:452
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:329
MoFEMErrorCode removeSkeletonField(const std::string &name)
Remove field form skeleton.
Definition: Simple.cpp:397
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: Simple.cpp:22
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:668
EntityHandle & getBoundaryMeshSet()
Get the BoundaryMeshSet object.
Definition: Simple.hpp:291
MoFEMErrorCode removeBoundaryField(const std::string &name)
Remove field form boundary.
Definition: Simple.cpp:382
std::vector< std::string > & getOtherFiniteElements()
Get the Other Finite Elements.
Definition: Simple.hpp:378
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:340
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:319
int dIm
dimension of problem
Definition: Simple.hpp:468
MoFEMErrorCode deleteDM()
Delete dm.
Definition: Simple.cpp:727
PetscLogEvent MOFEM_EVENT_SimpleSetUP
Definition: Simple.hpp:436
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:441
MoFEMErrorCode createBoundaryMeshset()
Definition: Simple.cpp:759
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:439
BitRefLevel bitLevelMask
BitRefLevel of the probelm.
Definition: Simple.hpp:434
std::string & getBoundaryFEName()
Get the Boundary FE Name.
Definition: Simple.hpp:354
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:440
MoFEMErrorCode deleteFiniteElements()
Delete finite elements.
Definition: Simple.cpp:738
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:200
BitRefLevel & getBitRefLevelMask()
Get the BitRefLevel.
Definition: Simple.hpp:312
EntityHandle skeletonMeshset
skeleton meshset with boundary
Definition: Simple.hpp:445
MoFEMErrorCode createSkeletonMeshset()
Definition: Simple.cpp:806
std::string & getProblemName()
Get the Problem Name.
Definition: Simple.hpp:368
bool & getAddSkeletonFE()
Get the Add Skeleton F E object.
Definition: Simple.hpp:403
base class for all interface classes