Logo  0.95.0-final
Finite Element Embedded Library and Language in C++
Feel++ Feel++ on Github Feel++ community
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
filterfromvtk.hpp
Go to the documentation of this file.
1 /* -*- mode: c++; coding: utf-8; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; show-trailing-whitespace: t -*- vim:fenc=utf-8:ft=tcl:et:sw=4:ts=4:sts=4
2 
3  This file is part of the Feel library
4 
5  Author(s): Gilles Steiner <gilles.steiner@epfl.ch>
6  Christophe Prud'homme <christophe.prudhomme@feelpp.org>
7  Date: 2004-11-18
8 
9  Copyright (C) 2004 EPFL
10  Copyright (C) 2012 Universite de Strasbbourg
11  Copyright (C) 2006-2012 Feel++ Consortium
12 
13  This library is free software; you can redistribute it and/or
14  modify it under the terms of the GNU Lesser General Public
15  License as published by the Free Software Foundation; either
16  version 3.0 of the License, or (at your option) any later version.
17 
18  This library is distributed in the hope that it will be useful,
19  but WITHOUT ANY WARRANTY; without even the implied warranty of
20  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
21  Lesser General Public License for more details.
22 
23  You should have received a copy of the GNU Lesser General Public
24  License along with this library; if not, write to the Free Software
25  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
26 */
33 #ifndef __filter_H
34 #define __filter_H 1
35 
37 //#include <feel/feelmesh/mesh1d.hpp>
38 #include <feel/feelmesh/mesh2d.hpp>
39 #include <feel/feelmesh/mesh3d.hpp>
40 #include <feel/feeldiscr/mesh.hpp>
41 
42 #include <feel/feelcore/feel.hpp>
43 
44 
45 #if defined(FEELPP_HAS_VTK)
46 // Vtk header files
47 //#include "vtkPointSet.h"
48 #include "vtkPolyData.h"
49 #include "vtkUnstructuredGrid.h"
50 #include "vtkCell.h"
51 
52 #endif /* FEELPP_HAS_VTK */
53 
54 namespace Feel
55 {
56 
65 template<typename MeshType>
67  :
68 public VisitorBase,
69 public Visitor<MeshType>
70 {
71 
72 public:
73 
74  static const uint16_type nDim = MeshType::nDim;
75  BOOST_STATIC_ASSERT( nDim == 2 || nDim == 3 );
76 
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)
87 
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;
91 
92 
93 #endif /* FEELPP_HAS_VTK */
94 
95 
99 
100 #if defined(FEELPP_HAS_VTK)
101  FilterFromVtk( vtkmesh_type* __vtkmesh )
102  {
103  M_vtkmesh = vtkmesh_type::New();
104  M_vtkmesh->CopyStructure( __vtkmesh );
105  }
106 #else
107  FilterFromVtk()
108  {}
109 #endif /* FEELPP_HAS_VTK */
110 
111  ~FilterFromVtk()
112  {
113 #if defined(FEELPP_HAS_VTK)
114  M_vtkmesh->Delete();
115 #endif /* FEELPP_HAS_VTK */
116  }
117 
118 #if defined(FEELPP_HAS_VTK)
119  vtkmesh_type* getVtkMesh()
120  {
121  return M_vtkmesh;
122  }
123 #endif /* FEELPP_HAS_VTK */
124 
126 
130 
131  void visit( mesh_type* mesh )
132  {
133  visit( mesh, mpl::int_<nDim>() );
134  }
135 
136 
138 
139 protected :
140 
141 #if defined(FEELPP_HAS_VTK)
142  vtkmesh_type* M_vtkmesh;
143 #endif /* FEELPP_HAS_VTK */
144 
145 private:
146 
147  void visit( mesh_type* mesh, mpl::int_<2> );
148  //void visit( mesh_type* mesh, mpl::int_<3> );
149 };
150 
159 template<typename MeshType>
161  :
162 public VisitorBase,
163 public Visitor<MeshType>
164 {
165 
166 public:
167 
168  static const uint16_type nDim = MeshType::nDim;
169  BOOST_STATIC_ASSERT( nDim == 2 || nDim == 3 );
170 
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;
184 
185 #endif /* FEELPP_HAS_VTK */
186 
187 
191 
192 #if defined(FEELPP_HAS_VTK)
193  FilterFromVtk3D( vtkmesh_type* __vtkmesh )
194  {
195  M_vtkmesh = vtkmesh_type::New();
196  M_vtkmesh->CopyStructure( __vtkmesh );
197  }
198 #else
200  {}
201 #endif /* FEELPP_HAS_VTK */
202 
203  ~FilterFromVtk3D()
204  {
205 #if defined(FEELPP_HAS_VTK)
206  M_vtkmesh->Delete();
207 #endif /* FEELPP_HAS_VTK */
208  }
209 
210 #if defined(FEELPP_HAS_VTK)
211  vtkmesh_type* getVtkMesh()
212  {
213  return M_vtkmesh;
214  }
215 #endif /* FEELPP_HAS_VTK */
216 
218 
222 
223  void visit( mesh_type* mesh )
224  {
225  visit( mesh, mpl::int_<nDim>() );
226  }
227 
228 
230 
231 protected :
232 
233 #if defined(FEELPP_HAS_VTK)
234  vtkmesh_type* M_vtkmesh;
235 #endif /* FEELPP_HAS_VTK */
236 
237 private:
238 
239  //void visit( mesh_type* mesh, mpl::int_<2> );
240  void visit( mesh_type* mesh, mpl::int_<3> );
241 };
242 
243 template<typename MeshType>
244 void
245 FilterFromVtk<MeshType>::visit( mesh_type* mesh, mpl::int_<2> )
246 {
247  detail::ignore_unused_variable_warning( mesh );
248 #if defined(FEELPP_HAS_VTK)
249  // std::cout <<"Start of mesh conversion !" << std::endl;
250 
251  vtkPolyData * _vtkMesh = this->getVtkMesh();
252 
253  uint __n = _vtkMesh->GetNumberOfPoints(); // Number of nodes
254 
255  DVLOG( 2 ) <<"Number of points : "<< __n << "\n";
256 
257  uint __nele = _vtkMesh->GetNumberOfPolys(); // Number of elements
258 
259  DVLOG( 2 ) <<"Number of elements : "<< __nele << "\n";
260 
261  // add the points to the mesh
262 
263  for ( uint __i = 0; __i < __n; ++__i )
264  {
265  node_type __nd( 2 );
266  __nd[0] = _vtkMesh->GetPoint( __i )[0];
267  __nd[1] = _vtkMesh->GetPoint( __i )[1];
268  point_type __pt( __i,__nd, false );
269  __pt.marker() = 0;
270 
271  if ( __nd[0] == -1 || __nd[1] == -1 || __nd[0] + __nd[1] == 0 )
272 
273  {
274  __pt.setOnBoundary( true );
275 
276  }
277 
278  else
279  {
280  __pt.setOnBoundary( false );
281  }
282 
283  mesh->addPoint( __pt );
284  }
285 
286 #if 0
287  size_type n_faces = 0;
288 
289  // Add Boundary faces
290 
291  face_type* pf0 = new face_type;
292 
293  pf0->setMarker( 1 );
294  pf0->setPoint( 0, mesh->point( 1 ) );
295  pf0->setPoint( 1, mesh->point( 2 ) );
296 
297  pf0->setId( n_faces++ );
298  pf0->setOnBoundary( true );
299  mesh->addFace( *pf0 );
300 
301  delete pf0;
302 
303  face_type* pf1 = new face_type;
304 
305  pf1->setMarker( 1 );
306  pf1->setPoint( 0, mesh->point( 2 ) );
307  pf1->setPoint( 1, mesh->point( 0 ) );
308 
309  pf1->setId( n_faces++ );
310  pf1->setOnBoundary( true );
311  mesh->addFace( *pf1 );
312 
313  delete pf1;
314 
315  face_type* pf2 = new face_type;
316 
317  pf2->setMarker( 1 );
318  pf2->setPoint( 0, mesh->point( 0 ) );
319  pf2->setPoint( 1, mesh->point( 1 ) );
320 
321  pf2->setId( n_faces++ );
322  pf2->setOnBoundary( true );
323  mesh->addFace( *pf2 );
324 
325  delete pf2;
326  FEELPP_ASSERT( n_faces == mesh->numFaces() )( n_faces )( mesh->numFaces() ).error( "invalid face container size" );
327 
328 #endif
329  // add the elements to the mesh
330 
331  for ( uint __i = 0; __i < __nele; ++__i )
332  {
333  DVLOG( 2 ) << "[FilterFromVtk] element " << __i << "\n";
334  // Here we only have triangular elements of order 1
335 
336  element_type * pf = new element_type;
337 
338  pf->setId( __i );
339  pf->setMarker( 0 );
340 
341  // Warning : Vtk orientation is not the same as Feel orientation !
342 
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 );
350  delete pf;
351  }
352 
353 
354 
355  DVLOG( 2 ) <<"[FilterFromVtk] done with element accumulation !\n";
356 
357  mesh->setNumVertices( __n );
358 
359  DVLOG( 2 ) <<"[FilterFromVtk] Face Update !\n";
360 
361 
362  DVLOG( 2 ) <<"[FilterFromVtk] Face Update Successful !\n";
363 
364 #else
365  std::cerr << "The library was not compiled with vtk support\n";
366 #endif /* FEELPP_HAS_VTK */
367 
368 }
369 
370 template<typename MeshType>
371 void
372 FilterFromVtk3D<MeshType>::visit( mesh_type* mesh, mpl::int_<3> )
373 {
374  detail::ignore_unused_variable_warning( mesh );
375 #if defined(FEELPP_HAS_VTK)
376  // std::cout <<"Start of mesh conversion !" << std::endl;
377 
378  vtkUnstructuredGrid * _vtkMesh = this->getVtkMesh();
379 
380  uint __n = _vtkMesh->GetNumberOfPoints(); // Nbre of nodes
381 
382  DVLOG( 2 ) <<"Number of points : "<< __n << "\n";
383 
384  uint __nele = _vtkMesh->GetNumberOfCells(); // Nbre of elements
385 
386  DVLOG( 2 ) <<"Number of elements : "<< __nele << "\n";
387 
388  // add the points to the mesh
389 
390  for ( uint __i = 0; __i < __n; ++__i )
391  {
392  node_type __nd( 3 );
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 );
397 
398  if ( __nd[0] == -1 || __nd[1] == -1 || __nd[2] == -1 || __nd[0] + __nd[1] + __nd[2] == 0 )
399 
400  {
401  __pt.setOnBoundary( true );
402  __pt.marker() = 0;
403  }
404 
405  else
406  {
407  __pt.setOnBoundary( false );
408  __pt.marker() = 1;
409  }
410 
411  mesh->addPoint( __pt );
412  }
413 
414  DVLOG( 2 ) << "[FilterFromVtk3D] mesh np = " << mesh->numPoints() << "\n";
415  FEELPP_ASSERT( mesh->numPoints() == __n )( __n )( mesh->numPoints() ).error( "invalid number of points" );
416 
417 #if 0
418  size_type n_faces = 0;
419 
420 
421  // Add Boundary faces
422  face_type* pf0 = new face_type;
423 
424  pf0->setMarker( 1 );
425  pf0->setPoint( 0, mesh->point( 1 ) );
426  pf0->setPoint( 1, mesh->point( 3 ) );
427  pf0->setPoint( 2, mesh->point( 2 ) );
428 
429 
430  pf0->setId( n_faces++ );
431  pf0->setOnBoundary( true );
432  mesh->addFace( *pf0 );
433 
434  delete pf0;
435 
436  face_type* pf1 = new face_type;
437 
438  pf1->setMarker( 1 );
439  pf1->setPoint( 0, mesh->point( 2 ) );
440  pf1->setPoint( 1, mesh->point( 3 ) );
441  pf1->setPoint( 2, mesh->point( 0 ) );
442 
443 
444  pf1->setId( n_faces++ );
445  pf1->setOnBoundary( true );
446  mesh->addFace( *pf1 );
447 
448  delete pf1;
449 
450  face_type* pf2 = new face_type;
451 
452  pf2->setMarker( 1 );
453  pf2->setPoint( 0, mesh->point( 3 ) );
454  pf2->setPoint( 1, mesh->point( 1 ) );
455  pf2->setPoint( 2, mesh->point( 0 ) );
456 
457  pf2->setId( n_faces++ );
458  pf2->setOnBoundary( true );
459  mesh->addFace( *pf2 );
460 
461  delete pf2;
462 
463  face_type* pf3 = new face_type;
464 
465  pf3->setMarker( 1 );
466  pf3->setPoint( 0, mesh->point( 0 ) );
467  pf3->setPoint( 1, mesh->point( 1 ) );
468  pf3->setPoint( 2, mesh->point( 2 ) );
469 
470  pf3->setId( n_faces++ );
471  pf3->setOnBoundary( true );
472  mesh->addFace( *pf3 );
473 
474  delete pf3;
475  FEELPP_ASSERT( n_faces == mesh->numFaces() )( n_faces )( mesh->numFaces() ).error( "invalid face container size" );
476 #endif
477  // add the elements to the mesh
478 
479  for ( uint __i = 0; __i < __nele; ++__i )
480  {
481  // Here we only have triangular elements of order 1
482 
483  element_type * pf = new element_type;
484 
485  pf->setId( __i );
486  pf->setMarker( 0 );
487 
488  // Warning : Vtk orientation is not the same as Feel orientation !
489 
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 ) ) );
494 
495  mesh->addElement( *pf );
496  delete pf;
497  }
498 
499 
500 
501 
502  DVLOG( 2 ) <<"[FilterFromVtk] done with element accumulation !\n";
503 
504  mesh->setNumVertices( __n );
505 
506  DVLOG( 2 ) <<"[FilterFromVtk] Face Update !\n";
507 
508  // do not renumber the mesh entities
509  //mesh->updateForUse( MESH_ALL_COMPONENTS & (~MESH_RENUMBER) );
510 
511  DVLOG( 2 ) <<"[FilterFromVtk] Face Update Successful !\n";
512 
513 #else
514  std::cerr << "The library was not compiled with vtk support\n";
515 #endif /* FEELPP_HAS_VTK */
516 
517 }
518 
519 } // Feel
520 
521 #endif /* __filter_H */

Generated on Fri Oct 25 2013 14:24:11 for Feel++ by doxygen 1.8.4