ViennaGrid - The Vienna Grid Library: viennagrid/mesh/element

Go to the documentation of this file.
00001 #ifndef VIENNAGRID_ELEMENT_CREATION_HPP
00002 #define VIENNAGRID_ELEMENT_CREATION_HPP
00003 
00004 /* =======================================================================
00005    Copyright (c) 2011-2014, Institute for Microelectronics,
00006                             Institute for Analysis and Scientific Computing,
00007                             TU Wien.
00008 
00009                             -----------------
00010                      ViennaGrid - The Vienna Grid Library
00011                             -----------------
00012 
00013    License:      MIT (X11), see file LICENSE in the base directory
00014 ======================================================================= */
00015 
00016 #include "viennagrid/mesh/mesh.hpp"
00017 #include "viennagrid/mesh/segmentation.hpp"
00018 #include "viennagrid/topology/plc.hpp"
00019 #include "viennagrid/algorithm/norm.hpp"
00020 
00025 namespace viennagrid
00026 {
00027   // implementation for make_element and make_element_with_id
00028   template<typename ElementTagT>
00029   struct make_element_impl
00030   {
00032     template<typename MeshOrSegmentHandleTypeT, typename VertexHandleIteratorT>
00033     static typename result_of::handle<MeshOrSegmentHandleTypeT, ElementTagT>::type
00034     make(MeshOrSegmentHandleTypeT & mesh_obj,
00035          VertexHandleIteratorT vertices_begin,
00036          VertexHandleIteratorT const & vertices_end)
00037     {
00038       typedef typename viennagrid::result_of::element<MeshOrSegmentHandleTypeT, ElementTagT>::type ElementType;
00039       ElementType element = ElementType( detail::inserter(mesh_obj).get_physical_container_collection() );
00040 
00041       unsigned int element_index = 0;
00042       for ( ; vertices_begin != vertices_end; ++vertices_begin, ++element_index )
00043           viennagrid::set_vertex( element, *vertices_begin, element_index );
00044 
00045       return detail::push_element<true, true>(mesh_obj, element).first;
00046     }
00047 
00049     template<typename MeshOrSegmentHandleTypeT, typename VertexHandleIteratorT, typename IDT>
00050     static typename result_of::handle<MeshOrSegmentHandleTypeT, ElementTagT>::type
00051     make(MeshOrSegmentHandleTypeT & mesh_obj,
00052          VertexHandleIteratorT vertices_begin,
00053          VertexHandleIteratorT const & vertices_end,
00054          IDT id)
00055     {
00056       typedef typename viennagrid::result_of::element<MeshOrSegmentHandleTypeT, ElementTagT>::type ElementType;
00057       ElementType element = ElementType( detail::inserter(mesh_obj).get_physical_container_collection() );
00058 
00059       element.id( id );
00060 
00061       unsigned int element_index = 0;
00062       for ( ; vertices_begin != vertices_end; ++vertices_begin, ++element_index )
00063           viennagrid::set_vertex( element, *vertices_begin, element_index );
00064 
00065       return detail::push_element<false, true>(mesh_obj, element ).first;
00066     }
00067   };
00068 
00069   // specialization for PLC: no implementation! Use make_plc instead
00070   template<>
00071   struct make_element_impl<plc_tag>
00072   {
00074     template<typename MeshOrSegmentHandleTypeT, typename VertexHandleIteratorT>
00075     static typename result_of::handle<MeshOrSegmentHandleTypeT, plc_tag>::type
00076     make(MeshOrSegmentHandleTypeT & mesh_obj,
00077          VertexHandleIteratorT vertices_begin,
00078          VertexHandleIteratorT const & vertices_end);
00079 
00081     template<typename MeshOrSegmentHandleTypeT, typename VertexHandleIteratorT, typename IDT>
00082     static typename result_of::handle<MeshOrSegmentHandleTypeT, plc_tag>::type
00083     make(MeshOrSegmentHandleTypeT & mesh_obj,
00084          VertexHandleIteratorT vertices_begin,
00085          VertexHandleIteratorT const & vertices_end,
00086          IDT id);
00087   };
00088 
00089 
00090 
00091   // doxygen doku in forwards.hpp
00092   template<typename ElementTypeOrTagT, typename MeshOrSegmentHandleTypeT, typename VertexHandleIteratorT>
00093   typename result_of::handle<MeshOrSegmentHandleTypeT, ElementTypeOrTagT>::type
00094   make_element(
00095         MeshOrSegmentHandleTypeT & mesh_obj,
00096         VertexHandleIteratorT vertices_begin,
00097         VertexHandleIteratorT const & vertices_end)
00098   {
00099     typedef typename viennagrid::result_of::element_tag<ElementTypeOrTagT>::type ElementTagT;
00100     return make_element_impl<ElementTagT>::make( mesh_obj, vertices_begin, vertices_end );
00101   }
00102 
00103   // doxygen doku in forwards.hpp
00104   template<typename ElementTypeOrTagT, typename MeshOrSegmentHandleTypeT, typename VertexHandleIteratorT, typename IDT>
00105   typename result_of::handle<MeshOrSegmentHandleTypeT, ElementTypeOrTagT>::type
00106   make_element_with_id(
00107         MeshOrSegmentHandleTypeT & mesh_obj,
00108         VertexHandleIteratorT vertices_begin,
00109         VertexHandleIteratorT const & vertices_end,
00110         IDT id)
00111   {
00112     typedef typename viennagrid::result_of::element_tag<ElementTypeOrTagT>::type ElementTagT;
00113     return make_element_impl<ElementTagT>::make( mesh_obj, vertices_begin, vertices_end, id );
00114   }
00115 
00116 
00117 
00125   template<typename MeshOrSegmentHandleTypeT, typename VertexHandleIteratorT>
00126   typename result_of::cell_handle<MeshOrSegmentHandleTypeT>::type
00127   make_cell(
00128         MeshOrSegmentHandleTypeT & mesh_segment,
00129         VertexHandleIteratorT vertices_begin,
00130         VertexHandleIteratorT const & vertices_end)
00131   {
00132     typedef typename viennagrid::result_of::cell_tag<MeshOrSegmentHandleTypeT>::type CellTagType;
00133     return make_element<CellTagType>( mesh_segment, vertices_begin, vertices_end );
00134   }
00135 
00136 
00145   template<typename MeshOrSegmentHandleTypeT, typename VertexHandleIteratorT>
00146   typename result_of::cell_handle<MeshOrSegmentHandleTypeT>::type
00147   make_cell_with_id(
00148         MeshOrSegmentHandleTypeT & mesh_segment,
00149         VertexHandleIteratorT vertices_begin,
00150         VertexHandleIteratorT const & vertices_end,
00151         typename viennagrid::result_of::cell<MeshOrSegmentHandleTypeT>::type::id_type id)
00152   {
00153     typedef typename viennagrid::result_of::cell_tag<MeshOrSegmentHandleTypeT>::type CellTagType;
00154     return make_element_with_id<CellTagType>( mesh_segment, vertices_begin, vertices_end, id );
00155   }
00156 
00157 
00158 
00159 
00160 
00161 
00162   // doxygen doku in forwards.hpp
00163   template<typename MeshOrSegmentHandleTypeT>
00164   typename result_of::vertex_handle<MeshOrSegmentHandleTypeT>::type make_vertex(MeshOrSegmentHandleTypeT & mesh_segment)
00165   {
00166     typedef typename result_of::element<MeshOrSegmentHandleTypeT, vertex_tag>::type element_type;
00167     return detail::push_element<true, true>(mesh_segment, element_type() ).first;
00168   }
00169 
00170   // doxygen doku in forwards.hpp
00171   template<typename MeshOrSegmentHandleTypeT>
00172   typename result_of::vertex_handle<MeshOrSegmentHandleTypeT>::type make_vertex(
00173         MeshOrSegmentHandleTypeT & mesh_obj,
00174         typename result_of::point<MeshOrSegmentHandleTypeT>::type const & point)
00175   {
00176     typename viennagrid::result_of::vertex_handle<MeshOrSegmentHandleTypeT>::type vtx_handle = make_vertex(mesh_obj);
00177     viennagrid::point(mesh_obj, vtx_handle) = point;
00178     return vtx_handle;
00179   }
00180 
00181   // doxygen doku in forwards.hpp
00182   template<typename MeshOrSegmentHandleTypeT>
00183   typename result_of::vertex_handle<MeshOrSegmentHandleTypeT>::type make_vertex_with_id(
00184         MeshOrSegmentHandleTypeT & mesh_obj,
00185         typename viennagrid::result_of::element<MeshOrSegmentHandleTypeT, vertex_tag>::type::id_type id,
00186         typename result_of::point<MeshOrSegmentHandleTypeT>::type const & point)
00187   {
00188     typedef typename result_of::vertex<MeshOrSegmentHandleTypeT>::type VertexType;
00189     VertexType element;
00190     element.id( id );
00191 
00192     typename result_of::vertex_handle<MeshOrSegmentHandleTypeT>::type ret = detail::push_element<false, true>(mesh_obj, element ).first;
00193     viennagrid::point(mesh_obj, ret) = point;
00194 
00195     return ret;
00196   }
00197 
00198   // doxygen doku in forwards.hpp
00199   template<typename MeshOrSegmentHandleTypeT>
00200   typename result_of::vertex_handle<MeshOrSegmentHandleTypeT>::type make_unique_vertex(
00201         MeshOrSegmentHandleTypeT & mesh_obj,
00202         typename result_of::point<MeshOrSegmentHandleTypeT>::type const & point,
00203         typename result_of::coord<MeshOrSegmentHandleTypeT>::type tolerance)
00204   {
00205     typedef typename result_of::element_range<MeshOrSegmentHandleTypeT, vertex_tag>::type vertex_range_type;
00206     typedef typename result_of::iterator<vertex_range_type>::type vertex_range_iterator;
00207 
00208     if (tolerance > 0)
00209     {
00210       vertex_range_type vertices(mesh_obj);
00211       for (vertex_range_iterator hit = vertices.begin(); hit != vertices.end(); ++hit)
00212       {
00213           if (viennagrid::norm_2( point - viennagrid::point(mesh_obj, *hit) ) < tolerance )
00214               return hit.handle();
00215       }
00216     }
00217 
00218     return make_vertex(mesh_obj, point);
00219   }
00220 
00221   // doxygen doku in forwards.hpp
00222   template<typename MeshOrSegmentHandleTypeT>
00223   typename result_of::vertex_handle<MeshOrSegmentHandleTypeT>::type make_unique_vertex(
00224         MeshOrSegmentHandleTypeT & mesh_obj,
00225         typename result_of::point<MeshOrSegmentHandleTypeT>::type const & p)
00226   {
00227     return make_unique_vertex( mesh_obj, p, viennagrid::norm_2(p) / 1e6 );
00228   }
00229 
00230 
00231   // doxygen doku in forwards.hpp
00232   template<typename MeshOrSegmentHandleTypeT, typename VertexHandleT>
00233   typename result_of::line_handle<MeshOrSegmentHandleTypeT>::type make_line(
00234         MeshOrSegmentHandleTypeT & mesh_obj,
00235         VertexHandleT v0, VertexHandleT v1)
00236   {
00237     viennagrid::static_array<VertexHandleT, 2> handles;
00238     handles[0] = v0;
00239     handles[1] = v1;
00240 
00241     return make_element<viennagrid::line_tag>( mesh_obj, handles.begin(), handles.end() );
00242   }
00243 
00244   // doxygen doku in forwards.hpp
00245   template<typename MeshOrSegmentHandleTypeT, typename VertexHandleT>
00246   typename result_of::edge_handle<MeshOrSegmentHandleTypeT>::type make_edge(
00247         MeshOrSegmentHandleTypeT & mesh_obj,
00248         VertexHandleT v0, VertexHandleT v1)
00249   {
00250     viennagrid::static_array<VertexHandleT, 2> handles;
00251     handles[0] = v0;
00252     handles[1] = v1;
00253 
00254     return make_element<viennagrid::edge_tag>( mesh_obj, handles.begin(), handles.end() );
00255   }
00256 
00257   // doxygen doku in forwards.hpp
00258   template<typename MeshOrSegmentHandleTypeT, typename VertexHandleT>
00259   typename result_of::triangle_handle<MeshOrSegmentHandleTypeT>::type make_triangle(
00260         MeshOrSegmentHandleTypeT & mesh_obj,
00261         VertexHandleT v0, VertexHandleT v1, VertexHandleT v2)
00262   {
00263     viennagrid::static_array<VertexHandleT, 3> handles;
00264     handles[0] = v0;
00265     handles[1] = v1;
00266     handles[2] = v2;
00267 
00268     return make_element<viennagrid::triangle_tag>( mesh_obj, handles.begin(), handles.end() );
00269   }
00270 
00271   // doxygen doku in forwards.hpp
00272   template<typename MeshOrSegmentHandleTypeT, typename VertexHandleT>
00273   typename result_of::quadrilateral_handle<MeshOrSegmentHandleTypeT>::type make_quadrilateral(
00274         MeshOrSegmentHandleTypeT & mesh_obj,
00275         VertexHandleT v0, VertexHandleT v1, VertexHandleT v2, VertexHandleT v3)
00276   {
00277     viennagrid::static_array<VertexHandleT, 4> handles;
00278     handles[0] = v0;
00279     handles[1] = v1;
00280     handles[2] = v2;
00281     handles[3] = v3;
00282 
00283     return make_element<viennagrid::quadrilateral_tag>( mesh_obj, handles.begin(), handles.end() );
00284   }
00285 
00286 
00287   // doxygen doku in forwards.hpp
00288   template<typename MeshOrSegmentHandleTypeT, typename LineHandleIteratorT, typename VertexHandleIteratorT, typename PointIteratorT>
00289   typename result_of::plc_handle<MeshOrSegmentHandleTypeT>::type make_plc(
00290         MeshOrSegmentHandleTypeT & mesh_obj,
00291         LineHandleIteratorT    lines_begin,           LineHandleIteratorT     lines_end,
00292         VertexHandleIteratorT  loose_vertices_begin,  VertexHandleIteratorT   loose_vertices_end,
00293         PointIteratorT         hole_points_begin,     PointIteratorT          hole_points_end)
00294   {
00295     typedef typename viennagrid::result_of::element<MeshOrSegmentHandleTypeT, plc_tag>::type PLCType;;
00296     typedef typename result_of::handle<MeshOrSegmentHandleTypeT, plc_tag>::type PLCHandleType;
00297     PLCType plc( detail::inserter(mesh_obj).get_physical_container_collection() );
00298 
00299     for ( ; lines_begin != lines_end; ++lines_begin)
00300       plc.container( viennagrid::line_tag() ).insert_unique_handle( *lines_begin );
00301 
00302     for ( ; loose_vertices_begin != loose_vertices_end; ++loose_vertices_begin)
00303       plc.container( viennagrid::vertex_tag() ).insert_unique_handle( *loose_vertices_begin );
00304 
00305     PLCHandleType handle = viennagrid::detail::push_element<true, true>(mesh_obj, plc ).first;
00306 
00307     PLCType & inserted_plc = viennagrid::dereference_handle(mesh_obj, handle);
00308 
00309     std::copy(hole_points_begin, hole_points_end, std::back_inserter( inserted_plc.appendix() ) );
00310     return handle;
00311   }
00312 
00313   // doxygen doku in forwards.hpp
00314   template<typename MeshOrSegmentHandleTypeT, typename LineHandleIteratorT, typename VertexHandleIteratorT>
00315   typename result_of::plc_handle<MeshOrSegmentHandleTypeT>::type make_plc(
00316         MeshOrSegmentHandleTypeT & mesh_obj,
00317         LineHandleIteratorT    lines_begin,           LineHandleIteratorT     lines_end,
00318         VertexHandleIteratorT  loose_vertices_begin,  VertexHandleIteratorT   loose_vertices_end)
00319   {
00320     typedef typename viennagrid::result_of::point<MeshOrSegmentHandleTypeT>::type PointType;
00321     PointType tmp;
00322     return make_plc(mesh_obj, lines_begin, lines_end, loose_vertices_begin, loose_vertices_end, &tmp, &tmp);
00323   }
00324 
00325   template<typename MeshOrSegmentHandleTypeT, typename LineHandleIteratorT>
00326   typename result_of::plc_handle<MeshOrSegmentHandleTypeT>::type make_plc(
00327         MeshOrSegmentHandleTypeT & mesh_obj,
00328         LineHandleIteratorT    lines_begin,           LineHandleIteratorT     lines_end)
00329   {
00330     typedef typename viennagrid::result_of::vertex_handle<MeshOrSegmentHandleTypeT>::type VertexHandleType;
00331     VertexHandleType tmp;
00332     return make_plc(mesh_obj, lines_begin, lines_end, &tmp, &tmp);
00333   }
00334 
00335 
00336 
00337 
00338 
00339   // doxygen doku in forwards.hpp
00340   template<typename MeshOrSegmentHandleTypeT, typename VertexHandleT>
00341   typename result_of::tetrahedron_handle<MeshOrSegmentHandleTypeT>::type make_tetrahedron(
00342         MeshOrSegmentHandleTypeT & mesh_obj,
00343         VertexHandleT v0, VertexHandleT v1, VertexHandleT v2, VertexHandleT v3)
00344   {
00345     viennagrid::static_array<VertexHandleT, 4> handles;
00346     handles[0] = v0;
00347     handles[1] = v1;
00348     handles[2] = v2;
00349     handles[3] = v3;
00350 
00351     return make_element<viennagrid::tetrahedron_tag>( mesh_obj, handles.begin(), handles.end() );
00352   }
00353 
00354 
00355   // doxygen doku in forwards.hpp
00356   template<typename MeshOrSegmentHandleTypeT, typename VertexHandleT>
00357   typename result_of::hexahedron_handle<MeshOrSegmentHandleTypeT>::type make_hexahedron(
00358         MeshOrSegmentHandleTypeT & mesh_obj,
00359         VertexHandleT v0, VertexHandleT v1, VertexHandleT v2, VertexHandleT v3,
00360         VertexHandleT v4, VertexHandleT v5, VertexHandleT v6, VertexHandleT v7)
00361   {
00362     viennagrid::static_array<VertexHandleT, 8> handles;
00363     handles[0] = v0;
00364     handles[1] = v1;
00365     handles[2] = v2;
00366     handles[3] = v3;
00367     handles[4] = v4;
00368     handles[5] = v5;
00369     handles[6] = v6;
00370     handles[7] = v7;
00371 
00372     return make_element<viennagrid::hexahedron_tag>( mesh_obj, handles.begin(), handles.end() );
00373   }
00374 
00375 
00376   // doxygen doku in forwards.hpp
00377   template<typename ElementT, typename MeshOrSegmentHandleT>
00378   typename viennagrid::result_of::handle<
00379       MeshOrSegmentHandleT,
00380       typename viennagrid::result_of::element_tag<ElementT>::type
00381     >::type copy_element( ElementT const & element, MeshOrSegmentHandleT & mesh_segment,
00382                           typename viennagrid::result_of::coord<MeshOrSegmentHandleT>::type tolerance )
00383   {
00384     typedef typename viennagrid::result_of::element_tag<ElementT>::type             ElementTag;
00385     typedef typename viennagrid::result_of::vertex_handle<MeshOrSegmentHandleT>::type   VertexHandleType;
00386     std::vector<VertexHandleType> vertex_handles;
00387 
00388     typedef typename viennagrid::result_of::const_vertex_range<ElementT>::type      VertexRangeType;
00389     typedef typename viennagrid::result_of::iterator<VertexRangeType>::type         VertexRangeIterator;
00390 
00391     VertexRangeType vertices(element);
00392     for (VertexRangeIterator vit = vertices.begin(); vit != vertices.end(); ++vit)
00393       vertex_handles.push_back( make_unique_vertex(mesh_segment, viennagrid::point(*vit)), tolerance );
00394 
00395     return make_element<ElementTag>( mesh_segment, vertex_handles.begin(), vertex_handles.end() );
00396   }
00397 
00398   // doxygen doku in forwards.hpp
00399   template<typename ElementT, typename MeshOrSegmentHandleT>
00400   typename viennagrid::result_of::handle<
00401       MeshOrSegmentHandleT,
00402       typename viennagrid::result_of::element_tag<ElementT>::type
00403     >::type copy_element( ElementT const & element, MeshOrSegmentHandleT & mesh_segment )
00404   {
00405     copy_element(element, mesh_segment, -1.0);
00406   }
00407 
00408   namespace detail
00409   {
00410     template<typename ElementTagT>
00411     struct copy_elements_impl
00412     {
00414       template<typename ElementIteratorT, typename OutputMeshOrSegmentHandleT>
00415       static void copy_elements(ElementIteratorT const & begin, ElementIteratorT const & end,
00416                       OutputMeshOrSegmentHandleT & output_mesh,
00417                       typename viennagrid::result_of::coord<OutputMeshOrSegmentHandleT>::type tolerance )
00418       {
00419         typedef typename std::iterator_traits<ElementIteratorT>::value_type ElementType;
00420         typedef typename viennagrid::result_of::element_tag<ElementType>::type ElementTagType;
00421 
00422         typedef typename viennagrid::result_of::vertex_id<ElementType>::type VertexIDType;
00423         //typedef typename viennagrid::result_of::vertex_handle<ElementType>::type VertexHandleType;
00424         typedef typename viennagrid::result_of::vertex_handle<OutputMeshOrSegmentHandleT>::type OutputVertexHandleType;
00425 
00426         typedef typename viennagrid::result_of::const_vertex_range<ElementType>::type ConstVertexOnElementRangeType;
00427         typedef typename viennagrid::result_of::iterator<ConstVertexOnElementRangeType>::type ConstVertexOnElementIteratorType;
00428 
00429         std::map<VertexIDType, OutputVertexHandleType> vertex_map;
00430         for (ElementIteratorT eit = begin; eit != end; ++eit)
00431         {
00432           ConstVertexOnElementRangeType vertices_on_element( *eit );
00433           std::vector<OutputVertexHandleType> vtx_handles( vertices_on_element.size() );
00434 
00435           unsigned int index = 0;
00436           for (ConstVertexOnElementIteratorType vit = vertices_on_element.begin(); vit != vertices_on_element.end(); ++vit, ++index)
00437           {
00438             typename std::map<VertexIDType, OutputVertexHandleType>::iterator it = vertex_map.find( vit->id() );
00439             if (it == vertex_map.end())
00440             {
00441               vtx_handles[index] = viennagrid::make_unique_vertex( output_mesh, viennagrid::point(*vit), tolerance );
00442               vertex_map[vit->id()] = vtx_handles[index];
00443             }
00444             else
00445               vtx_handles[index] = it->second;
00446           }
00447 
00448           viennagrid::make_element<ElementTagType>( output_mesh, vtx_handles.begin(), vtx_handles.end() );
00449         }
00450       }
00451 
00453       template<typename InputMeshOrSegmentHandleT, typename ElementHandleIteratorT, typename OutputMeshOrSegmentHandleT>
00454       static void copy_elements_by_handle(InputMeshOrSegmentHandleT const & input_mesh,
00455                                 ElementHandleIteratorT const & begin, ElementHandleIteratorT const & end,
00456                                 OutputMeshOrSegmentHandleT & output_mesh,
00457                                 typename viennagrid::result_of::coord<OutputMeshOrSegmentHandleT>::type tolerance )
00458       {
00459         typedef typename std::iterator_traits<ElementHandleIteratorT>::value_type ElementHandleType;
00460         typedef typename viennagrid::detail::result_of::value_type<ElementHandleType>::type ElementType;
00461         typedef typename viennagrid::result_of::element_tag<ElementType>::type ElementTagType;
00462 
00463         typedef typename viennagrid::result_of::vertex_id<ElementType>::type VertexIDType;
00464         typedef typename viennagrid::result_of::vertex_handle<ElementType>::type VertexHandleType;
00465 
00466         typedef typename viennagrid::result_of::const_vertex_range<ElementType>::type ConstVertexOnElementRangeType;
00467         typedef typename viennagrid::result_of::iterator<ConstVertexOnElementRangeType>::type ConstVertexOnElementIteratorType;
00468 
00469         std::map<VertexIDType, VertexHandleType> vertex_map;
00470         for (ElementHandleIteratorT eit = begin; eit != end; ++eit)
00471         {
00472           ConstVertexOnElementRangeType vertices_on_element( viennagrid::dereference_handle(input_mesh, *eit) );
00473           std::vector<VertexHandleType> vtx_handles( vertices_on_element.size() );
00474 
00475           unsigned int index = 0;
00476           for (ConstVertexOnElementIteratorType vit = vertices_on_element.begin(); vit != vertices_on_element.end(); ++vit, ++index)
00477           {
00478             typename std::map<VertexIDType, VertexHandleType>::iterator it = vertex_map.find( vit->id() );
00479             if (it == vertex_map.end())
00480             {
00481               vtx_handles[index] = viennagrid::make_unique_vertex( output_mesh, viennagrid::point(*vit), tolerance );
00482               vertex_map[vit->id()] = vtx_handles[index];
00483             }
00484             else
00485               vtx_handles[index] = it->second;
00486           }
00487 
00488           viennagrid::make_element<ElementTagType>( output_mesh, vtx_handles.begin(), vtx_handles.end() );
00489         }
00490       }
00491     };
00492 
00493 
00494     template<>
00495     struct copy_elements_impl<plc_tag>
00496     {
00498       template<typename ElementIteratorT, typename OutputMeshOrSegmentHandleT>
00499       static void copy_elements(ElementIteratorT const & begin, ElementIteratorT const & end,
00500                       OutputMeshOrSegmentHandleT & output_mesh,
00501                       typename viennagrid::result_of::coord<OutputMeshOrSegmentHandleT>::type tolerance )
00502       {
00503         typedef typename std::iterator_traits<ElementIteratorT>::value_type ElementType;
00504 
00505         typedef typename viennagrid::result_of::vertex_id<ElementType>::type VertexIDType;
00506         typedef typename viennagrid::result_of::vertex_handle<ElementType>::type VertexHandleType;
00507 
00508         typedef typename viennagrid::result_of::line_id<ElementType>::type LineIDType;
00509         typedef typename viennagrid::result_of::line_handle<ElementType>::type LineHandleType;
00510 
00511         typedef typename viennagrid::result_of::const_line_range<ElementType>::type ConstLineOnElementRangeType;
00512         typedef typename viennagrid::result_of::iterator<ConstLineOnElementRangeType>::type ConstLineOnElementIteratorType;
00513 
00514         std::map<VertexIDType, VertexHandleType> vertex_map;
00515         std::map<LineIDType, LineHandleType> line_map;
00516 
00517 
00518         for (ElementIteratorT eit = begin; eit != end; ++eit)
00519         {
00520           ConstLineOnElementRangeType lines_on_element( *eit );
00521           std::vector<LineHandleType> line_handles( lines_on_element.size() );
00522 
00523           unsigned int line_index = 0;
00524           for (ConstLineOnElementIteratorType lit = lines_on_element.begin(); lit != lines_on_element.end(); ++lit, ++line_index)
00525           {
00526             typename std::map<LineIDType, LineHandleType>::iterator it1 = line_map.find( lit->id() );
00527             if (it1 == line_map.end())
00528             {
00529               VertexHandleType vtx_handles[2];
00530               for (unsigned int vertex_index = 0; vertex_index != 2; ++vertex_index)
00531               {
00532                 typename std::map<VertexIDType, VertexHandleType>::iterator it2 = vertex_map.find( viennagrid::vertices(*lit)[vertex_index].id() );
00533                 if (it2 == vertex_map.end())
00534                 {
00535                   vtx_handles[vertex_index] = viennagrid::make_unique_vertex( output_mesh, viennagrid::point( viennagrid::vertices(*lit)[vertex_index] ), tolerance );
00536                   vertex_map[ viennagrid::vertices(*lit)[vertex_index].id() ] = vtx_handles[vertex_index];
00537                 }
00538                 else
00539                   vtx_handles[vertex_index] = it2->second;
00540               }
00541 
00542               line_handles[line_index] = viennagrid::make_line( output_mesh, vtx_handles[0], vtx_handles[1] );
00543               line_map[lit->id()] = line_handles[line_index];
00544             }
00545             else
00546               line_handles[line_index] = it1->second;
00547           }
00548 
00549           viennagrid::make_plc( output_mesh, line_handles.begin(), line_handles.end() );
00550         }
00551       }
00552 
00554       template<typename InputMeshOrSegmentHandleT, typename ElementHandleIteratorT, typename OutputMeshOrSegmentHandleT>
00555       static void copy_elements_by_handle(InputMeshOrSegmentHandleT const & input_mesh,
00556                                 ElementHandleIteratorT const & begin, ElementHandleIteratorT const & end,
00557                                 OutputMeshOrSegmentHandleT & output_mesh,
00558                                 typename viennagrid::result_of::coord<OutputMeshOrSegmentHandleT>::type tolerance )
00559       {
00560         typedef typename std::iterator_traits<ElementHandleIteratorT>::value_type ElementHandleType;
00561         typedef typename viennagrid::detail::result_of::value_type<ElementHandleType>::type ElementType;
00562 
00563         typedef typename viennagrid::result_of::vertex_id<ElementType>::type VertexIDType;
00564         typedef typename viennagrid::result_of::vertex_handle<ElementType>::type VertexHandleType;
00565 
00566         typedef typename viennagrid::result_of::line_id<ElementType>::type LineIDType;
00567         typedef typename viennagrid::result_of::line_handle<ElementType>::type LineHandleType;
00568 
00569         typedef typename viennagrid::result_of::const_line_range<ElementType>::type ConstLineOnElementRangeType;
00570         typedef typename viennagrid::result_of::iterator<ConstLineOnElementRangeType>::type ConstLineOnElementIteratorType;
00571 
00572         std::map<VertexIDType, VertexHandleType> vertex_map;
00573         std::map<LineIDType, LineHandleType> line_map;
00574 
00575 
00576         for (ElementHandleIteratorT eit = begin; eit != end; ++eit)
00577         {
00578           ConstLineOnElementRangeType lines_on_element( viennagrid::dereference_handle(input_mesh, *eit) );
00579           std::vector<LineHandleType> line_handles( lines_on_element.size() );
00580 
00581           unsigned int line_index = 0;
00582           for (ConstLineOnElementIteratorType lit = lines_on_element.begin(); lit != lines_on_element.end(); ++lit, ++line_index)
00583           {
00584             typename std::map<LineIDType, LineHandleType>::iterator it1 = line_map.find( lit->id() );
00585             if (it1 == line_map.end())
00586             {
00587               VertexHandleType vtx_handles[2];
00588               for (unsigned int vertex_index = 0; vertex_index != 2; ++vertex_index)
00589               {
00590                 typename std::map<VertexIDType, VertexHandleType>::iterator it2 = vertex_map.find( viennagrid::vertices(*lit)[vertex_index].id() );
00591                 if (it2 == vertex_map.end())
00592                 {
00593                   vtx_handles[vertex_index] = viennagrid::make_unique_vertex( output_mesh, viennagrid::point( viennagrid::vertices(*lit)[vertex_index] ), tolerance );
00594                   vertex_map[ viennagrid::vertices(*lit)[vertex_index].id() ] = vtx_handles[vertex_index];
00595                 }
00596                 else
00597                   vtx_handles[vertex_index] = it2->second;
00598               }
00599 
00600               line_handles[line_index] = viennagrid::make_line( output_mesh, vtx_handles[0], vtx_handles[1] );
00601               line_map[lit->id()] = line_handles[line_index];
00602             }
00603             else
00604               line_handles[line_index] = it1->second;
00605           }
00606 
00607           viennagrid::make_plc( output_mesh, line_handles.begin(), line_handles.end() );
00608         }
00609 
00610       }
00611     };
00612   }
00613 
00614 
00615 
00616 
00617   // doxygen doku in forwards.hpp
00618   template<typename ElementIteratorT, typename OutputMeshOrSegmentHandleT>
00619   void copy_elements(ElementIteratorT const & begin, ElementIteratorT const & end,
00620                      OutputMeshOrSegmentHandleT & output_mesh,
00621                      typename viennagrid::result_of::coord<OutputMeshOrSegmentHandleT>::type tolerance )
00622   {
00623     typedef typename std::iterator_traits<ElementIteratorT>::value_type ElementType;
00624     typedef typename viennagrid::result_of::element_tag<ElementType>::type ElementTagType;
00625 
00626     detail::copy_elements_impl<ElementTagType>::copy_elements(begin, end, output_mesh, tolerance);
00627   }
00628 
00629   // doxygen doku in forwards.hpp
00630   template<typename ElementIteratorT, typename OutputMeshOrSegmentHandleT>
00631   void copy_elements(ElementIteratorT const & begin, ElementIteratorT const & end,
00632                      OutputMeshOrSegmentHandleT & output_mesh )
00633   {
00634     copy_elements(begin, end, output_mesh, -1.0);
00635   }
00636 
00637 
00638   // doxygen doku in forwards.hpp
00639   template<typename InputMeshOrSegmentHandleT, typename ElementHandleIteratorT, typename OutputMeshOrSegmentHandleT>
00640   void copy_elements_by_handle(InputMeshOrSegmentHandleT const & input_mesh,
00641                             ElementHandleIteratorT const & begin, ElementHandleIteratorT const & end,
00642                             OutputMeshOrSegmentHandleT & output_mesh,
00643                             typename viennagrid::result_of::coord<OutputMeshOrSegmentHandleT>::type tolerance )
00644   {
00645     typedef typename std::iterator_traits<ElementHandleIteratorT>::value_type ElementHandleType;
00646     typedef typename viennagrid::detail::result_of::value_type<ElementHandleType>::type ElementType;
00647     typedef typename viennagrid::result_of::element_tag<ElementType>::type ElementTagType;
00648 
00649     detail::copy_elements_impl<ElementTagType>::copy_elements_by_handle(input_mesh, begin, end, output_mesh, tolerance);
00650   }
00651 
00652   // doxygen doku in forwards.hpp
00653   template<typename InputMeshOrSegmentHandleT, typename ElementHandleIteratorT, typename OutputMeshOrSegmentHandleT>
00654   void copy_elements_by_handle(InputMeshOrSegmentHandleT const & input_mesh,
00655                             ElementHandleIteratorT const & begin, ElementHandleIteratorT const & end,
00656                             OutputMeshOrSegmentHandleT & output_mesh )
00657   {
00658     copy_elements_by_handle(input_mesh, begin, end, output_mesh, -1.0);
00659   }
00660 
00661 }
00662 
00663 
00664 
00665 #endif