v0.9.2
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  */
57  MoFEMErrorCode loadFile(const std::string options,
58  const std::string mesh_file_name);
59  /**
60  * \brief Load mesh file with parallel options if number of cores > 1
61  * @param mesh_file_name file name if not set default or set by command line
62  * is used.
63  * @return error code
64  */
65  MoFEMErrorCode loadFile(const std::string mesh_file_name = "");
66  /**
67  * \brief Add field on domain
68  * @param name name of the filed
69  * @param space space (L2,H1,Hdiv,Hcurl)
70  * @param base approximation base, see FieldApproximationBase
71  * @param nb_of_coefficients number of field coefficients
72  * @param tag_type type of the tag MB_TAG_DENSE or MB_TAG_SPARSE
73  * (DENSE is faster and uses less memory, SPARSE is more flexible if you
74  * define field on subdomains)
75  * @param bh if MF_EXCL throws error if field exits, MF_ZERO
76  * no error if field exist
77  * @param verb verbosity level
78  * @return error code
79  */
81  addDomainField(const std::string &name, const FieldSpace space,
82  const FieldApproximationBase base,
83  const FieldCoefficientsNumber nb_of_coefficients,
84  const TagType tag_type = MB_TAG_SPARSE,
85  const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
86 
87  /**
88  * \brief Add field on boundary
89  * @param name name of the filed
90  * @param space space (L2,H1,Hdiv,Hcurl)
91  * @param base approximation base, see FieldApproximationBase
92  * @param nb_of_coefficients number of field coefficients
93  * @param tag_type type of the tag MB_TAG_DENSE or MB_TAG_SPARSE
94  * (DENSE is faster and uses less memory, SPARSE is more flexible if you
95  * define field on subdomains)
96  * @param bh if MF_EXCL throws error if field exits, MF_ZERO
97  * no error if field exist
98  * @param verb verbosity level
99  * @return error code
100  */
102  addBoundaryField(const std::string &name, const FieldSpace space,
103  const FieldApproximationBase base,
104  const FieldCoefficientsNumber nb_of_coefficients,
105  const TagType tag_type = MB_TAG_SPARSE,
106  const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
107 
108  /**
109  * \brief Add field on skeleton
110  * @param name name of the filed
111  * @param space space (L2,H1,Hdiv,Hcurl)
112  * @param base approximation base, see FieldApproximationBase
113  * @param nb_of_coefficients number of field coefficients
114  * @param tag_type type of the tag MB_TAG_DENSE or MB_TAG_SPARSE
115  * (DENSE is faster and uses less memory, SPARSE is more flexible if you
116  * define field on subdomains)
117  * @param bh if MF_EXCL throws error if field exits, MF_ZERO
118  * no error if field exist
119  * @param verb verbosity level
120  * @return error code
121  */
123  addSkeletonField(const std::string &name, const FieldSpace space,
124  const FieldApproximationBase base,
125  const FieldCoefficientsNumber nb_of_coefficients,
126  const TagType tag_type = MB_TAG_SPARSE,
127  const enum MoFEMTypes bh = MF_ZERO, int verb = -1);
128 
129  /**
130  * \brief Add field on domain
131  * @param name name of the filed
132  * @param space space (L2,H1,Hdiv,Hcurl)
133  * @param base approximation base, see FieldApproximationBase
134  * @param nb_of_coefficients number of field coefficients
135  * @param tag_type type of the tag MB_TAG_DENSE or MB_TAG_SPARSE
136  * (DENSE is faster and uses less memory, SPARSE is more flexible if you
137  * define field on subdomains)
138  * @param bh if MF_EXCL throws error if field exits, MF_ZERO
139  * no error if field exist
140  * @param verb verbosity level
141  * @return error code
142  */
143  MoFEMErrorCode addDataField(const std::string &name, const FieldSpace space,
144  const FieldApproximationBase base,
145  const FieldCoefficientsNumber nb_of_coefficients,
146  const TagType tag_type = MB_TAG_SPARSE,
147  const enum MoFEMTypes bh = MF_ZERO,
148  int verb = -1);
149 
150 
151  /**
152  * @brief Remove field form domain
153  *
154  * @param name
155  * @return MoFEMErrorCode
156  */
157  MoFEMErrorCode removeDomainField(const std::string &name);
158 
159  /**
160  * @brief Remove field form boundary
161  *
162  * @param name
163  * @return MoFEMErrorCode
164  */
165  MoFEMErrorCode removeBoundaryField(const std::string &name);
166 
167  /**
168  * @brief Remove field form skeleton
169  *
170  * @param name
171  * @return MoFEMErrorCode
172  */
173  MoFEMErrorCode removeSkeletonField(const std::string &name);
174 
175  /**
176  * \brief Define finite elements
177  * @return Error code
178  */
180 
181  /**
182  * \brief define problem
183  * @return error code
184  */
185  MoFEMErrorCode defineProblem(const PetscBool is_partitioned = PETSC_TRUE);
186 
187  /**
188  * \brief Set field order
189  * @param std::field_name field name
190  * @param order order
191  * @param range of entities to which order is set (If null it sat to all
192  * entities)
193  * @return error code
194  */
195  MoFEMErrorCode setFieldOrder(const std::string field_name, const int order,
196  const Range *ents = NULL);
197 
198  /**
199  * \brief Build fields
200  * @return error code
201  */
203 
204  /**
205  * \brief Build finite elements
206  * @return error code
207  */
209 
210  /**
211  * @brief Set the skeleton adjacency object
212  *
213  * @param dim
214  * @return MoFEMErrorCode
215  */
217 
218  /**
219  * \brief Build problem
220  * @return error code
221  */
223 
224  /**
225  * \brief Setup problem
226  * @return error code
227  */
228  MoFEMErrorCode setUp(const PetscBool is_partitioned = PETSC_TRUE);
229 
230  /**
231  * @brief Rebuild internal MoFEM data structures
232  *
233  * Call this function after you add field or remove it.
234  *
235  * \note If you add field, or remove it, finite element and problem needs to
236  * be rebuild. However DM can remain the same.
237  *
238  * @return MoFEMErrorCode
239  */
241 
242  /**
243  * \brief Get DM
244  * @param dm discrete manager
245  * @return error code
246  */
247  MoFEMErrorCode getDM(DM *dm);
248 
249  /**
250  * @brief Return smart DM object
251  *
252  * \code
253  * {
254  * auto dm = simple_interface->getDM();
255  *
256  * // ...
257  *
258  * // dm is automatically destroyed when out of the scope
259  * }
260  * \endcode
261  *
262  * @return SmartPetscObj<DM>
263  */
264  inline SmartPetscObj<DM> getDM() { return dM; }
265 
266  /**
267  * @brief Get the problem dimensio
268  *
269  * Problem dimension is determined by highest dimension of entities on the
270  * mesh.
271  *
272  * @return int
273  */
274  inline int getDim() const { return dIm; }
275 
276  /**
277  * @brief Get the Domain FE Name
278  *
279  * @return const std::string
280  */
281  inline const std::string getDomainFEName() const { return domainFE; }
282 
283  /**
284  * @brief Get the Boundary FE Name
285  *
286  * @return const std::string
287  */
288  inline const std::string getBoundaryFEName() const { return boundaryFE; }
289 
290  /**
291  * @brief Get the Skeleton FE Name
292  *
293  * @return const std::string
294  */
295  inline const std::string getSkeletonFEName() const { return skeletonFE; }
296 
297  /**
298  * @brief Get the Problem Name
299  *
300  * @return const std::string
301  */
302  inline const std::string getProblemName() const { return nameOfProblem; }
303 
304  /**
305  * @brief Get the Domain FE Name
306  *
307  * @return std::string&
308  */
309  inline std::string &getDomainFEName() { return domainFE; }
310 
311  /**
312  * @brief Get the Boundary FE Name
313  *
314  * @return std::string&
315  */
316  inline std::string &getBoundaryFEName() { return boundaryFE; }
317 
318  /**
319  * @brief Get the Skeleton FE Name
320  *
321  * @return std::string&
322  */
323  inline std::string &getSkeletonFEName() { return skeletonFE; }
324 
325  /**
326  * @brief Get the Problem Name
327  *
328  * @return std::string&
329  */
330  inline std::string &getProblemName() { return nameOfProblem; }
331 
332  /**
333  * @brief Get the Other Finite Elements
334  *
335  * User can create finite elements using directly core interface and
336  * and them to simple problem by this function
337  *
338  * @return std::vector<std::string>&
339  */
340  inline std::vector<std::string> &getOtherFiniteElements() { return otherFEs; }
341 
342 private:
344 
345  const BitRefLevel bitLevel; ///< BitRefLevel of the probelm
346 
352 
353  EntityHandle meshSet; ///< domain meshset
354  EntityHandle boundaryMeshset; ///< meshset with boundary
355  std::vector<std::string> domainFields; ///< domain fields
356  std::vector<std::string> boundaryFields; ///< boundary fields
357  std::vector<std::string> skeletonFields; ///< fields on the skeleton
358  std::vector<std::string> dataFields; ///< Data fields
359  std::vector<std::string> noFieldFields; ///< NOFIELD field name
360  std::vector<std::string> noFieldDataFields; ///< NOFIELD field name
361 
362  std::multimap<std::string, std::pair<int, Range>> fieldsOrder; ///< fields order
363 
364  std::string nameOfProblem; ///< problem name
365  std::string domainFE; ///< domain finite element
366  std::string boundaryFE; ///< boundary finite element
367  std::string skeletonFE; ///< skeleton finite element
368 
369  std::vector<std::string> otherFEs; ///< Other finite elements
370 
371  char meshFileName[255]; ///< mesh file name
372  int dIm; ///< dimension of problem
373 
375  dM; ///< Discrete manager (interface to PETSc/MoFEM functions)
376 
377  template<int DIM = -1>
379 
380 };
381 
382 } // namespace MoFEM
383 
384 #endif // __SIMPLE_HPP__
385 
386 /**
387  * \defgroup mofem_simple_interface Simple interface
388  * \brief Implementation of simple interface for fast problem set-up.
389  *
390  * \ingroup mofem
391  **/
std::string boundaryFE
boundary finite element
Definition: Simple.hpp:366
std::string domainFE
domain finite element
Definition: Simple.hpp:365
MoFEM::Core & cOre
Definition: Simple.hpp:343
MoFEMErrorCode buildFiniteElements()
Build finite elements.
Definition: Simple.cpp:632
EntityHandle meshSet
domain meshset
Definition: Simple.hpp:353
MoFEMErrorCode setUp(const PetscBool is_partitioned=PETSC_TRUE)
Setup problem.
Definition: Simple.cpp:672
int getDim() const
Get the problem dimensio.
Definition: Simple.hpp:274
MoFEM interface unique ID.
const std::string getBoundaryFEName() const
Get the Boundary FE Name.
Definition: Simple.hpp:288
MoFEMErrorCode defineFiniteElements()
Define finite elements.
Definition: Simple.cpp:387
char meshFileName[255]
mesh file name
Definition: Simple.hpp:371
int dIm
dimension of problem
Definition: Simple.hpp:372
SmartPetscObj< DM > dM
Discrete manager (interface to PETSc/MoFEM functions)
Definition: Simple.hpp:375
base class for all interface classes
const std::string getSkeletonFEName() const
Get the Skeleton FE Name.
Definition: Simple.hpp:295
std::vector< std::string > boundaryFields
boundary fields
Definition: Simple.hpp:356
Core (interface) class.
Definition: Core.hpp:50
MoFEM interface.
std::vector< std::string > noFieldFields
NOFIELD field name.
Definition: Simple.hpp:359
std::string & getSkeletonFEName()
Get the Skeleton FE Name.
Definition: Simple.hpp:323
std::vector< std::string > otherFEs
Other finite elements.
Definition: Simple.hpp:369
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:55
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MoFEMErrorCode reSetUp()
Rebuild internal MoFEM data structures.
Definition: Simple.cpp:684
char mesh_file_name[255]
std::string skeletonFE
skeleton finite element
Definition: Simple.hpp:367
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:372
const int dim
MoFEMErrorCode removeBoundaryField(const std::string &name)
Remove field form boundary.
Definition: Simple.cpp:357
std::vector< std::string > domainFields
domain fields
Definition: Simple.hpp:355
std::vector< std::string > noFieldDataFields
NOFIELD field name.
Definition: Simple.hpp:360
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:324
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:287
MoFEMErrorCode buildProblem()
Build problem.
Definition: Simple.cpp:658
std::multimap< std::string, std::pair< int, Range > > fieldsOrder
fields order
Definition: Simple.hpp:362
FieldApproximationBase
approximation base
Definition: definitions.h:150
std::vector< std::string > dataFields
Data fields.
Definition: Simple.hpp:358
Simple(const MoFEM::Core &core)
Definition: Simple.cpp:185
MoFEMErrorCode loadFile(const std::string options, const std::string mesh_file_name)
Load mesh file.
Definition: Simple.cpp:212
MoFEMErrorCode buildFields()
Build fields.
Definition: Simple.cpp:490
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
std::string & getDomainFEName()
Get the Domain FE Name.
Definition: Simple.hpp:309
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:269
PetscLogEvent MOFEM_EVENT_SimpleLoadMesh
Definition: Simple.hpp:347
PetscLogEvent MOFEM_EVENT_SimpleBuildFields
Definition: Simple.hpp:348
std::string & getProblemName()
Get the Problem Name.
Definition: Simple.hpp:330
FieldSpace
approximation spaces
Definition: definitions.h:174
std::string nameOfProblem
problem name
Definition: Simple.hpp:364
const std::string getProblemName() const
Get the Problem Name.
Definition: Simple.hpp:302
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:305
const std::string getDomainFEName() const
Get the Domain FE Name.
Definition: Simple.hpp:281
constexpr int order
std::string & getBoundaryFEName()
Get the Boundary FE Name.
Definition: Simple.hpp:316
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:50
EntityHandle boundaryMeshset
meshset with boundary
Definition: Simple.hpp:354
PetscLogEvent MOFEM_EVENT_SimpleBuildFiniteElements
Definition: Simple.hpp:349
MoFEMErrorCode removeDomainField(const std::string &name)
Remove field form domain.
Definition: Simple.cpp:341
SmartPetscObj< DM > getDM()
Return smart DM object.
Definition: Simple.hpp:264
PetscLogEvent MOFEM_EVENT_SimpleBuildProblem
Definition: Simple.hpp:350
std::vector< std::string > skeletonFields
fields on the skeleton
Definition: Simple.hpp:357
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:340
MoFEMErrorCode setFieldOrder(const std::string field_name, const int order, const Range *ents=NULL)
Set field order.
Definition: Simple.cpp:482
const BitRefLevel bitLevel
BitRefLevel of the probelm.
Definition: Simple.hpp:345
MoFEMErrorCode defineProblem(const PetscBool is_partitioned=PETSC_TRUE)
define problem
Definition: Simple.cpp:459
MoFEMErrorCode setSkeletonAdjacency()
Definition: Simple.cpp:35
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:189
PetscLogEvent MOFEM_EVENT_SimpleKSPSolve
Definition: Simple.hpp:351
virtual ~Simple()=default