Line data Source code
1 : /*
2 : * Copyright (c) 2015: G-CSC, Goethe University Frankfurt
3 : * Author: Michael Lampe
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 : /*
34 : * Singular sources and sinks for the FV discretizations:
35 : */
36 :
37 : #ifndef __H__UG__LIB_DISC__SPATIAL_DISC__FV_SINGULAR_SOURCES_AND_SINKS__
38 : #define __H__UG__LIB_DISC__SPATIAL_DISC__FV_SINGULAR_SOURCES_AND_SINKS__
39 :
40 : #include <vector>
41 :
42 : // ug4 headers
43 : #include "lib_grid/algorithms/bounding_box_util.h"
44 : #include "lib_grid/algorithms/geom_obj_util/geom_obj_util.h"
45 :
46 : namespace ug {
47 :
48 : /// Base class for the point sources and sinks
49 : template <int dim, typename TData>
50 0 : class FVPointSourceOrSink : public TData
51 : {
52 : private:
53 : MathVector<dim> point; ///< the point of the source/sink
54 :
55 : public:
56 :
57 : /// class constructor (with the default constructor for data)
58 : FVPointSourceOrSink
59 : (
60 : const MathVector<dim>& _point ///< coordinates of the source/sink
61 : )
62 : : point(_point) {}
63 :
64 : /// class constructor (with the copy constructor for data)
65 : FVPointSourceOrSink
66 : (
67 : const MathVector<dim>& _point, ///< coordinates of the source/sink
68 : const TData& _data ///< the data to copy
69 : )
70 : : TData(_data), point(_point) {}
71 :
72 : /// class constructor (with the default constructor for data)
73 0 : FVPointSourceOrSink
74 : (
75 : const std::vector<number>& _point ///< coordinates of the source/sink
76 : )
77 : {
78 0 : if (_point.size () != dim)
79 0 : UG_THROW ("FVPointSourceOrSink: Number of coordinates does not match the dimension.");
80 0 : for (size_t i = 0; i < dim; i++) point[i] = _point[i];
81 0 : }
82 :
83 : /// returns the point of the source/sink
84 0 : const MathVector<dim>& position () {return point;}
85 :
86 : /// test whether a source/sink point corresponds to a given corner of the element
87 : template <typename TElem, typename TAAPos, typename TFVGeom>
88 : bool corresponds_to
89 : (
90 : TElem* elem, ///< the element
91 : Grid& grid, ///< the grid
92 : TAAPos& aaPos, ///< position of the vertices
93 : const TFVGeom& geo, ///< FV geometry (initialized for 'elem')
94 : size_t co ///< corner to get the contribution for
95 : );
96 : };
97 :
98 : /// Base class for line sources and sinks
99 : template<int dim, typename TData>
100 0 : class FVLineSourceOrSink : public TData
101 : {
102 : private:
103 : MathVector<dim> point1; ///< beginning of the line segment
104 : MathVector<dim> point2; ///< end of the line segment
105 :
106 : public:
107 :
108 : /// class constructor (with the default constructor for data)
109 : FVLineSourceOrSink
110 : (
111 : const MathVector<dim>& _point1, ///< beginning of the line segment
112 : const MathVector<dim>& _point2 ///< end of the line segment
113 : )
114 : : point1(_point1), point2(_point2) {}
115 :
116 : /// class constructor (with the copy constructor for data)
117 : FVLineSourceOrSink
118 : (
119 : const MathVector<dim>& _point1, ///< beginning of the line segment
120 : const MathVector<dim>& _point2, ///< end of the line segment
121 : const TData& _data ///< the data to copy
122 : )
123 : : TData(_data), point1(_point1), point2(_point2) {}
124 :
125 : /// class constructor (with the default constructor for data)
126 0 : FVLineSourceOrSink
127 : (
128 : const std::vector<number>& _point1, ///< beginning of the line segment
129 : const std::vector<number>& _point2 ///< end of the line segment
130 : )
131 : {
132 0 : if (_point1.size () != dim || _point2.size () != dim)
133 0 : UG_THROW ("FVLineSourceOrSink: Number of coordinates does not match the dimension.");
134 0 : for (size_t i = 0; i < dim; i++)
135 : {
136 0 : point1[i] = _point1[i]; point2[i] = _point2[i];
137 : }
138 0 : }
139 :
140 : /// returns the beginning of the line segment
141 0 : const MathVector<dim>& from_position () {return point1;}
142 :
143 : /// returns the end of the line segment
144 0 : const MathVector<dim>& to_position () {return point2;}
145 :
146 : /// test whether a source/sink point corresponds to a given corner of the element
147 : template <typename TElem, typename TAAPos, typename TFVGeom>
148 : bool corresponds_to
149 : (
150 : TElem* elem, ///< [in] the element
151 : Grid& grid, ///< [in] the grid
152 : TAAPos& aaPos, ///< [in] position of the vertices
153 : const TFVGeom& geo, ///< [in] FV geometry (initialized for 'elem')
154 : size_t co, ///< [in] corner to get the contribution for
155 : MathVector<dim>& ls, ///< [out] beginning of the subsegment
156 : MathVector<dim>& le ///< [out] end of the subsegment
157 : );
158 : };
159 :
160 : /// Partial specialization of for 1d: No line sources and sinks in 1d
161 : /// Base class for line sources and sinks
162 : template<typename TData>
163 0 : class FVLineSourceOrSink<1, TData> : public TData
164 : {
165 : private:
166 : static const int dim = 1;
167 : MathVector<dim> point1; ///< beginning of the line segment
168 : MathVector<dim> point2; ///< end of the line segment
169 :
170 : public:
171 :
172 : /// class constructor (with the default constructor for data)
173 : FVLineSourceOrSink
174 : (
175 : const MathVector<dim>& _point1, ///< beginning of the line segment
176 : const MathVector<dim>& _point2 ///< end of the line segment
177 : )
178 : : point1(_point1), point2(_point2) {}
179 :
180 : /// class constructor (with the copy constructor for data)
181 : FVLineSourceOrSink
182 : (
183 : const MathVector<dim>& _point1, ///< beginning of the line segment
184 : const MathVector<dim>& _point2, ///< end of the line segment
185 : const TData& _data ///< the data to copy
186 : )
187 : : TData(_data), point1(_point1), point2(_point2) {}
188 :
189 : /// class constructor (with the default constructor for data)
190 0 : FVLineSourceOrSink
191 : (
192 : const std::vector<number>& _point1, ///< beginning of the line segment
193 : const std::vector<number>& _point2 ///< end of the line segment
194 : )
195 : {
196 0 : UG_THROW ("FVLineSourceOrSink: Line sources and sinks are not supported in 1d.");
197 : }
198 :
199 : /// returns the beginning of the line segment
200 0 : const MathVector<dim>& from_position () {return point1;}
201 :
202 : /// returns the end of the line segment
203 0 : const MathVector<dim>& to_position () {return point2;}
204 :
205 : /// test whether a source/sink point corresponds to a given corner of the element
206 : template <typename TElem, typename TAAPos, typename TFVGeom>
207 0 : bool corresponds_to
208 : (
209 : TElem* elem, ///< [in] the element
210 : Grid& grid, ///< [in] the grid
211 : TAAPos& aaPos, ///< [in] position of the vertices
212 : const TFVGeom& geo, ///< [in] FV geometry (initialized for 'elem')
213 : size_t co, ///< [in] corner to get the contribution for
214 : MathVector<dim>& ls, ///< [out] beginning of the subsegment
215 : MathVector<dim>& le ///< [out] end of the subsegment
216 : )
217 : {
218 0 : UG_THROW ("FVLineSourceOrSink: Line sources and sinks are not supported in 1d.");
219 : return false;
220 : }
221 : };
222 :
223 : /// Manager class for point and line sources and sinks
224 : template <int dim, typename TPointData, typename TLineData = TPointData>
225 0 : class FVSingularSourcesAndSinks
226 : {
227 : public:
228 : typedef FVPointSourceOrSink<dim, TPointData> point_sss_type;
229 : typedef FVLineSourceOrSink<dim, TLineData> line_sss_type;
230 :
231 : private:
232 : typename std::vector<SmartPtr<point_sss_type> > ListP;
233 : typename std::vector<SmartPtr<line_sss_type> > ListL;
234 :
235 : public:
236 :
237 : /// class constructor
238 0 : FVSingularSourcesAndSinks () {}
239 :
240 : /// add a point source or sink
241 0 : void add_point
242 : (
243 : SmartPtr<point_sss_type> point_sss ///< the object to add
244 : )
245 : {
246 0 : ListP.push_back (point_sss);
247 0 : }
248 :
249 : /// add a line source or sink
250 0 : void add_line
251 : (
252 : SmartPtr<line_sss_type> line_sss ///< the object to add
253 : )
254 : {
255 0 : ListL.push_back (line_sss);
256 0 : }
257 :
258 : /// returns the number of the point sources and sinks
259 : size_t num_points () {return ListP.size ();}
260 :
261 : /// return a point sink with a given index
262 : point_sss_type * point (size_t i) {return ListP[i].get();}
263 :
264 : /// calls the init-functions of all point sinks (if there are any in TPointData)
265 : void init_all_point_sss ()
266 : {
267 0 : for (size_t i = 0; i < ListP.size (); i++)
268 : ((TPointData *) ListP[i].get())->init ();
269 : }
270 :
271 : /// class of point source and sink iterator for a given scv
272 : template <typename TElem, typename TAAPos, typename TFVGeom>
273 : class point_iterator
274 : {
275 : typedef point_iterator<TElem, TAAPos, TFVGeom> this_type;
276 : typedef FVSingularSourcesAndSinks<dim, TPointData, TLineData> master_type;
277 :
278 : master_type * m_sss;
279 : TElem * m_elem;
280 : Grid & m_grid;
281 : TAAPos & m_aaPos;
282 : const TFVGeom & m_geo;
283 : size_t m_co;
284 : AABox<MathVector<dim> > m_elem_bbox;
285 : size_t m_index; ///< index of the current source/sink
286 :
287 : public:
288 :
289 : /// class constructor (sets the iterator to the beginnging of the list)
290 0 : point_iterator
291 : (
292 : master_type* sss, /// the manager
293 : TElem* elem, ///< the element to find the sources and sinks for
294 : Grid& grid, ///< the grid
295 : TAAPos& aaPos, ///< position attachment for the grid
296 : const TFVGeom& geo, ///< FV geometry for the element
297 : size_t co ///< corner to check the sources and sinks for
298 : )
299 0 : : m_sss (sss), m_elem (elem), m_grid (grid), m_aaPos (aaPos), m_geo (geo), m_co (co)
300 : {
301 0 : m_elem_bbox = CalculateBoundingBox(m_elem, m_aaPos);
302 0 : next_sss (0);
303 0 : }
304 :
305 : /// copy constructor
306 : point_iterator (this_type& op)
307 : : m_elem (op.m_elem), m_grid (op.m_grid), m_aaPos (op.m_aaPos),
308 : m_geo (op.m_geo), m_co (op.m_co), m_elem_bbox (op.m_elem_bbox),
309 : m_index (op.m_index)
310 : {}
311 :
312 : /// moves the iterator to the beginning of the list
313 : this_type& rewind () {m_index = 0; return * this;}
314 :
315 : /// access operator
316 : point_sss_type * operator * () {return m_sss->ListP[m_index].get();}
317 :
318 : /// checks whether we are at the end
319 0 : bool is_over () {return m_index >= m_sss->num_points ();}
320 :
321 : /// moves to the next valid source or sink
322 : this_type & operator++ ()
323 : {
324 0 : next_sss (m_index + 1);
325 0 : return * this;
326 : }
327 :
328 : private:
329 : void next_sss (size_t index);
330 : };
331 :
332 : /// returns the number of the line sources and sinks
333 : size_t num_lines () {return ListL.size ();}
334 :
335 : /// return a line sink with a given index
336 : line_sss_type * line (size_t i) {return ListL[i].get();}
337 :
338 : /// calls the init-functions of all line sinks (if there are any in TLineData)
339 : void init_all_line_sss ()
340 : {
341 : for (size_t i = 0; i < ListL.size (); i++)
342 : ((TLineData *) ListL[i].get())->init ();
343 : }
344 :
345 : /// class of line source and sink iterator for a given scv
346 : template <typename TElem, typename TAAPos, typename TFVGeom>
347 : class line_iterator
348 : {
349 : typedef line_iterator<TElem, TAAPos, TFVGeom> this_type;
350 : typedef FVSingularSourcesAndSinks<dim, TPointData, TLineData> master_type;
351 :
352 : master_type * m_sss;
353 : TElem * m_elem;
354 : Grid & m_grid;
355 : TAAPos & m_aaPos;
356 : const TFVGeom & m_geo;
357 : size_t m_co;
358 : AABox<MathVector<dim> > m_elem_bbox;
359 :
360 : size_t m_index; ///< index of the current source/sink
361 : MathVector<dim> m_ls; ///< the 1st of the intersection points with the bnd of the scv
362 : MathVector<dim> m_le; ///< the 2nd of the intersection points with the bnd of the scv
363 :
364 : public:
365 :
366 : /// class constructor (sets the iterator to the beginnging of the list)
367 0 : line_iterator
368 : (
369 : master_type* sss, /// the manager
370 : TElem* elem, ///< the element to find the sources and sinks for
371 : Grid& grid, ///< the grid
372 : TAAPos& aaPos, ///< position attachment for the grid
373 : const TFVGeom& geo, ///< FV geometry for the element
374 : size_t co ///< corner to check the sources and sinks for
375 : )
376 0 : : m_sss (sss), m_elem (elem), m_grid (grid), m_aaPos (aaPos), m_geo (geo), m_co (co)
377 : {
378 0 : m_elem_bbox = CalculateBoundingBox(m_elem, m_aaPos);
379 0 : next_sss (0, m_ls, m_le);
380 0 : }
381 :
382 : /// copy constructor
383 : line_iterator (this_type& op)
384 : : m_elem (op.m_elem), m_grid (op.m_grid), m_aaPos (op.m_aaPos),
385 : m_geo (op.m_geo), m_co (op.m_co), m_elem_bbox (op.m_elem_bbox),
386 : m_index (op.m_index)
387 : {}
388 :
389 : /// moves the iterator to the beginning of the list
390 : this_type& rewind () {m_index = 0; return * this;}
391 :
392 : /// access operator
393 : line_sss_type * operator * () {return m_sss->ListL[m_index].get();}
394 :
395 : /// the 1st of the intersection points with the bnd of the scv
396 : const MathVector<dim>& seg_start () {return m_ls;}
397 :
398 : /// the 2nd of the intersection points with the bnd of the scv
399 : const MathVector<dim>& seg_end () {return m_le;}
400 :
401 : /// checks whether we are at the end
402 0 : bool is_over () {return m_index >= m_sss->num_lines ();}
403 :
404 : /// moves to the next valid source or sink
405 : this_type & operator++ ()
406 : {
407 0 : next_sss (m_index + 1, m_ls, m_le);
408 0 : return * this;
409 : }
410 :
411 : private:
412 : void next_sss (size_t index, MathVector<dim>& ls, MathVector<dim>& le);
413 : };
414 : };
415 :
416 : } // namespace ug
417 :
418 : #include "fv1_sss_impl.h"
419 :
420 : #endif // __H__UG__LIB_DISC__SPATIAL_DISC__FV_SINGULAR_SOURCES_AND_SINKS__
421 :
422 : /* End of File */
|