Mesh.h 26 KB
Newer Older
1
2
3
4
// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
5
// ==  http://www.amdis-fem.org                                              ==
6
7
// ==                                                                        ==
// ============================================================================
8
9
10
11
12
13
14
15
16
17
18
19
//
// Software License for AMDiS
//
// Copyright (c) 2010 Dresden University of Technology 
// All rights reserved.
// Authors: Simon Vey, Thomas Witkowski et al.
//
// 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
35
36

/** \file Mesh.h */

/** \defgroup Triangulation Triangulation module
 * @{ <img src="triangulation.png"> @}
 *
 * Example:
 *
 * @{ <img src="hierarchicalMesh.png"> @}
 *
 * \brief
 * Contains all triangulation classes.
 */

#ifndef AMDIS_MESH_H
#define AMDIS_MESH_H

Thomas Witkowski's avatar
Thomas Witkowski committed
37
38
39
40
#include <deque>
#include <set>
#include <stdio.h>
#include "AMDiS_fwd.h"
41
42
43
44
45
46
47
48
49
50
51
52
#include "DOFAdmin.h"
#include "Line.h"
#include "Triangle.h"
#include "Tetrahedron.h"
#include "Element.h"
#include "ElInfo.h"
#include "FixVec.h"
#include "Serializable.h"
#include "BoundaryCondition.h"

namespace AMDiS {

53
54
55
  using namespace std;


56
57
58
59
60
61
62
  /** \ingroup Triangulation 
   * \brief
   * A Mesh holds all information about a triangulation. 
   */
  class Mesh : public Serializable
  {
  public:
Thomas Witkowski's avatar
Thomas Witkowski committed
63
    /// Creates a mesh with the given name of dimension dim
64
    Mesh(string name, int dim);
65

Thomas Witkowski's avatar
Thomas Witkowski committed
66
    /// Destructor
67
68
    virtual ~Mesh();

Thomas Witkowski's avatar
Thomas Witkowski committed
69
    /// Reads macro triangulation.
70
71
    void initialize();

Thomas Witkowski's avatar
Thomas Witkowski committed
72
    /// Assignment operator
73
74
75
76
77
78
79
80
    Mesh& operator=(const Mesh&);

    /** \name getting methods
     * \{
     */

    /** \brief
     * Returns geometric information about this mesh. With GeoIndex p it is 
Backofen, Rainer's avatar
Backofen, Rainer committed
81
     * specified which information is requested.
82
     */
83
84
    inline int getGeo(GeoIndex p) const 
    { 
85
      return Global::getGeo(p, dim); 
86
    }
87

Thomas Witkowski's avatar
Thomas Witkowski committed
88
    /// Returns \ref name of the mesh
89
    inline string getName() const 
90
    { 
91
      return name; 
92
    }
93

Thomas Witkowski's avatar
Thomas Witkowski committed
94
    /// Returns \ref dim of the mesh
95
    inline int getDim() const
96
    { 
97
      return dim; 
98
    }
99

100
101
    /// Returns \ref nDofEl of the mesh
    inline const int getNumberOfAllDofs() const 
102
    { 
103
      return nDofEl; 
104
    }
105

Thomas Witkowski's avatar
Thomas Witkowski committed
106
    /// Returns \ref nNodeEl of the mesh
107
108
    inline const int getNumberOfNodes() const 
    { 
109
      return nNodeEl; 
110
    }
111

Thomas Witkowski's avatar
Thomas Witkowski committed
112
    /// Returns \ref nVertices of the mesh
113
114
    inline const int getNumberOfVertices() const 
    { 
115
      return nVertices; 
116
    }
117

Thomas Witkowski's avatar
Thomas Witkowski committed
118
    /// Returns \ref nEdges of the mesh 
119
120
    inline const int getNumberOfEdges() const 
    { 
121
      return nEdges; 
122
    }
123

Thomas Witkowski's avatar
Thomas Witkowski committed
124
    /// Returns \ref nFaces of the mesh 
125
126
    inline const int getNumberOfFaces() const 
    { 
127
      return nFaces; 
128
    }
129

Thomas Witkowski's avatar
Thomas Witkowski committed
130
    /// Returns \ref nLeaves of the mesh 
131
132
    inline const int getNumberOfLeaves() const 
    { 
133
      return nLeaves; 
134
    }
135

Thomas Witkowski's avatar
Thomas Witkowski committed
136
    /// Returns \ref nElements of the mesh
137
138
    inline const int getNumberOfElements() const 
    { 
139
      return nElements; 
140
    }
141

Thomas Witkowski's avatar
Thomas Witkowski committed
142
    /// Returns \ref maxEdgeNeigh of the mesh
143
144
    inline const int getMaxEdgeNeigh() const 
    { 
145
      return maxEdgeNeigh; 
146
    }
147

Thomas Witkowski's avatar
Thomas Witkowski committed
148
    /// Returns \ref parametric of the mesh
149
150
    inline Parametric *getParametric() const 
    { 
151
      return parametric; 
152
    }
153

Thomas Witkowski's avatar
Thomas Witkowski committed
154
    /// Returns \ref diam of the mesh
155
156
    inline const WorldVector<double>& getDiameter() const 
    { 
157
      return diam; 
158
    }
159

160
161
    /// Returns nDof[i] of the mesh
    inline const int getNumberOfDofs(int i) const 
162
    { 
163
      return nDof[i]; 
164
    }
165

Thomas Witkowski's avatar
Thomas Witkowski committed
166
    /// Returns \ref elementPrototype of the mesh
167
168
    inline Element* getElementPrototype() 
    { 
169
      return elementPrototype; 
170
    }
171

Thomas Witkowski's avatar
Thomas Witkowski committed
172
    /// Returns \ref leafDataPrototype of the mesh
173
174
    inline ElementData* getElementDataPrototype() 
    { 
175
      return elementDataPrototype; 
176
    }
177

Thomas Witkowski's avatar
Thomas Witkowski committed
178
    /// Returns node[i] of the mesh 
179
180
    inline int getNode(int i) const 
    { 
181
      return node[i]; 
182
    }
183
184
185
186
187
188
189

    /** \brief
     * Allocates the number of DOFs needed at position and registers the DOFs
     * at the DOFAdmins. The number of needed DOFs is the sum over the needed
     * DOFs of all DOFAdmin objects belonging to this mesh. 
     * The return value is a pointer to the first allocated DOF. 
     */
190
    DegreeOfFreedom *getDof(GeoIndex position);
191

Thomas Witkowski's avatar
Thomas Witkowski committed
192
    /// Returns *(\ref admin[i]) of the mesh
193
    inline const DOFAdmin& getDofAdmin(int i) const 
194
    {
195
      return *(admin[i]);
196
    }
197
198

    /** \brief
199
     * Creates a DOFAdmin with name lname. nDof specifies how many DOFs 
200
201
202
     * are needed at the different positions (see \ref DOFAdmin::nrDOF).
     * A pointer to the created DOFAdmin is returned.
     */
203
    const DOFAdmin* createDOFAdmin(string lname, DimVec<int> nDof);
204
205
206
207
208

    /** \brief
     * Returns the size of \ref admin which is the number of the DOFAdmins
     * belonging to this mesh
     */
209
210
    const int getNumberOfDOFAdmin() const 
    {
211
      return admin.size();
212
    }
213
214
215
216
217

    /** \brief
     * Returns the size of \ref macroElements which is the number of
     * of macro elements of this mesh
     */
218
219
    const int getNumberOfMacros() const 
    {
220
      return macroElements.size();
221
    }
222

Thomas Witkowski's avatar
Thomas Witkowski committed
223
    /// Returns a DOFAdmin which at least manages vertex DOFs
224
225
    const DOFAdmin* getVertexAdmin() const;

Thomas Witkowski's avatar
Thomas Witkowski committed
226
227
    /// Allocates an array of DOF pointers. The array holds one pointer for 
    /// each node.
Thomas Witkowski's avatar
Thomas Witkowski committed
228
    DegreeOfFreedom **createDofPtrs();
229

Thomas Witkowski's avatar
Thomas Witkowski committed
230
    /// Returns \ref preserveCoarseDOFs of the mesh
231
232
    inline bool queryCoarseDOFs() const 
    { 
233
      return preserveCoarseDOFs;
234
    }
235

Thomas Witkowski's avatar
Thomas Witkowski committed
236
    /// Returns an iterator to the begin of \ref macroElements
237
    inline deque<MacroElement*>::iterator firstMacroElement() 
238
    {
239
      return macroElements.begin();
240
    }
241

Thomas Witkowski's avatar
Thomas Witkowski committed
242
    /// Returns macroElements[i].
243
244
    inline MacroElement *getMacroElement(int i) 
    { 
245
      return macroElements[i]; 
246
    }
247

Thomas Witkowski's avatar
Thomas Witkowski committed
248
    /// Returns an iterator to the end of \ref macroElements
249
    inline deque<MacroElement*>::iterator endOfMacroElements() 
250
    {
251
      return macroElements.end();
252
    }
253

254
    /// Returns \ref macroElements, the list of all macro elements in the mesh.
255
    deque<MacroElement*>& getMacroElements()
256
257
258
259
    {
      return macroElements;
    }

260
261
262
263
264
265
    /** \} */

    /** \name setting methods
     * \{
     */

Thomas Witkowski's avatar
Thomas Witkowski committed
266
    /// Sets \ref name of the mesh
267
    inline void setName(string aName) 
268
    { 
269
      name = aName;
270
    }
271

Thomas Witkowski's avatar
Thomas Witkowski committed
272
    /// Sets \ref nVertices of the mesh
273
274
    inline void setNumberOfVertices(int n) 
    { 
275
      nVertices = n; 
276
    }
277

Thomas Witkowski's avatar
Thomas Witkowski committed
278
    /// Sets \ref nFaces of the mesh
279
280
    inline void setNumberOfFaces(int n) 
    { 
281
      nFaces = n; 
282
    }
283

Thomas Witkowski's avatar
Thomas Witkowski committed
284
    /// Increments \ref nVertices by inc
285
286
    inline void incrementNumberOfVertices(int inc) 
    { 
287
      nVertices += inc; 
288
    }
289
 
Thomas Witkowski's avatar
Thomas Witkowski committed
290
    /// Sets \ref nEdges of the mesh
291
292
    inline void setNumberOfEdges(int n) 
    { 
293
      nEdges = n; 
294
    }
295

Thomas Witkowski's avatar
Thomas Witkowski committed
296
    /// Increments \ref nEdges by inc
297
298
    inline void incrementNumberOfEdges(int inc) 
    { 
299
      nEdges += inc; 
300
    }
301

Thomas Witkowski's avatar
Thomas Witkowski committed
302
    /// Increments \ref nFaces by inc
303
304
    inline void incrementNumberOfFaces(int inc) 
    { 
305
      nFaces += inc; 
306
    }
307

Thomas Witkowski's avatar
Thomas Witkowski committed
308
    /// Sets \ref nLeaves of the mesh
309
310
    inline void setNumberOfLeaves(int n) 
    { 
311
      nLeaves = n; 
312
    }
313

Thomas Witkowski's avatar
Thomas Witkowski committed
314
    /// Increments \ref nLeaves by inc
315
316
    inline void incrementNumberOfLeaves(int inc) 
    { 
317
      nLeaves += inc; 
318
    }
319

Thomas Witkowski's avatar
Thomas Witkowski committed
320
    /// Sets \ref nElements of the mesh
321
322
    inline void setNumberOfElements(int n) 
    { 
323
      nElements = n; 
324
    }
325

Thomas Witkowski's avatar
Thomas Witkowski committed
326
    /// Increments \ref nElements by inc
327
328
    inline void incrementNumberOfElements(int inc) 
    { 
329
      nElements += inc; 
330
    }
331

Thomas Witkowski's avatar
Thomas Witkowski committed
332
    /// Sets *\ref diam to w
333
334
    void setDiameter(const WorldVector<double>& w);

Thomas Witkowski's avatar
Thomas Witkowski committed
335
    /// Sets (*\ref diam)[i] to d
336
337
    void setDiameter(int i, double d);

Thomas Witkowski's avatar
Thomas Witkowski committed
338
    /// Sets \ref preserveCoarseDOFs = true
339
340
    inline void retainCoarseDOFs() 
    {
341
      preserveCoarseDOFs = true;
342
    }
343

Thomas Witkowski's avatar
Thomas Witkowski committed
344
    /// Sets \ref preserveCoarseDOFs = b
345
346
    inline void setPreserveCoarseDOFs(bool b) 
    {
347
      preserveCoarseDOFs = b;
348
    }
349

Thomas Witkowski's avatar
Thomas Witkowski committed
350
    /// Sets \ref preserveCoarseDOFs = false
351
352
    inline void noCoarseDOFs() 
    {
353
      preserveCoarseDOFs = false;
354
    }
355

Thomas Witkowski's avatar
Thomas Witkowski committed
356
    /// Sets \ref elementPrototype of the mesh
357
358
    inline void setElementPrototype(Element* prototype) 
    {
359
      elementPrototype = prototype;
360
361
    }
    
Thomas Witkowski's avatar
Thomas Witkowski committed
362
    /// Sets \ref elementDataPrototype of the mesh
363
364
    inline void setElementDataPrototype(ElementData* prototype) 
    {
365
      elementDataPrototype = prototype;
366
    }
367

Thomas Witkowski's avatar
Thomas Witkowski committed
368
    ///
369
370
    inline void setParametric(Parametric *param) 
    {
371
      parametric = param;
372
    }
373

Thomas Witkowski's avatar
Thomas Witkowski committed
374
    ///
375
376
    inline void setMaxEdgeNeigh(int m) 
    { 
377
      maxEdgeNeigh = m; 
378
    }
379
380
381
  
    /** \} */

Thomas Witkowski's avatar
Thomas Witkowski committed
382
    /// Creates a new Element by cloning \ref elementPrototype
383
384
    Element* createNewElement(Element *parent = NULL);

Thomas Witkowski's avatar
Thomas Witkowski committed
385
    /// Creates a new ElInfo dependent of \ref dim of the mesh
386
387
    ElInfo* createNewElInfo();

Thomas Witkowski's avatar
Thomas Witkowski committed
388
    /// Frees DOFs at the given position pointed by dof 
389
    void freeDof(DegreeOfFreedom* dof, GeoIndex position);
390

Thomas Witkowski's avatar
Thomas Witkowski committed
391
    /// Frees memory for the given element el
392
393
    void freeElement(Element* el);

Thomas Witkowski's avatar
Thomas Witkowski committed
394
    /// Performs DOF compression for all DOFAdmins (see \ref DOFAdmin::compress)
395
396
    void dofCompress();

Thomas Witkowski's avatar
Thomas Witkowski committed
397
    /// Adds a DOFAdmin to the mesh
398
    virtual void addDOFAdmin(DOFAdmin *admin);
399

400
401
402
    /// Recalculates the number of leave elements.
    void updateNumberOfLeaves();

Thomas Witkowski's avatar
Thomas Witkowski committed
403
    /// Clears \ref macroElements
404
405
    inline void clearMacroElements() 
    { 
406
      macroElements.clear();
407
    }
408
  
Thomas Witkowski's avatar
Thomas Witkowski committed
409
    /// Adds a macro element to the mesh
410
411
    void addMacroElement(MacroElement* me);

412
    /* \brief
413
414
415
     * Removes a set of macro elements from the mesh. This works only for the 
     * case, that there are no global or local refinements, i.e., all macro 
     * elements have no children.
416
     */
417
    void removeMacroElements(std::set<MacroElement*>& macros,
418
			     vector<const FiniteElemSpace*>& feSpaces);
419

Thomas Witkowski's avatar
Thomas Witkowski committed
420
    /// Frees the array of DOF pointers (see \ref createDofPtrs)
421
    void freeDofPtrs(DegreeOfFreedom **ptrs);
422

Thomas Witkowski's avatar
Thomas Witkowski committed
423
    /// Used by \ref findElementAtPoint. 
424
425
426
    bool findElInfoAtPoint(const WorldVector<double>& xy,
			   ElInfo *el_info,
			   DimVec<double>& bary,
427
			   const MacroElement *start_mel,
428
429
			   const WorldVector<double> *xy0,
			   double *sp);
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465

    /** \brief
     * Access to an element at world coordinates xy. Some applications need the 
     * access to elements at a special location in world coordinates. Examples 
     * are characteristic methods for convection problems, or the implementation
     * of a special right hand side like point evaluations or curve integrals.
     * For such purposes, a routine is available which returns an element pointer
     * and corresponding barycentric coordinates.
     *
     * \param xy world coordinates of point
     * \param elp return address for a pointer to the element at xy
     * \param pary returns barycentric coordinates of xy
     * \param start_mel initial guess for the macro element containing xy or NULL
     * \param xy0 start point from a characteristic method, see below, or NULL
     * \param sp return address for relative distance to domain boundary in a 
     *        characteristic method, see below, or NULL
     * \return true is xy is inside the domain , false otherwise
     * 
     * For a characteristic method, where \f$ xy = xy_0 - V\tau \f$, it may be 
     * convenient to know the point on the domain's boundary which lies on the 
     * line segment between the old point xy0 and the new point xy, in case that 
     * xy is outside the domain. Such information is returned when xy0 and a 
     * pointer sp!=NULL are supplied: *sp is set to the value s such that 
     * \f$ xy_0 +s (xy -xy_0) \in \partial Domain \f$, and the element and local 
     * coordinates corresponding to that boundary point will be returned via elp 
     * and bary.
     *
     * The implementation of findElementAtPoint() is based on the transformation 
     * from world to local coordinates, available via the routine worldToCoord(),
     * At the moment, findElementAtPoint() works correctly only for domains with 
     * non-curved boundary. This is due to the fact that the implementation first
     * looks for the macro-element containing xy and then finds its path through 
     * the corresponding element tree based on the macro barycentric coordinates.
     * For non-convex domains, it is possible that in some cases a point inside
     * the domain is considered as external.
     */
466
467
468
    bool findElementAtPoint(const WorldVector<double>& xy,
			    Element **elp, 
			    DimVec<double>& bary,
469
			    const MacroElement *start_mel,
470
471
			    const WorldVector<double> *xy0,
			    double *sp);
472

473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
    /** \brief
     * Returns for a given dof its world coordinates in this mesh. Because we do
     * not have any direct connection between dofs and coordinates, this function
     * has to search for the element in this mesh, that contains the dof. Than the
     * coordinates can be computed. Therefore, this function is very costly and
     * should be used for debugging purpose only.
     *
     * @param[in]    dof       A pointer to the dof we have to search for.
     * @param[in]    feSpace   The fe space to be used for the search.
     * @param[out]   coords    World vector that stores the coordinates of the dof.
     *
     * The function returns true, if the dof was found, otherwise false.
     */
    bool getDofIndexCoords(const DegreeOfFreedom* dof, 
			   const FiniteElemSpace* feSpace,
488
489
490
491
			   WorldVector<double>& coords)
    {
      return getDofIndexCoords(*dof, feSpace, coords);
    }
492

493
494
495
496
497

    /** \brief
     * This function is equal to \ref getDofIndexCoords as defined above, but takes
     * a DOF index instead of a DOF pointer.
     */
498
499
500
    bool getDofIndexCoords(DegreeOfFreedom dof, 
			   const FiniteElemSpace* feSpace,
			   WorldVector<double>& coords);
501

502
    /** \brief
503
504
505
     * Traverse the whole mesh and stores to each DOF of the given finite
     * element space the coordinates in a given DOFVector. Works in the same
     * way as the function \ref getDofIndexCoords defined above.
506
     *
507
508
     * @param[in]   feSpace   The FE space to be used for the search.
     * @param[out]  coords    DOF vector that stores the coordinates to each DOF.
509
510
511
512
     */
    void getDofIndexCoords(const FiniteElemSpace* feSpace,
			   DOFVector<WorldVector<double> >& coords);

513
514
515
516
517
518
    /** \brief
     * Traverse the mesh and get all DOFs in this mesh for a given FE space.
     *
     * @param[in]   feSpace   The FE space to be used for collecting DOFs.
     * @param[out]  allDofs   The set which is filled with all DOFs.
     */
519
    void getAllDofs(const FiniteElemSpace *feSpace, 
520
521
		    std::set<const DegreeOfFreedom*>& allDofs);

Thomas Witkowski's avatar
Thomas Witkowski committed
522
    /// Returns FILL_ANY_?D
523
524
    inline static const Flag& getFillAnyFlag(int dim) 
    {
525
      switch (dim) {
526
527
528
529
530
531
532
533
534
535
536
537
538
      case 1:
	return FILL_ANY_1D;
	break;
      case 2:
	return FILL_ANY_2D;
	break;
      case 3:
	return FILL_ANY_3D;
	break;
      default:
	ERROR_EXIT("invalid dim\n");
	return FILL_ANY_1D;
      }
539
    }
540

Thomas Witkowski's avatar
Thomas Witkowski committed
541
    /// Serialize the mesh to a file.
542
    void serialize(ostream &out);
543

Thomas Witkowski's avatar
Thomas Witkowski committed
544
    /// Deserialize a mesh from a file.
545
    void deserialize(istream &in);
546

Thomas Witkowski's avatar
Thomas Witkowski committed
547
    /// Returns \ref elementIndex and increments it by 1.
548
549
    inline int getNextElementIndex() 
    { 
550
      return elementIndex++; 
551
    }
552

Thomas Witkowski's avatar
Thomas Witkowski committed
553
    /// Returns \ref initialized.
554
555
    inline bool isInitialized() 
    {
556
557
      return initialized; 
    }
558
  
Thomas Witkowski's avatar
Thomas Witkowski committed
559
    ///
560
    inline map<BoundaryType, VertexVector*>& getPeriodicAssociations() 
561
    {
562
      return periodicAssociations;
563
    }
564

565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
    /// Returns the periodic association for a specific boundary type.
    inline VertexVector& getPeriodicAssociations(BoundaryType b)
    {
      FUNCNAME("Mesh::getPeriodicAssociations()");

      TEST_EXIT_DBG(periodicAssociations.count(b) == 1)
	("There are no periodic assoications for boundary type %d!\n", b);

      return (*(periodicAssociations[b]));
    }

    
    /// Returns whether the given boundary type is periodic, i.e., if there is
    /// a periodic association for this boundary type.
    inline bool isPeriodicAssociation(BoundaryType b)
    {
      return (periodicAssociations.count(b) == 1 ? true : false);
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
584
    ///
585
586
    bool associated(DegreeOfFreedom dof1, DegreeOfFreedom dof2);

Thomas Witkowski's avatar
Thomas Witkowski committed
587
    ///
588
589
    bool indirectlyAssociated(DegreeOfFreedom dof1, DegreeOfFreedom dof2);

Thomas Witkowski's avatar
Thomas Witkowski committed
590
    /// Returns \macroFileInfo
591
592
    inline MacroInfo* getMacroFileInfo() 
    { 
593
      return macroFileInfo;
594
    }
595

596
597
598
599
600
601
602
603
604
605
606
607
    /// Increment the value of mesh change index, see \ref changeIndex.
    inline void incChangeIndex()
    {
      changeIndex++;
    }

    /// Returns the mesh change index, see \ref changeIndex.
    inline long getChangeIndex()
    {
      return changeIndex;
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
608
    ///
609
610
    void clearMacroFileInfo();

Thomas Witkowski's avatar
Thomas Witkowski committed
611
    ///
612
613
    int calcMemoryUsage();

614
615
616
    ///
    void deleteMeshStructure();

617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
    /// In parallel computations the level of all macro elements is equal to the number
    /// of global pre refinements, \ref nParallelPreRefinements.
    inline int getMacroElementLevel()
    {
      return nParallelPreRefinements;
    }
#else
    /// In sequentiel computations the level of all macro elements is always 0.
    inline int getMacroElementLevel()
    {
      return 0;
    }
#endif

632
633
634
635
    /// Creates a map for all elements in mesh that maps from element indices
    /// to the corresponding pointers.
    void getElementIndexMap(map<int, Element*> &elIndexMap);

636
  public:
Thomas Witkowski's avatar
Thomas Witkowski committed
637
    ///
638
639
    static const Flag FILL_NOTHING;

Thomas Witkowski's avatar
Thomas Witkowski committed
640
    ///
641
    static const Flag FILL_COORDS; 
642

Thomas Witkowski's avatar
Thomas Witkowski committed
643
    ///
644
    static const Flag FILL_BOUND; 
645

Thomas Witkowski's avatar
Thomas Witkowski committed
646
    ///
647
    static const Flag FILL_NEIGH; 
648

Thomas Witkowski's avatar
Thomas Witkowski committed
649
    ///
650
    static const Flag FILL_OPP_COORDS; 
651

Thomas Witkowski's avatar
Thomas Witkowski committed
652
    ///
653
654
    static const Flag FILL_ORIENTATION; 

Thomas Witkowski's avatar
Thomas Witkowski committed
655
    ///
656
    static const Flag FILL_ADD_ALL; 
657
  
Thomas Witkowski's avatar
Thomas Witkowski committed
658
    ///
659
    static const Flag FILL_ANY_1D; 
660

Thomas Witkowski's avatar
Thomas Witkowski committed
661
    ///
662
    static const Flag FILL_ANY_2D; 
663

Thomas Witkowski's avatar
Thomas Witkowski committed
664
    ///
665
    static const Flag FILL_ANY_3D; 
666

Thomas Witkowski's avatar
Thomas Witkowski committed
667
    ///
668
    static const Flag FILL_DET;
669

Thomas Witkowski's avatar
Thomas Witkowski committed
670
    ///
671
672
673
674
675
676
    static const Flag FILL_GRD_LAMBDA;

    //**************************************************************************
    //  flags for Mesh traversal                                                
    //**************************************************************************

Thomas Witkowski's avatar
Thomas Witkowski committed
677
    ///
678
    static const Flag CALL_EVERY_EL_PREORDER;
679

Thomas Witkowski's avatar
Thomas Witkowski committed
680
    ///
681
    static const Flag CALL_EVERY_EL_INORDER;
682

Thomas Witkowski's avatar
Thomas Witkowski committed
683
    ///
684
    static const Flag CALL_EVERY_EL_POSTORDER;
685

Thomas Witkowski's avatar
Thomas Witkowski committed
686
    ///
687
    static const Flag CALL_LEAF_EL;
688

Thomas Witkowski's avatar
Thomas Witkowski committed
689
    ///
690
    static const Flag CALL_LEAF_EL_LEVEL;
691

Thomas Witkowski's avatar
Thomas Witkowski committed
692
    ///
693
    static const Flag CALL_EL_LEVEL;
694

Thomas Witkowski's avatar
Thomas Witkowski committed
695
    ///
696
697
    static const Flag CALL_MG_LEVEL;

698
699
700
    /// If set, left and right children are swapped in traverse.
    static const Flag CALL_REVERSE_MODE;

701
  protected:
Thomas Witkowski's avatar
Thomas Witkowski committed
702
    ///
703
    bool findElementAtPointRecursive(ElInfo *elinfo,
704
				     const DimVec<double>& lambda,
705
				     int outside,
706
707
				     ElInfo *final_el_info);

708
709
710
711
712
713
714
715
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
    /** \brief
     * This functions is called in parallel computations by the function \ref
     * Mesh::initialize(). It checks that the macro file has enough macro elements
     * for the number of used processors and that all macro elements are of type 0.
     * If this is not the case, that macro mesh is globally refined in an
     * apropriate way and is written to a new macro file.
     *
716
717
     * The function overwrittes the macro and periodic filenames, if a new macro
     * fule was created for the current parallel usage.
718
     *
719
720
721
722
723
724
725
     * \param[in/out]  macroFilename      Name of the macro mesh file.
     * \param[in/out]  periodicFilename   If periodic boundaries are used, name of the
     *                                    periodicity file. Otherwise, the string must
     *                                    be empty.
     * \param[in]      check              If the mesh should be checked to be a correct
     *                                    AMDiS macro mesh, the value must be 1 and 0
     *                                    otherwise.
726
     */
727
728
    void checkParallelMacroFile(string &macroFilename, 
				string &periodicFilename,
729
730
731
				int check);
#endif

732
  protected:
Thomas Witkowski's avatar
Thomas Witkowski committed
733
    /// maximal number of DOFs at one position
734
735
    static const int MAX_DOF;

Thomas Witkowski's avatar
Thomas Witkowski committed
736
    /// Name of this Mesh
737
    string name;
738

Thomas Witkowski's avatar
Thomas Witkowski committed
739
    /// Dimension of this Mesh. Doesn't have to be equal to dimension of world.
740
741
    int dim;

Thomas Witkowski's avatar
Thomas Witkowski committed
742
    /// Number of vertices in this Mesh
743
744
    int nVertices;

Thomas Witkowski's avatar
Thomas Witkowski committed
745
    /// Number of Edges in this Mesh
746
747
    int nEdges;

Thomas Witkowski's avatar
Thomas Witkowski committed
748
    /// Number of leaf elements in this Mesh
749
750
    int nLeaves;

Thomas Witkowski's avatar
Thomas Witkowski committed
751
    /// Total number of elements in this Mesh
752
753
    int nElements;

Thomas Witkowski's avatar
Thomas Witkowski committed
754
    /// Number of faces in this Mesh
755
756
757
758
759
760
761
762
763
    int nFaces;

    /** \brief
     * Maximal number of elements that share one edge; used to allocate memory 
     * to store pointers to the neighbour at the refinement/coarsening edge 
     * (only 3d);
     */
    int maxEdgeNeigh;

Thomas Witkowski's avatar
Thomas Witkowski committed
764
    /// Diameter of the mesh in the DIM_OF_WORLD directions
765
766
767
768
769
770
771
772
773
774
    WorldVector<double> diam;

    /** \brief
     * Is pointer to NULL if mesh contains no parametric elements else pointer 
     * to a Parametric object containing coefficients of the parameterization 
     * and related information
     */
    Parametric *parametric;

    /** \brief
775
776
777
778
779
780
781
     * When an element is refined, not all dofs of the coarse element must be 
     * part of the new elements. An example are centered dofs when using higher
     * lagrange basis functions. The midpoint dof of the parents element is not
     * a dof of the both children elements. Therefore, the dof can be deleted. In
     * some situation, e.g., when using multigrid techniques, it can be necessary to
     * store this coarse dofs. Then this variable must be set to true. If false, the
     * not required coarse dofs will be deleted.
782
783
784
     */
    bool preserveCoarseDOFs;

Thomas Witkowski's avatar
Thomas Witkowski committed
785
    /// Number of all DOFs on a single element
786
    int nDofEl;
787
788
789
790
791

    /** \brief
     * Number of DOFs at the different positions VERTEX, EDGE, (FACE,) CENTER on
     * an element:
     *
792
     * - nDof[VERTEX]: number of DOFs at a vertex (>= 1)
793
     *
794
     * - nDof[EDGE]: number of DOFs at an edge; if no DOFs are associated to
795
796
     *   edges, then this value is 0
     *
797
     * - nDof[FACE]: number of DOFs at a face; if no DOFs are associated to
798
799
     *   faces, then this value is 0 (only 3d)
     *
800
     * - nDof[CENTER]: number of DOFs at the barycenter; if no DOFs are 
801
802
     *   associated to the barycenter, then this value is 0
     */
803
    DimVec<int> nDof;
804

805
    /// Number of nodes on a single element where DOFs are located. Needed for 
806
    /// the (de-) allocation of the DOF-vector on the element (\ref Element::dof).
807
    /// Here "node" is equivalent to the number of basis functions on the element.
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
    int nNodeEl;

    /** \brief
     * Gives the index of the first node at vertex, edge, face (only 3d), and 
     * barycenter:
     *
     * - node[VERTEX]: has always value 0; dof[0],...,dof[N_VERTICES-1] are 
     *   always DOFs at the vertices;
     *
     * - node[EDGE]: dof[node[EDGE]],..., dof[node[EDGE]+N_EDGES-1] are the DOFs
     *   at the N_EDGES edges, if DOFs are located at edges;
     *
     * - node[FACE]: dof[node[FACE]],..., dof[node[FACE]+N_FACES-1] are the DOFs
     *   at the N_FACES faces, if DOFs are located at faces (only 3d);
     *
     * - node[CENTER]: dof[node[CENTER]] are the DOFs at the barycenter, if DOFs
     *   are located at the barycenter;
     */
    DimVec<int> node;

Thomas Witkowski's avatar
Thomas Witkowski committed
828
    /// List of all DOFAdmins
829
    vector<DOFAdmin*> admin;
830

Thomas Witkowski's avatar
Thomas Witkowski committed
831
    /// List of all MacroElements of this Mesh
832
    deque<MacroElement*> macroElements;
833

Thomas Witkowski's avatar
Thomas Witkowski committed
834
    /// Used by check functions
835
    static vector<DegreeOfFreedom> dof_used;
836

837
838
839
840
841
842
843
844
845
846
    /** \brief
     * This map is used for serialization and deserialization of mesh elements. 
     * During the serialization process, all elements are visited and their dof indices
     * are written to the file. If a dof index at a position, i.e. vertex, line or face,
     * was written to file, the combination of dof index and position is inserted to
     * this map. That ensures that the same dof at the same position, but being part of
     * another element, is not written twice to the file. 
     * When a state should be deserialized, the information can be used to construct
     * exactly the same dof structure.
     */
847
    static map<pair<DegreeOfFreedom, int>, DegreeOfFreedom*> serializedDOFs;
848
849
850
851
852
853
854
855
856

    /** \brief
     * Used while mesh refinement. To create new elements 
     * elementPrototype->clone() is called, which returns a Element of the
     * same type as elementPrototype. So e.g. Elements of the different
     * dimensions can be created in a uniform way. 
     */
    Element* elementPrototype;

Thomas Witkowski's avatar
Thomas Witkowski committed
857
858
    /// Prototype for leaf data. Used for creation of new leaf data while 
    /// refinement.
859
860
    ElementData* elementDataPrototype;

Thomas Witkowski's avatar
Thomas Witkowski committed
861
    /// Used for enumeration of all mesh elements
862
863
    int elementIndex;

Thomas Witkowski's avatar
Thomas Witkowski committed
864
    /// True if the mesh is already initialized, false otherwise.
865
866
    bool initialized;

Thomas Witkowski's avatar
Thomas Witkowski committed
867
    /// Map of managed periodic vertex associations.
868
    map<BoundaryType, VertexVector*> periodicAssociations;
869

Thomas Witkowski's avatar
Thomas Witkowski committed
870
871
    /// If the mesh has been created by reading a macro file, here the information
    /// are stored about the content of the file.
872
    MacroInfo *macroFileInfo;
873

Thomas Witkowski's avatar
Thomas Witkowski committed
874
875
876
877
    /// This index is incremented every time the mesh is changed, e.g. by the 
    /// refinement or the coarsening manager. It can be used by other object if 
    /// the mesh has been changed by first copying this variable elsewhere and 
    /// comparing its values.
878
879
    long changeIndex;

880
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
Thomas Witkowski's avatar
Thomas Witkowski committed
881
882
883
    /// In parallel computations the mesh may be globally prerefined to achieve a
    /// fine enought starting mesh for the given number of ranks. The value of the
    /// variable will be defined in function \ref checkParallelMacroFile.
884
885
886
    int nParallelPreRefinements;
#endif

887
  protected:
Thomas Witkowski's avatar
Thomas Witkowski committed
888
    /// for findElement-Fcts
889
    DimVec<double> final_lambda;
890
891

    /** \brief
892
893
     * Temporary variables that are used in functions \ref fineElInfoatPoint and
     * \ref fineElementAtPointRecursive.
894
     */
895
    const WorldVector<double> *g_xy0, *g_xy;
896
897

    /** \brief
898
899
     * Temporary variable that is used in functions \ref fineElInfoatPoint and
     * \ref fineElementAtPointRecursive.
900
901
     */    
    double *g_sp;
902
903
904
905
906
907
908
909
910
911
912
913
914
   
    friend class MacroInfo;
    friend class MacroReader;
    friend class MacroWriter;
    friend class MacroElement;
    friend class Element;
  };

}

#endif  // AMDIS_MESH_H