Loading [MathJax]/extensions/tex2jax.js
ug4
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
data_import.h
Go to the documentation of this file.
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
38namespace ug{
39
41// Data Import
43
44enum DiscPart { NONE = 0,
45 MASS = 1 << 0,
46 STIFF = 1 << 1,
47 RHS = 1 << 2,
48 EXPL = 1 << 3,
50
51inline
52std::ostream& operator<< (std::ostream& outStream, DiscPart part)
53{
54 switch(part)
55 {
56 case NONE: outStream << "(none)"; break;
57 case MASS: outStream << "mass"; break;
58 case STIFF: outStream << "stiff"; break;
59 case RHS: outStream << "rhs"; break;
60 case EXPL: outStream << "expl"; break;
61 default: UG_THROW("Unknown DiscPart in operator<<");
62 }
63 return outStream;
64};
65
67
70template <int dim>
72{
73 public:
75 IDataImport(bool compLinDefect = true)
77 m_bCompLinDefect(compLinDefect)
78 {}
79
80 virtual ~IDataImport() {}
81
84
87
90
93
96
99
102
104 DiscPart part() const {return m_part;}
105
107 virtual bool data_given() const = 0;
108
110
113 virtual bool constant() const = 0;
114
116
119 bool zero_derivative() const
120 {
121 if(!m_spICplUserData.valid()) return true;
122 else if (m_spICplUserData->zero_derivative()) return true;
123 else return !m_bCompLinDefect;
124 }
125
128
131
133 const FunctionIndexMapping& map() const{return m_map;}
134
136 const FunctionIndexMapping& conn_map() const{return m_spICplUserData->map();}
137
139 size_t num_fct() const {return m_map.num_fct();}
140
142 virtual void set_roid(ReferenceObjectID id) = 0;
143
145 virtual void check_setup() = 0;
146
148 virtual void compute_lin_defect(LocalVector& u) = 0;
149
151 virtual void update_dof_sizes(const LocalIndices& ind) = 0;
152
154 virtual void add_jacobian(LocalMatrix& J, const number scale) = 0;
155
157 virtual void clear_ips() = 0;
158
159 protected:
162
165
166 protected:
169
172};
173
175
178template <typename TData, int dim>
179class DataImport : public IDataImport<dim>
180{
181 public:
183 DataImport(bool bLinDefect = true) : IDataImport<dim>(bLinDefect),
185 m_seriesID(-1), m_spUserData(NULL), m_vValue(NULL),
187 {clear_fct();
188 }
189
191 ~DataImport();
192
195
197 SmartPtr<ICplUserData<dim> > data() {return m_spUserData.template cast_dynamic<ICplUserData<dim> >();}
198
201
203 virtual bool data_given() const {return m_spUserData.valid();}
204
205
207 // Data
209
211 virtual bool constant() const
212 {
213 if(data_given()) return m_spUserData->constant();
214 else return false;
215 }
216
218 const TData& operator[](size_t ip) const{check_ip(ip); return m_vValue[ip];}
219
221 const TData* values() const {check_values(); return m_vValue;}
222
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
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
240 // Positions
242
244 size_t num_ip() const {return m_numIP;}
245
247 template <int ldim>
248 void set_local_ips(const MathVector<ldim>* vPos, size_t numIP, int timePointSpec,
249 bool bMayChange = true);
250
252 void set_local_ips(const MathVector<dim>* vPos, size_t numIP, int timePointSpec,
253 bool bMayChange = true);
254
256 template <int ldim>
257 void set_local_ips(const MathVector<ldim>* vPos, size_t numIP,
258 bool bMayChange = true);
259
261 void set_local_ips(const MathVector<dim>* vPos, size_t numIP,
262 bool bMayChange = true);
263
265 void set_time_point(int timePointSpec);
266
268 void set_global_ips(const MathVector<dim>* vPos, size_t numIP);
269
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
279 virtual void clear_ips();
280
282 // Linearization of Defect
284
286 size_t num_sh(size_t fct) const
287 {
288 UG_ASSERT(fct < m_vvNumDoFPerFct[fct], "Invalid index");
289 return m_vvNumDoFPerFct[fct];
290 }
291
293 TData* lin_defect(size_t ip, size_t fct)
294 {check_ip_fct(ip,fct);return &(m_vvvLinDefect[ip][fct][0]);}
295
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
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
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
309 void add_jacobian(LocalMatrix& J, const number scale);
310
312 virtual void update_dof_sizes(const LocalIndices& ind);
313
314 public:
316 typedef boost::function<void (const LocalVector& u,
317 std::vector<std::vector<TData> > vvvLinDefect[],
318 const size_t nip)> LinDefectFunc;
319
321 virtual void set_roid(ReferenceObjectID id);
322
324 virtual void check_setup();
325
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
336 void (*func)(
337 const LocalVector& u,
338 std::vector<std::vector<TData> > vvvLinDefect[],
339 const size_t nip));
340
342 void clear_fct();
343
346 {
348 if (! m_spUserData->at_current_time (m_seriesID))
349 return;
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());
356 }
357
358 protected:
360 inline void check_ip_fct(size_t ip, size_t fct) const;
361
363 inline void check_ip_fct_sh(size_t ip, size_t fct, size_t sh) const;
364
366 inline void check_ip(size_t ip) const;
367
369 inline void check_values() const;
370
372 void cache_data_access();
373
375 void resize_defect_array();
376
379
382
385
388
390 const TData* m_vValue;
391
393 size_t m_numIP;
394
397
399 std::vector<size_t> m_vvNumDoFPerFct;
400
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__ */
Definition smart_pointer.h:108
Type based UserData.
Definition user_data.h:501
Data import.
Definition data_import.h:180
virtual bool data_given() const
returns true if data given
Definition data_import.h:203
void add_jacobian(LocalMatrix &J, const number scale)
compute jacobian for derivative w.r.t. non-system owned unknowns
Definition data_import_impl.h:238
virtual void clear_ips()
removes the positions
Definition data_import_impl.h:228
void check_ip(size_t ip) const
checks in debug mode the correct index
Definition data_import_impl.h:324
virtual void update_dof_sizes(const LocalIndices &ind)
resize lin defect arrays
Definition data_import_impl.h:271
void clear_fct()
clear all evaluation functions
Definition data_import_impl.h:109
virtual void set_roid(ReferenceObjectID id)
sets the geometric object type
Definition data_import_impl.h:51
void set_data(SmartPtr< CplUserData< TData, dim > > spData)
set the user data
Definition data_import_impl.h:117
~DataImport()
Destructor.
Definition data_import_impl.h:45
const TData * lin_defect(size_t ip, size_t fct) const
returns the pointer to all linearized defects at one ip
Definition data_import.h:297
LinDefectFunc m_vLinDefectFunc[NUM_REFERENCE_OBJECTS]
function pointers for all elem types
Definition data_import.h:381
SmartPtr< DependentUserData< TData, dim > > m_spDependentUserData
connected export (if depended data)
Definition data_import.h:396
void set_global_ips(const MathVector< dim > *vPos, size_t numIP)
sets the global positions
Definition data_import_impl.h:217
void cache_data_access()
caches data access
Definition data_import_impl.h:130
const TData * deriv(size_t ip, size_t fct) const
return the derivative w.r.t to local function at ip
Definition data_import.h:224
virtual bool constant() const
returns if data is constant
Definition data_import.h:211
TData & lin_defect(size_t ip, size_t fct, size_t sh)
returns the linearized defect
Definition data_import.h:301
ReferenceObjectID m_id
current Geom Object
Definition data_import.h:378
const TData & lin_defect(size_t ip, size_t fct, size_t sh) const
const access to lin defect
Definition data_import.h:305
SmartPtr< ICplUserData< dim > > data()
returns the connected ICplUserData
Definition data_import.h:197
virtual void compute_lin_defect(LocalVector &u)
compute lin defect
Definition data_import.h:345
std::vector< size_t > m_vvNumDoFPerFct
number of functions and their dofs
Definition data_import.h:399
const MathVector< dim > & position(size_t i) const
position of ip
Definition data_import.h:271
DataImport(bool bLinDefect=true)
Constructor.
Definition data_import.h:183
std::vector< std::vector< std::vector< TData > > > m_vvvLinDefect
linearized defect (num_ip) x (num_fct) x (num_dofs(i))
Definition data_import.h:402
const TData & deriv(size_t ip, size_t fct, size_t dof) const
return the derivative w.r.t to local function and dof at ip
Definition data_import.h:232
void resize_defect_array()
resizes the lin defect arrays for current number of ips.
Definition data_import_impl.h:287
void check_ip_fct_sh(size_t ip, size_t fct, size_t sh) const
checks in debug mode the correct index
Definition data_import_impl.h:317
void set_fct(ReferenceObjectID id, TClass *obj, void(TClass::*func)(const LocalVector &u, std::vector< std::vector< TData > > vvvLinDefect[], const size_t nip))
register evaluation of linear defect for a element
Definition data_import_impl.h:83
SmartPtr< CplUserData< TData, dim > > user_data()
returns the connected ICplUserData
Definition data_import.h:200
SmartPtr< CplUserData< TData, dim > > m_spUserData
connected UserData
Definition data_import.h:387
TData * lin_defect(size_t ip, size_t fct)
returns the pointer to all linearized defects at one ip
Definition data_import.h:293
virtual void check_setup()
checks if ready for evaluation
Definition data_import_impl.h:60
void set_local_ips(const MathVector< ldim > *vPos, size_t numIP, int timePointSpec, bool bMayChange=true)
set the local integration points
Definition data_import_impl.h:141
boost::function< void(const LocalVector &u, std::vector< std::vector< TData > > vvvLinDefect[], const size_t nip)> LinDefectFunc
type of evaluation function
Definition data_import.h:318
void check_values() const
checks in debug mode the correct index
Definition data_import_impl.h:330
const TData * m_vValue
cached access to the UserData field
Definition data_import.h:390
size_t num_sh(size_t fct) const
number of shapes for local function
Definition data_import.h:286
int m_seriesID
series number provided by export
Definition data_import.h:384
const TData * values() const
returns the data value at ip
Definition data_import.h:221
void check_ip_fct(size_t ip, size_t fct) const
checks in debug mode the correct index
Definition data_import_impl.h:309
size_t m_numIP
number of ips
Definition data_import.h:393
size_t num_ip() const
number of integration points
Definition data_import.h:244
const TData & operator[](size_t ip) const
returns the data value at ip
Definition data_import.h:218
void set_time_point(int timePointSpec)
set the time point specification
Definition data_import_impl.h:211
describes a mapping between two local index sets
Definition function_group.h:186
size_t num_fct() const
returns the number of indices that are mapped
Definition function_group.h:195
Base class for data import.
Definition data_import.h:72
DiscPart m_part
flag to indicate where import is located
Definition data_import.h:168
void set_mass_part()
sets if import is located in mass part (for time dependent problems)
Definition data_import.h:92
void set_part(DiscPart part)
sets part of import
Definition data_import.h:83
void set_map(const FunctionIndexMapping &map)
set function group for linearization of defect
Definition data_import.h:130
const FunctionIndexMapping & map() const
get function mapping
Definition data_import.h:133
bool m_bCompLinDefect
indicates iff lin defect should be computed
Definition data_import.h:171
DiscPart part() const
returns if import is located in mass part (for time dependent problems)
Definition data_import.h:104
virtual ~IDataImport()
Definition data_import.h:80
void set_stiff_part()
sets if import is located in stiff part (default)
Definition data_import.h:98
FunctionIndexMapping m_map
Mapping for import fct.
Definition data_import.h:164
virtual void add_jacobian(LocalMatrix &J, const number scale)=0
add jacobian entries introduced by this import
void set_rhs_part()
sets if import is located in rhs part
Definition data_import.h:95
bool zero_derivative() const
returns if data depends on unknown functions
Definition data_import.h:119
virtual void check_setup()=0
checks if ready for evaluation
virtual void set_roid(ReferenceObjectID id)=0
sets the geometric object type
virtual bool constant() const =0
returns if data is constant
void set_none_part()
set to no part
Definition data_import.h:86
virtual void clear_ips()=0
removes the positions
virtual bool data_given() const =0
returns if data is set
void set_comp_lin_defect(bool b)
sets if lin defect is to be computed
Definition data_import.h:101
virtual void update_dof_sizes(const LocalIndices &ind)=0
resize arrays
size_t num_fct() const
number of functions
Definition data_import.h:139
IDataImport(bool compLinDefect=true)
Constructor.
Definition data_import.h:75
virtual void compute_lin_defect(LocalVector &u)=0
compute lin defect
SmartPtr< ICplUserData< dim > > m_spICplUserData
connected iexport
Definition data_import.h:161
void set_expl_part()
set to explicit part
Definition data_import.h:89
const FunctionIndexMapping & conn_map() const
get function mapping of dependent data
Definition data_import.h:136
virtual SmartPtr< ICplUserData< dim > > data()=0
returns the connected user data
Definition local_algebra.h:50
Definition local_algebra.h:422
Definition local_algebra.h:198
void access_by_map(const FunctionIndexMapping &funcMap)
access only part of the functions using mapping (restrict functions)
Definition local_algebra.h:306
a mathematical Vector with N entries.
Definition math_vector.h:97
std::ostream & operator<<(std::ostream &outStream, const ug::MathMatrix< 2, 2 > &m)
Definition math_matrix.cpp:38
#define UG_ASSERT(expr, msg)
Definition assert.h:70
#define UG_THROW(msg)
Definition error.h:57
double number
Definition types.h:124
Definition smart_pointer.h:814
the ug namespace
ReferenceObjectID
these ids are used to identify the shape of a geometric object.
Definition grid_base_objects.h:74
@ NUM_REFERENCE_OBJECTS
Definition grid_base_objects.h:85
@ ROID_UNKNOWN
Definition grid_base_objects.h:75
DiscPart
Definition data_import.h:44
@ NONE
Definition data_import.h:44
@ STIFF
Definition data_import.h:46
@ MAX_PART
Definition data_import.h:49
@ RHS
Definition data_import.h:47
@ EXPL
Definition data_import.h:48
@ MASS
Definition data_import.h:45
function func(x, y, z, t, si)