1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 /*! \file gm.h
4 * \ingroup gm
5 */
6
7 /******************************************************************************/
8 /* */
9 /* File: gm.h */
10 /* */
11 /* Purpose: grid manager header file (the heart of ug) */
12 /* */
13 /* Author: Peter Bastian, Klaus Johannsen */
14 /* */
15 /* Institut fuer Computeranwendungen III */
16 /* Universitaet Stuttgart */
17 /* Pfaffenwaldring 27 */
18 /* 70569 Stuttgart */
19 /* email: ug@ica3.uni-stuttgart.de */
20 /* */
21 /* blockvector data structure: */
22 /* Christian Wrobel */
23 /* Institut fuer Computeranwendungen III */
24 /* Universitaet Stuttgart */
25 /* Pfaffenwaldring 27 */
26 /* 70569 Stuttgart */
27 /* email: ug@ica3.uni-stuttgart.de */
28 /* */
29 /* History: 09.03.92 begin, ug version 2.0 (as ugtypes2.h) */
30 /* 13.12.94 begin, ug version 3.0 */
31 /* 27.09.95 blockvector implemented (Christian Wrobel) */
32 /* */
33 /* Remarks: */
34 /* */
35 /* */
36 /******************************************************************************/
37
38 /****************************************************************************/
39 /* */
40 /* auto include mechanism and other include files */
41 /* */
42 /****************************************************************************/
43
44 #ifndef __GM__
45 #define __GM__
46
47 #include <climits>
48 #include <cassert>
49 #include <cstdlib>
50 #include <ctime>
51 #include <cmath>
52 #include <memory>
53
54 #include <unordered_map>
55 #include <array>
56 #include <numeric>
57
58 #include <dune/uggrid/domain/domain.h>
59 #include <dune/uggrid/low/debug.h>
60 #include <dune/uggrid/low/dimension.h>
61 #include <dune/uggrid/low/heaps.h>
62 #include <dune/uggrid/low/namespace.h>
63 #include <dune/uggrid/low/misc.h>
64 #include <dune/uggrid/low/ugenv.h>
65 #include <dune/uggrid/low/ugtypes.h>
66 #include "pargm.h"
67 #include "cw.h"
68
69 #include <dune/uggrid/parallel/ppif/ppiftypes.hh>
70 #ifdef ModelP
71 # include <dune/uggrid/parallel/ddd/dddcontext.hh>
72 #endif
73
74 /****************************************************************************/
75 /* */
76 /* consistency of commandline defines */
77 /* */
78 /****************************************************************************/
79
80 #ifdef __TWODIM__
81 #ifdef Sideon
82 #error **** two dimensional case cannot have side data ****
83 #endif
84 #endif
85
86 /****************************************************************************/
87 /* */
88 /* derive additional switches from commandline specified basic switches */
89 /* */
90 /****************************************************************************/
91
92 #ifdef Debug
93 #define DEBUG_MODE "ON"
94 #else
95 #define DEBUG_MODE "OFF"
96 #endif
97
98 START_UGDIM_NAMESPACE
99
100 /****************************************************************************/
101 /* */
102 /* "hard" switches for interpolation matrix and block-vectors */
103 /* */
104 /****************************************************************************/
105
106 /** \brief Define to have matrices > 4KB (control word too small, adds integer to matrix struct) */
107 #define __XXL_MSIZE__
108
109 /** \brief If pointer between element/centernode is stored */
110 #undef __CENTERNODE__
111
112 #ifdef ModelP
113 /* This ensures that for each master node-vector all matrix-neighbors in link depth 2 are
114 at leat as a copy on the same processor and all connections are copied (even for ghosts) */
115 /*#define __OVERLAP2__*/
116 #endif
117
118 /****************************************************************************/
119 /* */
120 /* defines in the following order */
121 /* */
122 /* compile time constants defining static data size (i.e. arrays) */
123 /* other constants */
124 /* macros */
125 /* */
126 /****************************************************************************/
127
128 /* Necessary for most C runtime libraries */
129 #undef DOMAIN
130
131 /** @name Some size parameters */
132 /*@{*/
133 /** \brief maximum depth of triangulation */
134 #define MAXLEVEL 32
135 /** \brief use 5 bits for object identification */
136 #define MAXOBJECTS 32
137 /*@}*/
138
139 /** @name Some size macros for allocation purposes */
140 /*@{*/
141 /** \brief max number of sides of an element */
142 #define MAX_SIDES_OF_ELEM 6
143 /** \brief max number of edges of an element */
144 #define MAX_EDGES_OF_ELEM 12
145 /** \brief max number of corners of an element */
146 #define MAX_CORNERS_OF_ELEM 8
147 /** \brief max number of edges of a side */
148 #define MAX_EDGES_OF_SIDE 4
149 /** \brief max number of edges meeting in a corner */
150 #define MAX_EDGES_OF_CORNER 4
151 /** \brief max number of corners of a side */
152 #define MAX_CORNERS_OF_SIDE 4
153 /** \brief an edge has always two corners.. */
154 #define MAX_CORNERS_OF_EDGE 2
155 /** \brief two sides have one edge in common */
156 #define MAX_SIDES_OF_EDGE 2
157 /** \brief max number of sons of an element */
158 enum {MAX_SONS = 30};
159 /** \brief max number of nodes on elem side */
160 #define MAX_SIDE_NODES 9
161 /** \brief max number of son edges of edge */
162 #define MAX_SON_EDGES 2
163 /** \brief max #fine sides touching a coarse*/
164 #define MAX_SIDES_TOUCHING 10
165
166 /** \todo Please doc me! */
167 #define MAX_ELEM_VECTORS (MAX_CORNERS_OF_ELEM+MAX_EDGES_OF_ELEM+1+MAX_SIDES_OF_ELEM)
168 /** \brief max number of doubles in a vector or matrix mod 32 */
169 #define MAX_NDOF_MOD_32 256
170 /** \brief max number of doubles in a vector or matrix */
171 #define MAX_NDOF 32*MAX_NDOF_MOD_32
172 /*@}*/
173
174
175 /****************************************************************************/
176 /* */
177 /* defines for algebra */
178 /* */
179 /****************************************************************************/
180
181 /** \brief Number of different data types */
182 #define MAXVOBJECTS 4
183 /** \brief max number of abstract vector types */
184 #define MAXVECTORS 4
185 #if (MAXVECTORS<MAXVOBJECTS)
186 #error *** MAXVECTORS must not be smaller than MAXVOBJECTS ***
187 #endif
188
189 /** \brief max number of geometric domain parts */
190 #define MAXDOMPARTS 4
191
192 /** \brief transforms type into bitpattern */
193 #define BITWISE_TYPE(t) (1<<(t))
194
195 /* derived sizes for algebra */
196 /** \brief max number of diff. matrix types */
197 #define MAXMATRICES MAXVECTORS*MAXVECTORS
198 /** \brief max number of diff. connections */
199 #define MAXCONNECTIONS (MAXMATRICES + MAXVECTORS)
200
201 /** \todo Please doc me! */
202 #define MATRIXTYPE(rt,ct) ((rt)*MAXVECTORS+(ct))
203 /** \todo Please doc me! */
204 #define DIAGMATRIXTYPE(rt) (MAXMATRICES+rt)
205
206 /** \brief Type of geometric entity which a certain vector is attached to */
207 enum VectorType {NOVTYPE=-1, //** Undefined */
208 NODEVEC, /**< Vector associated to a node */
209 EDGEVEC, /**< Vector associated to an edge */
210 ELEMVEC, /**< Vector associated to an element */
211 SIDEVEC /**< Vector associated to an element side */
212 };
213
214 /** @name Some constants for abstract vector type names */
215 /*@{*/
216 /** \todo Please doc me! */
217 #define FROM_VTNAME '0'
218 /** \todo Please doc me! */
219 #define TO_VTNAME 'z'
220 /** \todo Please doc me! */
221 #define MAXVTNAMES (1+TO_VTNAME-FROM_VTNAME)
222 /*@}*/
223
224 /****************************************************************************/
225 /* */
226 /* various defines */
227 /* */
228 /****************************************************************************/
229
230 /** 0 = OK as usual */
231 /** @name result codes of user supplied functions*/
232 /*@{*/
233 /** \brief coordinate out of range */
234 #define OUT_OF_RANGE 1
235 /** \brief configProblem could not init problem */
236 #define CANNOT_INIT_PROBLEM 1
237
238 /** \brief Use of GSTATUS (for grids), use power of 2 */
239 enum {GSTATUS_BDF = 1,
240 GSTATUS_INTERPOLATE = 2,
241 GSTATUS_ASSEMBLED = 4,
242 GSTATUS_ORDERED = 8};
243 /*@}*/
244
245 /** \brief Possible values for rule in MarkForRefinement */
246 enum RefinementRule
247 {NO_REFINEMENT = 0,
248 COPY = 1,
249 RED = 2,
250 #ifdef __TWODIM__
251 BLUE = 3, // For quadrilaterals
252 #endif
253 COARSE = 4,
254 #ifdef __TWODIM__
255 // The BISECTION* rules are all triangle rules
256 BISECTION_1 = 5,
257 BISECTION_2_Q = 6,
258 BISECTION_2_T1 = 7,
259 BISECTION_2_T2 = 8,
260 BISECTION_3 = 9
261 #endif
262 #ifdef __THREEDIM__
263
264 TETRA_RED_HEX = 5,
265
266 PRISM_BISECT_1_2 = 9,
267 PRISM_QUADSECT = 7,
268 PRISM_BISECT_HEX0 = 5,
269 PRISM_BISECT_HEX1 = 8,
270 PRISM_BISECT_HEX2 = 6,
271 PRISM_ROTATE_LEFT = 10,
272 PRISM_ROTATE_RGHT = 11,
273 PRISM_QUADSECT_HEXPRI0 = 14,
274 PRISM_RED_HEX = 15,
275 PRISM_BISECT_0_1, /* No explicit numbers. Maybe they are not necessary? */
276 PRISM_BISECT_0_2,
277 PRISM_BISECT_0_3,
278
279 HEX_BISECT_0_1 = 5,
280 HEX_BISECT_0_2 = 6,
281 HEX_BISECT_0_3 = 7,
282 HEX_TRISECT_0 = 8,
283 HEX_TRISECT_5 = 9,
284 HEX_QUADSECT_0 = 12,
285 HEX_QUADSECT_1 = 13,
286 HEX_QUADSECT_2 = 14,
287 HEX_BISECT_HEXPRI0 = 15,
288 HEX_BISECT_HEXPRI1 = 16
289
290 #endif
291 };
292
293 /** \brief Values for element class */
294 enum MarkClass {NO_CLASS,
295 YELLOW_CLASS,
296 GREEN_CLASS,
297 RED_CLASS,
298 SWITCH_CLASS};
299
300 /** \brief Values for node types (relative to the father element of the vertex) */
301 enum NodeType {CORNER_NODE,
302 MID_NODE,
303 SIDE_NODE,
304 CENTER_NODE,
305 LEVEL_0_NODE};
306
307 /* REMARK: TOPNODE no more available since 970411
308 because of problems in parallelisation
309 to use it in serial version uncomment define
310 #define TOPNODE(p) ((p)->iv.topnode)
311 */
312
313 /****************************************************************************/
314 /* */
315 /* format definition data structure */
316 /* */
317 /****************************************************************************/
318
319 /*----------- general typedefs ---------------------------------------------*/
320
321 /** @name General typedefs */
322 /*@{*/
323 typedef DOUBLE DOUBLE_VECTOR[DIM];
324 typedef DOUBLE DOUBLE_VECTOR_2D[2];
325 typedef DOUBLE DOUBLE_VECTOR_3D[3];
326 /*@}*/
327
328 /*----------- typedef for functions ----------------------------------------*/
329
330 /** \brief Print user data --> string
331 * @param pointer to user data
332 * @param Prefix for each line
333 * @param resulting string
334 */
335 typedef INT (*ConversionProcPtr)(void *, const char *, char *);
336
337 /** \brief Tagged print user data --> string
338 *
339 */
340 typedef INT (*TaggedConversionProcPtr)(INT, /**< Tag for data identification */
341 void *, /**< Pointer to user data */
342 const char *, /**< Prefix for each line */
343 char * /**< Resulting string */
344 );
345
346
347 /*----------- definition of structs ----------------------------------------*/
348
349 /* struct documentation is in gm.doc */
350 struct format {
351
352 /* variables of format */
353 /** \brief number of doubles in vectors */
354 INT VectorSizes[MAXVECTORS];
355
356 /** \brief a single char for abstract type name */
357 char VTypeNames[MAXVECTORS];
358
359 /** \brief number of doubles in matrices */
360 INT MatrixSizes[MAXCONNECTIONS];
361
362 /** \brief depth of connection for matrices */
363 INT ConnectionDepth[MAXCONNECTIONS];
364
365 /* table connecting parts, objects and types */
366
367 /** \brief (part,obj) --> vtype, -1 if not defined */
368 INT po2t[MAXDOMPARTS][MAXVOBJECTS];
369
370 /* derived components */
371 /** \brief maximal connection depth */
372 INT MaxConnectionDepth;
373 /** \brief geometrical depth corresponding */
374 INT NeighborhoodDepth;
375
376 /* algebraic con with depth 1 */
377 /* both derived from ConnectionDepth */
378 /** \brief type --> part, bitwise, not unique */
379 INT t2p[MAXVECTORS];
380 /** \brief type --> object, bitwise, not unique */
381 INT t2o[MAXVECTORS];
382
383 /* both derived from po2t */
384 /** \brief type --> type name */
385 char t2n[MAXVECTORS];
386
387 /** \brief type name --> type */
388 INT n2t[MAXVTNAMES];
389
390 /** \brief 0 if vector not needed for geom object */
391 INT OTypeUsed[MAXVOBJECTS];
392
393 /** \brief largest part used */
394 INT MaxPart;
395 /* both derived from po2t */
396
397 /** \brief largest type used */
398 INT MaxType;
399 /* derived from VectorSizes */
400 };
401
402 typedef struct {
403
404 /** \brief Abstract type is described here
405 *
406 * description only complete with po2t info
407 */
408 int tp;
409
410 /** \brief A single char as name of abstract type */
411 char name;
412
413 /* \brief Data size in bytes */
414 int size;
415
416 } VectorDescriptor ;
417
418 typedef struct {
419
420 /** \brief This connection goes from position 'from' */
421 int from;
422
423 /** \brief to position 'to' */
424 int to;
425
426 /** \brief 1 if diagonal, 0 if not */
427 int diag;
428
429 /** \brief Number of bytes per connection */
430 int size;
431
432 /** \brief Size of interpolation matrices */
433 int isize;
434
435 /** \brief Connect with depth in dual graph */
436 int depth;
437 } MatrixDescriptor ;
438
439 /****************************************************************************/
440 /* */
441 /* matrix/vector/blockvector data structure */
442 /* */
443 /****************************************************************************/
444
445 /* Struct documentation in gm.doc */
446 struct vector {
447
448 /** \brief object identification, various flags */
449 UINT control;
450
451 /** \brief associated object */
452 union geom_object *object;
453
454 #ifdef ModelP
455 /** \todo Please doc me! */
456 DDD_HEADER ddd;
457 #endif
458
459 /** \brief double linked list of vectors */
460 struct vector *pred,*succ;
461
462 /** \brief ordering of unknowns */
463 UINT index;
464
465 /** \brief Index if the vector is part of the leaf grid */
466 UINT leafIndex;
467
468 #ifndef ModelP // Dune uses ddd.gid for ids in parallel
469 /** \brief A unique and persistent, but not necessarily consecutive index
470
471 Used to implement face ids for Dune.
472 */
473 INT id;
474 #endif
475
476 /** \brief implements matrix */
477 struct matrix *start;
478
479 /** \brief User data */
480 DOUBLE value[1];
481 };
482 typedef struct vector VECTOR;
483
484
485 /* Struct documentation in gm.doc */
486 struct matrix {
487
488 /** \brief object identification, various flags */
489 UINT control;
490
491 #ifdef __XXL_MSIZE__
492 /** \brief for people needing large matrices */
493 UINT xxl_msize;
494 #endif
495
496 /** \brief row list */
497 struct matrix *next;
498
499 /** \brief destination vector */
500 struct vector *vect;
501
502 /** \brief User data */
503 DOUBLE value[1];
504 };
505
506 typedef struct matrix MATRIX;
507 typedef struct matrix CONNECTION;
508
509 /****************************************************************************/
510 /* */
511 /* unstructured grid data structures */
512 /* */
513 /****************************************************************************/
514
515 /*----------- typedef for functions ----------------------------------------*/
516
517
518 /*----------- definition of structs ----------------------------------------*/
519
520 /** \brief Inner vertex data structure */
521 struct ivertex {
522
523 /** \brief Object identification, various flags */
524 UINT control;
525
526 /** \brief Unique id used for load/store */
527 INT id;
528
529 /** \brief Vertex position */
530 DOUBLE x[DIM];
531
532 /** \brief Local coordinates in father element */
533 DOUBLE xi[DIM];
534
535 /* When UG is used as part of the DUNE numerics system we need
536 a few more indices per node */
537
538 /** \brief An index hat is unique and consecutive per level.
539 Controlled by DUNE */
540 int leafIndex;
541
542 #ifdef ModelP
543 /** \todo Please doc me! */
544 DDD_HEADER ddd;
545 #endif
546
547 /* pointers */
548 /** \brief Doubly linked list of vertices */
549 union vertex *pred,*succ;
550
551 /** \brief Associated user data structure */
552 void *data;
553
554 /** \brief Father element */
555 union element *father;
556
557 #ifdef TOPNODE
558 /** \brief Highest node where defect is valid
559 \todo REMARK: TOPNODE no more available since 970411
560 because of problems in parallelisation
561 to use it in serial version uncomment define TOPNODE */
562 struct node *topnode;
563 #endif
564 };
565
566 /** \brief Boundary vertex data structure */
567 struct bvertex {
568
569 /* variables */
570 /** \brief Object identification, various flags */
571 UINT control;
572
573 /** \brief Unique id used for load/store */
574 INT id;
575
576 /** \brief Vertex position */
577 DOUBLE x[DIM];
578
579 /** \brief Local coordinates in father element */
580 DOUBLE xi[DIM];
581
582 /* When UG is used as part of the DUNE numerics system we need
583 a few more indices per node */
584
585 /** \brief An index hat is unique and consecutive per level.
586 Controlled by DUNE */
587 int leafIndex;
588
589 #ifdef ModelP
590 /** \brief Information about the parallelization of this object */
591 DDD_HEADER ddd;
592 #endif
593
594 /* pointers */
595 /** \brief Doubly linked list of vertices */
596 union vertex *pred,*succ;
597
598 /** \brief Associated user data structure */
599 void *data;
600
601 /** \brief Father element */
602 union element *father;
603
604 #ifdef TOPNODE
605 /** \brief Highest node where defect is valid
606 \todo REMARK: TOPNODE no more available since 970411
607 because of problems in parallelisation
608 to use it in serial version uncomment define TOPNODE */
609 struct node *topnode;
610 #endif
611
612 /** \brief Pointer to boundary point decriptor */
613 BNDP *bndp;
614 };
615
616 /** \brief Only used to define pointer to vertex */
617 union vertex {
618 struct ivertex iv;
619 struct bvertex bv;
620 };
621
622 /** \brief A simply linked list of elements */
623 struct elementlist {
624 union element *el;
625 struct elementlist *next;
626 };
627
628 /** \brief Level-dependent part of a vertex */
629 struct node {
630
631 /** \brief Object identification, various flags */
632 UINT control;
633
634 /** \brief Unique id used for load/store */
635 INT id;
636
637 /* When UG is used as part of the DUNE numerics system we need
638 a few more indices per node */
639
640 /** \brief An index hat is unique and consecutive per level.
641 Controlled by DUNE */
642 int levelIndex;
643
644 /** \brief Information if this node is on the leaf. */
645 bool isLeaf;
646
647 #ifdef ModelP
648 /** \brief Per-node message buffer used by Dune for dynamic load-balancing */
649 char* message_buffer_;
650
651 /** \brief Size of the `message_buffer` */
652 std::size_t message_buffer_size_;
653 #endif
654
655 #ifdef ModelP
656 /** \brief Information about the parallelization of this object */
657 DDD_HEADER ddd;
658 #endif
659
660 /* pointers */
661 /** \brief Doubly linked list of nodes per level*/
662 struct node *pred,*succ;
663
664 /** \brief List of links */
665 struct link *start;
666
667 /** \brief Node or edge on coarser level (NULL if none) */
668 union geom_object *father;
669
670 /** \brief Node on finer level (NULL if none) */
671 struct node *son;
672
673 /** \brief Corresponding vertex structure */
674 union vertex *myvertex;
675
676 /** \brief Associated vector
677 *
678 * WARNING: the allocation of the vector pointer depends on the format */
679 VECTOR *vector;
680
681 /** \brief Associated data pointer
682 *
683 * WARNING: The allocation of the data pointer depends on the format */
684 void *data;
685
686 #ifdef ModelP
message_buffernode687 const char* message_buffer() const
688 { return message_buffer_; }
689
message_buffer_sizenode690 std::size_t message_buffer_size() const
691 { return message_buffer_size_; }
692
message_buffernode693 void message_buffer(char* p, std::size_t size)
694 {
695 message_buffer_ = p;
696 message_buffer_size_ = size;
697 }
698
message_buffer_freenode699 void message_buffer_free()
700 {
701 std::free(message_buffer_);
702 message_buffer(nullptr, 0);
703 }
704 #endif
705 };
706
707 /** \todo Please doc me! */
708 struct link {
709
710 /** \brief object identification, various flags */
711 UINT control;
712
713 /** \brief ptr to next link */
714 struct link *next;
715
716 /** \brief ptr to neighbor node */
717 struct node *nbnode;
718
719 /** \brief ptr to neighboring elem */
720 #if defined(__TWODIM__)
721 union element *elem;
722 #endif
723
724 };
725
726 /** \brief Undirected edge of the grid graph */
727 struct edge {
728
729 /* variables */
730 /* two links */
731 struct link links[2];
732
733 /* When UG is used as part of the DUNE numerics system we need
734 place to store codim dim-1 indices */
735
736 /** \brief An index hat is unique and consecutive per level.
737 Controlled by DUNE */
738 int levelIndex;
739
740 /** \brief An index hat is unique and consecutive on the grid surface.
741 Controlled by DUNE */
742 int leafIndex;
743
744 /** \brief A unique and persistent, but not necessarily consecutive index */
745 INT id;
746
747 #ifdef ModelP
748 /** Bookkeeping information for DDD */
749 DDD_HEADER ddd;
750 #endif
751
752 /** \brief Pointer to mid node on next finer grid */
753 struct node *midnode;
754
755 /** \brief associated vector
756 *
757 * WARNING: the allocation of the vector pointer depends on the format */
758 VECTOR *vector;
759 };
760
761 /** \brief A generic grid element
762
763 No difference between inner and boundary elements
764 */
765 struct generic_element {
766
767 /** \brief object identification, various flags */
768 UINT control;
769
770 /** \brief unique id used for load/store */
771 INT id;
772
773 /** \brief additional flags for elements */
774 UINT flag;
775
776 /** \brief to store NodeOrder for hexahedrons */
777 INT property;
778
779 /* When UG is used as part of the DUNE numerics system we need
780 a few more indices per element */
781
782 /** \brief An index hat is unique and consecutive per level.
783 Controlled by DUNE */
784 int levelIndex;
785
786 /** \brief An index hat is unique and consecutive on the grid surface.
787 Controlled by DUNE */
788 int leafIndex;
789
790 #ifdef ModelP
791 /** \brief Per-node message buffer used by Dune for dynamic load-balancing */
792 char* message_buffer;
793
794 /** \brief Size of the `message_buffer` */
795 std::size_t message_buffer_size;
796 #endif
797
798 #ifdef ModelP
799 /** \brief Information about the parallelization of this object */
800 DDD_HEADER ddd;
801
802 /** \brief Stores partition information */
803 INT lb1;
804 #endif
805
806 /** \brief double linked list of elements */
807 union element *pred, *succ;
808
809 #ifdef __CENTERNODE__
810 /** \brief Pointer to center node */
811 struct node *centernode;
812 #endif
813
814 /** \brief Element specific part of variable length array managed by ug */
815 void *refs[1];
816
817 };
818
819 /** \brief A triangle element in a 2d grid
820 \implements generic_element
821 */
822 struct triangle {
823
824 /** \brief Object identification, various flags */
825 UINT control;
826
827 /** \brief Unique id used for load/store */
828 INT id;
829
830 /** \brief Additional flags for elements */
831 UINT flag;
832
833 /** \brief Even more property bits */
834 INT property;
835
836 /* When UG is used as part of the DUNE numerics system we need
837 a few more indices per element */
838
839 /** \brief An index hat is unique and consecutive per level.
840 Controlled by DUNE */
841 int levelIndex;
842
843 /** \brief An index hat is unique and consecutive on the grid surface.
844 Controlled by DUNE */
845 int leafIndex;
846
847 #ifdef ModelP
848 /** \brief Per-node message buffer used by Dune for dynamic load-balancing */
849 char* message_buffer;
850
851 /** \brief Size of the `message_buffer` */
852 std::size_t message_buffer_size;
853 #endif
854
855 #ifdef ModelP
856 /** \brief Information about the parallelization of this object */
857 DDD_HEADER ddd;
858
859 /** \brief stores partition information */
860 INT lb1;
861 #endif
862
863 /** \brief Realize a doubly linked list of elements */
864 union element *pred, *succ;
865
866 #ifdef __CENTERNODE__
867 /** \brief Pointer to the center node */
868 struct node *centernode;
869 #endif
870
871 /** \brief Corners of this element */
872 struct node *n[3];
873
874 /** \brief Father element on next-coarser grid */
875 union element *father;
876
877 #ifdef ModelP
878 /** \brief Element tree */
879 union element *sons[2];
880 #else
881 /** \brief Element tree */
882 union element *sons[1];
883 #endif
884
885 /** \brief The neighboring elements */
886 union element *nb[3];
887
888 /** \brief Associated vector
889
890 WARNING: the allocation of the vector pointer depends on the format
891 void *ptr[4] would be possible too:
892 if there are no element vectors, the sides will be ptr[0],ptr[1],ptr[2]
893 Use the macros to find the correct address! */
894 VECTOR *vector;
895
896 /** \brief Only on the boundary, NULL if interior side */
897 BNDS *bnds[3];
898
899 /* \brief Associated data pointer
900
901 WARNING: the allocation of the data pointer depends on the format */
902 void *data;
903 };
904
905 /** \brief A quadrilateral element in a 2d grid
906 \implements generic_element
907 */
908 struct quadrilateral {
909
910 /** \brief object identification, various flags */
911 UINT control;
912
913 /** \brief unique id used for load/store */
914 INT id;
915
916 /** \brief additional flags for elements */
917 UINT flag;
918
919 /** \brief Even more flags */
920 INT property;
921
922 /* When UG is used as part of the DUNE numerics system we need
923 a few more indices per element */
924
925 /** \brief An index hat is unique and consecutive per level.
926 Controlled by DUNE */
927 int levelIndex;
928
929 /** \brief An index hat is unique and consecutive on the grid surface.
930 Controlled by DUNE */
931 int leafIndex;
932
933 #ifdef ModelP
934 /** \brief Per-node message buffer used by Dune for dynamic load-balancing */
935 char* message_buffer;
936
937 /** \brief Size of the `message_buffer` */
938 std::size_t message_buffer_size;
939 #endif
940
941 #ifdef ModelP
942 /** \brief Information about the parallelization of this object */
943 DDD_HEADER ddd;
944
945 /** \brief stores partition information */
946 INT lb1;
947 #endif
948
949 /** \brief doubly linked list of elements */
950 union element *pred, *succ;
951
952 #ifdef __CENTERNODE__
953 /** \brief pointer to center node */
954 struct node *centernode;
955 #endif
956
957 /** \brief corners of that element */
958 struct node *n[4];
959
960 /** \brief father element on coarser grid */
961 union element *father;
962
963 #ifdef ModelP
964 /** \brief Element tree */
965 union element *sons[2];
966 #else
967 /** \brief Element tree */
968 union element *sons[1];
969 #endif
970
971 /** \brief The neighbor elements */
972 union element *nb[4];
973
974 /** \brief Associated vector
975
976 WARNING: the allocation of the vector pointer depends on the format
977 void *ptr[5] would be possible too:
978 if there are no element vectors, the sides will be ptr[0],ptr[1], ..
979 Use the macros to find the correct address!
980 */
981 VECTOR *vector;
982
983 /** \brief only on bnd, NULL if interior side */
984 BNDS *bnds[4];
985
986 /** \brief Associated data pointer
987
988 WARNING: the allocation of the data pointer depends on the format */
989 void *data;
990 };
991
992 /** \brief A tetrahedral element in a 3d grid
993 \implements generic_element
994 */
995 struct tetrahedron {
996
997 /** \brief object identification, various flags */
998 UINT control;
999
1000 /** \brief Unique id used for load/store */
1001 INT id;
1002
1003 /** \brief Additional flags */
1004 UINT flag;
1005
1006 /** \brief Even more flags */
1007 INT property;
1008
1009 /* When UG is used as part of the DUNE numerics system we need
1010 a few more indices per element */
1011
1012 /** \brief An index hat is unique and consecutive per level.
1013 Controlled by DUNE */
1014 int levelIndex;
1015
1016 /** \brief An index hat is unique and consecutive on the grid surface.
1017 Controlled by DUNE */
1018 int leafIndex;
1019
1020 #ifdef ModelP
1021 /** \brief Per-node message buffer used by Dune for dynamic load-balancing */
1022 char* message_buffer;
1023
1024 /** \brief Size of the `message_buffer` */
1025 std::size_t message_buffer_size;
1026 #endif
1027
1028 #ifdef ModelP
1029
1030 /** \brief Information about the parallelization of this element */
1031 DDD_HEADER ddd;
1032
1033 /** \brief Stores partition information */
1034 INT lb1;
1035 #endif
1036
1037 /* pointers */
1038
1039 /** \brief Realize a double linked list of elements */
1040 union element *pred, *succ;
1041
1042 #ifdef __CENTERNODE__
1043 /** \brief Pointer to the center node */
1044 struct node *centernode;
1045 #endif
1046
1047 /** \brief Corners of this element */
1048 struct node *n[4];
1049
1050 /** \brief Father element on coarser grid */
1051 union element *father;
1052
1053 #ifdef ModelP
1054 /** \brief Element tree */
1055 union element *sons[2]; /* element tree */
1056 #else
1057 /** \brief Element tree */
1058 union element *sons[1];
1059 #endif
1060
1061 /** \brief The neighboring elements */
1062 union element *nb[4];
1063
1064 /** \brief Associated vector
1065
1066 WARNING: the allocation of the vector pointer depends on the format
1067 void *ptr[9] would be possible too:
1068 if there are no element vectors, the sides will be ptr[0],ptr[1], ..
1069 Use the macros to find the correct address! */
1070 VECTOR *vector;
1071
1072 /** \brief Associated vector for sides */
1073 VECTOR *sidevector[4];
1074
1075 /** \brief The boundary segments, NULL if interior side */
1076 BNDS *bnds[4];
1077
1078 /** \brief Associated data pointer
1079
1080 WARNING: the allocation of the data pointer depends on the format */
1081 void *data;
1082 };
1083
1084 /** \brief A pyramid element in a 3d grid
1085 \implements generic_element
1086 */
1087 struct pyramid {
1088
1089 /** \brief object identification, various flags */
1090 UINT control;
1091
1092 /** \brief Unique id used for load/store */
1093 INT id;
1094
1095 /** \brief Additional flags for elements */
1096 UINT flag;
1097
1098 /** \brief Even more flags */
1099 INT property;
1100
1101 /* When UG is used as part of the DUNE numerics system we need
1102 a few more indices per element */
1103
1104 /** \brief An index hat is unique and consecutive per level.
1105 Controlled by DUNE */
1106 int levelIndex;
1107
1108 /** \brief An index hat is unique and consecutive on the grid surface.
1109 Controlled by DUNE */
1110 int leafIndex;
1111
1112 #ifdef ModelP
1113 /** \brief Per-node message buffer used by Dune for dynamic load-balancing */
1114 char* message_buffer;
1115
1116 /** \brief Size of the `message_buffer` */
1117 std::size_t message_buffer_size;
1118 #endif
1119
1120 #ifdef ModelP
1121 /** \brief Information about the parallelization of this element */
1122 DDD_HEADER ddd;
1123
1124 /** \brief Stores partition information */
1125 INT lb1;
1126 #endif
1127
1128 /* pointers */
1129 /** \brief Realize a doubly linked list of elements */
1130 union element *pred, *succ;
1131
1132 #ifdef __CENTERNODE__
1133 /** \brief Pointer to center node */
1134 struct node *centernode;
1135 #endif
1136
1137 /** \brief Corners of this element */
1138 struct node *n[5];
1139
1140 /** \brief Father element on coarser grid */
1141 union element *father;
1142
1143 #ifdef ModelP
1144 /** \brief Element tree */
1145 union element *sons[2];
1146 #else
1147 /** \brief Element tree */
1148 union element *sons[1];
1149 #endif
1150
1151 /** \brief The neighbor elements */
1152 union element *nb[5];
1153
1154 /** \brief Associated vector
1155
1156 WARNING: the allocation of the vector pointer depends on the format
1157 void *ptr[11] would be possible too:
1158 if there are no element vectors, the sides will be ptr[0],ptr[1], ..
1159 Use the macros to find the correct address! */
1160 VECTOR *vector;
1161
1162 /** \brief Associated vector for sides */
1163 VECTOR *sidevector[5];
1164
1165 /** \brief The boundary segments, NULL if interior side */
1166 BNDS *bnds[5];
1167
1168 /** \brief Associated data pointer
1169
1170 WARNING: the allocation of the data pointer depends on the format */
1171 void *data;
1172 };
1173
1174 /** \brief A prism element in a 3d grid
1175 \implements generic_element
1176 */
1177 struct prism {
1178
1179 /** \brief object identification, various flags */
1180 UINT control;
1181
1182 /** \brief Unique id used for load/store */
1183 INT id;
1184
1185 /** \brief Additional flags for this element */
1186 UINT flag;
1187
1188 /** \brief Even more flags */
1189 INT property;
1190
1191 /* When UG is used as part of the DUNE numerics system we need
1192 a few more indices per element */
1193
1194 /** \brief An index hat is unique and consecutive per level.
1195 Controlled by DUNE */
1196 int levelIndex;
1197
1198 /** \brief An index hat is unique and consecutive on the grid surface.
1199 Controlled by DUNE */
1200 int leafIndex;
1201
1202 #ifdef ModelP
1203 /** \brief Per-node message buffer used by Dune for dynamic load-balancing */
1204 char* message_buffer;
1205
1206 /** \brief Size of the `message_buffer` */
1207 std::size_t message_buffer_size;
1208 #endif
1209
1210 #ifdef ModelP
1211 /** \brief Information about the parallelization of this element */
1212 DDD_HEADER ddd;
1213
1214 /** \brief Stores partition information */
1215 INT lb1;
1216 #endif
1217
1218 /* pointers */
1219 /** \brief Realize doubly linked list */
1220 union element *pred, *succ;
1221
1222 #ifdef __CENTERNODE__
1223 /** \brief Pointer to center node */
1224 struct node *centernode;
1225 #endif
1226
1227 /** \brief Corners of this element */
1228 struct node *n[6];
1229
1230 /** \brief Father element on next coarser grid */
1231 union element *father;
1232
1233 #ifdef ModelP
1234 /** \brief Element tree */
1235 union element *sons[2];
1236 #else
1237 /** \brief Element tree */
1238 union element *sons[1];
1239 #endif
1240
1241 /** \brief Neighbor elements */
1242 union element *nb[5];
1243
1244 /** \brief Associated vector
1245
1246 WARNING: the allocation of the vector pointer depends on the format
1247 void *ptr[11] would be possible too:
1248 if there are no element vectors, the sides will be ptr[0],ptr[1], ...
1249 Use the macros to find the correct address!
1250 */
1251 VECTOR *vector;
1252
1253 /** \brief Associated vectors for sides */
1254 VECTOR *sidevector[5];
1255
1256 /** \brief Boundary segments, NULL if interior side */
1257 BNDS *bnds[5];
1258
1259 /** \brief Associated data pointer
1260
1261 WARNING: the allocation of the data pointer depends on the format */
1262 void *data;
1263 };
1264
1265 /** \brief A hexahedral element in a 3d grid
1266 \implements generic_element
1267 */
1268 struct hexahedron {
1269
1270 /** \brief object identification, various flags */
1271 UINT control;
1272
1273 /** \brief Unique id used for load/store */
1274 INT id;
1275
1276 /** \brief Additional flags for this element */
1277 UINT flag;
1278
1279 /** \brief Even more flags */
1280 INT property;
1281
1282 /* When UG is used as part of the DUNE numerics system we need
1283 a few more indices per element */
1284
1285 /** \brief An index hat is unique and consecutive per level.
1286 Controlled by DUNE */
1287 int levelIndex;
1288
1289 /** \brief An index hat is unique and consecutive on the grid surface.
1290 Controlled by DUNE */
1291 int leafIndex;
1292
1293 #ifdef ModelP
1294 /** \brief Per-node message buffer used by Dune for dynamic load-balancing */
1295 char* message_buffer;
1296
1297 /** \brief Size of the `message_buffer` */
1298 std::size_t message_buffer_size;
1299 #endif
1300
1301 #ifdef ModelP
1302 /** \brief Information about the parallelization of this element */
1303 DDD_HEADER ddd;
1304
1305 /** \brief Stores partition information */
1306 INT lb1;
1307 #endif
1308
1309 /** \brief Realize doubly linked list of elements on one grid level */
1310 union element *pred, *succ;
1311
1312 #ifdef __CENTERNODE__
1313 /** \brief Pointer to center node */
1314 struct node *centernode;
1315 #endif
1316
1317 /** \brief Corners of this element */
1318 struct node *n[8];
1319
1320 /** \brief Father element on coarser grid */
1321 union element *father;
1322
1323 #ifdef ModelP
1324 /** \brief Element tree */
1325 union element *sons[2];
1326 #else
1327 /** \brief Element tree */
1328 union element *sons[1];
1329 #endif
1330
1331 /** \brief The neighboring elements */
1332 union element *nb[6];
1333
1334 /** \brief Associated vector
1335
1336 WARNING: the allocation of the vector pointer depends on the format
1337 void *ptr[13] would be possible too:
1338 if there are no element vectors, the sides will be ptr[0],ptr[1], ...
1339 Use the macros to find the correct address!
1340 */
1341 VECTOR *vector;
1342
1343 /** \brief Associated vectors for sides */
1344 VECTOR *sidevector[6];
1345
1346 /** \brief Boundary segments, NULL if interior side */
1347 BNDS *bnds[6];
1348
1349 /** \brief Associated data pointer
1350
1351 WARNING: the allocation of the data pointer depends on the format */
1352 void *data;
1353 } ;
1354
1355 /** \brief Objects that can hold an element */
1356 union element {
1357 struct generic_element ge;
1358 #ifdef __TWODIM__
1359 struct triangle tr;
1360 struct quadrilateral qu;
1361 #endif
1362 #ifdef __THREEDIM__
1363 struct tetrahedron te;
1364 struct pyramid py;
1365 struct prism pr;
1366 struct hexahedron he;
1367 #endif
1368
1369 #ifdef ModelP
message_buffer()1370 const char* message_buffer() const
1371 { return ge.message_buffer; }
1372
message_buffer_size()1373 std::size_t message_buffer_size() const
1374 { return ge.message_buffer_size; }
1375
message_buffer(char * p,std::size_t size)1376 void message_buffer(char* p, std::size_t size)
1377 {
1378 ge.message_buffer = p;
1379 ge.message_buffer_size = size;
1380 }
1381
message_buffer_free()1382 void message_buffer_free()
1383 {
1384 std::free(ge.message_buffer);
1385 message_buffer(nullptr, 0);
1386 }
1387 #endif
1388 };
1389
1390 /** \brief Objects that can hold a vector */
1391 union geom_object {
1392 struct node nd;
1393 struct edge ed;
1394 union element el;
1395 };
1396
1397 /** \brief Objects that can have a key */
1398 union object_with_key {
1399 struct node nd;
1400 union element el;
1401 struct vector ve;
1402 union vertex vertex;
1403 struct edge edge;
1404 };
1405
1406 typedef struct
1407 {
1408 UINT VecReserv[MAXVECTORS][MAX_NDOF_MOD_32];
1409 UINT MatReserv[MAXCONNECTIONS][MAX_NDOF_MOD_32];
1410 /** \todo Is this used anywhere? */
1411 UINT VecConsistentStatus[MAXMATRICES][MAX_NDOF_MOD_32];
1412 /** \todo Is this used anywhere? */
1413 UINT VecCollectStatus[MAXMATRICES][MAX_NDOF_MOD_32];
1414 } DATA_STATUS;
1415
1416 struct grid {
1417
1418 /** \brief Object identification, various flags */
1419 UINT control;
1420
1421 /** \brief level + 32; needed for control word check not detecting HEAPFAULT */
1422 INT attribut;
1423
1424 /** \brief A word storing status information. This can be used also by the
1425 problem class, e.g. to store whether the grid level is assembled or not. */
1426 INT status;
1427
1428 /** \brief Level within the multigrid structure */
1429 INT level;
1430
1431 /** \brief Number of vertices */
1432 INT nVert[NS_DIM_PREFIX MAX_PRIOS];
1433
1434 /** \brief Number of nodes on this grid level */
1435 INT nNode[NS_DIM_PREFIX MAX_PRIOS];
1436
1437 /** \brief Number of elements on this grid level */
1438 INT nElem[NS_DIM_PREFIX MAX_PRIOS];
1439
1440 /** \brief Number of edges on this grid level */
1441 INT nEdge;
1442
1443 /** \brief Number of vectors on this grid level */
1444 INT nVector[NS_DIM_PREFIX MAX_PRIOS];
1445
1446 /** \brief Number of connections on this grid level */
1447 INT nCon;
1448
1449 DATA_STATUS data_status; /* memory management for vectors|matrix */
1450 /* status for consistent and collect */
1451 /* pointers */
1452 union element *elements[NS_DIM_PREFIX ELEMENT_LISTPARTS]; /* pointer to first element*/
1453 union element *lastelement[NS_DIM_PREFIX ELEMENT_LISTPARTS]; /*pointer to last element*/
1454 union vertex *vertices[NS_DIM_PREFIX VERTEX_LISTPARTS]; /* pointer to first vertex */
1455 union vertex *lastvertex[NS_DIM_PREFIX VERTEX_LISTPARTS]; /* pointer to last vertex */
1456 struct node *firstNode[NS_DIM_PREFIX NODE_LISTPARTS]; /* pointer to first node */
1457 struct node *lastNode[NS_DIM_PREFIX NODE_LISTPARTS]; /* pointer to last node */
1458 VECTOR *firstVector[NS_DIM_PREFIX VECTOR_LISTPARTS]; /* pointer to first vector */
1459 VECTOR *lastVector[NS_DIM_PREFIX VECTOR_LISTPARTS]; /* pointer to last vector */
1460
1461 struct grid *coarser, *finer; /* coarser and finer grids */
1462 struct multigrid *mg; /* corresponding multigrid structure */
1463
1464 const PPIF::PPIFContext& ppifContext() const;
1465
1466 #ifdef ModelP
1467 const DDD::DDDContext& dddContext() const;
1468 DDD::DDDContext& dddContext();
1469 #endif
1470 };
1471
1472 struct multigrid {
1473
1474 /** \brief env item */
1475 NS_PREFIX ENVDIR v;
1476
1477 /** \brief Multigrid status word */
1478 INT status;
1479
1480 /** \brief used for identification */
1481 INT magic_cookie;
1482
1483 /** \brief count objects in that multigrid */
1484 INT vertIdCounter;
1485
1486 /** \brief count objects in that multigrid */
1487 INT nodeIdCounter;
1488
1489 /** \brief count objects in that multigrid */
1490 INT elemIdCounter;
1491
1492 /** \brief count objects in that multigrid */
1493 INT edgeIdCounter;
1494
1495 #ifndef ModelP
1496 /** \brief Count vector objects in that multigrid */
1497 INT vectorIdCounter;
1498 #endif
1499
1500 /** \brief depth of the element tree */
1501 INT topLevel;
1502
1503 /** \brief level we are working on */
1504 INT currentLevel;
1505
1506 /** \brief last level with complete surface */
1507 INT fullrefineLevel;
1508
1509 /** \brief pointer to BndValProblem */
1510 BVP *theBVP;
1511
1512 /** \brief description of BVP-properties */
1513 BVP_DESC theBVPD;
1514
1515 /** \brief pointer to format definitions */
1516 std::unique_ptr<format> theFormat;
1517
1518 /** \brief associated heap structure */
1519 NS_PREFIX HEAP *theHeap;
1520
1521 /** \brief max nb of properties used in elements*/
1522 INT nProperty;
1523
1524
1525 /** \brief memory management for vectors|matrix
1526 * status for consistent and collect */
1527 DATA_STATUS data_status;
1528
1529 /* pointers */
1530 /** \brief pointers to the grids */
1531 struct grid *grids[MAXLEVEL];
1532
1533 /** @{ Helper structures to for an O(n) InsertElement */
1534 /** \brief List of pointers to face nodes,
1535 used as an Id of a face
1536 */
1537 struct FaceNodes : public std::array<node*,MAX_CORNERS_OF_SIDE>
1538 {
1539 using array<node*,MAX_CORNERS_OF_SIDE>::array;
1540 };
1541 /** \brief Hasher for FaceNodes */
1542 struct FaceHasher {
1543 std::hash<node*> hasher;
operatormultigrid::FaceHasher1544 std::size_t operator() (const FaceNodes& key) const {
1545 return std::accumulate(key.begin(), key.end(),
1546 144451, [&](auto a, auto b){
1547 return hasher(a+b);
1548 });
1549 }
1550 };
1551 /** \brief hash-map used for an O(1) search of the neighboring element
1552 during InsertElement
1553 */
1554 std::unordered_map<FaceNodes,
1555 std::pair<element *,int>,
1556 FaceHasher> facemap;
1557 /** @} */
1558
1559 /* i/o handling */
1560 /** \brief 1 if multigrid saved */
1561 INT saved;
1562
1563 /** \brief filename if saved */
1564 char filename[NS_PREFIX NAMESIZE];
1565
1566 /** \brief coarse grid complete */
1567 INT CoarseGridFixed;
1568
1569 /** \brief coarse grid MarkKey for SIMPLE_HEAP Mark/Release */
1570 INT MarkKey;
1571
ppifContextmultigrid1572 const PPIF::PPIFContext& ppifContext() const
1573 { return *ppifContext_; }
1574
1575 std::shared_ptr<PPIF::PPIFContext> ppifContext_;
1576
1577 #ifdef ModelP
dddContextmultigrid1578 const DDD::DDDContext& dddContext() const
1579 { return *dddContext_; }
1580
dddContextmultigrid1581 DDD::DDDContext& dddContext()
1582 { return *dddContext_; }
1583
1584 std::shared_ptr<DDD::DDDContext> dddContext_;
1585 #endif
1586 };
1587
1588 /****************************************************************************/
1589 /* */
1590 /* typedef for structs */
1591 /* */
1592 /****************************************************************************/
1593
1594 /* geometrical part */
1595 typedef struct format FORMAT;
1596
1597 typedef union vertex VERTEX;
1598 typedef struct elementlist ELEMENTLIST;
1599 typedef struct node NODE;
1600 typedef union element ELEMENT;
1601 typedef struct link LINK;
1602 typedef struct edge EDGE;
1603 typedef union geom_object GEOM_OBJECT;
1604 typedef struct grid GRID;
1605 typedef struct multigrid MULTIGRID;
1606 typedef union object_with_key KEY_OBJECT;
1607
1608 /****************************************************************************/
1609 /* */
1610 /* algebraic dependency for vector ordering */
1611 /* */
1612 /****************************************************************************/
1613
1614 typedef INT (*DependencyProcPtr)(GRID *, const char *);
1615
1616 struct AlgebraicDependency {
1617
1618 /* fields for enironment list variable */
1619 NS_PREFIX ENVVAR v;
1620
1621 DependencyProcPtr DependencyProc; /* pointer to dependency function */
1622 };
1623
1624 typedef struct AlgebraicDependency ALG_DEP;
1625
1626 /****************************************************************************/
1627 /* */
1628 /* find cut for vector ordering */
1629 /* */
1630 /****************************************************************************/
1631
1632 typedef VECTOR *(*FindCutProcPtr)(GRID *, VECTOR *, INT *);
1633
1634 typedef struct {
1635
1636 /* fields for enironment list variable */
1637 NS_PREFIX ENVVAR v;
1638
1639 FindCutProcPtr FindCutProc; /* pointer to find cut function */
1640
1641 } FIND_CUT;
1642
1643 /****************************************************************************/
1644 /* */
1645 /* dynamic management of control words */
1646 /* */
1647 /****************************************************************************/
1648
1649 /** @name status of control word */
1650 /*@{*/
1651 #define CW_FREE 0
1652 #define CW_USED 1
1653 /*@}*/
1654
1655
1656 /** @name Status of control entry */
1657 /*@{*/
1658 #define CE_FREE 0
1659 #define CE_USED 1
1660 #define CE_LOCKED 1
1661 /*@}*/
1662
1663 /** @name Initializer macros for control entry and word predefines */
1664 /*@{*/
1665 #define CW_INIT(used,cw,objs) {used, STR(cw), cw ## CW, cw ## OFFSET,objs}
1666 #define CW_INIT_UNUSED {CW_FREE,0,0,0}
1667 #define CE_INIT(mode,cw,ce,objs) {mode, STR(ce), cw ## CW, ce ## CE, ce ## SHIFT, ce ## LEN, objs}
1668 #define CE_INIT_UNUSED {CE_FREE, 0, 0, 0, 0, 0, 0}
1669 /*@}*/
1670
1671 /* general query macros */
1672
1673 /* dynamic control words */
1674 /*#define _DEBUG_CW_*/
1675 #if (defined _DEBUG_CW_) && \
1676 !(defined __COMPILE_CW__) /* to avoid infinite recursion during ReadCW */
1677
1678 /* map cw read/write to functions */
1679 #define CW_READ(p,ce) ReadCW(p,ce)
1680 #define CW_READ_STATIC(p,s,t) ReadCW(p,s ## CE)
1681 #define CW_WRITE(p,ce,n) WriteCW(p,ce,n)
1682 #define CW_WRITE_STATIC(p,s,t,n) WriteCW(p,s ## CE,n)
1683
1684 #else /* _DEBUG_CW_ */
1685
1686 #define ControlWord(p,ce) (((UINT *)(p))[control_entries[ce].offset_in_object])
1687
1688 #ifndef __T3E__
1689 #define CW_READ(p,ce) ((ControlWord(p,ce) & control_entries[ce].mask)>>control_entries[ce].offset_in_word)
1690 #endif
1691
1692 /* very special hack */
1693 #ifdef __T3E__
1694 #define CW_READ(p,ce) ((int)((ControlWord(p,ce) & control_entries[ce].mask)>>control_entries[ce].offset_in_word) )
1695 #endif
1696
1697 #define CW_WRITE(p,ce,n) ControlWord(p,ce) = (ControlWord(p,ce)&control_entries[ce].xor_mask)|(((n)<<control_entries[ce].offset_in_word)&control_entries[ce].mask)
1698
1699 /* static control words */
1700 #define StaticControlWord(p,t) (((UINT *)(p))[t ## OFFSET])
1701 #define StaticControlWordMask(s) ((POW2(s ## LEN) - 1) << s ## SHIFT)
1702
1703 #ifndef __T3E__
1704 #define CW_READ_STATIC(p,s,t) \
1705 ((StaticControlWord(p,t) & StaticControlWordMask(s)) >> s ## SHIFT)
1706 #endif
1707
1708 /* very special hack */
1709 #ifdef __T3E__
1710 #define CW_READ_STATIC(p,s,t) \
1711 ((int) ((StaticControlWord(p,t) & StaticControlWordMask(s)) >> s ## SHIFT))
1712 #endif
1713
1714 #define CW_WRITE_STATIC(p,s,t,n) \
1715 StaticControlWord(p,t) = \
1716 (StaticControlWord(p,t) & (~StaticControlWordMask(s))) | \
1717 (((n) << s ## SHIFT) & StaticControlWordMask(s))
1718
1719 #endif /* _DEBUG_CW_ */
1720
1721 /** \brief Enumeration list of all control words of gm.h */
1722 enum GM_CW {
1723 VECTOR_CW,
1724 MATRIX_CW,
1725 VERTEX_CW,
1726 NODE_CW,
1727 LINK_CW,
1728 EDGE_CW,
1729 ELEMENT_CW,
1730 FLAG_CW,
1731 PROPERTY_CW,
1732 GRID_CW,
1733 GRID_STATUS_CW,
1734 MULTIGRID_STATUS_CW,
1735
1736 GM_N_CW
1737 };
1738
1739 /** \brief Enumeration list of all control entry of gm.h */
1740 enum GM_CE {
1741 VTYPE_CE,
1742 VOTYPE_CE,
1743 VPART_CE,
1744 VCOUNT_CE,
1745 VECTORSIDE_CE,
1746 VCLASS_CE,
1747 VDATATYPE_CE,
1748 VNCLASS_CE,
1749 VNEW_CE,
1750 VCCUT_CE,
1751 VCCOARSE_CE,
1752 NEW_DEFECT_CE,
1753 VACTIVE_CE,
1754 FINE_GRID_DOF_CE,
1755 MOFFSET_CE,
1756 MROOTTYPE_CE,
1757 MDESTTYPE_CE,
1758 MDIAG_CE,
1759 MSIZE_CE,
1760 MNEW_CE,
1761 CEXTRA_CE,
1762 MDOWN_CE,
1763 MUP_CE,
1764 MLOWER_CE,
1765 MUPPER_CE,
1766 MACTIVE_CE,
1767 OBJ_CE,
1768 USED_CE,
1769 TAG_CE,
1770 LEVEL_CE,
1771 THEFLAG_CE,
1772 MOVE_CE,
1773 MOVED_CE,
1774 ONEDGE_CE,
1775 ONSIDE_CE,
1776 ONNBSIDE_CE,
1777 NOOFNODE_CE,
1778 NSUBDOM_CE,
1779 NTYPE_CE,
1780 NPROP_CE,
1781 MODIFIED_CE,
1782 NCLASS_CE,
1783 NNCLASS_CE,
1784 LOFFSET_CE,
1785 NO_OF_ELEM_CE,
1786 AUXEDGE_CE,
1787 EDGENEW_CE,
1788 EDSUBDOM_CE,
1789 ECLASS_CE,
1790 NSONS_CE,
1791 NEWEL_CE,
1792 SUBDOMAIN_CE,
1793 NODEORD_CE,
1794 PROP_CE,
1795 #ifdef ModelP
1796 XFERVECTOR_CE,
1797 XFERMATX_CE,
1798 #endif
1799
1800 GM_N_CE
1801 };
1802
1803 enum LV_MODIFIERS {
1804
1805 LV_VO_INFO = (1<<1), /* vector object related info */
1806 LV_POS = (1<<2) /* position vector */
1807 };
1808
1809 enum LV_ID_TYPES {
1810 LV_ID,
1811 LV_GID,
1812 LV_KEY
1813 };
1814
1815 /****************************************************************************/
1816 /* */
1817 /* Macro definitions for algebra structures */
1818 /* */
1819 /* */
1820 /* Use of the control word: */
1821 /* */
1822 /* macro name|bits |V|M|use */
1823 /* */
1824 /* all objects: */
1825 /* */
1826 /* vectors: */
1827 /* VOTYPE |0 - 1 |*| | node-,edge-,side- or elemvector */
1828 /* VCFLAG |3 |*| | flag for general use */
1829 /* VCUSED |4 |*| | flag for general use */
1830 /* VCOUNT |5-6 |*| | The number of elements that reference */
1831 /* the vector (if it is a side vector) */
1832 /* VECTORSIDE|7 - 9 |*| | nb of side the side vect corr. to (in object elem)*/
1833 /* VCLASS |11-12 |*| | class of v. (3: if corr. to red/green elem) */
1834 /* (2: if corr. to first algebraic nb.) */
1835 /* (1: if corr. to second algebraic nb.) */
1836 /* VDATATYPE |13-16 |*| | data type used bitwise */
1837 /* VNCLASS |17-18 |*| | type of elem on finer grid the v. lies geom. in: */
1838 /* 0: no elem on finer grid */
1839 /* 1: elem of 'second alg. nbhood' only */
1840 /* 2: elem of 'first alg. nbhood' only */
1841 /* 3: red or green elem */
1842 /* VNEW |19 |*| | 1 if vector is new */
1843 /* VCNEW |20 |*| | 1 if vector has a new connection */
1844 /* VACTIVE |24 |*| | 1 if vector is active inside a smoother */
1845 /* VCCUT |26 |*| | */
1846 /* VTYPE |27-28 |*| | abstract vector type */
1847 /* VPART |29-30 |*| | domain part */
1848 /* VCCOARSE |31 |*| | indicate algebraic part of VECTOR-MATRIX graph */
1849 /* */
1850 /* matrices: */
1851 /* MOFFSET |0 | |*| 0 if first matrix in connection else 1 */
1852 /* MROOTTYPE |1 - 2 | |*| VTYPE of root vector */
1853 /* MDESTTYPE |3 - 4 | |*| VTYPE of destination vector */
1854 /* MDIAG |5 | |*| 1 if diagonal matrix element */
1855 /* MNEW |6 | |*| ??? */
1856 /* CEXTRA |7 | |*| 1 if is extra connection */
1857 /* MDOWN |8 | |*| ??? */
1858 /* MUP |9 | |*| ??? */
1859 /* MLOWER |10 | |*| 1 if matrix belongs to lower triangular part */
1860 /* MUPPER |11 | |*| 1 if matrix belongs to upper triangular part */
1861 /* UG_MSIZE |12-25 | |*| size of the matrix in bytes */
1862 /* MUSED |12 | |*| general purpose flag */
1863 /* MNEW |28 | |*| 1 if matrix/connection is new */
1864 /* */
1865 /****************************************************************************/
1866
1867 /****************************************************************************/
1868 /* */
1869 /* general macros */
1870 /* */
1871 /****************************************************************************/
1872
1873 /* macros to calculate from a coordinate (2D/3D) a hopefully unique ID */
1874 #define SIGNIFICANT_DIGITS(d,exp_ptr) (ceil(frexp((d),(exp_ptr))*1e5))
1875
1876 /* the idea to calculate from a 2d/3D position a (hopefully) unique key:
1877 add the weighted significant digits of the coordinates; the weights
1878 may not have a common divisor to ensure uniqueness of the result;
1879 take from this again the sigificant digits */
1880 #ifdef __TWODIM__
1881 #define COORDINATE_TO_KEY(coord,dummy_int_ptr) ((INT)(SIGNIFICANT_DIGITS((SIGNIFICANT_DIGITS((coord)[0],(dummy_int_ptr))*1.246509423749342 + \
1882 SIGNIFICANT_DIGITS((coord)[1],(dummy_int_ptr))*PI)\
1883 , (dummy_int_ptr))))
1884 #endif
1885
1886 #ifdef __THREEDIM__
1887 #define COORDINATE_TO_KEY(coord,dummy_int_ptr) ((INT)(SIGNIFICANT_DIGITS((SIGNIFICANT_DIGITS((coord)[0],(dummy_int_ptr))*1.246509423749342 + \
1888 SIGNIFICANT_DIGITS((coord)[1],(dummy_int_ptr))*PI + \
1889 SIGNIFICANT_DIGITS((coord)[2],(dummy_int_ptr))*0.76453456834568356936598)\
1890 , (dummy_int_ptr))))
1891 #endif
1892
1893 /****************************************************************************/
1894 /* */
1895 /* macros for VECTORs */
1896 /* */
1897 /****************************************************************************/
1898
1899 /* control word offset */
1900 #define VECTOR_OFFSET 0
1901
1902 /* predefined control word entries */
1903 #define VOTYPE_SHIFT 0
1904 #define VOTYPE_LEN 2
1905 #define VOTYPE(p) CW_READ_STATIC(p,VOTYPE_,VECTOR_)
1906 #define SETVOTYPE(p,n) CW_WRITE_STATIC(p,VOTYPE_,VECTOR_,n)
1907 #if (MAXVOBJECTS > POW2(VOTYPE_LEN))
1908 #error *** VOTYPE_LEN too small ***
1909 #endif
1910
1911 #define VTYPE_SHIFT 2
1912 #define VTYPE_LEN 2
1913 #define VTYPE(p) (enum VectorType)CW_READ_STATIC(p,VTYPE_,VECTOR_)
1914 #define SETVTYPE(p,n) CW_WRITE_STATIC(p,VTYPE_,VECTOR_,n)
1915 #if (MAXVTYPES > POW2(VTYPE_LEN))
1916 #error *** VTYPE_LEN too small ***
1917 #endif
1918
1919 #define VDATATYPE_SHIFT 4
1920 #define VDATATYPE_LEN 4
1921 #define VDATATYPE(p) CW_READ_STATIC(p,VDATATYPE_,VECTOR_)
1922 #define SETVDATATYPE(p,n) CW_WRITE_STATIC(p,VDATATYPE_,VECTOR_,n)
1923 #if (MAXVTYPES > VDATATYPE_LEN)
1924 #error *** VDATATYPE_LEN too small ***
1925 #endif
1926
1927 #define VCLASS_SHIFT 8
1928 #define VCLASS_LEN 2
1929 #define VCLASS(p) CW_READ_STATIC(p,VCLASS_,VECTOR_)
1930 #define SETVCLASS(p,n) CW_WRITE_STATIC(p,VCLASS_,VECTOR_,n)
1931
1932 #define VNCLASS_SHIFT 10
1933 #define VNCLASS_LEN 2
1934 #define VNCLASS(p) CW_READ_STATIC(p,VNCLASS_,VECTOR_)
1935 #define SETVNCLASS(p,n) CW_WRITE_STATIC(p,VNCLASS_,VECTOR_,n)
1936
1937 #define VNEW_SHIFT 12
1938 #define VNEW_LEN 1
1939 #define VNEW(p) CW_READ_STATIC(p,VNEW_,VECTOR_)
1940 #define SETVNEW(p,n) CW_WRITE_STATIC(p,VNEW_,VECTOR_,n)
1941
1942 #define VCCUT_SHIFT 13
1943 #define VCCUT_LEN 1
1944 #define VCCUT(p) CW_READ_STATIC(p,VCCUT_,VECTOR_)
1945 #define SETVCCUT(p,n) CW_WRITE_STATIC(p,VCCUT_,VECTOR_,n)
1946
1947 #define VCOUNT_SHIFT 14
1948 #define VCOUNT_LEN 2
1949 #define VCOUNT(p) CW_READ_STATIC(p,VCOUNT_,VECTOR_)
1950 #define SETVCOUNT(p,n) CW_WRITE_STATIC(p,VCOUNT_,VECTOR_,n)
1951
1952 #define VECTORSIDE_SHIFT 16
1953 #define VECTORSIDE_LEN 3
1954 #define VECTORSIDE(p) CW_READ_STATIC(p,VECTORSIDE_,VECTOR_)
1955 #define SETVECTORSIDE(p,n) CW_WRITE_STATIC(p,VECTORSIDE_,VECTOR_,n)
1956
1957 #define VCCOARSE_SHIFT 19
1958 #define VCCOARSE_LEN 1
1959 #define VCCOARSE(p) CW_READ_STATIC(p,VCCOARSE_,VECTOR_)
1960 #define SETVCCOARSE(p,n) CW_WRITE_STATIC(p,VCCOARSE_,VECTOR_,n)
1961
1962 #define FINE_GRID_DOF_SHIFT 20
1963 #define FINE_GRID_DOF_LEN 1
1964 #define FINE_GRID_DOF(p) CW_READ_STATIC(p,FINE_GRID_DOF_,VECTOR_)
1965 #define SETFINE_GRID_DOF(p,n) CW_WRITE_STATIC(p,FINE_GRID_DOF_,VECTOR_,n)
1966
1967 #define NEW_DEFECT_SHIFT 21
1968 #define NEW_DEFECT_LEN 1
1969 #define NEW_DEFECT(p) CW_READ_STATIC(p,NEW_DEFECT_,VECTOR_)
1970 #define SETNEW_DEFECT(p,n) CW_WRITE_STATIC(p,NEW_DEFECT_,VECTOR_,n)
1971
1972 #ifdef ModelP
1973 #define XFERVECTOR_SHIFT 20
1974 #define XFERVECTOR_LEN 2
1975 #define XFERVECTOR(p) CW_READ(p,XFERVECTOR_CE)
1976 #define SETXFERVECTOR(p,n) CW_WRITE(p,XFERVECTOR_CE,n)
1977 #endif /* ModelP */
1978
1979 #define VPART_SHIFT 22
1980 #define VPART_LEN 2
1981 #define VPART(p) CW_READ_STATIC(p,VPART_,VECTOR_)
1982 #define SETVPART(p,n) CW_WRITE_STATIC(p,VPART_,VECTOR_,n)
1983 #if (MAXDOMPARTS > POW2(VPART_LEN))
1984 #error *** VPART_LEN too small ***
1985 #endif
1986
1987 #define VACTIVE_SHIFT 24
1988 #define VACTIVE_LEN 1
1989 #define VACTIVE(p) CW_READ_STATIC(p,VACTIVE_,VECTOR_)
1990 #define SETVACTIVE(p,n) CW_WRITE_STATIC(p,VACTIVE_,VECTOR_,n)
1991
1992 #define VCFLAG(p) THEFLAG(p)
1993 #define SETVCFLAG(p,n) SETTHEFLAG(p,n)
1994
1995 #define VCUSED(p) USED(p)
1996 #define SETVCUSED(p,n) SETUSED(p,n)
1997
1998 #define VOBJECT(v) ((v)->object)
1999 #ifdef ModelP
2000 #define PPREDVC(p,v) (((v)==PRIO_FIRSTVECTOR(p,PrioMaster)) ? \
2001 PRIO_LASTVECTOR(p,PrioBorder) : (v)->pred)
2002 #else
2003 #define PPREDVC(p,v) ((v)->pred)
2004 #endif
2005 #define PREDVC(v) ((v)->pred)
2006 #define SUCCVC(v) ((v)->succ)
2007 #define VINDEX(v) ((v)->index)
2008 #define V_IN_DATATYPE(v,dt) (VDATATYPE(v) & (dt))
2009 #define VSTART(v) ((v)->start)
2010 #define VVALUE(v,n) ((v)->value[n])
2011 #define VVALUEPTR(v,n) (&((v)->value[n]))
2012 #define VMYNODE(v) ((NODE*)((v)->object))
2013 #define VMYEDGE(v) ((EDGE*)((v)->object))
2014 #define VMYELEMENT(v) ((ELEMENT*)((v)->object))
2015 #define VUP(p) LoWrd(VINDEX(p))
2016 #define SETVUP(p,n) SetLoWrd(VINDEX(p),n)
2017 #define VDOWN(p) HiWrd(VINDEX(p))
2018 #define SETVDOWN(p,n) SetHiWrd(VINDEX(p),n)
2019
2020 /****************************************************************************/
2021 /* */
2022 /* macros for MATRIXs */
2023 /* */
2024 /****************************************************************************/
2025
2026 /* control word offset */
2027 #define MATRIX_OFFSET 0
2028
2029 #define MOFFSET_SHIFT 0
2030 #define MOFFSET_LEN 1
2031 #define MOFFSET(p) CW_READ_STATIC(p,MOFFSET_,MATRIX_)
2032 #define SETMOFFSET(p,n) CW_WRITE_STATIC(p,MOFFSET_,MATRIX_,n)
2033
2034 #define MROOTTYPE_SHIFT 1
2035 #define MROOTTYPE_LEN 2
2036 #define MROOTTYPE(p) CW_READ_STATIC(p,MROOTTYPE_,MATRIX_)
2037 #define SETMROOTTYPE(p,n) CW_WRITE_STATIC(p,MROOTTYPE_,MATRIX_,n)
2038 #if (MAXVTYPES > POW2(MROOTTYPE_LEN))
2039 #error *** MROOTTYPE_LEN too small ***
2040 #endif
2041
2042 #define MDESTTYPE_SHIFT 3
2043 #define MDESTTYPE_LEN 2
2044 #define MDESTTYPE(p) CW_READ_STATIC(p,MDESTTYPE_,MATRIX_)
2045 #define SETMDESTTYPE(p,n) CW_WRITE_STATIC(p,MDESTTYPE_,MATRIX_,n)
2046 #if (MAXVTYPES > POW2(MDESTTYPE_LEN))
2047 #error *** MDESTTYPE_LEN too small ***
2048 #endif
2049
2050 #define MDIAG_SHIFT 5
2051 #define MDIAG_LEN 1
2052 #define MDIAG(p) CW_READ_STATIC(p,MDIAG_,MATRIX_)
2053 #define SETMDIAG(p,n) CW_WRITE_STATIC(p,MDIAG_,MATRIX_,n)
2054
2055 #define MNEW_SHIFT 6
2056 #define MNEW_LEN 1
2057 #define MNEW(p) CW_READ_STATIC(p,MNEW_,MATRIX_)
2058 #define SETMNEW(p,n) CW_WRITE_STATIC(p,MNEW_,MATRIX_,n)
2059
2060 #define CEXTRA_SHIFT 7
2061 #define CEXTRA_LEN 1
2062 #define CEXTRA(p) CW_READ_STATIC(p,CEXTRA_,MATRIX_)
2063 #define SETCEXTRA(p,n) CW_WRITE_STATIC(p,CEXTRA_,MATRIX_,n)
2064
2065 #define MDOWN_SHIFT 8
2066 #define MDOWN_LEN 1
2067 #define MDOWN(p) CW_READ_STATIC(p,MDOWN_,MATRIX_)
2068 #define SETMDOWN(p,n) CW_WRITE_STATIC(p,MDOWN_,MATRIX_,n)
2069
2070 #define MUP_SHIFT 9
2071 #define MUP_LEN 1
2072 #define MUP(p) CW_READ_STATIC(p,MUP_,MATRIX_)
2073 #define SETMUP(p,n) CW_WRITE_STATIC(p,MUP_,MATRIX_,n)
2074
2075 #define MLOWER_SHIFT 10
2076 #define MLOWER_LEN 1
2077 #define MLOWER(p) CW_READ_STATIC(p,MLOWER_,MATRIX_)
2078 #define SETMLOWER(p,n) CW_WRITE_STATIC(p,MLOWER_,MATRIX_,n)
2079
2080 #define MUPPER_SHIFT 11
2081 #define MUPPER_LEN 1
2082 #define MUPPER(p) CW_READ_STATIC(p,MUPPER_,MATRIX_)
2083 #define SETMUPPER(p,n) CW_WRITE_STATIC(p,MUPPER_,MATRIX_,n)
2084
2085 #define MACTIVE_SHIFT 12
2086 #define MACTIVE_LEN 1
2087 #define MACTIVE(p) CW_READ_STATIC(p,MACTIVE_,MATRIX_)
2088 #define SETMACTIVE(p,n) CW_WRITE_STATIC(p,MACTIVE_,MATRIX_,n)
2089
2090 #define MSIZE_SHIFT 13
2091 #define MSIZE_LEN 12
2092 #ifndef __XXL_MSIZE__
2093 #define MSIZEMAX (POW2(MSIZE_LEN)-1)
2094 #define UG_MSIZE(p) (CW_READ(p,MSIZE_CE)+sizeof(MATRIX)-sizeof(DOUBLE))
2095 #define SETMSIZE(p,n) CW_WRITE(p,MSIZE_CE,(n-sizeof(MATRIX)+sizeof(DOUBLE)))
2096 #else
2097 #define MSIZEMAX 10000000
2098 #define UG_MSIZE(p) ((p)->xxl_msize)
2099 #define SETMSIZE(p,n) (p)->xxl_msize = (n)
2100 #endif
2101
2102 #define MTYPE(p) (MDIAG(p) ? (MAXMATRICES+MROOTTYPE(p)) : (MROOTTYPE(p)*MAXVECTORS+MDESTTYPE(p)))
2103
2104 #define MUSED(p) USED(p)
2105 #define SETMUSED(p,n) SETUSED(p,n)
2106
2107 #ifdef ModelP
2108 #define XFERMATX_SHIFT 25
2109 #define XFERMATX_LEN 2
2110 #define XFERMATX(p) CW_READ(p,XFERMATX_CE)
2111 #define SETXFERMATX(p,n) CW_WRITE(p,XFERMATX_CE,n)
2112 #endif
2113
2114 #define MINC(m) ((MATRIX*)(((char *)(m))+UG_MSIZE(m)))
2115 #define MDEC(m) ((MATRIX*)(((char *)(m))-UG_MSIZE(m)))
2116 #define MNEXT(m) ((m)->next)
2117 #define MDEST(m) ((m)->vect)
2118 #define MADJ(m) ((MDIAG(m)) ? (m) : ((MOFFSET(m)) ? (MDEC(m)) : (MINC(m))))
2119 #define MROOT(m) MDEST(MADJ(m))
2120 #define MMYCON(m) ((MOFFSET(m)) ? (MDEC(m)) : (m))
2121 #define MVALUE(m,n) ((m)->value[n])
2122 #define MVALUEPTR(m,n) (&((m)->value[n]))
2123 #define MDESTINDEX(m) ((m)->vect->index)
2124 #define MSTRONG(p) (MDOWN(p) && MUP(p))
2125
2126 /****************************************************************************/
2127 /* */
2128 /* macros for CONNECTIONs */
2129 /* */
2130 /****************************************************************************/
2131
2132 #define CMATRIX0(m) (m)
2133 #define CMATRIX1(m) ((MDIAG(m)) ? (NULL) : (MINC(m)))
2134 #define SETCUSED(c,n) {SETMUSED(CMATRIX0(c),n); SETMUSED(MADJ(CMATRIX0(c)),n);}
2135
2136 /****************************************************************************/
2137 /* */
2138 /* Macro definitions for geometric objects */
2139 /* */
2140 /* */
2141 /* Use of the control word: */
2142 /* */
2143 /* macro name|bits |V|N|L|E|V|M| use */
2144 /* C */
2145 /* all objects: */
2146 /* TAG |18-20 | | | |*| | |general purpose tag field */
2147 /* LEVEL |21-25 |*|*| |*| | |level of a node/element (imp. for copies) */
2148 /* THEFLAG |26 |*|*|*|*| | |general purp., leave them as you found 'em*/
2149 /* USED |27 |*|*|*|*| | |object visited, leave them as you found 'em*/
2150 /* OBJT |28-31 |*|*|*|*| | |object type identification */
2151
2152 /* */
2153 /* vertices: */
2154 /* MOVED |0 |*| | | | | |boundary vertex not lying on edge midpoint */
2155 /* MOVE |1-2 |*| | | | | |vertex can be moved on a 0(1,2,3) dim subsp*/
2156 /* ONEDGE |3 - 6 |*| | | | | |no. of edge in father element */
2157 /* ONSIDE |3 - 5 |*| | | | | |no. of side in father element */
2158 /* ONNBSIDE |6 - 8 |*| | | | | |no. of side in the neigbor of the father */
2159 /* NOOFNODE |9 -13 |*| | | | | |??? */
2160 /* */
2161 /* nodes: */
2162 /* NSUBDOM |0-3 | |*| | | | |subdomain id */
2163 /* MODIFIED |6 | |*| | | | |1 if node must be assembled */
2164 /* N_OUTFLOW |0-7 | */
2165 /* N_INFLOW |8-15 | */
2166 /* */
2167 /* links and edges: */
2168 /* LOFFSET |0 | | |*| | | |position of link in links array */
2169 /* EDGENEW |1 | | |*| | | |status of edge */
2170 /* NOOFELEM |2-8 | | |*| | | |nb. of elem. the edge is part of */
2171 /* AUXEDGE |9 | */
2172 /* EDSUBDOM |12-17 | | | |*| | |subdomain of edge if inner edge, 0 else */
2173 /* */
2174 /* elements: */
2175 /* ECLASS |8-9 | | | |*| | |element class from enumeration type */
2176 /* NSONS |10-13 | | | |*| | |number of sons */
2177 /* NEWEL |14 | | | |*| | |element just created */
2178 /* VSIDES |11-14 | | | |*| | |viewable sides */
2179 /* NORDER |15-19 | | | |*| | |view position order of the nodes */
2180 /* CUTMODE |26-27 | | | |*| | |elem intersects cutplane or... */
2181 /* */
2182 /****************************************************************************/
2183
2184 /* object identification */
2185 enum GM_OBJECTS {
2186
2187 MGOBJ, /*!< Multigrid object */
2188 IVOBJ, /*!< Inner vertex */
2189 BVOBJ, /*!< Boundary vertex */
2190 IEOBJ, /*!< Inner element */
2191 BEOBJ, /*!< Boundary element */
2192 EDOBJ, /*!< Edge object */
2193 NDOBJ, /*!< Node object */
2194 GROBJ, /*!< Grid object */
2195
2196 /* object numbers for algebra */
2197 VEOBJ, /*!< Vector object */
2198 MAOBJ, /*!< Matrix object */
2199
2200 NPREDEFOBJ, /*!< Number of predefined objects */
2201
2202 NOOBJ = -1 /*!< No object */
2203 };
2204 #define LIOBJ EDOBJ /* link and edge are identified */
2205 #define COOBJ MAOBJ /* connection and matrix are identified */
2206
2207 /****************************************************************************/
2208 /* */
2209 /* general macros */
2210 /* */
2211 /****************************************************************************/
2212
2213 /* control word offset */
2214 #define GENERAL_CW NODE_CW /* any of the geom objects */
2215 #define GENERAL_OFFSET 0
2216
2217 #define OBJ_SHIFT 28
2218 #define OBJ_LEN 4
2219 #define OBJT(p) (enum GM_OBJECTS)CW_READ_STATIC(p,OBJ_,GENERAL_)
2220 #define SETOBJT(p,n) CW_WRITE_STATIC(p,OBJ_,GENERAL_,n)
2221 #define OBJT_MAX (POW2(OBJ_LEN)-1)
2222
2223 #define USED_SHIFT 27
2224 #define USED_LEN 1
2225 #define USED(p) CW_READ_STATIC(p,USED_,GENERAL_)
2226 #define SETUSED(p,n) CW_WRITE_STATIC(p,USED_,GENERAL_,n)
2227
2228 #define THEFLAG_SHIFT 26
2229 #define THEFLAG_LEN 1
2230 #define THEFLAG(p) CW_READ_STATIC(p,THEFLAG_,GENERAL_)
2231 #define SETTHEFLAG(p,n) CW_WRITE_STATIC(p,THEFLAG_,GENERAL_,n)
2232
2233 #define LEVEL_SHIFT 21
2234 #define LEVEL_LEN 5
2235 #define LEVEL(p) CW_READ_STATIC(p,LEVEL_,GENERAL_)
2236 #define SETLEVEL(p,n) CW_WRITE_STATIC(p,LEVEL_,GENERAL_,n)
2237
2238 #define TAG_SHIFT 18
2239 #define TAG_LEN 3
2240 #define TAG(p) CW_READ_STATIC(p,TAG_,GENERAL_)
2241 #define SETTAG(p,n) CW_WRITE_STATIC(p,TAG_,GENERAL_,n)
2242
2243 #define REF2TAG(n) (reference2tag[n])
2244
2245 #define CTRL(p) (*((UINT *)(p)))
2246 #define ID(p) (((INT *)(p))[1])
2247
2248 /****************************************************************************/
2249 /* */
2250 /* macros for vertices */
2251 /* */
2252 /****************************************************************************/
2253
2254 /* control word offset */
2255 #define VERTEX_OFFSET 0
2256
2257 #define MOVE_SHIFT 1
2258 #define MOVE_LEN 2
2259 #define MOVE(p) CW_READ_STATIC(p,MOVE_,VERTEX_)
2260 #define SETMOVE(p,n) CW_WRITE_STATIC(p,MOVE_,VERTEX_,n)
2261
2262 #define MOVED_SHIFT 0
2263 #define MOVED_LEN 1
2264 #define MOVED(p) CW_READ_STATIC(p,MOVED_,VERTEX_)
2265 #define SETMOVED(p,n) CW_WRITE_STATIC(p,MOVED_,VERTEX_,n)
2266
2267 #define ONEDGE_SHIFT 3
2268 #define ONEDGE_LEN 4
2269 #define ONEDGE(p) CW_READ_STATIC(p,ONEDGE_,VERTEX_)
2270 #define SETONEDGE(p,n) CW_WRITE_STATIC(p,ONEDGE_,VERTEX_,n)
2271
2272 /* the following two overlap with ONEDGE */
2273 #define ONSIDE_SHIFT 3
2274 #define ONSIDE_LEN 3
2275 #define ONSIDE(p) CW_READ_STATIC(p,ONSIDE_,VERTEX_)
2276 #define SETONSIDE(p,n) CW_WRITE_STATIC(p,ONSIDE_,VERTEX_,n)
2277
2278 #define ONNBSIDE_SHIFT 6
2279 #define ONNBSIDE_LEN 3
2280 #define ONNBSIDE(p) CW_READ_STATIC(p,ONNBSIDE_,VERTEX_)
2281 #define SETONNBSIDE(p,n) CW_WRITE_STATIC(p,ONNBSIDE_,VERTEX_,n)
2282
2283 #define NOOFNODE_SHIFT 9
2284 #define NOOFNODE_LEN 5
2285 #define NOOFNODEMAX POW2(NOOFNODE_LEN)
2286 #if (MAXLEVEL > NOOFNODEMAX)
2287 #error **** set NOOFNODEMAX/_LEN appropriate to MAXLEVEL: 2^NOOFNODE_LEN = NOOFNODEMAX >= MAXLEVEL ****
2288 #endif
2289 #define NOOFNODE(p) CW_READ_STATIC(p,NOOFNODE_,VERTEX_)
2290 #define SETNOOFNODE(p,n) CW_WRITE_STATIC(p,NOOFNODE_,VERTEX_,n)
2291 #define INCNOOFNODE(p) SETNOOFNODE(p,NOOFNODE(p)+1)
2292 #define DECNOOFNODE(p) SETNOOFNODE(p,NOOFNODE(p)-1)
2293
2294 #define PREDV(p) ((p)->iv.pred)
2295 #define SUCCV(p) ((p)->iv.succ)
2296 #define CVECT(p) ((p)->iv.x)
2297 #define XC(p) ((p)->iv.x[0])
2298 #define YC(p) ((p)->iv.x[1])
2299 #define ZC(p) ((p)->iv.x[2])
2300 #define LCVECT(p) ((p)->iv.xi)
2301 #define XI(p) ((p)->iv.xi[0])
2302 #define ETA(p) ((p)->iv.xi[1])
2303 #define NU(p) ((p)->iv.xi[2])
2304 #define VDATA(p) ((p)->iv.data)
2305 #define VFATHER(p) ((p)->iv.father)
2306
2307 /* for boundary vertices */
2308 #define V_BNDP(p) ((p)->bv.bndp)
2309
2310 /* parallel macros */
2311 #ifdef ModelP
2312 #define PARHDRV(p) (&((p)->iv.ddd))
2313 #endif /* ModelP */
2314
2315 /****************************************************************************/
2316 /* */
2317 /* macros for nodes */
2318 /* */
2319 /****************************************************************************/
2320
2321 /* control word offset */
2322 #define NODE_OFFSET 0
2323
2324 #define NTYPE_SHIFT 0
2325 #define NTYPE_LEN 3
2326 #define NTYPE(p) CW_READ_STATIC(p,NTYPE_,NODE_)
2327 #define SETNTYPE(p,n) CW_WRITE_STATIC(p,NTYPE_,NODE_,n)
2328
2329 #define NSUBDOM_SHIFT 3
2330 #define NSUBDOM_LEN 6
2331 #define NSUBDOM(p) CW_READ_STATIC(p,NSUBDOM_,NODE_)
2332 #define SETNSUBDOM(p,n) CW_WRITE_STATIC(p,NSUBDOM_,NODE_,n)
2333
2334 #define NPROP_SHIFT 11
2335 #define NPROP_LEN 4
2336 #define NPROP(p) CW_READ_STATIC(p,NPROP_,NODE_)
2337 #define SETNPROP(p,n) CW_WRITE_STATIC(p,NPROP_,NODE_,n)
2338
2339 #define MODIFIED_SHIFT 15
2340 #define MODIFIED_LEN 1
2341 #define MODIFIED(p) CW_READ_STATIC(p,MODIFIED_,NODE_)
2342 #define SETMODIFIED(p,n) CW_WRITE_STATIC(p,MODIFIED_,NODE_,n)
2343
2344 #define NCLASS_SHIFT 16
2345 #define NCLASS_LEN 2
2346 #define NCLASS(p) CW_READ_STATIC(p,NCLASS_,NODE_)
2347 #define SETNCLASS(p,n) CW_WRITE_STATIC(p,NCLASS_,NODE_,n)
2348
2349 #define NNCLASS_SHIFT 18
2350 #define NNCLASS_LEN 2
2351 #define NNCLASS(p) CW_READ_STATIC(p,NNCLASS_,NODE_)
2352 #define SETNNCLASS(p,n) CW_WRITE_STATIC(p,NNCLASS_,NODE_,n)
2353
2354 #if defined ModelP && defined __OVERLAP2__
2355 #define NO_DELETE_OVERLAP2_LEN 1
2356 #define NO_DELETE_OVERLAP2(p) CW_READ(p,ce_NO_DELETE_OVERLAP2)
2357 #define SETNO_DELETE_OVERLAP2(p,n) CW_WRITE(p,ce_NO_DELETE_OVERLAP2,n)
2358 #endif
2359
2360 #define PREDN(p) ((p)->pred)
2361 #define SUCCN(p) ((p)->succ)
2362 #define START(p) ((p)->start)
2363
2364 #define NFATHER(p) ((NODE*)(p)->father)
2365 #define SETNFATHER(p,n) ((p)->father = n)
2366 #define NFATHEREDGE(p) ((EDGE*)(p)->father)
2367 /*
2368 #define NFATHER(p) ((NTYPE(p) == CORNER_NODE) ? (p)->father : NULL)
2369 #define NFATHEREDGE(p) ((NTYPE(p) == MID_NODE) ? (EDGE *)(p)->father : NULL)
2370 #define SETNFATHEREDGE(p,e) ((p)->father = (NODE *) (e))
2371 */
2372 #define CORNERTYPE(p) (NTYPE(p) == CORNER_NODE)
2373 #define MIDTYPE(p) (NTYPE(p) == MID_NODE)
2374 #define SIDETYPE(p) (NTYPE(p) == SIDE_NODE)
2375 #define CENTERTYPE(p) (NTYPE(p) == CENTER_NODE)
2376
2377 #define SONNODE(p) ((p)->son)
2378 #define MYVERTEX(p) ((p)->myvertex)
2379 #define NDATA(p) ((p)->data)
2380 #define NVECTOR(p) ((p)->vector)
2381
2382 #define NODE_ELEMENT_LIST(p) ((ELEMENTLIST *)(p)->data)
2383 #define ELEMENT_PTR(p) ((p)->el)
2384
2385 /****************************************************************************/
2386 /* */
2387 /* macros for links */
2388 /* */
2389 /****************************************************************************/
2390
2391 /* CAUTION: the controlword of LINK0 and its edge are identical (AVOID overlapping of flags) */
2392
2393 /* control word offset */
2394 #define LINK_OFFSET 0
2395
2396 #define LOFFSET_SHIFT 0
2397 #define LOFFSET_LEN 1
2398 #define LOFFSET(p) CW_READ(p,LOFFSET_CE)
2399 #define SETLOFFSET(p,n) CW_WRITE(p,LOFFSET_CE,n)
2400
2401 #define NBNODE(p) ((p)->nbnode)
2402 #define NEXT(p) ((p)->next)
2403 #define LDATA(p) ((p)->matelem)
2404 #define MATELEM(p) ((p)->matelem) /* can be used for node and link */
2405
2406 #define MYEDGE(p) ((EDGE *)((p)-LOFFSET(p)))
2407 #define REVERSE(p) ((p)+(1-LOFFSET(p)*2))
2408
2409 #if defined(__TWODIM__)
2410 #define LELEM(p) ((p)->elem)
2411 #define SET_LELEM(p,e) ((p)->elem = (e))
2412 #endif
2413
2414 /****************************************************************************/
2415 /* */
2416 /* macros for edges */
2417 /* */
2418 /****************************************************************************/
2419
2420 /* control word offset */
2421 #define EDGE_OFFSET 0
2422
2423 #define NO_OF_ELEM_SHIFT 2
2424 #define NO_OF_ELEM_LEN 7
2425 #define NO_OF_ELEM_MAX 128
2426 #define NO_OF_ELEM(p) CW_READ(p,NO_OF_ELEM_CE)
2427 #define SET_NO_OF_ELEM(p,n) CW_WRITE(p,NO_OF_ELEM_CE,n)
2428 #define INC_NO_OF_ELEM(p) SET_NO_OF_ELEM(p,NO_OF_ELEM(p)+1)
2429 #define DEC_NO_OF_ELEM(p) SET_NO_OF_ELEM(p,NO_OF_ELEM(p)-1)
2430
2431 #define AUXEDGE_SHIFT 9
2432 #define AUXEDGE_LEN 1
2433 #define AUXEDGE(p) CW_READ(p,AUXEDGE_CE)
2434 #define SETAUXEDGE(p,n) CW_WRITE(p,AUXEDGE_CE,n)
2435
2436 #define EDGENEW_SHIFT 1
2437 #define EDGENEW_LEN 1
2438 #define EDGENEW(p) CW_READ(p,EDGENEW_CE)
2439 #define SETEDGENEW(p,n) CW_WRITE(p,EDGENEW_CE,n)
2440
2441 /* boundary edges will be indicated by a subdomain id of 0 */
2442 #define EDSUBDOM_SHIFT 12
2443 #define EDSUBDOM_LEN 6
2444 #define EDSUBDOM(p) CW_READ(p,EDSUBDOM_CE)
2445 #define SETEDSUBDOM(p,n) CW_WRITE(p,EDSUBDOM_CE,n)
2446
2447 #define LINK0(p) (&((p)->links[0]))
2448 #define LINK1(p) (&((p)->links[1]))
2449 #define MIDNODE(p) ((p)->midnode)
2450 #define EDDATA(p) ((p)->data)
2451 #define EDVECTOR(p) ((p)->vector)
2452
2453 /****************************************************************************/
2454 /* */
2455 /* macros for elements */
2456 /* */
2457 /****************************************************************************/
2458
2459 enum {TRIANGLE = 3,
2460 QUADRILATERAL = 4};
2461
2462 enum {TETRAHEDRON = 4,
2463 PYRAMID = 5,
2464 PRISM = 6,
2465 HEXAHEDRON = 7};
2466
2467 /* control word offsets */
2468 #define ELEMENT_OFFSET 0
2469 #define FLAG_OFFSET 2
2470 #define PROPERTY_OFFSET 3
2471
2472 /* macros for control word */
2473 #define ECLASS_SHIFT 8
2474 #define ECLASS_LEN 2
2475 #define ECLASS(p) CW_READ(p,ECLASS_CE)
2476 #define SETECLASS(p,n) CW_WRITE(p,ECLASS_CE,n)
2477
2478 #define NSONS_SHIFT 10
2479 #define NSONS_LEN 5
2480 #define NSONS(p) CW_READ(p,NSONS_CE)
2481 #define SETNSONS(p,n) CW_WRITE(p,NSONS_CE,n)
2482
2483 #define NEWEL_SHIFT 17
2484 #define NEWEL_LEN 1
2485 #define NEWEL(p) CW_READ(p,NEWEL_CE)
2486 #define SETNEWEL(p,n) CW_WRITE(p,NEWEL_CE,n)
2487
2488 /* macros for flag word */
2489 /* are obviously all for internal use */
2490
2491 /* the property field */
2492 #define SUBDOMAIN_SHIFT 24
2493 #define SUBDOMAIN_LEN 6
2494 #define SUBDOMAIN(p) CW_READ(p,SUBDOMAIN_CE)
2495 #define SETSUBDOMAIN(p,n) CW_WRITE(p,SUBDOMAIN_CE,n)
2496
2497 #define NODEORD_SHIFT 0
2498 #define NODEORD_LEN 24
2499 #define NODEORD(p) CW_READ(p,NODEORD_CE)
2500 #define SETNODEORD(p,n) CW_WRITE(p,NODEORD_CE,n)
2501
2502 #define PROP_SHIFT 30
2503 #define PROP_LEN 2
2504 #define PROP(p) CW_READ(p,PROP_CE)
2505 #define SETPROP(p,n) CW_WRITE(p,PROP_CE,n)
2506
2507 /* parallel macros */
2508 #ifdef ModelP
2509 #define PARTITION(p) ((p)->ge.lb1)
2510 #define PARHDRE(p) (&((p)->ge.ddd))
2511 #endif
2512
2513 /*******************************/
2514 /* the general element concept */
2515 /*******************************/
2516
2517 /** \brief This structure contains all topological properties
2518 of an element and more ..
2519 */
2520 struct GENERAL_ELEMENT {
2521 INT tag; /**< Element type to be defined */
2522
2523 /* the following parameters determine size of refs array in element */
2524 INT max_sons_of_elem; /**< Max number of sons for this type */
2525 INT sides_of_elem; /**< How many sides ? */
2526 INT corners_of_elem; /**< How many corners ? */
2527
2528 /* local geometric description of the element */
2529 DOUBLE_VECTOR local_corner[MAX_CORNERS_OF_ELEM]; /**< Local coordinates of the corners of the element */
2530
2531 /* more size parameters */
2532 INT edges_of_elem; /**< How many edges ? */
2533 INT edges_of_side[MAX_SIDES_OF_ELEM]; /**< Number of edges for each side */
2534 INT corners_of_side[MAX_SIDES_OF_ELEM]; /**< Number of corners for each side */
2535 INT corners_of_edge; /**< Is always 2 ! */
2536
2537 /* index computations */
2538 /* Within each element sides, edges, corners are numbered in some way. */
2539 /* Within each side the edges and corners are numbered, within the edge the */
2540 /* corners are numbered. The following arrays map the local numbers within */
2541 /* the side or edge to the numbering within the element. */
2542 INT edge_of_side[MAX_SIDES_OF_ELEM][MAX_EDGES_OF_SIDE];
2543 INT corner_of_side[MAX_SIDES_OF_ELEM][MAX_CORNERS_OF_SIDE];
2544 INT corner_of_edge[MAX_EDGES_OF_ELEM][MAX_CORNERS_OF_EDGE];
2545
2546 /* the following parameters are derived from data above */
2547 GM_OBJECTS mapped_inner_objt = NOOBJ; /* tag to objt mapping for free list*/
2548 GM_OBJECTS mapped_bnd_objt = NOOBJ; /* tag to objt mapping for free list*/
2549 INT inner_size, bnd_size; /* size in bytes used for alloc */
2550 INT edge_with_corners[MAX_CORNERS_OF_ELEM][MAX_CORNERS_OF_ELEM];
2551 INT side_with_edge[MAX_EDGES_OF_ELEM][MAX_SIDES_OF_EDGE];
2552 INT corner_of_side_inv[MAX_SIDES_OF_ELEM][MAX_CORNERS_OF_ELEM];
2553 INT edges_of_corner[MAX_CORNERS_OF_ELEM][MAX_EDGES_OF_ELEM];
2554 INT corner_of_oppedge[MAX_EDGES_OF_ELEM][MAX_CORNERS_OF_EDGE];
2555 INT corner_opp_to_side[MAX_SIDES_OF_ELEM];
2556 INT opposite_edge[MAX_EDGES_OF_ELEM];
2557 INT side_opp_to_corner[MAX_CORNERS_OF_ELEM];
2558 INT edge_of_corner[MAX_CORNERS_OF_ELEM][MAX_EDGES_OF_ELEM];
2559 INT edge_of_two_sides[MAX_SIDES_OF_ELEM][MAX_SIDES_OF_ELEM];
2560
2561 /* ... the refinement rules should be placed here later */
2562 };
2563
2564 END_UGDIM_NAMESPACE
2565
2566 /** \todo move this to include section, when other general element stuff is separated */
2567 #include "elements.h"
2568
2569 START_UGDIM_NAMESPACE
2570
2571 /****************************************************************************/
2572 /* */
2573 /* macros for element descriptors */
2574 /* */
2575 /****************************************************************************/
2576
2577 /** @name Macros to access element descriptors by element pointers */
2578 /*@{*/
2579 #define SIDES_OF_ELEM(p) (element_descriptors[TAG(p)]->sides_of_elem)
2580 #define EDGES_OF_ELEM(p) (element_descriptors[TAG(p)]->edges_of_elem)
2581 #define CORNERS_OF_ELEM(p) (element_descriptors[TAG(p)]->corners_of_elem)
2582 #define LOCAL_COORD_OF_ELEM(p,c) (element_descriptors[TAG(p)]->local_corner[(c)])
2583
2584
2585 #define SONS_OF_ELEM(p) (element_descriptors[TAG(p)]->max_sons_of_elem) /* this is the number of pointers ! */
2586
2587 #define EDGES_OF_SIDE(p,i) (element_descriptors[TAG(p)]->edges_of_side[(i)])
2588 #define CORNERS_OF_SIDE(p,i) (element_descriptors[TAG(p)]->corners_of_side[(i)])
2589
2590 #define CORNERS_OF_EDGE 2
2591
2592 #define EDGE_OF_SIDE(p,s,e) (element_descriptors[TAG(p)]->edge_of_side[(s)][(e)])
2593 #define EDGE_OF_TWO_SIDES(p,s,t) (element_descriptors[TAG(p)]->edge_of_two_sides[(s)][(t)])
2594 #define CORNER_OF_SIDE(p,s,c) (element_descriptors[TAG(p)]->corner_of_side[(s)][(c)])
2595 #define CORNER_OF_EDGE(p,e,c) (element_descriptors[TAG(p)]->corner_of_edge[(e)][(c)])
2596
2597 #define EDGE_WITH_CORNERS(p,c0,c1) (element_descriptors[TAG(p)]->edge_with_corners[(c0)][(c1)])
2598 #define SIDE_WITH_EDGE(p,e,k) (element_descriptors[TAG(p)]->side_with_edge[(e)][(k)])
2599 #define CORNER_OF_SIDE_INV(p,s,c) (element_descriptors[TAG(p)]->corner_of_side_inv[(s)][(c)])
2600 #define EDGES_OF_CORNER(p,c,k) (element_descriptors[TAG(p)]->edges_of_corner[(c)][(k)])
2601 #define CORNER_OF_OPPEDGE(p,e,c) (element_descriptors[TAG(p)]->corner_of_oppedge[(e)][(c)])
2602 #define CORNER_OPP_TO_SIDE(p,s) (element_descriptors[TAG(p)]->corner_opp_to_side[(s)])
2603 #define OPPOSITE_EDGE(p,e) (element_descriptors[TAG(p)]->opposite_edge[(e)])
2604 #define SIDE_OPP_TO_CORNER(p,c) (element_descriptors[TAG(p)]->side_opp_to_corner[(c)])
2605 #define EDGE_OF_CORNER(p,c,e) (element_descriptors[TAG(p)]->edge_of_corner[(c)][(e)])
2606
2607 #define CTRL2(p) ((p)->ge.flag)
2608 #define FLAG(p) ((p)->ge.flag)
2609 #define SUCCE(p) ((p)->ge.succ)
2610 #define PREDE(p) ((p)->ge.pred)
2611
2612 #ifdef __CENTERNODE__
2613 #define CENTERNODE(p) ((p)->ge.centernode)
2614 #endif
2615 #define CORNER(p,i) ((NODE *) (p)->ge.refs[n_offset[TAG(p)]+(i)])
2616 #define EFATHER(p) ((ELEMENT *) (p)->ge.refs[father_offset[TAG(p)]])
2617 #define SON(p,i) ((ELEMENT *) (p)->ge.refs[sons_offset[TAG(p)]+(i)])
2618 /** \todo NbElem is declared in ugm.h, but never defined.
2619 We need a clean solution. */
2620 /*
2621 #if defined(__TWODIM__)
2622 #define NBELEM(p,i) NbElem((p),(i))
2623 #else
2624 */
2625 #define NBELEM(p,i) ((ELEMENT *) (p)->ge.refs[nb_offset[TAG(p)]+(i)])
2626 /*
2627 #endif
2628 */
2629 #define ELEM_BNDS(p,i) ((BNDS *) (p)->ge.refs[side_offset[TAG(p)]+(i)])
2630 #define EVECTOR(p) ((VECTOR *) (p)->ge.refs[evector_offset[TAG(p)]])
2631 #define SVECTOR(p,i) ((VECTOR *) (p)->ge.refs[svector_offset[TAG(p)]+(i)])
2632 #define SIDE_ON_BND(p,i) (ELEM_BNDS(p,i) != NULL)
2633 #define INNER_SIDE(p,i) (ELEM_BNDS(p,i) == NULL)
2634 #define INNER_BOUNDARY(p,i) (InnerBoundary(p,i))
2635 /* TODO: replace by function call */
2636
2637 #ifdef __TWODIM__
2638 #define EDGE_ON_BND(p,i) (ELEM_BNDS(p,i) != NULL)
2639 #endif
2640
2641 #ifdef __THREEDIM__
2642 #define EDGE_ON_BND(p,i) (SIDE_ON_BND(p,SIDE_WITH_EDGE(p,i,0)) || \
2643 SIDE_ON_BND(p,SIDE_WITH_EDGE(p,i,1)))
2644 #endif
2645 /*@}*/
2646
2647 /* use the following macros to assign values, since definition */
2648 /* above is no proper lvalue. */
2649 #ifdef __CENTERNODE__
2650 #define SET_CENTERNODE(p,q) ((p)->ge.centernode = q)
2651 #endif
2652 #define SET_CORNER(p,i,q) ((p)->ge.refs[n_offset[TAG(p)]+(i)] = q)
2653 #define SET_EFATHER(p,q) ((p)->ge.refs[father_offset[TAG(p)]] = q)
2654 #define SET_SON(p,i,q) ((p)->ge.refs[sons_offset[TAG(p)]+(i)] = q)
2655 /** \todo Set_NbElem is declared in ugm.h, but never defined.
2656 We need a clean solution. */
2657 /*
2658 #if defined(__TWODIM__)
2659 #define SET_NBELEM(p,i,q) Set_NbElem((p),(i),(q))
2660 #else
2661 */
2662 #define SET_NBELEM(p,i,q) ((p)->ge.refs[nb_offset[TAG(p)]+(i)] = q)
2663 /*
2664 #endif
2665 */
2666 #if defined(__TWODIM__)
2667 #define VOID_NBELEM(p,i) NBELEM(p,i)
2668 #else
2669 #define VOID_NBELEM(p,i) ((p)->ge.refs[nb_offset[TAG(p)]+(i)])
2670 #endif
2671 #define SET_BNDS(p,i,q) ((p)->ge.refs[side_offset[TAG(p)]+(i)] = q)
2672 #define SET_EVECTOR(p,q) ((p)->ge.refs[evector_offset[TAG(p)]] = q)
2673 #define SET_SVECTOR(p,i,q) ((p)->ge.refs[svector_offset[TAG(p)]+(i)] = q)
2674
2675 /** @name Macros to access corner pointers directly */
2676 /*@{*/
2677 #define CORNER_OF_EDGE_PTR(e,i,j) (CORNER(e,CORNER_OF_EDGE(e,i,j)))
2678 #define CORNER_OF_SIDE_PTR(e,i,j) (CORNER(e,CORNER_OF_SIDE(e,i,j)))
2679 /*@}*/
2680
2681 /** @name Macros to access element descriptors by element tags */
2682 /*@{*/
2683 #define INNER_SIZE_TAG(t) (element_descriptors[t]->inner_size)
2684 #define BND_SIZE_TAG(t) (element_descriptors[t]->bnd_size)
2685 #define MAPPED_INNER_OBJT_TAG(t) (element_descriptors[t]->mapped_inner_objt)
2686 #define MAPPED_BND_OBJT_TAG(t) (element_descriptors[t]->mapped_bnd_objt)
2687
2688 #define SIDES_OF_TAG(t) (element_descriptors[t]->sides_of_elem)
2689 #define EDGES_OF_TAG(t) (element_descriptors[t]->edges_of_elem)
2690 #define CORNERS_OF_TAG(t) (element_descriptors[t]->corners_of_elem)
2691 #define LOCAL_COORD_OF_TAG(t,c) (element_descriptors[t]->local_corner[(c)])
2692
2693 #define SONS_OF_TAG(t) (element_descriptors[t]->max_sons_of_elem) /* this is the number of pointers ! */
2694
2695 #define EDGES_OF_SIDE_TAG(t,i) (element_descriptors[t]->edges_of_side[(i)])
2696 #define CORNERS_OF_SIDE_TAG(t,i) (element_descriptors[t]->corners_of_side[(i)])
2697
2698 #define EDGE_OF_SIDE_TAG(t,s,e) (element_descriptors[t]->edge_of_side[(s)][(e)])
2699 #define EDGE_OF_TWO_SIDES_TAG(t,s,u) (element_descriptors[t]->edge_of_two_sides[(s)][(u)])
2700 #define CORNER_OF_SIDE_TAG(t,s,c) (element_descriptors[t]->corner_of_side[(s)][(c)])
2701 #define CORNER_OF_EDGE_TAG(t,e,c) (element_descriptors[t]->corner_of_edge[(e)][(c)])
2702
2703 #define EDGE_WITH_CORNERS_TAG(t,c0,c1) (element_descriptors[t]->edge_with_corners[(c0)][(c1)])
2704 #define SIDE_WITH_EDGE_TAG(t,e,k) (element_descriptors[t]->side_with_edge[(e)][(k)])
2705 #define CORNER_OF_SIDE_INV_TAG(t,s,c) (element_descriptors[t]->corner_of_side_inv[(s)][(c)])
2706 #define EDGES_OF_CORNER_TAG(t,c,k) (element_descriptors[t]->edges_of_corner[(c)][(k)])
2707 #define CORNER_OF_OPPEDGE_TAG(t,e,c) (element_descriptors[t]->corner_of_oppedge[(e)][(c)])
2708 #define CORNER_OPP_TO_SIDE_TAG(t,s) (element_descriptors[t]->corner_opp_to_side[(s)])
2709 #define OPPOSITE_EDGE_TAG(t,e) (element_descriptors[t]->opposite_edge[(e)])
2710 #define SIDE_OPP_TO_CORNER_TAG(t,c) (element_descriptors[t]->side_opp_to_corner[(c)])
2711 #define EDGE_OF_CORNER_TAG(t,c,e) (element_descriptors[t]->edge_of_corner[(c)][(e)])
2712 /*@}*/
2713
2714 /** @name Macros to access reference descriptors by number of element corners */
2715 /*@{*/
2716 #define SIDES_OF_REF(n) (reference_descriptors[n]->sides_of_elem)
2717 #define EDGES_OF_REF(n) (reference_descriptors[n]->edges_of_elem)
2718 #define CORNERS_OF_REF(n) (reference_descriptors[n]->corners_of_elem)
2719 #define LOCAL_COORD_OF_REF(n,c) (reference_descriptors[n]->local_corner[(c)])
2720 #define EDGES_OF_SIDE_REF(n,i) (reference_descriptors[n]->edges_of_side[(i)])
2721 #define CORNERS_OF_SIDE_REF(n,i) (reference_descriptors[n]->corners_of_side[(i)])
2722 #define EDGE_OF_SIDE_REF(n,s,e) (reference_descriptors[n]->edge_of_side[(s)][(e)])
2723 #define EDGE_OF_TWO_SIDES_REF(n,s,t) (reference_descriptors[n]->edge_of_two_sides[(s)][(t)])
2724 #define CORNER_OF_SIDE_REF(n,s,c) (reference_descriptors[n]->corner_of_side[(s)][(c)])
2725 #define CORNER_OF_EDGE_REF(n,e,c) (reference_descriptors[n]->corner_of_edge[(e)][(c)])
2726 #define EDGE_WITH_CORNERS_REF(n,c0,c1) (reference_descriptors[n]->edge_with_corners[(c0)][(c1)])
2727 #define SIDE_WITH_EDGE_REF(n,e,k) (reference_descriptors[n]->side_with_edge[(e)][(k)])
2728 #define CORNER_OF_SIDE_INV_REF(n,s,c) (reference_descriptors[n]->corner_of_side_inv[(s)][(c)])
2729 #define EDGES_OF_CORNER_REF(n,c,k) (reference_descriptors[n]->edges_of_corner[(c)][(k)])
2730 #define CORNER_OF_OPPEDGE_REF(n,e,c) (reference_descriptors[n]->corner_of_oppedge[(e)][(c)])
2731 #define CORNER_OPP_TO_SIDE_REF(n,s) (reference_descriptors[n]->corner_opp_to_side[(s)])
2732 #define OPPOSITE_EDGE_REF(n,e) (reference_descriptors[n]->opposite_edge[(e)])
2733 #define SIDE_OPP_TO_CORNER_REF(n,c) (reference_descriptors[n]->side_opp_to_corner[(c)])
2734 #define EDGE_OF_CORNER_REF(n,c,e) (reference_descriptors[n]->edge_of_corner[(c)][(e)])
2735 /*@}*/
2736
2737 /****************************************************************************/
2738 /* */
2739 /* macros for grids */
2740 /* */
2741 /****************************************************************************/
2742
2743 /* control word offset */
2744 #define GRID_OFFSET 0
2745 #define GRID_STATUS_OFFSET 1
2746
2747 #define GLEVEL(p) ((p)->level)
2748 #define GATTR(p) ((p)->attribut)
2749 #define GFORMAT(p) MGFORMAT(MYMG(p))
2750 #define SETGLOBALGSTATUS(p) ((p)->status=~0)
2751 #define GSTATUS(p,n) ((p)->status&(n))
2752 #define RESETGSTATUS(p,n) ((p)->status&=~(n))
2753
2754 #ifdef ModelP
2755 #define PFIRSTELEMENT(p) ((LISTPART_FIRSTELEMENT(p,0)!=NULL) ?\
2756 (LISTPART_FIRSTELEMENT(p,0)) : (FIRSTELEMENT(p)))
2757 #define PRIO_FIRSTELEMENT(p,prio) ((p)->elements[PRIO2LISTPART(ELEMENT_LIST,prio)])
2758 #define LISTPART_FIRSTELEMENT(p,part) ((p)->elements[part])
2759 #define FIRSTELEMENT(p) ((p)->elements[PRIO2LISTPART(ELEMENT_LIST,PrioMaster)])
2760
2761 #define PLASTELEMENT(p) LASTELEMENT(p)
2762 #define PRIO_LASTELEMENT(p,prio) ((p)->lastelement[PRIO2LISTPART(ELEMENT_LIST,prio)])
2763 #define LISTPART_LASTELEMENT(p,part) ((p)->lastelement[part])
2764 #define LASTELEMENT(p) ((p)->lastelement[PRIO2LISTPART(ELEMENT_LIST,PrioMaster)])
2765 #else
2766 #define FIRSTELEMENT(p) ((p)->elements[0])
2767 #define PFIRSTELEMENT(p) FIRSTELEMENT(p)
2768 #define LASTELEMENT(p) ((p)->lastelement[0])
2769 #define PLASTELEMENT(p) LASTELEMENT(p)
2770 #endif
2771
2772 #ifdef ModelP
2773 #define PFIRSTVERTEX(p) ((LISTPART_FIRSTVERTEX(p,0)!=NULL) ?\
2774 (LISTPART_FIRSTVERTEX(p,0)) :\
2775 ((LISTPART_FIRSTVERTEX(p,1)!=NULL) ?\
2776 (LISTPART_FIRSTVERTEX(p,1)) : (FIRSTVERTEX(p))))
2777 #define PRIO_FIRSTVERTEX(p,prio) ((p)->vertices[PRIO2LISTPART(VERTEX_LIST,prio)])
2778 #define LISTPART_FIRSTVERTEX(p,part) ((p)->vertices[part])
2779 #define FIRSTVERTEX(p) (((p)->vertices[PRIO2LISTPART(VERTEX_LIST,\
2780 PrioBorder)]!=NULL) ?\
2781 (p)->vertices[PRIO2LISTPART(VERTEX_LIST,PrioBorder)] :\
2782 (p)->vertices[PRIO2LISTPART(VERTEX_LIST,PrioMaster)])
2783 #define SFIRSTVERTEX(p) (p)->vertices[PRIO2LISTPART(VERTEX_LIST,PrioMaster)]
2784
2785 #define PLASTVERTEX(p) LASTVERTEX(p)
2786 #define PRIO_LASTVERTEX(p,prio) ((p)->lastvertex[PRIO2LISTPART(VERTEX_LIST,prio)])
2787 #define LISTPART_LASTVERTEX(p,part) ((p)->lastvertex[part])
2788 #define LASTVERTEX(p) ((p)->lastvertex[PRIO2LISTPART(VERTEX_LIST,PrioMaster)])
2789 #else
2790 #define FIRSTVERTEX(p) ((p)->vertices[0])
2791 #define PFIRSTVERTEX(p) FIRSTVERTEX(p)
2792 #define SFIRSTVERTEX(p) FIRSTVERTEX(p)
2793 #define LASTVERTEX(p) ((p)->lastvertex[0])
2794 #define PLASTVERTEX(p) LASTVERTEX(p)
2795 #endif
2796
2797 #define FIRSTELEMSIDE(p) ((p)->sides)
2798
2799 #ifdef ModelP
2800 #define PFIRSTNODE(p) ((LISTPART_FIRSTNODE(p,0)!=NULL) ?\
2801 (LISTPART_FIRSTNODE(p,0)) :\
2802 ((LISTPART_FIRSTNODE(p,1)!=NULL) ?\
2803 (LISTPART_FIRSTNODE(p,1)) : (FIRSTNODE(p))))
2804 #define PRIO_FIRSTNODE(p,prio) ((p)->firstNode[PRIO2LISTPART(NODE_LIST,prio)])
2805 #define LISTPART_FIRSTNODE(p,part) ((p)->firstNode[part])
2806 #define FIRSTNODE(p) (((p)->firstNode[PRIO2LISTPART(NODE_LIST,\
2807 PrioBorder)]!=NULL) ?\
2808 (p)->firstNode[PRIO2LISTPART(NODE_LIST,PrioBorder)] :\
2809 (p)->firstNode[PRIO2LISTPART(NODE_LIST,PrioMaster)])
2810 #define SFIRSTNODE(p) (p)->firstNode[PRIO2LISTPART(NODE_LIST,PrioMaster)]
2811
2812 #define PLASTNODE(p) LASTNODE(p)
2813 #define PRIO_LASTNODE(p,prio) ((p)->lastNode[PRIO2LISTPART(NODE_LIST,prio)])
2814 #define LISTPART_LASTNODE(p,part) ((p)->lastNode[part])
2815 #define LASTNODE(p) ((p)->lastNode[PRIO2LISTPART(NODE_LIST,PrioMaster)])
2816 #else
2817 #define FIRSTNODE(p) ((p)->firstNode[0])
2818 #define PFIRSTNODE(p) FIRSTNODE(p)
2819 #define SFIRSTNODE(p) FIRSTNODE(p)
2820 #define LASTNODE(p) ((p)->lastNode[0])
2821 #define PLASTNODE(p) LASTNODE(p)
2822 #endif
2823
2824 #ifdef ModelP
2825 #define PFIRSTVECTOR(p) ((LISTPART_FIRSTVECTOR(p,0)!=NULL) ?\
2826 (LISTPART_FIRSTVECTOR(p,0)) :\
2827 ((LISTPART_FIRSTVECTOR(p,1)!=NULL) ?\
2828 (LISTPART_FIRSTVECTOR(p,1)) : (FIRSTVECTOR(p))))
2829 #define PRIO_FIRSTVECTOR(p,prio) ((p)->firstVector[PRIO2LISTPART(VECTOR_LIST,prio)])
2830 #define LISTPART_FIRSTVECTOR(p,part) ((p)->firstVector[part])
2831 #define FIRSTVECTOR(p) (((p)->firstVector[PRIO2LISTPART(VECTOR_LIST,\
2832 PrioBorder)]!=NULL) ?\
2833 (p)->firstVector[PRIO2LISTPART(VECTOR_LIST,PrioBorder)] :\
2834 (p)->firstVector[PRIO2LISTPART(VECTOR_LIST,PrioMaster)])
2835 #define SFIRSTVECTOR(p) (p)->firstVector[PRIO2LISTPART(VECTOR_LIST,PrioMaster)]
2836
2837 #define PLASTVECTOR(p) LASTVECTOR(p)
2838 #define PRIO_LASTVECTOR(p,prio) ((p)->lastVector[PRIO2LISTPART(VECTOR_LIST,prio)])
2839 #define LISTPART_LASTVECTOR(p,part) ((p)->lastVector[part])
2840 #define LASTVECTOR(p) ((p)->lastVector[PRIO2LISTPART(VECTOR_LIST,PrioMaster)])
2841 #else
2842 #define FIRSTVECTOR(p) ((p)->firstVector[0])
2843 #define PFIRSTVECTOR(p) FIRSTVECTOR(p)
2844 #define SFIRSTVECTOR(p) FIRSTVECTOR(p)
2845 #define LASTVECTOR(p) ((p)->lastVector[0])
2846 #define PLASTVECTOR(p) LASTVECTOR(p)
2847 #endif
2848
2849 #define UPGRID(p) ((p)->finer)
2850 #define DOWNGRID(p) ((p)->coarser)
2851 #define MYMG(p) ((p)->mg)
2852 #define NV(p) ((p)->nVert[0])
2853 #define NN(p) ((p)->nNode[0])
2854 #define NT(p) ((p)->nElem[0])
2855 #define NVEC(p) ((p)->nVector[0])
2856 #ifdef ModelP
2857 #define NV_PRIO(p,prio) ((p)->nVert[prio])
2858 #define NN_PRIO(p,prio) ((p)->nNode[prio])
2859 #define NT_PRIO(p,prio) ((p)->nElem[prio])
2860 #define NVEC_PRIO(p,prio) ((p)->nVector[prio])
2861 #endif
2862 #define NE(p) ((p)->nEdge)
2863 #define NS(p) ((p)->nSide)
2864 #define NC(p) ((p)->nCon)
2865 #define VEC_DEF_IN_OBJ_OF_GRID(p,tp) (GFORMAT(p)->OTypeUsed[(tp)]>0)
2866 #define NIMAT(p) ((p)->nIMat)
2867
2868 #define GRID_ATTR(g) ((unsigned char) (GLEVEL(g)+32))
2869 #define ATTR_TO_GLEVEL(i) (i-32)
2870
2871 inline const PPIF::PPIFContext&
ppifContext()2872 grid::ppifContext() const
2873 {
2874 return mg->ppifContext();
2875 }
2876
2877 #ifdef ModelP
2878 inline const DDD::DDDContext&
dddContext()2879 grid::dddContext() const
2880 {
2881 return mg->dddContext();
2882 }
2883
2884 inline DDD::DDDContext&
dddContext()2885 grid::dddContext()
2886 {
2887 return mg->dddContext();
2888 }
2889 #endif
2890
2891 /****************************************************************************/
2892 /* */
2893 /* macros for multigrids */
2894 /* */
2895 /****************************************************************************/
2896
2897 /* control word offset */
2898 #define MULTIGRID_STATUS_OFFSET ((sizeof(ENVDIR))/sizeof(UINT))
2899
2900 #define MGSTATUS(p) ((p)->status)
2901 #define RESETMGSTATUS(p) {(p)->status=0; (p)->magic_cookie = (int)time(NULL); (p)->saved=0;}
2902 #define MG_MAGIC_COOKIE(p) ((p)->magic_cookie)
2903 #define VIDCNT(p) ((p)->vertIdCounter)
2904 #define NIDCNT(p) ((p)->nodeIdCounter)
2905 #define EIDCNT(p) ((p)->elemIdCounter)
2906 #define TOPLEVEL(p) ((p)->topLevel)
2907 #define CURRENTLEVEL(p) ((p)->currentLevel)
2908 #define FULLREFINELEVEL(p) ((p)->fullrefineLevel)
2909 #define MGFORMAT(p) ((p)->theFormat)
2910 #define MG_BVP(p) ((p)->theBVP)
2911 #define MG_BVPD(p) (&((p)->theBVPD))
2912 #define MGBNDSEGDESC(p,i) (&((p)->segments[i]))
2913 #define MGVERTEX(p,k) ((p)->corners[k])
2914 #define MGNOOFCORNERS(p) ((p)->numOfCorners)
2915 #define MGHEAP(p) ((p)->theHeap)
2916 #define MG_NPROPERTY(p) ((p)->nProperty)
2917 #define GRID_ON_LEVEL(p,i) ((p)->grids[i])
2918 #define MGNAME(p) ((p)->v.name)
2919 #define VEC_DEF_IN_OBJ_OF_MG(p,tp) (MGFORMAT(p)->OTypeUsed[(tp)]>0)
2920 #define NELIST_DEF_IN_MG(p) (MGFORMAT(p)->nodeelementlist)
2921 #define NDATA_DEF_IN_MG(p) (MGFORMAT(p)->nodedata)
2922 #define MG_SAVED(p) ((p)->saved)
2923 #define MG_FILENAME(p) ((p)->filename)
2924 #define MG_COARSE_FIXED(p) ((p)->CoarseGridFixed)
2925 #define MG_MARK_KEY(p) ((p)->MarkKey)
2926
2927 /****************************************************************************/
2928 /* */
2929 /* macros for formats */
2930 /* */
2931 /****************************************************************************/
2932
2933 #define FMT_S_VEC_TP(f,t) ((f)->VectorSizes[t])
2934 #define FMT_VTYPE_NAME(f,t) ((f)->VTypeNames[t])
2935 #define FMT_S_MAT_TP(f,t) ((f)->MatrixSizes[t])
2936 #define FMT_S_MATPTR(f) ((f)->MatrixSizes)
2937 #define FMT_S_IMAT_TP(f,t) ((f)->IMatrixSizes[t])
2938 #define FMT_CONN_DEPTH_TP(f,t) ((f)->ConnectionDepth[t])
2939 #define FMT_CONN_DEPTH_PTR(f) ((f)->ConnectionDepth)
2940 #define FMT_CONN_DEPTH_MAX(f) ((f)->MaxConnectionDepth)
2941 #define FMT_NB_DEPTH(f) ((f)->NeighborhoodDepth)
2942 #define FMT_PO2T(f,p,o) ((f)->po2t[p][o])
2943 #define FMT_T2P(f,t) ((f)->t2p[t])
2944 #define FMT_TYPE_IN_PART(f,t,o) ((f)->t2p[o] & (1<<o))
2945 #define FMT_T2O(f,o) ((f)->t2o[o])
2946 #define FMT_TYPE_USES_OBJ(f,t,o) ((f)->t2o[o] & (1<<o))
2947 #define FMT_USES_OBJ(f,o) ((f)->OTypeUsed[o])
2948 #define FMT_MAX_PART(f) ((f)->MaxPart)
2949 #define FMT_MAX_TYPE(f) ((f)->MaxType)
2950
2951 #define FMT_N2T(f,c) (((c)<FROM_VTNAME) ? NOVTYPE : ((c)>TO_VTNAME) ? NOVTYPE : (f)->n2t[(c)-FROM_VTNAME])
2952 #define FMT_SET_N2T(f,c,t) ((f)->n2t[(c)-FROM_VTNAME] = t)
2953 #define FMT_T2N(f,t) (((f)->t2n[t]))
2954
2955 /** \brief Constants for USED flags of objects */
2956 enum {MG_ELEMUSED = 1,
2957 MG_NODEUSED = 2,
2958 MG_EDGEUSED = 4,
2959 MG_VERTEXUSED = 8,
2960 MG_VECTORUSED = 16,
2961 MG_MATRIXUSED = 32};
2962
2963
2964 /****************************************************************************/
2965 /* */
2966 /* declaration of exported global variables */
2967 /* */
2968 /****************************************************************************/
2969
2970 #if defined ModelP && defined __OVERLAP2__
2971 extern INT ce_NO_DELETE_OVERLAP2;
2972 #endif
2973
2974 /****************************************************************************/
2975 /* */
2976 /* interface functions for module grid manager */
2977 /* */
2978 /****************************************************************************/
2979
2980 /** \brief Return values for functions returning an INT. The usual rule is: 0 ok, >0 error */
2981 enum {GM_OK = 0,
2982 GM_ERROR = 1,
2983 GM_FILEOPEN_ERROR = 2,
2984 GM_RULE_WITH_ORIENTATION = 3,
2985 GM_RULE_WITHOUT_ORIENTATION = 4,
2986 GM_OUT_OF_MEM = 5,
2987 GM_OUT_OF_RANGE = 6,
2988 GM_NOT_FOUND = 7,
2989 GM_INCONSISTENCY = 8,
2990 GM_COARSE_NOT_FIXED = 9,
2991 GM_FATAL = 999};
2992
2993
2994 /** @name Some constants passed as parameters */
2995 /*@{*/
2996 enum {GM_KEEP_BOUNDARY_NODES,
2997 GM_MOVE_BOUNDARY_NODES,
2998 GM_REFINE_TRULY_LOCAL,
2999 GM_COPY_ALL,
3000 GM_REFINE_NOT_CLOSED};
3001
3002 enum {GM_REFINE_PARALLEL, GM_REFINE_SEQUENTIAL};
3003
3004 enum {GM_REFINE_NOHEAPTEST, GM_REFINE_HEAPTEST};
3005
3006 enum {GM_FCFCLL = 1,
3007 GM_FFCCLL = 2,
3008 GM_FFLLCC = 3,
3009 GM_FFLCLC = 4,
3010 GM_CCFFLL = 5};
3011
3012 enum {GM_LOV_BEGIN = 1,
3013 GM_LOV_END = 2};
3014
3015 enum {GM_GEN_FIRST, GM_GEN_LAST, GM_GEN_CUT};
3016
3017 enum {GM_ALL_LEVELS = 1,
3018 GM_CURRENT_LEVEL = 2};
3019
3020 enum {GM_ORDER_IN_COLS, GM_ORDER_IN_ROWS};
3021
3022 /*@}*/
3023
3024 /* get/set current multigrid, loop through multigrids */
3025 MULTIGRID *MakeMGItem (const char *name, std::shared_ptr<PPIF::PPIFContext> context);
3026 MULTIGRID *GetMultigrid (const char *name);
3027 MULTIGRID *GetFirstMultigrid (void);
3028 MULTIGRID *GetNextMultigrid (const MULTIGRID *theMG);
3029
3030 /* format definition */
3031 std::unique_ptr<FORMAT> CreateFormat ();
3032
3033 /* create, saving and disposing a multigrid structure */
3034 MULTIGRID *CreateMultiGrid (char *MultigridName, char *BndValProblem,
3035 const char *format,
3036 INT optimizedIE, INT insertMesh,
3037 std::shared_ptr<PPIF::PPIFContext> ppifContext = nullptr);
3038 MULTIGRID *OpenMGFromDataFile(MULTIGRID *theMG, INT number, char *type,
3039 char *DataFileName, NS_PREFIX MEM heapSize);
3040 MULTIGRID *LoadMultiGrid (const char *MultigridName, const char *name, const char *type,
3041 const char *BndValProblem, const char *format,
3042 unsigned long heapSize,INT force,INT optimizedIE, INT autosave,
3043 std::shared_ptr<PPIF::PPIFContext> ppifContext = nullptr);
3044 INT SaveMultiGrid (MULTIGRID *theMG, const char *name, const char *type, const char *comment, INT autosave, INT rename);
3045 INT DisposeGrid (GRID *theGrid);
3046 INT DisposeMultiGrid (MULTIGRID *theMG);
3047 INT Collapse (MULTIGRID *theMG);
3048
3049 /* coarse grid manipulations */
3050 NODE *InsertInnerNode (GRID *theGrid, const DOUBLE *pos);
3051 NODE *InsertBoundaryNode (GRID *theGrid, BNDP *bndp);
3052
3053 INT DeleteNode (GRID *theGrid, NODE *theNode);
3054 ELEMENT *InsertElement (GRID *theGrid, INT n, NODE **NodeList, ELEMENT **ElemList, INT *NbgSdList, INT *bnds_flag);
3055 INT InsertMesh (MULTIGRID *theMG, MESH *theMesh);
3056 INT DeleteElement (MULTIGRID *theMG, ELEMENT *theElement);
3057
3058 /* refinement */
3059 /** \todo !!! should be moved to rm.h [Thimo] */
3060 INT EstimateHere (const ELEMENT *theElement);
3061 INT MarkForRefinement (ELEMENT *theElement, enum RefinementRule rule, INT data);
3062 INT GetRefinementMark (ELEMENT *theElement, INT *rule, void *data);
3063 INT GetRefinementMarkType (ELEMENT *theElement);
3064 INT AdaptMultiGrid (MULTIGRID *theMG, INT flag, INT seq, INT mgtest);
3065 INT TestRefineInfo (MULTIGRID *theMG);
3066 INT SetRefineInfo (MULTIGRID *theMG);
3067
3068
3069 /* moving nodes */
3070 #ifdef __THREEDIM__
3071 INT GetSideIDFromScratch (ELEMENT *theElement, NODE *theNode);
3072 #endif
3073
3074 /* algebraic connections */
3075 CONNECTION *CreateExtraConnection (GRID *theGrid, VECTOR *from, VECTOR *to);
3076 INT DisposeExtraConnections (GRID *theGrid);
3077 INT DisposeConnectionsInGrid (GRID *theGrid);
3078 MATRIX *GetMatrix (const VECTOR *FromVector, const VECTOR *ToVector);
3079 CONNECTION *GetConnection (const VECTOR *FromVector, const VECTOR *ToVector);
3080 INT GetAllVectorsOfElement (GRID *theGrid, ELEMENT *theElement,
3081 VECTOR **vec);
3082
3083 /* searching */
3084 ELEMENT *FindElementOnSurface (MULTIGRID *theMG, DOUBLE *global);
3085 INT InnerBoundary (ELEMENT *t, INT side);
3086
3087 /* list */
3088 void ListMultiGridHeader (const INT longformat);
3089 void ListMultiGrid (const MULTIGRID *theMG, const INT isCurrent, const INT longformat);
3090 INT MultiGridStatus (const MULTIGRID *theMG, INT gridflag, INT greenflag, INT lbflag, INT verbose);
3091 void ListGrids (const MULTIGRID *theMG);
3092 void ListNode (const MULTIGRID *theMG, const NODE *theNode, INT dataopt, INT bopt, INT nbopt, INT vopt);
3093 void ListElement (const MULTIGRID *theMG, const ELEMENT *theElement, INT dataopt, INT bopt, INT nbopt, INT vopt);
3094 void ListVector (const MULTIGRID *theMG, const VECTOR *theVector, INT matrixopt, INT dataopt, INT modifiers);
3095
3096 /* query */
3097 LINK *GetLink (const NODE *from, const NODE *to);
3098 EDGE *GetSonEdge (const EDGE *theEdge);
3099 INT GetSonEdges (const EDGE *theEdge, EDGE *SonEdges[MAX_SON_EDGES]);
3100 EDGE *GetFatherEdge (const EDGE *theEdge);
3101 #ifdef __THREEDIM__
3102 EDGE *FatherEdge (NODE **SideNodes, INT ncorners, NODE **Nodes, EDGE *theEdge);
3103 #endif
3104 EDGE *GetEdge (const NODE *from, const NODE *to);
3105 INT GetSons (const ELEMENT *theElement, ELEMENT *SonList[MAX_SONS]);
3106 #ifdef ModelP
3107 INT GetAllSons (const ELEMENT *theElement, ELEMENT *SonList[MAX_SONS]);
3108 #endif
3109 INT VectorPosition (const VECTOR *theVector, DOUBLE *position);
3110 INT VectorInElement (ELEMENT *theElement, VECTOR *theVector);
3111
3112 /* check */
3113 #ifndef ModelP
3114 INT CheckGrid (GRID *theGrid, INT checkgeom, INT checkalgebra, INT checklists);
3115 #else
3116 INT CheckGrid (GRID *theGrid, INT checkgeom, INT checkalgebra, INT checklists, INT checkif);
3117 #endif
3118 INT CheckLists (GRID *theGrid);
3119 INT CheckSubdomains (MULTIGRID *theMG);
3120
3121 /* multigrid user data space management (using the heaps.c block heap management) */
3122 INT AllocateControlEntry (INT cw_id, INT length, INT *ce_id);
3123 INT FreeControlEntry (INT ce_id);
3124 void ListCWofObject (const void *obj, INT offset);
3125 void ListAllCWsOfObject (const void *obj);
3126 void ListAllCWsOfAllObjectTypes (PrintfProcPtr myprintf);
3127 UINT ReadCW (const void *obj, INT ce);
3128 void WriteCW (void *obj, INT ce, INT n);
3129
3130 /* ordering of degrees of freedom */
3131 ALG_DEP *CreateAlgebraicDependency (const char *name, DependencyProcPtr DependencyProc);
3132 FIND_CUT *CreateFindCutProc (const char *name, FindCutProcPtr FindCutProc);
3133
3134 /* miscellaneous */
3135 INT RenumberMultiGrid (MULTIGRID *theMG, INT *nboe, INT *nioe, INT *nbov, INT *niov, NODE ***vid_n, INT *foid, INT *non, INT MarkKey);
3136 INT MGSetVectorClasses (MULTIGRID *theMG);
3137 INT SetEdgeSubdomainFromElements (GRID *theGrid);
3138 INT SetSubdomainIDfromBndInfo (MULTIGRID *theMG);
3139 INT FixCoarseGrid (MULTIGRID *theMG);
3140 INT ClearMultiGridUsedFlags (MULTIGRID *theMG, INT FromLevel, INT ToLevel, INT mask);
3141 void CalculateCenterOfMass (ELEMENT *theElement, DOUBLE_VECTOR center_of_mass);
3142 INT KeyForObject (KEY_OBJECT *obj);
3143
3144 /** \todo remove the following functions after the code will never need any debugging */
3145 char *PrintElementInfo (ELEMENT *theElement,INT full);
3146
3147 END_UGDIM_NAMESPACE
3148
3149 #endif
3150