Line data Source code
1 : /*
2 : * Copyright (c) 2012-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__REFERENCE_ELEMENT__REFERENCE_MAPPING_PROVIDER__
34 : #define __H__UG__LIB_DISC__REFERENCE_ELEMENT__REFERENCE_MAPPING_PROVIDER__
35 :
36 :
37 : #include "common/common.h"
38 : #include "common/math/ugmath.h"
39 : #include "lib_grid/grid/grid_base_objects.h"
40 :
41 : namespace ug{
42 :
43 : /// virtual base class for reference mappings
44 : /**
45 : * This class is the base class for reference mappings in order to make them
46 : * selectable through a provider (on the price of virtual functions).
47 : *
48 : * \tparam TDim reference element dimension
49 : * \tparam TWorldDim (physical) world dimension
50 : */
51 : template <int TDim, int TWorldDim = TDim>
52 : class DimReferenceMapping
53 : {
54 : public:
55 : /// world dimension (range space dimension)
56 : static const int worldDim = TWorldDim;
57 :
58 : /// reference dimension (domain space dimension)
59 : static const int dim = TDim;
60 :
61 : public:
62 : /// returns if mapping is affine
63 : virtual bool is_linear() const = 0;
64 :
65 : /// refresh mapping for new set of corners
66 : virtual void update(const MathVector<worldDim>* vCornerCoord) = 0;
67 :
68 : /// refresh mapping for new set of corners
69 : virtual void update(const std::vector<MathVector<worldDim> >& vCornerCoord) = 0;
70 :
71 : /// map local coordinate to global coordinate
72 : virtual void local_to_global(MathVector<worldDim>& globPos,
73 : const MathVector<dim>& locPos) const = 0;
74 :
75 : /// map n local coordinate to global coordinate
76 : virtual void local_to_global(MathVector<worldDim>* vGlobPos,
77 : const MathVector<dim>* vLocPos, size_t n) const = 0;
78 :
79 : /// map local coordinate to global coordinate for a vector of local positions
80 : virtual void local_to_global(std::vector<MathVector<worldDim> >& vGlobPos,
81 : const std::vector<MathVector<dim> >& vLocPos) const = 0;
82 :
83 : /// map global coordinate to local coordinate
84 : virtual void global_to_local(MathVector<dim>& locPos,
85 : const MathVector<worldDim>& globPos,
86 : const size_t maxIter = 1000,
87 : const number tol = 1e-10) const = 0;
88 :
89 : /// map global coordinate to local coordinate for n local positions
90 : virtual void global_to_local(MathVector<dim>* vLocPos,
91 : const MathVector<worldDim>* vGlobPos, size_t n,
92 : const size_t maxIter = 1000,
93 : const number tol = 1e-10) const = 0;
94 :
95 : /// map global coordinate to local coordinate for a vector of local positions
96 : virtual void global_to_local(std::vector<MathVector<dim> >& vLocPos,
97 : const std::vector<MathVector<worldDim> >& vGlobPos,
98 : const size_t maxIter = 1000,
99 : const number tol = 1e-10) const = 0;
100 :
101 : /// returns jacobian
102 : virtual void jacobian(MathMatrix<worldDim, dim>& J,
103 : const MathVector<dim>& locPos) const = 0;
104 :
105 : /// returns jacobian for n local positions
106 : virtual void jacobian(MathMatrix<worldDim, dim>* vJ,
107 : const MathVector<dim>* vLocPos, size_t n) const = 0;
108 :
109 : /// returns jacobian for a vector of local positions
110 : virtual void jacobian(std::vector<MathMatrix<worldDim, dim> >& vJ,
111 : const std::vector<MathVector<dim> >& vLocPos) const = 0;
112 :
113 : /// returns transposed of jacobian
114 : virtual void jacobian_transposed(MathMatrix<dim, worldDim>& JT,
115 : const MathVector<dim>& locPos) const = 0;
116 :
117 : /// returns transposed of jacobian for n local positions
118 : virtual void jacobian_transposed(MathMatrix<dim, worldDim>* vJT,
119 : const MathVector<dim>* vLocPos, size_t n) const = 0;
120 :
121 : /// returns transposed of jacobian for a vector of positions
122 : virtual void jacobian_transposed(std::vector<MathMatrix<dim, worldDim> >& vJT,
123 : const std::vector<MathVector<dim> >& vLocPos) const = 0;
124 :
125 : /// returns transposed of the inverse of the jacobian and returns sqrt of gram determinante
126 : virtual number jacobian_transposed_inverse(MathMatrix<worldDim, dim>& JTInv,
127 : const MathVector<dim>& locPos) const = 0;
128 :
129 : /// returns transposed of the inverse of the jacobian for n local positions
130 : virtual void jacobian_transposed_inverse(MathMatrix<worldDim, dim>* vJTInv,
131 : const MathVector<dim>* vLocPos, size_t n) const = 0;
132 :
133 : /// returns transposed of the inverse of the jacobian for n local positions
134 : virtual void jacobian_transposed_inverse(MathMatrix<worldDim, dim>* vJTInv,
135 : number* vDet,
136 : const MathVector<dim>* vLocPos, size_t n) const = 0;
137 :
138 : /// returns transposed of the inverse of the jacobian for a vector of positions
139 : virtual void jacobian_transposed_inverse(std::vector<MathMatrix<worldDim, dim> >& vJTInv,
140 : const std::vector<MathVector<dim> >& vLocPos) const = 0;
141 :
142 : /// returns transposed of the inverse of the jacobian for a vector of positions
143 : virtual void jacobian_transposed_inverse(std::vector<MathMatrix<worldDim, dim> >& vJTInv,
144 : std::vector<number>& vDet,
145 : const std::vector<MathVector<dim> >& vLocPos) const = 0;
146 :
147 : /// returns the determinate of the jacobian
148 : virtual number sqrt_gram_det(const MathVector<dim>& locPos) const = 0;
149 :
150 : /// returns the determinate of the jacobian for n local positions
151 : virtual void sqrt_gram_det(number* vDet, const MathVector<dim>* vLocPos, size_t n) const = 0;
152 :
153 : /// returns the determinate of the jacobian for a vector of local positions
154 : virtual void sqrt_gram_det(std::vector<number>& vDet,
155 : const std::vector<MathVector<dim> >& vLocPos) const = 0;
156 :
157 : /// virtual destructor
158 : virtual ~DimReferenceMapping() {}
159 : };
160 :
161 : ////////////////////////////////////////////////////////////////////////////////
162 : ////////////////////////////////////////////////////////////////////////////////
163 :
164 : /// class to provide reference mappings
165 : /**
166 : * This class provides references mappings. It is implemented as a Singleton.
167 : */
168 : class ReferenceMappingProvider {
169 : private:
170 : // disallow private constructor
171 : ReferenceMappingProvider();
172 :
173 : // disallow copy and assignment (intentionally left unimplemented)
174 : ReferenceMappingProvider(const ReferenceMappingProvider&);
175 : ReferenceMappingProvider& operator=(const ReferenceMappingProvider&);
176 :
177 : // private destructor
178 0 : ~ReferenceMappingProvider(){};
179 :
180 : // Singleton provider
181 0 : static ReferenceMappingProvider& inst()
182 : {
183 0 : static ReferenceMappingProvider myInst;
184 0 : return myInst;
185 : };
186 :
187 : // This is very dirty implementation, since casting to void. But, it is
188 : // efficient, easy and typesafe. Maybe it should be changed to something
189 : // inherent typesafe (not using casts) some day
190 : // holding all mappings (worldDim x dim x roid)
191 : void* m_vvvMapping[4][4][NUM_REFERENCE_OBJECTS];
192 :
193 : // casts void to map
194 : template <int TDim, int TWorldDim>
195 : DimReferenceMapping<TDim, TWorldDim>* get_mapping(ReferenceObjectID roid)
196 : {
197 : UG_STATIC_ASSERT(TDim <= 3, only_implemented_for_ref_dim_smaller_equal_3);
198 : UG_STATIC_ASSERT(TWorldDim <= 3, only_implemented_for_ref_dim_smaller_equal_3);
199 : UG_ASSERT(roid < NUM_REFERENCE_OBJECTS, "Roid specified incorrectly.");
200 : UG_ASSERT(roid >= 0, "Roid specified incorrectly.");
201 : return reinterpret_cast<DimReferenceMapping<TDim, TWorldDim>*>
202 0 : (m_vvvMapping[TDim][TWorldDim][roid]);
203 : }
204 :
205 : // casts map to void
206 : template <int TDim, int TWorldDim>
207 : void set_mapping(ReferenceObjectID roid, DimReferenceMapping<TDim, TWorldDim>& map)
208 : {
209 0 : m_vvvMapping[TDim][TWorldDim][roid] = reinterpret_cast<void*>(&map);
210 : }
211 :
212 : public:
213 : /// returns a reference to a DimReferenceMapping
214 : /**
215 : * This class returns a reference mapping for a ReferenceObjectID. The
216 : * reference dimension and the world dimension must be chosen as
217 : * template arguments. An exception is throw if such an mapping does not
218 : * exist.
219 : *
220 : * \param[in] roid Reference Object ID
221 : * \tparam TDim reference element dimension
222 : * \tparam TWorldDim (physical) world dimension
223 : */
224 : template <int TDim, int TWorldDim>
225 0 : static DimReferenceMapping<TDim, TWorldDim>& get(ReferenceObjectID roid)
226 : {
227 0 : DimReferenceMapping<TDim, TWorldDim>* pMap = inst().get_mapping<TDim, TWorldDim>(roid);
228 0 : if(!pMap){
229 0 : UG_THROW("ReferenceMappingProvider: ReferenceMapping not found for "
230 : <<roid<<" from R^"<<TDim<<" to R^"<<TWorldDim);
231 : }
232 0 : else return *pMap;
233 : }
234 :
235 : /// returns a reference to a DimReferenceMapping with updated element corners
236 : /**
237 : * This class returns a reference mapping for a ReferenceObjectID. The
238 : * reference dimension and the world dimension must be chosen as
239 : * template arguments. An exception is throw if such an mapping does not
240 : * exist.
241 : *
242 : * \param[in] roid Reference Object ID
243 : * \param[in] vCornerCoord The corner coordinates of the element
244 : * \tparam TDim reference element dimension
245 : * \tparam TWorldDim (physical) world dimension
246 : */
247 : template <int TDim, int TWorldDim>
248 : static DimReferenceMapping<TDim, TWorldDim>& get(ReferenceObjectID roid,
249 : const std::vector<MathVector<TWorldDim> >& vCornerCoord)
250 : {
251 0 : DimReferenceMapping<TDim, TWorldDim>& map = get<TDim, TWorldDim>(roid);
252 0 : map.update(vCornerCoord);
253 : return map;
254 : }
255 :
256 : /// returns a reference to a DimReferenceMapping with updated element corners
257 : /**
258 : * This class returns a reference mapping for a ReferenceObjectID. The
259 : * reference dimension and the world dimension must be chosen as
260 : * template arguments. An exception is throw if such an mapping does not
261 : * exist.
262 : *
263 : * \param[in] roid Reference Object ID
264 : * \param[in] vCornerCoord The corner coordinates of the element
265 : * \tparam TDim reference element dimension
266 : * \tparam TWorldDim (physical) world dimension
267 : */
268 : template <int TDim, int TWorldDim>
269 : static DimReferenceMapping<TDim, TWorldDim>& get(ReferenceObjectID roid,
270 : const MathVector<TWorldDim>* vCornerCoord)
271 : {
272 0 : DimReferenceMapping<TDim, TWorldDim>& map = get<TDim, TWorldDim>(roid);
273 0 : map.update(vCornerCoord);
274 : return map;
275 : }
276 :
277 : };
278 :
279 : } // end namespace ug
280 :
281 : #endif /* __H__UG__LIB_DISC__REFERENCE_ELEMENT__REFERENCE_MAPPING_PROVIDER__ */
|
