]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
added new hierachic_to_lexicographic_numbering function
authorangela.klewinghaus <angela.klewinghaus@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 25 Feb 2013 13:21:31 +0000 (13:21 +0000)
committerangela.klewinghaus <angela.klewinghaus@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 25 Feb 2013 13:21:31 +0000 (13:21 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_face_elements@28549 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/fe/fe_poly_face.templates.h
deal.II/include/deal.II/fe/fe_tools.h
deal.II/source/fe/fe_tools.cc

index 513213aeed0be8ff035a3403825a5aa21e22928a..61c6bfdd2d5bfed4377332971559299be72206d3 100644 (file)
@@ -205,14 +205,23 @@ FE_PolyFace<POLY,dim,spacedim>::fill_fe_face_values (
 
   const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
 
-  const unsigned int foffset = fe_data.shape_values.size() * face;
   if (flags & update_values)
     for (unsigned int i=0; i<quadrature.size(); ++i)
       {
         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
           data.shape_values(k,i) = 0.;
+       // Fill data for vertex shape functions
+       
+       // Fill data for line shape functions
+
+       // in 3D, this should be a loop over all lines of the face,
+       // and the multiplier should not be the face itself
+       const unsigned int foffset = this->first_line_index + this->dofs_per_line * face;
+
         for (unsigned int k=0; k<fe_data.shape_values.size(); ++k)
-          data.shape_values(foffset+k,i) = fe_data.shape_values[k][i];
+          data.shape_values(foffset+k,i) = fe_data.shape_values[k+this->first_face_line_index][i];
+
+       // Fill data for quad shape functions
       }
 }
 
index d0cb7482c679dca5d9fa371b67e634c37c737f0f..b4b82667a18b1ab8d4d194284f2b8aca5d2d4ce0 100644 (file)
@@ -1148,6 +1148,12 @@ namespace FETools
    * of degrees of freedom in the
    * finite element.
    */
+
+  template <int dim>
+  void
+  hierarchic_to_lexicographic_numbering (unsigned int degree, 
+                                        std::vector<unsigned int> &h2l);
+
   template <int dim>
   void
   hierarchic_to_lexicographic_numbering (const FiniteElementData<dim> &fe_data,
@@ -1171,7 +1177,8 @@ namespace FETools
    * numbering. All the remarks
    * made about the above function
    * are also valid here.
-   */
+   */  
+
   template <int dim>
   void
   lexicographic_to_hierarchic_numbering (const FiniteElementData<dim> &fe_data,
index 26ed3a42fa02dd22feb85a5f149ee4766192e3c9..dba1b97ff305accc24f3f01e6efdb5b3f83b2ddd 100644 (file)
@@ -2602,30 +2602,21 @@ namespace FETools
 
   template <int dim>
   void
-  hierarchic_to_lexicographic_numbering (const FiniteElementData<dim> &fe,
-                                         std::vector<unsigned int> &h2l)
-  {
-    Assert (h2l.size() == fe.dofs_per_cell,
-            ExcDimensionMismatch (h2l.size(), fe.dofs_per_cell));
-    h2l = hierarchic_to_lexicographic_numbering (fe);
-  }
-
-
-
-  template <int dim>
-  std::vector<unsigned int>
-  hierarchic_to_lexicographic_numbering (const FiniteElementData<dim> &fe)
+  hierarchic_to_lexicographic_numbering (unsigned int degree, std::vector<unsigned int>& h2l)
   {
-    Assert (fe.n_components() == 1, ExcInvalidFE());
-
-    std::vector<unsigned int> h2l (fe.dofs_per_cell);
-
-    const unsigned int dofs_per_cell = fe.dofs_per_cell;
-    // polynomial degree
-    const unsigned int degree = fe.dofs_per_line+1;
-    // number of grid points in each
+    // number of support points in each
     // direction
     const unsigned int n = degree+1;
+    
+    unsigned int dofs_per_cell = n;
+    for (unsigned int i=1;i<dim;++i)
+      dofs_per_cell *= n;
+    
+    // Assert size maches degree
+    AssertDimension (h2l.size(), dofs_per_cell);
+
+    // polynomial degree
+    const unsigned int dofs_per_line = degree - 1;
 
     // the following lines of code are
     // somewhat odd, due to the way the
@@ -2661,29 +2652,29 @@ namespace FETools
         h2l[next_index++] = n*n-1;
 
         // left   line
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = (1+i)*n;
 
         // right  line
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = (2+i)*n-1;
 
         // bottom line
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = 1+i;
 
         // top    line
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = n*(n-1)+i+1;
 
         // inside quad
-        Assert (fe.dofs_per_quad == fe.dofs_per_line*fe.dofs_per_line,
-                ExcInternalError());
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
-          for (unsigned int j=0; j<fe.dofs_per_line; ++j)
+       // Assert (fe.dofs_per_quad == fe.dofs_per_line*fe.dofs_per_line,
+       //         ExcInternalError());
+        for (unsigned int i=0; i<dofs_per_line; ++i)
+          for (unsigned int j=0; j<dofs_per_line; ++j)
             h2l[next_index++] = n*(i+1)+j+1;
 
-        Assert (next_index == fe.dofs_per_cell, ExcInternalError());
+        Assert (next_index == dofs_per_cell, ExcInternalError());
 
         break;
       }
@@ -2702,82 +2693,82 @@ namespace FETools
         h2l[next_index++] = (n*n+n+1)*degree;  // 7
 
         // line 0
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = (i+1)*n;
         // line 1
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = n-1+(i+1)*n;
         // line 2
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = 1+i;
         // line 3
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = 1+i+n*(n-1);
 
         // line 4
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = (n-1)*n*n+(i+1)*n;
         // line 5
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = (n-1)*(n*n+1)+(i+1)*n;
         // line 6
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = n*n*(n-1)+i+1;
         // line 7
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = n*n*(n-1)+i+1+n*(n-1);
 
         // line 8
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = (i+1)*n*n;
         // line 9
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = n-1+(i+1)*n*n;
         // line 10
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = (i+1)*n*n+n*(n-1);
         // line 11
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
           h2l[next_index++] = n-1+(i+1)*n*n+n*(n-1);
 
 
         // inside quads
-        Assert (fe.dofs_per_quad == fe.dofs_per_line*fe.dofs_per_line,
-                ExcInternalError());
+       // Assert (fe.dofs_per_quad == fe.dofs_per_line*fe.dofs_per_line,
+        //        ExcInternalError());
         // face 0
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
-          for (unsigned int j=0; j<fe.dofs_per_line; ++j)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
+          for (unsigned int j=0; j<dofs_per_line; ++j)
             h2l[next_index++] = (i+1)*n*n+n*(j+1);
         // face 1
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
-          for (unsigned int j=0; j<fe.dofs_per_line; ++j)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
+          for (unsigned int j=0; j<dofs_per_line; ++j)
             h2l[next_index++] = (i+1)*n*n+n-1+n*(j+1);
         // face 2, note the orientation!
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
-          for (unsigned int j=0; j<fe.dofs_per_line; ++j)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
+          for (unsigned int j=0; j<dofs_per_line; ++j)
             h2l[next_index++] = (j+1)*n*n+i+1;
         // face 3, note the orientation!
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
-          for (unsigned int j=0; j<fe.dofs_per_line; ++j)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
+          for (unsigned int j=0; j<dofs_per_line; ++j)
             h2l[next_index++] = (j+1)*n*n+n*(n-1)+i+1;
         // face 4
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
-          for (unsigned int j=0; j<fe.dofs_per_line; ++j)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
+          for (unsigned int j=0; j<dofs_per_line; ++j)
             h2l[next_index++] = n*(i+1)+j+1;
         // face 5
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
-          for (unsigned int j=0; j<fe.dofs_per_line; ++j)
+        for (unsigned int i=0; i<dofs_per_line; ++i)
+          for (unsigned int j=0; j<dofs_per_line; ++j)
             h2l[next_index++] = (n-1)*n*n+n*(i+1)+j+1;
 
         // inside hex
-        Assert (fe.dofs_per_hex == fe.dofs_per_quad*fe.dofs_per_line,
-                ExcInternalError());
-        for (unsigned int i=0; i<fe.dofs_per_line; ++i)
-          for (unsigned int j=0; j<fe.dofs_per_line; ++j)
-            for (unsigned int k=0; k<fe.dofs_per_line; ++k)
+       // Assert (fe.dofs_per_hex == fe.dofs_per_quad*fe.dofs_per_line,
+        //        ExcInternalError());
+        for (unsigned int i=0; i<dofs_per_line; ++i)
+          for (unsigned int j=0; j<dofs_per_line; ++j)
+            for (unsigned int k=0; k<dofs_per_line; ++k)
               h2l[next_index++]       = n*n*(i+1)+n*(j+1)+k+1;
 
-        Assert (next_index == fe.dofs_per_cell, ExcInternalError());
+        Assert (next_index == dofs_per_cell, ExcInternalError());
 
         break;
       }
@@ -2786,11 +2777,33 @@ namespace FETools
         Assert (false, ExcNotImplemented());
       }
 
-    return h2l;
+    // return h2l;
   }
 
 
 
+  template <int dim>
+  void
+  hierarchic_to_lexicographic_numbering (const FiniteElementData<dim> &fe,
+                                         std::vector<unsigned int> &h2l)
+  {
+    Assert (h2l.size() == fe.dofs_per_cell,
+            ExcDimensionMismatch (h2l.size(), fe.dofs_per_cell));
+    hierarchic_to_lexicographic_numbering (h21, fe.dofs_per_line+1);
+  }
+
+
+
+  template <int dim>
+  std::vector<unsigned int>
+  hierarchic_to_lexicographic_numbering (const FiniteElementData<dim> &fe)
+  {
+    Assert (fe.n_components() == 1, ExcInvalidFE());
+    std::vector<unsigned int> h2l(fe.dofs_per_cell);
+    hierarchic_to_lexicographic_numbering (h2l, fe.dofs_per_line+1);
+    return (h2l);
+  }
+
   template <int dim>
   void
   lexicographic_to_hierarchic_numbering (const FiniteElementData<dim> &fe,

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.