Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #pragma once
      11             : 
      12             : #include "libmesh/libmesh_config.h"
      13             : 
      14             : #include "NonlinearSystemBase.h"
      15             : #include "SlepcEigenSolverConfiguration.h"
      16             : 
      17             : #include "libmesh/condensed_eigen_system.h"
      18             : 
      19             : // forward declarations
      20             : class EigenProblem;
      21             : class KernelBase;
      22             : class ResidualObject;
      23             : 
      24             : #ifdef LIBMESH_HAVE_SLEPC
      25             : 
      26             : #include <slepceps.h>
      27             : 
      28             : namespace Moose
      29             : {
      30             : void assemble_matrix(EquationSystems & es, const std::string & system_name);
      31             : }
      32             : 
      33             : /**
      34             :  * Nonlinear eigenvalue system to be solved
      35             :  */
      36             : class NonlinearEigenSystem : public NonlinearSystemBase
      37             : {
      38             : public:
      39             :   NonlinearEigenSystem(EigenProblem & problem, const std::string & name);
      40             : 
      41             :   virtual void solve() override;
      42             : 
      43             :   /**
      44             :    * Quit the current solve as soon as possible.
      45             :    */
      46             :   virtual void stopSolve(const ExecFlagType & exec_flag,
      47             :                          const std::set<TagID> & vector_tags_to_close) override;
      48             : 
      49             :   /**
      50             :    * Returns the current nonlinear iteration number.  In libmesh, this is
      51             :    * updated during the nonlinear solve, so it should be up-to-date.
      52             :    */
      53             :   virtual unsigned int getCurrentNonlinearIterationNumber() override;
      54             : 
      55             :   virtual void setupFiniteDifferencedPreconditioner() override;
      56             : 
      57             :   /**
      58             :    * Returns the convergence state
      59             :    *
      60             :    * @return true if converged, otherwise false
      61             :    */
      62             :   virtual bool converged() override;
      63             : 
      64             :   virtual NumericVector<Number> & RHS() override;
      65             : 
      66             :   /*
      67             :    *  A residual vector at MOOSE side for AX
      68             :    */
      69             :   NumericVector<Number> & residualVectorAX();
      70             : 
      71             :   /*
      72             :    * A residual vector at MOOSE side for BX
      73             :    */
      74             :   NumericVector<Number> & residualVectorBX();
      75             : 
      76             :   void attachSLEPcCallbacks();
      77             : 
      78             :   /**
      79             :    * Get the number of converged eigenvalues
      80             :    *
      81             :    * @return The number of converged eigenvalues
      82             :    */
      83        4878 :   unsigned int getNumConvergedEigenvalues() const { return _eigen_sys.get_n_converged(); };
      84             : 
      85             :   virtual libMesh::NonlinearSolver<Number> * nonlinearSolver() override;
      86             : 
      87             :   /**
      88             :    * Retrieve snes from slepc eigen solver. It is valid for only nonlinear eigen solver.
      89             :    * You should see a big error if you do this for linear solver.
      90             :    */
      91             :   virtual SNES getSNES() override;
      92             : 
      93             :   /**
      94             :    * Retrieve EPS (SLEPc eigen solver)
      95             :    */
      96             :   virtual EPS getEPS();
      97             : 
      98       42532 :   libMesh::CondensedEigenSystem & sys() { return _eigen_sys; }
      99             : 
     100             :   /**
     101             :    * For eigenvalue problems (including standard and generalized), inhomogeneous (Dirichlet or
     102             :    * Neumann)
     103             :    * boundary conditions are not allowed.
     104             :    */
     105             :   void checkIntegrity();
     106             : 
     107             :   /**
     108             :    * Return the Nth converged eigenvalue.
     109             :    *
     110             :    * @return The Nth converged eigenvalue as a complex number, i.e. the first and the second number
     111             :    * is the real and the imaginary part of
     112             :    * the eigenvalue, respectively.
     113             :    */
     114             :   std::pair<Real, Real> getConvergedEigenvalue(dof_id_type n) const;
     115             : 
     116             :   /**
     117             :    * Return the Nth converged eigenvalue and copies the respective eigen vector to the solution
     118             :    * vector.
     119             :    *
     120             :    * @return The Nth converged eigenvalue as a complex number, i.e. the first and the second number
     121             :    * is the real and the imaginary part of
     122             :    * the eigenvalue, respectively.
     123             :    */
     124             :   std::pair<Real, Real> getConvergedEigenpair(dof_id_type n) const;
     125             : 
     126             :   /**
     127             :    * Get the number of converged eigenvalues
     128             :    *
     129             :    * @return all converged eigenvalues as complex numbers
     130             :    */
     131         673 :   const std::vector<std::pair<Real, Real>> & getAllConvergedEigenvalues() const
     132             :   {
     133         673 :     return _eigen_values;
     134             :   }
     135             : 
     136             :   /**
     137             :    * Vector tag ID of right hand side
     138             :    */
     139       47081 :   TagID eigenVectorTag() const { return _Bx_tag; }
     140             : 
     141             :   /**
     142             :    * Vector tag ID of left hand side
     143             :    */
     144       56936 :   TagID nonEigenVectorTag() const { return _Ax_tag; }
     145             : 
     146             :   /**
     147             :    * Matrix tag ID of right hand side
     148             :    */
     149        4225 :   TagID eigenMatrixTag() const { return _B_tag; }
     150             : 
     151             :   /**
     152             :    * Matrix tag ID of left hand side
     153             :    */
     154        4258 :   TagID nonEigenMatrixTag() const { return _A_tag; }
     155             : 
     156             :   std::set<TagID> defaultVectorTags() const override;
     157             :   std::set<TagID> defaultMatrixTags() const override;
     158             : 
     159             :   /**
     160             :    * If the preconditioning matrix includes eigen kernels
     161             :    */
     162         558 :   void precondMatrixIncludesEigenKernels(bool precond_matrix_includes_eigen)
     163             :   {
     164         558 :     _precond_matrix_includes_eigen = precond_matrix_includes_eigen;
     165         558 :   }
     166             : 
     167             :   bool precondMatrixIncludesEigenKernels() const { return _precond_matrix_includes_eigen; }
     168             : 
     169        3808 :   TagID precondMatrixTag() const { return _precond_tag; }
     170             : 
     171             :   virtual void attachPreconditioner(libMesh::Preconditioner<Number> * preconditioner) override;
     172             : 
     173        1590 :   libMesh::Preconditioner<Number> * preconditioner() const { return _preconditioner; }
     174             : 
     175             :   virtual void turnOffJacobian() override;
     176             : 
     177             :   void residualAndJacobianTogether() override;
     178             : 
     179             :   /**
     180             :    * Initialize the condensed matrices. This is a no-op if there are no
     181             :    * constraints in the DofMap
     182             :    */
     183             :   void initializeCondensedMatrices();
     184             : 
     185             :   virtual void postInit() override;
     186             :   virtual void reinit() override;
     187             : 
     188             : protected:
     189             :   virtual void postAddResidualObject(ResidualObject & object) override;
     190             : 
     191             :   void computeScalingJacobian() override;
     192             :   void computeScalingResidual() override;
     193             : 
     194             :   libMesh::CondensedEigenSystem & _eigen_sys;
     195             :   EigenProblem & _eigen_problem;
     196             :   std::unique_ptr<SlepcEigenSolverConfiguration> _solver_configuration;
     197             :   std::vector<std::pair<Real, Real>> _eigen_values;
     198             :   unsigned int _n_eigen_pairs_required;
     199             :   NumericVector<Number> & _work_rhs_vector_AX;
     200             :   NumericVector<Number> & _work_rhs_vector_BX;
     201             :   TagID _Ax_tag;
     202             :   TagID _Bx_tag;
     203             :   TagID _A_tag;
     204             :   TagID _B_tag;
     205             :   TagID _precond_tag;
     206             :   bool _precond_matrix_includes_eigen;
     207             :   // Libmesh preconditioner
     208             :   libMesh::Preconditioner<Number> * _preconditioner;
     209             : 
     210             :   /// The number of degrees of freedom constrained at the libMesh level, e.g. via hanging node or
     211             :   /// periodic boundary constraints
     212             :   dof_id_type _num_constrained_dofs;
     213             : };
     214             : 
     215             : #else
     216             : 
     217             : class NonlinearEigenSystem : public libMesh::ParallelObject
     218             : {
     219             : public:
     220             :   NonlinearEigenSystem(EigenProblem & problem, const std::string & name);
     221             : 
     222             :   /**
     223             :    * Returns the convergence state
     224             :    * @return true if converged, otherwise false
     225             :    */
     226             :   bool converged() { return false; }
     227             : 
     228             :   void checkIntegrity() {}
     229             : };
     230             : 
     231             : #endif