Line data Source code
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__FINITE_VOLUME_OUTPUT__
34 : #define __H__UG__LIB_DISC__SPATIAL_DISC__DISC_HELPER__FINITE_VOLUME_OUTPUT__
35 :
36 : // other ug4 modules
37 : #include "common/common.h"
38 : #include "lib_disc/domain_util.h"
39 :
40 : // finite volume geometry
41 : #include "fv1_geom.h"
42 : #include "hfv1_geom.h"
43 :
44 : namespace ug{
45 :
46 : ////////////////////////////////////////////////////////////////////////////////
47 : // ConstructGridOfSCVF
48 : ////////////////////////////////////////////////////////////////////////////////
49 :
50 : template <typename TElem, template <class, int> class TFVGeom, int TWorldDim>
51 0 : void CreateSCVF(const TElem& elem, TFVGeom<TElem, TWorldDim>& geo, ISubsetHandler& shOut,
52 : Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPosOut)
53 : {
54 : // extract dimensions
55 : static const int refDim = TFVGeom<TElem, TWorldDim>::dim;
56 :
57 : // extract grid
58 0 : Grid& grid = *shOut.grid();
59 :
60 : // tmp vector for vertices
61 : std::vector<Vertex*> vVert;
62 :
63 : // loop all scv of the element
64 0 : for(size_t i = 0; i < geo.num_scvf(); ++i)
65 : {
66 : const typename TFVGeom<TElem, TWorldDim>::SCVF& scvf = geo.scvf(i);
67 :
68 : // clear vertices
69 : vVert.clear();
70 :
71 : // loop corners of scv
72 0 : for(size_t co = 0; co < scvf.num_corners(); ++co)
73 : {
74 : // create a new vertex
75 0 : RegularVertex* vrt = *(grid.create<RegularVertex>());
76 0 : vVert.push_back(vrt);
77 :
78 : // set the coordinates
79 : aaPosOut[vrt] = scvf.global_corner(co);
80 : }
81 :
82 : // edge
83 : if(refDim == 2)
84 : {
85 : if (scvf.num_corners() == 2) {
86 0 : grid.template create<RegularEdge>(EdgeDescriptor(vVert[0], vVert[1]));
87 : }
88 : else {
89 : UG_THROW("SCVF has a number of nodes, that is not drawable.");
90 : }
91 : }
92 : // face
93 : else if(refDim == 3)
94 : {
95 : if (scvf.num_corners() == 4) {
96 0 : grid.template create<Quadrilateral>(QuadrilateralDescriptor(vVert[0], vVert[1],
97 : vVert[2], vVert[3]));
98 : }
99 : else if (scvf.num_corners() == 3) {
100 0 : grid.template create<Triangle>(TriangleDescriptor(vVert[0], vVert[1],
101 : vVert[2]));
102 : }
103 : else{
104 : UG_THROW("SCVF has a number of nodes, that is not drawable.");
105 : }
106 : }
107 : }
108 0 : }
109 :
110 : template <typename TElem, template <class, int> class TFVGeom, int TWorldDim>
111 0 : void ConstructGridOfSCVF(ISubsetHandler& shOut,
112 : const SurfaceView& surfView,
113 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPos,
114 : Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPosOut,
115 : int si)
116 : {
117 : // Create Geometry
118 0 : TFVGeom<TElem, TWorldDim> geo;
119 :
120 : // iterators for primary grid
121 : typename SurfaceView::traits<TElem>::const_iterator iter, iterBegin, iterEnd;
122 0 : iterBegin = surfView.begin<TElem>(GridLevel(GridLevel::TOP, GridLevel::SURFACE, false), SurfaceView::SURFACE);
123 : iterEnd = surfView.end<TElem>(GridLevel(GridLevel::TOP, GridLevel::SURFACE, false), SurfaceView::SURFACE);
124 :
125 : // corners of element
126 : std::vector<MathVector<TWorldDim> > vCornerCoords;
127 :
128 : // iterate over primary grid
129 0 : for(iter = iterBegin; iter != iterEnd; ++iter)
130 : {
131 : // get element
132 : TElem* elem = *iter;
133 :
134 : // get corner coordinates
135 : CollectCornerCoordinates(vCornerCoords, *elem, aaPos);
136 :
137 : // update finite volume geometry
138 0 : geo.update(elem, &vCornerCoords[0], &(*surfView.subset_handler()));
139 :
140 : // Create dual grid
141 0 : CreateSCVF<TElem, TFVGeom, TWorldDim>(*elem, geo, shOut, aaPosOut);
142 : }
143 0 : }
144 :
145 :
146 : template <template <class, int> class TFVGeom, int TWorldDim>
147 : struct ConstructGridOfSCVFWrapper{};
148 :
149 : template <template <class, int> class TFVGeom>
150 : struct ConstructGridOfSCVFWrapper<TFVGeom, 1>
151 : {
152 0 : static void apply(ISubsetHandler& shOut, const SurfaceView& surfView,
153 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<1> > >& aaPos,
154 : Grid::VertexAttachmentAccessor<Attachment<MathVector<1> > >& aaPosOut,
155 : int si, int siDim)
156 : {
157 0 : switch(siDim)
158 : {
159 0 : case 1: ConstructGridOfSCVF<RegularEdge, TFVGeom, 1>(shOut, surfView, aaPos, aaPosOut, si);
160 : break;
161 0 : default: UG_THROW("CreateDualGrid: Dimension " << siDim << " not supported. World dimension is " << 1 <<".");
162 : }
163 0 : }
164 : };
165 :
166 : template <template <class, int> class TFVGeom>
167 : struct ConstructGridOfSCVFWrapper<TFVGeom, 2>
168 : {
169 0 : static void apply(ISubsetHandler& shOut, const SurfaceView& surfView,
170 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<2> > >& aaPos,
171 : Grid::VertexAttachmentAccessor<Attachment<MathVector<2> > >& aaPosOut,
172 : int si, int siDim)
173 : {
174 0 : switch(siDim)
175 : {
176 0 : case 1: ConstructGridOfSCVF<RegularEdge, TFVGeom, 2>(shOut, surfView, aaPos, aaPosOut, si);
177 0 : break;
178 0 : case 2: ConstructGridOfSCVF<Triangle, TFVGeom, 2>(shOut, surfView, aaPos, aaPosOut, si);
179 0 : ConstructGridOfSCVF<Quadrilateral, TFVGeom, 2>(shOut, surfView, aaPos, aaPosOut, si);
180 0 : break;
181 0 : default: UG_THROW("CreateDualGrid: Dimension " << siDim << " not supported. World dimension is " << 2);
182 : }
183 0 : }
184 : };
185 :
186 : template <template <class, int> class TFVGeom>
187 : struct ConstructGridOfSCVFWrapper<TFVGeom, 3>
188 : {
189 0 : static void apply(ISubsetHandler& shOut, const SurfaceView& surfView,
190 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<3> > >& aaPos,
191 : Grid::VertexAttachmentAccessor<Attachment<MathVector<3> > >& aaPosOut,
192 : int si, int siDim)
193 : {
194 0 : switch(siDim)
195 : {
196 0 : case 1: ConstructGridOfSCVF<RegularEdge, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
197 0 : break;
198 0 : case 2: ConstructGridOfSCVF<Triangle, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
199 0 : ConstructGridOfSCVF<Quadrilateral, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
200 0 : break;
201 0 : case 3: ConstructGridOfSCVF<Tetrahedron, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
202 0 : ConstructGridOfSCVF<Hexahedron, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
203 0 : ConstructGridOfSCVF<Prism, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
204 0 : ConstructGridOfSCVF<Pyramid, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
205 0 : ConstructGridOfSCVF<Octahedron, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
206 0 : break;
207 0 : default: UG_THROW("CreateDualGrid: Dimension " << siDim << " not supported. World dimension is " << 3);
208 : }
209 0 : }
210 : };
211 :
212 : template <template <class, int> class TFVGeom, int TWorldDim>
213 : void ConstructGridOfSCVF(ISubsetHandler& shOut, const SurfaceView& surfView,
214 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPos,
215 : Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPosOut,
216 : int si, int siDim)
217 : {
218 0 : ConstructGridOfSCVFWrapper<TFVGeom, TWorldDim>::apply(shOut, surfView, aaPos, aaPosOut, si, siDim);
219 0 : }
220 :
221 :
222 : ////////////////////////////////////////////////////////////////////////////////
223 : // ConstructGridOfSCV
224 : ////////////////////////////////////////////////////////////////////////////////
225 :
226 : template <typename TElem, template <class, int> class TFVGeom, int TWorldDim>
227 0 : void CreateSCV(const TElem& elem, TFVGeom<TElem, TWorldDim>& geo, ISubsetHandler& shOut,
228 : Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPosOut)
229 : {
230 : // extract dimensions
231 : static const int refDim = TFVGeom<TElem, TWorldDim>::dim;
232 :
233 : // extract grid
234 0 : Grid& grid = *shOut.grid();
235 :
236 : // tmp vector for vertices
237 : std::vector<Vertex*> vVert;
238 :
239 : // loop all scv of the element
240 0 : for(size_t i = 0; i < geo.num_scv(); ++i)
241 : {
242 : const typename TFVGeom<TElem, TWorldDim>::SCV& scv = geo.scv(i);
243 :
244 : // clear vertices
245 : vVert.clear();
246 :
247 : // loop corners of scv
248 0 : for(size_t co = 0; co < scv.num_corners(); ++co)
249 : {
250 : // create a new vertex
251 0 : RegularVertex* vrt = *(grid.create<RegularVertex>());
252 0 : vVert.push_back(vrt);
253 :
254 : // set the coordinates
255 : aaPosOut[vrt] = scv.global_corner(co);
256 :
257 : // if co == 0, it is a vertex of the primary grid
258 : // We assign all of those vertices to subset 0
259 : // The other vertices remain in subset -1
260 0 : if(co == 0) shOut.assign_subset(vrt, 0);
261 0 : else shOut.assign_subset(vrt, -1);
262 :
263 : }
264 :
265 : // edge
266 : if(refDim == 1)
267 : {
268 0 : grid.template create<RegularEdge>(EdgeDescriptor(vVert[0], vVert[1]));
269 : }
270 : // face
271 : else if(refDim == 2)
272 : {
273 0 : if (scv.num_corners() == 4) {
274 0 : grid.template create<Quadrilateral>(QuadrilateralDescriptor(vVert[0], vVert[1],
275 : vVert[2], vVert[3]));
276 : }
277 : else {
278 0 : UG_THROW("SCV has a number of nodes, that is not drawable.");
279 : }
280 : }
281 : // volume
282 : else if(refDim == 3)
283 : {
284 0 : if (scv.num_corners() == 8) {
285 0 : grid.template create<Hexahedron>(HexahedronDescriptor( vVert[0], vVert[1], vVert[2], vVert[3],
286 : vVert[4], vVert[5], vVert[6], vVert[7]));
287 : }
288 0 : else if (scv.num_corners() == 4) {
289 0 : grid.template create<Tetrahedron>(TetrahedronDescriptor( vVert[0], vVert[1], vVert[2], vVert[3]));
290 : }
291 : else {
292 0 : UG_THROW("SCV has a number of nodes, that is not drawable.");
293 : }
294 : }
295 : }
296 0 : }
297 :
298 : template <typename TElem, template <class, int> class TFVGeom, int TWorldDim>
299 0 : void ConstructGridOfSCV(ISubsetHandler& shOut, const SurfaceView& surfView,
300 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPos,
301 : Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPosOut,
302 : int si)
303 : {
304 : // Create Geometry
305 0 : TFVGeom<TElem, TWorldDim> geo;
306 :
307 : // iterators for primary grid
308 : typename SurfaceView::traits<TElem>::const_iterator iter, iterBegin, iterEnd;
309 0 : iterBegin = surfView.begin<TElem>(GridLevel(GridLevel::TOP, GridLevel::SURFACE, false), SurfaceView::SURFACE);
310 0 : iterEnd = surfView.end<TElem>(GridLevel(GridLevel::TOP, GridLevel::SURFACE, false), SurfaceView::SURFACE);
311 :
312 : // corners of element
313 : std::vector<MathVector<TWorldDim> > vCornerCoords;
314 :
315 : // iterate over primary grid
316 0 : for(iter = iterBegin; iter != iterEnd; ++iter)
317 : {
318 : // get element
319 : TElem* elem = *iter;
320 :
321 : // get corner coordinates
322 : CollectCornerCoordinates(vCornerCoords, *elem, aaPos);
323 :
324 : // update finite volume geometry
325 0 : geo.update(elem, &vCornerCoords[0], &(*surfView.subset_handler()));
326 :
327 : // Create dual grid
328 0 : CreateSCV<TElem, TFVGeom, TWorldDim>(*elem, geo, shOut, aaPosOut);
329 : }
330 0 : }
331 :
332 :
333 : template <template <class, int> class TFVGeom, int TWorldDim>
334 : struct ConstructGridOfSCVWrapper{};
335 :
336 : template <template <class TElem, int TWorldDim> class TFVGeom>
337 : struct ConstructGridOfSCVWrapper<TFVGeom, 1>
338 : {
339 0 : static void apply(ISubsetHandler& shOut, const SurfaceView& surfView,
340 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<1> > >& aaPos,
341 : Grid::VertexAttachmentAccessor<Attachment<MathVector<1> > >& aaPosOut,
342 : int si, int siDim)
343 : {
344 0 : switch(siDim)
345 : {
346 0 : case 1: ConstructGridOfSCV<RegularEdge, TFVGeom, 1>(shOut, surfView, aaPos, aaPosOut, si);
347 : break;
348 0 : default: UG_THROW("CreateDualGrid: Dimension " << siDim << " not supported. World dimension is " << 1);
349 : }
350 0 : }
351 : };
352 :
353 : template <template <class, int> class TFVGeom>
354 : struct ConstructGridOfSCVWrapper<TFVGeom, 2>
355 : {
356 0 : static void apply(ISubsetHandler& shOut, const SurfaceView& surfView,
357 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<2> > >& aaPos,
358 : Grid::VertexAttachmentAccessor<Attachment<MathVector<2> > >& aaPosOut,
359 : int si, int siDim)
360 : {
361 0 : switch(siDim)
362 : {
363 0 : case 1: ConstructGridOfSCV<RegularEdge, TFVGeom, 2>(shOut, surfView, aaPos, aaPosOut, si);
364 0 : break;
365 0 : case 2: ConstructGridOfSCV<Triangle, TFVGeom, 2>(shOut, surfView, aaPos, aaPosOut, si);
366 0 : ConstructGridOfSCV<Quadrilateral, TFVGeom, 2>(shOut, surfView, aaPos, aaPosOut, si);
367 0 : break;
368 0 : default: UG_THROW("CreateDualGrid: Dimension " << siDim << " not supported. World dimension is " << 2);
369 : }
370 0 : }
371 : };
372 :
373 : template <template <class, int> class TFVGeom>
374 : struct ConstructGridOfSCVWrapper<TFVGeom, 3>
375 : {
376 0 : static void apply(ISubsetHandler& shOut, const SurfaceView& surfView,
377 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<3> > >& aaPos,
378 : Grid::VertexAttachmentAccessor<Attachment<MathVector<3> > >& aaPosOut,
379 : int si, int siDim)
380 : {
381 0 : switch(siDim)
382 : {
383 0 : case 1: ConstructGridOfSCV<RegularEdge, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
384 0 : break;
385 0 : case 2: ConstructGridOfSCV<Triangle, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
386 0 : ConstructGridOfSCV<Quadrilateral, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
387 0 : break;
388 0 : case 3: ConstructGridOfSCV<Tetrahedron, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
389 0 : ConstructGridOfSCV<Hexahedron, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
390 0 : ConstructGridOfSCV<Prism, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
391 0 : ConstructGridOfSCV<Pyramid, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
392 0 : ConstructGridOfSCV<Octahedron, TFVGeom, 3>(shOut, surfView, aaPos, aaPosOut, si);
393 0 : break;
394 0 : default: UG_THROW("CreateDualGrid: Dimension " << siDim << " not supported. World dimension is " << 3);
395 : }
396 0 : }
397 : };
398 :
399 : template <template <class, int> class TFVGeom, int TWorldDim>
400 : void ConstructGridOfSCV(ISubsetHandler& shOut, const SurfaceView& surfView,
401 : const Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPos,
402 : Grid::VertexAttachmentAccessor<Attachment<MathVector<TWorldDim> > >& aaPosOut,
403 : int si, int siDim)
404 : {
405 0 : ConstructGridOfSCVWrapper<TFVGeom, TWorldDim>::apply(shOut, surfView, aaPos, aaPosOut, si, siDim);
406 0 : }
407 :
408 :
409 : ////////////////////////////////////////////////////////////////////////////////
410 : // Assignement of Subsets
411 : ////////////////////////////////////////////////////////////////////////////////
412 :
413 : template <typename TElem>
414 0 : void ColorSubControlVolumeFaces(ISubsetHandler& shOut)
415 : {
416 : // extract grid
417 0 : Grid& grid = *shOut.grid();
418 :
419 : // iterators for primary grid
420 : typename geometry_traits<TElem>::iterator iter, iterBegin, iterEnd;
421 : iterBegin = grid.begin<TElem>();
422 : iterEnd = grid.end<TElem>();
423 :
424 : // iterate over primary grid
425 : int si = 0;
426 0 : for(iter = iterBegin; iter != iterEnd; ++iter)
427 : {
428 0 : shOut.assign_subset(*iter, si++);
429 : }
430 0 : }
431 :
432 : template <typename TElem>
433 0 : void ColorSubControlVolume(ISubsetHandler& shOut)
434 : {
435 : // extract grid
436 0 : Grid& grid = *shOut.grid();
437 :
438 : // iterators for primary grid
439 : typename geometry_traits<TElem>::iterator iter, iterBegin, iterEnd;
440 : iterBegin = grid.begin<TElem>();
441 : iterEnd = grid.end<TElem>();
442 :
443 : // iterate over primary grid
444 : int si = 0;
445 0 : for(iter = iterBegin; iter != iterEnd; ++iter)
446 : {
447 0 : shOut.assign_subset(*iter, si++);
448 : }
449 0 : }
450 :
451 :
452 : template <int TRefDim>
453 0 : void ColorControlVolume(ISubsetHandler& shOut)
454 : {
455 : // extract grid
456 0 : Grid& grid = *shOut.grid();
457 :
458 : std::vector<Volume*> vVols;
459 : std::vector<Face*> vFaces;
460 : std::vector<Edge*> vEdges;
461 :
462 : int si = 0;
463 0 : for(VertexIterator iter = shOut.grid()->begin<Vertex>();
464 0 : iter != shOut.grid()->end<Vertex>(); ++iter, ++si)
465 : {
466 0 : if(shOut.get_subset_index(*iter) != 0) continue;
467 :
468 : switch(TRefDim)
469 : {
470 0 : case 1: CollectEdges(vEdges, grid, *iter);
471 : shOut.assign_subset(vEdges.begin(), vEdges.end(), si);
472 : break;
473 0 : case 2: CollectFaces(vFaces, grid, *iter);
474 : shOut.assign_subset(vFaces.begin(), vFaces.end(), si);
475 : break;
476 0 : case 3: CollectVolumes(vVols, grid, *iter);
477 : shOut.assign_subset(vVols.begin(), vVols.end(), si);
478 : break;
479 : default: UG_THROW("Dimension " << TRefDim << " is not supported.");
480 : }
481 : }
482 0 : }
483 :
484 :
485 : template <template <class, int> class TFVGeom, typename TAAPosition>
486 0 : void CreateGridOfSubControlVolumes(ISubsetHandler& shOut, TAAPosition& aaPosOut, const ISubsetHandler& sh, const TAAPosition& aaPos, const SurfaceView& surfView, int si = -1)
487 : {
488 : static const int dim = TAAPosition::ValueType::Size;
489 :
490 : // Construct dual domain for scv and given subset
491 0 : if(si >= 0)
492 : {
493 0 : const int siDim = DimensionOfSubset(sh, si);
494 : ConstructGridOfSCV<TFVGeom, dim>(shOut, surfView, aaPos, aaPosOut, si, siDim);
495 : }
496 : // if no subset selected, construct dual grid for all subsets with dim == worldDim
497 : else
498 : {
499 0 : for(si = 0; si < sh.num_subsets(); ++si)
500 : {
501 0 : const int siDim = DimensionOfSubset(sh, si);
502 0 : if(siDim != dim) continue;
503 :
504 : ConstructGridOfSCV<TFVGeom, dim>(shOut, surfView, aaPos, aaPosOut, si, siDim);
505 : }
506 : }
507 :
508 : // Let each SubControlVolume be one subset
509 : switch(dim)
510 : {
511 0 : case 1: ColorSubControlVolume<Edge>(shOut); break;
512 0 : case 2: ColorSubControlVolume<Face>(shOut); break;
513 0 : case 3: ColorSubControlVolume<Volume>(shOut); break;
514 : default: UG_THROW("WriteDualGridToFile: Dimension "<<dim<<" not supported.");
515 : }
516 0 : }
517 :
518 : template <template <class, int> class TFVGeom, typename TAAPosition, typename TAPosition>
519 0 : void CreateGridOfControlVolumes(ISubsetHandler& shOut, TAAPosition& aaPosOut, TAPosition& aPosOut, const ISubsetHandler& sh, const TAAPosition& aaPos, const SurfaceView& surfView, int si = -1)
520 : {
521 : static const int dim = TAAPosition::ValueType::Size;
522 :
523 : // Construct dual domain for scv and given subset
524 0 : if(si >= 0)
525 : {
526 0 : const int siDim = DimensionOfSubset(sh, si);
527 : ConstructGridOfSCV<TFVGeom, dim>(shOut, surfView, aaPos, aaPosOut, si, siDim);
528 : }
529 : // if no subset selected, construct dual grid for all subsets with dim == worldDim
530 : else
531 : {
532 0 : for(si = 0; si < sh.num_subsets(); ++si)
533 : {
534 0 : const int siDim = DimensionOfSubset(sh, si);
535 0 : if(siDim != dim) continue;
536 :
537 : ConstructGridOfSCV<TFVGeom, dim>(shOut, surfView, aaPos, aaPosOut, si, siDim);
538 : }
539 : }
540 :
541 : // remove doubles
542 0 : RemoveDoubles<dim>(*shOut.grid(),
543 : shOut.grid()->begin<Vertex>(), shOut.grid()->end<Vertex>(),
544 : aPosOut, 1e-5);
545 :
546 : // Let each SubControlVolume be one subset
547 : switch(dim)
548 : {
549 0 : case 1: ColorControlVolume<1>(shOut); break;
550 0 : case 2: ColorControlVolume<2>(shOut); break;
551 0 : case 3: ColorControlVolume<3>(shOut); break;
552 : default: UG_THROW("WriteDualGridToFile: Dimension "<<dim<<" not supported.");
553 : }
554 0 : }
555 :
556 : template <template <class, int> class TFVGeom, typename TAAPosition>
557 0 : void CreateGridOfSubControlVolumeFaces(ISubsetHandler& shOut, TAAPosition& aaPosOut, const ISubsetHandler& sh, const TAAPosition& aaPos, const SurfaceView& surfView, int si = -1)
558 : {
559 : static const int dim = TAAPosition::ValueType::Size;
560 :
561 : // Construct dual domain for scv and given subset
562 0 : if(si >= 0)
563 : {
564 0 : const int siDim = DimensionOfSubset(sh, si);
565 :
566 : ConstructGridOfSCVF<TFVGeom, dim>(shOut, surfView, aaPos, aaPosOut, si, siDim);
567 : }
568 : // if no subset selected, construct dual grid for all subsets with dim == worldDim
569 : else
570 : {
571 0 : for(si = 0; si < sh.num_subsets(); ++si)
572 : {
573 0 : const int siDim = DimensionOfSubset(sh, si);
574 0 : if(siDim != dim) continue;
575 :
576 : ConstructGridOfSCVF<TFVGeom, dim>(shOut, surfView, aaPos, aaPosOut, si, siDim);
577 : }
578 : }
579 :
580 : // Let each SubControlVolume be one subset
581 : switch(dim)
582 : {
583 0 : case 1: ColorSubControlVolumeFaces<Vertex>(shOut); break;
584 0 : case 2: ColorSubControlVolumeFaces<Edge>(shOut); break;
585 0 : case 3: ColorSubControlVolumeFaces<Face>(shOut); break;
586 : default: UG_THROW("WriteDualGridToFile: Dimension " << dim << " not supported.");
587 : }
588 0 : }
589 :
590 : /**
591 : * Creates a grid consisting of the sub-control-volume faces w.r.t to a given
592 : * domain.
593 : *
594 : * @param domOut domain to be filled
595 : * @param domIn original domain
596 : * @param si subset used (-1 for whole domain)
597 : */
598 : template <template <class, int> class TFVGeom, typename TDomain>
599 0 : void CreateSubControlVolumeFaceDomain(TDomain& domOut, const TDomain& domIn, const SurfaceView& surfView, int si = -1)
600 : {
601 0 : if(&domOut == &domIn)
602 0 : UG_THROW("CreateSubControlVolumeFaceDomain: Domains must be different.");
603 :
604 : CreateGridOfSubControlVolumeFaces<TFVGeom, typename TDomain::position_accessor_type>
605 0 : ( *domOut.subset_handler(),
606 : domOut.position_accessor(),
607 0 : *domIn.subset_handler(),
608 : domIn.position_accessor(),
609 : surfView, si);
610 0 : }
611 :
612 : /**
613 : * Creates a grid consisting of the sub-control-volume w.r.t to a given
614 : * domain.
615 : *
616 : * @param domOut domain to be filled
617 : * @param domIn original domain
618 : * @param si subset used (-1 for whole domain)
619 : */
620 : template <template <class, int> class TFVGeom, typename TDomain>
621 0 : void CreateSubControlVolumeDomain(TDomain& domOut, const TDomain& domIn, const SurfaceView& surfView, int si = -1)
622 : {
623 0 : if(&domOut == &domIn)
624 0 : UG_THROW("CreateSubControlVolumeDomain: Domains must be different.");
625 :
626 : CreateGridOfSubControlVolumes<TFVGeom, typename TDomain::position_accessor_type>
627 0 : ( *domOut.subset_handler(),
628 : domOut.position_accessor(),
629 0 : *domIn.subset_handler(),
630 : domIn.position_accessor(),
631 : surfView, si);
632 0 : }
633 :
634 : /**
635 : * Creates a grid consisting of the sub-control-volume faces w.r.t to a given
636 : * domain.
637 : *
638 : * @param domOut domain to be filled
639 : * @param domIn original domain
640 : * @param si subset used (-1 for whole domain)
641 : */
642 : template <template <class, int> class TFVGeom, typename TDomain>
643 0 : void CreateControlVolumeDomain(TDomain& domOut, const TDomain& domIn, const SurfaceView& surfView, int si = -1)
644 : {
645 0 : if(&domOut == &domIn)
646 0 : UG_THROW("CreateControlVolumeDomain: Domains must be different.");
647 :
648 : CreateGridOfControlVolumes<TFVGeom, typename TDomain::position_accessor_type, typename TDomain::position_attachment_type>
649 0 : ( *domOut.subset_handler(),
650 : domOut.position_accessor(),
651 : domOut.position_attachment(),
652 0 : *domIn.subset_handler(),
653 : domIn.position_accessor(),
654 : surfView, si);
655 0 : }
656 :
657 : } // end namespace ug
658 :
659 :
660 : #endif /* __H__UG__LIB_DISC__SPATIAL_DISC__DISC_HELPER__FINITE_VOLUME_OUTPUT__ */
|
