libmesh_common.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_LIBMESH_COMMON_H
21 #define LIBMESH_LIBMESH_COMMON_H
22 
23 
24 // The library configuration options
25 #include "libmesh/libmesh_config.h"
26 
27 // Include the MPI definition
28 #ifdef LIBMESH_HAVE_MPI
29 # include "libmesh/ignore_warnings.h"
30 # include <mpi.h>
32 #endif
33 
34 // C/C++ includes everyone should know about
35 #include <complex>
36 #ifdef LIBMESH_HAVE_STDLIB_H
37 # include <cstdlib>
38 #endif
39 #include <typeinfo> // std::bad_cast
40 
41 // _basic_ library functionality
42 #include "libmesh/libmesh_base.h"
44 extern "C" {
46 }
47 
48 // Proxy class for libMesh::out/err output
49 #include "libmesh/ostream_proxy.h"
50 
51 // Here we add missing types to the standard namespace. For example,
52 // std::max(double, float) etc... are well behaved but not defined
53 // by the standard. This also includes workarounds for super-strict
54 // implementations, for example Sun Studio and PGI C++. However,
55 // this necessarily requires breaking the ISO-C++ standard, and is
56 // really just a hack. As such, only do it if we are building the
57 // libmesh library itself. Specifially, *DO NOT* export this to
58 // user code or install this header.
59 #ifdef LIBMESH_IS_COMPILING_ITSELF
61 #endif
62 
63 
64 namespace libMesh
65 {
66 
67 
68 // A namespace for functions used in the bodies of the macros below.
69 // The macros generally call these functions with __FILE__, __LINE__,
70 // __DATE__, and __TIME__ in the appropriate order. These should not
71 // be called by users directly! The implementations can be found in
72 // libmesh_common.C.
73 namespace MacroFunctions
74 {
75  void here(const char* file, int line, const char* date, const char* time);
76  void stop(const char* file, int line, const char* date, const char* time);
77  void report_error(const char* file, int line, const char* date, const char* time);
78 }
79 
80 // Undefine any existing macros
81 #ifdef Real
82 # undef Real
83 #endif
84 
85 //#ifdef REAL
86 //# undef REAL
87 //#endif
88 
89 #ifdef Complex
90 # undef Complex
91 #endif
92 
93 #ifdef COMPLEX
94 # undef COMPLEX
95 #endif
96 
97 #ifdef MPI_REAL
98 # undef MPI_REAL
99 #endif
100 
101 // Check to see if TOLERANCE has been defined by another
102 // package, if so we might want to change the name...
103 #ifdef TOLERANCE
104  DIE A HORRIBLE DEATH HERE...
105 # undef TOLERANCE
106 #endif
107 
108 
109 
110 // Define the type to use for real numbers
111 
112 typedef LIBMESH_DEFAULT_SCALAR_TYPE Real;
113 
114 // Define a corresponding tolerance. This is what should be
115 // considered "good enough" when doing floating point comparisons.
116 // For example, v == 0 is changed to std::abs(v) < TOLERANCE.
117 
118 #ifndef LIBMESH_DEFAULT_SINGLE_PRECISION
119  #ifdef LIBMESH_DEFAULT_TRIPLE_PRECISION
120 static const Real TOLERANCE = 1.e-8;
121  # define MPI_REAL MPI_LONG_DOUBLE
122  #else
123 static const Real TOLERANCE = 1.e-6;
124  # define MPI_REAL MPI_DOUBLE
125  #endif
126 #else
127 static const Real TOLERANCE = 2.5e-3;
128  # define MPI_REAL MPI_FLOAT
129 #endif
130 
131 // Define the type to use for complex numbers
132 // Always use std::complex<double>, as required by Petsc?
133 // If your version of Petsc doesn't support
134 // std::complex<other_precision>, then you'd better just leave
135 // Real==double
136 typedef std::complex<Real> Complex;
137 typedef std::complex<Real> COMPLEX;
138 
139 
140 // Helper functions for complex/real numbers
141 // to clean up #ifdef LIBMESH_USE_COMPLEX_NUMBERS elsewhere
142 template<typename T> inline T libmesh_real(T a) { return a; }
143 template<typename T> inline T libmesh_conj(T a) { return a; }
144 
145 template<typename T>
146 inline T libmesh_real(std::complex<T> a) { return std::real(a); }
147 
148 template<typename T>
149 inline std::complex<T> libmesh_conj(std::complex<T> a) { return std::conj(a); }
150 
151 // isnan isn't actually C++ standard yet; in contexts where it's not defined in
152 // cmath, libmesh_isnan will just end up returning false.
153 inline bool libmesh_isnan(float a) { return libmesh_C_isnan_float(a); }
154 inline bool libmesh_isnan(double a) { return libmesh_C_isnan_double(a); }
155 inline bool libmesh_isnan(long double a) { return libmesh_C_isnan_longdouble(a); }
156 
157 template <typename T>
158 inline bool libmesh_isnan(std::complex<T> a) { return (libmesh_isnan(std::real(a)) || libmesh_isnan(std::imag(a))); }
159 
160 
161 // Define the value type for unknowns in simulations.
162 // This is either Real or Complex, depending on how
163 // the library was configures
164 #if defined (LIBMESH_USE_REAL_NUMBERS)
165  typedef Real Number;
166 #elif defined (LIBMESH_USE_COMPLEX_NUMBERS)
167  typedef Complex Number;
168 #else
169  DIE A HORRIBLE DEATH HERE...
170 #endif
171 
172 
173 // Define the value type for error estimates.
174 // Since AMR/C decisions don't have to be precise,
175 // we default to float for memory efficiency.
176 typedef float ErrorVectorReal;
177 #define MPI_ERRORVECTORREAL MPI_FLOAT
178 
179 
180 #ifdef LIBMESH_HAVE_MPI
181 
184  extern MPI_Comm COMM_WORLD;
185 #else
186 
190  extern int COMM_WORLD;
191 #endif
192 
193 // Let's define a couple output streams - these will default
194 // to cout/cerr, but LibMeshInit (or the user) can also set them to
195 // something more sophisticated.
196 //
197 // We use a proxy class rather than references so they can be
198 // reseated at runtime.
199 
202 
203 // These are useful macros that behave like functions in the code.
204 // If you want to make sure you are accessing a section of code just
205 // stick a libmesh_here(); in it, for example
206 #define libmesh_here() \
207  do { \
208  libMesh::MacroFunctions::here(__FILE__, __LINE__, __DATE__, __TIME__); \
209  } while(0)
210 
211 // the libmesh_stop() macro will stop the code until a SIGCONT signal
212 // is recieved. This is useful, for example, when determining the
213 // memory used by a given operation. A libmesh_stop() could be
214 // instered before and after a questionable operation and the delta
215 // memory can be obtained from a ps or top. This macro only works for
216 // serial cases.
217 #define libmesh_stop() \
218  do { \
219  libMesh::MacroFunctions::stop(__FILE__, __LINE__, __DATE__, __TIME__); \
220  } while(0)
221 
222 // The libmesh_dbg_var() macro indicates that an argument to a function
223 // is used only in debug mode (i.e., when NDEBUG is not defined).
224 #ifndef NDEBUG
225 #define libmesh_dbg_var(var) var
226 #else
227 #define libmesh_dbg_var(var)
228 #endif
229 
230 // The libmesh_dbg_var() macro indicates that an argument to a function
231 // is used only in debug mode (i.e., when NDEBUG is not defined).
232 #ifndef NDEBUG
233 #define libmesh_dbg_var(var) var
234 #else
235 #define libmesh_dbg_var(var)
236 #endif
237 
238 // The libmesh_assert() macro acts like C's assert(), but throws a
239 // libmesh_error() (including stack trace, etc) instead of just exiting
240 #ifdef NDEBUG
241 
242 #define libmesh_assert(asserted) ((void) 0)
243 #define libmesh_assert_msg(asserted, msg) ((void) 0)
244 #define libmesh_assert_equal_to(expr1,expr2) ((void) 0)
245 #define libmesh_assert_not_equal_to(expr1,expr2) ((void) 0)
246 #define libmesh_assert_less(expr1,expr2) ((void) 0)
247 #define libmesh_assert_greater(expr1,expr2) ((void) 0)
248 #define libmesh_assert_less_equal(expr1,expr2) ((void) 0)
249 #define libmesh_assert_greater_equal(expr1,expr2) ((void) 0)
250 
251 #else
252 
253 #define libmesh_assert(asserted) \
254  do { \
255  if (!(asserted)) { \
256  libMesh::err << "Assertion `" #asserted "' failed." << std::endl; \
257  libmesh_error(); \
258  } } while(0)
259 
260 #define libmesh_assert_msg(asserted, msg) \
261  do { \
262  if (!(asserted)) { \
263  libMesh::err << "Assertion `" #asserted "' failed." << std::endl; \
264  libmesh_error_msg(msg); \
265  } } while(0)
266 
267 #define libmesh_assert_equal_to(expr1,expr2) \
268  do { \
269  if (!(expr1 == expr2)) { \
270  libMesh::err << "Assertion `" #expr1 " == " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << std::endl; \
271  libmesh_error(); \
272  } } while(0)
273 
274 #define libmesh_assert_not_equal_to(expr1,expr2) \
275  do { \
276  if (!(expr1 != expr2)) { \
277  libMesh::err << "Assertion `" #expr1 " != " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << std::endl; \
278  libmesh_error(); \
279  } } while(0)
280 
281 #define libmesh_assert_less(expr1,expr2) \
282  do { \
283  if (!(expr1 < expr2)) { \
284  libMesh::err << "Assertion `" #expr1 " < " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << std::endl; \
285  libmesh_error(); \
286  } } while(0)
287 
288 #define libmesh_assert_greater(expr1,expr2) \
289  do { \
290  if (!(expr1 > expr2)) { \
291  libMesh::err << "Assertion `" #expr1 " > " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << std::endl; \
292  libmesh_error(); \
293  } } while(0)
294 
295 #define libmesh_assert_less_equal(expr1,expr2) \
296  do { \
297  if (!(expr1 <= expr2)) { \
298  libMesh::err << "Assertion `" #expr1 " <= " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << std::endl; \
299  libmesh_error(); \
300  } } while(0)
301 
302 #define libmesh_assert_greater_equal(expr1,expr2) \
303  do { \
304  if (!(expr1 >= expr2)) { \
305  libMesh::err << "Assertion `" #expr1 " >= " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << std::endl; \
306  libmesh_error(); \
307  } } while(0)
308 
309 #endif
310 
311 // The libmesh_error() macro prints a message and throws a LogicError
312 // exception
313 //
314 // The libmesh_not_implemented() macro prints a message and throws a
315 // NotImplemented exception
316 //
317 // The libmesh_file_error(const std::string& filename) macro prints a message
318 // and throws a FileError exception
319 //
320 // The libmesh_convergence_failure() macro
321 // throws a ConvergenceFailure exception
322 #define libmesh_error() \
323  do { \
324  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, __DATE__, __TIME__); \
325  LIBMESH_THROW(libMesh::LogicError()); \
326  } while(0)
327 
328 #define libmesh_error_msg(msg) \
329  do { \
330  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, __DATE__, __TIME__); \
331  libMesh::err << msg << std::endl; \
332  LIBMESH_THROW(libMesh::LogicError()); \
333  } while(0)
334 
335 #define libmesh_not_implemented() \
336  do { \
337  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, __DATE__, __TIME__); \
338  LIBMESH_THROW(libMesh::NotImplemented()); \
339  } while(0)
340 
341 #define libmesh_not_implemented_msg(msg) \
342  do { \
343  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, __DATE__, __TIME__); \
344  libMesh::err << msg << std::endl; \
345  LIBMESH_THROW(libMesh::NotImplemented()); \
346  } while(0)
347 
348 #define libmesh_file_error(filename) \
349  do { \
350  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, __DATE__, __TIME__); \
351  LIBMESH_THROW(libMesh::FileError(filename)); \
352  } while(0)
353 
354 #define libmesh_file_error_msg(filename, msg) \
355  do { \
356  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, __DATE__, __TIME__); \
357  libMesh:err << msg << std::endl; \
358  LIBMESH_THROW(libMesh::FileError(filename)); \
359  } while(0)
360 
361 #define libmesh_convergence_failure() \
362  do { \
363  LIBMESH_THROW(libMesh::ConvergenceFailure()); \
364  } while(0)
365 
366 // The libmesh_example_assert() macro prints a message and calls
367 // "return 0;" if the assertion specified by the macro is not true. This
368 // macro is used in the example executables, which should run when the
369 // configure-time libMesh options support them but which should exit
370 // without failure otherwise.
371 //
372 // This macro only works in main(), because we have no better way than
373 // "return" from main to immediately exit successfully - std::exit()
374 // gets seen by at least some MPI stacks as failure.
375 //
376 // We now return 77, the automake code for a skipped test.
377 
378 #define libmesh_example_assert(asserted, requirement) \
379  do { \
380  if (!(asserted)) { \
381  libMesh::out << "Assertion `" #asserted "' failed. Configuring libMesh with " requirement " may be required to run this code." << std::endl; \
382  return 77; \
383  } } while(0)
384 
385 // libmesh_cast_ref and libmesh_cast_ptr do a dynamic cast and assert
386 // the result, if we have RTTI enabled and we're in debug or
387 // development modes, but they just do a faster static cast if we're
388 // in optimized mode.
389 //
390 // Use these casts when you're certain that a cast will succeed in
391 // correct code but you want to be able to double-check.
392 template <typename Tnew, typename Told>
393 inline Tnew libmesh_cast_ref(Told& oldvar)
394 {
395 #if !defined(NDEBUG) && defined(LIBMESH_HAVE_RTTI)
396  try
397  {
398  Tnew newvar = dynamic_cast<Tnew>(oldvar);
399  return newvar;
400  }
401  catch (std::bad_cast)
402  {
403  libMesh::err << "Failed to convert " << typeid(Told).name()
404  << " reference to " << typeid(Tnew).name()
405  << std::endl;
406  libMesh::err << "The " << typeid(Told).name()
407  << " appears to be a "
408  << typeid(*(&oldvar)).name() << std::endl;
409  libmesh_error();
410  }
411 #else
412  return(static_cast<Tnew>(oldvar));
413 #endif
414 }
415 
416 // We use two different function names to avoid an odd overloading
417 // ambiguity bug with icc 10.1.008
418 template <typename Tnew, typename Told>
419 inline Tnew libmesh_cast_ptr (Told* oldvar)
420 {
421 #if !defined(NDEBUG) && defined(LIBMESH_HAVE_RTTI)
422  Tnew newvar = dynamic_cast<Tnew>(oldvar);
423  if (!newvar)
424  {
425  libMesh::err << "Failed to convert " << typeid(Told).name()
426  << " pointer to " << typeid(Tnew).name()
427  << std::endl;
428  libMesh::err << "The " << typeid(Told).name()
429  << " appears to be a "
430  << typeid(*oldvar).name() << std::endl;
431  libmesh_error();
432  }
433  return newvar;
434 #else
435  return(static_cast<Tnew>(oldvar));
436 #endif
437 }
438 
439 
440 // libmesh_cast_int asserts that the value of the castee is within the
441 // bounds which are exactly representable by the output type, if we're
442 // in debug or development modes, but it just does a faster static
443 // cast if we're in optimized mode.
444 //
445 // Use these casts when you're certain that a cast will succeed in
446 // correct code but you want to be able to double-check.
447 template <typename Tnew, typename Told>
448 inline Tnew libmesh_cast_int (Told oldvar)
449 {
450  libmesh_assert_equal_to
451  (oldvar, static_cast<Told>(static_cast<Tnew>(oldvar)));
452 
453  return(static_cast<Tnew>(oldvar));
454 }
455 
456 
457 // The libmesh_do_once macro helps us avoid redundant repeated
458 // repetitions of the same warning messages
459 #undef libmesh_do_once
460 #define libmesh_do_once(do_this) \
461  do { \
462  static bool did_this_already = false; \
463  if (!did_this_already) { \
464  did_this_already = true; \
465  do_this; \
466  } } while (0)
467 
468 
469 // The libmesh_experimental macro warns that you are using
470 // bleeding-edge code
471 #undef libmesh_experimental
472 #ifdef LIBMESH_ENABLE_WARNINGS
473 #define libmesh_experimental() \
474  libmesh_do_once(libMesh::out << "*** Warning, This code is untested, experimental, or likely to see future API changes: " \
475  << __FILE__ << ", line " << __LINE__ << ", compiled " << __DATE__ << " at " << __TIME__ << " ***" << std::endl;)
476 #else
477 #define libmesh_experimental() ((void) 0)
478 #endif
479 
480 
481 // The libmesh_deprecated macro warns that you are using obsoleted code
482 #undef libmesh_deprecated
483 #ifdef LIBMESH_ENABLE_WARNINGS
484 #define libmesh_deprecated() \
485  libmesh_do_once(libMesh::out << "*** Warning, This code is deprecated, and likely to be removed in future library versions! " \
486  << __FILE__ << ", line " << __LINE__ << ", compiled " << __DATE__ << " at " << __TIME__ << " ***" << std::endl;)
487 #else
488 #define libmesh_deprecated() ((void) 0)
489 #endif
490 
491 
492 // A function template for ignoring unused variables. This is a way
493 // to shut up unused variable compiler warnings on a case by case
494 // basis.
495 template<class T> inline void libmesh_ignore( const T& ) { }
496 
497 
498 // build a integer representation of version
499 #define LIBMESH_VERSION_ID(major,minor,patch) (((major) << 16) | ((minor) << 8) | ((patch) & 0xFF))
500 
501 } // namespace libMesh
502 
503 
504 #endif // LIBMESH_LIBMESH_COMMON_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:05 UTC

Hosted By:
SourceForge.net Logo