From e8209d94ea8bca5ea2bdf7f6dd1b063c7fe67ebc Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Tue, 20 Jan 2026 12:02:45 -0600 Subject: [PATCH 1/9] Draft --- ...lidMechanicsAugmentedLagrangianContact.cpp | 24 +++++++++++++++---- ...lidMechanicsAugmentedLagrangianContact.hpp | 2 ++ 2 files changed, 22 insertions(+), 4 deletions(-) diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp index df469f26ea1..d4c1fb9b4c1 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp @@ -35,6 +35,8 @@ #include "constitutive/contact/FrictionSelector.hpp" #include "constitutive/solid/PorousSolid.hpp" #include "constitutive/solid/SolidFields.hpp" +#include "finiteElement/FiniteElementDiscretization.hpp" +//#include "mesh/ElementType.hpp" #include "mesh/DomainPartition.hpp" namespace geos @@ -48,10 +50,6 @@ SolidMechanicsAugmentedLagrangianContact::SolidMechanicsAugmentedLagrangianConta Group * const parent ): ContactSolverBase( name, parent ) { - - m_faceTypeToFiniteElements.insert( {"Quadrilateral", std::make_unique< finiteElement::H1_QuadrilateralFace_Lagrange1_GaussLegendre2 >()} ); - m_faceTypeToFiniteElements.insert( {"Triangle", std::make_unique< finiteElement::H1_TriangleFace_Lagrange1_Gauss1 >()} ); - registerWrapper( viewKeyStruct::simultaneousString(), &m_simultaneous ). setInputFlag( InputFlags::OPTIONAL ). setApplyDefaultValue( 1 ). @@ -219,6 +217,24 @@ void SolidMechanicsAugmentedLagrangianContact::setupSystem( DomainPartition & do PhysicsSolverBase::setupSystem( domain, dofManager, localMatrix, rhs, solution, setSparsity ); } +void SolidMechanicsAugmentedLagrangianContact::postInputInitialization() +{ + ContactSolverBase::postInputInitialization(); + + DomainPartition & domain = this->getGroupByPath< DomainPartition >( "/Problem/domain" ); + NumericalMethodsManager const & numericalMethodManager = domain.getNumericalMethodManager(); + FiniteElementDiscretizationManager const & feDiscretizationManager = + numericalMethodManager.getFiniteElementDiscretizationManager(); + FiniteElementDiscretization const & feDiscretization = + feDiscretizationManager.getGroup< FiniteElementDiscretization >( getDiscretizationName() ); + + m_faceTypeToFiniteElements["Quadrilateral"] = + feDiscretization.factory( ElementType::Quadrilateral ); + + m_faceTypeToFiniteElements["Triangle"] = + feDiscretization.factory( ElementType::Triangle ); +} + void SolidMechanicsAugmentedLagrangianContact::setSparsityPattern( DomainPartition & domain, DofManager & dofManager, CRSMatrix< real64, globalIndex > & GEOS_UNUSED_PARAM( localMatrix ), diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp index af1dc42caf8..2ebed39073e 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp @@ -47,6 +47,8 @@ class SolidMechanicsAugmentedLagrangianContact : public ContactSolverBase */ string getCatalogName() const override { return catalogName(); } + virtual void postInputInitialization() override; + virtual void registerDataOnMesh( dataRepository::Group & meshBodies ) override final; virtual void setupDofs( DomainPartition const & domain, From 403b4b1a92c86e8ae2445667434a80445dd2a4fd Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Tue, 20 Jan 2026 12:11:21 -0600 Subject: [PATCH 2/9] Removed includ --- .../contact/SolidMechanicsAugmentedLagrangianContact.cpp | 1 - 1 file changed, 1 deletion(-) diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp index d4c1fb9b4c1..6cbb712624f 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp @@ -36,7 +36,6 @@ #include "constitutive/solid/PorousSolid.hpp" #include "constitutive/solid/SolidFields.hpp" #include "finiteElement/FiniteElementDiscretization.hpp" -//#include "mesh/ElementType.hpp" #include "mesh/DomainPartition.hpp" namespace geos From 8a3adca4f9eb5030702f7e83bffb2800c07a1512 Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Tue, 20 Jan 2026 12:20:11 -0600 Subject: [PATCH 3/9] Added SolidMechanicsLagrangeContactBubbleStab --- ...olidMechanicsLagrangeContactBubbleStab.cpp | 23 +++++++++++++++---- ...olidMechanicsLagrangeContactBubbleStab.hpp | 2 ++ 2 files changed, 21 insertions(+), 4 deletions(-) diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp index ff06db81bfb..2271e0775f9 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp @@ -27,7 +27,7 @@ #include "physicsSolvers/solidMechanics/contact/kernels/SolidMechanicsContactFaceBubbleKernels.hpp" #include "physicsSolvers/solidMechanics/contact/LogLevelsInfo.hpp" #include "physicsSolvers/LogLevelsInfo.hpp" - +#include "finiteElement/FiniteElementDiscretization.hpp" #include "constitutive/contact/FrictionSelector.hpp" #include "fieldSpecification/FieldSpecificationManager.hpp" @@ -43,9 +43,6 @@ SolidMechanicsLagrangeContactBubbleStab::SolidMechanicsLagrangeContactBubbleStab Group * const parent ): ContactSolverBase( name, parent ) { - m_faceTypeToFiniteElements.insert( {"Quadrilateral", std::make_unique< finiteElement::H1_QuadrilateralFace_Lagrange1_GaussLegendre2 >()} ); - m_faceTypeToFiniteElements.insert( {"Triangle", std::make_unique< finiteElement::H1_TriangleFace_Lagrange1_Gauss1 >()} ); - LinearSolverParameters & linSolParams = m_linearSolverParameters.get(); linSolParams.mgr.strategy = LinearSolverParameters::MGR::StrategyType::lagrangianContactMechanicsBubbleStab; linSolParams.mgr.separateComponents = true; @@ -192,6 +189,24 @@ void SolidMechanicsLagrangeContactBubbleStab::setupDofs( DomainPartition const & meshTargets ); } +void SolidMechanicsLagrangeContactBubbleStab::postInputInitialization() +{ + ContactSolverBase::postInputInitialization(); + + DomainPartition & domain = this->getGroupByPath< DomainPartition >( "/Problem/domain" ); + NumericalMethodsManager const & numericalMethodManager = domain.getNumericalMethodManager(); + FiniteElementDiscretizationManager const & feDiscretizationManager = + numericalMethodManager.getFiniteElementDiscretizationManager(); + FiniteElementDiscretization const & feDiscretization = + feDiscretizationManager.getGroup< FiniteElementDiscretization >( getDiscretizationName() ); + + m_faceTypeToFiniteElements["Quadrilateral"] = + feDiscretization.factory( ElementType::Quadrilateral ); + + m_faceTypeToFiniteElements["Triangle"] = + feDiscretization.factory( ElementType::Triangle ); +} + void SolidMechanicsLagrangeContactBubbleStab::setupSystem( DomainPartition & domain, DofManager & dofManager, CRSMatrix< real64, globalIndex > & localMatrix, diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp index 6c8562044e4..aaff54e8fff 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp @@ -50,6 +50,8 @@ class SolidMechanicsLagrangeContactBubbleStab : public ContactSolverBase */ string getCatalogName() const override { return catalogName(); } + virtual void postInputInitialization() override; + virtual void registerDataOnMesh( Group & MeshBodies ) override final; real64 solverStep( real64 const & time_n, From bfa6b023ab4eea009590e5ee1b763cbea7770f99 Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Tue, 20 Jan 2026 20:31:37 -0600 Subject: [PATCH 4/9] switched to insert --- .../contact/SolidMechanicsAugmentedLagrangianContact.cpp | 7 ++----- .../contact/SolidMechanicsLagrangeContactBubbleStab.cpp | 7 ++----- 2 files changed, 4 insertions(+), 10 deletions(-) diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp index 6cbb712624f..54acc173fa4 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp @@ -227,11 +227,8 @@ void SolidMechanicsAugmentedLagrangianContact::postInputInitialization() FiniteElementDiscretization const & feDiscretization = feDiscretizationManager.getGroup< FiniteElementDiscretization >( getDiscretizationName() ); - m_faceTypeToFiniteElements["Quadrilateral"] = - feDiscretization.factory( ElementType::Quadrilateral ); - - m_faceTypeToFiniteElements["Triangle"] = - feDiscretization.factory( ElementType::Triangle ); + m_faceTypeToFiniteElements.insert( {"Quadrilateral", feDiscretization.factory( ElementType::Quadrilateral )} ); + m_faceTypeToFiniteElements.insert( {"Triangle", feDiscretization.factory( ElementType::Triangle )} ); } void SolidMechanicsAugmentedLagrangianContact::setSparsityPattern( DomainPartition & domain, diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp index 2271e0775f9..5b854b00b4d 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp @@ -200,11 +200,8 @@ void SolidMechanicsLagrangeContactBubbleStab::postInputInitialization() FiniteElementDiscretization const & feDiscretization = feDiscretizationManager.getGroup< FiniteElementDiscretization >( getDiscretizationName() ); - m_faceTypeToFiniteElements["Quadrilateral"] = - feDiscretization.factory( ElementType::Quadrilateral ); - - m_faceTypeToFiniteElements["Triangle"] = - feDiscretization.factory( ElementType::Triangle ); + m_faceTypeToFiniteElements.insert( {"Quadrilateral", feDiscretization.factory( ElementType::Quadrilateral )} ); + m_faceTypeToFiniteElements.insert( {"Triangle", feDiscretization.factory( ElementType::Triangle )} ); } void SolidMechanicsLagrangeContactBubbleStab::setupSystem( DomainPartition & domain, From 56f33db16d26e07544adebde829ad33cda1656d5 Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Wed, 21 Jan 2026 14:01:16 -0600 Subject: [PATCH 5/9] Added log info --- .../contact/SolidMechanicsAugmentedLagrangianContact.cpp | 8 ++++++++ .../contact/SolidMechanicsLagrangeContactBubbleStab.cpp | 8 ++++++++ 2 files changed, 16 insertions(+) diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp index 54acc173fa4..93fb332c2d0 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp @@ -229,6 +229,14 @@ void SolidMechanicsAugmentedLagrangianContact::postInputInitialization() m_faceTypeToFiniteElements.insert( {"Quadrilateral", feDiscretization.factory( ElementType::Quadrilateral )} ); m_faceTypeToFiniteElements.insert( {"Triangle", feDiscretization.factory( ElementType::Triangle )} ); + + GEOS_LOG_RANK_0( GEOS_FMT( "{} using finite element discretization {}:", + getName(), getDiscretizationName() ) ); + for( auto const & [name, fePtr] : m_faceTypeToFiniteElements ) + { + GEOS_LOG_RANK_0( GEOS_FMT( " {} face elements: {} quadrature points", + name, fePtr->getNumQuadraturePoints() ) ); + } } void SolidMechanicsAugmentedLagrangianContact::setSparsityPattern( DomainPartition & domain, diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp index 5b854b00b4d..05bd49ab6bb 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp @@ -202,6 +202,14 @@ void SolidMechanicsLagrangeContactBubbleStab::postInputInitialization() m_faceTypeToFiniteElements.insert( {"Quadrilateral", feDiscretization.factory( ElementType::Quadrilateral )} ); m_faceTypeToFiniteElements.insert( {"Triangle", feDiscretization.factory( ElementType::Triangle )} ); + + GEOS_LOG_RANK_0( GEOS_FMT( "{} using finite element discretization {}:", + getName(), getDiscretizationName() ) ); + for( auto const & [name, fePtr] : m_faceTypeToFiniteElements ) + { + GEOS_LOG_RANK_0( GEOS_FMT( " {} face elements: {} quadrature points", + name, fePtr->getNumQuadraturePoints() ) ); + } } void SolidMechanicsLagrangeContactBubbleStab::setupSystem( DomainPartition & domain, From e672590cbe76ad8b380c4d94e2975a76dc2ce12f Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Wed, 21 Jan 2026 18:27:24 -0600 Subject: [PATCH 6/9] Checked order for tet meshes --- ...lidMechanicsAugmentedLagrangianContact.cpp | 36 +++++++++++++++++++ 1 file changed, 36 insertions(+) diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp index 93fb332c2d0..06a79511a65 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp @@ -109,6 +109,42 @@ void SolidMechanicsAugmentedLagrangianContact::registerDataOnMesh( dataRepositor { ContactSolverBase::registerDataOnMesh( meshBodies ); + // Check for tetrahedral elements and validate high-order quadrature requirement + NumericalMethodsManager const & numericalMethodManager = + this->getGroupByPath< DomainPartition >( "/Problem/domain" ).getNumericalMethodManager(); + FiniteElementDiscretizationManager const & feDiscretizationManager = + numericalMethodManager.getFiniteElementDiscretizationManager(); + FiniteElementDiscretization const & feDiscretization = + feDiscretizationManager.getGroup< FiniteElementDiscretization >( getDiscretizationName() ); + + integer const useHighOrderQuadrature = + feDiscretization.getReference< integer >( "useHighOrderQuadratureRule" ); + + bool hasTetrahedra = false; + forDiscretizationOnMeshTargets( meshBodies, [&]( string const &, + MeshLevel const & mesh, + string_array const & regionNames ) + { + ElementRegionManager const & elemManager = mesh.getElemManager(); + elemManager.forElementRegions< CellElementRegion >( regionNames, [&]( localIndex const, + CellElementRegion const & region ) + { + region.forElementSubRegions< CellElementSubRegion >( [&]( CellElementSubRegion const & subRegion ) + { + if( subRegion.getElementType() == ElementType::Tetrahedron ) + { + hasTetrahedra = true; + } + } ); + } ); + } ); + + GEOS_ERROR_IF( hasTetrahedra && useHighOrderQuadrature != 1, + GEOS_FMT( "{}: Tetrahedral meshes require useHighOrderQuadratureRule=\"1\" for correct integration of bubble contributions. " + "Please add this attribute to your FiniteElements/{} XML block.", + getName(), getDiscretizationName() ) ); + + forDiscretizationOnMeshTargets( meshBodies, [&] ( string const &, MeshLevel & meshLevel, string_array const & ) From 8873133ee4cdc5eac710d942933a51e7bf1d00cb Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Wed, 21 Jan 2026 18:44:28 -0600 Subject: [PATCH 7/9] Checked order for tet meshes in bubbleStab solver --- ...olidMechanicsLagrangeContactBubbleStab.cpp | 36 +++++++++++++++++++ 1 file changed, 36 insertions(+) diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp index 05bd49ab6bb..394aa7f5ce9 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp @@ -99,6 +99,42 @@ void SolidMechanicsLagrangeContactBubbleStab::registerDataOnMesh( Group & meshBo { ContactSolverBase::registerDataOnMesh( meshBodies ); + // Check for tetrahedral elements and validate high-order quadrature requirement + NumericalMethodsManager const & numericalMethodManager = + this->getGroupByPath< DomainPartition >( "/Problem/domain" ).getNumericalMethodManager(); + FiniteElementDiscretizationManager const & feDiscretizationManager = + numericalMethodManager.getFiniteElementDiscretizationManager(); + FiniteElementDiscretization const & feDiscretization = + feDiscretizationManager.getGroup< FiniteElementDiscretization >( getDiscretizationName() ); + + integer const useHighOrderQuadrature = + feDiscretization.getReference< integer >( "useHighOrderQuadratureRule" ); + + bool hasTetrahedra = false; + forDiscretizationOnMeshTargets( meshBodies, [&]( string const &, + MeshLevel const & mesh, + string_array const & regionNames ) + { + ElementRegionManager const & elemManager = mesh.getElemManager(); + elemManager.forElementRegions< CellElementRegion >( regionNames, [&]( localIndex const, + CellElementRegion const & region ) + { + region.forElementSubRegions< CellElementSubRegion >( [&]( CellElementSubRegion const & subRegion ) + { + if( subRegion.getElementType() == ElementType::Tetrahedron ) + { + hasTetrahedra = true; + } + } ); + } ); + } ); + + GEOS_ERROR_IF( hasTetrahedra && useHighOrderQuadrature != 1, + GEOS_FMT( "{}: Tetrahedral meshes require useHighOrderQuadratureRule=\"1\" for correct integration of bubble contributions. " + "Please add this attribute to your FiniteElements/{} XML block.", + getName(), getDiscretizationName() ) ); + + forDiscretizationOnMeshTargets( meshBodies, [&] ( string const &, MeshLevel & meshLevel, string_array const & ) From d63a5cb43e9aac10d616bca88a9c5fed3eae53de Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Thu, 22 Jan 2026 13:16:46 -0600 Subject: [PATCH 8/9] Fixed schema generation --- ...lidMechanicsAugmentedLagrangianContact.cpp | 84 +++++++++++-------- ...lidMechanicsAugmentedLagrangianContact.hpp | 8 ++ ...olidMechanicsLagrangeContactBubbleStab.cpp | 82 ++++++++++-------- ...olidMechanicsLagrangeContactBubbleStab.hpp | 8 ++ 4 files changed, 111 insertions(+), 71 deletions(-) diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp index 06a79511a65..dbf13d324c8 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp @@ -109,42 +109,6 @@ void SolidMechanicsAugmentedLagrangianContact::registerDataOnMesh( dataRepositor { ContactSolverBase::registerDataOnMesh( meshBodies ); - // Check for tetrahedral elements and validate high-order quadrature requirement - NumericalMethodsManager const & numericalMethodManager = - this->getGroupByPath< DomainPartition >( "/Problem/domain" ).getNumericalMethodManager(); - FiniteElementDiscretizationManager const & feDiscretizationManager = - numericalMethodManager.getFiniteElementDiscretizationManager(); - FiniteElementDiscretization const & feDiscretization = - feDiscretizationManager.getGroup< FiniteElementDiscretization >( getDiscretizationName() ); - - integer const useHighOrderQuadrature = - feDiscretization.getReference< integer >( "useHighOrderQuadratureRule" ); - - bool hasTetrahedra = false; - forDiscretizationOnMeshTargets( meshBodies, [&]( string const &, - MeshLevel const & mesh, - string_array const & regionNames ) - { - ElementRegionManager const & elemManager = mesh.getElemManager(); - elemManager.forElementRegions< CellElementRegion >( regionNames, [&]( localIndex const, - CellElementRegion const & region ) - { - region.forElementSubRegions< CellElementSubRegion >( [&]( CellElementSubRegion const & subRegion ) - { - if( subRegion.getElementType() == ElementType::Tetrahedron ) - { - hasTetrahedra = true; - } - } ); - } ); - } ); - - GEOS_ERROR_IF( hasTetrahedra && useHighOrderQuadrature != 1, - GEOS_FMT( "{}: Tetrahedral meshes require useHighOrderQuadratureRule=\"1\" for correct integration of bubble contributions. " - "Please add this attribute to your FiniteElements/{} XML block.", - getName(), getDiscretizationName() ) ); - - forDiscretizationOnMeshTargets( meshBodies, [&] ( string const &, MeshLevel & meshLevel, string_array const & ) @@ -200,6 +164,54 @@ void SolidMechanicsAugmentedLagrangianContact::registerDataOnMesh( dataRepositor } +void SolidMechanicsAugmentedLagrangianContact::initializePostInitialConditionsPreSubGroups() +{ + ContactSolverBase::initializePostInitialConditionsPreSubGroups(); + + DomainPartition & domain = this->getGroupByPath< DomainPartition >( "/Problem/domain" ); + validateTetrahedralQuadrature( domain.getMeshBodies() ); +} + +void SolidMechanicsAugmentedLagrangianContact::validateTetrahedralQuadrature( Group & meshBodies ) +{ + string const discretizationName = getDiscretizationName(); + + NumericalMethodsManager const & numericalMethodManager = + this->getGroupByPath< DomainPartition >( "/Problem/domain" ).getNumericalMethodManager(); + FiniteElementDiscretizationManager const & feDiscretizationManager = + numericalMethodManager.getFiniteElementDiscretizationManager(); + FiniteElementDiscretization const & feDiscretization = + feDiscretizationManager.getGroup< FiniteElementDiscretization >( discretizationName ); + + + integer const useHighOrderQuadrature = + feDiscretization.getReference< integer >( "useHighOrderQuadratureRule" ); + + bool hasTetrahedra = false; + forDiscretizationOnMeshTargets( meshBodies, [&]( string const &, + MeshLevel const & mesh, + string_array const & regionNames ) + { + ElementRegionManager const & elemManager = mesh.getElemManager(); + elemManager.forElementRegions< CellElementRegion >( regionNames, [&]( localIndex const, + CellElementRegion const & region ) + { + region.forElementSubRegions< CellElementSubRegion >( [&]( CellElementSubRegion const & subRegion ) + { + if( subRegion.getElementType() == ElementType::Tetrahedron ) + { + hasTetrahedra = true; + } + } ); + } ); + } ); + + GEOS_ERROR_IF( hasTetrahedra && useHighOrderQuadrature != 1, + GEOS_FMT( "{}: Tetrahedral meshes require useHighOrderQuadratureRule=\"1\" for correct integration of bubble contributions. " + "Please add this attribute to your FiniteElements/{} XML block.", + getName(), discretizationName ) ); +} + void SolidMechanicsAugmentedLagrangianContact::setupDofs( DomainPartition const & domain, DofManager & dofManager ) const { diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp index 2ebed39073e..c3758556deb 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp @@ -51,6 +51,8 @@ class SolidMechanicsAugmentedLagrangianContact : public ContactSolverBase virtual void registerDataOnMesh( dataRepository::Group & meshBodies ) override final; + virtual void initializePostInitialConditionsPreSubGroups() override; + virtual void setupDofs( DomainPartition const & domain, DofManager & dofManager ) const override; @@ -215,6 +217,12 @@ class SolidMechanicsAugmentedLagrangianContact : public ContactSolverBase private: + /** + * @brief Validate that tetrahedral meshes use high-order quadrature rules + * @param meshBodies the group containing the mesh bodies + */ + void validateTetrahedralQuadrature( Group & meshBodies ); + /** * @brief add the number of non-zero elements induced by the coupling between * nodal and bubble displacement. diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp index 394aa7f5ce9..d73fb84619a 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp @@ -99,13 +99,57 @@ void SolidMechanicsLagrangeContactBubbleStab::registerDataOnMesh( Group & meshBo { ContactSolverBase::registerDataOnMesh( meshBodies ); - // Check for tetrahedral elements and validate high-order quadrature requirement + forDiscretizationOnMeshTargets( meshBodies, [&] ( string const &, + MeshLevel & meshLevel, + string_array const & ) + { + FaceManager & faceManager = meshLevel.getFaceManager(); + + // Register the total bubble displacement + faceManager.registerField< contact::totalBubbleDisplacement >( this->getName() ). + reference().resizeDimension< 1 >( 3 ); + + // Register the incremental bubble displacement + faceManager.registerField< contact::incrementalBubbleDisplacement >( this->getName() ). + reference().resizeDimension< 1 >( 3 ); + } ); + + forFractureRegionOnMeshTargets( meshBodies, [&] ( SurfaceElementRegion & fractureRegion ) + { + fractureRegion.forElementSubRegions< SurfaceElementSubRegion >( [&]( SurfaceElementSubRegion & subRegion ) + { + // Register the rotation matrix + subRegion.registerField< contact::rotationMatrix >( this->getName() ). + reference().resizeDimension< 1, 2 >( 3, 3 ); + + subRegion.registerField< contact::deltaTraction >( getName() ). + reference().resizeDimension< 1 >( 3 ); + + subRegion.registerField< contact::targetIncrementalJump >( getName() ). + reference().resizeDimension< 1 >( 3 ); + } ); + } ); +} + +void SolidMechanicsLagrangeContactBubbleStab::initializePostInitialConditionsPreSubGroups() +{ + ContactSolverBase::initializePostInitialConditionsPreSubGroups(); + + DomainPartition & domain = this->getGroupByPath< DomainPartition >( "/Problem/domain" ); + validateTetrahedralQuadrature( domain.getMeshBodies() ); +} + +void SolidMechanicsLagrangeContactBubbleStab::validateTetrahedralQuadrature( Group & meshBodies ) +{ + string const discretizationName = getDiscretizationName(); + NumericalMethodsManager const & numericalMethodManager = this->getGroupByPath< DomainPartition >( "/Problem/domain" ).getNumericalMethodManager(); FiniteElementDiscretizationManager const & feDiscretizationManager = numericalMethodManager.getFiniteElementDiscretizationManager(); FiniteElementDiscretization const & feDiscretization = - feDiscretizationManager.getGroup< FiniteElementDiscretization >( getDiscretizationName() ); + feDiscretizationManager.getGroup< FiniteElementDiscretization >( discretizationName ); + integer const useHighOrderQuadrature = feDiscretization.getReference< integer >( "useHighOrderQuadratureRule" ); @@ -132,39 +176,7 @@ void SolidMechanicsLagrangeContactBubbleStab::registerDataOnMesh( Group & meshBo GEOS_ERROR_IF( hasTetrahedra && useHighOrderQuadrature != 1, GEOS_FMT( "{}: Tetrahedral meshes require useHighOrderQuadratureRule=\"1\" for correct integration of bubble contributions. " "Please add this attribute to your FiniteElements/{} XML block.", - getName(), getDiscretizationName() ) ); - - - forDiscretizationOnMeshTargets( meshBodies, [&] ( string const &, - MeshLevel & meshLevel, - string_array const & ) - { - FaceManager & faceManager = meshLevel.getFaceManager(); - - // Register the total bubble displacement - faceManager.registerField< contact::totalBubbleDisplacement >( this->getName() ). - reference().resizeDimension< 1 >( 3 ); - - // Register the incremental bubble displacement - faceManager.registerField< contact::incrementalBubbleDisplacement >( this->getName() ). - reference().resizeDimension< 1 >( 3 ); - } ); - - forFractureRegionOnMeshTargets( meshBodies, [&] ( SurfaceElementRegion & fractureRegion ) - { - fractureRegion.forElementSubRegions< SurfaceElementSubRegion >( [&]( SurfaceElementSubRegion & subRegion ) - { - // Register the rotation matrix - subRegion.registerField< contact::rotationMatrix >( this->getName() ). - reference().resizeDimension< 1, 2 >( 3, 3 ); - - subRegion.registerField< contact::deltaTraction >( getName() ). - reference().resizeDimension< 1 >( 3 ); - - subRegion.registerField< contact::targetIncrementalJump >( getName() ). - reference().resizeDimension< 1 >( 3 ); - } ); - } ); + getName(), discretizationName ) ); } void SolidMechanicsLagrangeContactBubbleStab::setupDofs( DomainPartition const & domain, diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp index aaff54e8fff..d17b43a3c2d 100644 --- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp +++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp @@ -52,6 +52,8 @@ class SolidMechanicsLagrangeContactBubbleStab : public ContactSolverBase virtual void postInputInitialization() override; + virtual void initializePostInitialConditionsPreSubGroups() override; + virtual void registerDataOnMesh( Group & MeshBodies ) override final; real64 solverStep( real64 const & time_n, @@ -202,6 +204,12 @@ class SolidMechanicsLagrangeContactBubbleStab : public ContactSolverBase private: + /** + * @brief Validate that tetrahedral meshes use high-order quadrature rules + * @param meshBodies the group containing the mesh bodies + */ + void validateTetrahedralQuadrature( Group & meshBodies ); + /** * @brief add the number of non-zero elements induced by the coupling between * nodal and bubble displacement. From 15625fcc7bccc8e7a499eb7fff35eea0e7e8c2b7 Mon Sep 17 00:00:00 2001 From: DENEL Bertrand Date: Thu, 22 Jan 2026 16:45:28 -0600 Subject: [PATCH 9/9] Draft --- .github/workflows/ci_tests.yml | 15 +++++++++++++++ 1 file changed, 15 insertions(+) diff --git a/.github/workflows/ci_tests.yml b/.github/workflows/ci_tests.yml index 63b2924cfc8..e5cd9019559 100644 --- a/.github/workflows/ci_tests.yml +++ b/.github/workflows/ci_tests.yml @@ -415,6 +415,21 @@ jobs: DOCKER_CERTS_UPDATE_COMMAND: "update-ca-trust" HOST_CONFIG: /spack-generated-wave-solver-only.cmake + - name: Rockylinux CUDA on T4 (8, clang 17.0.6, cuda 12.9.1) + # GitHub-hosted runner with 1x T4 GPU (16GB VRAM) + BUILD_AND_TEST_CLI_ARGS: "--no-install-schema" + CMAKE_BUILD_TYPE: Release + BUILD_GENERATOR: "--ninja" + ENABLE_HYPRE_DEVICE: CUDA + ENABLE_HYPRE: ON + ENABLE_TRILINOS: OFF + GEOS_ENABLE_BOUNDS_CHECK: OFF + DOCKER_REPOSITORY: geosx/rockylinux8-clang17-cuda12.9.1 + RUNS_ON: GPU + NPROC: 8 + DOCKER_RUN_ARGS: "--cpus=8 --memory=128g --runtime=nvidia --gpus all" + HOST_CONFIG: /spack-generated.cmake + #- name: Sherlock GPU (centos 7.9.2009, gcc 10.1.0, open-mpi 4.1.2, openblas 0.3.10, cuda 12.4.0,) # BUILD_AND_TEST_CLI_ARGS: "--no-run-unit-tests --no-install-schema" # BUILD_GENERATOR: "--ninja"