Adding PYRAMID13 support in miscellaneous places:

* string_to_enum() * Various functions in fe_lagrange.C * Quadrature initialization routines * FEAbstract::on_reference_element() * lagrange_nodal_soln() helper function * FEAbstract::get_refspace_nodes() * FEInterface::max_order() * l2_lagrange_n_dofs() * monomial_n_dofs() * xyz_n_dofs() * Quality::valid() * ElementTypes::basic_name() * ElementTypes::name() * QGrid::init_3D()
colibius · Feb 24, 2014 · 29aae59 · 29aae59
1 parent b6a276e
commit 29aae59
Show file tree

Hide file tree

Showing 12 changed files with 93 additions and 0 deletions.
diff --git a/src/fe/fe_abstract.C b/src/fe/fe_abstract.C
@@ -644,6 +644,33 @@ void FEAbstract::get_refspace_nodes(const ElemType itemType, std::vector<Point>&
        nodes[4] = Point (0.,0.,1.);
        return;
     }
+    case PYRAMID13:
+    {
+       nodes.resize(13);
+
+       // base corners
+       nodes[0] = Point (-1.,-1.,0.);
+       nodes[1] = Point (1.,-1.,0.);
+       nodes[2] = Point (1.,1.,0.);
+       nodes[3] = Point (-1.,1.,0.);
+
+       // apex
+       nodes[4] = Point (0.,0.,1.);
+
+       // base midedge
+       nodes[5] = Point (0.,-1.,0.);
+       nodes[6] = Point (1.,0.,0.);
+       nodes[7] = Point (0.,1.,0.);
+       nodes[8] = Point (-1,0.,0.);
+
+       // lateral midedge
+       nodes[9] = Point (-.5,-.5,.5);
+       nodes[10] = Point (.5,-.5,.5);
+       nodes[11] = Point (.5,.5,.5);
+       nodes[12] = Point (-.5,.5,.5);
+
+       return;
+    }
     case PYRAMID14:
     {
        nodes.resize(14);
@@ -811,6 +838,7 @@ bool FEAbstract::on_reference_element(const Point& p, const ElemType t, const Re
 
 
     case PYRAMID5:
+    case PYRAMID13:
     case PYRAMID14:
       {
         // Check that the point is on the same side of all the faces

diff --git a/src/fe/fe_interface.C b/src/fe/fe_interface.C
@@ -1038,6 +1038,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	      return 2;
 	    case PYRAMID5:
 	      return 1;
+            case PYRAMID13:
             case PYRAMID14:
               return 2;
 	    default:
@@ -1064,6 +1065,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	    case PRISM15:
 	    case PRISM18:
 	    case PYRAMID5:
+            case PYRAMID13:
             case PYRAMID14:
 	      return unlimited;
 	    default:
@@ -1101,6 +1103,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	    case PRISM15:
 	    case PRISM18:
 	    case PYRAMID5:
+            case PYRAMID13:
             case PYRAMID14:
 	      return 0;
 	    default:
@@ -1132,6 +1135,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	    case PRISM15:
 	    case PRISM18:
 	    case PYRAMID5:
+            case PYRAMID13:
             case PYRAMID14:
 	      return 0;
 	    default:
@@ -1159,6 +1163,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	    case PRISM15:
 	    case PRISM18:
 	    case PYRAMID5:
+            case PYRAMID13:
             case PYRAMID14:
 	      return unlimited;
 	    default:
@@ -1188,6 +1193,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	    case PRISM15:
 	    case PRISM18:
 	    case PYRAMID5:
+            case PYRAMID13:
             case PYRAMID14:
 	      return 0;
 	    default:
@@ -1221,6 +1227,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	    case PRISM15:
 	    case PRISM18:
 	    case PYRAMID5:
+            case PYRAMID13:
             case PYRAMID14:
 	      return 0;
 	    default:
@@ -1255,6 +1262,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	    case PRISM15:
 	    case PRISM18:
 	    case PYRAMID5:
+            case PYRAMID13:
 	    case PYRAMID14:
 	      return 0;
 	    default:
@@ -1289,6 +1297,7 @@ unsigned int FEInterface::max_order(const FEType& fe_t,
 	    case PRISM15:
 	    case PRISM18:
 	    case PYRAMID5:
+            case PYRAMID13:
 	    case PYRAMID14:
 	      return 0;
 	    default:

diff --git a/src/fe/fe_l2_lagrange.C b/src/fe/fe_l2_lagrange.C
@@ -329,6 +329,7 @@ namespace libMesh
 		return 6;
 
 	      case PYRAMID5:
+              case PYRAMID13:
               case PYRAMID14:
 		return 5;
 

diff --git a/src/fe/fe_lagrange.C b/src/fe/fe_lagrange.C
@@ -227,6 +227,28 @@ namespace libMesh
 		  return;
 		}
 
+              case PYRAMID13:
+                {
+                  libmesh_assert_equal_to (elem_soln.size(), 5);
+                  libmesh_assert_equal_to (nodal_soln.size(), 13);
+
+                  nodal_soln[0]  = elem_soln[0];
+                  nodal_soln[1]  = elem_soln[1];
+                  nodal_soln[2]  = elem_soln[2];
+                  nodal_soln[3]  = elem_soln[3];
+                  nodal_soln[4]  = elem_soln[4];
+                  nodal_soln[5]  = .5*(elem_soln[0] + elem_soln[1]);
+                  nodal_soln[6]  = .5*(elem_soln[1] + elem_soln[2]);
+                  nodal_soln[7]  = .5*(elem_soln[2] + elem_soln[3]);
+                  nodal_soln[8]  = .5*(elem_soln[3] + elem_soln[0]);
+                  nodal_soln[9]  = .5*(elem_soln[0] + elem_soln[4]);
+                  nodal_soln[10] = .5*(elem_soln[1] + elem_soln[4]);
+                  nodal_soln[11] = .5*(elem_soln[2] + elem_soln[4]);
+                  nodal_soln[12] = .5*(elem_soln[3] + elem_soln[4]);
+
+                  return;
+                }
+
               case PYRAMID14:
                 {
                   libmesh_assert_equal_to (elem_soln.size(), 5);
@@ -349,6 +371,7 @@ namespace libMesh
 		return 6;
 
 	      case PYRAMID5:
+              case PYRAMID13:
               case PYRAMID14:
 		return 5;
 
@@ -400,6 +423,9 @@ namespace libMesh
 	      case PRISM18:
 		return 18;
 
+              case PYRAMID13:
+                return 13;
+
               case PYRAMID14:
                 return 14;
 
@@ -568,6 +594,7 @@ namespace libMesh
 		}
 
 	      case PYRAMID5:
+              case PYRAMID13:
               case PYRAMID14:
 		{
 		  switch (n)
@@ -612,6 +639,7 @@ namespace libMesh
 	      case HEX27:
 	      case PRISM15:
 	      case PRISM18:
+              case PYRAMID13:
               case PYRAMID14:
 		return 1;
 

diff --git a/src/fe/fe_monomial.C b/src/fe/fe_monomial.C
@@ -140,6 +140,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
               case PYRAMID14:
 		return 4;
 
@@ -186,6 +187,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 10;
 
@@ -232,6 +234,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 20;
 
@@ -277,6 +280,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 35;
 
@@ -321,6 +325,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return (order+1)*(order+2)*(order+3)/6;
 

diff --git a/src/fe/fe_xyz.C b/src/fe/fe_xyz.C
@@ -137,6 +137,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 4;
 
@@ -183,6 +184,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 10;
 
@@ -229,6 +231,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 20;
 
@@ -274,6 +277,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 35;
 
@@ -318,6 +322,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return (order+1)*(order+2)*(order+3)/6;
 
@@ -385,6 +390,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 4;
 
@@ -434,6 +440,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 10;
 
@@ -480,6 +487,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 20;
 
@@ -526,6 +534,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return 35;
 
@@ -569,6 +578,7 @@ namespace libMesh
 	      case PRISM15:
 	      case PRISM18:
 	      case PYRAMID5:
+              case PYRAMID13:
 	      case PYRAMID14:
 		return (order+1)*(order+2)*(order+3)/6;
 

diff --git a/src/geom/elem_quality.C b/src/geom/elem_quality.C
@@ -408,6 +408,7 @@ std::vector<ElemQuality> Quality::valid(const ElemType t)
       }
 
     case PYRAMID5:
+    case PYRAMID13:
     case PYRAMID14:
       {
 	// None yet

diff --git a/src/geom/elem_type.C b/src/geom/elem_type.C
@@ -79,6 +79,7 @@ std::string ElementTypes::basic_name (const ElemType t)
       }
 
     case PYRAMID5:
+    case PYRAMID13:
     case PYRAMID14:
       {
 	its_name = "Pyramid";
@@ -232,6 +233,12 @@ std::string ElementTypes::name(const ElemType t)
 	break;
       }
 
+    case PYRAMID13:
+      {
+        its_name = "Pyramid 13";
+        break;
+      }
+
     case PYRAMID14:
       {
 	its_name = "Pyramid 14";

diff --git a/src/quadrature/quadrature_conical_3D.C b/src/quadrature/quadrature_conical_3D.C
@@ -38,6 +38,7 @@ void QConical::init_3D(const ElemType type_in,
       } // end case TET4, TET10
 
     case PYRAMID5:
+    case PYRAMID13:
     case PYRAMID14:
       {
 	this->conical_product_pyramid(p);

diff --git a/src/quadrature/quadrature_gauss_3D.C b/src/quadrature/quadrature_gauss_3D.C
@@ -525,6 +525,7 @@ void QGauss::init_3D(const ElemType type_in,
       //---------------------------------------------
       // Pyramid
     case PYRAMID5:
+    case PYRAMID13:
     case PYRAMID14:
       {
 	// We compute the Pyramid rule as a conical product of a

diff --git a/src/quadrature/quadrature_grid_3D.C b/src/quadrature/quadrature_grid_3D.C
@@ -107,6 +107,7 @@ void QGrid::init_3D(const ElemType type_in,
       //---------------------------------------------
       // Pyramid
     case PYRAMID5:
+    case PYRAMID13:
     case PYRAMID14:
       {
         _points.resize((_order+1)*(_order+2)*(_order+3)/6);

diff --git a/src/utils/string_to_enum.C b/src/utils/string_to_enum.C
@@ -134,6 +134,7 @@ INSTANTIATE_ENUM_MAPS(ElemType, elem_type)
 
 	elem_type_to_enum["PYRAMID"   ]=PYRAMID5;
 	elem_type_to_enum["PYRAMID5"  ]=PYRAMID5;
+        elem_type_to_enum["PYRAMID13" ]=PYRAMID13;
 	elem_type_to_enum["PYRAMID14" ]=PYRAMID14;
 
 	elem_type_to_enum["INFEDGE"   ]=INFEDGE2;