Linear algebra updates #2427

klevzoff · 2023-05-09T09:26:42Z

This is a collection of LAI-related changes extracted from my multiscale branch to be reviewed separately. I will post some comments later highlighting the most important changes.

Rebaseline needed due to addition of new input parameters.

TODO:

Update TPL tag to include Build Trilinos with ParMETIS TPL thirdPartyLibs#196

klevzoff · 2023-05-09T22:53:48Z

src/coreComponents/linearAlgebra/solvers/BlockPreconditioner.hpp

+  std::array< Matrix const *, 2 > m_prolongators{};

  /// Prolongation operators for each sub-block
-  std::array< Matrix, 2 > m_prolongators;
+  std::array< Matrix, 2 > m_prolongatorsOwned{};

  /// Matrix blocks
  BlockOperator< Vector, Matrix > m_matBlocks;

-  /// Individual block preconditioners
-  std::array< std::unique_ptr< PreconditionerBase< LAI > >, 2 > m_solvers;
+  /// Individual block preconditioner pointers
+  std::array< PreconditionerBase< LAI > *, 2 > m_solvers{};
+
+  /// Individual block preconditioner operators (when owned)
+  std::array< std::unique_ptr< PreconditionerBase< LAI > >, 2 > m_solversOwned{};


These changes add support for constructing a BlockPreconditioner out of sub-block solvers that exist elsewhere, which is occasionally useful (e.g. coupled multiscale smoother).

klevzoff · 2023-05-09T22:56:11Z

src/coreComponents/linearAlgebra/solvers/KrylovSolver.cpp

+    constexpr char const * headFormat = "{:>4}   {:>12}   {:>9}   {:>12}";
+    constexpr char const * lineFormat = "{:>4}   {:>12.6e}   {:>9.6f}   {:>12.6e}";
+    integer const iter = m_result.numIterations;
+    if( iter == 0 )
+    {
+      GEOS_LOG_RANK_0( GEOS_FMT( "[{}] start iteration", methodName() ) );
+      GEOS_LOG_RANK_0( GEOS_FMT( headFormat, "iter", "resid.norm", "conv.rate", "rel.res.norm" ) );
+    }
+    real64 const norm = m_residualNorms[iter];
+    real64 const relNorm = m_residualNorms[0] > 0.0 ? norm / m_residualNorms[0] : 0.0;
+    real64 const convRate = iter > 0 ? norm / m_residualNorms[iter - 1] : 1.0;
+    GEOS_LOG_RANK_0( GEOS_FMT( lineFormat, iter, norm, convRate, relNorm ) );


Iteration progress output improved to provide more information and match that of hypre

src/coreComponents/linearAlgebra/unitTests/testKrylovSolvers.cpp

src/coreComponents/linearAlgebra/CMakeLists.txt

klevzoff · 2023-05-09T23:41:31Z

src/coreComponents/physicsSolvers/LinearSolverParameters.cpp

+template< typename CHILD, typename T >
+void registerInputBlock( Group * parent, char const * const key, T & data )
+{
+  parent->registerGroup( key, std::make_unique< CHILD >( key, parent, data ) ).setInputFlags( InputFlags::OPTIONAL );
+}
+
+class BlockParametersInput final : public dataRepository::Group
+{
+public:
+
+  /// Constructor
+  BlockParametersInput( string const & name,
+                        Group * const parent,
+                        LinearSolverParameters::Block & params );
+
+  virtual Group * createChild( string const & childKey, string const & childName ) override
+  {
+    GEOS_UNUSED_VAR( childKey, childName );
+    return nullptr;
+  }
+
+  /// Keys appearing in XML
+  struct viewKeyStruct
+  {
+    static constexpr char const * shapeString()  { return "shape"; }
+    static constexpr char const * schurTypeString() { return "schurType"; }
+    static constexpr char const * scalingString()  { return "scaling"; }
+  };
+
+private:
+
+  LinearSolverParameters::Block & m_parameters;
+};
+
+BlockParametersInput::BlockParametersInput( string const & name,
+                                            Group * const parent,
+                                            LinearSolverParameters::Block & params )
+  :
+  Group( name, parent ),
+  m_parameters( params )
+{
+  registerWrapper( viewKeyStruct::shapeString(), &m_parameters.shape ).
+    setDefaultValue( m_parameters.shape ).
+    setInputFlag( InputFlags::OPTIONAL ).
+    setDescription( "Block preconditioner shape, options: "
+                    "``" + EnumStrings< LinearSolverParameters::Block::Shape >::concat( "``, ``" ) + "``" );
+
+  registerWrapper( viewKeyStruct::schurTypeString(), &m_parameters.schurType ).
+    setDefaultValue( m_parameters.schurType ).
+    setInputFlag( InputFlags::OPTIONAL ).
+    setDescription( "Schur complement type, options: "
+                    "``" + EnumStrings< LinearSolverParameters::Block::SchurType >::concat( "``, ``" ) + "``" );
+
+  registerWrapper( viewKeyStruct::scalingString(), &m_parameters.scaling ).
+    setDefaultValue( m_parameters.scaling ).
+    setInputFlag( InputFlags::OPTIONAL ).
+    setDescription( "Block scaling type, options: "
+                    "``" + EnumStrings< LinearSolverParameters::Block::Scaling >::concat( "``, ``" ) + "``" );
+}


This is a new way of structuring linear solver parameter inputs: instead of an ever-growing list of prefixed attributes like

<LinearSolverParameters blockShape="..." blockSchurType="..." blockScaling="..."/>

we will start adding nested structure as follows:

<LinearSolverParameters> <Block shape="..." schurType="..." scaling="..."/> </LinearSolverParameters>

This improves readability by indenting parameter blocks according to their logical structure and avoiding repetition of the prefixes. I use this extensively in multiscale work, but will also propose a similar refactoring for AMG parameters.

This is very nice, specially if we want to expose some MGR options per level...

klevzoff · 2023-05-09T23:47:49Z

src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseBase.cpp

-    FieldIdentifiers fieldsToBeSync;
-    fieldsToBeSync.addElementFields( { fields::flow::pressure::key() },
-                                     regionNames );
-
-    CommunicationTools::getInstance().synchronizeFields( fieldsToBeSync, mesh, domain.getNeighbors(), false );
-


This synchronization was useless because consistency of values is already guaranteed by virtue of ghosted cells being included in sets to which initial conditions are applied.

how about other flow solvers?

klevzoff · 2023-05-09T23:49:01Z

src/coreComponents/physicsSolvers/fluidFlow/SinglePhaseFVM.cpp

+  switch( linParams.preconditionerType )
+  {
+    default:
+    {
+      return LAInterface::createPreconditioner( linParams );
+    }


The switch statement may seem pointless at the moment, but it's left as a placeholder for adding custom preconditioners (e.g. multiscale)

klevzoff · 2023-05-09T23:51:53Z

src/coreComponents/physicsSolvers/solidMechanics/SolidMechanicsLagrangianFEM.hpp

+  /// Rigid body modes; TODO remove mutable hack
+  mutable array1d< ParallelVector > m_rigidBodyModes;


This is now mutable because it's populated on demand in createPreconditioner() function which is const. This is isn't great, but I don't have a better solution for now.

src/coreComponents/physicsSolvers/solidMechanics/SolidMechanicsLagrangianFEM.cpp

klevzoff added the flag: ready for review label May 9, 2023

klevzoff requested a review from castelletto1 May 9, 2023 09:26

klevzoff self-assigned this May 9, 2023

klevzoff commented May 9, 2023

View reviewed changes

klevzoff added the flag: requires rebaseline Requires rebaseline branch in integratedTests label May 10, 2023

klevzoff force-pushed the feature/klevzoff/linear-algebra-updates branch 2 times, most recently from 6aa101b to a167fdc Compare May 12, 2023 06:12

klevzoff force-pushed the feature/klevzoff/linear-algebra-updates branch 2 times, most recently from dfbc927 to 742c689 Compare July 18, 2023 09:09

klevzoff force-pushed the feature/klevzoff/linear-algebra-updates branch from 742c689 to 9004d43 Compare October 8, 2023 01:45

Linear algebra changes

2bb59cd

klevzoff force-pushed the feature/klevzoff/linear-algebra-updates branch from 9004d43 to 2bb59cd Compare January 7, 2024 07:58

rrsettgast requested review from rrsettgast, CusiniM, jhuang2601, frankfeifan, ryar9534, matteofrigo5, Guotong-Ren, wrtobin, untereiner, MelReyCG and corbett5 as code owners December 3, 2024 20:15

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Linear algebra updates #2427

Linear algebra updates #2427

klevzoff commented May 9, 2023 •

edited

Loading

klevzoff May 9, 2023

klevzoff May 9, 2023

klevzoff May 9, 2023 •

edited

Loading

victorapm May 10, 2023

klevzoff May 9, 2023

paveltomin May 10, 2023

klevzoff May 9, 2023

klevzoff May 9, 2023

		/// Rigid body modes; TODO remove mutable hack
		mutable array1d< ParallelVector > m_rigidBodyModes;

Linear algebra updates #2427

Are you sure you want to change the base?

Linear algebra updates #2427

Conversation

klevzoff commented May 9, 2023 • edited Loading

klevzoff May 9, 2023

Choose a reason for hiding this comment

klevzoff May 9, 2023

Choose a reason for hiding this comment

klevzoff May 9, 2023 • edited Loading

Choose a reason for hiding this comment

victorapm May 10, 2023

Choose a reason for hiding this comment

klevzoff May 9, 2023

Choose a reason for hiding this comment

paveltomin May 10, 2023

Choose a reason for hiding this comment

klevzoff May 9, 2023

Choose a reason for hiding this comment

klevzoff May 9, 2023

Choose a reason for hiding this comment

klevzoff commented May 9, 2023 •

edited

Loading

klevzoff May 9, 2023 •

edited

Loading