Loading [MathJax]/extensions/tex2jax.js
ug4
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
hfv1_geom.h
Go to the documentation of this file.
1/*
2 * Copyright (c) 2010-2015: G-CSC, Goethe University Frankfurt
3 * Author: Andreas Vogel
4 *
5 * This file is part of UG4.
6 *
7 * UG4 is free software: you can redistribute it and/or modify it under the
8 * terms of the GNU Lesser General Public License version 3 (as published by the
9 * Free Software Foundation) with the following additional attribution
10 * requirements (according to LGPL/GPL v3 §7):
11 *
12 * (1) The following notice must be displayed in the Appropriate Legal Notices
13 * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 *
15 * (2) The following notice must be displayed at a prominent place in the
16 * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 *
18 * (3) The following bibliography is recommended for citation and must be
19 * preserved in all covered files:
20 * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 * parallel geometric multigrid solver on hierarchically distributed grids.
22 * Computing and visualization in science 16, 4 (2013), 151-164"
23 * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 * flexible software system for simulating pde based models on high performance
25 * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 *
27 * This program is distributed in the hope that it will be useful,
28 * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 * GNU Lesser General Public License for more details.
31 */
32
33#ifndef __H__UG__LIB_DISC__SPATIAL_DISC__DISC_HELPER__HANGING_FINITE_VOLUME_GEOMETRY__
34#define __H__UG__LIB_DISC__SPATIAL_DISC__DISC_HELPER__HANGING_FINITE_VOLUME_GEOMETRY__
35
36// extern libraries
37#include <cmath>
38#include <vector>
39
40// other ug4 modules
41#include "common/common.h"
42
43// library intern includes
46#include "fv_util.h"
47#include "fv1_geom.h"
48
49namespace ug{
50
52// Hanging node FV1 Geometry for Reference Element Type
54
56
60template < typename TElem, int TWorldDim>
62 private:
65
67 static const size_t m_numNaturalSCV = ref_elem_type::numCorners;
68
70 static const size_t m_numNaturalSCVF = ref_elem_type::numEdges;
71
72 public:
74 static const int dim = ref_elem_type::dim;
75
77 static const int worldDim = TWorldDim;
78
80 static const bool usesHangingNodes = true;
81
83 static const bool staticLocalData = true;
84
85 protected:
86 struct MidID
87 {
88 MidID() : dim(0), id(0) {};
89 MidID(size_t dim_, size_t id_) : dim(dim_), id(id_) {};
90 size_t dim;
91 size_t id;
92 };
93
94 public:
95 class SCVF
96 {
97 public:
99 static const int dim = ref_elem_type::dim;
100
102 static const int worldDim = TWorldDim;
103
104 private:
106 friend class HFV1Geometry<TElem, TWorldDim>;
107
109 static const size_t m_numIP = 1;
110
113
114 public:
115 SCVF() {};
116
119 inline size_t from() const {return m_from;}
120
123 inline size_t to() const {return m_to;}
124
126 inline size_t num_ip() const {return m_numIP;}
127
129 inline const MathVector<dim>& local_ip() const {return localIP;}
130
132 inline const MathVector<worldDim>& global_ip() const {return globalIP;}
133
135 inline const MathVector<worldDim>& normal() const {return Normal;} // includes area
136
138 inline const MathMatrix<worldDim,dim>& JTInv() const {return JtInv;}
139
141 inline number detJ() const {return detj;}
142
144 inline size_t num_sh() const {return vShape.size();}
145
147 inline number shape(size_t sh) const
148 {UG_ASSERT(sh < num_sh(), "Invalid ip index"); return vShape[sh];}
149
151 inline const MathVector<dim>& local_grad(size_t sh) const
152 {UG_ASSERT(sh < num_sh(), "Invalid ip index"); return localGrad[sh];}
153
155 inline const MathVector<dim>* local_grad_vector() const {return &localGrad[0];}
156
158 inline const MathVector<worldDim>& global_grad(size_t sh) const
159 {UG_ASSERT(sh < num_sh(), "Invalid ip index"); return globalGrad[sh];}
160
162 inline const MathVector<worldDim>* global_grad_vector() const {return &globalGrad[0];}
163
165 inline size_t num_corners() const {return m_numCorners;}
166
168 inline const MathVector<dim>& local_corner(size_t co) const
169 {UG_ASSERT(co < num_corners(), "Invalid corner index."); return m_vLocPos[co];}
170
172 inline const MathVector<worldDim>& global_corner(size_t co) const
173 {UG_ASSERT(co < num_corners(), "Invalid corner index."); return m_vGloPos[co];}
174
175 private:
176 // This scvf separates the scv with the ids given in "from" and "to"
177 // The computed normal points in direction from->to
178 size_t m_from, m_to;
179
180 // ordering is:
181 // 2D: edgeMidPoint, CenterOfElement
182 // 3D: edgeMidPoint, Side one, CenterOfElement, Side two
183 MathVector<dim> m_vLocPos[m_numCorners]; // local corners of scvf
184 MathVector<worldDim> m_vGloPos[m_numCorners]; // global corners of scvf
185 MidID m_midId[m_numCorners]; // dimension and id of object, that's midpoint bounds the scvf
186
187 // scvf part
188 MathVector<dim> localIP; // local integration point
189 MathVector<worldDim> globalIP; // global intergration point
190 MathVector<worldDim> Normal; // normal (incl. area)
191
192 // shapes and derivatives
193 std::vector<number> vShape; // shapes at ip
194 std::vector<MathVector<dim> > localGrad; // local grad at ip
195 std::vector<MathVector<worldDim> > globalGrad; // global grad at ip
196 MathMatrix<worldDim,dim> JtInv; // Jacobian transposed at ip
197 number detj; // Jacobian det at ip
198 };
199
200 class SCV
201 {
202 private:
203 // let outer class access private members
204 friend class HFV1Geometry<TElem, TWorldDim>;
205
207 static const size_t m_numIP = 1;
208
211
214
215 public:
217
219 inline size_t node_id() const {return nodeId;}
220
222 inline size_t num_ip() const {return m_numIP;}
223
225 inline const MathVector<dim>& local_ip() const {return m_vLocPos[0];}
226
228 inline const MathVector<worldDim>& global_ip() const {return m_vGloPos[0];}
229
231 inline number volume() const {return vol;}
232
234 inline size_t num_corners() const {return m_numCorners;}
235
237 inline const MathVector<dim>& local_corner(size_t co) const
238 {UG_ASSERT(co < num_corners(), "Invalid corner index."); return m_vLocPos[co];}
239
241 inline const MathVector<worldDim>& global_corner(size_t co) const
242 {UG_ASSERT(co < num_corners(), "Invalid corner index."); return m_vGloPos[co];}
243
244 private:
245 size_t nodeId; // node id of associated node
247
248 // The ordering is: Corner, ...
249 MathVector<dim> m_vLocPos[m_maxNumCorners]; // local position of node
250 MathVector<worldDim> m_vGloPos[m_maxNumCorners]; // global position of node
251 MidID m_midId[m_maxNumCorners]; // dimension and id of object, that's midpoint bounds the scv
253 };
254
255
256 public:
258 HFV1Geometry();
259
261 void update(GridObject* pElem, const MathVector<worldDim>* vCornerCoords,
262 const ISubsetHandler* ish = NULL);
263
264 // debug output
265 void print();
266
268 GridObject* elem() const {return m_pElem;}
269
270 public:
272 inline size_t num_scvf() const {return m_vSCVF.size();}
273
275 inline const SCVF& scvf(size_t i) const
276 {UG_ASSERT(i < num_scvf(), "Invalid Index."); return m_vSCVF[i];}
277
279 size_t num_scv() const {return m_vSCV.size();}
280
282 inline const SCV& scv(size_t i) const
283 {UG_ASSERT(i < num_scv(), "Invalid Index."); return m_vSCV[i];}
284
286 inline size_t num_sh() const
287 {
288 return m_numSh;
289 };
290
291 public:
293 size_t num_scvf_ips() const {return m_vGlobSCVFIP.size();}
294
297
299 const MathVector<dim>* scvf_local_ips() const {return &m_vLocSCVFIP[0];}
300
302 size_t num_scv_ips() const {return m_vGlobSCVIP.size();}
303
306
308 const MathVector<dim>* scv_local_ips() const {return &m_vLocSCVIP[0];}
309
311 const MathVector<dim>& local_node_position(size_t nodeID) const
312 {
313 UG_ASSERT(nodeID < m_locMid[0].size(), "Invalid node id.");
314 return m_locMid[0][nodeID];
315 }
316
318 const MathVector<worldDim>& global_node_position(size_t nodeID) const
319 {
320 UG_ASSERT(nodeID < m_gloMid[0].size(), "Invalid node id.");
321 return m_gloMid[0][nodeID];
322 }
323
324 protected:
325 std::vector<MathVector<worldDim> > m_vGlobSCVFIP;
326 std::vector<MathVector<dim> > m_vLocSCVFIP;
327 std::vector<MathVector<worldDim> > m_vGlobSCVIP;
328 std::vector<MathVector<dim> > m_vLocSCVIP;
329
330 protected:
332 {
333 for(size_t i = 0; i < scvf.num_corners(); ++i)
334 {
335 const size_t dim = scvf.m_midId[i].dim;
336 const size_t id = scvf.m_midId[i].id;
337 scvf.m_vLocPos[i] = m_locMid[dim][id];
338 }
339 }
340
342 {
343 for(size_t i = 0; i < scvf.num_corners(); ++i)
344 {
345 const size_t dim = scvf.m_midId[i].dim;
346 const size_t id = scvf.m_midId[i].id;
347 scvf.m_vGloPos[i] = m_gloMid[dim][id];
348 }
349 }
350
352 {
353 for(size_t i = 0; i < scv.num_corners(); ++i)
354 {
355 const size_t dim_ = scv.m_midId[i].dim;
356 const size_t id = scv.m_midId[i].id;
357 UG_ASSERT(dim_ <= (size_t)dim, "Dimension wrong");
358 UG_ASSERT(id < m_locMid[dim_].size(), "id " << id << " in dim="
359 <<dim_<<" wrong. (size is "<< m_locMid[dim_].size()<<")\n");
360 scv.m_vLocPos[i] = m_locMid[dim_][id];
361 }
362 }
363
365 {
366 for(size_t i = 0; i < scv.num_corners(); ++i)
367 {
368 const size_t dim_ = scv.m_midId[i].dim;
369 const size_t id = scv.m_midId[i].id;
370 UG_ASSERT(dim_ <= (size_t)dim, "Dimension wrong");
371 UG_ASSERT(id < m_gloMid[dim_].size(), "id " << id << " in dim="
372 <<dim_<<" wrong. (size is "<< m_gloMid[dim_].size()<<")\n");
373 scv.m_vGloPos[i] = m_gloMid[dim_][id];
374 }
375 }
376
377 // i = number of side
379 MathVector<dim>& locSideMid,
380 MathVector<worldDim>& gloSideMid)
381 {
382 if(dim>1)
383 {
384 const size_t coID0 = m_rRefElem.id(2, i, 0, 0);
385 locSideMid = m_locMid[0][coID0];
386 gloSideMid = m_gloMid[0][coID0];
387
388 // add corner coordinates of the corners of the geometric object
389 for(size_t j = 1; j < m_rRefElem.num(2, i, 0); ++j)
390 {
391 const size_t coID = m_rRefElem.id(2, i, 0, j);
392 locSideMid += m_locMid[0][coID];
393 gloSideMid += m_gloMid[0][coID];
394 }
395
396 // scale for correct averaging
397 locSideMid *= 1./(m_rRefElem.num(2, i, 0));
398 gloSideMid *= 1./(m_rRefElem.num(2, i, 0));
399 }
400 }
401
402 // i, j, k, l = number nodes
403 void compute_side_midpoints(size_t i, size_t j, size_t k, size_t l,
404 MathVector<dim>& locSideMid,
405 MathVector<worldDim>& gloSideMid)
406 {
407 VecAdd(locSideMid, m_locMid[0][i], m_locMid[0][j], m_locMid[0][k], m_locMid[0][l]);
408 VecAdd(gloSideMid, m_gloMid[0][i], m_gloMid[0][j], m_gloMid[0][k], m_gloMid[0][l]);
409
410 // scale for correct averaging
411 locSideMid *= 0.25;
412 gloSideMid *= 0.25;
413 }
414
415 // i, j, k = number nodes
416 void compute_side_midpoints(size_t i, size_t j, size_t k,
417 MathVector<dim>& locSideMid, MathVector<worldDim>& gloSideMid)
418 {
419 VecAdd(locSideMid, m_locMid[0][i], m_locMid[0][j], m_locMid[0][k]);
420 VecAdd(gloSideMid, m_gloMid[0][i], m_gloMid[0][j], m_gloMid[0][k]);
421
422 // scale for correct averaging
423 locSideMid *= 1./3.;
424 gloSideMid *= 1./3.;
425 }
426
427 // returns edgeID of child edge, that has corner co. i is the natural edge index
428 size_t get_child_edge_of_corner(size_t i, size_t co)
429 {
430 for(size_t e = 0; e < m_vNatEdgeInfo[i].num_child_edges(); ++e)
431 {
432 const size_t childId = m_vNatEdgeInfo[i].child_edge(e);
433 if(m_vNewEdgeInfo[childId].from() == co)
434 return childId;
435 if(m_vNewEdgeInfo[childId].to() == co)
436 return childId;
437 }
438 return -1;
439 }
440
441 private:
443 {
444 friend class HFV1Geometry<TElem, TWorldDim>;
445 public:
446 NewEdgeInfo() : m_from(-1), m_to(-1){}
447 size_t from() const {return m_from;}
448 size_t to() const {return m_to;}
449 private:
450 size_t m_from;
451 size_t m_to;
452 };
453
455 {
456 friend class HFV1Geometry<TElem, TWorldDim>;
457 public:
458 NatEdgeInfo() : nodeId(-1), numChildEdges(0) {for(size_t i=0; i<2; ++i) childEdge[i] = -1;}
459 bool is_hanging() const {return numChildEdges == 2;}
460 size_t node_id() const {UG_ASSERT(is_hanging(), "Should only be called, if edge is hanging."); return nodeId;}
461 size_t num_child_edges() const {return numChildEdges;}
462 size_t child_edge(size_t i) const {UG_ASSERT(i < num_child_edges(), "Wrong id."); return childEdge[i];}
463
464 private:
465 size_t nodeId;
467 size_t childEdge[2];
468 };
469
470 // help array: maps NaturalEdge -> NodeOnEdge (-1 if no node on edge)
471 std::vector<NatEdgeInfo> m_vNatEdgeInfo;
472 std::vector<NewEdgeInfo> m_vNewEdgeInfo;
473
474 private:
475 // pointer to current element
477
478 std::vector<MathVector<dim> > m_locMid[dim+1];
479 std::vector<MathVector<worldDim> > m_gloMid[dim+1];
480
481 // SubControlVolumeFaces
482 std::vector<SCVF> m_vSCVF;
483
484 // SubControlVolumes
485 std::vector<SCV> m_vSCV;
486
487 // number of shape functions
488 size_t m_numSh;
489
490 // Reference Mapping
492
493 // Reference Element
495};
496
497template < int TDim, int TWorldDim = TDim>
499 private:
502
505
508
509 public:
511 static const int dim = TDim;
512
514 static const int worldDim = TWorldDim;
515
517 static const bool usesHangingNodes = true;
518
520 static const bool staticLocalData = true;
521
524
526 typedef typename traits::elem_type_0 elem_type_0;
527 typedef typename traits::elem_type_1 elem_type_1;
528 typedef typename traits::elem_type_2 elem_type_2;
529 typedef typename traits::elem_type_3 elem_type_3;
530 typedef typename traits::elem_type_4 elem_type_4;
531
532 protected:
533 struct MidID
534 {
535 MidID() : dim(0), id(0) {};
536 MidID(size_t dim_, size_t id_) : dim(dim_), id(id_) {};
537 size_t dim;
538 size_t id;
539 };
540
541 public:
542 class SCVF
543 {
544 public:
546 static const int dim = TDim;
547
549 static const int worldDim = TWorldDim;
550
551 private:
553 friend class DimHFV1Geometry<TDim, TWorldDim>;
554
556 static const size_t m_numIP = 1;
557
559 static const size_t m_numCorners = traits::NumCornersOfSCVF;
560
561 public:
562 SCVF() {};
563
565 inline size_t from() const {return m_from;}
566
568 inline size_t to() const {return m_to;}
569
571 inline size_t num_ip() const {return m_numIP;}
572
574 inline const MathVector<dim>& local_ip() const {return localIP;}
575
577 inline const MathVector<worldDim>& global_ip() const {return globalIP;}
578
580 inline const MathVector<worldDim>& normal() const {return Normal;} // includes area
581
583 inline const MathMatrix<worldDim,dim>& JTInv() const {return JtInv;}
584
586 inline number detJ() const {return detj;}
587
589 inline size_t num_sh() const {return vShape.size();}
590
592 inline number shape(size_t sh) const
593 {UG_ASSERT(sh < num_sh(), "Invalid ip index"); return vShape[sh];}
594
596 inline const MathVector<dim>& local_grad(size_t sh) const
597 {UG_ASSERT(sh < num_sh(), "Invalid ip index"); return localGrad[sh];}
598
600 inline const MathVector<dim>* local_grad_vector() const {return &localGrad[0];}
601
603 inline const MathVector<worldDim>& global_grad(size_t sh) const
604 {UG_ASSERT(sh < num_sh(), "Invalid ip index"); return globalGrad[sh];}
605
607 inline const MathVector<worldDim>* global_grad_vector() const {return &globalGrad[0];}
608
610 inline size_t num_corners() const {return m_numCorners;}
611
613 inline const MathVector<dim>& local_corner(size_t co) const
614 {UG_ASSERT(co < num_corners(), "Invalid corner index."); return m_vLocPos[co];}
615
617 inline const MathVector<worldDim>& global_corner(size_t co) const
618 {UG_ASSERT(co < num_corners(), "Invalid corner index."); return m_vGloPos[co];}
619
620 private:
621 // This scvf separates the scv with the ids given in "from" and "to"
622 // The computed normal points in direction from->to
623 size_t m_from, m_to;
624
625 // ordering is:
626 // 2D: edgeMidPoint, CenterOfElement
627 // 3D: edgeMidPoint, Side one, CenterOfElement, Side two
628 MathVector<dim> m_vLocPos[m_numCorners]; // local corners of scvf
629 MathVector<worldDim> m_vGloPos[m_numCorners]; // global corners of scvf
630 MidID m_midId[m_numCorners]; // dimension and id of object, that's midpoint bounds the scvf
631
632 // scvf part
633 MathVector<dim> localIP; // local integration point
634 MathVector<worldDim> globalIP; // global intergration point
635 MathVector<worldDim> Normal; // normal (incl. area)
636
637 // shapes and derivatives
638 std::vector<number> vShape; // shapes at ip
639 std::vector<MathVector<dim> > localGrad; // local grad at ip
640 std::vector<MathVector<worldDim> > globalGrad; // global grad at ip
641 MathMatrix<worldDim,dim> JtInv; // Jacobian transposed at ip
642 number detj; // Jacobian det at ip
643 };
644
645 class SCV
646 {
647 private:
648 // let outer class access private members
649 friend class DimHFV1Geometry<TDim, TWorldDim>;
650
652 static const size_t m_numIP = 1;
653
655 static const size_t m_maxNumCorners = traits::MaxNumCornersOfSCV;
656
658 typedef typename traits::scv_type scv_type;
659
660 public:
662
664 inline size_t node_id() const {return nodeId;}
665
667 inline size_t num_ip() const {return m_numIP;}
668
670 inline const MathVector<dim>& local_ip() const {return m_vLocPos[0];}
671
673 inline const MathVector<worldDim>& global_ip() const {return m_vGloPos[0];}
674
676 inline number volume() const {return vol;}
677
679 inline size_t num_corners() const {return m_numCorners;}
680
682 inline const MathVector<dim>& local_corner(size_t co) const
683 {UG_ASSERT(co < num_corners(), "Invalid corner index."); return m_vLocPos[co];}
684
686 inline const MathVector<worldDim>& global_corner(size_t co) const
687 {UG_ASSERT(co < num_corners(), "Invalid corner index."); return m_vGloPos[co];}
688
690 inline size_t num_sh() const {return numSH;}
691
693 inline number shape(size_t sh) const {return vShape[sh];}
694
696 inline const number* shape_vector() const {return vShape;}
697
699 inline const MathVector<dim>& local_grad(size_t sh) const
700 {UG_ASSERT(sh < num_sh(), "Invalid index"); return localGrad[sh];}
701
703 inline const MathVector<dim>* local_grad_vector() const {return localGrad;}
704
706 inline const MathVector<worldDim>& global_grad(size_t sh) const
707 {UG_ASSERT(sh < num_sh(), "Invalid index"); return globalGrad[sh];}
708
710 inline const MathVector<worldDim>* global_grad_vector() const {return globalGrad;}
711
713 inline const MathMatrix<worldDim,dim>& JTInv() const {return JtInv;}
714
716 inline number detJ() const {return detj;}
717
718 private:
719 size_t nodeId; // node id of associated node
721
722 // The ordering is: Corner, ...
723 MathVector<dim> m_vLocPos[m_maxNumCorners]; // local position of node
724 MathVector<worldDim> m_vGloPos[m_maxNumCorners]; // global position of node
725 MidID m_midId[m_maxNumCorners]; // dimension and id of object, that's midpoint bounds the scv
727
728 // shapes and derivatives
729 size_t numSH;
730 number vShape[m_maxNSH]; // shapes at ip
731 MathVector<dim> localGrad[m_maxNSH]; // local grad at ip
733 MathMatrix<worldDim,dim> JtInv; // Jacobian transposed at ip
734 number detj; // Jacobian det at ip
735 };
736
737
738 public:
741
743 void update_local_data();
744
746 void update(GridObject* pElem, const MathVector<worldDim>* vCornerCoords,
747 const ISubsetHandler* ish = NULL);
748
749 // debug output
750 // void print();
751
752 public:
754 inline size_t num_scvf() const {return m_vSCVF.size();}
755
757 inline const SCVF& scvf(size_t i) const
758 {UG_ASSERT(i < num_scvf(), "Invalid Index."); return m_vSCVF[i];}
759
761 size_t num_scv() const {return m_vSCV.size();}
762
764 inline const SCV& scv(size_t i) const
765 {UG_ASSERT(i < num_scv(), "Invalid Index."); return m_vSCV[i];}
766
768 inline size_t num_sh() const
769 {
770 return m_numSh;
771 };
772
773 public:
775 size_t num_scvf_ips() const {return m_vGlobSCVFIP.size();}
776
779
781 const MathVector<dim>* scvf_local_ips() const {return &m_vLocSCVFIP[0];}
782
784 size_t num_scv_ips() const {return m_vGlobSCVIP.size();}
785
788
790 const MathVector<dim>* scv_local_ips() const {return &m_vLocSCVIP[0];}
791
792 protected:
793 std::vector<MathVector<worldDim> > m_vGlobSCVFIP;
794 std::vector<MathVector<dim> > m_vLocSCVFIP;
795 std::vector<MathVector<worldDim> > m_vGlobSCVIP;
796 std::vector<MathVector<dim> > m_vLocSCVIP;
797
798 protected:
800 {
801 for(size_t i = 0; i < scvf.num_corners(); ++i)
802 {
803 const size_t dim = scvf.m_midId[i].dim;
804 const size_t id = scvf.m_midId[i].id;
805 scvf.m_vLocPos[i] = m_locMid[dim][id];
806 }
807 }
808
810 {
811 for(size_t i = 0; i < scvf.num_corners(); ++i)
812 {
813 const size_t dim = scvf.m_midId[i].dim;
814 const size_t id = scvf.m_midId[i].id;
815 scvf.m_vGloPos[i] = m_gloMid[dim][id];
816 }
817 }
818
820 {
821 for(size_t i = 0; i < scv.num_corners(); ++i)
822 {
823 const size_t dim_ = scv.m_midId[i].dim;
824 const size_t id = scv.m_midId[i].id;
825 UG_ASSERT(dim_ <= (size_t)dim, "Dimension wrong");
826 UG_ASSERT(id < m_locMid[dim_].size(), "id " << id << " in dim="
827 <<dim_<<" wrong. (size is "<< m_locMid[dim_].size()<<")\n");
828 scv.m_vLocPos[i] = m_locMid[dim_][id];
829 }
830 }
831
833 {
834 for(size_t i = 0; i < scv.num_corners(); ++i)
835 {
836 const size_t dim_ = scv.m_midId[i].dim;
837 const size_t id = scv.m_midId[i].id;
838 UG_ASSERT(dim_ <= (size_t)dim, "Dimension wrong");
839 UG_ASSERT(id < m_gloMid[dim_].size(), "id " << id << " in dim="
840 <<dim_<<" wrong. (size is "<< m_gloMid[dim_].size()<<")\n");
841 scv.m_vGloPos[i] = m_gloMid[dim_][id];
842 }
843 }
844
845 // i = number of side
847 MathVector<dim>& locSideMid,
848 MathVector<worldDim>& gloSideMid)
849 {
850 if(dim>1)
851 {
852 const size_t coID0 = m_rRefElem.id(2, i, 0, 0);
853 locSideMid = m_locMid[0][coID0];
854 gloSideMid = m_gloMid[0][coID0];
855
856 // add corner coordinates of the corners of the geometric object
857 for(size_t j = 1; j < m_rRefElem.num(2, i, 0); ++j)
858 {
859 const size_t coID = m_rRefElem.id(2, i, 0, j);
860 locSideMid += m_locMid[0][coID];
861 gloSideMid += m_gloMid[0][coID];
862 }
863
864 // scale for correct averaging
865 locSideMid *= 1./(m_rRefElem.num(2, i, 0));
866 gloSideMid *= 1./(m_rRefElem.num(2, i, 0));
867 }
868 }
869
870 // i, j, k, l = number nodes
871 void compute_side_midpoints(size_t i, size_t j, size_t k, size_t l,
872 MathVector<dim>& locSideMid,
873 MathVector<worldDim>& gloSideMid)
874 {
875 VecAdd(locSideMid, m_locMid[0][i], m_locMid[0][j], m_locMid[0][k], m_locMid[0][l]);
876 VecAdd(gloSideMid, m_gloMid[0][i], m_gloMid[0][j], m_gloMid[0][k], m_gloMid[0][l]);
877
878 // scale for correct averaging
879 locSideMid *= 0.25;
880 gloSideMid *= 0.25;
881 }
882
883 // i, j, k = number nodes
884 void compute_side_midpoints(size_t i, size_t j, size_t k,
885 MathVector<dim>& locSideMid, MathVector<worldDim>& gloSideMid)
886 {
887 VecAdd(locSideMid, m_locMid[0][i], m_locMid[0][j], m_locMid[0][k]);
888 VecAdd(gloSideMid, m_gloMid[0][i], m_gloMid[0][j], m_gloMid[0][k]);
889
890 // scale for correct averaging
891 locSideMid *= 1./3.;
892 gloSideMid *= 1./3.;
893 }
894
895 // returns edgeID of child edge, that has corner co. i is the natural edge index
896 size_t get_child_edge_of_corner(size_t i, size_t co)
897 {
898 for(size_t e = 0; e < m_vNatEdgeInfo[i].num_child_edges(); ++e)
899 {
900 const size_t childId = m_vNatEdgeInfo[i].child_edge(e);
901 if(m_vNewEdgeInfo[childId].from() == co)
902 return childId;
903 if(m_vNewEdgeInfo[childId].to() == co)
904 return childId;
905 }
906 return -1;
907 }
908
909 private:
911 {
912 friend class DimHFV1Geometry<TDim, TWorldDim>;
913 public:
914 NewEdgeInfo() : m_from(-1), m_to(-1){}
915 size_t from() const {return m_from;}
916 size_t to() const {return m_to;}
917 private:
918 size_t m_from;
919 size_t m_to;
920 };
921
923 {
924 friend class DimHFV1Geometry<TDim, TWorldDim>;
925 public:
926 NatEdgeInfo() : nodeId(-1), numChildEdges(0) {for(size_t i=0; i<2; ++i) childEdge[i] = -1;}
927 bool is_hanging() const {return numChildEdges == 2;}
928 size_t node_id() const {UG_ASSERT(is_hanging(), "Should only be called, if edge is hanging."); return nodeId;}
929 size_t num_child_edges() const {return numChildEdges;}
930 size_t child_edge(size_t i) const {UG_ASSERT(i < num_child_edges(), "Wrong id."); return childEdge[i];}
931
932 private:
933 size_t nodeId;
935 size_t childEdge[2];
936 };
937
938 // help array: maps NaturalEdge -> NodeOnEdge (-1 if no node on edge)
939 std::vector<NatEdgeInfo> m_vNatEdgeInfo;
940 std::vector<NewEdgeInfo> m_vNewEdgeInfo;
941
942 private:
945
948
949 std::vector<MathVector<dim> > m_locMid[dim+1];
950 std::vector<MathVector<worldDim> > m_gloMid[dim+1];
951
952 // SubControlVolumeFaces
953 std::vector<SCVF> m_vSCVF;
954
955 // SubControlVolumes
956 std::vector<SCV> m_vSCV;
957
958 // number of shape functions
959 size_t m_numSh;
960
963
966
967};
968
969
971// Hanging FV1 Manifold Geometry
973
974template <typename TElem, int TWorldDim>
976{
977 public:
978 // type of element
979 typedef TElem elem_type;
980
981 // type of reference element
983
984 public:
985 // order
986 static const int order = 1;
987
988 // number of natural bfs (eq. of SCV)
989 static const size_t m_numNaturalBF = ref_elem_type::numCorners;
990
991 // number of natural boundary face sides (eq. of SCVF)
992 static const size_t m_numNaturalBFS = ref_elem_type::numEdges;
993
994 // dimension of reference element
995 static const int dim = ref_elem_type::dim;
996
997 // dimension of world
998 static const int worldDim = TWorldDim;
999
1000 // type of BoundaryFaces
1002
1003 // Hanging node flag: this geometry supports hanging nodes
1004 static const bool usesHangingNodes = true;
1005
1007 static const bool staticLocalData = true;
1008
1009 protected:
1010 struct MidID
1011 {
1012 MidID() : dim(0), id(0) {};
1013 MidID(size_t dim_, size_t id_) : dim(dim_), id(id_) {};
1014 size_t dim;
1015 size_t id;
1016 };
1017
1018 public:
1019 class BF
1020 {
1021 private:
1022 // let outer class access private members
1023 friend class HFV1ManifoldGeometry<TElem, TWorldDim>;
1024
1025 // number of integration points
1026 static const size_t m_numIP = 1;
1027
1028 // max number of corners of bf
1030
1031 public:
1032 BF() {};
1033
1035 inline size_t node_id() const {return nodeId;}
1036
1038 inline size_t num_ip() const {return m_numIP;}
1039
1041 inline const MathVector<dim>& local_ip() const //{return localIP;}
1042 {return m_vLocPos[0];} // <-- always the vertex
1043
1045 inline const MathVector<worldDim>& global_ip() const //{return globalIP;}
1046 {return m_vGloPos[0];} // <-- here too
1047
1049 inline number volume() const {return vol;}
1050
1052 inline size_t num_corners() const {return numCorners;}
1053
1055 inline const MathVector<dim>& local_corner(size_t i) const
1056 {UG_ASSERT(i < num_corners(), "Invalid index."); return m_vLocPos[i];}
1057
1059 inline const MathVector<worldDim>& global_corner(size_t i) const
1060 {UG_ASSERT(i < num_corners(), "Invalid index."); return m_vGloPos[i];}
1061
1063 inline size_t num_sh() const {return vShape.size();}
1064
1066 inline number shape(size_t i, size_t ip) const
1067 {UG_ASSERT(ip < num_ip(), "Invalid index"); return vShape[i];}
1068
1069 private:
1070 size_t nodeId; // id of associated node
1071
1072 // CORNERS: ordering is:
1073 // 1D: corner, centerOfElement
1074 // 2D: corner, side one, centerOfElement, side two
1075 MathVector<dim> m_vLocPos[numCorners]; // local position of node
1076 MathVector<worldDim> m_vGloPos[numCorners]; // global position of node
1077
1078 //IPs & shapes
1079 MathVector<dim> localIP; // local integration point
1080 MathVector<worldDim> globalIP; // global integration point
1081 std::vector<number> vShape; // shapes at ip
1082
1084 };
1085
1086 protected:
1087 std::vector<MathVector<dim> > m_vLocBFIP;
1088 std::vector<MathVector<worldDim> > m_vGlobBFIP;
1089
1090 public:
1093
1095 void update(GridObject* elem, const MathVector<worldDim>* vCornerCoords,
1096 const ISubsetHandler* ish = NULL);
1097
1099 void print();
1100/*
1102 const MathVector<worldDim>* corners() const {return m_gloMid[0];}
1103*/
1105 inline size_t num_bf() const {return m_vBF.size();}
1106
1108 inline const BF& bf(size_t i) const
1109 {UG_ASSERT(i < num_bf(), "Invalid Index."); return m_vBF[i];}
1110
1113
1115 size_t num_bf_global_ips() const {return m_vGlobBFIP.size();}
1116
1118 const MathVector<dim>* bf_local_ips() const {return &m_vLocBFIP[0];}
1119
1121 size_t num_bf_local_ips() const {return m_vLocBFIP.size();}
1122
1123 protected:
1125 MathVector<worldDim>& gloSideMid)
1126 {
1127 locSideMid = m_locMid[0][0];
1128 gloSideMid = m_gloMid[0][0];
1129
1130 // add corner coordinates of the corners of the geometric object
1131 for (size_t j = 1; j < m_numNaturalBF; ++j)
1132 {
1133 locSideMid += m_locMid[0][j];
1134 gloSideMid += m_gloMid[0][j];
1135 }
1136
1137 // scale for correct averaging
1138 locSideMid /= m_numNaturalBF;
1139 gloSideMid /= m_numNaturalBF;
1140 }
1141
1142 // i, j, k, l = number nodes
1143 void compute_side_midpoints(size_t i, size_t j, size_t k, size_t l,
1144 MathVector<dim>& locSideMid,
1145 MathVector<worldDim>& gloSideMid)
1146 {
1147 VecAdd(locSideMid, m_locMid[0][i], m_locMid[0][j], m_locMid[0][k], m_locMid[0][l]);
1148 VecAdd(gloSideMid, m_gloMid[0][i], m_gloMid[0][j], m_gloMid[0][k], m_gloMid[0][l]);
1149
1150 // scale for correct averaging
1151 locSideMid *= 0.25;
1152 gloSideMid *= 0.25;
1153 }
1154
1155 // i, j, k = number nodes
1156 void compute_side_midpoints(size_t i, size_t j, size_t k,
1157 MathVector<dim>& locSideMid, MathVector<worldDim>& gloSideMid)
1158 {
1159 VecAdd(locSideMid, m_locMid[0][i], m_locMid[0][j], m_locMid[0][k]);
1160 VecAdd(gloSideMid, m_gloMid[0][i], m_gloMid[0][j], m_gloMid[0][k]);
1161
1162 // scale for correct averaging
1163 locSideMid *= 1./3.;
1164 gloSideMid *= 1./3.;
1165 }
1166
1167 private:
1169 {
1170 friend class HFV1ManifoldGeometry<TElem, TWorldDim>;
1171
1172 public:
1173 NewEdgeInfo() : m_from(-1), m_to(-1){}
1174 size_t from() const {return m_from;}
1175 size_t to() const {return m_to;}
1176 private:
1177 size_t m_from;
1178 size_t m_to;
1179 };
1180
1182 {
1183 friend class HFV1ManifoldGeometry<TElem, TWorldDim>;
1184 public:
1185 NatEdgeInfo() : nodeId(-1), numChildEdges(0) {for(size_t i=0; i<2; ++i) childEdge[i] = -1;}
1186 bool is_hanging() const {return numChildEdges == 2;}
1187 size_t node_id() const {UG_ASSERT(is_hanging(), "Should only be called, if edge is hanging."); return nodeId;}
1188 size_t num_child_edges() const {return numChildEdges;}
1189 size_t child_edge(size_t i) const {UG_ASSERT(i < num_child_edges(), "Wrong id."); return childEdge[i];}
1190
1191 private:
1192 size_t nodeId;
1194 size_t childEdge[2];
1195 };
1196
1197 // help array: maps NaturalEdge -> NodeOnEdge (-1 if no node on edge)
1198 std::vector<NatEdgeInfo> m_vNatEdgeInfo;
1199 std::vector<NewEdgeInfo> m_vNewEdgeInfo;
1200
1201 private:
1202 // pointer to current element
1204
1205 // local and global geom object midpoints for each dimension
1206 std::vector<MathVector<dim> > m_locMid[dim+1];
1207 std::vector<MathVector<worldDim> > m_gloMid[dim+1];
1208
1209 // BndFaces
1210 std::vector<BF> m_vBF;
1211
1212 // Reference Mapping
1214
1215 // Reference Element
1217};
1218
1219
1220}
1221
1222#endif /* __H__UG__LIB_DISC__SPATIAL_DISC__DISC_HELPER__HANGING_FINITE_VOLUME_GEOMETRY__ */
Definition hfv1_geom.h:543
MidID m_midId[m_numCorners]
Definition hfv1_geom.h:630
const MathVector< worldDim > * global_grad_vector() const
vector of global gradients in ip point
Definition hfv1_geom.h:607
SCVF()
Definition hfv1_geom.h:562
size_t m_from
Definition hfv1_geom.h:623
size_t num_sh() const
number of shape functions
Definition hfv1_geom.h:589
size_t num_corners() const
number of corners, that bound the scvf
Definition hfv1_geom.h:610
MathMatrix< worldDim, dim > JtInv
Definition hfv1_geom.h:641
std::vector< MathVector< worldDim > > globalGrad
Definition hfv1_geom.h:640
const MathVector< worldDim > & normal() const
normal on scvf (points direction "from"->"to"). Norm is equal to area
Definition hfv1_geom.h:580
number detj
Definition hfv1_geom.h:642
size_t from() const
index of SubControlVolume on one side of the scvf
Definition hfv1_geom.h:565
size_t to() const
index of SubControlVolume on one side of the scvf
Definition hfv1_geom.h:568
const MathVector< worldDim > & global_grad(size_t sh) const
value of global gradient of shape function i in integration point
Definition hfv1_geom.h:603
MathVector< worldDim > globalIP
Definition hfv1_geom.h:634
number shape(size_t sh) const
value of shape function i in integration point
Definition hfv1_geom.h:592
static const int dim
dimension of reference element
Definition hfv1_geom.h:546
static const size_t m_numCorners
Number of corners of scvf.
Definition hfv1_geom.h:559
const MathVector< dim > * local_grad_vector() const
vector of local gradients in ip point
Definition hfv1_geom.h:600
MathVector< worldDim > m_vGloPos[m_numCorners]
Definition hfv1_geom.h:629
size_t m_to
Definition hfv1_geom.h:623
const MathVector< dim > & local_grad(size_t sh) const
value of local gradient of shape function i in integration point
Definition hfv1_geom.h:596
MathVector< dim > localIP
Definition hfv1_geom.h:633
size_t num_ip() const
number of integration points on scvf
Definition hfv1_geom.h:571
std::vector< number > vShape
Definition hfv1_geom.h:638
const MathVector< dim > & local_ip() const
local integration point of scvf
Definition hfv1_geom.h:574
std::vector< MathVector< dim > > localGrad
Definition hfv1_geom.h:639
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition hfv1_geom.h:613
const MathMatrix< worldDim, dim > & JTInv() const
Transposed Inverse of Jacobian in integration point.
Definition hfv1_geom.h:583
static const int worldDim
dimension of world
Definition hfv1_geom.h:549
const MathVector< worldDim > & global_corner(size_t co) const
return glbal corner number i
Definition hfv1_geom.h:617
MathVector< worldDim > Normal
Definition hfv1_geom.h:635
static const size_t m_numIP
number of integration points
Definition hfv1_geom.h:556
number detJ() const
Determinante of Jacobian in integration point.
Definition hfv1_geom.h:586
MathVector< dim > m_vLocPos[m_numCorners]
Definition hfv1_geom.h:628
const MathVector< worldDim > & global_ip() const
global integration point of scvf
Definition hfv1_geom.h:577
Definition hfv1_geom.h:646
size_t node_id() const
node id that this scv is associated to
Definition hfv1_geom.h:664
number shape(size_t sh) const
value of shape function i in integration point
Definition hfv1_geom.h:693
size_t num_sh() const
number of shape functions
Definition hfv1_geom.h:690
number volume() const
volume of scv
Definition hfv1_geom.h:676
size_t numSH
Definition hfv1_geom.h:729
MathMatrix< worldDim, dim > JtInv
Definition hfv1_geom.h:733
number detj
Definition hfv1_geom.h:734
SCV()
Definition hfv1_geom.h:661
static const size_t m_maxNumCorners
Number of corners of scvf.
Definition hfv1_geom.h:655
const MathVector< dim > & local_grad(size_t sh) const
value of local gradient of shape function i in integration point
Definition hfv1_geom.h:699
MathVector< dim > m_vLocPos[m_maxNumCorners]
Definition hfv1_geom.h:723
MidID m_midId[m_maxNumCorners]
Definition hfv1_geom.h:725
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition hfv1_geom.h:682
const number * shape_vector() const
vector of shape functions in ip point
Definition hfv1_geom.h:696
number detJ() const
Determinant of Jacobian in integration point.
Definition hfv1_geom.h:716
MathVector< worldDim > globalGrad[m_maxNSH]
Definition hfv1_geom.h:732
const MathVector< dim > & local_ip() const
local integration point of scv
Definition hfv1_geom.h:670
MathVector< worldDim > m_vGloPos[m_maxNumCorners]
Definition hfv1_geom.h:724
size_t m_numCorners
Definition hfv1_geom.h:720
size_t num_corners() const
number of corners, that bound the scvf
Definition hfv1_geom.h:679
MathVector< dim > localGrad[m_maxNSH]
Definition hfv1_geom.h:731
const MathVector< worldDim > & global_ip() const
global integration point
Definition hfv1_geom.h:673
number vShape[m_maxNSH]
Definition hfv1_geom.h:730
const MathVector< worldDim > & global_corner(size_t co) const
return global corner number i
Definition hfv1_geom.h:686
size_t num_ip() const
number of integration points
Definition hfv1_geom.h:667
const MathVector< worldDim > & global_grad(size_t sh) const
value of global gradient of shape function i in integration point
Definition hfv1_geom.h:706
const MathVector< dim > * local_grad_vector() const
vector of local gradients in ip point
Definition hfv1_geom.h:703
traits::scv_type scv_type
type of element the subcontrol volume represents
Definition hfv1_geom.h:658
const MathVector< worldDim > * global_grad_vector() const
vector of global gradients in ip point
Definition hfv1_geom.h:710
number vol
Definition hfv1_geom.h:726
static const size_t m_numIP
number of integration points
Definition hfv1_geom.h:652
const MathMatrix< worldDim, dim > & JTInv() const
Transposed Inverse of Jacobian in integration point.
Definition hfv1_geom.h:713
size_t nodeId
Definition hfv1_geom.h:719
Definition hfv1_geom.h:498
std::vector< MathVector< worldDim > > m_gloMid[dim+1]
Definition hfv1_geom.h:950
const SCV & scv(size_t i) const
const access to SubControlVolume number i
Definition hfv1_geom.h:764
size_t num_sh() const
number of shape functions
Definition hfv1_geom.h:768
size_t num_scvf_ips() const
returns number of all scv ips
Definition hfv1_geom.h:775
static const bool staticLocalData
flag indicating if local data may change
Definition hfv1_geom.h:520
std::vector< SCVF > m_vSCVF
Definition hfv1_geom.h:953
const SCVF & scvf(size_t i) const
const access to SubControlVolumeFace number i
Definition hfv1_geom.h:757
traits::elem_type_3 elem_type_3
Definition hfv1_geom.h:529
void copy_global_corners(SCVF &scvf)
Definition hfv1_geom.h:809
std::vector< MathVector< worldDim > > m_vGlobSCVFIP
Definition hfv1_geom.h:793
void copy_local_corners(SCVF &scvf)
Definition hfv1_geom.h:799
std::vector< MathVector< worldDim > > m_vGlobSCVIP
Definition hfv1_geom.h:795
void update_local_data()
update local data for an element
Definition hfv1_geom.cpp:688
size_t m_numNaturalSCVF
number of SubControlVolumeFaces
Definition hfv1_geom.h:504
void update(GridObject *pElem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
update values for an element
Definition hfv1_geom.cpp:742
size_t get_child_edge_of_corner(size_t i, size_t co)
Definition hfv1_geom.h:896
void copy_local_corners(SCV &scv)
Definition hfv1_geom.h:819
traits::elem_type_0 elem_type_0
element type names
Definition hfv1_geom.h:526
ReferenceObjectID m_roid
current reference object id
Definition hfv1_geom.h:947
static const bool usesHangingNodes
Hanging node flag: this Geometry does support hanging nodes.
Definition hfv1_geom.h:517
DimHFV1Geometry()
constructor
Definition hfv1_geom.h:740
traits::elem_type_4 elem_type_4
Definition hfv1_geom.h:530
std::vector< NatEdgeInfo > m_vNatEdgeInfo
Definition hfv1_geom.h:939
size_t num_scvf() const
number of SubControlVolumeFaces
Definition hfv1_geom.h:754
std::vector< SCV > m_vSCV
Definition hfv1_geom.h:956
size_t num_scv() const
number of SubControlVolumes
Definition hfv1_geom.h:761
size_t m_numSh
Definition hfv1_geom.h:959
const MathVector< dim > * scvf_local_ips() const
returns all ips of scv as they appear in scv loop
Definition hfv1_geom.h:781
std::vector< MathVector< dim > > m_vLocSCVFIP
Definition hfv1_geom.h:794
std::vector< MathVector< dim > > m_vLocSCVIP
Definition hfv1_geom.h:796
void copy_global_corners(SCV &scv)
Definition hfv1_geom.h:832
static const int dim
dimension of reference element
Definition hfv1_geom.h:511
const MathVector< worldDim > * scvf_global_ips() const
returns all ips of scv as they appear in scv loop
Definition hfv1_geom.h:778
static const size_t m_maxNSH
max number of shapes
Definition hfv1_geom.h:507
static const int worldDim
dimension of world
Definition hfv1_geom.h:514
const MathVector< dim > * scv_local_ips() const
returns all ips of scv as they appear in scv loop
Definition hfv1_geom.h:790
traits::elem_type_1 elem_type_1
Definition hfv1_geom.h:527
std::vector< NewEdgeInfo > m_vNewEdgeInfo
Definition hfv1_geom.h:940
traits::elem_type_2 elem_type_2
Definition hfv1_geom.h:528
DimReferenceElement< dim > m_rRefElem
reference element
Definition hfv1_geom.h:962
void compute_side_midpoints(size_t i, size_t j, size_t k, MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:884
size_t num_scv_ips() const
returns number of all scv ips
Definition hfv1_geom.h:784
void compute_side_midpoints(size_t i, MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:846
DimReferenceMapping< dim, worldDim > * m_rMapping
reference mapping
Definition hfv1_geom.h:965
void compute_side_midpoints(size_t i, size_t j, size_t k, size_t l, MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:871
std::vector< MathVector< dim > > m_locMid[dim+1]
Definition hfv1_geom.h:949
size_t m_numNaturalSCV
number of SubControlVolumes
Definition hfv1_geom.h:501
const MathVector< worldDim > * scv_global_ips() const
returns all ips of scv as they appear in scv loop
Definition hfv1_geom.h:787
hdimfv1_traits< TDim > traits
traits
Definition hfv1_geom.h:523
GridObject * m_pElem
pointer to current element
Definition hfv1_geom.h:944
dimension dependent base class for reference elements
Definition reference_element.h:183
virtual base class for reference mappings
Definition reference_mapping_provider.h:53
base class for all FVGeometries
Definition fv_geom_base.h:76
The base class for all geometric objects, such as vertices, edges, faces, volumes,...
Definition grid_base_objects.h:157
Definition hfv1_geom.h:96
const MathVector< worldDim > & global_corner(size_t co) const
return glbal corner number i
Definition hfv1_geom.h:172
MathMatrix< worldDim, dim > JtInv
Definition hfv1_geom.h:196
number detJ() const
Determinante of Jacobian in integration point.
Definition hfv1_geom.h:141
const MathVector< worldDim > & global_grad(size_t sh) const
value of global gradient of shape function i in integration point
Definition hfv1_geom.h:158
const MathVector< worldDim > & normal() const
normal on scvf (points direction "from"->"to"). Norm is equal to area
Definition hfv1_geom.h:135
const MathMatrix< worldDim, dim > & JTInv() const
Transposed Inverse of Jacobian in integration point.
Definition hfv1_geom.h:138
size_t num_ip() const
number of integration points on scvf
Definition hfv1_geom.h:126
static const int worldDim
dimension of world
Definition hfv1_geom.h:102
MathVector< worldDim > m_vGloPos[m_numCorners]
Definition hfv1_geom.h:184
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition hfv1_geom.h:168
std::vector< MathVector< dim > > localGrad
Definition hfv1_geom.h:194
const MathVector< dim > & local_ip() const
local integration point of scvf
Definition hfv1_geom.h:129
MathVector< worldDim > globalIP
Definition hfv1_geom.h:189
size_t from() const
Definition hfv1_geom.h:119
std::vector< MathVector< worldDim > > globalGrad
Definition hfv1_geom.h:195
static const size_t m_numCorners
Number of corners of scvf.
Definition hfv1_geom.h:112
number shape(size_t sh) const
value of shape function i in integration point
Definition hfv1_geom.h:147
size_t m_from
Definition hfv1_geom.h:178
const MathVector< dim > & local_grad(size_t sh) const
value of local gradient of shape function i in integration point
Definition hfv1_geom.h:151
size_t to() const
Definition hfv1_geom.h:123
const MathVector< worldDim > & global_ip() const
global integration point of scvf
Definition hfv1_geom.h:132
size_t m_to
Definition hfv1_geom.h:178
size_t num_corners() const
number of corners, that bound the scvf
Definition hfv1_geom.h:165
size_t num_sh() const
number of shape functions
Definition hfv1_geom.h:144
const MathVector< dim > * local_grad_vector() const
vector of local gradients in ip point
Definition hfv1_geom.h:155
MathVector< dim > m_vLocPos[m_numCorners]
Definition hfv1_geom.h:183
MidID m_midId[m_numCorners]
Definition hfv1_geom.h:185
number detj
Definition hfv1_geom.h:197
MathVector< dim > localIP
Definition hfv1_geom.h:188
static const size_t m_numIP
number of integration points
Definition hfv1_geom.h:109
const MathVector< worldDim > * global_grad_vector() const
vector of global gradients in ip point
Definition hfv1_geom.h:162
static const int dim
dimension of reference element
Definition hfv1_geom.h:99
std::vector< number > vShape
Definition hfv1_geom.h:193
MathVector< worldDim > Normal
Definition hfv1_geom.h:190
SCVF()
Definition hfv1_geom.h:115
Definition hfv1_geom.h:201
static const size_t m_maxNumCorners
Number of corners of scvf.
Definition hfv1_geom.h:210
MathVector< dim > m_vLocPos[m_maxNumCorners]
Definition hfv1_geom.h:249
MidID m_midId[m_maxNumCorners]
Definition hfv1_geom.h:251
const MathVector< worldDim > & global_corner(size_t co) const
return global corner number i
Definition hfv1_geom.h:241
const MathVector< dim > & local_ip() const
local integration point of scv
Definition hfv1_geom.h:225
const MathVector< worldDim > & global_ip() const
global integration point
Definition hfv1_geom.h:228
number vol
Definition hfv1_geom.h:252
size_t m_numCorners
Definition hfv1_geom.h:246
MathVector< worldDim > m_vGloPos[m_maxNumCorners]
Definition hfv1_geom.h:250
size_t num_ip() const
number of integration points
Definition hfv1_geom.h:222
size_t node_id() const
node id that this scv is associated to
Definition hfv1_geom.h:219
size_t num_corners() const
number of corners, that bound the scvf
Definition hfv1_geom.h:234
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition hfv1_geom.h:237
SCV()
Definition hfv1_geom.h:216
static const size_t m_numIP
number of integration points
Definition hfv1_geom.h:207
number volume() const
volume of scv
Definition hfv1_geom.h:231
size_t nodeId
Definition hfv1_geom.h:245
hfv1_traits< ref_elem_type, TWorldDim >::scv_type scv_type
type of element the subcontrol volume represents
Definition hfv1_geom.h:213
Geometry and shape functions for 1st order Vertex-Centered Finite Volume with the hanging nodes.
Definition hfv1_geom.h:61
size_t num_scvf() const
number of SubControlVolumeFaces
Definition hfv1_geom.h:272
ReferenceMapping< ref_elem_type, worldDim > m_rMapping
Definition hfv1_geom.h:491
void copy_local_corners(SCV &scv)
Definition hfv1_geom.h:351
reference_element_traits< TElem >::reference_element_type ref_elem_type
type of reference element
Definition hfv1_geom.h:64
const MathVector< dim > & local_node_position(size_t nodeID) const
return local coords for node ID
Definition hfv1_geom.h:311
void compute_side_midpoints(size_t i, size_t j, size_t k, MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:416
std::vector< MathVector< worldDim > > m_vGlobSCVFIP
Definition hfv1_geom.h:325
void update(GridObject *pElem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
update values for an element
Definition hfv1_geom.cpp:67
const MathVector< worldDim > * scvf_global_ips() const
returns all ips of scv as they appear in scv loop
Definition hfv1_geom.h:296
const ref_elem_type & m_rRefElem
Definition hfv1_geom.h:494
std::vector< NatEdgeInfo > m_vNatEdgeInfo
Definition hfv1_geom.h:471
const MathVector< dim > * scv_local_ips() const
returns all ips of scv as they appear in scv loop
Definition hfv1_geom.h:308
static const int dim
dimension of reference element
Definition hfv1_geom.h:74
std::vector< MathVector< dim > > m_vLocSCVIP
Definition hfv1_geom.h:328
HFV1Geometry()
constructor
Definition hfv1_geom.cpp:40
std::vector< MathVector< dim > > m_vLocSCVFIP
Definition hfv1_geom.h:326
std::vector< NewEdgeInfo > m_vNewEdgeInfo
Definition hfv1_geom.h:472
const MathVector< dim > * scvf_local_ips() const
returns all ips of scv as they appear in scv loop
Definition hfv1_geom.h:299
void compute_side_midpoints(size_t i, MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:378
size_t num_sh() const
number of shape functions
Definition hfv1_geom.h:286
void print()
Definition hfv1_geom.cpp:644
const SCVF & scvf(size_t i) const
const access to SubControlVolumeFace number i
Definition hfv1_geom.h:275
void compute_side_midpoints(size_t i, size_t j, size_t k, size_t l, MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:403
std::vector< MathVector< worldDim > > m_vGlobSCVIP
Definition hfv1_geom.h:327
size_t m_numSh
Definition hfv1_geom.h:488
void copy_global_corners(SCVF &scvf)
Definition hfv1_geom.h:341
std::vector< MathVector< dim > > m_locMid[dim+1]
Definition hfv1_geom.h:478
const SCV & scv(size_t i) const
const access to SubControlVolume number i
Definition hfv1_geom.h:282
size_t num_scv() const
number of SubControlVolumes
Definition hfv1_geom.h:279
size_t num_scvf_ips() const
returns number of all scv ips
Definition hfv1_geom.h:293
static const int worldDim
dimension of world
Definition hfv1_geom.h:77
const MathVector< worldDim > * scv_global_ips() const
returns all ips of scv as they appear in scv loop
Definition hfv1_geom.h:305
void copy_local_corners(SCVF &scvf)
Definition hfv1_geom.h:331
size_t num_scv_ips() const
returns number of all scv ips
Definition hfv1_geom.h:302
static const bool usesHangingNodes
Hanging node flag: this Geometry does support hanging nodes.
Definition hfv1_geom.h:80
GridObject * m_pElem
Definition hfv1_geom.h:476
void copy_global_corners(SCV &scv)
Definition hfv1_geom.h:364
std::vector< SCV > m_vSCV
Definition hfv1_geom.h:485
static const bool staticLocalData
flag indicating if local data may change
Definition hfv1_geom.h:83
std::vector< MathVector< worldDim > > m_gloMid[dim+1]
Definition hfv1_geom.h:479
std::vector< SCVF > m_vSCVF
Definition hfv1_geom.h:482
static const size_t m_numNaturalSCV
number of SubControlVolumes
Definition hfv1_geom.h:67
const MathVector< worldDim > & global_node_position(size_t nodeID) const
return global coords for node ID
Definition hfv1_geom.h:318
GridObject * elem() const
get the element
Definition hfv1_geom.h:268
size_t get_child_edge_of_corner(size_t i, size_t co)
Definition hfv1_geom.h:428
static const size_t m_numNaturalSCVF
number of SubControlVolumeFaces
Definition hfv1_geom.h:70
Definition hfv1_geom.h:1020
size_t nodeId
Definition hfv1_geom.h:1070
MathVector< dim > localIP
Definition hfv1_geom.h:1079
number volume() const
volume of bf
Definition hfv1_geom.h:1049
BF()
Definition hfv1_geom.h:1032
static const size_t numCorners
Definition hfv1_geom.h:1029
const MathVector< worldDim > & global_corner(size_t i) const
return global position of corner number i
Definition hfv1_geom.h:1059
const MathVector< worldDim > & global_ip() const
global integration point of scvf
Definition hfv1_geom.h:1045
size_t node_id() const
node id that this bf is associated to
Definition hfv1_geom.h:1035
const MathVector< dim > & local_corner(size_t i) const
return local position of corner number i
Definition hfv1_geom.h:1055
number shape(size_t i, size_t ip) const
value of shape function i in integration point
Definition hfv1_geom.h:1066
MathVector< worldDim > m_vGloPos[numCorners]
Definition hfv1_geom.h:1076
size_t num_sh() const
number of shape functions
Definition hfv1_geom.h:1063
size_t num_ip() const
number of integration points
Definition hfv1_geom.h:1038
MathVector< worldDim > globalIP
Definition hfv1_geom.h:1080
size_t num_corners() const
number of corners, that bound the bf
Definition hfv1_geom.h:1052
static const size_t m_numIP
Definition hfv1_geom.h:1026
std::vector< number > vShape
Definition hfv1_geom.h:1081
number vol
Definition hfv1_geom.h:1083
MathVector< dim > m_vLocPos[numCorners]
Definition hfv1_geom.h:1075
const MathVector< dim > & local_ip() const
local integration point of scvf
Definition hfv1_geom.h:1041
Definition hfv1_geom.h:976
size_t num_bf_local_ips() const
returns number of all scvf ips
Definition hfv1_geom.h:1121
std::vector< MathVector< dim > > m_locMid[dim+1]
Definition hfv1_geom.h:1206
std::vector< MathVector< dim > > m_vLocBFIP
Definition hfv1_geom.h:1087
size_t num_bf_global_ips() const
returns number of all scvf ips
Definition hfv1_geom.h:1115
static const int dim
Definition hfv1_geom.h:995
static const bool usesHangingNodes
Definition hfv1_geom.h:1004
const BF & bf(size_t i) const
const access to Boundary Face number i
Definition hfv1_geom.h:1108
void print()
print information about hfvg to log
Definition hfv1_geom.cpp:1848
GridObject * m_pElem
Definition hfv1_geom.h:1203
ReferenceMapping< ref_elem_type, worldDim > m_rMapping
Definition hfv1_geom.h:1213
void compute_side_midpoints(size_t i, size_t j, size_t k, MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:1156
const MathVector< dim > * bf_local_ips() const
returns all ips of scvf as they appear in scv loop
Definition hfv1_geom.h:1118
size_t num_bf() const
number of BoundaryFaces
Definition hfv1_geom.h:1105
std::vector< NewEdgeInfo > m_vNewEdgeInfo
Definition hfv1_geom.h:1199
TElem elem_type
Definition hfv1_geom.h:979
hfv1_traits< ref_elem_type, dim >::scv_type bf_type
Definition hfv1_geom.h:1001
static const size_t m_numNaturalBFS
Definition hfv1_geom.h:992
HFV1ManifoldGeometry()
constructor
Definition hfv1_geom.cpp:1249
std::vector< MathVector< worldDim > > m_gloMid[dim+1]
Definition hfv1_geom.h:1207
void compute_side_midpoints(MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:1124
std::vector< MathVector< worldDim > > m_vGlobBFIP
Definition hfv1_geom.h:1088
const MathVector< worldDim > * bf_global_ips() const
returns all ips of scvf as they appear in scv loop
Definition hfv1_geom.h:1112
std::vector< NatEdgeInfo > m_vNatEdgeInfo
Definition hfv1_geom.h:1198
static const int order
Definition hfv1_geom.h:986
static const size_t m_numNaturalBF
Definition hfv1_geom.h:989
std::vector< BF > m_vBF
Definition hfv1_geom.h:1210
reference_element_traits< TElem >::reference_element_type ref_elem_type
Definition hfv1_geom.h:982
static const int worldDim
Definition hfv1_geom.h:998
static const bool staticLocalData
flag indicating if local data may change
Definition hfv1_geom.h:1007
void compute_side_midpoints(size_t i, size_t j, size_t k, size_t l, MathVector< dim > &locSideMid, MathVector< worldDim > &gloSideMid)
Definition hfv1_geom.h:1143
const ref_elem_type & m_rRefElem
Definition hfv1_geom.h:1216
void update(GridObject *elem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
update data for given element
Definition hfv1_geom.cpp:1270
Definition subset_handler_interface.h:223
A class for fixed size, dense matrices.
Definition math_matrix.h:63
a mathematical Vector with N entries.
Definition math_vector.h:97
size_t num(int dim) const
returns the number of geometric objects of dim
Definition reference_element.h:95
int id(int dim_i, size_t i, int dim_j, size_t j) const
id of object j in dimension dim_j of obj i in dimension dim_i
Definition reference_element.h:127
Definition reference_mapping.h:65
#define UG_ASSERT(expr, msg)
Definition assert.h:70
double number
Definition types.h:124
void VecAdd(vector_t &vOut, const vector_t &v1, const vector_t &v2)
adds two MathVector<N>s and stores the result in a third one
Definition math_vector_functions_common_impl.hpp:185
the ug namespace
ReferenceObjectID
these ids are used to identify the shape of a geometric object.
Definition grid_base_objects.h:74
@ ROID_UNKNOWN
Definition grid_base_objects.h:75
Definition hfv1_geom.h:534
size_t dim
Definition hfv1_geom.h:537
size_t id
Definition hfv1_geom.h:538
MidID(size_t dim_, size_t id_)
Definition hfv1_geom.h:536
MidID()
Definition hfv1_geom.h:535
Definition hfv1_geom.h:923
bool is_hanging() const
Definition hfv1_geom.h:927
size_t nodeId
Definition hfv1_geom.h:933
size_t numChildEdges
Definition hfv1_geom.h:934
size_t num_child_edges() const
Definition hfv1_geom.h:929
size_t node_id() const
Definition hfv1_geom.h:928
size_t childEdge[2]
Definition hfv1_geom.h:935
size_t child_edge(size_t i) const
Definition hfv1_geom.h:930
NatEdgeInfo()
Definition hfv1_geom.h:926
Definition hfv1_geom.h:911
size_t m_from
Definition hfv1_geom.h:918
NewEdgeInfo()
Definition hfv1_geom.h:914
size_t from() const
Definition hfv1_geom.h:915
size_t to() const
Definition hfv1_geom.h:916
size_t m_to
Definition hfv1_geom.h:919
Definition hfv1_geom.h:87
MidID(size_t dim_, size_t id_)
Definition hfv1_geom.h:89
MidID()
Definition hfv1_geom.h:88
size_t id
Definition hfv1_geom.h:91
size_t dim
Definition hfv1_geom.h:90
Definition hfv1_geom.h:455
size_t num_child_edges() const
Definition hfv1_geom.h:461
size_t numChildEdges
Definition hfv1_geom.h:466
NatEdgeInfo()
Definition hfv1_geom.h:458
size_t childEdge[2]
Definition hfv1_geom.h:467
bool is_hanging() const
Definition hfv1_geom.h:459
size_t nodeId
Definition hfv1_geom.h:465
size_t child_edge(size_t i) const
Definition hfv1_geom.h:462
size_t node_id() const
Definition hfv1_geom.h:460
Definition hfv1_geom.h:443
NewEdgeInfo()
Definition hfv1_geom.h:446
size_t to() const
Definition hfv1_geom.h:448
size_t m_to
Definition hfv1_geom.h:451
size_t m_from
Definition hfv1_geom.h:450
size_t from() const
Definition hfv1_geom.h:447
Definition hfv1_geom.h:1011
size_t id
Definition hfv1_geom.h:1015
size_t dim
Definition hfv1_geom.h:1014
MidID()
Definition hfv1_geom.h:1012
MidID(size_t dim_, size_t id_)
Definition hfv1_geom.h:1013
Definition hfv1_geom.h:1182
size_t node_id() const
Definition hfv1_geom.h:1187
size_t child_edge(size_t i) const
Definition hfv1_geom.h:1189
size_t childEdge[2]
Definition hfv1_geom.h:1194
bool is_hanging() const
Definition hfv1_geom.h:1186
size_t numChildEdges
Definition hfv1_geom.h:1193
size_t nodeId
Definition hfv1_geom.h:1192
size_t num_child_edges() const
Definition hfv1_geom.h:1188
NatEdgeInfo()
Definition hfv1_geom.h:1185
Definition hfv1_geom.h:1169
size_t to() const
Definition hfv1_geom.h:1175
size_t m_to
Definition hfv1_geom.h:1178
size_t m_from
Definition hfv1_geom.h:1177
size_t from() const
Definition hfv1_geom.h:1174
NewEdgeInfo()
Definition hfv1_geom.h:1173
Traits for Finite Volumes for a generic element of the fixed dimensionalities.
Definition fv_util.h:1325
Definition fv_util.h:1472
Traits for hanging finite volume (dummy implementation)
Definition fv_util.h:1356
void scv_type
Definition fv_util.h:1362
traits for reference elements
Definition reference_element_traits.h:48