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 "FEProblem.h"
      13             : #include "libmesh/libmesh_config.h"
      14             : #include <petscsnes.h>
      15             : 
      16             : class NonlinearSystem;
      17             : 
      18             : /**
      19             :  * A problem that handles Schur complement preconditioning of the incompressible Navier-Stokes
      20             :  * equations
      21             :  */
      22             : class NavierStokesProblem : public FEProblem
      23             : {
      24             : public:
      25             :   static InputParameters validParams();
      26             : 
      27             :   NavierStokesProblem(const InputParameters & parameters);
      28             : 
      29             : #if PETSC_RELEASE_GREATER_EQUALS(3, 20, 0)
      30             :   /**
      31             :    * @returns the mass matrix tag ID
      32             :    */
      33        1906 :   TagID massMatrixTagID() const { return getMatrixTagID(_mass_matrix); }
      34             : 
      35             :   /**
      36             :    * @returns the poisson operator matrix tag ID
      37             :    */
      38          42 :   TagID LMatrixTagID() const { return getMatrixTagID(_L_matrix); }
      39             : 
      40             :   /**
      41             :    * Clear the field split index sets
      42             :    */
      43        1906 :   void clearIndexSets() { _index_sets.clear(); }
      44             : 
      45             :   /*
      46             :    * Given a \p KSP \p node and where we are in the field split tree, given by \p tree_position,
      47             :    * return the next \p KSP object in the tree on the way to the Schur complement \p KSP object.
      48             :    * Each invocation of this method moves through one level of our \p _index_sets data member. This
      49             :    * method will call itself recursively until it reaches the Schur complement \p KSP
      50             :    */
      51             :   KSP findSchurKSP(KSP node, unsigned int tree_position);
      52             : 
      53             :   /**
      54             :    * Setup the Least Squares Commutator (LSC) preconditioner given the Schur complement \p KSP
      55             :    * object
      56             :    */
      57             :   void setupLSCMatrices(KSP schur_ksp);
      58             : 
      59             :   /**
      60             :    * Will destroy any matrices we allocated
      61             :    */
      62             :   virtual ~NavierStokesProblem();
      63             : 
      64             : protected:
      65             :   /**
      66             :    * Reinitialize PETSc output for proper linear/nonlinear iteration display
      67             :    */
      68             :   virtual void initPetscOutputAndSomeSolverSettings() override;
      69             : 
      70             : private:
      71             :   /// Whether to commute operators in the style of Olshanskii. If this is true, then the user must
      72             :   /// provide both (pressure) mass matrices and a Poisson operator for the velocity
      73             :   const bool _commute_lsc;
      74             :   /// The tag name of the mass matrix
      75             :   const TagName & _mass_matrix;
      76             :   /// The tag name of the Poisson operator
      77             :   const TagName & _L_matrix;
      78             :   /// Whether the user attached a mass matrix
      79             :   const bool _have_mass_matrix;
      80             :   /// Whether the user attached a Poisson operator matrix
      81             :   const bool _have_L_matrix;
      82             : 
      83             :   /// Whether to directly use the pressure mass matrix to form the Schur complement
      84             :   /// preconditioner. This is only appropriate for Stokes flow in which the pressure mass matrix is
      85             :   /// spectrally equivalent to the Schur complement
      86             :   const bool _pressure_mass_matrix_as_pre;
      87             : 
      88             :   /// The length of this vector should correspond to the number of split nesting levels there are in
      89             :   /// the field split. Then the integers should indicate the path one shold take in the nesting tree
      90             :   /// to get to the location of the Schur complement field split
      91             :   const std::vector<unsigned int> & _schur_fs_index;
      92             : 
      93             :   /// The mass matrix used for scaling
      94             :   Mat _Q_scale = nullptr;
      95             :   /// The Poisson operator
      96             :   Mat _L = nullptr;
      97             : 
      98             :   /// This will end up being the same length as \p _schur_fs_index. Let's give an example of what
      99             :   /// this data member means. If the user sets "schur_fs_index = '1'", then this means the Schur
     100             :   /// complement field split is nested within another field split, and the Schur complement field
     101             :   /// split is at the 1st index of the top split (some other set of degrees of freedom take up the
     102             :   /// 0th index of the top split). So in this example \p _index_sets will be of length 1, and the
     103             :   /// Index Set (IS) held by this container will hold all the Schur complement field split degrees
     104             :   /// of freedom (e.g. all the system degrees of freedom minus the degrees of freedom held in the
     105             :   /// 0th index of the top split). An example of this example is if we split out all the velocity
     106             :   /// Dirichlet degrees of freedom into the 0th index of the top split, and then our Schur
     107             :   /// complement at index 1 of the top split handles all non-Dirichlet velocity degrees of freedom
     108             :   /// and all pressure degrees of freedom
     109             :   std::vector<IS> _index_sets;
     110             : #endif
     111             : };