33#ifndef __H__UG__LIB_DISC__FUNCTION_SPACES__INTERPOLATE__
34#define __H__UG__LIB_DISC__FUNCTION_SPACES__INTERPOLATE__
71template <
typename TGr
idFunction>
81 typedef typename TGridFunction::domain_type domain_type;
82 typedef typename domain_type::position_type position_type;
87 const typename domain_type::position_accessor_type& aaPos
88 = spGridFct->domain()->position_accessor();
91 std::vector<DoFIndex> ind;
92 typename TGridFunction::template dim_traits<0>::const_iterator iterEnd, iter;
94 for(
size_t i = 0; i < ssGrp.
size(); ++i)
97 const int si = ssGrp[i];
100 if(!spGridFct->is_def_in_subset(fct, si))
continue;
103 iterEnd = spGridFct->template end<Vertex>(si);
104 iter = spGridFct->template begin<Vertex>(si);
105 for(; iter != iterEnd; ++iter)
111 position_type glob_pos = aaPos[vrt];
112 position_type rel_pos=glob_pos;
117 (*spInterpolFunction)(val, rel_pos, time, si);
120 spGridFct->dof_indices(vrt, fct, ind);
123 for(
size_t i = 0; i < ind.size(); ++i)
126 DoFRef(*spGridFct, ind[i]) = val;
134template <
typename TGr
idFunction>
136 typename TGridFunction::domain_type::position_type pos,
141 (*spInterpolFunction)(val, pos, time, si);
158template <
typename TGr
idFunction>
166 const int dim = TGridFunction::dim;
169 InterpolateOnDiffVertices<TGridFunction>(spInterpolFunction, spGridFct, fct, time, ssGrp, diff_pos);
189template <
typename TElem,
typename TGr
idFunction>
200 const int dim = ref_elem_type::dim;
204 typedef typename TGridFunction::domain_type domain_type;
205 typedef typename domain_type::position_type position_type;
208 typename TGridFunction::template traits<TElem>::const_iterator iterEnd, iter;
209 iterEnd = spGridFct->template end<TElem>(si);
210 iter = spGridFct->template begin<TElem>(si);
213 if(iter == iterEnd)
return;
216 LFEID id = spGridFct->local_finite_element_id(fct);
220 LocalFiniteElementProvider::get<dim>(roid,
id);
223 const size_t nsh = trialSpace.
num_sh();
226 std::vector<MathVector<dim> > loc_pos(nsh);
227 for(
size_t i = 0; i < nsh; ++i)
228 if(!trialSpace.
position(i, loc_pos[i]))
229 UG_THROW(
"InterpolateOnElem: Cannot find meaningful"
230 " local positions of dofs.");
236 for( ; iter != iterEnd; ++iter)
242 std::vector<position_type> vCorner;
246 mapping.
update(&vCorner[0]);
249 std::vector<DoFIndex> ind;
250 spGridFct->dof_indices(elem, fct, ind);
253 if(ind.size() != nsh)
254 UG_THROW(
"InterpolateOnElem: On subset "<<si<<
": Number of shapes is "
255 <<nsh<<
", but got "<<ind.size()<<
" multi indices.");
258 for(
size_t i = 0; i < nsh; ++i)
261 position_type glob_pos;
267 position_type rel_pos=glob_pos;
272 (*spInterpolFunction)(val, rel_pos, time, si);
275 DoFRef(*spGridFct, ind[i]) = val;
291template <
typename TGr
idFunction>
299 for(
size_t i = 0; i < ssGrp.
size(); ++i)
302 const int si = ssGrp[i];
305 if(!spGridFct->is_def_in_subset(fct, si))
continue;
310 const int dim = ssGrp.
dim(i);
312 if(dim > TGridFunction::dim)
313 UG_THROW(
"InterpolateOnElements: Dimension of subset is "<<dim<<
", but "
314 " World Dimension is "<<TGridFunction::dim<<
". Cannot interpolate this.");
330 InterpolateOnDiffElements<RegularEdge, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time,diff_pos);
333 InterpolateOnDiffElements<Triangle, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
334 InterpolateOnDiffElements<Quadrilateral, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
337 InterpolateOnDiffElements<Tetrahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
338 InterpolateOnDiffElements<Hexahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
339 InterpolateOnDiffElements<Prism, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
340 InterpolateOnDiffElements<Pyramid, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
341 InterpolateOnDiffElements<Octahedron, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time, diff_pos);
343 default:
UG_THROW(
"InterpolateOnElements: Dimension " <<dim<<
344 " not possible for world dim "<<3<<
".");
347 UG_CATCH_THROW(
"InterpolateOnElements: Cannot interpolate on elements.");
362template <
typename TElem,
typename TGr
idFunction>
368 const int dim = TGridFunction::dim;
371 InterpolateOnDiffElements<TElem,TGridFunction>(spInterpolFunction, spGridFct, fct, si,time, diff_pos);
385template <
typename TGr
idFunction>
393 for(
size_t i = 0; i < ssGrp.
size(); ++i)
396 const int si = ssGrp[i];
399 if(!spGridFct->is_def_in_subset(fct, si))
continue;
404 const int dim = ssGrp.
dim(i);
406 if(dim > TGridFunction::dim)
407 UG_THROW(
"InterpolateOnElements: Dimension of subset is "<<dim<<
", but "
408 " World Dimension is "<<TGridFunction::dim<<
". Cannot interpolate this.");
424 InterpolateOnElements<RegularEdge, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
427 InterpolateOnElements<Triangle, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
428 InterpolateOnElements<Quadrilateral, TGridFunction>(spInterpolFunction, spGridFct, fct, si, time);
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);
437 default:
UG_THROW(
"InterpolateOnElements: Dimension " <<dim<<
438 " not possible for world dim "<<3<<
".");
441 UG_CATCH_THROW(
"InterpolateOnElements: Cannot interpolate on elements.");
460template <
typename TGr
idFunction>
463 const char* subsets,
number time)
466 const int dim = TGridFunction::dim;
468 Interpolate(spInterpolFunction, spGridFct, cmp, subsets, time, diff_pos);
482template <
typename TGr
idFunction>
489 if(spInterpolFunction->requires_grid_fct())
490 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).");
496 bool bUseP1Interpolation =
false;
497 if(spGridFct->local_finite_element_id(fct).type() ==
LFEID::LAGRANGE &&
498 spGridFct->local_finite_element_id(fct).order() == 1)
499 bUseP1Interpolation =
true;
502 if(bUseP1Interpolation){
503 InterpolateOnDiffVertices<TGridFunction>(spInterpolFunction, spGridFct, fct, time, ssGrp, diff_pos);
505 InterpolateOnDiffElements<TGridFunction>(spInterpolFunction, spGridFct, fct, time, ssGrp, diff_pos);
524template <
typename TGr
idFunction>
532 const size_t fct = spGridFct->fct_id_by_name(cmp);
535 if(fct > spGridFct->num_fct())
536 UG_THROW(
"Interpolate: Name of component '"<<cmp<<
"' not found.");
538 Interpolate(spInterpolFunction, spGridFct, fct, subsets, time, diff_pos);
553template <
typename TGr
idFunction>
558 const bool bAllowManyfoldInterpolation =
562 SubsetGroup ssGrp(spGridFct->domain()->subset_handler());
566 if(!bAllowManyfoldInterpolation)
568 UG_THROW(
"Interpolate: Subsets '"<<subsets<<
"' do not have same dimension."
569 "Can not integrate on subsets of different dimensions.");
575 if(!bAllowManyfoldInterpolation)
579 Interpolate(spInterpolFunction, spGridFct, fct, ssGrp, time, diff_pos);
583template <
typename TGr
idFunction>
588 Interpolate(spInterpolFunction, spGridFct, cmp, NULL, time);
590template <
typename TGr
idFunction>
595 Interpolate(spInterpolFunction, spGridFct, cmp, subsets, 0.0);
597template <
typename TGr
idFunction>
601 Interpolate(spInterpolFunction, spGridFct, cmp, NULL, 0.0);
604template <
typename TGr
idFunction>
608 Interpolate(spInterpolFunction, spGridFct, cmp, NULL, 0.0, m_diff_pos);
615template <
typename TGr
idFunction>
618 const char* subsets,
number time)
620 static const int dim = TGridFunction::dim;
625template <
typename TGr
idFunction>
630template <
typename TGr
idFunction>
635template <
typename TGr
idFunction>
641template <
typename TGr
idFunction>
646 static const int dim = TGridFunction::dim;
650 InterpolateDiff(sp, spGridFct, cmp, subsets, time,m_diff_pos);
661template <
typename TGr
idFunction>
664 const char* subsets,
number time)
666 static const int dim = TGridFunction::dim;
671template <
typename TGr
idFunction>
675{
Interpolate(LuaFunction, spGridFct, cmp, NULL, time);}
676template <
typename TGr
idFunction>
680{
Interpolate(LuaFunction, spGridFct, cmp, subsets, 0.0);}
681template <
typename TGr
idFunction>
684{
Interpolate(LuaFunction, spGridFct, cmp, NULL, 0.0);}
687template <
typename TGr
idFunction>
690 const char* subsets,
number time)
692 static const int dim = TGridFunction::dim;
694 make_sp(
new LuaUserData<number, dim>(LuaFunction));
697template <
typename TGr
idFunction>
701{
Interpolate(LuaFunction, spGridFct, cmp, NULL, time);}
702template <
typename TGr
idFunction>
706{
Interpolate(LuaFunction, spGridFct, cmp, subsets, 0.0);}
707template <
typename TGr
idFunction>
710{
Interpolate(LuaFunction, spGridFct, cmp, NULL, 0.0);}
713template <
typename TGr
idFunction>
714void InterpolateDiff(
const char* LuaFunction,
716 const char* subsets,
number time,
SmartPtr<CplUserData<MathVector<TGridFunction::dim>, TGridFunction::dim> > m_diff_pos )
718 static const int dim = TGridFunction::dim;
721 InterpolateDiff(sp, spGridFct, cmp,subsets, time, m_diff_pos);
724template <
typename TGr
idFunction>
725void InterpolateDiff(LuaFunctionHandle LuaFunction,
727 const char* subsets,
number time,
SmartPtr<CplUserData<MathVector<TGridFunction::dim>, TGridFunction::dim> > m_diff_pos)
729 static const int dim = TGridFunction::dim;
731 make_sp(
new LuaUserData<number, dim>(LuaFunction));
732 InterpolateDiff(sp, spGridFct, cmp, subsets, time, m_diff_pos);
734template <
typename TGr
idFunction>
737{InterpolateDiff(LuaFunction, spGridFct, cmp, NULL, 0.0, m_diff_pos);}
Definition smart_pointer.h:108
constant scalar user data
Definition const_user_data.h:153
Type based UserData.
Definition user_data.h:501
virtual bool position(size_t i, MathVector< TDim > &pos) const =0
local position of DoF i
Identifier for Local Finite Elements.
Definition local_finite_element_id.h:98
@ LAGRANGE
Definition local_finite_element_id.h:104
virtual size_t num_sh() const
Definition local_dof_set.cpp:46
virtual base class for local shape function sets
Definition local_shape_function_set.h:70
static SmartPtr< LuaUserData< TData, dim, TRet > > create(const std::string &name)
Definition lua_user_data.h:223
a mathematical Vector with N entries.
Definition math_vector.h:97
Definition reference_mapping.h:65
void local_to_global(MathVector< worldDim > &globPos, const MathVector< dim > &locPos) const
map local coordinate to global coordinate
void update(const MathVector< worldDim > *vCornerCoord)
refresh mapping for new set of corners
Group of subsets.
Definition subset_group.h:51
size_t size() const
number of subsets in this group
Definition subset_group.h:122
void add_all()
select all subsets of underlying subset handler
Definition subset_group.cpp:133
int dim(size_t i) const
dimension of subset
Definition subset_group.cpp:237
void add(int si)
adds a subset by number to this group
Definition subset_group.cpp:64
Type based UserData.
Definition user_data.h:143
Base-class for all vertex-types.
Definition grid_base_objects.h:231
@ PST_CONSISTENT
Definition parallel_storage_type.h:68
void CollectCornerCoordinates(std::vector< typename TAAPos::ValueType > &vCornerCoordsOut, const TElem &elem, const TAAPos &aaPos, bool clearContainer=true)
returns the corner coordinates of a geometric object
Definition domain_util_impl.h:75
#define UG_CATCH_THROW(msg)
Definition error.h:64
#define UG_THROW(msg)
Definition error.h:57
double number
Definition types.h:124
void InterpolateOnDiffVertices(SmartPtr< UserData< number, TGridFunction::dim > > spInterpolFunction, SmartPtr< TGridFunction > spGridFct, size_t fct, number time, const SubsetGroup &ssGrp, const MathVector< TGridFunction::dim > diff_pos)
Definition interpolate.h:72
ReferenceObjectID
these ids are used to identify the shape of a geometric object.
Definition grid_base_objects.h:74
void InterpolateOnElements(SmartPtr< UserData< number, TGridFunction::dim > > spInterpolFunction, SmartPtr< TGridFunction > spGridFct, size_t fct, int si, number time)
Definition interpolate.h:363
void Interpolate(SmartPtr< UserData< number, TGridFunction::dim > > spInterpolFunction, SmartPtr< TGridFunction > spGridFct, const char *cmp, const char *subsets, number time)
interpolates a function on a subset
Definition interpolate.h:461
void TokenizeString(const string &str, vector< string > &vToken, const char delimiter)
Definition string_util.cpp:56
bool SameDimensionsInAllSubsets(const SubsetGroup &subsetGroup)
Definition subset_group.cpp:299
void InterpolateOnVertices(SmartPtr< UserData< number, TGridFunction::dim > > spInterpolFunction, SmartPtr< TGridFunction > spGridFct, size_t fct, number time, const SubsetGroup &ssGrp)
Definition interpolate.h:159
void RemoveLowerDimSubsets(SubsetGroup &subsetGroup)
Definition subset_group.cpp:315
void InterpolateOnDiffElements(SmartPtr< UserData< number, TGridFunction::dim > > spInterpolFunction, SmartPtr< TGridFunction > spGridFct, size_t fct, int si, number time, const MathVector< TGridFunction::dim > diff_pos)
Definition interpolate.h:190
number get_number_on_coords(SmartPtr< UserData< number, TGridFunction::dim > > spInterpolFunction, typename TGridFunction::domain_type::position_type pos, number time, const int si)
Definition interpolate.h:135
number & DoFRef(TMatrix &mat, const DoFIndex &iInd, const DoFIndex &jInd)
Definition multi_index.h:276
@ DIM_SUBSET_EMPTY_GRID
Definition subset_dim_util.h:92
SmartPtr< T, FreePolicy > make_sp(T *inst)
returns a SmartPtr for the passed raw pointer
Definition smart_pointer.h:836
traits for reference elements
Definition reference_element_traits.h:48