Line data    Source code

       1              : // Copyright Amazon.com, Inc. or its affiliates. All Rights Reserved.
       2              : // SPDX-License-Identifier: Apache-2.0
       3              : 
       4              : #ifndef PALACE_LINALG_DIST_RELAXATION_SMOOTHER_HPP
       5              : #define PALACE_LINALG_DIST_RELAXATION_SMOOTHER_HPP
       6              : 
       7              : #include <memory>
       8              : #include "linalg/operator.hpp"
       9              : #include "linalg/solver.hpp"
      10              : #include "linalg/vector.hpp"
      11              : 
      12              : namespace mfem
      13              : {
      14              : 
      15              : template <typename T>
      16              : class Array;
      17              : 
      18              : }  // namespace mfem
      19              : 
      20              : namespace palace
      21              : {
      22              : 
      23              : //
      24              : // Hiptmair distributive relaxation smoother applying smoothers to both the operator in the
      25              : // primary space as well as its projection into an auxiliary space.
      26              : // Reference: Hiptmair, Multigrid method for Maxwell's equations, SIAM J. Numer. Anal.
      27              : //            (1998).
      28              : //
      29              : template <typename OperType>
      30              : class DistRelaxationSmoother : public Solver<OperType>
      31              : {
      32              :   using VecType = typename Solver<OperType>::VecType;
      33              : 
      34              : private:
      35              :   // Number of smoother iterations.
      36              :   const int pc_it;
      37              : 
      38              :   // Discrete gradient matrix (not owned).
      39              :   const Operator *G;
      40              : 
      41              :   // System matrix and its projection GᵀAG (not owned).
      42              :   const OperType *A, *A_G;
      43              :   const mfem::Array<int> *dbc_tdof_list_G;
      44              : 
      45              :   // Point smoother objects for each matrix.
      46              :   mutable std::unique_ptr<Solver<OperType>> B;
      47              :   std::unique_ptr<Solver<OperType>> B_G;
      48              : 
      49              :   // Temporary vectors for smoother application.
      50              :   mutable VecType x_G, y_G, r_G, r;
      51              : 
      52              : public:
      53              :   DistRelaxationSmoother(MPI_Comm comm, const Operator &G, int smooth_it,
      54              :                          int cheby_smooth_it, int cheby_order, double cheby_sf_max,
      55              :                          double cheby_sf_min, bool cheby_4th_kind);
      56              : 
      57            0 :   void SetOperator(const OperType &op) override
      58              :   {
      59            0 :     MFEM_ABORT("SetOperator with a single operator is not implemented for "
      60              :                "DistRelaxationSmoother, use the two argument signature instead!");
      61              :   }
      62              : 
      63              :   void SetOperators(const OperType &op, const OperType &op_G);
      64              : 
      65            0 :   void Mult(const VecType &x, VecType &y) const override
      66              :   {
      67            0 :     if (r.Size() != y.Size())
      68              :     {
      69            0 :       r.SetSize(y.Size());
      70            0 :       r.UseDevice(true);
      71              :     }
      72            0 :     Mult2(x, y, r);
      73            0 :   }
      74              : 
      75            0 :   void MultTranspose(const VecType &x, VecType &y) const override
      76              :   {
      77            0 :     if (r.Size() != y.Size())
      78              :     {
      79            0 :       r.SetSize(y.Size());
      80            0 :       r.UseDevice(true);
      81              :     }
      82            0 :     MultTranspose2(x, y, r);
      83            0 :   }
      84              : 
      85              :   void Mult2(const VecType &x, VecType &y, VecType &r) const override;
      86              : 
      87              :   void MultTranspose2(const VecType &x, VecType &y, VecType &r) const override;
      88              : };
      89              : 
      90              : }  // namespace palace
      91              : 
      92              : #endif  // PALACE_LINALG_DIST_RELAXATION_SMOOTHER_HPP