1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2020 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 // The library configuration options
24 #include "libmesh/libmesh_config.h"
25
26 // Use actual timestamps or constant dummies (to aid ccache)
27 #ifdef LIBMESH_ENABLE_TIMESTAMPS
28 # define LIBMESH_TIME __TIME__
29 # define LIBMESH_DATE __DATE__
30 #else
31 # define LIBMESH_TIME "notime"
32 # define LIBMESH_DATE "nodate"
33 #endif
34
35 // C/C++ includes everyone should know about
36 #include <cstdlib>
37 #ifdef __PGI
38 // BSK, Thu Feb 20 08:32:06 CST 2014 - For some reason, unless PGI gets
39 // <cmath> early this nonsense shows up:
40 // "/software/x86_64/pgi/12.9/linux86-64/12.9/include/CC/cmath", line 57: error:
41 // the global scope has no "abs"
42 // using _STLP_VENDOR_CSTD::abs;
43 // So include <cmath> as early as possible under the PGI compilers.
44 # include <cmath>
45 #endif
46 #include <complex>
47 #include <typeinfo> // std::bad_cast
48 #include <type_traits> // std::decay
49 #include <functional> // std::less, etc
50
51 // Include the MPI definition
52 #ifdef LIBMESH_HAVE_MPI
53 # include "libmesh/ignore_warnings.h"
54 # include <mpi.h>
55 # include "libmesh/restore_warnings.h"
56 #endif
57
58 // Quad precision if we need it
59 #ifdef LIBMESH_DEFAULT_QUADRUPLE_PRECISION
60 #include "libmesh/float128_shims.h"
61 #endif
62
63 // _basic_ library functionality
64 #include "libmesh/libmesh_base.h"
65 #include "libmesh/libmesh_exceptions.h"
66
67 // Proxy class for libMesh::out/err output
68 #include "libmesh/ostream_proxy.h"
69
70 // Make sure the libmesh_nullptr define is available for backwards
71 // compatibility, although we no longer use it in the library.
72 #include "libmesh/libmesh_nullptr.h"
73
74 namespace libMesh
75 {
76
77 // A namespace for functions used in the bodies of the macros below.
78 // The macros generally call these functions with __FILE__, __LINE__,
79 // __DATE__, and __TIME__ in the appropriate order. These should not
80 // be called by users directly! The implementations can be found in
81 // libmesh_common.C.
82 namespace MacroFunctions
83 {
84 void here(const char * file, int line, const char * date, const char * time);
85 void stop(const char * file, int line, const char * date, const char * time);
86 void report_error(const char * file, int line, const char * date, const char * time);
87 }
88
89 // Undefine any existing macros
90 #ifdef Real
91 # undef Real
92 #endif
93
94 //#ifdef REAL
95 //# undef REAL
96 //#endif
97
98 #ifdef Complex
99 # undef Complex
100 #endif
101
102 #ifdef COMPLEX
103 # undef COMPLEX
104 #endif
105
106 // Check to see if TOLERANCE has been defined by another
107 // package, if so we might want to change the name...
108 #ifdef TOLERANCE
109 DIE A HORRIBLE DEATH HERE...
110 # undef TOLERANCE
111 #endif
112
113
114
115 // Define the type to use for real numbers
116
117 typedef LIBMESH_DEFAULT_SCALAR_TYPE Real;
118
119 // Define a corresponding tolerance. This is what should be
120 // considered "good enough" when doing floating point comparisons.
121 // For example, v == 0 is changed to std::abs(v) < TOLERANCE.
122
123 #ifdef LIBMESH_DEFAULT_SINGLE_PRECISION
124 static const Real TOLERANCE = 2.5e-3;
125 # if defined (LIBMESH_DEFAULT_TRIPLE_PRECISION) || \
126 defined (LIBMESH_DEFAULT_QUADRUPLE_PRECISION)
127 # error Cannot define multiple precision levels
128 # endif
129 #endif
130
131 #ifdef LIBMESH_DEFAULT_TRIPLE_PRECISION
132 static const Real TOLERANCE = 1.e-8;
133 # if defined (LIBMESH_DEFAULT_QUADRUPLE_PRECISION)
134 # error Cannot define multiple precision levels
135 # endif
136 #endif
137
138 #ifdef LIBMESH_DEFAULT_QUADRUPLE_PRECISION
139 static const Real TOLERANCE = 1.e-11;
140 #endif
141
142 #if !defined (LIBMESH_DEFAULT_SINGLE_PRECISION) && \
143 !defined (LIBMESH_DEFAULT_TRIPLE_PRECISION) && \
144 !defined (LIBMESH_DEFAULT_QUADRUPLE_PRECISION)
145 static const Real TOLERANCE = 1.e-6;
146 #endif
147
148 // Define the type to use for complex numbers
149 // Always use std::complex<double>, as required by Petsc?
150 // If your version of Petsc doesn't support
151 // std::complex<other_precision>, then you'd better just leave
152 // Real==double
153 typedef std::complex<Real> Complex;
154 typedef std::complex<Real> COMPLEX;
155
156
157 // Helper functions for complex/real numbers
158 // to clean up #ifdef LIBMESH_USE_COMPLEX_NUMBERS elsewhere
libmesh_real(T a)159 template<typename T> inline T libmesh_real(T a) { return a; }
libmesh_conj(T a)160 template<typename T> inline T libmesh_conj(T a) { return a; }
161
162 template<typename T>
libmesh_real(std::complex<T> a)163 inline T libmesh_real(std::complex<T> a) { return std::real(a); }
164
165 template<typename T>
libmesh_conj(std::complex<T> a)166 inline std::complex<T> libmesh_conj(std::complex<T> a) { return std::conj(a); }
167
168 // std::isnan() is in <cmath> as of C++11.
169 template <typename T>
libmesh_isnan(T x)170 inline bool libmesh_isnan(T x) { return std::isnan(x); }
171
172 template <typename T>
libmesh_isnan(std::complex<T> a)173 inline bool libmesh_isnan(std::complex<T> a)
174 { return (std::isnan(std::real(a)) || std::isnan(std::imag(a))); }
175
176 // std::isinf() is in <cmath> as of C++11.
177 template <typename T>
libmesh_isinf(T x)178 inline bool libmesh_isinf(T x) { return std::isinf(x); }
179
180 template <typename T>
libmesh_isinf(std::complex<T> a)181 inline bool libmesh_isinf(std::complex<T> a)
182 { return (std::isinf(std::real(a)) || std::isinf(std::imag(a))); }
183
184 // Define the value type for unknowns in simulations.
185 // This is either Real or Complex, depending on how
186 // the library was configures
187 #if defined (LIBMESH_USE_REAL_NUMBERS)
188 typedef Real Number;
189 #elif defined (LIBMESH_USE_COMPLEX_NUMBERS)
190 typedef Complex Number;
191 #else
192 DIE A HORRIBLE DEATH HERE...
193 #endif
194
195
196 // Define the value type for error estimates.
197 // Since AMR/C decisions don't have to be precise,
198 // we default to float for memory efficiency.
199 typedef float ErrorVectorReal;
200 #define MPI_ERRORVECTORREAL MPI_FLOAT
201
202
203 #ifdef LIBMESH_HAVE_MPI
204
205 /**
206 * MPI Communicator used to initialize libMesh.
207 */
208 extern MPI_Comm GLOBAL_COMM_WORLD;
209 #else
210
211 /**
212 * Something to use with CHKERRABORT if we're just using PETSc's MPI
213 * "uni" stub.
214 */
215 extern int GLOBAL_COMM_WORLD;
216 #endif
217
218 // Let's define a couple output streams - these will default
219 // to cout/cerr, but LibMeshInit (or the user) can also set them to
220 // something more sophisticated.
221 //
222 // We use a proxy class rather than references so they can be
223 // reseated at runtime.
224
225 extern OStreamProxy out;
226 extern OStreamProxy err;
227
228 // This global variable is to help us deprecate AutoPtr. We can't
229 // just use libmesh_deprecated() because then you get one print out
230 // per template instantiation, instead of one total print out.
231 extern bool warned_about_auto_ptr;
232
233 // These are useful macros that behave like functions in the code.
234 // If you want to make sure you are accessing a section of code just
235 // stick a libmesh_here(); in it, for example
236 #define libmesh_here() \
237 do { \
238 libMesh::MacroFunctions::here(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
239 } while (0)
240
241 // the libmesh_stop() macro will stop the code until a SIGCONT signal
242 // is received. This is useful, for example, when determining the
243 // memory used by a given operation. A libmesh_stop() could be
244 // inserted before and after a questionable operation and the delta
245 // memory can be obtained from a ps or top. This macro only works for
246 // serial cases.
247 #define libmesh_stop() \
248 do { \
249 libMesh::MacroFunctions::stop(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
250 } while (0)
251
252 // The libmesh_dbg_var() macro indicates that an argument to a function
253 // is used only in debug mode (i.e., when NDEBUG is not defined).
254 #ifndef NDEBUG
255 #define libmesh_dbg_var(var) var
256 #else
257 #define libmesh_dbg_var(var)
258 #endif
259
260 // The libmesh_inf_var() macro indicates that an argument to a function
261 // is used only when infinite elements are enabled
262 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
263 #define libmesh_inf_var(var) var
264 #else
265 #define libmesh_inf_var(var)
266 #endif
267
268 // The libmesh_assert() macro acts like C's assert(), but throws a
269 // libmesh_error() (including stack trace, etc) instead of just exiting
270 #ifdef NDEBUG
271
272 #define libmesh_assert_msg(asserted, msg) ((void) 0)
273 #define libmesh_exceptionless_assert_msg(asserted, msg) ((void) 0)
274 #define libmesh_assert_equal_to_msg(expr1,expr2, msg) ((void) 0)
275 #define libmesh_assert_not_equal_to_msg(expr1,expr2, msg) ((void) 0)
276 #define libmesh_assert_less_msg(expr1,expr2, msg) ((void) 0)
277 #define libmesh_assert_greater_msg(expr1,expr2, msg) ((void) 0)
278 #define libmesh_assert_less_equal_msg(expr1,expr2, msg) ((void) 0)
279 #define libmesh_assert_greater_equal_msg(expr1,expr2, msg) ((void) 0)
280
281 #else
282
283 #define libmesh_assertion_types(expr1,expr2) \
284 typedef typename std::decay<decltype(expr1)>::type libmesh_type1; \
285 typedef typename std::decay<decltype(expr2)>::type libmesh_type2
286
287 #define libmesh_assert_msg(asserted, msg) \
288 do { \
289 if (!(asserted)) { \
290 libMesh::err << "Assertion `" #asserted "' failed." << std::endl; \
291 libmesh_error_msg(msg); \
292 } } while (0)
293
294 #define libmesh_exceptionless_assert_msg(asserted, msg) \
295 do { \
296 if (!(asserted)) { \
297 libMesh::err << "Assertion `" #asserted "' failed." << std::endl; \
298 libmesh_exceptionless_error(); \
299 } } while (0)
300
301 #define libmesh_assert_equal_to_msg(expr1,expr2, msg) \
302 do { \
303 if (!((expr1) == (expr2))) { \
304 libMesh::err << "Assertion `" #expr1 " == " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
305 libmesh_error(); \
306 } } while (0)
307
308 #define libmesh_assert_not_equal_to_msg(expr1,expr2, msg) \
309 do { \
310 if (!((expr1) != (expr2))) { \
311 libMesh::err << "Assertion `" #expr1 " != " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
312 libmesh_error(); \
313 } } while (0)
314
315 // Older Clang has a parsing bug that can't seem to handle the handle the decay statement when
316 // expanding these macros. We'll just fall back to the previous behavior with those older compilers
317 #if defined(__clang__) && (__clang_major__ <= 3)
318
319 #define libmesh_assert_less_msg(expr1,expr2, msg) \
320 do { \
321 if (!((expr1) < (expr2))) { \
322 libMesh::err << "Assertion `" #expr1 " < " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
323 libmesh_error(); \
324 } } while (0)
325
326 #define libmesh_assert_greater_msg(expr1,expr2, msg) \
327 do { \
328 if (!((expr1) > (expr2))) { \
329 libMesh::err << "Assertion `" #expr1 " > " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
330 libmesh_error(); \
331 } } while (0)
332
333 #define libmesh_assert_less_equal_msg(expr1,expr2, msg) \
334 do { \
335 if (!((expr1) <= (expr2))) { \
336 libMesh::err << "Assertion `" #expr1 " <= " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
337 libmesh_error(); \
338 } } while (0)
339
340 #define libmesh_assert_greater_equal_msg(expr1,expr2, msg) \
341 do { \
342 if (!((expr1) >= (expr2))) { \
343 libMesh::err << "Assertion `" #expr1 " >= " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
344 libmesh_error(); \
345 } } while (0)
346
347 // Newer clangs and all other supported compilers
348 #else
349
350 template <template <class> class Comp>
351 struct casting_compare {
352
353 template <typename T1, typename T2>
operatorcasting_compare354 bool operator()(const T1 & e1, const T2 & e2) const
355 {
356 typedef typename std::decay<T1>::type DT1;
357 typedef typename std::decay<T2>::type DT2;
358 return (Comp<DT2>()(static_cast<DT2>(e1), e2) &&
359 Comp<DT1>()(e1, static_cast<DT1>(e2)));
360 }
361
362 template <typename T1>
operatorcasting_compare363 bool operator()(const T1 & e1, const T1 & e2) const
364 {
365 return Comp<T1>()(e1, e2);
366 }
367 };
368
369 #define libmesh_assert_less_msg(expr1,expr2, msg) \
370 do { \
371 if (!libMesh::casting_compare<std::less>()(expr1, expr2)) { \
372 libMesh::err << "Assertion `" #expr1 " < " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
373 libmesh_error(); \
374 } } while (0)
375
376 #define libmesh_assert_greater_msg(expr1,expr2, msg) \
377 do { \
378 if (!libMesh::casting_compare<std::greater>()(expr1, expr2)) { \
379 libMesh::err << "Assertion `" #expr1 " > " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
380 libmesh_error(); \
381 } } while (0)
382
383 #define libmesh_assert_less_equal_msg(expr1,expr2, msg) \
384 do { \
385 if (!libMesh::casting_compare<std::less_equal>()(expr1, expr2)) { \
386 libMesh::err << "Assertion `" #expr1 " <= " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
387 libmesh_error(); \
388 } } while (0)
389
390 #define libmesh_assert_greater_equal_msg(expr1,expr2, msg) \
391 do { \
392 if (!libMesh::casting_compare<std::greater_equal>()(expr1, expr2)) { \
393 libMesh::err << "Assertion `" #expr1 " >= " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
394 libmesh_error(); \
395 } } while (0)
396
397 #endif // clang <4.0
398
399 #endif
400
401
402 #define libmesh_assert(asserted) libmesh_assert_msg(asserted, "")
403 #define libmesh_exceptionless_assert(asserted) libmesh_exceptionless_assert_msg(asserted, "")
404 #define libmesh_assert_equal_to(expr1,expr2) libmesh_assert_equal_to_msg(expr1,expr2, "")
405 #define libmesh_assert_not_equal_to(expr1,expr2) libmesh_assert_not_equal_to_msg(expr1,expr2, "")
406 #define libmesh_assert_less(expr1,expr2) libmesh_assert_less_msg(expr1,expr2, "")
407 #define libmesh_assert_greater(expr1,expr2) libmesh_assert_greater_msg(expr1,expr2, "")
408 #define libmesh_assert_less_equal(expr1,expr2) libmesh_assert_less_equal_msg(expr1,expr2, "")
409 #define libmesh_assert_greater_equal(expr1,expr2) libmesh_assert_greater_equal_msg(expr1,expr2, "")
410
411 // The libmesh_error() macro prints a message and throws a LogicError
412 // exception
413 //
414 // The libmesh_not_implemented() macro prints a message and throws a
415 // NotImplemented exception
416 //
417 // The libmesh_file_error(const std::string & filename) macro prints a message
418 // and throws a FileError exception
419 //
420 // The libmesh_convergence_failure() macro
421 // throws a ConvergenceFailure exception
422 #define libmesh_error_msg(msg) \
423 do { \
424 libMesh::err << msg << std::endl; \
425 std::stringstream msg_stream; \
426 msg_stream << msg; \
427 libMesh::MacroFunctions::report_error(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
428 LIBMESH_THROW(libMesh::LogicError(msg_stream.str())); \
429 } while (0)
430
431 #define libmesh_error() libmesh_error_msg("")
432
433 #define libmesh_error_msg_if(cond, msg) \
434 do { \
435 if (cond) \
436 libmesh_error_msg(msg); \
437 } while (0)
438
439 #define libmesh_exceptionless_error_msg(msg) \
440 do { \
441 libMesh::err << msg << std::endl; \
442 libMesh::MacroFunctions::report_error(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
443 std::terminate(); \
444 } while (0)
445
446 #define libmesh_exceptionless_error() libmesh_exceptionless_error_msg("")
447
448 #define libmesh_not_implemented_msg(msg) \
449 do { \
450 libMesh::err << msg << std::endl; \
451 libMesh::MacroFunctions::report_error(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
452 LIBMESH_THROW(libMesh::NotImplemented()); \
453 } while (0)
454
455 #define libmesh_not_implemented() libmesh_not_implemented_msg("")
456
457 #define libmesh_file_error_msg(filename, msg) \
458 do { \
459 libMesh::err << "Error with file `" << filename << "'" << std::endl; \
460 libMesh::MacroFunctions::report_error(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
461 libMesh::err << msg << std::endl; \
462 LIBMESH_THROW(libMesh::FileError(filename)); \
463 } while (0)
464
465 #define libmesh_file_error(filename) libmesh_file_error_msg(filename,"")
466
467 #define libmesh_convergence_failure() \
468 do { \
469 LIBMESH_THROW(libMesh::ConvergenceFailure()); \
470 } while (0)
471
472 // The libmesh_example_requires() macro prints a message and calls
473 // "return 77;" if the condition specified by the macro is not true. This
474 // macro is used in the example executables, which should run when the
475 // configure-time libMesh options support them but which should exit
476 // without failure otherwise.
477 //
478 // This macro only works in main(), because we have no better way than
479 // "return" from main to immediately exit successfully - std::exit()
480 // gets seen by at least some MPI stacks as failure.
481 //
482 // 77 is the automake code for a skipped test.
483
484 #define libmesh_example_requires(condition, option) \
485 do { \
486 if (!(condition)) { \
487 libMesh::out << "Configuring libMesh with " << option << " is required to run this example." << std::endl; \
488 return 77; \
489 } } while (0)
490
491 // The libmesh_do_once macro helps us avoid redundant repeated
492 // repetitions of the same warning messages
493 #undef libmesh_do_once
494 #define libmesh_do_once(do_this) \
495 do { \
496 static bool did_this_already = false; \
497 if (!did_this_already) { \
498 did_this_already = true; \
499 do_this; \
500 } } while (0)
501
502
503 // The libmesh_warning macro outputs a file/line/time stamped warning
504 // message, if warnings are enabled.
505 #ifdef LIBMESH_ENABLE_WARNINGS
506 #define libmesh_warning(message) \
507 libmesh_do_once(libMesh::out << message \
508 << __FILE__ << ", line " << __LINE__ << ", compiled " << LIBMESH_DATE << " at " << LIBMESH_TIME << " ***" << std::endl;)
509 #else
510 #define libmesh_warning(message) ((void) 0)
511 #endif
512
513 // The libmesh_experimental macro warns that you are using
514 // bleeding-edge code
515 #undef libmesh_experimental
516 #define libmesh_experimental() \
517 libmesh_warning("*** Warning, This code is untested, experimental, or likely to see future API changes: ");
518
519
520 // The libmesh_deprecated macro warns that you are using obsoleted code
521 #undef libmesh_deprecated
522 #ifndef LIBMESH_ENABLE_DEPRECATED
523 #define libmesh_deprecated() \
524 libmesh_error_msg("*** Error, This code is deprecated, and likely to be removed in future library versions! ");
525 #else
526 #define libmesh_deprecated() \
527 libmesh_warning("*** Warning, This code is deprecated, and likely to be removed in future library versions! ");
528 #endif
529
530 // A function template for ignoring unused variables. This is a way
531 // to shut up unused variable compiler warnings on a case by case
532 // basis.
libmesh_ignore(const Args &...)533 template<class ...Args> inline void libmesh_ignore( const Args&... ) { }
534
535
536 // cast_ref and cast_ptr do a dynamic cast and assert
537 // the result, if we have RTTI enabled and we're in debug or
538 // development modes, but they just do a faster static cast if we're
539 // in optimized mode.
540 //
541 // Use these casts when you're certain that a cast will succeed in
542 // correct code but you want to be able to double-check.
543 template <typename Tnew, typename Told>
cast_ref(Told & oldvar)544 inline Tnew cast_ref(Told & oldvar)
545 {
546 #if !defined(NDEBUG) && defined(LIBMESH_HAVE_RTTI) && defined(LIBMESH_ENABLE_EXCEPTIONS)
547 try
548 {
549 Tnew newvar = dynamic_cast<Tnew>(oldvar);
550 return newvar;
551 }
552 catch (std::bad_cast &)
553 {
554 libMesh::err << "Failed to convert " << typeid(Told).name()
555 << " reference to " << typeid(Tnew).name()
556 << std::endl;
557 libMesh::err << "The " << typeid(Told).name()
558 << " appears to be a "
559 << typeid(*(&oldvar)).name() << std::endl;
560 libmesh_error();
561 }
562 #else
563 return(static_cast<Tnew>(oldvar));
564 #endif
565 }
566
567 // We use two different function names to avoid an odd overloading
568 // ambiguity bug with icc 10.1.008
569 template <typename Tnew, typename Told>
cast_ptr(Told * oldvar)570 inline Tnew cast_ptr (Told * oldvar)
571 {
572 #if !defined(NDEBUG) && defined(LIBMESH_HAVE_RTTI)
573 Tnew newvar = dynamic_cast<Tnew>(oldvar);
574 if (!newvar)
575 {
576 libMesh::err << "Failed to convert " << typeid(Told).name()
577 << " pointer to " << typeid(Tnew).name()
578 << std::endl;
579 libMesh::err << "The " << typeid(Told).name()
580 << " appears to be a "
581 << typeid(*oldvar).name() << std::endl;
582 libmesh_error();
583 }
584 return newvar;
585 #else
586 return(static_cast<Tnew>(oldvar));
587 #endif
588 }
589
590
591 #ifdef LIBMESH_ENABLE_DEPRECATED
592 template <typename Tnew, typename Told>
libmesh_cast_ptr(Told * oldvar)593 inline Tnew libmesh_cast_ptr (Told * oldvar)
594 {
595 libmesh_deprecated();
596
597 // we use the less redundantly named libMesh::cast_ptr now
598 return cast_ptr<Tnew>(oldvar);
599 }
600 #endif // LIBMESH_ENABLE_DEPRECATED
601
602
603 // cast_int asserts that the value of the castee is within the
604 // bounds which are exactly representable by the output type, if we're
605 // in debug or development modes, but it just does a faster static
606 // cast if we're in optimized mode.
607 //
608 // Use these casts when you're certain that a cast will succeed in
609 // correct code but you want to be able to double-check.
610 template <typename Tnew, typename Told>
cast_int(Told oldvar)611 inline Tnew cast_int (Told oldvar)
612 {
613 libmesh_assert_equal_to
614 (oldvar, static_cast<Told>(static_cast<Tnew>(oldvar)));
615
616 return(static_cast<Tnew>(oldvar));
617 }
618
619
620 template <typename Tnew, typename Told>
libmesh_cast_int(Told oldvar)621 inline Tnew libmesh_cast_int (Told oldvar)
622 {
623 // we use the less redundantly named libMesh::cast_int now
624 return cast_int<Tnew>(oldvar);
625 }
626
627
628 // build a integer representation of version
629 #define LIBMESH_VERSION_ID(major,minor,patch) (((major) << 16) | ((minor) << 8) | ((patch) & 0xFF))
630
631
632 // libmesh_override is simply a synonym for override as we now require
633 // a C++11 compiler that supports this keyword.
634 #define libmesh_override override
635
636 // libmesh_delete is simply a synonym for '=delete' as we now require
637 // a C++11 compiler that supports this keyword.
638 #define libmesh_delete =delete
639
640 // libmesh_final is simply a synonym for 'final' as we now require
641 // a C++11 compiler that supports this keyword.
642 #define libmesh_final final
643
644 // Define backwards-compatible fallthrough attribute. We could
645 // eventually also add support for other compiler-specific fallthrough
646 // attributes.
647 #ifdef LIBMESH_HAVE_CXX17_FALLTHROUGH_ATTRIBUTE
648 #define libmesh_fallthrough() [[fallthrough]]
649 #elif defined(LIBMESH_HAVE_DOUBLE_UNDERSCORE_ATTRIBUTE_FALLTHROUGH)
650 #define libmesh_fallthrough() __attribute__((fallthrough))
651 #else
652 #define libmesh_fallthrough() ((void) 0)
653 #endif
654
655 } // namespace libMesh
656
657
658 // Backwards compatibility
659 namespace libMeshEnums
660 {
661 using namespace libMesh;
662 }
663
664 // Backwards compatibility with pre-TIMPI reference
665 namespace TIMPI {}
666
667 namespace libMesh {
668 namespace Parallel {
669 using namespace TIMPI;
670 }
671 }
672
673
674 // Here we add missing types to the standard namespace. For example,
675 // std::max(double, float) etc... are well behaved but not defined
676 // by the standard. This also includes workarounds for super-strict
677 // implementations, for example Sun Studio and PGI C++. However,
678 // this necessarily requires breaking the ISO-C++ standard, and is
679 // really just a hack. As such, only do it if we are building the
680 // libmesh library itself. Specifically, *DO NOT* export this to
681 // user code or install this header.
682 //
683 // We put this at the end of libmesh_common.h so we can make use of
684 // any exotic definitions of Real above.
685 #ifdef LIBMESH_IS_COMPILING_ITSELF
686 # include "libmesh/libmesh_augment_std_namespace.h"
687 #endif
688
689
690 #endif // LIBMESH_LIBMESH_COMMON_H
691