core/matvec/Matrix.h

// -*- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature -*-
//-----------------------------------------------------------------------------
// Copyright 2000-2026 CEA (www.cea.fr) IFPEN (www.ifpenergiesnouvelles.com)
// See the top-level COPYRIGHT file for details.
// SPDX-License-Identifier: Apache-2.0
//-----------------------------------------------------------------------------
/*---------------------------------------------------------------------------*/
/* Matrix.h                                                    (C) 2000-20265 */
/*                                                                           */
/* Matrix d'algèbre linéraire.                                               */
/*---------------------------------------------------------------------------*/
#ifndef ARCANE_CORE_MATVEC_MATRIX_H
#define ARCANE_CORE_MATVEC_MATRIX_H
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
#include "arcane/utils/Array.h"
#include "arcane/utils/Numeric.h"
 
#include "arcane/core/matvec/Vector.h"
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
namespace Arcane::MatVec
{
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
class Vector;
class MatrixImpl;
class IPreconditioner;
class AMG;
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
class ARCANE_CORE_EXPORT Matrix
{
 public:
 
  Matrix() = default;
  Matrix(Integer nb_row, Integer nb_column);
  Matrix(const Matrix& rhs);
  ~Matrix();
  void operator=(const Matrix& rhs);
 
 private:
 
  explicit Matrix(MatrixImpl* impl);
 
 public:

  Matrix clone() const;
 
 public:
 
  // Nombre de lignes de la matrice
  Integer nbRow() const;
  // Nombre de colonnes de la matrice
  Integer nbColumn() const;
  void setRowsSize(IntegerConstArrayView rows_size);
  void setValues(IntegerConstArrayView columns, RealConstArrayView values);
  void dump(std::ostream& o) const;
  RealConstArrayView values() const;
  RealArrayView values();
  IntegerConstArrayView rowsIndex() const;
  IntegerConstArrayView columns() const;
  IntegerArrayView rowsIndex();
  IntegerArrayView columns();
  void setValue(Integer row, Integer column, Real value);
  Real value(Integer row, Integer column) const;
  void sortDiagonale();
  static Matrix read(const String& filename);
  static Matrix readHypre(const String& filename);
 
 private:

  MatrixImpl* m_impl = nullptr;
};
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
class ARCANE_CORE_EXPORT IPreconditioner
{
 public:
 
  virtual ~IPreconditioner() {}
 
 public:
 
  virtual void apply(Vector& out_vec, const Vector& vec) = 0;
};
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
class ARCANE_CORE_EXPORT ConjugateGradientSolver
{
 public:
 
  ConjugateGradientSolver() = default;
  bool solve(const Matrix& a, const Vector& b, Vector& x, Real epsilon, IPreconditioner* p = 0);
  Integer nbIteration() const { return m_nb_iteration; }
  Real residualNorm() const { return m_residual_norm; }
  void setMaxIteration(Integer max_iteration)
  {
    m_max_iteration = max_iteration;
  }
 
 private:
 
  Integer m_nb_iteration = 0;
  Real m_residual_norm = 0.0;
  Integer m_max_iteration = 5000;
 
 private:
 
  void _applySolver(const Matrix& a, const Vector& b, Vector& x, Real epsilon, IPreconditioner* p);
  void _applySolver2(const Matrix& a, const Vector& b, Vector& x, Real epsilon, IPreconditioner* precond);
  void _applySolverAsHypre(const Matrix& a, const Vector& b, Vector& x, Real tol, IPreconditioner* precond);
};
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
class ARCANE_CORE_EXPORT DiagonalPreconditioner
: public IPreconditioner
{
 public:
 
  explicit DiagonalPreconditioner(const Matrix& matrix);
 
 public:
 
  virtual void apply(Vector& out_vec, const Vector& vec);
 
 private:
 
  Vector m_inverse_diagonal;
};
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
class ARCANE_CORE_EXPORT MatrixOperation
{
 public:
 
  void matrixVectorProduct(const Matrix& mat, const Vector& vec, Vector& out_vec);
  void matrixVectorProduct2(const Matrix& mat, const Vector& vec, Vector& out_vec);
  Real dot(const Vector& vec);
  Real dot(const Vector& vec1, const Vector& vec2);
  void negateVector(Vector& vec);
  void scaleVector(Vector& vec, Real mul);
  void addVector(Vector& out_vec, const Vector& vec);
};
 
class ARCANE_CORE_EXPORT MatrixOperation2
{
 public:
 
  Matrix matrixMatrixProduct(const Matrix& left_matrix, const Matrix& right_matrix);
  Matrix matrixMatrixProductFast(const Matrix& left_matrix, const Matrix& right_matrix);
  Matrix transpose(const Matrix& matrix);
  Matrix transposeFast(const Matrix& matrix);
  // Applique le produit de matrice L * M * R avec L = transpose(R)
  Matrix applyGalerkinOperator(const Matrix& left_matrix, const Matrix& matrix, const Matrix& right_matrix);
  Matrix applyGalerkinOperator2(const Matrix& left_matrix, const Matrix& matrix,
                                const Matrix& right_matrix);
 
 private:
 
  void _dumpColumnMatrix(std::ostream& o, IntegerConstArrayView columns_index,
                         IntegerConstArrayView rows,
                         RealConstArrayView values);
};
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
class ARCANE_CORE_EXPORT AMGPreconditioner
: public IPreconditioner
{
 public:
 
  explicit AMGPreconditioner(ITraceMng* tm)
  : m_trace_mng(tm)
  {}
  ~AMGPreconditioner() override;
 
 public:
 
  virtual void build(const Matrix& matrix);
 
 public:
 
  void apply(Vector& out_vec, const Vector& vec) override;
 
 private:
 
  ITraceMng* m_trace_mng = nullptr;
  AMG* m_amg = nullptr;
};
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
class ARCANE_CORE_EXPORT AMGSolver
{
 public:
 
  AMGSolver(ITraceMng* tm)
  : m_trace_mng(tm)
  , m_amg(0)
  {}
  virtual ~AMGSolver();
 
 public:
 
  virtual void build(const Matrix& matrix);
 
 public:
 
  virtual void solve(const Vector& vector_b, Vector& vector_x);
 
 private:
 
  ITraceMng* m_trace_mng = nullptr;
  AMG* m_amg = nullptr;
};
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
class ARCANE_CORE_EXPORT DirectSolver
{
 public:
 
  void solve(const Matrix& matrix, const Vector& vector_b, Vector& vector_x);
 
 private:
 
  void _solve(RealArrayView mat_values, RealArrayView vec_values, Integer size);
};
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
} // namespace Arcane::MatVec
 
/*---------------------------------------------------------------------------*/
/*---------------------------------------------------------------------------*/
 
#endif