v0.9.1
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 static const MOFEMuuid IDD_MOFEMSimple =
31 
32 /**
33  * \brief Simple interface for fast problem set-up
34  * \ingroup mofem_simple_interface
35  */
36 struct Simple : public UnknownInterface {
37 
39  UnknownInterface **iface) const;
40 
41  Simple(const MoFEM::Core &core);
42  virtual ~Simple() = default;
43 
44  /**
45  * \brief get options
46  * @return error code
47  */
49 
50  /**
51  * \brief Load mesh file
52  * @param options file load options
53  * @param mesh_file_name file name if not set default or set by command line
54  * is used.
55  * @return error code
56  */
58  loadFile(const std::string options = "PARALLEL=READ_PART;"
59  "PARALLEL_RESOLVE_SHARED_ENTS;"
60  "PARTITION=PARALLEL_PARTITION;",
61  const std::string mesh_file_name = "");
62 
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  /**
149  * @brief Remove field form domain
150  *
151  * @param name
152  * @return MoFEMErrorCode
153  */
154  MoFEMErrorCode removeDomainField(const std::string &name);
155 
156  /**
157  * @brief Remove field form boundary
158  *
159  * @param name
160  * @return MoFEMErrorCode
161  */
162  MoFEMErrorCode removeBoundaryField(const std::string &name);
163 
164  /**
165  * @brief Remove field form skeleton
166  *
167  * @param name
168  * @return MoFEMErrorCode
169  */
170  MoFEMErrorCode removeSkeletonField(const std::string &name);
171 
172  /**
173  * \brief Define finite elements
174  * @return Error code
175  */
177 
178  /**
179  * \brief define problem
180  * @return error code
181  */
182  MoFEMErrorCode defineProblem(const PetscBool is_partitioned = PETSC_TRUE);
183 
184  /**
185  * \brief Set field order
186  * @param std::field_name field name
187  * @param order order
188  * @param range of entities to which order is set (If null it sat to all
189  * entities)
190  * @return error code
191  */
192  MoFEMErrorCode setFieldOrder(const std::string field_name, const int order,
193  const Range *ents = NULL);
194 
195  /**
196  * \brief Build fields
197  * @return error code
198  */
200 
201  /**
202  * \brief Build finite elements
203  * @return error code
204  */
206 
207  /**
208  * @brief Set the skeleton adjacency object
209  *
210  * @param dim
211  * @return MoFEMErrorCode
212  */
213  MoFEMErrorCode setSkeletonAdjacency(int dim = -1);
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  */
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 dimensio
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 Get the Domain FE Name
275  *
276  * @return const std::string
277  */
278  inline const std::string getDomainFEName() const { return domainFE; }
279 
280  /**
281  * @brief Get the Boundary FE Name
282  *
283  * @return const std::string
284  */
285  inline const std::string getBoundaryFEName() const { return boundaryFE; }
286 
287  /**
288  * @brief Get the Skeleton FE Name
289  *
290  * @return const std::string
291  */
292  inline const std::string getSkeletonFEName() const { return skeletonFE; }
293 
294  /**
295  * @brief Get the Problem Name
296  *
297  * @return const std::string
298  */
299  inline const std::string getProblemName() const { return nameOfProblem; }
300 
301  /**
302  * @brief Get the Domain FE Name
303  *
304  * @return std::string&
305  */
306  inline std::string &getDomainFEName() { return domainFE; }
307 
308  /**
309  * @brief Get the Boundary FE Name
310  *
311  * @return std::string&
312  */
313  inline std::string &getBoundaryFEName() { return boundaryFE; }
314 
315  /**
316  * @brief Get the Skeleton FE Name
317  *
318  * @return std::string&
319  */
320  inline std::string &getSkeletonFEName() { return skeletonFE; }
321 
322  /**
323  * @brief Get the Problem Name
324  *
325  * @return std::string&
326  */
327  inline std::string &getProblemName() { return nameOfProblem; }
328 
329  /**
330  * @brief Get the Other Finite Elements
331  *
332  * User can create finite elements using directly core interface and
333  * and them to simple problem by this function
334  *
335  * @return std::vector<std::string>&
336  */
337  inline std::vector<std::string> &getOtherFiniteElements() { return otherFEs; }
338 
339 private:
341 
342  const BitRefLevel bitLevel; ///< BitRefLevel of the probelm
343 
349 
350  EntityHandle meshSet; ///< domain meshset
351  EntityHandle boundaryMeshset; ///< meshset with boundary
352  std::vector<std::string> domainFields; ///< domain fields
353  std::vector<std::string> boundaryFields; ///< boundary fields
354  std::vector<std::string> skeletonFields; ///< fields on the skeleton
355  std::vector<std::string> dataFields; ///< Data fields
356  std::vector<std::string> noFieldFields; ///< NOFIELD field name
357  std::vector<std::string> noFieldDataFields; ///< NOFIELD field name
358 
359  std::multimap<std::string, std::pair<int, Range>> fieldsOrder; ///< fields order
360 
361  std::string nameOfProblem; ///< problem name
362  std::string domainFE; ///< domain finite element
363  std::string boundaryFE; ///< boundary finite element
364  std::string skeletonFE; ///< skeleton finite element
365 
366  std::vector<std::string> otherFEs; ///< Other finite elements
367 
368  char meshFileName[255]; ///< mesh file name
369  int dIm; ///< dimension of problem
370 
372  dM; ///< Discrete manager (interface to PETSc/MoFEM functions)
373 
374  template<int DIM = -1>
376 
377 };
378 
379 } // namespace MoFEM
380 
381 #endif // __SIMPLE_HPP__
382 
383 /**
384  * \defgroup mofem_simple_interface Simple interface
385  * \brief Implementation of simple interface for fast problem set-up.
386  *
387  * \ingroup mofem
388  **/
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:363
std::string domainFE
domain finite element
Definition: Simple.hpp:362
MoFEM::Core & cOre
Definition: Simple.hpp:340
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:620
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:350
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:660
int getDim() const
Get the problem dimensio.
Definition: Simple.hpp:271
MoFEM interface unique ID.
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
Definition: Simple.hpp:285
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:374
char meshFileName[255]
mesh file name
Definition: Simple.hpp:368
int dIm
dimension of problem
Definition: Simple.hpp:369
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:372
base class for all interface classes
const std::string getSkeletonFEName() const
Get the Skeleton FE Name.
Definition: Simple.hpp:292
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:353
Core (interface) class.
Definition: Core.hpp:50
MoFEM interface.
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:356
std::string & getSkeletonFEName()
Get the Skeleton FE Name.
Definition: Simple.hpp:320
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:366
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MoFEMErrorCode reSetUp()
Rebuild internal MoFEM data structures.
Definition: Simple.cpp:672
char mesh_file_name[255]
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:364
Simple interface for fast problem set-up.
Definition: Simple.hpp:36
static const MOFEMuuid IDD_MOFEMSimple
Definition: Simple.hpp:29
MoFEMErrorCode removeSkeletonField(const std::string &name)
Remove field form skeleton.
Definition: Simple.cpp:359
MoFEMErrorCode removeBoundaryField(const std::string &name)
Remove field form boundary.
Definition: Simple.cpp:344
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:352
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:357
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:311
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:274
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:646
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:359
FieldApproximationBase
approximation base
Definition: definitions.h:149
MoFEMErrorCode loadFile(const std::string options="PARALLEL=READ_PART;" "PARALLEL_RESOLVE_SHARED_ENTS;" "PARTITION=PARALLEL_PARTITION;", const std::string mesh_file_name="")
Load mesh file.
Definition: Simple.cpp:212
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:355
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:185
constexpr int order
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:478
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
std::string & getDomainFEName()
Get the Domain FE Name.
Definition: Simple.hpp:306
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Definition: Simple.cpp:22
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:256
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:344
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:345
std::string & getProblemName()
Get the Problem Name.
Definition: Simple.hpp:327
FieldSpace
approximation spaces
Definition: definitions.h:173
std::string nameOfProblem
problem name
Definition: Simple.hpp:361
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:299
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:292
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:278
std::string & getBoundaryFEName()
Get the Boundary FE Name.
Definition: Simple.hpp:313
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:351
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:346
MoFEMErrorCode removeDomainField(const std::string &name)
Remove field form domain.
Definition: Simple.cpp:328
SmartPetscObj< DM > getDM()
Return smart DM object.
Definition: Simple.hpp:261
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:347
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:354
int FieldCoefficientsNumber
Number of field coefficients.
Definition: Types.hpp:38
MoFEMErrorCode getOptions()
get options
Definition: Simple.cpp:198
std::vector< std::string > & getOtherFiniteElements()
Get the Other Finite Elements.
Definition: Simple.hpp:337
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:470
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:342
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:447
MoFEMErrorCode setSkeletonAdjacency()
Definition: Simple.cpp:35
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:188
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:348
virtual ~Simple()=default