Line data Source code
1 : /*
2 : * Copyright (c) 2009-2015: G-CSC, Goethe University Frankfurt
3 : * Author: Sebastian Reiter
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 : #include "common/static_assert.h"
34 : #include "multi_grid.h"
35 :
36 : #ifndef __H__LIB_GRID__MULTI_GRID_IMPL__
37 : #define __H__LIB_GRID__MULTI_GRID_IMPL__
38 :
39 : namespace ug
40 : {
41 : /*
42 : template <class TChild, class TElem>
43 : int MultiGrid::num_children(TElem* elem)
44 : {
45 :
46 : int sharedPipeSec = elem->container_section();
47 : assert(sharedPipeSec != -1 && "bad shared pipe section!");
48 :
49 : return get_elem_info(elem)->m_children.num_elements(sharedPipeSec);
50 : }
51 :
52 : // child access
53 : template <class TChild, class TElem>
54 : int MultiGrid::get_children(std::vector<TChild*>& vChildrenOut, TElem* elem)
55 : {
56 :
57 : int sharedPipeSec = geometry_traits<TChild>::CONTAINER_SECTION;
58 : assert(sharedPipeSec != -1 && "bad shared pipe section!");
59 :
60 : MGElementInfo* elemInfo = get_elem_info(elem);
61 :
62 : typename SectionContainer::iterator iter =
63 : elemInfo->m_children.section_begin(sharedPipeSec);
64 :
65 : typename SectionContainer::iterator iterEnd =
66 : elemInfo->m_children.section_end(sharedPipeSec);
67 :
68 : int numChildren = elemInfo->m_children.num_elements(sharedPipeSec);
69 :
70 : if(vChildrenOut.capacity() < numChildren)
71 : vChildrenOut.reserve(numChildren);
72 :
73 : vChildrenOut.clear();
74 :
75 : for(; iter != iterEnd; ++iter)
76 : {
77 : vChildrenOut.push_back((TChild*) *iter);
78 : }
79 :
80 : return vChildrenOut.size();
81 :
82 : }
83 : */
84 :
85 : inline size_t MultiGrid::
86 : top_level() const
87 : {
88 0 : if(m_hierarchy.num_subsets() <= 0) return 0;
89 0 : return size_t(m_hierarchy.num_subsets() - 1);
90 : }
91 :
92 : template <class TElem>
93 : size_t MultiGrid::
94 : num_children_total(TElem* elem) const
95 : {
96 : size_t numChildren = num_children<TElem>(elem);
97 : size_t numChildrenTotal = numChildren;
98 :
99 : for(size_t i = 0; i < numChildren; ++i)
100 : numChildrenTotal += num_children_total(get_child<TElem>(elem, i));
101 :
102 : return numChildrenTotal;
103 : }
104 :
105 :
106 : template<class TGeomObj>
107 : typename geometry_traits<TGeomObj>::iterator
108 0 : MultiGrid::create(size_t level)
109 : {
110 0 : typename geometry_traits<TGeomObj>::iterator iter =
111 : Grid::create<TGeomObj>();
112 : // put the element into the hierarchy
113 : // (by default it already was assigned to level 0)
114 0 : if(level > 0){
115 0 : level_required(level);
116 0 : m_hierarchy.assign_subset(*iter, level);
117 : }
118 0 : return iter;
119 : }
120 :
121 : template <class TGeomObj>
122 : typename geometry_traits<TGeomObj>::iterator
123 0 : MultiGrid::create(const typename geometry_traits<TGeomObj>::Descriptor& descriptor,
124 : size_t level)
125 : {
126 0 : typename geometry_traits<TGeomObj>::iterator iter =
127 : Grid::create<TGeomObj>(descriptor);
128 : // put the element into the hierarchy
129 : // (by default it already was assigned to level 0)
130 0 : if(level > 0){
131 0 : level_required(level);
132 0 : m_hierarchy.assign_subset(*iter, level);
133 : }
134 0 : return iter;
135 : }
136 :
137 : inline void MultiGrid::level_required(int lvl)
138 : {
139 0 : if(m_hierarchy.num_subsets() <= lvl){
140 0 : create_levels(lvl - m_hierarchy.num_subsets() + 1);
141 : }
142 : }
143 :
144 :
145 : template <class TChild>
146 0 : size_t MultiGrid::num_children(GridObject* elem) const
147 : {
148 0 : switch(elem->base_object_id()){
149 : case VERTEX: return num_children<TChild>(static_cast<Vertex*>(elem));
150 0 : case EDGE: return num_children<TChild>(static_cast<Edge*>(elem));
151 0 : case FACE: return num_children<TChild>(static_cast<Face*>(elem));
152 0 : case VOLUME: return num_children<TChild>(static_cast<Volume*>(elem));
153 : }
154 : return 0;
155 : }
156 :
157 : template <class TChild>
158 0 : TChild* MultiGrid::get_child(GridObject* elem, size_t ind) const
159 : {
160 0 : switch(elem->base_object_id()){
161 0 : case VERTEX: return get_child<TChild>(static_cast<Vertex*>(elem), ind);
162 0 : case EDGE: return get_child<TChild>(static_cast<Edge*>(elem), ind);
163 0 : case FACE: return get_child<TChild>(static_cast<Face*>(elem), ind);
164 0 : case VOLUME: return get_child<TChild>(static_cast<Volume*>(elem), ind);
165 : }
166 : return NULL;
167 : }
168 :
169 : template <class TElem>
170 : void MultiGrid::
171 : clear_child_connections(TElem* parent)
172 : {
173 : if(has_children(parent))
174 : get_info(parent).unregister_from_children(*this);
175 : }
176 :
177 : template <class TElem>
178 0 : void MultiGrid::
179 : associate_parent(TElem* elem, GridObject* parent)
180 : {
181 0 : if(elem->base_object_id() > parent->base_object_id()){
182 0 : UG_THROW("Dimension of parent too low.");
183 : }
184 :
185 : GridObject* oldParent = get_parent(elem);
186 0 : if(oldParent == parent)
187 : return;
188 :
189 0 : if(oldParent)
190 0 : remove_child(oldParent, elem);
191 :
192 : if(parent){
193 0 : add_child(parent, elem);
194 0 : set_parent_type(elem, (char)parent->base_object_id());
195 : }
196 :
197 : set_parent(elem, parent);
198 : }
199 :
200 : template <class TElem>
201 : char MultiGrid::
202 : parent_type(TElem* elem) const
203 : {
204 0 : return m_aaParentType[elem];
205 : }
206 :
207 : template <class TElem>
208 : void MultiGrid::
209 : set_parent_type(TElem* elem, char type)
210 : {
211 0 : m_aaParentType[elem] = type;
212 0 : }
213 :
214 : // info-access
215 : inline MultiGrid::VertexInfo& MultiGrid::get_info(Vertex* v)
216 : {
217 : return m_aaVrtInf[v];
218 : }
219 :
220 : inline MultiGrid::EdgeInfo& MultiGrid::get_info(Edge* e)
221 : {
222 : return m_aaEdgeInf[e];
223 : }
224 :
225 0 : inline MultiGrid::FaceInfo& MultiGrid::get_info(Face* f)
226 : {
227 0 : if(FaceInfo* info = m_aaFaceInf[f])
228 0 : return *info;
229 0 : UG_THROW("MultiGrid::get_info(...): No face info available!");
230 : }
231 :
232 0 : inline MultiGrid::VolumeInfo& MultiGrid::get_info(Volume* v)
233 : {
234 0 : if(VolumeInfo* info = m_aaVolInf[v])
235 0 : return *info;
236 0 : UG_THROW("MultiGrid::get_info(...): No vertex info available!");
237 : }
238 :
239 : // const info-access
240 : inline const MultiGrid::VertexInfo& MultiGrid::get_info(Vertex* v) const
241 : {
242 : return m_aaVrtInf[v];
243 : }
244 :
245 : inline const MultiGrid::EdgeInfo& MultiGrid::get_info(Edge* e) const
246 : {
247 : return m_aaEdgeInf[e];
248 : }
249 :
250 0 : inline const MultiGrid::FaceInfo& MultiGrid::get_info(Face* f) const
251 : {
252 0 : static FaceInfo emptyInfo;
253 0 : if(FaceInfo* info = m_aaFaceInf[f])
254 0 : return *info;
255 : return emptyInfo;
256 : }
257 :
258 0 : inline const MultiGrid::VolumeInfo& MultiGrid::get_info(Volume* v) const
259 : {
260 0 : static VolumeInfo emptyInfo;
261 0 : if(VolumeInfo* info = m_aaVolInf[v])
262 0 : return *info;
263 : return emptyInfo;
264 : }
265 :
266 : template <class TParent, class TChild>
267 0 : void MultiGrid::add_child(TParent* p, TChild* c)
268 : {
269 0 : create_child_info(p);
270 0 : get_info(p).add_child(c);
271 0 : }
272 :
273 : template <class TChild>
274 0 : void MultiGrid::add_child(GridObject* p, TChild* c)
275 : {
276 0 : switch(p->base_object_id()){
277 0 : case VERTEX: add_child(static_cast<Vertex*>(p), c); break;
278 0 : case EDGE: add_child(static_cast<Edge*>(p), c); break;
279 0 : case FACE: add_child(static_cast<Face*>(p), c); break;
280 0 : case VOLUME: add_child(static_cast<Volume*>(p), c); break;
281 : }
282 0 : }
283 :
284 : template <class TParent, class TChild>
285 0 : void MultiGrid::remove_child(TParent* p, TChild* c)
286 : {
287 0 : get_info(p).remove_child(c);
288 0 : }
289 :
290 : template <class TChild>
291 0 : void MultiGrid::remove_child(GridObject* p, TChild* c)
292 : {
293 0 : switch(p->base_object_id()){
294 0 : case VERTEX: remove_child(static_cast<Vertex*>(p), c); break;
295 0 : case EDGE: remove_child(static_cast<Edge*>(p), c); break;
296 0 : case FACE: remove_child(static_cast<Face*>(p), c); break;
297 0 : case VOLUME: remove_child(static_cast<Volume*>(p), c); break;
298 : }
299 0 : }
300 :
301 : template <class TElem, class TParent>
302 0 : void MultiGrid::element_created(TElem* elem, TParent* pParent)
303 : {
304 : // if hierarchical_insertion is enabled, the element will be put
305 : // into the next higher level of pParents level.
306 :
307 : int level = 0;
308 0 : if(pParent)
309 : {
310 : // the element is inserted into a new layer.
311 0 : level = get_level(pParent) + 1;
312 0 : set_parent_type(elem, pParent->base_object_id());
313 : }
314 : else
315 : set_parent_type(elem, -1);
316 :
317 : // register parent and child
318 : //typename mginfo_traits<TElem>::info_type& info = get_info(elem);
319 : //info.m_pParent = pParent;
320 : set_parent(elem, pParent);
321 0 : if(pParent)
322 : {
323 : // make sure that the parent has an info object
324 0 : create_child_info(pParent);
325 :
326 : // add the element to the parents children list
327 0 : typename mginfo_traits<TParent>::info_type& parentInfo = get_info(pParent);
328 0 : parentInfo.add_child(elem);
329 : }
330 :
331 : // put the element into the hierarchy
332 : level_required(level);
333 0 : m_hierarchy.assign_subset(elem, level);
334 0 : }
335 :
336 : template <class TElem, class TParent>
337 0 : void MultiGrid::element_created(TElem* elem, TParent* pParent,
338 : TElem* pReplaceMe)
339 : {
340 : UG_ASSERT(pReplaceMe, "Only call this method with a valid element which shall be replaced.");
341 : int level = get_level(pReplaceMe);
342 :
343 : // register parent and child
344 : set_parent(elem, pParent);
345 :
346 0 : if(pParent)
347 : {
348 : // add the element to the parents children list
349 : // pParent should have an info object at this time!
350 0 : typename mginfo_traits<TParent>::info_type& parentInfo = get_info(pParent);
351 : parentInfo.replace_child(elem, pReplaceMe);
352 : }
353 :
354 : // put the element into the hierarchy
355 : level_required(level);
356 0 : m_hierarchy.assign_subset(elem, level);
357 :
358 : // explicitly copy the parent-type from pReplaceMe to the new vrt.
359 : // This has to be done explicitly since a parent may not exist locally in
360 : // a parallel environment.
361 : set_parent_type(elem, parent_type(pReplaceMe));
362 0 : }
363 :
364 : template <class TElem>
365 0 : void MultiGrid::element_to_be_erased(TElem* elem)
366 : {
367 : // we have to remove the elements children as well.
368 0 : if(has_children(elem)){
369 0 : get_info(elem).unregister_from_children(*this);
370 : }
371 : // we have to remove the associated info object
372 0 : release_child_info(elem);
373 0 : }
374 :
375 : template <class TElem, class TParent>
376 0 : void MultiGrid::element_to_be_erased(TElem* elem, TParent* pParent)
377 : {
378 : // unregister the element from its parent.
379 : // parents always have an info object
380 0 : typename mginfo_traits<TParent>::info_type& parentInfo = get_info(pParent);
381 0 : parentInfo.remove_child(elem);
382 0 : element_to_be_erased(elem);
383 : if(!parentInfo.has_children()){
384 0 : release_child_info(pParent);
385 : }
386 0 : }
387 :
388 : /*
389 : template <class TElem>
390 : void MultiGrid::element_to_be_replaced(TElem* elemOld, TElem* elemNew)
391 : {
392 : }
393 : */
394 :
395 :
396 :
397 : ////////////////////////////////////////////////////////////////////////
398 : ////////////////////////////////////////////////////////////////////////
399 : // specialization of wrapper classes
400 :
401 : ////////////////////////////////////////////////////////////////////////
402 : // specialization for Grid
403 : template <>
404 : class MGWrapper<Grid>
405 : {
406 : public:
407 : MGWrapper(Grid& grid) : m_grid(grid) {}
408 :
409 : inline uint num_levels() const
410 : {return 1;}
411 :
412 : template <class TElem> inline
413 : uint num(int level) const
414 : {return m_grid.num<TElem>();}
415 :
416 : template <class TElem> inline
417 : typename geometry_traits<TElem>::iterator
418 : begin(int level)
419 : {return m_grid.begin<TElem>();}
420 :
421 : template <class TElem> inline
422 : typename geometry_traits<TElem>::iterator
423 : end(int level)
424 : {return m_grid.end<TElem>();}
425 :
426 : protected:
427 : Grid& m_grid;
428 : };
429 :
430 : ////////////////////////////////////////////////////////////////////////
431 : // specialization for MultiGrid
432 : template <>
433 : class MGWrapper<MultiGrid>
434 : {
435 : public:
436 : MGWrapper(MultiGrid& grid) : m_grid(grid) {}
437 :
438 : inline uint num_levels() const
439 : {return (uint)m_grid.num_levels();}
440 :
441 : template <class TElem> inline
442 : uint num(int level) const
443 : {return m_grid.num<TElem>(level);}
444 :
445 : template <class TElem> inline
446 : typename geometry_traits<TElem>::iterator
447 : begin(int level)
448 : {return m_grid.begin<TElem>(level);}
449 :
450 : template <class TElem> inline
451 : typename geometry_traits<TElem>::iterator
452 : end(int level)
453 : {return m_grid.end<TElem>(level);}
454 :
455 : protected:
456 : MultiGrid& m_grid;
457 : };
458 :
459 : }// end of namespace
460 :
461 : #endif
|
