42 #include <feel/feelcore/feel.hpp>
45 #if defined(FEELPP_HAS_VTK)
48 #include "vtkPolyData.h"
49 #include "vtkUnstructuredGrid.h"
65 template<
typename MeshType>
74 static const uint16_type nDim = MeshType::nDim;
75 BOOST_STATIC_ASSERT( nDim == 2 || nDim == 3 );
80 typedef MeshType mesh_type;
81 typedef typename mesh_type::point_type point_type;
82 typedef typename point_type::node_type node_type;
83 typedef typename mesh_type::edge_type edge_type;
84 typedef typename mesh_type::face_type face_type;
85 typedef typename mesh_type::element_type element_type;
86 #if defined(FEELPP_HAS_VTK)
88 typedef typename mpl::if_<mpl::equal_to<mpl::int_<MeshType::nDim>,mpl::int_<2> >,
89 mpl::identity<vtkPolyData>,
90 mpl::identity<vtkUnstructuredGrid> >::type::type vtkmesh_type;
100 #if defined(FEELPP_HAS_VTK)
103 M_vtkmesh = vtkmesh_type::New();
104 M_vtkmesh->CopyStructure( __vtkmesh );
113 #if defined(FEELPP_HAS_VTK)
118 #if defined(FEELPP_HAS_VTK)
119 vtkmesh_type* getVtkMesh()
133 visit( mesh, mpl::int_<nDim>() );
141 #if defined(FEELPP_HAS_VTK)
142 vtkmesh_type* M_vtkmesh;
147 void visit( mesh_type* mesh, mpl::int_<2> );
159 template<
typename MeshType>
168 static const uint16_type nDim = MeshType::nDim;
169 BOOST_STATIC_ASSERT( nDim == 2 || nDim == 3 );
174 typedef MeshType mesh_type;
175 typedef typename mesh_type::point_type point_type;
176 typedef typename point_type::node_type node_type;
177 typedef typename mesh_type::edge_type edge_type;
178 typedef typename mesh_type::face_type face_type;
179 typedef typename mesh_type::element_type element_type;
180 #if defined(FEELPP_HAS_VTK)
181 typedef typename mpl::if_<mpl::equal_to<mpl::int_<MeshType::nDim>,mpl::int_<2> >,
182 mpl::identity<vtkPolyData>,
183 mpl::identity<vtkUnstructuredGrid> >::type::type vtkmesh_type;
192 #if defined(FEELPP_HAS_VTK)
195 M_vtkmesh = vtkmesh_type::New();
196 M_vtkmesh->CopyStructure( __vtkmesh );
205 #if defined(FEELPP_HAS_VTK)
210 #if defined(FEELPP_HAS_VTK)
211 vtkmesh_type* getVtkMesh()
225 visit( mesh, mpl::int_<nDim>() );
233 #if defined(FEELPP_HAS_VTK)
234 vtkmesh_type* M_vtkmesh;
240 void visit( mesh_type* mesh, mpl::int_<3> );
243 template<
typename MeshType>
247 detail::ignore_unused_variable_warning( mesh );
248 #if defined(FEELPP_HAS_VTK)
251 vtkPolyData * _vtkMesh = this->getVtkMesh();
253 uint __n = _vtkMesh->GetNumberOfPoints();
255 DVLOG( 2 ) <<
"Number of points : "<< __n <<
"\n";
257 uint __nele = _vtkMesh->GetNumberOfPolys();
259 DVLOG( 2 ) <<
"Number of elements : "<< __nele <<
"\n";
263 for ( uint __i = 0; __i < __n; ++__i )
266 __nd[0] = _vtkMesh->GetPoint( __i )[0];
267 __nd[1] = _vtkMesh->GetPoint( __i )[1];
268 point_type __pt( __i,__nd,
false );
271 if ( __nd[0] == -1 || __nd[1] == -1 || __nd[0] + __nd[1] == 0 )
274 __pt.setOnBoundary(
true );
280 __pt.setOnBoundary(
false );
283 mesh->addPoint( __pt );
291 face_type* pf0 =
new face_type;
294 pf0->setPoint( 0, mesh->point( 1 ) );
295 pf0->setPoint( 1, mesh->point( 2 ) );
297 pf0->setId( n_faces++ );
298 pf0->setOnBoundary(
true );
299 mesh->addFace( *pf0 );
303 face_type* pf1 =
new face_type;
306 pf1->setPoint( 0, mesh->point( 2 ) );
307 pf1->setPoint( 1, mesh->point( 0 ) );
309 pf1->setId( n_faces++ );
310 pf1->setOnBoundary(
true );
311 mesh->addFace( *pf1 );
315 face_type* pf2 =
new face_type;
318 pf2->setPoint( 0, mesh->point( 0 ) );
319 pf2->setPoint( 1, mesh->point( 1 ) );
321 pf2->setId( n_faces++ );
322 pf2->setOnBoundary(
true );
323 mesh->addFace( *pf2 );
326 FEELPP_ASSERT( n_faces == mesh->numFaces() )( n_faces )( mesh->numFaces() ).error(
"invalid face container size" );
331 for ( uint __i = 0; __i < __nele; ++__i )
333 DVLOG( 2 ) <<
"[FilterFromVtk] element " << __i <<
"\n";
336 element_type * pf =
new element_type;
343 pf->setPoint( 0, mesh->point( _vtkMesh->GetCell( __i )->GetPointId( 0 ) ) );
344 pf->setPoint( 1, mesh->point( _vtkMesh->GetCell( __i )->GetPointId( 1 ) ) );
345 pf->setPoint( 2, mesh->point( _vtkMesh->GetCell( __i )->GetPointId( 2 ) ) );
346 DVLOG( 2 ) <<
"[FilterFromVtk] point 0 " << pf->point( 0 ).node() <<
" global id: " << pf->point( 0 ).id() <<
"\n"
347 <<
"[FilterFromVtk] point 1 " << pf->point( 1 ).node() <<
" global id: " << pf->point( 1 ).id() <<
"\n"
348 <<
"[FilterFromVtk] point 2 " << pf->point( 2 ).node() <<
" global id: " << pf->point( 2 ).id() <<
"\n";
349 mesh->addElement( *pf );
355 DVLOG( 2 ) <<
"[FilterFromVtk] done with element accumulation !\n";
357 mesh->setNumVertices( __n );
359 DVLOG( 2 ) <<
"[FilterFromVtk] Face Update !\n";
362 DVLOG( 2 ) <<
"[FilterFromVtk] Face Update Successful !\n";
365 std::cerr <<
"The library was not compiled with vtk support\n";
370 template<
typename MeshType>
374 detail::ignore_unused_variable_warning( mesh );
375 #if defined(FEELPP_HAS_VTK)
378 vtkUnstructuredGrid * _vtkMesh = this->getVtkMesh();
380 uint __n = _vtkMesh->GetNumberOfPoints();
382 DVLOG( 2 ) <<
"Number of points : "<< __n <<
"\n";
384 uint __nele = _vtkMesh->GetNumberOfCells();
386 DVLOG( 2 ) <<
"Number of elements : "<< __nele <<
"\n";
390 for ( uint __i = 0; __i < __n; ++__i )
393 __nd[0] = _vtkMesh->GetPoint( __i )[0];
394 __nd[1] = _vtkMesh->GetPoint( __i )[1];
395 __nd[2] = _vtkMesh->GetPoint( __i )[2];
396 point_type __pt( __i,__nd,
false );
398 if ( __nd[0] == -1 || __nd[1] == -1 || __nd[2] == -1 || __nd[0] + __nd[1] + __nd[2] == 0 )
401 __pt.setOnBoundary(
true );
407 __pt.setOnBoundary(
false );
411 mesh->addPoint( __pt );
414 DVLOG( 2 ) <<
"[FilterFromVtk3D] mesh np = " << mesh->numPoints() <<
"\n";
415 FEELPP_ASSERT( mesh->numPoints() == __n )( __n )( mesh->numPoints() ).error(
"invalid number of points" );
422 face_type* pf0 =
new face_type;
425 pf0->setPoint( 0, mesh->point( 1 ) );
426 pf0->setPoint( 1, mesh->point( 3 ) );
427 pf0->setPoint( 2, mesh->point( 2 ) );
430 pf0->setId( n_faces++ );
431 pf0->setOnBoundary(
true );
432 mesh->addFace( *pf0 );
436 face_type* pf1 =
new face_type;
439 pf1->setPoint( 0, mesh->point( 2 ) );
440 pf1->setPoint( 1, mesh->point( 3 ) );
441 pf1->setPoint( 2, mesh->point( 0 ) );
444 pf1->setId( n_faces++ );
445 pf1->setOnBoundary(
true );
446 mesh->addFace( *pf1 );
450 face_type* pf2 =
new face_type;
453 pf2->setPoint( 0, mesh->point( 3 ) );
454 pf2->setPoint( 1, mesh->point( 1 ) );
455 pf2->setPoint( 2, mesh->point( 0 ) );
457 pf2->setId( n_faces++ );
458 pf2->setOnBoundary(
true );
459 mesh->addFace( *pf2 );
463 face_type* pf3 =
new face_type;
466 pf3->setPoint( 0, mesh->point( 0 ) );
467 pf3->setPoint( 1, mesh->point( 1 ) );
468 pf3->setPoint( 2, mesh->point( 2 ) );
470 pf3->setId( n_faces++ );
471 pf3->setOnBoundary(
true );
472 mesh->addFace( *pf3 );
475 FEELPP_ASSERT( n_faces == mesh->numFaces() )( n_faces )( mesh->numFaces() ).error(
"invalid face container size" );
479 for ( uint __i = 0; __i < __nele; ++__i )
483 element_type * pf =
new element_type;
490 pf->setPoint( 0, mesh->point( _vtkMesh->GetCell( __i )->GetPointId( 0 ) ) );
491 pf->setPoint( 1, mesh->point( _vtkMesh->GetCell( __i )->GetPointId( 1 ) ) );
492 pf->setPoint( 2, mesh->point( _vtkMesh->GetCell( __i )->GetPointId( 2 ) ) );
493 pf->setPoint( 3, mesh->point( _vtkMesh->GetCell( __i )->GetPointId( 3 ) ) );
495 mesh->addElement( *pf );
502 DVLOG( 2 ) <<
"[FilterFromVtk] done with element accumulation !\n";
504 mesh->setNumVertices( __n );
506 DVLOG( 2 ) <<
"[FilterFromVtk] Face Update !\n";
511 DVLOG( 2 ) <<
"[FilterFromVtk] Face Update Successful !\n";
514 std::cerr <<
"The library was not compiled with vtk support\n";