[Getfem-commits] r5105 - in /trunk/getfem/src: getfem/getfem_export.h getfem_export.cc

[logari81]

Author: logari81
Date: Tue Oct 20 23:52:15 2015
New Revision: 5105

URL: http://svn.gna.org/viewcvs/getfem?rev=5105&view=rev
Log:
support 9-node quads and 27-node hexas in export_to_vtk (can be read with VTK 
5.1 and above)

Modified:
    trunk/getfem/src/getfem/getfem_export.h
    trunk/getfem/src/getfem_export.cc

Modified: trunk/getfem/src/getfem/getfem_export.h
URL: 
http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem/getfem_export.h?rev=5105&r1=5104&r2=5105&view=diff
==============================================================================
--- trunk/getfem/src/getfem/getfem_export.h     (original)
+++ trunk/getfem/src/getfem/getfem_export.h     Tue Oct 20 23:52:15 2015
@@ -79,13 +79,23 @@
     enum { EMPTY, HEADER_WRITTEN, STRUCTURE_WRITTEN, IN_CELL_DATA,
            IN_POINT_DATA } state;
   public:
-    typedef enum { VTK_VERTEX = 1, VTK_LINE = 3, VTK_QUADRATIC_EDGE = 21,
-                   VTK_TRIANGLE = 5, VTK_QUADRATIC_TRIANGLE = 22,
-                   VTK_PIXEL = 8, VTK_QUAD = 9, VTK_QUADRATIC_QUAD = 23,
-                   VTK_TETRA = 10, VTK_QUADRATIC_TETRA = 24,
-                   VTK_WEDGE = 13, /*VTK_QUADRATIC_WEDGE = 26,*/
-                   VTK_VOXEL = 11, VTK_HEXAHEDRON = 12,
-                   VTK_QUADRATIC_HEXAHEDRON = 25 } vtk_cell_type;
+    typedef enum { VTK_VERTEX = 1,
+                   VTK_LINE = 3,
+                   VTK_TRIANGLE = 5,
+                   VTK_PIXEL = 8,
+                   VTK_QUAD = 9,
+                   VTK_TETRA = 10,
+                   VTK_VOXEL = 11,
+                   VTK_HEXAHEDRON = 12,
+                   VTK_WEDGE = 13,
+                   VTK_QUADRATIC_EDGE = 21,
+                   VTK_QUADRATIC_TRIANGLE = 22,
+                   VTK_QUADRATIC_QUAD = 23,
+                   VTK_QUADRATIC_TETRA = 24,
+                   VTK_QUADRATIC_HEXAHEDRON = 25,
+                   /*VTK_QUADRATIC_WEDGE = 26,*/
+                   VTK_BIQUADRATIC_QUAD = 28,
+                   VTK_TRIQUADRATIC_HEXAHEDRON = 29 } vtk_cell_type;
     vtk_export(const std::string& fname, bool ascii_ = false);
     vtk_export(std::ostream &os_, bool ascii_ = false);
 

Modified: trunk/getfem/src/getfem_export.cc
URL: 
http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem_export.cc?rev=5105&r1=5104&r2=5105&view=diff
==============================================================================
--- trunk/getfem/src/getfem_export.cc   (original)
+++ trunk/getfem/src/getfem_export.cc   Tue Oct 20 23:52:15 2015
@@ -43,14 +43,16 @@
       bgeot::sc(dm[vtk_export::VTK_QUADRATIC_TRIANGLE]) = 0, 2, 5, 1, 4, 3;
       bgeot::sc(dm[vtk_export::VTK_QUAD]) = 0, 1, 3, 2;
       bgeot::sc(dm[vtk_export::VTK_PIXEL]) = 0, 1, 2, 3;
-      bgeot::sc(dm[vtk_export::VTK_QUADRATIC_QUAD]) = 0, 2, 8, 6, 1, 5, 7, 3;
+      bgeot::sc(dm[vtk_export::VTK_QUADRATIC_QUAD]) = 0, 2, 7, 5, 1, 4, 6, 3;
       bgeot::sc(dm[vtk_export::VTK_TETRA]) = 0, 1, 2, 3;
       bgeot::sc(dm[vtk_export::VTK_QUADRATIC_TETRA]) = 0, 2, 5, 9, 1, 4, 3, 6, 
7, 8;
       bgeot::sc(dm[vtk_export::VTK_WEDGE]) = 0, 1, 2, 3, 4, 5;
       //bgeot::sc(dm[vtk_export::VTK_QUADRATIC_WEDGE]) = 0, 6, 1, 7, 2, 8, 3, 
4, 5;
       bgeot::sc(dm[vtk_export::VTK_VOXEL]) = 0, 1, 2, 3, 4, 5, 6, 7;
       bgeot::sc(dm[vtk_export::VTK_HEXAHEDRON]) = 0, 1, 3, 2, 4, 5, 7, 6;
-      bgeot::sc(dm[vtk_export::VTK_QUADRATIC_HEXAHEDRON]) = 0, 2, 8, 6, 18, 
20, 26, 24, 1, 5, 7, 3, 19, 23, 25, 21, 9, 11, 17, 15;
+      bgeot::sc(dm[vtk_export::VTK_QUADRATIC_HEXAHEDRON]) = 0, 2, 7, 5, 12, 
14, 19, 17, 1, 4, 6, 3, 13, 16, 18, 15, 8, 9, 11, 10;
+      bgeot::sc(dm[vtk_export::VTK_BIQUADRATIC_QUAD]) = 0, 2, 8, 6, 1, 5, 7, 
3, 4;
+      bgeot::sc(dm[vtk_export::VTK_TRIQUADRATIC_HEXAHEDRON]) = 0, 2, 8, 6, 18, 
20, 26, 24, 1, 5, 7, 3, 19, 23, 25, 21, 9, 11, 17, 15, 12, 14, 10, 16, 4, 22;
     }
     return dm[t];
   }
@@ -144,7 +146,7 @@
     if (&mf != pmf.get()) pmf.reset(new mesh_fem(mf.linked_mesh(),1));
 
     /* initialize pmf with finite elements suitable for VTK (which only knows
-       isoparametric FEMs of order 1 and 2 (with inner points missing)) */
+       isoparametric FEMs of order 1 and 2) */
     for (dal::bv_visitor cv(mf.convex_index()); !cv.finished(); ++cv) {
       bgeot::pgeometric_trans pgt = mf.linked_mesh().trans_of_convex(cv);
       pfem pf = mf.fem_of_element(cv);
@@ -154,15 +156,21 @@
         /* could be a better test for discontinuity .. */
         if (!dof_linkable(pf->dof_types()[i])) { discontinuous = true; break; }
       }
-      pfem classical_pf1 = discontinuous ? classical_discontinuous_fem(pgt, 1)
-                                        : classical_fem(pgt, 1);
-      short_type degree = short_type(((pf != classical_pf1 &&
-                                   pf->estimated_degree() > 1) ||
-                                  pgt->structure() !=
-                                  pgt->basic_structure()) ? 2 : 1);
-      pmf->set_finite_element(cv, discontinuous ?
-                              classical_discontinuous_fem(pgt, degree) :
-                              classical_fem(pgt, degree));
+
+      if (pf == fem_descriptor("FEM_Q2_INCOMPLETE(2)") ||
+          pf == fem_descriptor("FEM_Q2_INCOMPLETE(3)"))
+        pmf->set_finite_element(cv, pf);
+      else {
+        pfem classical_pf1 = discontinuous ? classical_discontinuous_fem(pgt, 
1)
+                                           : classical_fem(pgt, 1);
+        short_type degree = short_type(((pf != classical_pf1 &&
+                                         pf->estimated_degree() > 1) ||
+                                        pgt->structure() !=
+                                        pgt->basic_structure()) ? 2 : 1);
+        pmf->set_finite_element(cv, discontinuous ?
+                                classical_discontinuous_fem(pgt, degree) :
+                                classical_fem(pgt, degree));
+      }
     }
     /* find out which dof will be exported to VTK */
 
@@ -181,12 +189,14 @@
           if (nbd == 3) t = VTK_TRIANGLE;
           else if (nbd == 4) t = check_voxel(m.points_of_convex(cv)) ? 
VTK_PIXEL : VTK_QUAD;
           else if (nbd == 6) t = VTK_QUADRATIC_TRIANGLE;
-          else if (nbd == 9) t = VTK_QUADRATIC_QUAD; break;
+          else if (nbd == 8) t = VTK_QUADRATIC_QUAD;
+          else if (nbd == 9) t = VTK_BIQUADRATIC_QUAD; break;
         case 3:
           if (nbd == 4) t = VTK_TETRA;
           else if (nbd == 10) t = VTK_QUADRATIC_TETRA;
           else if (nbd == 8) t = check_voxel(m.points_of_convex(cv)) ? 
VTK_VOXEL : VTK_HEXAHEDRON;
-          else if (nbd == 27) t = VTK_QUADRATIC_HEXAHEDRON;
+          else if (nbd == 20) t = VTK_QUADRATIC_HEXAHEDRON;
+          else if (nbd == 27) t = VTK_TRIQUADRATIC_HEXAHEDRON;
           else if (nbd == 6) t = VTK_WEDGE; break;
       }
       GMM_ASSERT1(t != -1, "semi internal error.. could not map " <<