Increased Intrepid2::MaxOrder from 8 to 10.

Modified StructuredIntegrationPerformance to allow saving calibrations to file, and reading calibrations from file.
trilinos · May 20, 2024 · 1f22c13 · 1f22c13
1 parent 573fbbd
commit 1f22c13
Show file tree

Hide file tree

Showing 3 changed files with 32 additions and 8 deletions.
diff --git a/packages/intrepid2/src/Shared/Intrepid2_Polylib.hpp b/packages/intrepid2/src/Shared/Intrepid2_Polylib.hpp
@@ -208,7 +208,8 @@ namespace Intrepid2 {
     static constexpr ordinal_type  MaxPolylibOrder =
         (Parameters::MaxOrder > Parameters::MaxCubatureDegreeEdge) ? Parameters::MaxOrder :
                                                                      Parameters::MaxCubatureDegreeEdge;
-    static constexpr ordinal_type  MaxPolylibPoint = MaxPolylibOrder/2+2;
+    // NVR: HVOL bases on tri/tet use Polylib with order + spaceDim + 2 points; in 3D this can be up to Parameters::MaxOrder + 5.
+    static constexpr ordinal_type  MaxPolylibPoint = (MaxPolylibOrder/2+2 > Parameters::MaxOrder + 5) ? MaxPolylibOrder/2+2 : Parameters::MaxOrder + 5;
 
     struct Serial {
 

diff --git a/packages/intrepid2/src/Shared/Intrepid2_Types.hpp b/packages/intrepid2/src/Shared/Intrepid2_Types.hpp
@@ -120,7 +120,7 @@ namespace Intrepid2 {
     /// The maximum number of points to eval in serial mode.
     static constexpr ordinal_type MaxNumPtsPerBasisEval= 1;
     /// The maximum reconstruction order.
-    static constexpr ordinal_type MaxOrder             = 8;
+    static constexpr ordinal_type MaxOrder             = 10;
     /// The maximum number of integration points for direct cubature rules.
     static constexpr ordinal_type MaxIntegrationPoints = 4893;    
     /// The maximum degree of the polynomial that can be integrated exactly by a direct edge rule.

diff --git a/...ntrepid2/unit-test/performance/StructuredIntegration/StructuredIntegrationPerformance.cpp b/...ntrepid2/unit-test/performance/StructuredIntegration/StructuredIntegrationPerformance.cpp
@@ -1257,6 +1257,14 @@ int main( int argc, char* argv[] )
       return -1;
     }
 
+    if (mode == Test)
+    {
+      // then overwrite poly orders to match the test cases we define in getWorksetSizeMap
+      polyOrderFixed = -1;
+      polyOrderMin = 1;
+      polyOrderMax = 3;
+    }
+
     using Scalar = double;
     using ExecutionSpace = Kokkos::DefaultExecutionSpace;
     using DeviceType = Kokkos::DefaultExecutionSpace::device_type;
@@ -1284,16 +1292,31 @@ int main( int argc, char* argv[] )
     map<int, Kokkos::Array<int,spaceDim> > gridCellCountsForPolyOrder;
     for (int p=polyOrderMin; p<=polyOrderMax; p++)
     {
-      int maxBasisCardinality = 1;
-      for (auto formulationChoice : formulationChoices)
+      if (mode != Test)
+      {
+        int maxBasisCardinality = 1;
+        for (auto formulationChoice : formulationChoices)
+        {
+          for (auto basisFamilyChoice : basisFamilyChoices)
+          {
+            auto basis = getHypercubeBasisForFormulation<Scalar, Scalar, DeviceType, spaceDim>(formulationChoice, basisFamilyChoice, p);
+            maxBasisCardinality = std::max(basis->getCardinality(), maxBasisCardinality);
+          }
+        }
+        gridCellCountsForPolyOrder[p] = getMeshWidths<spaceDim>(maxBasisCardinality, maxEntryCount, maxCellCount);
+      }
+      else // mode == Test; set the gridCellCounts according to mesh widths indicated in getWorksetSizes above.
       {
-        for (auto basisFamilyChoice : basisFamilyChoices)
+        Kokkos::Array<int,spaceDim> meshWidths_8, meshWidths_4;
+        for (int d=0; d<spaceDim; d++)
         {
-          auto basis = getHypercubeBasisForFormulation<Scalar, Scalar, DeviceType, spaceDim>(formulationChoice, basisFamilyChoice, p);
-          maxBasisCardinality = std::max(basis->getCardinality(), maxBasisCardinality);
+          meshWidths_4[d] = 4;
+          meshWidths_8[d] = 8;
         }
+        gridCellCountsForPolyOrder[1] = meshWidths_8;
+        gridCellCountsForPolyOrder[2] = meshWidths_8;
+        gridCellCountsForPolyOrder[3] = meshWidths_4;
       }
-      gridCellCountsForPolyOrder[p] = getMeshWidths<spaceDim>(maxBasisCardinality, maxEntryCount, maxCellCount);
 
       int cellCount = 1;
       for (int d=0; d<spaceDim; d++)