30 #ifndef __quadptlocalization_H
31 #define __quadptlocalization_H 1
36 #define FEELPP_EXPORT_QUADLOCALIZATION 0
37 #if FEELPP_EXPORT_QUADLOCALIZATION
38 #include <feel/feelfilters/exporter.hpp>
46 template <
typename IteratorRange,
typename Im,
typename Expr>
47 class QuadPtLocalization
51 static const size_type iDim = boost::tuples::template element<0, IteratorRange>::type::value;
54 static const size_type context = mpl::if_< boost::is_same< mpl::int_<iDim>, mpl::int_<MESH_FACES> >,
55 mpl::int_<Expr::context|vm::JACOBIAN|vm::KB|vm::NORMAL|vm::POINT>,
56 mpl::int_<Expr::context|vm::POINT> >::type::value;
61 typedef IteratorRange range_iterator;
64 typedef typename boost::tuples::template element<1, IteratorRange>::type element_iterator_type;
65 typedef typename element_iterator_type::value_type::GeoShape GeoShape;
68 typedef typename mesh_type::gm_type gm_type;
69 typedef typename mesh_type::element_type geoelement_type;
75 typedef typename boost::remove_reference<typename element_iterator_type::reference>::type const_t;
76 typedef typename boost::remove_const<const_t>::type the_face_element_type;
77 typedef typename the_face_element_type::super2::template Element<the_face_element_type>::type the_element_type;
78 typedef the_element_type geoelement_type;
79 typedef typename geoelement_type::gm_type gm_type;
86 typedef typename gm_type::template Context<context, geoelement_type> gmc_type;
87 typedef boost::shared_ptr<gmc_type> gmc_ptrtype;
88 typedef typename gm_type::precompute_type pc_type;
89 typedef typename gm_type::precompute_ptrtype pc_ptrtype;
91 typedef typename mesh_type::value_type value_type;
92 typedef typename mesh_type::node_type node_type;
93 typedef typename matrix_node<value_type>::type matrix_node_type;
95 typedef typename geoelement_type::value_type value_type;
96 typedef typename geoelement_type::node_type node_type;
97 typedef typename geoelement_type::matrix_node_type matrix_node_type;
101 typedef typename mpl::if_<mpl::bool_<geoelement_type::is_simplex>,
102 mpl::identity<typename Im::template apply<geoelement_type::nDim, value_type, Simplex>::type >,
103 mpl::identity<typename Im::template apply<geoelement_type::nDim, value_type, Hypercube>::type >
104 >::type::type im_type;
106 typedef typename im_type::face_quadrature_type im_face_type;
108 typedef typename QuadMapped<im_type>::permutation_type permutation_type;
109 typedef typename QuadMapped<im_type>::permutation_points_type permutation_points_type;
114 typedef std::map< size_type,std::list<boost::tuple<size_type,size_type,node_type,node_type> > > result_temp1_type;
116 typedef std::map< size_type,boost::tuple< std::vector< boost::tuple<size_type,size_type > >,matrix_node_type,matrix_node_type> > result_temp2_type;
118 typedef std::list< boost::tuple< size_type, result_temp2_type> > result_container_type;
121 typedef std::list< boost::tuple<size_type,std::vector< boost::tuple<size_type,size_type > >,matrix_node_type > > result_container_linear_type;
126 QuadPtLocalization( IteratorRange
const& elts )
131 M_ppts( M_qm( this->im() ) ),
132 M_hasPrecompute(false)
134 M_listRange.push_back( elts );
137 QuadPtLocalization( std::list<IteratorRange>
const& elts )
142 M_ppts( M_qm( this->im() ) ),
143 M_hasPrecompute(false)
150 im_type
const& im()
const
158 im_face_type im( uint16_type f )
const
160 return M_im.face( f );
166 element_iterator_type beginElement()
const
168 return M_listRange.front().template get<1>();
174 element_iterator_type endElement()
const
176 return M_listRange.front().template get<2>();
182 result_container_linear_type
const & resultLinear()
const
190 result_container_type
const & result()
const
192 return M_resBilinear;
196 bool hasPrecompute()
const {
return M_hasPrecompute; }
203 return theIdq/this->im().nPoints();
207 return theIdq/this->im().nPointsOnFace();
212 return theIdq%this->im().nPoints();
216 return theIdq%this->im().nPointsOnFace();
221 gmc_ptrtype gmcForThisElt(
size_type theIdElt,
222 std::vector<boost::tuple<size_type,size_type> >
const& indexLocalToQuad,
223 mpl::int_<MESH_ELEMENTS> )
226 element_iterator_type elt_it = this->beginElement(), elt_en = this->beginElement();
227 auto itListRange = M_listRange.begin();
228 auto const enListRange = M_listRange.end();
230 for (
size_type ide = 0 ; itListRange!=enListRange && !findElt; ++itListRange)
232 boost::tie( boost::tuples::ignore, elt_it, elt_en ) = *itListRange;
233 const size_type distRange = std::distance(elt_it, elt_en);
234 if ( (ide+distRange-1) < theIdElt) { ide+=distRange;
continue; }
236 for ( ; ide<theIdElt ; ++ide ) ++elt_it;
241 auto elt_it = this->beginElement();
242 for (
size_type i=0; i<theIdElt; ++i ) ++elt_it;
245 uint16_type nContextPt = indexLocalToQuad.size();
246 matrix_node_type
const& quadPtsRef = this->im().points();
247 matrix_node_type newquadPtsRef( quadPtsRef.size1() , nContextPt );
249 for ( uint16_type i=0; i<nContextPt; ++i )
251 ublas::column( newquadPtsRef, i ) = ublas::column( quadPtsRef, indexLocalToQuad[i].
template get<0>() );
255 pc_ptrtype geopc(
new pc_type( elt_it->gm(), newquadPtsRef ) );
256 gmc_ptrtype gmc(
new gmc_type( elt_it->gm(),*elt_it, geopc ) );
263 gmc_ptrtype gmcForThisElt(
size_type theIdElt,
264 std::vector<boost::tuple<size_type,size_type> >
const& indexLocalToQuad,
265 mpl::int_<MESH_FACES> )
268 element_iterator_type elt_it=this->beginElement(), elt_en=this->beginElement();
269 auto itListRange = M_listRange.begin();
270 auto const enListRange = M_listRange.end();
272 for (
size_type ide = 0 ; itListRange!=enListRange && !findElt; ++itListRange)
274 boost::tie( boost::tuples::ignore, elt_it, elt_en ) = *itListRange;
275 const size_type distRange = std::distance(elt_it, elt_en);
276 if ( (ide+distRange-1) < theIdElt) { ide+=distRange;
continue; }
278 for ( ; ide<theIdElt ; ++ide ) ++elt_it;
283 auto elt_it = this->beginElement();
284 for (
size_type i=0; i<theIdElt; ++i ) ++elt_it;
287 const uint16_type nContextPt = indexLocalToQuad.size();
288 const uint16_type __face_id_in_elt_0 = elt_it->pos_first();
290 auto const __perm = elt_it->element( 0 ).permutation( __face_id_in_elt_0 );
291 matrix_node_type
const& quadPtsRef = M_ppts[ __face_id_in_elt_0].find( __perm )->second;
293 matrix_node_type newquadPtsRef( quadPtsRef.size1() , nContextPt );
295 for ( uint16_type i=0; i<nContextPt; ++i )
297 ublas::column( newquadPtsRef, i ) = ublas::column( quadPtsRef, indexLocalToQuad[i].
template get<0>() );
301 std::vector<std::map<permutation_type, pc_ptrtype> > __geopc( this->im().nFaces() );
303 for ( uint16_type __f = 0; __f < this->im().nFaces(); ++__f )
305 for ( permutation_type __p( permutation_type::IDENTITY );
306 __p < permutation_type( permutation_type::N_PERMUTATIONS ); ++__p )
308 __geopc[__f][__p] = pc_ptrtype(
new pc_type( elt_it->element( 0 ).gm(),
313 gmc_ptrtype gmc(
new gmc_type( elt_it->element( 0 ).gm(),
314 elt_it->element( 0 ),
316 __face_id_in_elt_0 ) );
323 std::vector< boost::tuple< std::vector<boost::tuple<size_type,size_type> >,gmc_ptrtype,matrix_node_type > >
324 getUsableDataInFormContext( std::vector<boost::tuple<size_type,size_type> >
const& indexLocalToQuad,
325 matrix_node_type
const & ptsRefTest )
328 uint16_type nContextPt = indexLocalToQuad.size();
329 std::map<size_type,std::list<size_type> > mapEltId;
331 for ( uint16_type i=0; i<nContextPt; ++i )
333 size_type eltId = this->eltForThisQuadPt( indexLocalToQuad[i].
template get<1>(),mpl::int_<iDim>() );
334 mapEltId[eltId].push_back( i );
338 auto nEltInContext = mapEltId.size();
341 std::vector< boost::tuple< std::vector<boost::tuple<size_type,size_type> >, gmc_ptrtype, matrix_node_type > > vec_res( nEltInContext );
343 auto map_it = mapEltId.begin();
344 auto map_en = mapEltId.end();
347 for ( ; map_it!= map_en ; ++map_it, ++cptId )
350 std::vector<boost::tuple<size_type,size_type> > newindexLocalToQuad( map_it->second.size() );
352 matrix_node_type newptsRefTest( ptsRefTest.size1(), map_it->second.size() );
353 auto sublist_it = map_it->second.begin();
354 auto sublist_en = map_it->second.end();
357 for ( ; sublist_it != sublist_en ; ++sublist_it,++index )
359 newindexLocalToQuad[index] = indexLocalToQuad[*sublist_it];
360 ublas::column( newptsRefTest, index ) = ublas::column( ptsRefTest,*sublist_it );
363 auto thegmc = gmcForThisElt( map_it->first,newindexLocalToQuad,mpl::int_<iDim>() );
365 vec_res[cptId] = boost::make_tuple( newindexLocalToQuad, thegmc, newptsRefTest );
373 std::vector< boost::tuple< std::vector<boost::tuple<size_type,size_type> >,gmc_ptrtype,matrix_node_type, matrix_node_type > >
374 getUsableDataInFormContext( std::vector<boost::tuple<size_type,size_type> >
const& indexLocalToQuad,
375 matrix_node_type
const & ptsRefTest,
376 matrix_node_type
const & ptsRefTrial )
379 uint16_type nContextPt = indexLocalToQuad.size();
380 std::map<size_type,std::list<size_type> > mapEltId;
382 for ( uint16_type i=0; i<nContextPt; ++i )
384 size_type eltId = this->eltForThisQuadPt( indexLocalToQuad[i].
template get<1>(),mpl::int_<iDim>() );
385 mapEltId[eltId].push_back( i );
389 auto nEltInContext = mapEltId.size();
392 std::vector< boost::tuple< std::vector<boost::tuple<size_type,size_type> >, gmc_ptrtype, matrix_node_type, matrix_node_type > > vec_res( nEltInContext );
394 auto map_it = mapEltId.begin();
395 auto map_en = mapEltId.end();
398 for ( ; map_it!= map_en ; ++map_it, ++cptId )
401 std::vector<boost::tuple<size_type,size_type> > newindexLocalToQuad( map_it->second.size() );
403 matrix_node_type newptsRefTest( ptsRefTest.size1(), map_it->second.size() );
404 matrix_node_type newptsRefTrial( ptsRefTrial.size1(), map_it->second.size() );
406 auto sublist_it = map_it->second.begin();
407 auto sublist_en = map_it->second.end();
410 for ( ; sublist_it != sublist_en ; ++sublist_it,++index )
412 newindexLocalToQuad[index] = indexLocalToQuad[*sublist_it];
413 ublas::column( newptsRefTest, index ) = ublas::column( ptsRefTest,*sublist_it );
414 ublas::column( newptsRefTrial, index ) = ublas::column( ptsRefTrial,*sublist_it );
418 auto thegmc = gmcForThisElt( map_it->first,newindexLocalToQuad,mpl::int_<iDim>() );
420 vec_res[cptId] = boost::make_tuple( newindexLocalToQuad, thegmc, newptsRefTest,newptsRefTrial );
428 template <
typename Mesh1Type>
430 localization( mpl::int_<MESH_FACES> ,
431 boost::shared_ptr<Mesh1Type> meshTest,
432 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > & testEltToPtsQuad )
435 auto begin_elt_it = this->beginElement();
437 std::vector<std::map<permutation_type, pc_ptrtype> > __geopc( this->im().nFaces() );
442 for ( uint16_type __f = 0; __f < this->im().nFaces(); ++__f )
446 for ( permutation_type __p( permutation_type::IDENTITY );
447 __p < permutation_type( permutation_type::N_PERMUTATIONS ); ++__p )
450 __geopc[__f][__p] = pc_ptrtype(
new pc_type( begin_elt_it->element( 0 ).gm(), M_ppts[__f].find( __p )->second ) );
455 uint16_type __face_id_in_elt_0 = begin_elt_it->pos_first();
457 gmc_ptrtype gmc(
new gmc_type( begin_elt_it->element( 0 ).gm(),
458 begin_elt_it->element( 0 ),
460 __face_id_in_elt_0 ) );
462 auto meshTestLocalization = meshTest->tool_localization();
463 meshTestLocalization->updateForUse();
464 const auto nbNearNeighborAtStartTest = meshTestLocalization->kdtree()->nPtMaxNearNeighbor();
465 bool notUseOptLocTest = Mesh1Type::nDim!=Mesh1Type::nRealDim;
466 if (notUseOptLocTest) meshTestLocalization->kdtree()->nbNearNeighbor(Mesh1Type::element_type::numPoints);
468 matrix_node_type ptsReal( begin_elt_it->vertices().size1(), 1 );
470 node_type testNodeRef;
472 element_iterator_type elt_it, elt_en;
474 auto itListRange = M_listRange.begin();
475 auto const enListRange = M_listRange.end();
476 for (
size_type ide = 0 ; itListRange!=enListRange ; ++itListRange)
478 boost::tie( boost::tuples::ignore, elt_it, elt_en ) = *itListRange;
479 for ( ; elt_it != elt_en; ++elt_it, ++ide )
481 __face_id_in_elt_0 = elt_it->pos_first();
484 gmc->update( elt_it->element( 0 ), __face_id_in_elt_0 );
487 for (
int q = 0; q < gmc->nPoints(); ++ q )
491 auto testAnalysis = meshTestLocalization->searchElement( gmc->xReal( q ) );
492 testIdElt = testAnalysis.template get<1>();
493 testNodeRef = testAnalysis.template get<2>();
495 ublas::column(ptsReal,0 ) = gmc->xReal( q );
497 auto resLocalisationTest = meshTestLocalization->run_analysis(ptsReal,testIdElt,elt_it->vertices(),mpl::int_<0>());
498 testIdElt = resLocalisationTest.template get<1>();
499 testNodeRef = meshTestLocalization->result_analysis().begin()->second.begin()->template get<1>();
501 testEltToPtsQuad[testIdElt].push_back( boost::make_tuple( idq,q,testNodeRef ) );
506 if (notUseOptLocTest) meshTestLocalization->kdtree()->nbNearNeighbor(nbNearNeighborAtStartTest);
512 template <
typename Mesh1Type>
514 localization( mpl::int_<MESH_ELEMENTS> ,
515 boost::shared_ptr<Mesh1Type> meshTest,
516 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > & testEltToPtsQuad )
519 auto begin_elt_it = this->beginElement();
522 pc_ptrtype geopc(
new pc_type( begin_elt_it->gm(), this->im().points() ) );
523 gmc_ptrtype gmc(
new gmc_type( begin_elt_it->gm(),*begin_elt_it, geopc ) );
525 auto meshTestLocalization = meshTest->tool_localization();
526 meshTestLocalization->updateForUse();
527 const auto nbNearNeighborAtStartTest = meshTestLocalization->kdtree()->nPtMaxNearNeighbor();
528 bool notUseOptLocTest = Mesh1Type::nDim!=Mesh1Type::nRealDim;
529 if (notUseOptLocTest) meshTestLocalization->kdtree()->nbNearNeighbor(Mesh1Type::element_type::numPoints);
531 matrix_node_type ptsReal( begin_elt_it->vertices().size1(), 1 );
533 node_type testNodeRef;
535 element_iterator_type elt_it, elt_en;
537 auto itListRange = M_listRange.begin();
538 auto const enListRange = M_listRange.end();
539 for (
size_type ide = 0 ; itListRange!=enListRange ; ++itListRange)
541 boost::tie( boost::tuples::ignore, elt_it, elt_en ) = *itListRange;
542 for ( ; elt_it != elt_en; ++elt_it, ++ide )
544 gmc->update( *elt_it );
546 for (
int q = 0; q < gmc->nPoints(); ++ q )
552 auto testAnalysis = meshTestLocalization->searchElement( gmc->xReal( q ) );
553 testIdElt = testAnalysis.template get<1>();
554 testNodeRef = testAnalysis.template get<2>();
556 ublas::column(ptsReal,0 ) = gmc->xReal( q );
558 auto resLocalisationTest = meshTestLocalization->run_analysis(ptsReal,testIdElt,elt_it->vertices(),mpl::int_<0>());
559 testIdElt = resLocalisationTest.template get<1>();
560 testNodeRef = meshTestLocalization->result_analysis().begin()->second.begin()->template get<1>();
563 testEltToPtsQuad[testIdElt].push_back( boost::make_tuple( idq,q,testNodeRef ) );
568 if (notUseOptLocTest) meshTestLocalization->kdtree()->nbNearNeighbor(nbNearNeighborAtStartTest);
574 template <
typename Mesh1Type,
typename Mesh2Type>
576 localization( mpl::int_<MESH_FACES> ,
577 boost::shared_ptr<Mesh1Type> meshTest,
578 boost::shared_ptr<Mesh2Type> meshTrial,
579 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > & testEltToPtsQuad,
580 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > & trialEltToPtsQuad,
581 std::vector<std::list<size_type> > & EltCoupled )
585 auto begin_elt_it = this->beginElement();
587 std::vector<std::map<permutation_type, pc_ptrtype> > __geopc( this->im().nFaces() );
589 for ( uint16_type __f = 0; __f < this->im().nFaces(); ++__f )
591 for ( permutation_type __p( permutation_type::IDENTITY );
592 __p < permutation_type( permutation_type::N_PERMUTATIONS ); ++__p )
595 __geopc[__f][__p] = pc_ptrtype(
new pc_type( begin_elt_it->element( 0 ).gm(), M_ppts[__f].find( __p )->second ) );
599 uint16_type __face_id_in_elt_0 = begin_elt_it->pos_first();
601 gmc_ptrtype gmc(
new gmc_type( begin_elt_it->element( 0 ).gm(),
602 begin_elt_it->element( 0 ),
604 __face_id_in_elt_0 ) );
607 auto meshTrialLocalization = meshTrial->tool_localization();
608 meshTrialLocalization->updateForUse();
609 const auto nbNearNeighborAtStartTrial = meshTrialLocalization->kdtree()->nPtMaxNearNeighbor();
610 bool notUseOptLocTrial = Mesh2Type::nDim!=Mesh2Type::nRealDim;
611 if (notUseOptLocTrial) meshTrialLocalization->kdtree()->nbNearNeighbor(Mesh2Type::element_type::numPoints);
613 auto meshTestLocalization = meshTest->tool_localization();
614 meshTestLocalization->updateForUse();
615 const auto nbNearNeighborAtStartTest = meshTestLocalization->kdtree()->nPtMaxNearNeighbor();
616 bool notUseOptLocTest = Mesh1Type::nDim!=Mesh1Type::nRealDim;
617 if (notUseOptLocTest) meshTestLocalization->kdtree()->nbNearNeighbor(Mesh1Type::element_type::numPoints);
620 matrix_node_type ptsReal( begin_elt_it->vertices().size1(), 1 );
622 node_type trialNodeRef,testNodeRef;
624 bool quadMeshIsSameThatTrialMesh=
false,quadMeshIsSameThatTestMesh=
false;
625 if ( dynamic_cast<void*>( const_cast<MeshBase*>( begin_elt_it->mesh() ) ) ==
dynamic_cast<void*
>( meshTrial.get() ) )
626 quadMeshIsSameThatTrialMesh=
true;
627 if ( dynamic_cast<void*>( const_cast<MeshBase*>( begin_elt_it->mesh() ) ) ==
dynamic_cast<void*
>( meshTest.get() ) )
628 quadMeshIsSameThatTestMesh=
true;
630 #if FEELPP_EXPORT_QUADLOCALIZATION
631 std::map<size_type,std::list<size_type> > mapBetweenMeshes_test;
632 std::map<size_type,std::list<size_type> > mapBetweenMeshes_trial;
635 element_iterator_type elt_it, elt_en;
637 auto itListRange = M_listRange.begin();
638 auto const enListRange = M_listRange.end();
639 for (
size_type ide = 0 ; itListRange!=enListRange ; ++itListRange)
641 boost::tie( boost::tuples::ignore, elt_it, elt_en ) = *itListRange;
642 for ( ; elt_it != elt_en; ++elt_it, ++ide )
644 __face_id_in_elt_0 = elt_it->pos_first();
647 gmc->update( elt_it->element( 0 ), __face_id_in_elt_0 );
651 for (
int q = 0; q < gmc->nPoints(); ++ q )
657 ublas::column(ptsReal,0 ) = gmc->xReal( q );
658 if (!quadMeshIsSameThatTrialMesh)
661 auto resLocalisationTrial = meshTrialLocalization->run_analysis(ptsReal,trialIdElt,elt_it->vertices(),mpl::int_<0>());
662 trialIdElt = resLocalisationTrial.template get<1>();
663 trialNodeRef = meshTrialLocalization->result_analysis().begin()->second.begin()->template get<1>();
667 trialIdElt = gmc->id();
668 trialNodeRef = gmc->xRef(q);
670 trialEltToPtsQuad[trialIdElt].push_back( boost::make_tuple( idq,q,trialNodeRef ) );
673 if (!quadMeshIsSameThatTestMesh)
676 auto resLocalisationTest = meshTestLocalization->run_analysis(ptsReal,testIdElt,elt_it->vertices(),mpl::int_<0>());
677 testIdElt = resLocalisationTest.template get<1>();
678 testNodeRef = meshTestLocalization->result_analysis().begin()->second.begin()->template get<1>();
682 testIdElt = gmc->id();
683 testNodeRef = gmc->xRef(q);
685 testEltToPtsQuad[testIdElt].push_back( boost::make_tuple( idq,q,testNodeRef ) );
688 if ( std::find( EltCoupled[testIdElt].begin(),EltCoupled[testIdElt].end(),trialIdElt )==EltCoupled[testIdElt].end() )
690 EltCoupled[testIdElt].push_back( trialIdElt );
695 #if FEELPP_EXPORT_QUADLOCALIZATION
696 mapBetweenMeshes_test[testIdElt].push_back( trialIdElt );
697 mapBetweenMeshes_trial[trialIdElt].push_back( testIdElt );
703 #if FEELPP_EXPORT_QUADLOCALIZATION
704 this->exportQuadLocalization(meshTest,meshTrial,mapBetweenMeshes_test,mapBetweenMeshes_trial);
707 if (notUseOptLocTest) meshTestLocalization->kdtree()->nbNearNeighbor(nbNearNeighborAtStartTest);
708 if (notUseOptLocTrial) meshTrialLocalization->kdtree()->nbNearNeighbor(nbNearNeighborAtStartTrial);
716 template <
typename Mesh1Type,
typename Mesh2Type>
718 localization( mpl::int_<MESH_ELEMENTS> ,
719 boost::shared_ptr<Mesh1Type> meshTest,
720 boost::shared_ptr<Mesh2Type> meshTrial,
721 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > & testEltToPtsQuad,
722 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > & trialEltToPtsQuad,
723 std::vector<std::list<size_type> > & EltCoupled )
727 auto begin_elt_it = this->beginElement();
731 pc_ptrtype geopc(
new pc_type( begin_elt_it->gm(), this->im().points() ) );
732 gmc_ptrtype gmc(
new gmc_type( begin_elt_it->gm(),*begin_elt_it, geopc ) );
734 auto meshTrialLocalization = meshTrial->tool_localization();
735 meshTrialLocalization->updateForUse();
736 const auto nbNearNeighborAtStartTrial = meshTrialLocalization->kdtree()->nPtMaxNearNeighbor();
737 bool notUseOptLocTrial = Mesh2Type::nDim!=Mesh2Type::nRealDim;
738 if (notUseOptLocTrial) meshTrialLocalization->kdtree()->nbNearNeighbor(Mesh2Type::element_type::numPoints);
740 auto meshTestLocalization = meshTest->tool_localization();
741 meshTestLocalization->updateForUse();
742 const auto nbNearNeighborAtStartTest = meshTestLocalization->kdtree()->nPtMaxNearNeighbor();
743 bool notUseOptLocTest = Mesh1Type::nDim!=Mesh1Type::nRealDim;
744 if (notUseOptLocTest) meshTestLocalization->kdtree()->nbNearNeighbor(Mesh1Type::element_type::numPoints);
746 #if FEELPP_EXPORT_QUADLOCALIZATION
747 std::map<size_type,std::list<size_type> > mapBetweenMeshes_test;
748 std::map<size_type,std::list<size_type> > mapBetweenMeshes_trial;
758 matrix_node_type ptsReal( begin_elt_it->vertices().size1(), 1 );
760 node_type trialNodeRef,testNodeRef;
762 bool quadMeshIsSameThatTrialMesh=
false,quadMeshIsSameThatTestMesh=
false;
763 if ( dynamic_cast<void*>( const_cast<MeshBase*>( begin_elt_it->mesh() ) ) ==
dynamic_cast<void*
>( meshTrial.get() ) )
764 quadMeshIsSameThatTrialMesh=
true;
765 if ( dynamic_cast<void*>( const_cast<MeshBase*>( begin_elt_it->mesh() ) ) ==
dynamic_cast<void*
>( meshTest.get() ) )
766 quadMeshIsSameThatTestMesh=
true;
768 element_iterator_type elt_it, elt_en;
770 auto itListRange = M_listRange.begin();
771 auto const enListRange = M_listRange.end();
772 for (
size_type ide = 0 ; itListRange!=enListRange ; ++itListRange)
774 boost::tie( boost::tuples::ignore, elt_it, elt_en ) = *itListRange;
775 for ( ; elt_it != elt_en; ++elt_it, ++ide )
777 gmc->update( *elt_it );
779 for (
int q = 0; q < gmc->nPoints(); ++ q )
785 ublas::column(ptsReal,0 ) = gmc->xReal( q );
786 if (!quadMeshIsSameThatTrialMesh)
789 auto resLocalisationTrial = meshTrialLocalization->run_analysis(ptsReal,trialIdElt,elt_it->vertices(),mpl::int_<0>());
790 trialIdElt = resLocalisationTrial.template get<1>();
791 trialNodeRef = meshTrialLocalization->result_analysis().begin()->second.begin()->template get<1>();
795 trialIdElt = gmc->id();
796 trialNodeRef = gmc->xRef(q);
799 trialEltToPtsQuad[trialIdElt].push_back( boost::make_tuple( idq,q,trialNodeRef ) );
802 if (!quadMeshIsSameThatTestMesh)
805 auto resLocalisationTest = meshTestLocalization->run_analysis(ptsReal,testIdElt,elt_it->vertices(),mpl::int_<0>());
806 testIdElt = resLocalisationTest.template get<1>();
807 testNodeRef = meshTestLocalization->result_analysis().begin()->second.begin()->template get<1>();
811 testIdElt = gmc->id();
812 testNodeRef = gmc->xRef(q);
814 testEltToPtsQuad[testIdElt].push_back( boost::make_tuple( idq,q,testNodeRef ) );
819 if ( std::find( EltCoupled[testIdElt].begin(),EltCoupled[testIdElt].end(),trialIdElt )==EltCoupled[testIdElt].end() )
821 EltCoupled[testIdElt].push_back( trialIdElt );
823 #if FEELPP_EXPORT_QUADLOCALIZATION
824 mapBetweenMeshes_test[testIdElt].push_back( trialIdElt );
825 mapBetweenMeshes_trial[trialIdElt].push_back( testIdElt );
830 #if 0 // try with neighbors
831 auto const& geoeltTrial = meshTrial->element( trialIdElt );
832 for ( uint16_type ms=0; ms < geoeltTrial.nNeighbors(); ms++ )
834 const size_type neighborTrial_id = geoeltTrial.neighbor( ms ).first;
837 auto resNeighboor = meshTrialLocalization->isIn( neighborTrial_id,gmc->xReal( q ),elt_it->vertices(),mpl::int_<1>() );
838 if (resNeighboor.get<0>())
840 trialEltToPtsQuad[neighborTrial_id].push_back( boost::make_tuple( idq,q,resNeighboor.get<1>() ) );
842 if ( std::find( EltCoupled[testIdElt].begin(),EltCoupled[testIdElt].end(),neighborTrial_id )==EltCoupled[testIdElt].end() )
844 EltCoupled[testIdElt].push_back( neighborTrial_id );
846 #if FEELPP_EXPORT_QUADLOCALIZATION
847 mapBetweenMeshes_trial[neighborTrial_id].push_back( testIdElt );
852 auto const& geoeltTest = meshTest->element( testIdElt );
853 for ( uint16_type ms=0; ms < geoeltTest.nNeighbors(); ms++ )
855 const size_type neighborTest_id = geoeltTest.neighbor( ms ).first;
858 auto resNeighboor = meshTestLocalization->isIn( neighborTest_id,gmc->xReal( q ),elt_it->vertices(),mpl::int_<1>() );
859 if (resNeighboor.get<0>())
861 testEltToPtsQuad[neighborTest_id].push_back( boost::make_tuple( idq,q,resNeighboor.get<1>() ) );
863 if ( std::find( EltCoupled[neighborTest_id].begin(),EltCoupled[neighborTest_id].end(),trialIdElt )==EltCoupled[neighborTest_id].end() )
865 EltCoupled[neighborTest_id].push_back( trialIdElt );
867 #if FEELPP_EXPORT_QUADLOCALIZATION
868 mapBetweenMeshes_test[neighborTest_id].push_back( trialIdElt );
876 #if FEELPP_EXPORT_QUADLOCALIZATION
877 this->exportQuadLocalization(meshTest,meshTrial,mapBetweenMeshes_test,mapBetweenMeshes_trial);
880 if (notUseOptLocTest) meshTestLocalization->kdtree()->nbNearNeighbor(nbNearNeighborAtStartTest);
881 if (notUseOptLocTrial) meshTrialLocalization->kdtree()->nbNearNeighbor(nbNearNeighborAtStartTrial);
888 template <
typename Mesh1Type>
890 update( boost::shared_ptr<Mesh1Type> meshTest )
892 typedef Mesh1Type mesh_test_type;
898 auto nEltTest= meshTest->numElements();
899 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > testEltToPtsQuad( nEltTest );
900 this->localization( mpl::int_<iDim>(),meshTest,testEltToPtsQuad );
904 auto eltTest_it = testEltToPtsQuad.begin();
905 auto eltTest_en = testEltToPtsQuad.end();
907 for ( ; eltTest_it != eltTest_en ; ++eltTest_it, ++theIdEltTest )
909 if ( eltTest_it->size()>0 )
911 auto nPtsRef = eltTest_it->size();
912 matrix_node_type ptsRefTest( mesh_test_type::nDim,nPtsRef );
913 std::vector<boost::tuple<size_type,size_type> > indexLocalToQuad( nPtsRef );
915 auto idq_it = eltTest_it->begin();
916 auto idq_en = eltTest_it->end();
917 uint16_type cptIdq = 0;
919 for ( ; idq_it != idq_en ; ++idq_it, ++cptIdq )
921 indexLocalToQuad[cptIdq] = boost::make_tuple( idq_it->template get<1>(),idq_it->template get<0>() );
922 ublas::column( ptsRefTest, cptIdq ) = idq_it->template get<2>();
925 M_resLinear.push_back( boost::make_tuple( theIdEltTest, indexLocalToQuad, ptsRefTest ) );
933 template <
typename Mesh1Type,
typename Mesh2Type>
935 update( boost::shared_ptr<Mesh1Type> meshTest, boost::shared_ptr<Mesh2Type> meshTrial )
937 typedef Mesh1Type mesh_test_type;
938 typedef Mesh2Type mesh_trial_type;
941 M_resBilinear.clear();
943 auto nEltTrial= meshTrial->numElements();
944 auto nEltTest= meshTest->numElements();
945 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > trialEltToPtsQuad( nEltTrial ) ;
946 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > > testEltToPtsQuad( nEltTest );
947 std::vector<std::list<size_type> > EltCoupled( nEltTest );
950 this->localization( mpl::int_<iDim>(),
958 auto eltCoupled_it = EltCoupled.begin();
959 auto eltCoupled_en = EltCoupled.end();
962 for ( ; eltCoupled_it != eltCoupled_en ; ++eltCoupled_it, ++theIdEltTest )
964 if (eltCoupled_it->size()==0)
continue;
967 result_temp1_type mapTrial2idq;
970 auto quadPtTest_it = testEltToPtsQuad[theIdEltTest].begin();
971 auto const quadPtTest_en = testEltToPtsQuad[theIdEltTest].end();
972 for ( ; quadPtTest_it != quadPtTest_en ; ++quadPtTest_it )
975 auto theIdq = quadPtTest_it->template get<0>();
976 auto theq = quadPtTest_it->template get<1>();
978 auto theNodeRefTest = quadPtTest_it->template get<2>();
980 auto theRes = this->findQuadPt( *eltCoupled_it,trialEltToPtsQuad,theIdq );
982 auto const theIdEltTrial = theRes.template get<0>();
983 auto const theNodeRefTrial = theRes.template get<1>();
987 mapTrial2idq[theIdEltTrial].push_back( boost::make_tuple( theq,theIdq,theNodeRefTest,theNodeRefTrial ) );
990 result_temp2_type mapiIdTrial2qAndPtRef;
993 auto map_it = mapTrial2idq.begin();
994 auto map_en = mapTrial2idq.end();
995 for ( ; map_it != map_en ; ++map_it )
998 auto nPtsRef = map_it->second.size();
1000 matrix_node_type ptsRefTest( mesh_test_type::nDim,nPtsRef );
1001 matrix_node_type ptsRefTrial( mesh_trial_type::nDim,nPtsRef );
1003 std::vector<boost::tuple<size_type,size_type> > indexLocalToQuad( nPtsRef );
1005 auto idq_it = map_it->second.begin();
1006 auto idq_en = map_it->second.end();
1007 uint16_type cptIdq = 0;
1009 for ( ; idq_it != idq_en ; ++idq_it, ++cptIdq )
1011 indexLocalToQuad[cptIdq] = boost::make_tuple( idq_it->template get<0>(),idq_it->template get<1>() );
1012 ublas::column( ptsRefTest, cptIdq ) = idq_it->template get<2>();
1013 ublas::column( ptsRefTrial, cptIdq ) = idq_it->template get<3>();
1016 mapiIdTrial2qAndPtRef[map_it->first] = boost::make_tuple( indexLocalToQuad,ptsRefTest,ptsRefTrial );
1020 M_resBilinear.push_back( boost::make_tuple( theIdEltTest,mapiIdTrial2qAndPtRef ) );
1027 boost::tuple<size_type,node_type >
1028 findQuadPt( std::list<size_type>
const & eltCoupled_it,
1029 std::vector<std::list<boost::tuple< size_type,size_type,node_type> > >
const & trialPtQuadToElt,
1033 ublas::vector<double> resNodeRef;
1036 auto eltTrial_it = eltCoupled_it.begin();
1037 auto eltTrial_en = eltCoupled_it.end();
1039 while ( !find && eltTrial_it != eltTrial_en )
1041 auto quadPtTrial_it = trialPtQuadToElt[*eltTrial_it].begin();
1042 auto quadPtTrial_en = trialPtQuadToElt[*eltTrial_it].end();
1044 while ( !find && quadPtTrial_it != quadPtTrial_en )
1046 if ( quadPtTrial_it->template get<0>() == idq )
1049 resNodeRef=quadPtTrial_it->template get<2>();
1052 else ++quadPtTrial_it;
1055 if ( !find ) ++eltTrial_it;
1057 else resId=*eltTrial_it;
1060 if ( !find ) std::cout <<
"BUG in findQuadPt " << std::endl;
1062 return boost::make_tuple( resId,resNodeRef );
1066 #if FEELPP_EXPORT_QUADLOCALIZATION
1067 template <
typename Mesh1Type,
typename Mesh2Type>
1069 exportQuadLocalization( boost::shared_ptr<Mesh1Type> meshTest,
1070 boost::shared_ptr<Mesh2Type> meshTrial,
1071 std::map<
size_type,std::list<size_type> >
const& mapBetweenMeshes_test,
1072 std::map<
size_type,std::list<size_type> >
const& mapBetweenMeshes_trial )
const
1074 typedef FunctionSpace<Mesh1Type, bases<Lagrange<0, Scalar,Discontinuous> > > space_test_disc_type;
1075 auto spaceGraphProjTest = space_test_disc_type::New( meshTest );
1076 auto elemTest_itt = meshTest->beginElementWithProcessId( 0 );
1077 auto elemTest_ent = meshTest->endElementWithProcessId( 0 );
1078 auto graphProjTest = spaceGraphProjTest->element();graphProjTest.zero();
1079 for ( ; elemTest_itt != elemTest_ent; ++elemTest_itt )
1081 if ( mapBetweenMeshes_test.find( elemTest_itt->id() ) != mapBetweenMeshes_test.end() )
1083 if ( mapBetweenMeshes_test.find( elemTest_itt->id() )->second.size()>0 )
1085 auto element_dof1 = spaceGraphProjTest->dof()->getIndices( elemTest_itt->id() );
1086 graphProjTest( element_dof1[0] ) = 1;
1092 exporterTest->step( 0 )->setMesh( graphProjTest.mesh() );
1093 exporterTest->step( 0 )->add(
"quadLocalizationProjTest", graphProjTest );
1094 exporterTest->save();
1096 typedef FunctionSpace<Mesh2Type, bases<Lagrange<0, Scalar,Discontinuous> > > space_trial_disc_type;
1097 auto spaceGraphProjTrial = space_trial_disc_type::New( meshTrial );
1098 auto elemTrial_itt = meshTrial->beginElementWithProcessId( 0 );
1099 auto elemTrial_ent = meshTrial->endElementWithProcessId( 0 );
1100 auto graphProjTrial = spaceGraphProjTrial->element();graphProjTrial.zero();
1101 for ( ; elemTrial_itt != elemTrial_ent; ++elemTrial_itt )
1103 if ( mapBetweenMeshes_trial.find( elemTrial_itt->id() ) != mapBetweenMeshes_trial.end() )
1105 if ( mapBetweenMeshes_trial.find( elemTrial_itt->id() )->second.size()>0 )
1107 auto element_dof2 = spaceGraphProjTrial->dof()->getIndices( elemTrial_itt->id() );
1108 graphProjTrial( element_dof2[0] ) = 1;
1114 exporterTrial->step( 0 )->setMesh( graphProjTrial.mesh() );
1115 exporterTrial->step( 0 )->add(
"quadLocalizationProjTrial", graphProjTrial );
1116 exporterTrial->save();
1125 template<
typename FE1,
typename FE2,
typename ElemContType>
1126 struct bilinearformContext
1128 typedef Expr expression_type;
1131 typedef vf::detail::BilinearForm<FE1,FE2,ElemContType> FormType;
1133 typedef typename FormType::gm_1_type gm_formTest_type;
1134 typedef typename FormType::mesh_element_1_type geoelement_formTest_type;
1135 typedef typename gm_formTest_type::template Context<expression_type::context|vm::POINT,geoelement_formTest_type> gmc_formTest_type;
1136 typedef boost::shared_ptr<gmc_formTest_type> gmc_formTest_ptrtype;
1137 typedef typename gm_formTest_type::precompute_type pc_formTest_type;
1138 typedef typename gm_formTest_type::precompute_ptrtype pc_formTest_ptrtype;
1141 typedef typename FormType::gm_2_type gm_formTrial_type;
1142 typedef typename FormType::mesh_element_2_type geoelement_formTrial_type;
1143 typedef typename gm_formTrial_type::template Context<expression_type::context|vm::POINT,geoelement_formTrial_type> gmc_formTrial_type;
1144 typedef boost::shared_ptr<gmc_formTrial_type> gmc_formTrial_ptrtype;
1145 typedef typename gm_formTrial_type::precompute_type pc_formTrial_type;
1146 typedef typename gm_formTrial_type::precompute_ptrtype pc_formTrial_ptrtype;
1149 typedef std::list<boost::tuple< std::vector<boost::tuple<size_type,size_type> >,
1151 gmc_formTest_ptrtype,
1152 gmc_formTrial_ptrtype > > return_loc_type;
1153 typedef std::list<return_loc_type> return_type;
1158 template<
typename FE1,
typename FE2,
typename ElemContType>
1160 precompute(vf::detail::BilinearForm<FE1,FE2,ElemContType>
const& __form)
1162 typedef typename bilinearformContext<FE1,FE2,ElemContType>::gmc_formTest_type gmc_formTest_type;
1163 typedef typename bilinearformContext<FE1,FE2,ElemContType>::gmc_formTest_ptrtype gmc_formTest_ptrtype;
1164 typedef typename bilinearformContext<FE1,FE2,ElemContType>::pc_formTest_type pc_formTest_type;
1165 typedef typename bilinearformContext<FE1,FE2,ElemContType>::pc_formTest_ptrtype pc_formTest_ptrtype;
1167 typedef typename bilinearformContext<FE1,FE2,ElemContType>::gmc_formTrial_type gmc_formTrial_type;
1168 typedef typename bilinearformContext<FE1,FE2,ElemContType>::gmc_formTrial_ptrtype gmc_formTrial_ptrtype;
1169 typedef typename bilinearformContext<FE1,FE2,ElemContType>::pc_formTrial_type pc_formTrial_type;
1170 typedef typename bilinearformContext<FE1,FE2,ElemContType>::pc_formTrial_ptrtype pc_formTrial_ptrtype;
1172 auto meshTrial = __form.trialSpace()->mesh();
1173 auto meshTest = __form.testSpace()->mesh();
1175 typename bilinearformContext<FE1,FE2,ElemContType>::return_type theres;
1177 auto res_it = this->result().begin();
1178 auto const res_en = this->result().end();
1179 for ( ; res_it != res_en ; ++res_it )
1181 auto const idEltTest = res_it->template get<0>();
1182 auto const& map = res_it->template get<1>();
1183 auto map_it = map.begin();
1184 auto const map_en = map.end();
1185 for ( ; map_it != map_en ; ++map_it )
1187 auto const idEltTrial = map_it->first;
1188 auto const& eltTrial = meshTrial->element( idEltTrial );
1189 auto const& eltTest = meshTest->element( idEltTest );
1191 auto const& ptRefTest = map_it->second.template get<1>();
1192 auto const& ptRefTrial = map_it->second.template get<2>();
1193 auto const& themapQuad = map_it->second.template get<0>();
1194 auto vec_gmcExpr = this->getUsableDataInFormContext( themapQuad,ptRefTest,ptRefTrial );
1195 auto gmcExpr_it = vec_gmcExpr.begin();
1196 auto const gmcExpr_en = vec_gmcExpr.end();
1197 typename bilinearformContext<FE1,FE2,ElemContType>::return_loc_type theresloc;
1198 for ( ; gmcExpr_it != gmcExpr_en ; ++gmcExpr_it )
1200 pc_formTest_ptrtype geopcFormTest(
new pc_formTest_type( __form.gm(), gmcExpr_it->template get<2>() ) );
1201 gmc_formTest_ptrtype gmcFormTest(
new gmc_formTest_type( __form.gm(), eltTest , geopcFormTest ) );
1203 pc_formTrial_ptrtype geopcFormTrial(
new pc_formTrial_type( __form.gmTrial(), gmcExpr_it->template get<3>() ) );
1204 gmc_formTrial_ptrtype gmcFormTrial(
new gmc_formTrial_type( __form.gmTrial(), eltTrial, geopcFormTrial ) );
1206 theresloc.push_back(boost::make_tuple(gmcExpr_it->template get<0>(),gmcExpr_it->template get<1>(),gmcFormTest,gmcFormTrial));
1208 theres.push_back(theresloc);
1212 M_precompute = theres;
1213 M_hasPrecompute=
true;
1218 template<
typename FE,
typename VectorType,
typename ElemContType>
1219 struct linearformContext
1221 typedef Expr expression_type;
1224 typedef vf::detail::LinearForm<FE,VectorType,ElemContType> FormType;
1225 typedef typename FormType::gm_type gm_form_type;
1226 typedef typename FormType::mesh_test_element_type geoelement_form_type;
1227 typedef typename gm_form_type::template Context<expression_type::context|vm::POINT,geoelement_form_type> gmc_form_type;
1228 typedef boost::shared_ptr<gmc_form_type> gmc_form_ptrtype;
1229 typedef typename gm_form_type::precompute_type pc_form_type;
1230 typedef typename gm_form_type::precompute_ptrtype pc_form_ptrtype;
1232 typedef std::list<boost::tuple< std::vector<boost::tuple<size_type,size_type> >,
1234 gmc_form_ptrtype > > return_loc_type;
1235 typedef std::list<return_loc_type> return_type;
1239 template<
typename FE,
typename VectorType,
typename ElemContType>
1241 precompute(vf::detail::LinearForm<FE,VectorType,ElemContType>
const& __form)
1244 typedef typename linearformContext<FE,VectorType,ElemContType>::gmc_form_type gmc_form_type;
1245 typedef typename linearformContext<FE,VectorType,ElemContType>::gmc_form_ptrtype gmc_form_ptrtype;
1246 typedef typename linearformContext<FE,VectorType,ElemContType>::pc_form_type pc_form_type;
1247 typedef typename linearformContext<FE,VectorType,ElemContType>::pc_form_ptrtype pc_form_ptrtype;
1249 auto meshTest = __form.testSpace()->mesh();
1251 typename linearformContext<FE,VectorType,ElemContType>::return_type theres;
1253 auto res_it = this->resultLinear().begin();
1254 auto const res_en = this->resultLinear().end();
1255 for ( ; res_it != res_en ; ++res_it )
1257 auto const idEltTest = res_it->template get<0>();
1258 auto const& eltTest = meshTest->element( idEltTest );
1260 auto ptRefTest = res_it->template get<2>();
1261 auto themapQuad = res_it->template get<1>();
1263 auto vec_gmcExpr = this->getUsableDataInFormContext( themapQuad,ptRefTest );
1264 auto gmcExpr_it = vec_gmcExpr.begin();
1265 auto const gmcExpr_en = vec_gmcExpr.end();
1266 typename linearformContext<FE,VectorType,ElemContType>::return_loc_type theresloc;
1267 for ( ; gmcExpr_it != gmcExpr_en ; ++gmcExpr_it )
1270 pc_form_ptrtype geopcForm(
new pc_form_type( __form.gm(), gmcExpr_it->template get<2>() ) );
1271 gmc_form_ptrtype gmcForm(
new gmc_form_type( __form.gm(), eltTest, geopcForm ) );
1272 theresloc.push_back(boost::make_tuple(gmcExpr_it->template get<0>(),gmcExpr_it->template get<1>(),gmcForm));
1274 theres.push_back(theresloc);
1278 M_precompute = theres;
1279 M_hasPrecompute=
true;
1284 template<
typename FE1,
typename FE2,
typename ElemContType>
1285 typename bilinearformContext<FE1,FE2,ElemContType>::return_type
const&
1286 getPrecompute(vf::detail::BilinearForm<FE1,FE2,ElemContType>
const& __form)
const
1288 return boost::any_cast<
typename bilinearformContext<FE1,FE2,ElemContType>::return_type
const&>( M_precompute);
1291 template<
typename FE,
typename VectorType,
typename ElemContType>
1292 typename linearformContext<FE,VectorType,ElemContType>::return_type
const&
1293 getPrecompute(vf::detail::LinearForm<FE,VectorType,ElemContType>
const& __form)
const
1295 return boost::any_cast<
typename linearformContext<FE,VectorType,ElemContType>::return_type
const&>( M_precompute);
1300 std::list<range_iterator> M_listRange;
1302 element_iterator_type M_eltbegin;
1303 element_iterator_type M_eltend;
1304 mutable im_type M_im;
1305 QuadMapped<im_type> M_qm;
1306 permutation_points_type M_ppts;
1308 result_container_type M_resBilinear;
1310 result_container_linear_type M_resLinear;
1312 bool M_hasPrecompute;
1313 boost::any M_precompute;
1323 template <
typename RangeType>
1324 struct quadptlocrangetype
1326 typedef typename mpl::if_< boost::is_std_list<RangeType>,
1327 mpl::identity<RangeType>,
1328 mpl::identity<std::list<RangeType> > >::type::type::value_type type;
1333 template<
typename MeshTestType,
typename MeshTrialType,
typename IteratorRange,
typename Im,
typename Expr>
1334 boost::shared_ptr<QuadPtLocalization<typename Feel::detail::quadptlocrangetype<IteratorRange>::type,Im,Expr> >
1335 quadPtLocPtr( boost::shared_ptr<MeshTestType> meshTest, boost::shared_ptr<MeshTrialType> meshTrial, IteratorRange
const& elts,Im
const& im,Expr
const& expr )
1337 typedef QuadPtLocalization<typename Feel::detail::quadptlocrangetype<IteratorRange>::type,Im,Expr> quadptloc_type;
1338 boost::shared_ptr<quadptloc_type> res(
new quadptloc_type(elts) );
1339 res->update( meshTest,meshTrial );
1343 template<
typename MeshTestType,
typename IteratorRange,
typename Im,
typename Expr>
1344 boost::shared_ptr<QuadPtLocalization<typename Feel::detail::quadptlocrangetype<IteratorRange>::type,Im,Expr> >
1345 quadPtLocPtr( boost::shared_ptr<MeshTestType> meshTest, IteratorRange
const& elts,Im
const& im,Expr
const& expr )
1347 typedef QuadPtLocalization<typename Feel::detail::quadptlocrangetype<IteratorRange>::type,Im,Expr> quadptloc_type;
1348 boost::shared_ptr<quadptloc_type> res(
new quadptloc_type(elts) );
1349 res->update( meshTest );
const size_type invalid_size_type_value
Definition: feelcore/feel.hpp:360
static boost::shared_ptr< Exporter< MeshType, N > > New(std::string const &exportername, std::string prefix=Environment::about().appName(), WorldComm const &worldComm=Environment::worldComm())
Definition: exporterimpl.hpp:127
size_t size_type
Indices (starting from 0)
Definition: feelcore/feel.hpp:319