libmesh_common.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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 
49 
50 // Include the MPI definition
51 #ifdef LIBMESH_HAVE_MPI
52 # include "libmesh/ignore_warnings.h"
53 # include <mpi.h>
55 #endif
56 
57 // _basic_ library functionality
58 #include "libmesh/libmesh_base.h"
60 extern "C" {
62 }
63 
64 // Proxy class for libMesh::out/err output
65 #include "libmesh/ostream_proxy.h"
66 
67 // Here we add missing types to the standard namespace. For example,
68 // std::max(double, float) etc... are well behaved but not defined
69 // by the standard. This also includes workarounds for super-strict
70 // implementations, for example Sun Studio and PGI C++. However,
71 // this necessarily requires breaking the ISO-C++ standard, and is
72 // really just a hack. As such, only do it if we are building the
73 // libmesh library itself. Specifically, *DO NOT* export this to
74 // user code or install this header.
75 #ifdef LIBMESH_IS_COMPILING_ITSELF
77 #endif
78 
79 // Make sure the C++03 compatible libmesh_nullptr is available
80 // throughout the library.
82 
83 namespace libMesh
84 {
85 
86 
87 // A namespace for functions used in the bodies of the macros below.
88 // The macros generally call these functions with __FILE__, __LINE__,
89 // __DATE__, and __TIME__ in the appropriate order. These should not
90 // be called by users directly! The implementations can be found in
91 // libmesh_common.C.
92 namespace MacroFunctions
93 {
94 void here(const char * file, int line, const char * date, const char * time);
95 void stop(const char * file, int line, const char * date, const char * time);
96 void report_error(const char * file, int line, const char * date, const char * time);
97 }
98 
99 // Undefine any existing macros
100 #ifdef Real
101 # undef Real
102 #endif
103 
104 //#ifdef REAL
105 //# undef REAL
106 //#endif
107 
108 #ifdef Complex
109 # undef Complex
110 #endif
111 
112 #ifdef COMPLEX
113 # undef COMPLEX
114 #endif
115 
116 #ifdef MPI_REAL
117 # undef MPI_REAL
118 #endif
119 
120 // Check to see if TOLERANCE has been defined by another
121 // package, if so we might want to change the name...
122 #ifdef TOLERANCE
123 DIE A HORRIBLE DEATH HERE...
124 # undef TOLERANCE
125 #endif
126 
127 
128 
129 // Define the type to use for real numbers
130 
131 typedef LIBMESH_DEFAULT_SCALAR_TYPE Real;
132 
133 // Define a corresponding tolerance. This is what should be
134 // considered "good enough" when doing floating point comparisons.
135 // For example, v == 0 is changed to std::abs(v) < TOLERANCE.
136 
137 #ifndef LIBMESH_DEFAULT_SINGLE_PRECISION
138 #ifdef LIBMESH_DEFAULT_TRIPLE_PRECISION
139 static const Real TOLERANCE = 1.e-8;
140 # define MPI_REAL MPI_LONG_DOUBLE
141 #else
142 static const Real TOLERANCE = 1.e-6;
143 # define MPI_REAL MPI_DOUBLE
144 #endif
145 #else
146 static const Real TOLERANCE = 2.5e-3;
147 # define MPI_REAL MPI_FLOAT
148 #endif
149 
150 // Define the type to use for complex numbers
151 // Always use std::complex<double>, as required by Petsc?
152 // If your version of Petsc doesn't support
153 // std::complex<other_precision>, then you'd better just leave
154 // Real==double
155 typedef std::complex<Real> Complex;
156 typedef std::complex<Real> COMPLEX;
157 
158 
159 // Helper functions for complex/real numbers
160 // to clean up #ifdef LIBMESH_USE_COMPLEX_NUMBERS elsewhere
161 template<typename T> inline T libmesh_real(T a) { return a; }
162 template<typename T> inline T libmesh_conj(T a) { return a; }
163 
164 template<typename T>
165 inline T libmesh_real(std::complex<T> a) { return std::real(a); }
166 
167 template<typename T>
168 inline std::complex<T> libmesh_conj(std::complex<T> a) { return std::conj(a); }
169 
170 // isnan isn't actually C++ standard yet; in contexts where it's not defined in
171 // cmath, libmesh_isnan will just end up returning false.
172 inline bool libmesh_isnan(float a) { return libmesh_C_isnan_float(a); }
173 inline bool libmesh_isnan(double a) { return libmesh_C_isnan_double(a); }
174 inline bool libmesh_isnan(long double a) { return libmesh_C_isnan_longdouble(a); }
175 
176 template <typename T>
177 inline bool libmesh_isnan(std::complex<T> a) { return (libmesh_isnan(std::real(a)) || libmesh_isnan(std::imag(a))); }
178 
179 // Same goes for isinf, which can be implemented in terms of isnan.
180 // http://stackoverflow.com/a/2249173/659433
181 template <typename T>
182 inline bool libmesh_isinf(T x) { return !libmesh_isnan(x) && libmesh_isnan(x - x); }
183 
184 template <typename T>
185 inline bool libmesh_isinf(std::complex<T> a) { return (libmesh_isinf(std::real(a)) || libmesh_isinf(std::imag(a))); }
186 
187 // Define the value type for unknowns in simulations.
188 // This is either Real or Complex, depending on how
189 // the library was configures
190 #if defined (LIBMESH_USE_REAL_NUMBERS)
191 typedef Real Number;
192 #elif defined (LIBMESH_USE_COMPLEX_NUMBERS)
193 typedef Complex Number;
194 #else
195 DIE A HORRIBLE DEATH HERE...
196 #endif
197 
198 
199 // Define the value type for error estimates.
200 // Since AMR/C decisions don't have to be precise,
201 // we default to float for memory efficiency.
202 typedef float ErrorVectorReal;
203 #define MPI_ERRORVECTORREAL MPI_FLOAT
204 
205 
206 #ifdef LIBMESH_HAVE_MPI
207 #ifndef LIBMESH_DISABLE_COMMWORLD
208 
211 extern MPI_Comm COMM_WORLD;
212 #endif
213 
217 extern MPI_Comm GLOBAL_COMM_WORLD;
218 #else
219 #ifndef LIBMESH_DISABLE_COMMWORLD
220 
224 extern int COMM_WORLD;
225 #endif
226 
231 extern int GLOBAL_COMM_WORLD;
232 #endif
233 
234 // Let's define a couple output streams - these will default
235 // to cout/cerr, but LibMeshInit (or the user) can also set them to
236 // something more sophisticated.
237 //
238 // We use a proxy class rather than references so they can be
239 // reseated at runtime.
240 
243 
244 // This global variable is to help us deprecate AutoPtr. We can't
245 // just use libmesh_deprecated() because then you get one print out
246 // per template instantiation, instead of one total print out.
247 extern bool warned_about_auto_ptr;
248 
249 // These are useful macros that behave like functions in the code.
250 // If you want to make sure you are accessing a section of code just
251 // stick a libmesh_here(); in it, for example
252 #define libmesh_here() \
253  do { \
254  libMesh::MacroFunctions::here(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
255  } while (0)
256 
257 // the libmesh_stop() macro will stop the code until a SIGCONT signal
258 // is received. This is useful, for example, when determining the
259 // memory used by a given operation. A libmesh_stop() could be
260 // inserted before and after a questionable operation and the delta
261 // memory can be obtained from a ps or top. This macro only works for
262 // serial cases.
263 #define libmesh_stop() \
264  do { \
265  libMesh::MacroFunctions::stop(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
266  } while (0)
267 
268 // The libmesh_dbg_var() macro indicates that an argument to a function
269 // is used only in debug mode (i.e., when NDEBUG is not defined).
270 #ifndef NDEBUG
271 #define libmesh_dbg_var(var) var
272 #else
273 #define libmesh_dbg_var(var)
274 #endif
275 
276 // The libmesh_assert() macro acts like C's assert(), but throws a
277 // libmesh_error() (including stack trace, etc) instead of just exiting
278 #ifdef NDEBUG
279 
280 #define libmesh_assert_msg(asserted, msg) ((void) 0)
281 #define libmesh_exceptionless_assert_msg(asserted, msg) ((void) 0)
282 #define libmesh_assert_equal_to_msg(expr1,expr2, msg) ((void) 0)
283 #define libmesh_assert_not_equal_to_msg(expr1,expr2, msg) ((void) 0)
284 #define libmesh_assert_less_msg(expr1,expr2, msg) ((void) 0)
285 #define libmesh_assert_greater_msg(expr1,expr2, msg) ((void) 0)
286 #define libmesh_assert_less_equal_msg(expr1,expr2, msg) ((void) 0)
287 #define libmesh_assert_greater_equal_msg(expr1,expr2, msg) ((void) 0)
288 
289 #else
290 
291 #define libmesh_assert_msg(asserted, msg) \
292  do { \
293  if (!(asserted)) { \
294  libMesh::err << "Assertion `" #asserted "' failed." << std::endl; \
295  libmesh_error_msg(msg); \
296  } } while (0)
297 
298 #define libmesh_exceptionless_assert_msg(asserted, msg) \
299  do { \
300  if (!(asserted)) { \
301  libMesh::err << "Assertion `" #asserted "' failed." << std::endl; \
302  libmesh_exceptionless_error(); \
303  } } while (0)
304 
305 #define libmesh_assert_equal_to_msg(expr1,expr2, msg) \
306  do { \
307  if (!(expr1 == expr2)) { \
308  libMesh::err << "Assertion `" #expr1 " == " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
309  libmesh_error(); \
310  } } while (0)
311 
312 #define libmesh_assert_not_equal_to_msg(expr1,expr2, msg) \
313  do { \
314  if (!(expr1 != expr2)) { \
315  libMesh::err << "Assertion `" #expr1 " != " #expr2 "' failed.\n" #expr1 " = " << (expr1) << "\n" #expr2 " = " << (expr2) << '\n' << msg << std::endl; \
316  libmesh_error(); \
317  } } while (0)
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 #endif
347 
348 
349 #define libmesh_assert(asserted) libmesh_assert_msg(asserted, "")
350 #define libmesh_exceptionless_assert(asserted) libmesh_exceptionless_assert_msg(asserted, "")
351 #define libmesh_assert_equal_to(expr1,expr2) libmesh_assert_equal_to_msg(expr1,expr2, "")
352 #define libmesh_assert_not_equal_to(expr1,expr2) libmesh_assert_not_equal_to_msg(expr1,expr2, "")
353 #define libmesh_assert_less(expr1,expr2) libmesh_assert_less_msg(expr1,expr2, "")
354 #define libmesh_assert_greater(expr1,expr2) libmesh_assert_greater_msg(expr1,expr2, "")
355 #define libmesh_assert_less_equal(expr1,expr2) libmesh_assert_less_equal_msg(expr1,expr2, "")
356 #define libmesh_assert_greater_equal(expr1,expr2) libmesh_assert_greater_equal_msg(expr1,expr2, "")
357 
358 // The libmesh_error() macro prints a message and throws a LogicError
359 // exception
360 //
361 // The libmesh_not_implemented() macro prints a message and throws a
362 // NotImplemented exception
363 //
364 // The libmesh_file_error(const std::string & filename) macro prints a message
365 // and throws a FileError exception
366 //
367 // The libmesh_convergence_failure() macro
368 // throws a ConvergenceFailure exception
369 #define libmesh_error_msg(msg) \
370  do { \
371  libMesh::err << msg << std::endl; \
372  std::stringstream msg_stream; \
373  msg_stream << msg; \
374  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
375  LIBMESH_THROW(libMesh::LogicError(msg_stream.str())); \
376  } while (0)
377 
378 #define libmesh_error() libmesh_error_msg("")
379 
380 #define libmesh_exceptionless_error_msg(msg) \
381  do { \
382  libMesh::err << msg << std::endl; \
383  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
384  std::terminate(); \
385  } while (0)
386 
387 #define libmesh_exceptionless_error() libmesh_exceptionless_error_msg("")
388 
389 #define libmesh_not_implemented_msg(msg) \
390  do { \
391  libMesh::err << msg << std::endl; \
392  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
393  LIBMESH_THROW(libMesh::NotImplemented()); \
394  } while (0)
395 
396 #define libmesh_not_implemented() libmesh_not_implemented_msg("")
397 
398 #define libmesh_file_error_msg(filename, msg) \
399  do { \
400  libMesh::err << "Error with file `" << filename << "'" << std::endl;\
401  libMesh::MacroFunctions::report_error(__FILE__, __LINE__, LIBMESH_DATE, LIBMESH_TIME); \
402  libMesh::err << msg << std::endl; \
403  LIBMESH_THROW(libMesh::FileError(filename)); \
404  } while (0)
405 
406 #define libmesh_file_error(filename) libmesh_file_error_msg(filename,"")
407 
408 #define libmesh_convergence_failure() \
409  do { \
410  LIBMESH_THROW(libMesh::ConvergenceFailure()); \
411  } while (0)
412 
413 // The libmesh_example_requires() macro prints a message and calls
414 // "return 77;" if the condition specified by the macro is not true. This
415 // macro is used in the example executables, which should run when the
416 // configure-time libMesh options support them but which should exit
417 // without failure otherwise.
418 //
419 // This macro only works in main(), because we have no better way than
420 // "return" from main to immediately exit successfully - std::exit()
421 // gets seen by at least some MPI stacks as failure.
422 //
423 // 77 is the automake code for a skipped test.
424 
425 #define libmesh_example_requires(condition, option) \
426  do { \
427  if (!(condition)) { \
428  libMesh::out << "Configuring libMesh with " #option " is required to run this example." << std::endl; \
429  return 77; \
430  } } while (0)
431 
432 // The libmesh_do_once macro helps us avoid redundant repeated
433 // repetitions of the same warning messages
434 #undef libmesh_do_once
435 #define libmesh_do_once(do_this) \
436  do { \
437  static bool did_this_already = false; \
438  if (!did_this_already) { \
439  did_this_already = true; \
440  do_this; \
441  } } while (0)
442 
443 
444 // The libmesh_warning macro outputs a file/line/time stamped warning
445 // message, if warnings are enabled.
446 #ifdef LIBMESH_ENABLE_WARNINGS
447 #define libmesh_warning(message) \
448  libmesh_do_once(libMesh::out << message \
449  << __FILE__ << ", line " << __LINE__ << ", compiled " << LIBMESH_DATE << " at " << LIBMESH_TIME << " ***" << std::endl;)
450 #else
451 #define libmesh_warning(message) ((void) 0)
452 #endif
453 
454 // The libmesh_experimental macro warns that you are using
455 // bleeding-edge code
456 #undef libmesh_experimental
457 #define libmesh_experimental() \
458  libmesh_warning("*** Warning, This code is untested, experimental, or likely to see future API changes: ");
459 
460 
461 // The libmesh_deprecated macro warns that you are using obsoleted code
462 #undef libmesh_deprecated
463 #define libmesh_deprecated() \
464  libmesh_warning("*** Warning, This code is deprecated, and likely to be removed in future library versions! ");
465 
466 // A function template for ignoring unused variables. This is a way
467 // to shut up unused variable compiler warnings on a case by case
468 // basis.
469 template<class T> inline void libmesh_ignore( const T & ) { }
470 
471 
472 // cast_ref and cast_ptr do a dynamic cast and assert
473 // the result, if we have RTTI enabled and we're in debug or
474 // development modes, but they just do a faster static cast if we're
475 // in optimized mode.
476 //
477 // Use these casts when you're certain that a cast will succeed in
478 // correct code but you want to be able to double-check.
479 template <typename Tnew, typename Told>
480 inline Tnew cast_ref(Told & oldvar)
481 {
482 #if !defined(NDEBUG) && defined(LIBMESH_HAVE_RTTI) && defined(LIBMESH_ENABLE_EXCEPTIONS)
483  try
484  {
485  Tnew newvar = dynamic_cast<Tnew>(oldvar);
486  return newvar;
487  }
488  catch (std::bad_cast)
489  {
490  libMesh::err << "Failed to convert " << typeid(Told).name()
491  << " reference to " << typeid(Tnew).name()
492  << std::endl;
493  libMesh::err << "The " << typeid(Told).name()
494  << " appears to be a "
495  << typeid(*(&oldvar)).name() << std::endl;
496  libmesh_error();
497  }
498 #else
499  return(static_cast<Tnew>(oldvar));
500 #endif
501 }
502 
503 template <typename Tnew, typename Told>
504 inline Tnew libmesh_cast_ref(Told & oldvar)
505 {
506  // we use the less redundantly named libMesh::cast_ref now
507  libmesh_deprecated();
508  return cast_ref<Tnew>(oldvar);
509 }
510 
511 // We use two different function names to avoid an odd overloading
512 // ambiguity bug with icc 10.1.008
513 template <typename Tnew, typename Told>
514 inline Tnew cast_ptr (Told * oldvar)
515 {
516 #if !defined(NDEBUG) && defined(LIBMESH_HAVE_RTTI)
517  Tnew newvar = dynamic_cast<Tnew>(oldvar);
518  if (!newvar)
519  {
520  libMesh::err << "Failed to convert " << typeid(Told).name()
521  << " pointer to " << typeid(Tnew).name()
522  << std::endl;
523  libMesh::err << "The " << typeid(Told).name()
524  << " appears to be a "
525  << typeid(*oldvar).name() << std::endl;
526  libmesh_error();
527  }
528  return newvar;
529 #else
530  return(static_cast<Tnew>(oldvar));
531 #endif
532 }
533 
534 
535 template <typename Tnew, typename Told>
536 inline Tnew libmesh_cast_ptr (Told * oldvar)
537 {
538  // we use the less redundantly named libMesh::cast_ptr now
539  return cast_ptr<Tnew>(oldvar);
540 }
541 
542 
543 // cast_int asserts that the value of the castee is within the
544 // bounds which are exactly representable by the output type, if we're
545 // in debug or development modes, but it just does a faster static
546 // cast if we're in optimized mode.
547 //
548 // Use these casts when you're certain that a cast will succeed in
549 // correct code but you want to be able to double-check.
550 template <typename Tnew, typename Told>
551 inline Tnew cast_int (Told oldvar)
552 {
553  libmesh_assert_equal_to
554  (oldvar, static_cast<Told>(static_cast<Tnew>(oldvar)));
555 
556  return(static_cast<Tnew>(oldvar));
557 }
558 
559 
560 template <typename Tnew, typename Told>
561 inline Tnew libmesh_cast_int (Told oldvar)
562 {
563  // we use the less redundantly named libMesh::cast_int now
564  return cast_int<Tnew>(oldvar);
565 }
566 
567 
568 // build a integer representation of version
569 #define LIBMESH_VERSION_ID(major,minor,patch) (((major) << 16) | ((minor) << 8) | ((patch) & 0xFF))
570 
571 
579 #ifdef LIBMESH_HAVE_CXX11_OVERRIDE
580 #define libmesh_override override
581 #else
582 #define libmesh_override
583 #endif
584 
585 // Define C++03 backwards-compatible function deletion keyword.
586 #ifdef LIBMESH_HAVE_CXX11_DELETED_FUNCTIONS
587 #define libmesh_delete =delete
588 #else
589 #define libmesh_delete
590 #endif
591 
592 // Define C++03 backwards-compatible final keyword.
593 #ifdef LIBMESH_HAVE_CXX11_FINAL
594 #define libmesh_final final
595 #else
596 #define libmesh_final
597 #endif
598 
599 // Define backwards-compatible fallthrough attribute. We could
600 // eventually also add support for other compiler-specific fallthrough
601 // attributes.
602 #ifdef LIBMESH_HAVE_CXX17_FALLTHROUGH_ATTRIBUTE
603 #define libmesh_fallthrough() [[fallthrough]]
604 #elif defined(LIBMESH_HAVE_DOUBLE_UNDERSCORE_ATTRIBUTE_FALLTHROUGH)
605 #define libmesh_fallthrough() __attribute__((fallthrough))
606 #else
607 #define libmesh_fallthrough() ((void) 0)
608 #endif
609 
610 } // namespace libMesh
611 
612 
613 // Backwards compatibility
614 namespace libMeshEnums
615 {
616 using namespace libMesh;
617 }
618 
619 #endif // LIBMESH_LIBMESH_COMMON_H
T libmesh_real(T a)
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
T libmesh_conj(T a)
Tnew cast_ref(Told &oldvar)
int libmesh_C_isnan_double(double a)
Definition: libmesh_isnan.c:26
Tnew cast_ptr(Told *oldvar)
bool warned_about_auto_ptr
static const Real TOLERANCE
MPI_Comm GLOBAL_COMM_WORLD
Definition: libmesh.C:199
Tnew libmesh_cast_ptr(Told *oldvar)
MPI_Comm COMM_WORLD
Definition: libmesh.C:197
DIE A HORRIBLE DEATH HERE typedef float ErrorVectorReal
bool libmesh_isinf(T x)
Tnew cast_int(Told oldvar)
PetscErrorCode Vec x
std::complex< Real > COMPLEX
void stop(const char *file, int line, const char *date, const char *time)
int libmesh_C_isnan_float(float a)
Definition: libmesh_isnan.c:25
Tnew libmesh_cast_ref(Told &oldvar)
OStreamProxy err(std::cerr)
void libmesh_ignore(const T &)
std::complex< Real > Complex
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool libmesh_isnan(float a)
static PetscErrorCode Mat * A
void here(const char *file, int line, const char *date, const char *time)
void report_error(const char *file, int line, const char *date, const char *time)
OStreamProxy out(std::cout)
Tnew libmesh_cast_int(Told oldvar)
int libmesh_C_isnan_longdouble(long double a)
Definition: libmesh_isnan.c:27