Global.h 14.5 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
/******************************************************************************
 *
 * AMDiS - Adaptive multidimensional simulations
 *
 * Copyright (C) 2013 Dresden University of Technology. All Rights Reserved.
 * Web: https://fusionforge.zih.tu-dresden.de/projects/amdis
 *
 * Authors: 
 * Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al.
 *
 * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
 * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
 *
 *
 * This file is part of AMDiS
 *
 * See also license.opensource.txt in the distribution.
 * 
 ******************************************************************************/
20

21
22
23
24
25
26
27
28
29
30
31
32
33
34

/** \file Global.h */

/** \mainpage AMDiS
 * @{ <img src="vis.png"> @}
 */

/** \defgroup Common Common
 */

#ifndef AMDIS_GLOBAL_H
#define AMDIS_GLOBAL_H

#if (__GNUC__) && (__GNUC__ > 2) 
35
#define OPENMODE std::ios::openmode
36
#else
37
#define OPENMODE std::ios::open_mode
38
39
40
#endif

/** \brief current AMDiS version */
41
42
43
#ifndef AMDIS_VERSION
#define AMDIS_VERSION  "AMDiS: Version 0.9.1"
#endif
44
45
46

#include <string>
#include <vector>
47
#include <set>
48
49
50
51
52
53
#include <fstream>
#include <math.h>
#include <iostream>
#include <stdio.h>
#include <functional>
#include <float.h>
54
#include <time.h>
55

56
57
58
59
60
61
#ifdef _WIN32
#include <io.h>
#else
#include <unistd.h>
#endif

62
#if HAVE_PARALLEL_DOMAIN_AMDIS
63
#include <mpi.h>
64
65
#endif

66
#include <boost/algorithm/string.hpp>
Praetorius, Simon's avatar
Praetorius, Simon committed
67
#include <boost/algorithm/string/trim.hpp>
Praetorius, Simon's avatar
Praetorius, Simon committed
68
#include <boost/math/special_functions/fpclassify.hpp>
69
#include "boost/tuple/tuple.hpp"
70
#include "AMDiS_fwd.h"
Thomas Witkowski's avatar
Blub  
Thomas Witkowski committed
71
#include "OpenMP.h"
72

73
namespace AMDiS {
74

75
76
  using namespace std;

77
78
  extern const char *funcName;

79
  const int amdisRevisionNumber = 1700;
80

81
  /// Used by matrix vector multiplication
82
83
84
85
  typedef enum { NoTranspose,
		 Transpose,
		 ConjugateTranspose } MatrixTranspose;

86
  /// Speciefies the norm used by Estimator.
87
88
89
  typedef enum { NO_NORM = 0, 
		 H1_NORM = 1, 
		 L2_NORM = 2 } Norm;
90

91
  /// Datatype for degrees of freedom 
92
93
  typedef signed int DegreeOfFreedom;

94
  /// Defines type for a vector of DOF pointers.
95
  typedef vector<const DegreeOfFreedom*> DofContainer;
96

97
98
  typedef std::set<const DegreeOfFreedom*> DofContainerSet;

99
  /// Defines a type for global edge identification via its DOFs.
100
  typedef pair<DegreeOfFreedom, DegreeOfFreedom> DofEdge;
101
102
103

  /// Defines a tzpe for global face identiication via its DOFs.
  typedef boost::tuple<DegreeOfFreedom, DegreeOfFreedom, DegreeOfFreedom> DofFace;
104

105

106
  /// Returns the GeoIndex of d for dimension dim.
107
108
#define INDEX_OF_DIM(d, dim) (static_cast<GeoIndex>((d == dim) ? CENTER : d + 1))

109
  /// Returns the dimension of GeoIndex ind for dimension dim
110
111
#define DIM_OF_INDEX(ind, dim) ((static_cast<int>(ind) == 0) ? dim : static_cast<int>(ind) - 1)

Thomas Witkowski's avatar
Blub  
Thomas Witkowski committed
112
113
114
115

#if SUPPRESS_OUTPUT
#define PRINT_LINE(stream, line)
#else
116
#if HAVE_PARALLEL_DOMAIN_AMDIS
117
#define PRINT_LINE(stream, line) stream << "[" << MPI::COMM_WORLD.Get_rank() << "] " << line
118
#else
119
#define PRINT_LINE(stream, line) stream << line
Thomas Witkowski's avatar
Blub  
Thomas Witkowski committed
120
#endif
121
122
#endif

123
  /// Calculates factorial of i
124
125
  int fac(int i);

126
127
  void waitSec(int seconds);

128
  void processMemUsage(double& vm_usage, double& resident_set, bool inMegaByte = true);
129

130
  /// Content comparision of two pointers. Used e.g. for find_if
131
  template<typename T>
132
  struct comparePtrContents : public binary_function<T*, T*, bool>
133
  {
134
    /// Overrides binary_function::operator()
135
136
    bool operator()(T* a, T* b) const 
    {
137
      return (*a == *b);
138
    }
139
140
  };

141
  
Praetorius, Simon's avatar
Praetorius, Simon committed
142
143
144
  /// check for file existence
  inline bool file_exists(const std::string filename)
  {
145
146
147
#ifdef _WIN32
    return _access(filename.c_str(), 0) == 0;
#else
Praetorius, Simon's avatar
Praetorius, Simon committed
148
    return access(filename.c_str(), F_OK) == 0;
149
#endif
150
  }
Praetorius, Simon's avatar
Praetorius, Simon committed
151
152
153
154
155
156
  
  /// check for inf and nan values
  inline bool isNumber(double val)
  {
    return !boost::math::isnan(val) && !boost::math::isinf(val);
  }
Praetorius, Simon's avatar
Praetorius, Simon committed
157

158
  
Praetorius, Simon's avatar
Praetorius, Simon committed
159
160
161
162
163
164
  /// trim std::string
  inline std::string trim(const std::string& oldStr)
  {
    std::string swap(oldStr);
    boost::algorithm::trim(swap);
    return swap;
165
  }
Praetorius, Simon's avatar
Praetorius, Simon committed
166

167
168
169
170
171
172
173
174
175
176
177
178
  // ===== some simple template functions ====================================== 

  template<typename T> inline T abs(T a) 
  {
    return  ((a) >= 0 ? (a) : -(a));
  }

  template<typename T> inline T sqr(T a) 
  {
    return  ((a)*(a));
  }

179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
  template<typename T> inline void nullify(T &a)
  {
    a = 0.0;
  }

  template<typename T> inline void nullify(std::vector<T> &a)
  {
    typename std::vector<T>::iterator it;
    for (it = a.begin(); it != a.end(); it++)
      nullify(*it);
  }

  template<typename T> inline void nullify(mtl::dense_vector<T> &a)
  {
    T null;
    nullify(null);
    a = null;
  }

  template<typename T> inline void nullify(WorldVector<T> &a)
  {
200
201
202
203
204
205
206
    T null; nullify(null);
    a = null;
  }

  template<typename T> inline void nullify(WorldMatrix<T> &a)
  {
    T null; nullify(null);
207
208
209
    a = null;
  }

210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
  // ===== some predefined values ===============================================
  const double m_e = 		2.7182818284590452354;
  const double m_log2e = 		1.4426950408889634074;
  const double m_log10e = 	0.43429448190325182765;
  const double m_ln2 = 		0.69314718055994530942;
  const double m_ln10 = 		2.30258509299404568402;
  const double m_pi = 		3.14159265358979323846;
  const double m_pi_2 = 		1.57079632679489661923;
  const double m_pi_4 = 		0.78539816339744830962;
  const double m_1_pi = 		0.31830988618379067154;
  const double m_2_pi = 		0.63661977236758134308;
  const double m_2_sqrtpi = 	1.12837916709551257390;
  const double m_sqrt2 = 		1.41421356237309504880;
  const double m_sqrt1_2 = 	0.70710678118654752440;

  // ===== tolerance for floating point comparison ==============================
#define DBL_TOL DBL_EPSILON
#define FLT_TOL FLT_EPSILON


  /** \brief
   * Manages the output of messages, warnings, errors, ...
   * Used by the macros FUNCNAME, ERROR, ERROR_EXIT, WARNING, TEST, MSG, INFO,
   * PRINT_INFO, WAIT, WAIT_REALLY.
   * Don't use this class directly but only via these macros!
   */
  class Msg
  {
  public:
239
    /// Prints a formated message to the message stream
240
241
    static void print(const char *format, ...);

242
    /// Prints a formated message to the error stream
243
244
    static void print_error(const char *format, ...);

245
    /// Prints a formated message to the error stream and exits 
246
247
    static void print_error_exit(const char *format, ...);

248
249
    ///
    static void catch_error_exit(const char *format, ...) {}
250

251
    /// Prints an error message with funcname, file, and line to the error stream
252
253
254
255
    static void print_error_funcname(const char *funcname,
				     const char *file, 
				     int line);

256
    /// Prints a warning to the message stream
257
258
    static void print_warn(const char *format, ...);

259
    /// Prints a warning with funcname, file, and line to the message stream
260
261
262
263
    static void print_warn_funcname(const char *funcname,
				    const char *file, 
				    int line);

264
    /// Prints the funcname to the message stream
265
266
    static void print_funcname(const char *funcname);

267
    /// Changes the message stream
268
    static void change_out(ostream*);
269

270
    /// Changes the error stream 
271
    static void change_error_out(ofstream *fp);
272

273
    /// Creates a filestream and sets the error stream to this filestream
274
275
    static void open_error_file(const char *filename, OPENMODE);

276
    /// Sets \ref msgInfo
277
278
279
280
    static void setMsgInfo(int info) 
    { 
      msgInfo = info; 
    }
281

282
    /// Returns \ref msgInfo
283
284
285
286
    static int  getMsgInfo() 
    { 
      return msgInfo; 
    }
287

288
    /// Sets \ref msgWait
289
290
291
292
    static void setMsgWait(bool wait) 
    { 
      msgWait = wait; 
    }
293

294
    /// Returns \ref msgWait
295
296
297
298
    static bool getMsgWait() 
    { 
      return msgWait; 
    }
299

300
    /// Waits for enter if w is true
301
302
    static void wait(bool w);

303
    /// Returns \ref out
304
    static ostream *getOutStream() 
305
306
307
    { 
      return out; 
    }
308

309
    /// Returns \ref error
310
    static ostream *getErrorStream() 
311
312
313
    { 
      return error; 
    }
314

315
316
317
318
319
320
321
322
  public:
#if HAVE_PARALLEL_DOMAIN_AMDIS
    /// In parallel computations, when this variable is true, only the 0 rank will
    /// print messages to the output stream. Error messages and warnings are always
    /// printed from all ranks.
    static bool outputMainRank;
#endif

323
  protected:
324
    /// Message stram
325
    static ostream *out;
326

327
    /// Error stream
328
    static ostream *error;
329

Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
330
331
    /// Remember funcName to avoid multiple output of funcName within the same
    /// function call
Thomas Witkowski's avatar
Blub  
Thomas Witkowski committed
332
    static ThreadPrivate<const char*> oldFuncName;
333

334
    /// Global info level
335
336
    static int msgInfo;

337
    /// Spezifies whether to wait when WAIT is called
338
339
340
    static bool msgWait;
  };

Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
341
342
343
  // ===========================================================================
  // ===== message macros ======================================================
  // ===========================================================================
344

Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
345
346
  /// Should be the first call in every functions. It defines the current 
  /// function name nn for message output via MSG, WARNING, ...
347
#define FUNCNAME(nn) const char *funcName; funcName = nn;
348
349
350
351
352
353
  
#if (DEBUG == 0) 
  #define FUNCNAME_DBG(nn)
#else
  #define FUNCNAME_DBG(nn) const char *funcName; funcName = nn;
#endif
354

355
  /// prints an error message 
356
357
358
#define ERROR Msg::print_error_funcname(funcName,__FILE__, __LINE__),	\
    Msg::print_error

359
  /// prints an error message and exits 
360
361
362
#define ERROR_EXIT Msg::print_error_funcname(funcName,__FILE__, __LINE__), \
    Msg::print_error_exit

363
  /// prints a warning
364
365
366
#define WARNING Msg::print_warn_funcname(funcName,__FILE__, __LINE__),	\
    Msg::print_warn

367
  /// if test is false, an error message is printed
368
369
#define TEST(test) if ((test));else ERROR

370
  /// if test is false, an error message is printed and the program exits
Thomas Witkowski's avatar
Thomas Witkowski committed
371
#define TEST_EXIT(test) if ((test));else ERROR_EXIT
372

373
  /// In debug mode, it corresponds to ERROR_EXIT, otherwise it is noop.
374
375
#if (DEBUG == 0) 
  #define TEST_EXIT_DBG(test) if (false) Msg::catch_error_exit
376
  #define DBG_VAR(var)
377
378
#else
#define TEST_EXIT_DBG(test) if ((test));else ERROR_EXIT
379
  #define DBG_VAR(var) var
380
#endif
381

382
  /// prints a message
383
384
#define MSG Msg::print_funcname(funcName), Msg::print

385
386
387
388
389
390
#if (DEBUG == 0) 
  #define MSG_DBG
#else
  #define MSG_DBG Msg::print_funcname(funcName), Msg::print
#endif

391
  /// prints a message, if min(Msg::msgInfo, info) >= noinfo
392
#define INFO(info,noinfo)						\
393
  if (Msg::getMsgInfo() && (std::min(Msg::getMsgInfo(), (info)) >= (noinfo))) MSG
394

395
  /// prints a message, if min(Msg::msgInfo, info) >= noinfo
396
#define PRINT_INFO(info,noinfo)						\
397
  if (Msg::getMsgInfo() && (std::min(Msg::getMsgInfo(), (info)) >= (noinfo))) Msg::print
398
399


Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
400
401
402
  /// If the value of Msg::wait is not zero the macro will produce the message 
  /// 'wait for <enter> ...' and will continue after pressing the return or enter
  /// key. Otherwise the program continues without a message.
403
404
#define WAIT Msg::wait(Msg::getMsgWait())

Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
405
406
  /// produces the message 'wait for <enter> ...' and will continue after
  /// pressing the return or enter key.
407
408
409
410
#define WAIT_REALLY Msg::wait(true)

#define TIME_USED(f,s) ((double)((s)-(f))/(double)CLOCKS_PER_SEC)

Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
411
412
  /// internal used indices to represent the different geometrical objects.
  /// Used as parameter for getGeo() and as template parameter for FixVec. 
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
  typedef enum
    {
      CENTER   = 0, /**< in 1d the center is at the edge, in 2d at the face, in 3d 
		     * at the interior of an element. So a Line has no edge but
		     * only a center, a Triangle has no face but only a center.
		     */
      VERTEX   = 1, /**< index for element vertices.
		     * number of vertices is equal to number of parts and 
		     * neighbours.
		     */
      EDGE     = 2, /**< index for element edges */
      FACE     = 3, /**< index for element faces */
      DIMEN    =-1, /**< index for problem dimension */
      PARTS    =-2, /**< index for parts of an element (vertices in 1d, edges in 2d
		     * , faces in 3d). Number of element parts is equal to number
		     * of vertices and neighbours. 
		     */
      NEIGH    =-3, /**< index for neighbours of an element.
		     * Number of element neighbours is equal to number of 
		     * vertices and parts.
		     */
      WORLD    =-4, /**< index for world dimension */
      BOUNDARY =-5, /**< index for boundary nodes of an element. This could be
		     * vertices, edges or faces.
		     */
      PROJECTION=-6 /**< index for element and boundary projections */
    } GeoIndex;

#define MAXPART FACE
#define MINPART PROJECTION


  /** \ingroup Common
   * \brief
   * Static class wich holds global information about the world and all types of
   * elements.
   */
  class Global
  {
  public:
Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
453
454
    /// returns a pointer to \ref referenceElement [dim]. With this pointer you
    /// can get information about the element via Element's getGeo method.
455
456
    static const Element *getReferenceElement(int dim) 
    {
457
      FUNCNAME("Global::getReferenceElement()");
458
      TEST_EXIT(dim > 0 && dim < 4)("invalid dim: %d\n", dim);
459
      return referenceElement[dim];
460
    }
461

462
    /// returns geometrical information. Currently this is only dimOfWorld.
463
464
    static inline int getGeo(GeoIndex p) 
    {
465
466
467
      if (WORLD == p) 
	return dimOfWorld; 

468
      ERROR_EXIT("Illegal request for geometry data: part = %d!\n", p);
469
      return 0;
470
    }
471

Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
472
473
    /// Returns geometrical information about elements of the dimension dim.
    /// getGeo(VERTEX, 3) returns 4 because a Tetrahedron has 4 vertices.
474
475
    static inline int getGeo(GeoIndex p, int dim) 
    {
476
      TEST_EXIT_DBG(p >= MINPART && p <= MAXPART)
477
	("Calling for invalid geometry value %d\n",p);
478
      TEST_EXIT_DBG(dim >= 0 && dim < 4)
479
480
481
482
483
484
	("invalid dim: %d\n", dim);
      TEST_EXIT_DBG((dim != 0) || (p == PARTS || p == VERTEX || p == EDGE || p == FACE))
	("dim = 0\n");

      return geoIndexTable[dim][p - MINPART];
    }
485

486
    /// Inits the Global class with the help of Parameters.
487
488
    static void init();

489
490
    static void clear();

Thomas Witkowski's avatar
Thomas Witkowski committed
491
  private:
Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
492
493
    /// Global is a pure static class. So the constructor is private to avoid
    /// instantiation.
Thomas Witkowski's avatar
Thomas Witkowski committed
494
    Global();
495
496

  private:
497
    /// Dimension of the simulated world
498
499
    static int dimOfWorld;

Thomas Witkowski's avatar
bla  
Thomas Witkowski committed
500
501
502
503
504
    /// contains a pointer to a Line, a Triangle, and a Tetrahedron.
    /// This allows the access to information of the concrete elements via
    /// the dimension index.
    /// referenceElement[3]->getGeo(VERTEX) gives the number of vertices of a
    /// Tetrahedron wich is 4 => no switch statement is necessary.
505
    static Element *referenceElement[4];
506

507
    /// Stores the precalculated results that should be returned by Global::getGeo.
508
    static vector<vector<int> > geoIndexTable;
509
510
511
512
513
514
  };

#define COMMA ,

  const int RescheduleErrorCode = 23;

515
516
517
518
519
520
521
522
523
  /**
   * \ingroup Assembler
   * \brief
   * Specifies the type of a FirstOrderTerm 
   */
  enum FirstOrderType {
    GRD_PSI,
    GRD_PHI
  };
524

525
526
527
528
}

#endif // AMDIS_GLOBAL_H