feat: split FiniteElementBase into nonvirtual and virtual interfaces

.github/workflows/ci_tests.yml

@@ -374,8 +374,8 @@ jobs:
             GEOS_ENABLE_BOUNDS_CHECK: OFF
             DOCKER_REPOSITORY: geosx/rockylinux8-clang17-cuda12.9.1
             RUNS_ON: streak
-            NPROC: 8
-            DOCKER_RUN_ARGS: "--cpus=8 --memory=256g --runtime=nvidia --gpus all -v /etc/pki/ca-trust/source/anchors/:/usr/local/share/ca-certificates/llnl:ro"
+            NPROC: 16
+            DOCKER_RUN_ARGS: "--cpus=16 --memory=448g --runtime=nvidia --gpus all -v /etc/pki/ca-trust/source/anchors/:/usr/local/share/ca-certificates/llnl:ro"
             DOCKER_CERTS_DIR: "/usr/local/share/ca-certificates"
             DOCKER_CERTS_UPDATE_COMMAND: "update-ca-trust"
             HOST_CONFIG: /spack-generated.cmake

scripts/setupPythonEnvironment.bash

@@ -71,7 +71,7 @@ case $key in
     echo "-p/--python-target \"Target parent python bin\""
     echo "-b/--bin-dir \"Directory to link new scripts\""
     echo "-d/--pkg-dir \"Directory containing target python packages\""
-    echo "-t/--tool-branch \"Target branch for geosPythonPackages (default=main) \""
+    echo "-r/--python-pkg-branch \"Target branch for geosPythonPackages (default=main) \""
     echo "-v/--verbose \"Increase verbosity level\""
     echo ""
     exit

src/CMakeLists.txt

@@ -139,6 +139,11 @@ blt_add_executable( NAME geosx
                                ${extraComponentsLinkList}
                                ${externalComponentsLinkList} )
 
+if(APPLE)
+  target_link_options(geosx PRIVATE "-Wl,-export_dynamic")
+endif()
+
+
 # Seems to be required on some CMake versions (e.g. 3.16) to get enforce device linking
 if( ${ENABLE_HYPRE_DEVICE} STREQUAL "CUDA"  )
   set_target_properties( geosx PROPERTIES CUDA_RESOLVE_DEVICE_SYMBOLS TRUE )
@@ -220,8 +225,8 @@ install( FILES ${CMAKE_BINARY_DIR}/schema.xsd
 ################################
 # Add python environment setup
 ################################
-#message(WARNING "Temporarily changing the geosPythonBranch to cusini/remove-relative-paths")
-#set(GEOS_PYTHON_PACKAGES_BRANCH "cusini/remove-relative-paths" CACHE STRING "" FORCE)
+message(WARNING "Temporarily changing the geosPythonBranch to feature/remove-dndx-detj")
+set(GEOS_PYTHON_PACKAGES_BRANCH "feature/remove-dndx-detj" CACHE STRING "" FORCE)
 
 
 if ( Python3_EXECUTABLE )

src/cmake/thirdparty/SetupGeosxThirdParty.cmake

@@ -887,7 +887,7 @@ if( ${CMAKE_VERSION} VERSION_LESS "3.19" )
     set( PYTHON_AND_VERSION Python3 )
     set( PYTHON_OPTIONAL_COMPONENTS)
 else()
-    set( PYTHON_AND_VERSION Python3 3.6.0...<4 )
+    set( PYTHON_AND_VERSION Python3 3.6.0...<4.0.0 )
     set( PYTHON_OPTIONAL_COMPONENTS OPTIONAL_COMPONENTS Development NumPy)
 endif()
 if(ENABLE_PYGEOSX)

src/coreComponents/finiteElement/CMakeLists.txt

@@ -29,6 +29,7 @@ set( finiteElement_headers
      LinearFormUtilities.hpp
      PDEUtilities.hpp
      elementFormulations/FiniteElementBase.hpp
+     elementFormulations/FiniteElementOperators.hpp
      elementFormulations/H1_Hexahedron_Lagrange1_GaussLegendre2.hpp
      elementFormulations/H1_Pyramid_Lagrange1_Gauss5.hpp
      elementFormulations/H1_QuadrilateralFace_Lagrange1_GaussLegendre2.hpp

src/coreComponents/finiteElement/FiniteElementDiscretization.hpp

@@ -27,23 +27,9 @@
 #include "LvArray/src/tensorOps.hpp"
 #include "FiniteElementDispatch.hpp"
 
-
 namespace geos
 {
 
-// TODO remove when these quantities are placed inside the FiniteElementBase
-// class.
-namespace dataRepository
-{
-namespace keys
-{
-string const dNdX = "dNdX";
-string const detJ = "detJ";
-}
-}
-
-
-
 class FiniteElementDiscretization : public dataRepository::Group
 {
 public:
@@ -71,14 +57,6 @@ class FiniteElementDiscretization : public dataRepository::Group
 
   ///@}
 
-  template< typename SUBREGION_TYPE,
-            typename FE_TYPE >
-  void calculateShapeFunctionGradients( arrayView2d< real64 const, nodes::REFERENCE_POSITION_USD > const & X,
-                                        SUBREGION_TYPE * const elementSubRegion,
-                                        typename FE_TYPE::template MeshData< SUBREGION_TYPE > meshData,
-                                        FE_TYPE & fe ) const;
-
-
   /**
    * @brief Factory method to instantiate a type of finite element formulation.
    * @param parentElementShape String key that indicates the type of
@@ -125,65 +103,6 @@ ENUM_STRINGS( FiniteElementDiscretization::Formulation,
               "SEM",
               "DG" );
 
-template< typename SUBREGION_TYPE,
-          typename FE_TYPE >
-void
-FiniteElementDiscretization::
-  calculateShapeFunctionGradients( arrayView2d< real64 const, nodes::REFERENCE_POSITION_USD > const & X,
-                                   SUBREGION_TYPE * const elementSubRegion,
-                                   typename FE_TYPE::template MeshData< SUBREGION_TYPE > meshData,
-                                   FE_TYPE & finiteElement ) const
-{
-  GEOS_MARK_FUNCTION;
-
-  // do not precompute shape functions in case of SEM or DG formulation (not needed)
-  if( m_formulation == Formulation::SEM || m_formulation == Formulation::DG )
-    return;
-
-  array4d< real64 > & dNdX = elementSubRegion->dNdX();
-  array2d< real64 > & detJ = elementSubRegion->detJ();
-  auto const & elemsToNodes = elementSubRegion->nodeList().toViewConst();
-
-  constexpr localIndex numNodesPerElem = FE_TYPE::maxSupportPoints;
-  constexpr localIndex numQuadraturePointsPerElem = FE_TYPE::numQuadraturePoints;
-  dNdX.resizeWithoutInitializationOrDestruction( elementSubRegion->size(), numQuadraturePointsPerElem, numNodesPerElem, 3 );
-  detJ.resize( elementSubRegion->size(), numQuadraturePointsPerElem );
-
-  finiteElement.setGradNView( dNdX.toViewConst() );
-  finiteElement.setDetJView( detJ.toViewConst() );
-
-  for( localIndex k = 0; k < elementSubRegion->size(); ++k )
-  {
-    typename FE_TYPE::StackVariables feStack;
-    finiteElement.template setup< FE_TYPE >( k, meshData, feStack );
-    real64 xLocal[numNodesPerElem][3];
-    localIndex numSupportPoints = finiteElement.template numSupportPoints< FE_TYPE >( feStack );
-    for( localIndex a=0; a< numSupportPoints; ++a )
-    {
-      localIndex const nodeIndex = elemsToNodes[ k][ a ];
-      for( int i=0; i<3; ++i )
-      {
-        xLocal[ a ][ i ] = X[ nodeIndex ][ i ];
-      }
-    }
-
-
-
-    for( localIndex q = 0; q < numQuadraturePointsPerElem; ++q )
-    {
-      real64 dNdXLocal[numNodesPerElem][3];
-      detJ( k, q ) = finiteElement.calcGradN( q, xLocal, feStack, dNdXLocal );
-
-      for( localIndex b = 0; b < numSupportPoints; ++b )
-      {
-        LvArray::tensorOps::copy< 3 >( dNdX[ k ][ q ][ b ], dNdXLocal[b] );
-      }
-    }
-  }
-
-}
-
-
 } /* namespace geos */
 
 #endif /* GEOS_FINITEELEMENT_FINITEELEMENTDISCRETIZATION_HPP_ */

src/coreComponents/finiteElement/FiniteElementDispatch.hpp

@@ -159,44 +159,41 @@ struct FiniteElementDispatchHandler< FE_TYPE, FE_TYPES... >
   dispatch3D( FiniteElementBase const & input,
               LAMBDA && lambda )
   {
-    if( auto const * const ptr = dynamic_cast< FE_TYPE const * >(&input) )
+    if( auto const * const wrapper = dynamic_cast< FE_TYPE const * >(&input) )
     {
-      lambda( *ptr );
-    }
-    else
-    {
-      FiniteElementDispatchHandler< FE_TYPES... >::dispatch3D( input, lambda );
+      lambda( *wrapper->getImpl() );
+      return;
     }
+
+    FiniteElementDispatchHandler< FE_TYPES... >::dispatch3D( input, lambda );
   }
 
   template< typename LAMBDA >
   static void
   dispatch3D( FiniteElementBase & input,
               LAMBDA && lambda )
   {
-    if( auto * const ptr = dynamic_cast< FE_TYPE * >(&input) )
+    if( auto * const wrapper = dynamic_cast< FE_TYPE * >(&input) )
     {
-      lambda( *ptr );
-    }
-    else
-    {
-      FiniteElementDispatchHandler< FE_TYPES... >::dispatch3D( input, lambda );
+      lambda( *wrapper->getImpl() );
+      return;
     }
+
+    FiniteElementDispatchHandler< FE_TYPES... >::dispatch3D( input, lambda );
   }
 
   template< typename LAMBDA >
   static void
   dispatch2D( FiniteElementBase const & input,
               LAMBDA && lambda )
   {
-    if( auto const * const ptr = dynamic_cast< FE_TYPE const * >(&input) )
+    if( auto const * const wrapper = dynamic_cast< FE_TYPE const * >(&input) )
     {
-      lambda( *ptr );
-    }
-    else
-    {
-      FiniteElementDispatchHandler< FE_TYPES... >::dispatch2D( input, lambda );
+      lambda( *wrapper->getImpl() );
+      return;
     }
+
+    FiniteElementDispatchHandler< FE_TYPES... >::dispatch2D( input, lambda );
   }
 };
 
@@ -207,21 +204,21 @@ void
 dispatchlowOrder3D( FiniteElementBase const & input,
                     LAMBDA && lambda )
 {
-  if( auto const * const ptr1 = dynamic_cast< H1_Hexahedron_Lagrange1_GaussLegendre2 const * >(&input) )
+  if( auto const * const wrapper1 = dynamic_cast< H1_Hexahedron_Lagrange1_GaussLegendre2 const * >(&input) )
   {
-    lambda( *ptr1 );
+    lambda( *wrapper1->getImpl() );
   }
-  else if( auto const * const ptr2 = dynamic_cast< H1_Wedge_Lagrange1_Gauss6 const * >(&input) )
+  else if( auto const * const wrapper2 = dynamic_cast< H1_Wedge_Lagrange1_Gauss6 const * >(&input) )
   {
-    lambda( *ptr2 );
+    lambda( *wrapper2->getImpl() );
   }
-  else if( auto const * const ptr3 = dynamic_cast< H1_Tetrahedron_Lagrange1_Gauss1 const * >(&input) )
+  else if( auto const * const wrapper3 = dynamic_cast< H1_Tetrahedron_Lagrange1_Gauss1 const * >(&input) )
   {
-    lambda( *ptr3 );
+    lambda( *wrapper3->getImpl() );
   }
-  else if( auto const * const ptr4 = dynamic_cast< H1_Pyramid_Lagrange1_Gauss5 const * >(&input) )
+  else if( auto const * const wrapper4 = dynamic_cast< H1_Pyramid_Lagrange1_Gauss5 const * >(&input) )
   {
-    lambda( *ptr4 );
+    lambda( *wrapper4->getImpl() );
   }
   else
   {
@@ -233,6 +230,4 @@ dispatchlowOrder3D( FiniteElementBase const & input,
 
 } // namespace geos
 
-
-
 #endif /* GEOS_FINITEELEMENT_FINITEELEMENTDISPATCH_HPP_ */