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__DATA_IMPORT__
34 : #define __H__UG__LIB_DISC__SPATIAL_DISC__DATA_IMPORT__
35 :
36 : #include "user_data.h"
37 :
38 : namespace ug{
39 :
40 : ////////////////////////////////////////////////////////////////////////////////
41 : // Data Import
42 : ////////////////////////////////////////////////////////////////////////////////
43 :
44 : enum DiscPart { NONE = 0,
45 : MASS = 1 << 0,
46 : STIFF = 1 << 1,
47 : RHS = 1 << 2,
48 : EXPL = 1 << 3,
49 : MAX_PART};
50 :
51 : inline
52 0 : std::ostream& operator<< (std::ostream& outStream, DiscPart part)
53 : {
54 0 : switch(part)
55 : {
56 0 : case NONE: outStream << "(none)"; break;
57 0 : case MASS: outStream << "mass"; break;
58 0 : case STIFF: outStream << "stiff"; break;
59 0 : case RHS: outStream << "rhs"; break;
60 0 : case EXPL: outStream << "expl"; break;
61 0 : default: UG_THROW("Unknown DiscPart in operator<<");
62 : }
63 0 : return outStream;
64 : };
65 :
66 : /// Base class for data import
67 : /**
68 : * An IDataImport is the base class for importing data to ElemDiscs
69 : */
70 : template <int dim>
71 : class IDataImport
72 : {
73 : public:
74 : /// Constructor
75 0 : IDataImport(bool compLinDefect = true)
76 0 : : m_spICplUserData(NULL), m_part(STIFF),
77 0 : m_bCompLinDefect(compLinDefect)
78 : {}
79 :
80 0 : virtual ~IDataImport() {}
81 :
82 : /// sets part of import
83 : void set_part(DiscPart part) {m_part = part;}
84 :
85 : /// set to no part
86 : void set_none_part() {m_part = NONE;}
87 :
88 : /// set to explicit part
89 0 : void set_expl_part() {m_part = EXPL;}
90 :
91 : /// sets if import is located in mass part (for time dependent problems)
92 0 : void set_mass_part() {m_part = MASS;}
93 :
94 : /// sets if import is located in rhs part
95 0 : void set_rhs_part() {m_part = RHS;}
96 :
97 : /// sets if import is located in stiff part (default)
98 : void set_stiff_part() {m_part = STIFF;}
99 :
100 : /// sets if lin defect is to be computed
101 0 : void set_comp_lin_defect(bool b) {m_bCompLinDefect = b;}
102 :
103 : /// returns if import is located in mass part (for time dependent problems)
104 0 : DiscPart part() const {return m_part;}
105 :
106 : /// returns if data is set
107 : virtual bool data_given() const = 0;
108 :
109 : /// returns if data is constant
110 : /**
111 : * This method, returns if the connected data is constant.
112 : */
113 : virtual bool constant() const = 0;
114 :
115 : /// returns if data depends on unknown functions
116 : /**
117 : * This method returns if the data depends on the unknown functions.
118 : */
119 : bool zero_derivative() const
120 : {
121 0 : if(!m_spICplUserData.valid()) return true;
122 0 : else if (m_spICplUserData->zero_derivative()) return true;
123 0 : else return !m_bCompLinDefect;
124 : }
125 :
126 : /// returns the connected user data
127 : virtual SmartPtr<ICplUserData<dim> > data() = 0;
128 :
129 : /// set function group for linearization of defect
130 : void set_map(const FunctionIndexMapping& map){m_map = map;}
131 :
132 : /// get function mapping
133 0 : const FunctionIndexMapping& map() const{return m_map;}
134 :
135 : /// get function mapping of dependent data
136 : const FunctionIndexMapping& conn_map() const{return m_spICplUserData->map();}
137 :
138 : /// number of functions
139 : size_t num_fct() const {return m_map.num_fct();}
140 :
141 : /// sets the geometric object type
142 : virtual void set_roid(ReferenceObjectID id) = 0;
143 :
144 : /// checks if ready for evaluation
145 : virtual void check_setup() = 0;
146 :
147 : /// compute lin defect
148 : virtual void compute_lin_defect(LocalVector& u) = 0;
149 :
150 : /// resize arrays
151 : virtual void update_dof_sizes(const LocalIndices& ind) = 0;
152 :
153 : /// add jacobian entries introduced by this import
154 : virtual void add_jacobian(LocalMatrix& J, const number scale) = 0;
155 :
156 : /// removes the positions
157 : virtual void clear_ips() = 0;
158 :
159 : protected:
160 : /// connected iexport
161 : SmartPtr<ICplUserData<dim> > m_spICplUserData;
162 :
163 : /// Mapping for import fct
164 : FunctionIndexMapping m_map;
165 :
166 : protected:
167 : /// flag to indicate where import is located
168 : DiscPart m_part;
169 :
170 : /// indicates iff lin defect should be computed
171 : bool m_bCompLinDefect;
172 : };
173 :
174 : /// Data import
175 : /**
176 : * A DataImport is used to import data into an ElemDisc.
177 : */
178 : template <typename TData, int dim>
179 : class DataImport : public IDataImport<dim>
180 : {
181 : public:
182 : /// Constructor
183 0 : DataImport(bool bLinDefect = true) : IDataImport<dim>(bLinDefect),
184 0 : m_id(ROID_UNKNOWN),
185 0 : m_seriesID(-1), m_spUserData(NULL), m_vValue(NULL),
186 0 : m_numIP(0), m_spDependentUserData(NULL)
187 0 : {clear_fct();
188 0 : }
189 :
190 : /// Destructor
191 : ~DataImport();
192 :
193 : /// set the user data
194 : void set_data(SmartPtr<CplUserData<TData, dim> > spData);
195 :
196 : /// returns the connected ICplUserData
197 0 : SmartPtr<ICplUserData<dim> > data() {return m_spUserData.template cast_dynamic<ICplUserData<dim> >();}
198 :
199 : /// returns the connected ICplUserData
200 : SmartPtr<CplUserData<TData, dim> > user_data(){return m_spUserData;}
201 :
202 : /// returns true if data given
203 0 : virtual bool data_given() const {return m_spUserData.valid();}
204 :
205 :
206 : /////////////////////////////////////////
207 : // Data
208 : /////////////////////////////////////////
209 :
210 : /// \copydoc IDataImport::constant()
211 0 : virtual bool constant() const
212 : {
213 0 : if(data_given()) return m_spUserData->constant();
214 : else return false;
215 : }
216 :
217 : /// returns the data value at ip
218 0 : const TData& operator[](size_t ip) const{check_ip(ip); return m_vValue[ip];}
219 :
220 : /// returns the data value at ip
221 0 : const TData* values() const {check_values(); return m_vValue;}
222 :
223 : /// return the derivative w.r.t to local function at ip
224 : const TData* deriv(size_t ip, size_t fct) const
225 : {
226 : UG_ASSERT(m_spDependentUserData.valid(), "No Dependent Data set");
227 : UG_ASSERT(m_seriesID >= 0, "No series ticket set");
228 : return m_spDependentUserData->deriv(m_seriesID, ip, fct);
229 : }
230 :
231 : /// return the derivative w.r.t to local function and dof at ip
232 : const TData& deriv(size_t ip, size_t fct, size_t dof) const
233 : {
234 : UG_ASSERT(m_spDependentUserData.valid(), "No Dependent Data set");
235 : UG_ASSERT(m_seriesID >= 0, "No series ticket set");
236 : return m_spDependentUserData->deriv(m_seriesID, ip, fct, dof);
237 : }
238 :
239 : /////////////////////////////////////////
240 : // Positions
241 : /////////////////////////////////////////
242 :
243 : /// number of integration points
244 0 : size_t num_ip() const {return m_numIP;}
245 :
246 : /// set the local integration points
247 : template <int ldim>
248 : void set_local_ips(const MathVector<ldim>* vPos, size_t numIP, int timePointSpec,
249 : bool bMayChange = true);
250 :
251 : /// set the local integration points
252 : void set_local_ips(const MathVector<dim>* vPos, size_t numIP, int timePointSpec,
253 : bool bMayChange = true);
254 :
255 : /// set the local integration points
256 : template <int ldim>
257 : void set_local_ips(const MathVector<ldim>* vPos, size_t numIP,
258 : bool bMayChange = true);
259 :
260 : /// set the local integration points
261 : void set_local_ips(const MathVector<dim>* vPos, size_t numIP,
262 : bool bMayChange = true);
263 :
264 : /// set the time point specification
265 : void set_time_point(int timePointSpec);
266 :
267 : /// sets the global positions
268 : void set_global_ips(const MathVector<dim>* vPos, size_t numIP);
269 :
270 : /// position of ip
271 : const MathVector<dim>& position(size_t i) const
272 : {
273 : if(data_given()) return m_spUserData->ip(m_seriesID, i);
274 : UG_THROW("DataImport::position: "
275 : "No Data set, but positions requested.");
276 : }
277 :
278 : /// removes the positions
279 : virtual void clear_ips();
280 :
281 : /////////////////////////////////////////
282 : // Linearization of Defect
283 : /////////////////////////////////////////
284 :
285 : /// number of shapes for local function
286 : size_t num_sh(size_t fct) const
287 : {
288 : UG_ASSERT(fct < m_vvNumDoFPerFct[fct], "Invalid index");
289 0 : return m_vvNumDoFPerFct[fct];
290 : }
291 :
292 : /// returns the pointer to all linearized defects at one ip
293 : TData* lin_defect(size_t ip, size_t fct)
294 : {check_ip_fct(ip,fct);return &(m_vvvLinDefect[ip][fct][0]);}
295 :
296 : /// returns the pointer to all linearized defects at one ip
297 : const TData* lin_defect(size_t ip, size_t fct) const
298 : {check_ip_fct(ip,fct);return &(m_vvvLinDefect[ip][fct][0]);}
299 :
300 : /// returns the linearized defect
301 : TData& lin_defect(size_t ip, size_t fct, size_t sh)
302 : {check_ip_fct_sh(ip,fct,sh);return m_vvvLinDefect[ip][fct][sh];}
303 :
304 : /// const access to lin defect
305 : const TData& lin_defect(size_t ip, size_t fct, size_t sh) const
306 : {check_ip_fct_sh(ip,fct,sh);return m_vvvLinDefect[ip][fct][sh];}
307 :
308 : /// compute jacobian for derivative w.r.t. non-system owned unknowns
309 : void add_jacobian(LocalMatrix& J, const number scale);
310 :
311 : /// resize lin defect arrays
312 : virtual void update_dof_sizes(const LocalIndices& ind);
313 :
314 : public:
315 : /// type of evaluation function
316 : typedef boost::function<void (const LocalVector& u,
317 : std::vector<std::vector<TData> > vvvLinDefect[],
318 : const size_t nip)> LinDefectFunc;
319 :
320 : /// sets the geometric object type
321 : virtual void set_roid(ReferenceObjectID id);
322 :
323 : /// checks if ready for evaluation
324 : virtual void check_setup();
325 :
326 : /// register evaluation of linear defect for a element
327 : template <typename TClass>
328 : void set_fct(ReferenceObjectID id, TClass* obj,
329 : void (TClass::*func)(
330 : const LocalVector& u,
331 : std::vector<std::vector<TData> > vvvLinDefect[],
332 : const size_t nip));
333 :
334 : /// register evaluation of linear defect for a element
335 : void set_fct(ReferenceObjectID id,
336 : void (*func)(
337 : const LocalVector& u,
338 : std::vector<std::vector<TData> > vvvLinDefect[],
339 : const size_t nip));
340 :
341 : /// clear all evaluation functions
342 : void clear_fct();
343 :
344 : /// compute lin defect
345 0 : virtual void compute_lin_defect(LocalVector& u)
346 : {
347 : /// compute the linearization only if the export parameter is 'at current time'
348 0 : if (! m_spUserData->at_current_time (m_seriesID))
349 : return;
350 : /// compute the linearization
351 : UG_ASSERT(m_vLinDefectFunc[m_id] != NULL, "No evaluation function.");
352 : UG_ASSERT(num_ip() == 0 || m_vvvLinDefect.size() >= num_ip(),
353 : "DataImport: Num ip "<<num_ip()<<", but memory: "<<m_vvvLinDefect.size());
354 : u.access_by_map(this->map());
355 0 : (m_vLinDefectFunc[m_id])(u, &m_vvvLinDefect[0], m_numIP);
356 : }
357 :
358 : protected:
359 : /// checks in debug mode the correct index
360 : inline void check_ip_fct(size_t ip, size_t fct) const;
361 :
362 : /// checks in debug mode the correct index
363 : inline void check_ip_fct_sh(size_t ip, size_t fct, size_t sh) const;
364 :
365 : /// checks in debug mode the correct index
366 : inline void check_ip(size_t ip) const;
367 :
368 : /// checks in debug mode the correct index
369 : inline void check_values() const;
370 :
371 : /// caches data access
372 : void cache_data_access();
373 :
374 : /// resizes the lin defect arrays for current number of ips.
375 : void resize_defect_array();
376 :
377 : /// current Geom Object
378 : ReferenceObjectID m_id;
379 :
380 : /// function pointers for all elem types
381 : LinDefectFunc m_vLinDefectFunc[NUM_REFERENCE_OBJECTS];
382 :
383 : /// series number provided by export
384 : int m_seriesID;
385 :
386 : /// connected UserData
387 : SmartPtr<CplUserData<TData, dim> > m_spUserData;
388 :
389 : /// cached access to the UserData field
390 : const TData* m_vValue;
391 :
392 : /// number of ips
393 : size_t m_numIP;
394 :
395 : /// connected export (if depended data)
396 : SmartPtr<DependentUserData<TData, dim> > m_spDependentUserData;
397 :
398 : /// number of functions and their dofs
399 : std::vector<size_t> m_vvNumDoFPerFct;
400 :
401 : /// linearized defect (num_ip) x (num_fct) x (num_dofs(i))
402 : std::vector<std::vector<std::vector<TData> > > m_vvvLinDefect;
403 : };
404 :
405 : } // end namespace ug
406 :
407 : // include implementation
408 : #include "data_import_impl.h"
409 :
410 : #endif /* __H__UG__LIB_DISC__SPATIAL_DISC__DATA_IMPORT__ */
|
