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__FUNCTION_SPACES__INTERPOLATE__
34 : #define __H__UG__LIB_DISC__FUNCTION_SPACES__INTERPOLATE__
35 :
36 : #include "common/common.h"
37 : #include "common/util/smart_pointer.h"
38 :
39 : #include "lib_grid/tools/subset_group.h"
40 :
41 : #include "lib_disc/domain_util.h"
42 : #include "lib_disc/common/groups_util.h"
43 : #include "lib_disc/local_finite_element/local_finite_element_provider.h"
44 : #include "lib_disc/spatial_disc/user_data/const_user_data.h"
45 : #include "lib_disc/reference_element/reference_mapping.h"
46 :
47 : #ifdef UG_FOR_LUA
48 : #include "bindings/lua/lua_user_data.h"
49 : #endif
50 :
51 : namespace ug{
52 :
53 : ////////////////////////////////////////////////////////////////////////////////
54 : // Interpolate on Vertices only
55 : ////////////////////////////////////////////////////////////////////////////////
56 :
57 : /**
58 : * This function interpolates a grid function on a vertex loop. Thus, it can only
59 : * be used if all degrees of freedom are located in the vertices only (e.g. P1
60 : * finite elements). In those cases it is faster than the element by element
61 : * loop.
62 : *
63 : * @param[in] spInterpolFunction data providing interpolation values
64 : * @param[out] spGridFct interpolated grid function
65 : * @param[in] fct symbolic name of function component
66 : * @param[in] ssGrp subsets, where to interpolate
67 : * @param[in] time time point
68 : * @param[in] diff_pos different vector between Interpolatin values and spGridFct.
69 : */
70 :
71 : template <typename TGridFunction>
72 0 : void InterpolateOnDiffVertices(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
73 : SmartPtr<TGridFunction> spGridFct,
74 : size_t fct,
75 : number time,
76 : const SubsetGroup& ssGrp,
77 : const MathVector<TGridFunction::dim> diff_pos)
78 : {
79 :
80 : // domain type and position_type
81 : typedef typename TGridFunction::domain_type domain_type;
82 : typedef typename domain_type::position_type position_type;
83 :
84 : //std::cout<<"Interpolate diff_vector: "<<diff_pos;
85 :
86 : // get position accessor (of interpolated grid function)
87 : const typename domain_type::position_accessor_type& aaPos
88 0 : = spGridFct->domain()->position_accessor();
89 :
90 :
91 : std::vector<DoFIndex> ind;
92 : typename TGridFunction::template dim_traits<0>::const_iterator iterEnd, iter;
93 :
94 0 : for(size_t i = 0; i < ssGrp.size(); ++i)
95 : {
96 : // get subset index
97 : const int si = ssGrp[i];
98 :
99 : // skip if function is not defined in subset
100 0 : if(!spGridFct->is_def_in_subset(fct, si)) continue;
101 :
102 : // iterate over all elements
103 0 : iterEnd = spGridFct->template end<Vertex>(si);
104 0 : iter = spGridFct->template begin<Vertex>(si);
105 0 : for(; iter != iterEnd; ++iter)
106 : {
107 : // get element
108 : Vertex* vrt = *iter;
109 :
110 : // global position
111 : position_type glob_pos = aaPos[vrt]; // position (of interpolated grid function)
112 : position_type rel_pos=glob_pos;
113 : rel_pos -=diff_pos;
114 :
115 : // value at position
116 : number val;
117 0 : (*spInterpolFunction)(val, rel_pos, time, si);
118 :
119 : // get multiindices of element
120 : spGridFct->dof_indices(vrt, fct, ind);
121 :
122 : // loop all dofs
123 0 : for(size_t i = 0; i < ind.size(); ++i)
124 : {
125 : // set value
126 0 : DoFRef(*spGridFct, ind[i]) = val;
127 : }
128 : }
129 : }
130 0 : }
131 : //getting value of spInterpolFunction at position
132 :
133 :
134 : template <typename TGridFunction>
135 : number get_number_on_coords(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
136 : typename TGridFunction::domain_type::position_type pos,
137 : number time,
138 : const int si
139 : ){
140 : number val;
141 : (*spInterpolFunction)(val, pos, time, si);
142 :
143 : return val;
144 : }
145 :
146 : /**
147 : * This function interpolates a grid function on a vertex loop. Thus, it can only
148 : * be used if all degrees of freedom are located in the vertices only (e.g. P1
149 : * finite elements). In those cases it is faster than the element by element
150 : * loop.
151 : *
152 : * @param[in] spInterpolFunction data providing interpolation values
153 : * @param[out] spGridFct interpolated grid function
154 : * @param[in] fct symbolic name of function component
155 : * @param[in] ssGrp subsets, where to interpolate
156 : * @param[in] time time point
157 : */
158 : template <typename TGridFunction>
159 0 : void InterpolateOnVertices(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
160 : SmartPtr<TGridFunction> spGridFct,
161 : size_t fct,
162 : number time,
163 : const SubsetGroup& ssGrp)
164 : {
165 : // dimension of reference element
166 : const int dim = TGridFunction::dim;
167 :
168 : MathVector<dim> diff_pos(0.0);
169 0 : InterpolateOnDiffVertices<TGridFunction>(spInterpolFunction, spGridFct, fct, time, ssGrp, diff_pos);
170 0 : }
171 :
172 :
173 :
174 : ////////////////////////////////////////////////////////////////////////////////
175 : // Interpolate on Elements
176 : ////////////////////////////////////////////////////////////////////////////////
177 :
178 : /**
179 : * This function interpolates a grid function on an element by element loop. On
180 : * each element the all associated (up to the boundary of the element) are
181 : * interpolated and the values are stored in the grid function.
182 : *
183 : * @param[in] spInterpolFunction data providing interpolation values
184 : * @param[out] spGridFct interpolated grid function
185 : * @param[in] fct symbolic name of function component
186 : * @param[in] si subset, where to interpolate
187 : * @param[in] time time point
188 : */
189 : template <typename TElem, typename TGridFunction>
190 0 : void InterpolateOnDiffElements(
191 : SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
192 : SmartPtr<TGridFunction> spGridFct, size_t fct, int si, number time,
193 : const MathVector<TGridFunction::dim> diff_pos)
194 : {
195 : // get reference element type
196 : typedef typename reference_element_traits<TElem>::reference_element_type ref_elem_type;
197 : const ReferenceObjectID roid = ref_elem_type::REFERENCE_OBJECT_ID;
198 :
199 : // dimension of reference element
200 : const int dim = ref_elem_type::dim;
201 :
202 :
203 : // domain type and position_type
204 : typedef typename TGridFunction::domain_type domain_type;
205 : typedef typename domain_type::position_type position_type;
206 :
207 : // get iterators
208 : typename TGridFunction::template traits<TElem>::const_iterator iterEnd, iter;
209 0 : iterEnd = spGridFct->template end<TElem>(si);
210 0 : iter = spGridFct->template begin<TElem>(si);
211 :
212 : // check if something to do:
213 0 : if(iter == iterEnd) return;
214 :
215 : // id of shape functions used
216 0 : LFEID id = spGridFct->local_finite_element_id(fct);
217 :
218 : // get trial space
219 : const LocalShapeFunctionSet<dim>& trialSpace =
220 : LocalFiniteElementProvider::get<dim>(roid, id);
221 :
222 : // number of dofs on element
223 0 : const size_t nsh = trialSpace.num_sh();
224 :
225 : // load local positions of dofs for the trial space on element
226 0 : std::vector<MathVector<dim> > loc_pos(nsh);
227 0 : for(size_t i = 0; i < nsh; ++i)
228 0 : if(!trialSpace.position(i, loc_pos[i]))
229 0 : UG_THROW("InterpolateOnElem: Cannot find meaningful"
230 : " local positions of dofs.");
231 :
232 : // create a reference mapping
233 : ReferenceMapping<ref_elem_type, domain_type::dim> mapping;
234 :
235 : // iterate over all elements
236 0 : for( ; iter != iterEnd; ++iter)
237 : {
238 : // get element
239 : TElem* elem = *iter;
240 :
241 : // get all corner coordinates
242 : std::vector<position_type> vCorner;
243 0 : CollectCornerCoordinates(vCorner, *elem, *spGridFct->domain());
244 :
245 : // update the reference mapping for the corners
246 0 : mapping.update(&vCorner[0]);
247 :
248 : // get multiindices of element
249 : std::vector<DoFIndex> ind;
250 : spGridFct->dof_indices(elem, fct, ind);
251 :
252 : // check multi indices
253 0 : if(ind.size() != nsh)
254 0 : UG_THROW("InterpolateOnElem: On subset "<<si<<": Number of shapes is "
255 : <<nsh<<", but got "<<ind.size()<<" multi indices.");
256 :
257 : // loop all dofs
258 0 : for(size_t i = 0; i < nsh; ++i)
259 : {
260 : // global position
261 : position_type glob_pos;
262 :
263 :
264 : // map local dof position to global position
265 0 : mapping.local_to_global(glob_pos, loc_pos[i]);
266 :
267 : position_type rel_pos=glob_pos;
268 : rel_pos -=diff_pos;
269 :
270 : // value at position
271 : number val;
272 0 : (*spInterpolFunction)(val, rel_pos, time, si);
273 :
274 : // set value
275 0 : DoFRef(*spGridFct, ind[i]) = val;
276 : }
277 : }
278 0 : }
279 :
280 : /**
281 : * This function interpolates a grid function on an element by element loop. On
282 : * each element the all associated (up to the boundary of the element) are
283 : * interpolated and the values are stored in the grid function.
284 : *
285 : * @param[in] spInterpolFunction data providing interpolation values
286 : * @param[out] spGridFct interpolated grid function
287 : * @param[in] fct symbolic name of function component
288 : * @param[in] ssGrp subsets, where to interpolate
289 : * @param[in] time time point
290 : */
291 : template <typename TGridFunction>
292 0 : void InterpolateOnDiffElements(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
293 : SmartPtr<TGridFunction> spGridFct,
294 : size_t fct,
295 : number time,
296 : const SubsetGroup& ssGrp,const MathVector<TGridFunction::dim> diff_pos)
297 : {
298 : // loop subsets
299 0 : for(size_t i = 0; i < ssGrp.size(); ++i)
300 : {
301 : // get subset index
302 : const int si = ssGrp[i];
303 :
304 : // skip if function is not defined in subset
305 0 : if(!spGridFct->is_def_in_subset(fct, si)) continue;
306 :
307 : // switch dimensions
308 : try
309 : {
310 0 : const int dim = ssGrp.dim(i);
311 :
312 0 : if(dim > TGridFunction::dim)
313 0 : UG_THROW("InterpolateOnElements: Dimension of subset is "<<dim<<", but "
314 : " World Dimension is "<<TGridFunction::dim<<". Cannot interpolate this.");
315 :
316 : // FIXME (at least for Lagrange, order > 1, parallel)
317 : // In a parallel scenario, the distribution CAN cause elements of of lower
318 : // dimension than the rest of their subset to be located disconnected from
319 : // the rest of the subset on a processor. For example, in 2D, think of a
320 : // (1D) boundary subset and a distribution where the boundary of a proc's
321 : // domain only touches the boundary subset in a vertex, but intersects with
322 : // the boundary subset in another place.
323 : // This vertex will not be considered during interpolation even though it
324 : // should be!
325 0 : switch(dim)
326 : {
327 : case DIM_SUBSET_EMPTY_GRID: break;
328 : case 0: /* \TODO: do nothing may be wrong */ break;
329 : case 1:
330 0 : InterpolateOnDiffElements<RegularEdge, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time,diff_pos);
331 0 : break;
332 : case 2:
333 0 : InterpolateOnDiffElements<Triangle, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
334 0 : InterpolateOnDiffElements<Quadrilateral, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
335 0 : break;
336 : case 3:
337 0 : InterpolateOnDiffElements<Tetrahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
338 0 : InterpolateOnDiffElements<Hexahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
339 0 : InterpolateOnDiffElements<Prism, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
340 0 : InterpolateOnDiffElements<Pyramid, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
341 0 : InterpolateOnDiffElements<Octahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
342 0 : break;
343 0 : default: UG_THROW("InterpolateOnElements: Dimension " <<dim<<
344 : " not possible for world dim "<<3<<".");
345 : }
346 : }
347 0 : UG_CATCH_THROW("InterpolateOnElements: Cannot interpolate on elements.");
348 : }
349 0 : }
350 :
351 : /**
352 : * This function interpolates a grid function on an element by element loop. On
353 : * each element the all associated (up to the boundary of the element) are
354 : * interpolated and the values are stored in the grid function.
355 : *
356 : * @param[in] spInterpolFunction data providing interpolation values
357 : * @param[out] spGridFct interpolated grid function
358 : * @param[in] fct symbolic name of function component
359 : * @param[in] si subset, where to interpolate
360 : * @param[in] time time point
361 : */
362 : template <typename TElem, typename TGridFunction>
363 : void InterpolateOnElements(
364 : SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
365 : SmartPtr<TGridFunction> spGridFct, size_t fct, int si, number time)
366 : {
367 : // dimension of reference element
368 : const int dim = TGridFunction::dim;
369 :
370 : MathVector<dim> diff_pos(0.0);
371 : InterpolateOnDiffElements<TElem,TGridFunction>(spInterpolFunction, spGridFct, fct, si,time, diff_pos);
372 : }
373 :
374 : /**
375 : * This function interpolates a grid function on an element by element loop. On
376 : * each element the all associated (up to the boundary of the element) are
377 : * interpolated and the values are stored in the grid function.
378 : *
379 : * @param[in] spInterpolFunction data providing interpolation values
380 : * @param[out] spGridFct interpolated grid function
381 : * @param[in] fct symbolic name of function component
382 : * @param[in] ssGrp subsets, where to interpolate
383 : * @param[in] time time point
384 : */
385 : template <typename TGridFunction>
386 : void InterpolateOnElements(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
387 : SmartPtr<TGridFunction> spGridFct,
388 : size_t fct,
389 : number time,
390 : const SubsetGroup& ssGrp)
391 : {
392 : // loop subsets
393 : for(size_t i = 0; i < ssGrp.size(); ++i)
394 : {
395 : // get subset index
396 : const int si = ssGrp[i];
397 :
398 : // skip if function is not defined in subset
399 : if(!spGridFct->is_def_in_subset(fct, si)) continue;
400 :
401 : // switch dimensions
402 : try
403 : {
404 : const int dim = ssGrp.dim(i);
405 :
406 : if(dim > TGridFunction::dim)
407 : UG_THROW("InterpolateOnElements: Dimension of subset is "<<dim<<", but "
408 : " World Dimension is "<<TGridFunction::dim<<". Cannot interpolate this.");
409 :
410 : // FIXME (at least for Lagrange, order > 1, parallel)
411 : // In a parallel scenario, the distribution CAN cause elements of of lower
412 : // dimension than the rest of their subset to be located disconnected from
413 : // the rest of the subset on a processor. For example, in 2D, think of a
414 : // (1D) boundary subset and a distribution where the boundary of a proc's
415 : // domain only touches the boundary subset in a vertex, but intersects with
416 : // the boundary subset in another place.
417 : // This vertex will not be considered during interpolation even though it
418 : // should be!
419 : switch(dim)
420 : {
421 : case DIM_SUBSET_EMPTY_GRID: break;
422 : case 0: /* \TODO: do nothing may be wrong */ break;
423 : case 1:
424 : InterpolateOnElements<RegularEdge, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
425 : break;
426 : case 2:
427 : InterpolateOnElements<Triangle, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
428 : InterpolateOnElements<Quadrilateral, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
429 : break;
430 : case 3:
431 : InterpolateOnElements<Tetrahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
432 : InterpolateOnElements<Hexahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
433 : InterpolateOnElements<Prism, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
434 : InterpolateOnElements<Pyramid, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
435 : InterpolateOnElements<Octahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
436 : break;
437 : default: UG_THROW("InterpolateOnElements: Dimension " <<dim<<
438 : " not possible for world dim "<<3<<".");
439 : }
440 : }
441 : UG_CATCH_THROW("InterpolateOnElements: Cannot interpolate on elements.");
442 : }
443 : }
444 :
445 : ////////////////////////////////////////////////////////////////////////////////
446 : // Interpolate routine
447 : ////////////////////////////////////////////////////////////////////////////////
448 :
449 : /// interpolates a function on a subset
450 : /**
451 : * This function interpolates a GridFunction. To evaluate the function on every
452 : * point a functor must be passed.
453 : *
454 : * @param[in] spInterpolFunction data providing interpolation values
455 : * @param[out] spGridFct interpolated grid function
456 : * @param[in] cmp symbolic name of function component
457 : * @param[in] subsets subsets, where to interpolate (NULL = everywhere)
458 : * @param[in] time time point
459 : */
460 : template <typename TGridFunction>
461 0 : void Interpolate(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
462 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
463 : const char* subsets, number time)
464 : {
465 : // dimension of reference element
466 : const int dim = TGridFunction::dim;
467 : MathVector<dim> diff_pos(0.0);
468 0 : Interpolate(spInterpolFunction, spGridFct, cmp, subsets, time, diff_pos);
469 0 : }
470 :
471 : /// interpolates a function on a subset
472 : /**
473 : * This function interpolates a GridFunction. To evaluate the function on every
474 : * point a functor must be passed.
475 : *
476 : * @param[in] spInterpolFunction data providing interpolation values
477 : * @param[out] spGridFct interpolated grid function
478 : * @param[in] fct id of function component
479 : * @param[in] subsets subsets, where to interpolate (NULL = everywhere)
480 : * @param[in] time time point
481 : */
482 : template <typename TGridFunction>
483 0 : void Interpolate(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
484 : SmartPtr<TGridFunction> spGridFct, const size_t fct,
485 : const SubsetGroup& ssGrp, number time, const MathVector<TGridFunction::dim> diff_pos)
486 : {
487 :
488 : // check, that values do not depend on a solution
489 0 : if(spInterpolFunction->requires_grid_fct())
490 0 : UG_THROW("Interpolate: The interpolation values depend on a grid function."
491 : " This is not allowed in the current implementation. Use constant,"
492 : " lua-callback or vrl-callback user data only (even within linkers).");
493 :
494 : // check if fast P1 interpolation can be used
495 : // \TODO: This should be improved. Manifold admissible if space continuous
496 : bool bUseP1Interpolation = false;
497 0 : if(spGridFct->local_finite_element_id(fct).type() == LFEID::LAGRANGE &&
498 : spGridFct->local_finite_element_id(fct).order() == 1)
499 : bUseP1Interpolation = true;
500 :
501 : //forward
502 : if(bUseP1Interpolation){
503 0 : InterpolateOnDiffVertices<TGridFunction>(spInterpolFunction, spGridFct, fct, time, ssGrp, diff_pos);
504 : }else{
505 0 : InterpolateOnDiffElements<TGridFunction>(spInterpolFunction, spGridFct, fct, time, ssGrp, diff_pos);
506 : }
507 : // adjust parallel storage state
508 : #ifdef UG_PARALLEL
509 : spGridFct->set_storage_type(PST_CONSISTENT);
510 : #endif
511 0 : }
512 :
513 : /// interpolates a function on a subset
514 : /**
515 : * This function interpolates a GridFunction. To evaluate the function on every
516 : * point a functor must be passed.
517 : *
518 : * @param[in] spInterpolFunction data providing interpolation values
519 : * @param[out] spGridFct interpolated grid function
520 : * @param[in] cmp symbolic name of function component
521 : * @param[in] subsets subsets, where to interpolate (NULL = everywhere)
522 : * @param[in] time time point
523 : */
524 : template <typename TGridFunction>
525 0 : void Interpolate(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
526 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
527 : const char* subsets, number time, const MathVector<TGridFunction::dim> diff_pos)
528 : {
529 :
530 :
531 : // get function id of name
532 : const size_t fct = spGridFct->fct_id_by_name(cmp);
533 :
534 : // check that function found
535 0 : if(fct > spGridFct->num_fct())
536 0 : UG_THROW("Interpolate: Name of component '"<<cmp<<"' not found.");
537 :
538 0 : Interpolate(spInterpolFunction, spGridFct, fct, subsets, time, diff_pos);
539 0 : }
540 :
541 :
542 : /// interpolates a function on a subset
543 : /**
544 : * This function interpolates a GridFunction. To evaluate the function on every
545 : * point a functor must be passed.
546 : *
547 : * @param[in] spInterpolFunction data providing interpolation values
548 : * @param[out] spGridFct interpolated grid function
549 : * @param[in] cmp symbolic name of function component
550 : * @param[in] subsets subsets, where to interpolate (NULL = everywhere)
551 : * @param[in] time time point
552 : */
553 : template <typename TGridFunction>
554 0 : void Interpolate(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
555 : SmartPtr<TGridFunction> spGridFct, const size_t fct,
556 : const char* subsets, number time, const MathVector<TGridFunction::dim> diff_pos)
557 : {
558 : const bool bAllowManyfoldInterpolation =
559 : (spGridFct->local_finite_element_id(fct).type() == LFEID::LAGRANGE);
560 :
561 : // create subset group
562 0 : SubsetGroup ssGrp(spGridFct->domain()->subset_handler());
563 0 : if(subsets != NULL)
564 : {
565 0 : ssGrp.add(TokenizeString(subsets));
566 0 : if(!bAllowManyfoldInterpolation)
567 0 : if(!SameDimensionsInAllSubsets(ssGrp))
568 0 : UG_THROW("Interpolate: Subsets '"<<subsets<<"' do not have same dimension."
569 : "Can not integrate on subsets of different dimensions.");
570 : }
571 : else
572 : {
573 : // add all subsets and remove lower dim subsets afterwards
574 0 : ssGrp.add_all();
575 0 : if(!bAllowManyfoldInterpolation)
576 0 : RemoveLowerDimSubsets(ssGrp);
577 : }
578 :
579 0 : Interpolate(spInterpolFunction, spGridFct, fct, ssGrp, time, diff_pos);
580 0 : }
581 :
582 :
583 : template <typename TGridFunction>
584 0 : void Interpolate(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
585 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
586 : number time)
587 : {
588 0 : Interpolate(spInterpolFunction, spGridFct, cmp, NULL, time);
589 0 : }
590 : template <typename TGridFunction>
591 0 : void Interpolate(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
592 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
593 : const char* subsets)
594 : {
595 0 : Interpolate(spInterpolFunction, spGridFct, cmp, subsets, 0.0);
596 0 : }
597 : template <typename TGridFunction>
598 0 : void Interpolate(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
599 : SmartPtr<TGridFunction> spGridFct, const char* cmp)
600 : {
601 0 : Interpolate(spInterpolFunction, spGridFct, cmp, NULL, 0.0);
602 0 : }
603 : //interpolate with diff_vector
604 : template <typename TGridFunction>
605 0 : void Interpolate(SmartPtr<UserData<number, TGridFunction::dim> > spInterpolFunction,
606 : SmartPtr<TGridFunction> spGridFct, const char* cmp, const MathVector<TGridFunction::dim>& m_diff_pos)
607 : {
608 0 : Interpolate(spInterpolFunction, spGridFct, cmp, NULL, 0.0, m_diff_pos);
609 0 : }
610 :
611 : ///////////////
612 : // const data
613 : ///////////////
614 :
615 : template <typename TGridFunction>
616 0 : void Interpolate(number val,
617 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
618 : const char* subsets, number time)
619 : {
620 : static const int dim = TGridFunction::dim;
621 0 : SmartPtr<UserData<number, dim> > sp =
622 0 : make_sp(new ConstUserNumber<dim>(val));
623 0 : Interpolate(sp, spGridFct, cmp, subsets, time);
624 0 : }
625 : template <typename TGridFunction>
626 0 : void Interpolate(number val,
627 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
628 : number time)
629 0 : {Interpolate(val, spGridFct, cmp, NULL, time);}
630 : template <typename TGridFunction>
631 0 : void Interpolate(number val,
632 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
633 : const char* subsets)
634 0 : {Interpolate(val, spGridFct, cmp, subsets, 0.0);}
635 : template <typename TGridFunction>
636 0 : void Interpolate(number val,
637 : SmartPtr<TGridFunction> spGridFct, const char* cmp)
638 0 : {Interpolate(val, spGridFct, cmp, NULL, 0.0);}
639 :
640 : //interpolate with diff_vector
641 : template <typename TGridFunction>
642 : void Interpolate(number val,
643 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
644 : const char* subsets, number time,const SmartPtr<CplUserData<MathVector<TGridFunction::dim>, TGridFunction::dim> > m_diff_pos)
645 : {
646 : static const int dim = TGridFunction::dim;
647 : SmartPtr<UserData<number, dim> > sp =
648 : make_sp(new ConstUserNumber<dim>(val));
649 :
650 : InterpolateDiff(sp, spGridFct, cmp, subsets, time,m_diff_pos);
651 : }
652 :
653 :
654 :
655 : ///////////////
656 : // lua data
657 : ///////////////
658 :
659 : #ifdef UG_FOR_LUA
660 : // function-name as string
661 : template <typename TGridFunction>
662 0 : void Interpolate(const char* LuaFunction,
663 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
664 : const char* subsets, number time)
665 : {
666 : static const int dim = TGridFunction::dim;
667 0 : SmartPtr<UserData<number, dim> > sp =
668 : LuaUserDataFactory<number, dim>::create(LuaFunction);
669 0 : Interpolate(sp, spGridFct, cmp, subsets, time);
670 0 : }
671 : template <typename TGridFunction>
672 0 : void Interpolate(const char* LuaFunction,
673 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
674 : number time)
675 0 : {Interpolate(LuaFunction, spGridFct, cmp, NULL, time);}
676 : template <typename TGridFunction>
677 0 : void Interpolate(const char* LuaFunction,
678 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
679 : const char* subsets)
680 0 : {Interpolate(LuaFunction, spGridFct, cmp, subsets, 0.0);}
681 : template <typename TGridFunction>
682 0 : void Interpolate(const char* LuaFunction,
683 : SmartPtr<TGridFunction> spGridFct, const char* cmp)
684 0 : {Interpolate(LuaFunction, spGridFct, cmp, NULL, 0.0);}
685 :
686 : // function as function handle
687 : template <typename TGridFunction>
688 0 : void Interpolate(LuaFunctionHandle LuaFunction,
689 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
690 : const char* subsets, number time)
691 : {
692 : static const int dim = TGridFunction::dim;
693 0 : SmartPtr<UserData<number, dim> > sp =
694 0 : make_sp(new LuaUserData<number, dim>(LuaFunction));
695 0 : Interpolate(sp, spGridFct, cmp, subsets, time);
696 0 : }
697 : template <typename TGridFunction>
698 0 : void Interpolate(LuaFunctionHandle LuaFunction,
699 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
700 : number time)
701 0 : {Interpolate(LuaFunction, spGridFct, cmp, NULL, time);}
702 : template <typename TGridFunction>
703 0 : void Interpolate(LuaFunctionHandle LuaFunction,
704 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
705 : const char* subsets)
706 0 : {Interpolate(LuaFunction, spGridFct, cmp, subsets, 0.0);}
707 : template <typename TGridFunction>
708 0 : void Interpolate(LuaFunctionHandle LuaFunction,
709 : SmartPtr<TGridFunction> spGridFct, const char* cmp)
710 0 : {Interpolate(LuaFunction, spGridFct, cmp, NULL, 0.0);}
711 :
712 : //interpolate with Diff-vector:
713 : template <typename TGridFunction>
714 : void InterpolateDiff(const char* LuaFunction,
715 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
716 : const char* subsets, number time, SmartPtr<CplUserData<MathVector<TGridFunction::dim>, TGridFunction::dim> > m_diff_pos )
717 : {
718 : static const int dim = TGridFunction::dim;
719 : SmartPtr<UserData<number, dim> > sp =
720 : LuaUserDataFactory<number, dim>::create(LuaFunction);
721 : InterpolateDiff(sp, spGridFct, cmp,subsets, time, m_diff_pos);
722 : }
723 :
724 : template <typename TGridFunction>
725 : void InterpolateDiff(LuaFunctionHandle LuaFunction,
726 : SmartPtr<TGridFunction> spGridFct, const char* cmp,
727 : const char* subsets, number time,SmartPtr<CplUserData<MathVector<TGridFunction::dim>, TGridFunction::dim> > m_diff_pos)
728 : {
729 : static const int dim = TGridFunction::dim;
730 : SmartPtr<UserData<number, dim> > sp =
731 : make_sp(new LuaUserData<number, dim>(LuaFunction));
732 : InterpolateDiff(sp, spGridFct, cmp, subsets, time, m_diff_pos);
733 : }
734 : template <typename TGridFunction>
735 : void Interpolate(LuaFunctionHandle LuaFunction,
736 : SmartPtr<TGridFunction> spGridFct, const char* cmp,const SmartPtr<CplUserData<MathVector<TGridFunction::dim>, TGridFunction::dim> > m_diff_pos)
737 : {InterpolateDiff(LuaFunction, spGridFct, cmp, NULL, 0.0, m_diff_pos);}
738 :
739 : #endif
740 :
741 :
742 :
743 :
744 : } // namespace ug
745 :
746 : #endif /*__H__UG__LIB_DISC__FUNCTION_SPACES__INTERPOLATE__*/
|