solid_preconditioners.h

Go to the documentation of this file.

// LIC// ====================================================================
// LIC// This file forms part of oomph-lib, the object-oriented,
// LIC// multi-physics finite-element library, available
// LIC// at http://www.oomph-lib.org.
// LIC//
// LIC// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
// LIC//
// LIC// This library is free software; you can redistribute it and/or
// LIC// modify it under the terms of the GNU Lesser General Public
// LIC// License as published by the Free Software Foundation; either
// LIC// version 2.1 of the License, or (at your option) any later version.
// LIC//
// LIC// This library is distributed in the hope that it will be useful,
// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// LIC// Lesser General Public License for more details.
// LIC//
// LIC// You should have received a copy of the GNU Lesser General Public
// LIC// License along with this library; if not, write to the Free Software
// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
// LIC// 02110-1301  USA.
// LIC//
// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
// LIC//
// LIC//====================================================================
#ifndef OOMPH_SOLID_PRECONDITIONERS_HEADER
#define OOMPH_SOLID_PRECONDITIONERS_HEADER
 
 
// Config header generated by autoconfig
#ifdef HAVE_CONFIG_H
#include <oomph-lib-config.h>
#endif
 
 
// oomphlib headers
#include "../generic/matrices.h"
#include "../generic/assembly_handler.h"
#include "../generic/problem.h"
#include "../generic/block_preconditioner.h"
#include "../generic/preconditioner.h"
#include "../generic/SuperLU_preconditioner.h"
#include "../generic/matrix_vector_product.h"
 
 
namespace oomph
{
  //===========================================================================
  //===========================================================================
  class PressureBasedSolidLSCPreconditioner
    : public BlockPreconditioner<CRDoubleMatrix>
  {
  public:
    PressureBasedSolidLSCPreconditioner()
      : BlockPreconditioner<CRDoubleMatrix>()
    {
      // Flag to indicate that the preconditioner has been setup
      // previously -- if setup() is called again, data can
      // be wiped.
      Preconditioner_has_been_setup = false;
 
      // By default we use SuperLU for both p and f blocks
      Using_default_p_preconditioner = true;
      Using_default_f_preconditioner = true;
 
      // resize the mesh pt
      // note: meaningless if subsidiary preconditioner
      this->set_nmesh(1);
      Solid_mesh_pt = 0;
 
      // Set default preconditioners (inexact solvers) -- they are
      // members of this class!
      P_preconditioner_pt = 0;
      F_preconditioner_pt = 0;
 
      // Flag to determine if mass matrix diagonal Q^{-1}
      // is used for scaling.
      P_matrix_using_scaling = true;
 
      // set Doc_time to false
      Doc_time = false;
 
      // null the off diagonal Block matrix pt
      Bt_mat_vec_pt = 0;
 
      // null the F matrix vector product helper
      F_mat_vec_pt = 0;
 
      // null the QBt matrix vector product pt
      QBt_mat_vec_pt = 0;
 
      // null the E matrix vector product helper
      E_mat_vec_pt = 0;
 
      // by default we do not form the E matrix (BQFQBt)
      Form_BFBt_product = false;
    }
 
    ~PressureBasedSolidLSCPreconditioner()
    {
      clean_up_memory();
    }
 
    PressureBasedSolidLSCPreconditioner(
      const PressureBasedSolidLSCPreconditioner&) = delete;
 
    // Commented out broken assignment operator because this can lead to a
    // conflict warning when used in the virtual inheritence hierarchy.
    // Essentially the compiler doesn't realise that two separate
    // implementations of the broken function are the same and so, quite
    // rightly, it shouts.
    /*void operator=(const PressureBasedSolidLSCPreconditioner&) = delete;*/
 
    void setup();
 
    void preconditioner_solve(const DoubleVector& r, DoubleVector& z);
 
    void set_solid_mesh(Mesh* mesh_pt)
    {
      Solid_mesh_pt = mesh_pt;
    }
 
    void enable_p_matrix_scaling()
    {
      P_matrix_using_scaling = true;
    }
 
    void disable_p_matrix_scaling()
    {
      P_matrix_using_scaling = false;
    }
 
    bool is_p_matrix_using_scaling() const
    {
      return P_matrix_using_scaling;
    }
 
    void set_p_preconditioner(Preconditioner* new_p_preconditioner_pt)
    {
      // If the default preconditioner has been used
      // clean it up now...
      if (Using_default_p_preconditioner)
      {
        delete P_preconditioner_pt;
      }
      P_preconditioner_pt = new_p_preconditioner_pt;
      Using_default_p_preconditioner = false;
    }
 
    void set_p_superlu_preconditioner()
    {
      if (!Using_default_p_preconditioner)
      {
        P_preconditioner_pt = new SuperLUPreconditioner;
        Using_default_p_preconditioner = true;
      }
    }
 
    void set_f_preconditioner(Preconditioner* new_f_preconditioner_pt)
    {
      // If the default preconditioner has been used
      // clean it up now...
      if (Using_default_f_preconditioner)
      {
        delete F_preconditioner_pt;
      }
      F_preconditioner_pt = new_f_preconditioner_pt;
      Using_default_f_preconditioner = false;
    }
 
    void set_f_superlu_preconditioner()
    {
      if (!Using_default_f_preconditioner)
      {
        F_preconditioner_pt = new SuperLUPreconditioner;
        Using_default_f_preconditioner = true;
      }
    }
 
 
    void enable_doc_time()
    {
      Doc_time = true;
    }
 
    void disable_doc_time()
    {
      Doc_time = false;
    }
 
    void enable_form_BFBt_product()
    {
      Form_BFBt_product = true;
    }
 
    void disable_form_BFBt_product()
    {
      Form_BFBt_product = false;
    }
 
    void clean_up_memory();
 
  private:
    // oomph-lib objects
    // -----------------
 
    // Pointers to preconditioner (=inexact solver) objects
    // -----------------------------------------------------
    Preconditioner* P_preconditioner_pt;
 
    Preconditioner* F_preconditioner_pt;
 
    bool Using_default_f_preconditioner;
 
    bool Using_default_p_preconditioner;
 
    bool Preconditioner_has_been_setup;
 
    bool P_matrix_using_scaling;
 
    CRDoubleMatrix* assemble_mass_matrix_diagonal();
 
    bool F_preconditioner_is_block_preconditioner;
 
    bool Doc_time;
 
    MatrixVectorProduct* F_mat_vec_pt;
 
    MatrixVectorProduct* QBt_mat_vec_pt;
 
    MatrixVectorProduct* Bt_mat_vec_pt;
 
    MatrixVectorProduct* E_mat_vec_pt;
 
    bool Form_BFBt_product;
 
    Mesh* Solid_mesh_pt;
  };
 
 
 
 
  //============================================================================
  //=============================================================================
  template<typename MATRIX>
  class PressureBasedSolidExactPreconditioner
    : public BlockPreconditioner<MATRIX>
  {
  public:
    PressureBasedSolidExactPreconditioner() : BlockPreconditioner<MATRIX>() {}
 
 
    ~PressureBasedSolidExactPreconditioner() {}
 
 
    PressureBasedSolidExactPreconditioner(
      const PressureBasedSolidExactPreconditioner&) = delete;
 
 
    /*void operator=(const PressureBasedSolidExactPreconditioner&) = delete;*/
 
 
    void setup();
 
    void preconditioner_solve(const Vector<double>& r, Vector<double>& z);
 
  protected:
    MATRIX P_matrix;
  };
 
} // namespace oomph
#endif