diff --git a/src/coreComponents/physicsSolvers/PhysicsSolverBase.cpp b/src/coreComponents/physicsSolvers/PhysicsSolverBase.cpp
index 7448c5ec46d..77abdd58dab 100644
--- a/src/coreComponents/physicsSolvers/PhysicsSolverBase.cpp
+++ b/src/coreComponents/physicsSolvers/PhysicsSolverBase.cpp
@@ -144,9 +144,11 @@ void PhysicsSolverBase::postInputInitialization()
                                            m_writeStatisticsCSV == StatsOutputType::all );
   getConvergenceStats().setCSVOutputRequest( m_writeStatisticsCSV == StatsOutputType::convergence ||
                                              m_writeStatisticsCSV == StatsOutputType::all );
-}
 
-PhysicsSolverBase::~PhysicsSolverBase() = default;
+  LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
+  if( linearSolverParameters.preconditionerType == LinearSolverParameters::PreconditionerType::mgr )
+    setMGRStrategy();
+}
 
 void PhysicsSolverBase::initialize_postMeshGeneration()
 {
@@ -1388,7 +1390,7 @@ void PhysicsSolverBase::solveLinearSystem( DofManager const & dofManager,
 
   matrix.setDofManager( &dofManager );
 
-  GEOS_WARNING_IF( isDirectSolver && dofManager.numGlobalDofs() > 100000,
+  GEOS_WARNING_IF( logger::internal::rank == 0 && isDirectSolver && dofManager.numGlobalDofs() > 100000,
                    "Direct solver used for large system ( > 100,000 DOFs ). "
                    "This may lead to high memory consumption and long computation times. "
                    "Consider using an iterative solver for better performance." );
@@ -1696,6 +1698,11 @@ bool PhysicsSolverBase::detectOscillations() const
   return f > oscillationFraction;
 }
 
+void PhysicsSolverBase::setMGRStrategy()
+{
+  GEOS_ERROR( GEOS_FMT( "{}: MGR strategy is not implemented for {}", this->getName(), this->getCatalogName()));
+}
+
 #if defined(GEOS_USE_PYGEOSX)
 PyTypeObject * PhysicsSolverBase::getPythonType() const
 { return python::getPySolverType(); }
diff --git a/src/coreComponents/physicsSolvers/PhysicsSolverBase.hpp b/src/coreComponents/physicsSolvers/PhysicsSolverBase.hpp
index b02c106a751..10e59d610ee 100644
--- a/src/coreComponents/physicsSolvers/PhysicsSolverBase.hpp
+++ b/src/coreComponents/physicsSolvers/PhysicsSolverBase.hpp
@@ -73,11 +73,6 @@ class PhysicsSolverBase : public ExecutableGroup
    */
   PhysicsSolverBase( PhysicsSolverBase && ) = default;
 
-  /**
-   * @brief Destructor for PhysicsSolverBase
-   */
-  virtual ~PhysicsSolverBase() override;
-
   /**
    * @brief Deleted constructor
    */
@@ -1005,6 +1000,11 @@ class PhysicsSolverBase : public ExecutableGroup
 
   virtual void postInputInitialization() override;
 
+  /**
+   * @brief Set MGR strategy
+   */
+  virtual void setMGRStrategy();
+
   /**
    * @brief Eisenstat-Walker adaptive tolerance
    *
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp
index 4f2f9fe43c1..33bf8b5e411 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp
@@ -194,10 +194,6 @@ void CompositionalMultiphaseFVM::initializePreSubGroups()
                           EnumStrings< ScalingType >::toString( ScalingType::Local )) );
   }
 
-  m_linearSolverParameters.get().mgr.strategy = m_isThermal
-                                                ? LinearSolverParameters::MGR::StrategyType::thermalCompositionalMultiphaseFVM
-                                                : LinearSolverParameters::MGR::StrategyType::compositionalMultiphaseFVM;
-
   checkDiscretizationName();
 
   DomainPartition & domain = this->getGroupByPath< DomainPartition >( "/Problem/domain" );
@@ -209,6 +205,13 @@ void CompositionalMultiphaseFVM::initializePreSubGroups()
                            getName(), EnumStrings< UpwindingScheme >::toString( UpwindingScheme::HU2PH )));
 }
 
+void CompositionalMultiphaseFVM::setMGRStrategy()
+{
+  m_linearSolverParameters.get().mgr.strategy = m_isThermal
+                                                    ? LinearSolverParameters::MGR::StrategyType::thermalCompositionalMultiphaseFVM
+                                                    : LinearSolverParameters::MGR::StrategyType::compositionalMultiphaseFVM;
+}
+
 void CompositionalMultiphaseFVM::setupDofs( DomainPartition const & domain,
                                             DofManager & dofManager ) const
 {
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.hpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.hpp
index b08b523a353..bf38a9a8ec1 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.hpp
@@ -206,6 +206,8 @@ class CompositionalMultiphaseFVM : public CompositionalMultiphaseBase
 
   virtual void initializePreSubGroups() override;
 
+  virtual void setMGRStrategy() override;
+
   real64 setNextDtBasedOnCFL( real64 const & currentDt,
                               DomainPartition & domain );
 
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.cpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.cpp
index 11b4fead408..eaad4625e02 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.cpp
@@ -53,6 +53,9 @@ CompositionalMultiphaseHybridFVM::CompositionalMultiphaseHybridFVM( const std::s
                                                                     Group * const parent ):
   CompositionalMultiphaseBase( name, parent ),
   m_lengthTolerance( 0 )
+{}
+
+void CompositionalMultiphaseHybridFVM::setMGRStrategy()
 {
   m_linearSolverParameters.get().mgr.strategy = LinearSolverParameters::MGR::StrategyType::compositionalMultiphaseHybridFVM;
 }
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.hpp b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.hpp
index e0be3309446..5f64eb0f542 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseHybridFVM.hpp
@@ -161,7 +161,7 @@ class CompositionalMultiphaseHybridFVM : public CompositionalMultiphaseBase
 
   struct viewKeyStruct : CompositionalMultiphaseBase::viewKeyStruct
   {
-    static constexpr char const * faceDofFieldString() { return "faceCenteredVariables"; }
+    static constexpr char const *faceDofFieldString() { return "faceCenteredVariables"; }
   };
 
   virtual void initializePostInitialConditionsPreSubGroups() override;
@@ -170,6 +170,8 @@ class CompositionalMultiphaseHybridFVM : public CompositionalMultiphaseBase
 
 protected:
 
+  virtual void setMGRStrategy() override;
+
   /// precompute the minGravityCoefficient for the buoyancy term
   void precomputeData( MeshLevel & mesh, string_array const & regionNames ) override;
 
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/ImmiscibleMultiphaseFlow.cpp b/src/coreComponents/physicsSolvers/fluidFlow/ImmiscibleMultiphaseFlow.cpp
index 389a33357b7..549ad565900 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/ImmiscibleMultiphaseFlow.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/ImmiscibleMultiphaseFlow.cpp
@@ -93,11 +93,6 @@ ImmiscibleMultiphaseFlow::ImmiscibleMultiphaseFlow( const string & name,
     setDescription( "Target (absolute) change in phase volume fraction in a time step" );
 }
 
-void ImmiscibleMultiphaseFlow::postInputInitialization()
-{
-  FlowSolverBase::postInputInitialization();
-}
-
 void ImmiscibleMultiphaseFlow::registerDataOnMesh( Group & meshBodies )
 {
   FlowSolverBase::registerDataOnMesh( meshBodies );
@@ -185,8 +180,6 @@ void ImmiscibleMultiphaseFlow::setConstitutiveNames( ElementSubRegionBase & subR
 
 void ImmiscibleMultiphaseFlow::initializePreSubGroups()
 {
-  m_linearSolverParameters.get().mgr.strategy = LinearSolverParameters::MGR::StrategyType::immiscibleMultiphaseFVM;
-
   FlowSolverBase::initializePreSubGroups();
 
   DomainPartition & domain = this->getGroupByPath< DomainPartition >( "/Problem/domain" );
@@ -205,6 +198,10 @@ void ImmiscibleMultiphaseFlow::initializePreSubGroups()
   } );
 }
 
+void ImmiscibleMultiphaseFlow::setMGRStrategy()
+{
+  m_linearSolverParameters.get().mgr.strategy = LinearSolverParameters::MGR::StrategyType::immiscibleMultiphaseFVM;
+}
 
 void ImmiscibleMultiphaseFlow::updateFluidModel( ObjectManagerBase & dataGroup ) const
 {
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/ImmiscibleMultiphaseFlow.hpp b/src/coreComponents/physicsSolvers/fluidFlow/ImmiscibleMultiphaseFlow.hpp
index 7526130316c..53e34198a6f 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/ImmiscibleMultiphaseFlow.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/ImmiscibleMultiphaseFlow.hpp
@@ -233,10 +233,11 @@ class ImmiscibleMultiphaseFlow : public FlowSolverBase
     static constexpr char const * elemDofFieldString() { return "elemDofField"; }
   };
 
+protected:
 
-private:
+  virtual void setMGRStrategy() override;
 
-  virtual void postInputInitialization() override;
+private:
 
   /**
    * @brief Update all relevant fluid models using current values of pressure and phase volume fraction
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.cpp b/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.cpp
index 345fe31f8fa..448e2f6e23f 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.cpp
@@ -119,7 +119,10 @@ ReactiveCompositionalMultiphaseOBL::ReactiveCompositionalMultiphaseOBL( const st
     setRTTypeName( rtTypes::CustomTypes::groupNameRefArray ).
     setInputFlag( InputFlags::OPTIONAL ).
     setDescription( "List of fluid phases" );
+}
 
+void ReactiveCompositionalMultiphaseOBL::setMGRStrategy()
+{
   m_linearSolverParameters.get().mgr.strategy = LinearSolverParameters::MGR::StrategyType::reactiveCompositionalMultiphaseOBL;
 }
 
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.hpp b/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.hpp
index c50d111d2e0..32ee709c24c 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/ReactiveCompositionalMultiphaseOBL.hpp
@@ -308,6 +308,9 @@ class ReactiveCompositionalMultiphaseOBL : public FlowSolverBase
 
   virtual void initializePreSubGroups() override;
 
+protected:
+
+  virtual void setMGRStrategy() override;
 
 private:
 
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseHybridFVM.cpp b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseHybridFVM.cpp
index 4b28b4a2780..6009715de8c 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseHybridFVM.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseHybridFVM.cpp
@@ -51,9 +51,12 @@ SinglePhaseHybridFVM::SinglePhaseHybridFVM( const string & name,
   SinglePhaseBase( name, parent ),
   m_areaRelTol( 1e-8 )
 {
-
   // one cell-centered dof per cell
   m_numDofPerCell = 1;
+}
+
+void SinglePhaseHybridFVM::setMGRStrategy()
+{
   m_linearSolverParameters.get().mgr.strategy = LinearSolverParameters::MGR::StrategyType::singlePhaseHybridFVM;
 }
 
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseHybridFVM.hpp b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseHybridFVM.hpp
index 90b08cd2271..06cf318a96d 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseHybridFVM.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseHybridFVM.hpp
@@ -44,7 +44,6 @@ class SinglePhaseHybridFVM : public SinglePhaseBase
   SinglePhaseHybridFVM( const string & name,
                         Group * const parent );
 
-
   /// deleted default constructor
   SinglePhaseHybridFVM() = delete;
 
@@ -70,7 +69,9 @@ class SinglePhaseHybridFVM : public SinglePhaseBase
    * @return string that contains the catalog name to generate a new NodeManager object through the object catalog.
    */
   static string catalogName()
-  { return "SinglePhaseHybridFVM"; }
+  {
+    return "SinglePhaseHybridFVM";
+  }
   /**
    * @copydoc PhysicsSolverBase::getCatalogName()
    */
@@ -185,6 +186,10 @@ class SinglePhaseHybridFVM : public SinglePhaseBase
 
   virtual void initializePostInitialConditionsPreSubGroups() override;
 
+protected:
+
+  virtual void setMGRStrategy() override;
+
 private:
 
   /// relative tolerance (redundant with FluxApproximationBase)
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/proppantTransport/ProppantTransport.cpp b/src/coreComponents/physicsSolvers/fluidFlow/proppantTransport/ProppantTransport.cpp
index e44833d58aa..fdc8a5d661c 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/proppantTransport/ProppantTransport.cpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/proppantTransport/ProppantTransport.cpp
@@ -83,11 +83,6 @@ ProppantTransport::ProppantTransport( const string & name,
     setDescription( "Flag that enables/disables proppant-packing update" );
 }
 
-void ProppantTransport::postInputInitialization()
-{
-  FlowSolverBase::postInputInitialization();
-}
-
 void ProppantTransport::registerDataOnMesh( Group & meshBodies )
 {
   FlowSolverBase::registerDataOnMesh( meshBodies );
diff --git a/src/coreComponents/physicsSolvers/fluidFlow/proppantTransport/ProppantTransport.hpp b/src/coreComponents/physicsSolvers/fluidFlow/proppantTransport/ProppantTransport.hpp
index 4b53d67075c..6ae89310673 100644
--- a/src/coreComponents/physicsSolvers/fluidFlow/proppantTransport/ProppantTransport.hpp
+++ b/src/coreComponents/physicsSolvers/fluidFlow/proppantTransport/ProppantTransport.hpp
@@ -215,10 +215,6 @@ class ProppantTransport : public FlowSolverBase
    */
   void updateState( ObjectManagerBase & dataGroup );
 
-protected:
-
-  virtual void postInputInitialization() override;
-
 private:
 
   /**
diff --git a/src/coreComponents/physicsSolvers/inducedSeismicity/QuasiDynamicEarthQuake.cpp b/src/coreComponents/physicsSolvers/inducedSeismicity/QuasiDynamicEarthQuake.cpp
index 17be5d463b5..8c6f36e1b27 100644
--- a/src/coreComponents/physicsSolvers/inducedSeismicity/QuasiDynamicEarthQuake.cpp
+++ b/src/coreComponents/physicsSolvers/inducedSeismicity/QuasiDynamicEarthQuake.cpp
@@ -49,11 +49,10 @@ QuasiDynamicEarthQuake< RSSOLVER_TYPE >::QuasiDynamicEarthQuake( const string &
 template< typename RSSOLVER_TYPE >
 void QuasiDynamicEarthQuake< RSSOLVER_TYPE >::postInputInitialization()
 {
+  PhysicsSolverBase::postInputInitialization();
 
   // Initialize member stress solver as specified in XML input
   m_stressSolver = &this->getParent().template getGroup< PhysicsSolverBase >( m_stressSolverName );
-
-  PhysicsSolverBase::postInputInitialization();
 }
 
 template< typename RSSOLVER_TYPE >
diff --git a/src/coreComponents/physicsSolvers/inducedSeismicity/SeismicityRate.cpp b/src/coreComponents/physicsSolvers/inducedSeismicity/SeismicityRate.cpp
index f614a96c67b..f075ebf6686 100644
--- a/src/coreComponents/physicsSolvers/inducedSeismicity/SeismicityRate.cpp
+++ b/src/coreComponents/physicsSolvers/inducedSeismicity/SeismicityRate.cpp
@@ -70,6 +70,8 @@ SeismicityRate::SeismicityRate( const string & name,
 
 void SeismicityRate::postInputInitialization()
 {
+  PhysicsSolverBase::postInputInitialization();
+
   // Check orthogonality of user-specified faults
   if( std::abs( LvArray::tensorOps::AiBi< 3 >( m_faultNormalDirection, m_faultShearDirection )) > 1e-8 )
   {
@@ -87,8 +89,6 @@ void SeismicityRate::postInputInitialization()
   {
     m_stressSolver = &this->getParent().getGroup< PhysicsSolverBase >( m_stressSolverName );
   }
-
-  PhysicsSolverBase::postInputInitialization();
 }
 
 SeismicityRate::~SeismicityRate()
diff --git a/src/coreComponents/physicsSolvers/multiphysics/CompositionalMultiphaseReservoirAndWells.cpp b/src/coreComponents/physicsSolvers/multiphysics/CompositionalMultiphaseReservoirAndWells.cpp
index b95f3493187..671bb4cb1ab 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/CompositionalMultiphaseReservoirAndWells.cpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/CompositionalMultiphaseReservoirAndWells.cpp
@@ -77,9 +77,6 @@ setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
@@ -118,9 +115,6 @@ setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/CoupledReservoirAndWellsBase.hpp b/src/coreComponents/physicsSolvers/multiphysics/CoupledReservoirAndWellsBase.hpp
index daf2baa05bb..b91b68e5914 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/CoupledReservoirAndWellsBase.hpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/CoupledReservoirAndWellsBase.hpp
@@ -183,8 +183,6 @@ class CoupledReservoirAndWellsBase : public CoupledSolver< RESERVOIR_SOLVER, WEL
 
     // assume that reservoir solver discretization is the primary one
     this->m_discretizationName = reservoirSolver()->getDiscretizationName();
-
-    setMGRStrategy();
   }
 
   virtual void
@@ -282,12 +280,6 @@ class CoupledReservoirAndWellsBase : public CoupledSolver< RESERVOIR_SOLVER, WEL
                               DofManager const & dofManager,
                               SparsityPatternView< globalIndex > const & pattern ) const = 0;
 
-  virtual void setMGRStrategy()
-  {
-    if( this->m_linearSolverParameters.get().preconditionerType == LinearSolverParameters::PreconditionerType::mgr )
-      GEOS_ERROR( GEOS_FMT( "{}: MGR strategy is not implemented for {}", this->getName(), this->getCatalogName()));
-  }
-
   /// Flag to determine whether the well transmissibility needs to be computed
   bool m_isWellTransmissibilityComputed;
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/CoupledSolver.hpp b/src/coreComponents/physicsSolvers/multiphysics/CoupledSolver.hpp
index bb485b2fade..5bbecd0ea54 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/CoupledSolver.hpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/CoupledSolver.hpp
@@ -721,7 +721,7 @@ class CoupledSolver : public PhysicsSolverBase
   virtual void
   postInputInitialization() override
   {
-    setSubSolvers();
+    setSubSolvers(); // has to go first
 
     PhysicsSolverBase::postInputInitialization();
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/HydrofractureSolver.cpp b/src/coreComponents/physicsSolvers/multiphysics/HydrofractureSolver.cpp
index e36249bbc76..1708575a29d 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/HydrofractureSolver.cpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/HydrofractureSolver.cpp
@@ -92,9 +92,6 @@ void HydrofractureSolver< POROMECHANICS_SOLVER >::setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = this->m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
@@ -223,8 +220,6 @@ void HydrofractureSolver< POROMECHANICS_SOLVER >::postInputInitialization()
 {
   Base::postInputInitialization();
 
-  setMGRStrategy();
-
   static const std::set< integer > binaryOptions = { 0, 1 };
   GEOS_ERROR_IF( binaryOptions.count( m_isMatrixPoroelastic ) == 0, viewKeyStruct::isMatrixPoroelasticString() << " option can be either 0 (false) or 1 (true)" );
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/MultiphasePoromechanics.cpp b/src/coreComponents/physicsSolvers/multiphysics/MultiphasePoromechanics.cpp
index 075df069e75..f97cb9d28fe 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/MultiphasePoromechanics.cpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/MultiphasePoromechanics.cpp
@@ -176,9 +176,6 @@ void MultiphasePoromechanics<>::setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = this->m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
@@ -210,9 +207,6 @@ void MultiphasePoromechanics< CompositionalMultiphaseReservoirAndWells<>, SolidM
 
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/MultiphasePoromechanics.hpp b/src/coreComponents/physicsSolvers/multiphysics/MultiphasePoromechanics.hpp
index de5de39883a..2020dde6d60 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/MultiphasePoromechanics.hpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/MultiphasePoromechanics.hpp
@@ -103,8 +103,7 @@ class MultiphasePoromechanics : public PoromechanicsSolver< FLOW_SOLVER, MECHANI
 
   virtual void setMGRStrategy() override
   {
-    if( this->m_linearSolverParameters.get().preconditionerType == LinearSolverParameters::PreconditionerType::mgr )
-      GEOS_ERROR( GEOS_FMT( "{}: MGR strategy is not implemented for {}", this->getName(), this->getCatalogName()));
+    Base::setMGRStrategy();
   }
 
   /**
diff --git a/src/coreComponents/physicsSolvers/multiphysics/OneWayCoupledFractureFlowContactMechanics.cpp b/src/coreComponents/physicsSolvers/multiphysics/OneWayCoupledFractureFlowContactMechanics.cpp
index e92dea5c400..e04399177da 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/OneWayCoupledFractureFlowContactMechanics.cpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/OneWayCoupledFractureFlowContactMechanics.cpp
@@ -38,12 +38,12 @@ OneWayCoupledFractureFlowContactMechanics< FLOW_SOLVER >::OneWayCoupledFractureF
 template< typename FLOW_SOLVER >
 void OneWayCoupledFractureFlowContactMechanics< FLOW_SOLVER >::postInputInitialization()
 {
+  Base::postInputInitialization();
+
   bool const isSequential = this->getNonlinearSolverParameters().couplingType() == NonlinearSolverParameters::CouplingType::Sequential;
   GEOS_THROW_IF( !isSequential,
                  "Only sequential coupling is allowed for this solver.",
                  InputError );
-
-  Base::postInputInitialization();
 }
 
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/PoromechanicsSolver.hpp b/src/coreComponents/physicsSolvers/multiphysics/PoromechanicsSolver.hpp
index f553759e327..94ad1c85609 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/PoromechanicsSolver.hpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/PoromechanicsSolver.hpp
@@ -126,8 +126,6 @@ class PoromechanicsSolver : public CoupledSolver< FLOW_SOLVER, MECHANICS_SOLVER
                              this->solidMechanicsSolver()->getName(),
                              EnumStrings< SolidMechanicsLagrangianFEM::TimeIntegrationOption >::toString( SolidMechanicsLagrangianFEM::TimeIntegrationOption::QuasiStatic ) ),
                    InputError );
-
-    setMGRStrategy();
   }
 
   virtual void setConstitutiveNamesCallSuper( ElementSubRegionBase & subRegion ) const override final
@@ -737,8 +735,6 @@ class PoromechanicsSolver : public CoupledSolver< FLOW_SOLVER, MECHANICS_SOLVER
     } );
   }
 
-  virtual void setMGRStrategy() = 0;
-
   virtual void updateBulkDensity( ElementSubRegionBase & subRegion ) = 0;
 
   virtual string getFlowDofKey() const = 0;
diff --git a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanics.cpp b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanics.cpp
index 7a3fdc6386a..38326f249ee 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanics.cpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanics.cpp
@@ -101,9 +101,6 @@ void SinglePhasePoromechanics<>::setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = this->m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
@@ -142,9 +139,6 @@ void SinglePhasePoromechanics< SinglePhaseReservoirAndWells<>, SolidMechanicsLag
 
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanics.hpp b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanics.hpp
index d7e40171f79..3986d21b155 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanics.hpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanics.hpp
@@ -118,8 +118,7 @@ class SinglePhasePoromechanics : public PoromechanicsSolver< FLOW_SOLVER, MECHAN
 
   virtual void setMGRStrategy() override
   {
-    if( this->m_linearSolverParameters.get().preconditionerType == LinearSolverParameters::PreconditionerType::mgr )
-      GEOS_ERROR( GEOS_FMT( "{}: MGR strategy is not implemented for {}", this->getName(), this->getCatalogName()));
+    Base::setMGRStrategy();
   }
 
   virtual void mapSolutionBetweenSolvers( DomainPartition & domain, integer const solverType ) override
diff --git a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsConformingFractures.cpp b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsConformingFractures.cpp
index 947a0a6ef4b..74af1be7547 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsConformingFractures.cpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsConformingFractures.cpp
@@ -37,9 +37,6 @@ void SinglePhasePoromechanicsConformingFractures<>::setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = this->m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsConformingFractures.hpp b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsConformingFractures.hpp
index f0ddd9bc8f1..250bca3bc89 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsConformingFractures.hpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsConformingFractures.hpp
@@ -84,8 +84,7 @@ class SinglePhasePoromechanicsConformingFractures : public PoromechanicsConformi
 
   virtual void setMGRStrategy() override
   {
-    if( this->m_linearSolverParameters.get().preconditionerType == LinearSolverParameters::PreconditionerType::mgr )
-      GEOS_ERROR( GEOS_FMT( "{}: MGR strategy is not implemented for {}", this->getName(), this->getCatalogName()));
+    Base::setMGRStrategy();
   }
 
   virtual void assembleSystem( real64 const time_n,
diff --git a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsEmbeddedFractures.cpp b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsEmbeddedFractures.cpp
index a327d48e569..9bf6db904be 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsEmbeddedFractures.cpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsEmbeddedFractures.cpp
@@ -42,15 +42,10 @@ SinglePhasePoromechanicsEmbeddedFractures::SinglePhasePoromechanicsEmbeddedFract
   SinglePhasePoromechanics( name, parent )
 { }
 
-SinglePhasePoromechanicsEmbeddedFractures::~SinglePhasePoromechanicsEmbeddedFractures() = default;
-
 void SinglePhasePoromechanicsEmbeddedFractures::setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsEmbeddedFractures.hpp b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsEmbeddedFractures.hpp
index 4b7d6fa3437..9a60e1b9577 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsEmbeddedFractures.hpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/SinglePhasePoromechanicsEmbeddedFractures.hpp
@@ -37,7 +37,6 @@ class SinglePhasePoromechanicsEmbeddedFractures : public SinglePhasePoromechanic
 
   SinglePhasePoromechanicsEmbeddedFractures( const std::string & name,
                                              Group * const parent );
-  ~SinglePhasePoromechanicsEmbeddedFractures() override;
 
   /**
    * @brief name of the node manager in the object catalog
diff --git a/src/coreComponents/physicsSolvers/multiphysics/SinglePhaseReservoirAndWells.cpp b/src/coreComponents/physicsSolvers/multiphysics/SinglePhaseReservoirAndWells.cpp
index 46cf3257088..bddc93c452b 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/SinglePhaseReservoirAndWells.cpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/SinglePhaseReservoirAndWells.cpp
@@ -72,9 +72,6 @@ setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
@@ -98,9 +95,6 @@ setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
diff --git a/src/coreComponents/physicsSolvers/multiphysics/SinglePhaseReservoirAndWells.hpp b/src/coreComponents/physicsSolvers/multiphysics/SinglePhaseReservoirAndWells.hpp
index 0adeee58026..040d3da4a9e 100644
--- a/src/coreComponents/physicsSolvers/multiphysics/SinglePhaseReservoirAndWells.hpp
+++ b/src/coreComponents/physicsSolvers/multiphysics/SinglePhaseReservoirAndWells.hpp
@@ -116,8 +116,7 @@ class SinglePhaseReservoirAndWells : public CoupledReservoirAndWellsBase< RESERV
 
   virtual void setMGRStrategy() override
   {
-    if( this->m_linearSolverParameters.get().preconditionerType == LinearSolverParameters::PreconditionerType::mgr )
-      GEOS_ERROR( GEOS_FMT( "{}: MGR strategy is not implemented for {}", this->getName(), this->getCatalogName()));
+    Base::setMGRStrategy();
   }
 
 private:
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/ContactSolverBase.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/ContactSolverBase.hpp
index ccfc68f6ea7..df94b2c5046 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/ContactSolverBase.hpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/ContactSolverBase.hpp
@@ -33,8 +33,6 @@ class ContactSolverBase : public SolidMechanicsLagrangianFEM
   ContactSolverBase( const string & name,
                      Group * const parent );
 
-  ~ContactSolverBase() override = default;
-
   virtual void registerDataOnMesh( dataRepository::Group & meshBodies ) override;
 
   virtual real64
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp
index 9a2291600ab..6173cde4d59 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.cpp
@@ -86,7 +86,10 @@ SolidMechanicsAugmentedLagrangianContact::SolidMechanicsAugmentedLagrangianConta
     setInputFlag( InputFlags::OPTIONAL ).
     setApplyDefaultValue( 5.e-02 ).
     setDescription( "Tolerance for the sliding check" );
+}
 
+void SolidMechanicsAugmentedLagrangianContact::setMGRStrategy()
+{
   // Set the default linear solver parameters
   LinearSolverParameters & linSolParams = m_linearSolverParameters.get();
 
@@ -100,9 +103,6 @@ SolidMechanicsAugmentedLagrangianContact::SolidMechanicsAugmentedLagrangianConta
   linSolParams.mgr.separateComponents = true;
 }
 
-SolidMechanicsAugmentedLagrangianContact::~SolidMechanicsAugmentedLagrangianContact()
-{}
-
 void SolidMechanicsAugmentedLagrangianContact::registerDataOnMesh( dataRepository::Group & meshBodies )
 {
   ContactSolverBase::registerDataOnMesh( meshBodies );
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp
index 390e8cb3802..d9af264f2c7 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsAugmentedLagrangianContact.hpp
@@ -32,8 +32,6 @@ class SolidMechanicsAugmentedLagrangianContact : public ContactSolverBase
   SolidMechanicsAugmentedLagrangianContact( const string & name,
                                             Group * const parent );
 
-  ~SolidMechanicsAugmentedLagrangianContact() override;
-
   /**
    * @brief name of the node manager in the object catalog
    * @return string that contains the catalog name to generate a new NodeManager object through the object catalog.
@@ -198,6 +196,10 @@ class SolidMechanicsAugmentedLagrangianContact : public ContactSolverBase
    */
   void createBubbleCellList( DomainPartition & domain ) const;
 
+protected:
+
+  virtual void setMGRStrategy() override;
+
 private:
 
   /**
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsEmbeddedFractures.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsEmbeddedFractures.cpp
index 0fe6810a514..2ac7625a721 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsEmbeddedFractures.cpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsEmbeddedFractures.cpp
@@ -71,19 +71,12 @@ void SolidMechanicsEmbeddedFractures::postInputInitialization()
     linearSolverParameters.amg.separateComponents = true;
     linearSolverParameters.dofsPerNode = 3;
   }
-  else
-  {
-    setMGRStrategy();
-  }
 }
 
 void SolidMechanicsEmbeddedFractures::setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsEmbeddedFractures.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsEmbeddedFractures.hpp
index 2834db10804..ef2dcd0a576 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsEmbeddedFractures.hpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsEmbeddedFractures.hpp
@@ -132,7 +132,7 @@ class SolidMechanicsEmbeddedFractures : public ContactSolverBase
 
   virtual void postInputInitialization() override final;
 
-  void setMGRStrategy();
+  virtual void setMGRStrategy() override;
 
 private:
 
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContact.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContact.cpp
index 5b67c829b29..05221ae628a 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContact.cpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContact.cpp
@@ -82,20 +82,10 @@ SolidMechanicsLagrangeContact::SolidMechanicsLagrangeContact( const string & nam
     setDescription( "Buffer parameter for local yield acceleration." );
 }
 
-void SolidMechanicsLagrangeContact::postInputInitialization()
-{
-  ContactSolverBase::postInputInitialization();
-
-  setMGRStrategy();
-}
-
 void SolidMechanicsLagrangeContact::setMGRStrategy()
 {
   LinearSolverParameters & linearSolverParameters = m_linearSolverParameters.get();
 
-  if( linearSolverParameters.preconditionerType != LinearSolverParameters::PreconditionerType::mgr )
-    return;
-
   linearSolverParameters.mgr.separateComponents = true;
   linearSolverParameters.dofsPerNode = 3;
 
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContact.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContact.hpp
index 674a08876c1..07a4792dfc1 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContact.hpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContact.hpp
@@ -175,9 +175,7 @@ class SolidMechanicsLagrangeContact : public ContactSolverBase
                                        DofManager const & dofManager,
                                        arrayView1d< real64 const > const & localRhs );
 
-  virtual void postInputInitialization() override final;
-
-  void setMGRStrategy();
+  virtual void setMGRStrategy() override;
 
 private:
   string m_stabilizationName;
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp
index ff06db81bfb..b3c30d81f04 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.cpp
@@ -46,17 +46,15 @@ SolidMechanicsLagrangeContactBubbleStab::SolidMechanicsLagrangeContactBubbleStab
   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;
-  linSolParams.dofsPerNode = 3;
-
   addLogLevel< logInfo::ResidualNorm >();
 }
 
-SolidMechanicsLagrangeContactBubbleStab::~SolidMechanicsLagrangeContactBubbleStab()
+void SolidMechanicsLagrangeContactBubbleStab::setMGRStrategy()
 {
-  // TODO Auto-generated destructor stub
+  LinearSolverParameters & linSolParams = m_linearSolverParameters.get();
+  linSolParams.mgr.strategy = LinearSolverParameters::MGR::StrategyType::lagrangianContactMechanicsBubbleStab;
+  linSolParams.mgr.separateComponents = true;
+  linSolParams.dofsPerNode = 3;
 }
 
 real64 SolidMechanicsLagrangeContactBubbleStab::solverStep( real64 const & time_n,
diff --git a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp
index 6c8562044e4..b7b7ea7168a 100644
--- a/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp
+++ b/src/coreComponents/physicsSolvers/solidMechanics/contact/SolidMechanicsLagrangeContactBubbleStab.hpp
@@ -31,12 +31,9 @@ class NumericalMethodsManager;
 class SolidMechanicsLagrangeContactBubbleStab : public ContactSolverBase
 {
 public:
-
   SolidMechanicsLagrangeContactBubbleStab( const string & name,
                                            Group * const parent );
 
-  ~SolidMechanicsLagrangeContactBubbleStab() override;
-
   /**
    * @brief name of the node manager in the object catalog
    * @return string that contains the catalog name to generate a new NodeManager object through the object catalog.
@@ -45,6 +42,7 @@ class SolidMechanicsLagrangeContactBubbleStab : public ContactSolverBase
   {
     return "SolidMechanicsLagrangeContactBubbleStab";
   }
+
   /**
    * @copydoc SolverBase::getCatalogName()
    */
@@ -198,6 +196,9 @@ class SolidMechanicsLagrangeContactBubbleStab : public ContactSolverBase
    */
   void computeRotationMatrices( DomainPartition & domain ) const;
 
+protected:
+
+  virtual void setMGRStrategy() override;
 
 private:
   /**