schwarz_precond.hpp

// FEAT3: Finite Element Analysis Toolbox, Version 3
// Copyright (C) 2010 by Stefan Turek & the FEAT group
// FEAT3 is released under the GNU General Public License version 3,
// see the file 'copyright.txt' in the top level directory for details.
 
#pragma once
 
#include <kernel/solver/base.hpp>
#include <kernel/global/vector.hpp>
#include <kernel/global/filter.hpp>
 
namespace FEAT
{
  namespace Solver
  {
    template<typename Vector_, typename Filter_>
    class SchwarzPrecond;
 
    template<typename LocalVector_, typename LocalFilter_, typename Mirror_>
    class SchwarzPrecond<Global::Vector<LocalVector_, Mirror_>, Global::Filter<LocalFilter_, Mirror_>> :
      public SolverBase<Global::Vector<LocalVector_, Mirror_>>
    {
    public:
      typedef LocalVector_ LocalVectorType;
      typedef Global::Vector<LocalVector_, Mirror_> GlobalVectorType;
      typedef Global::Filter<LocalFilter_, Mirror_> GlobalFilterType;
      typedef SolverBase<Global::Vector<LocalVector_, Mirror_>> BaseClass;
 
      typedef SolverBase<LocalVector_> LocalSolverType;
 
    protected:
      std::shared_ptr<LocalSolverType> _local_solver;
      const GlobalFilterType& _filter;
      bool _ignore_status;
 
    public:
      explicit SchwarzPrecond(std::shared_ptr<LocalSolverType> local_solver, const GlobalFilterType& filter, bool ignore_status = false) :
        _local_solver(local_solver),
        _filter(filter),
        _ignore_status(ignore_status)
      {
      }
 
      explicit SchwarzPrecond(const String& section_name, const PropertyMap* section,
        std::shared_ptr<LocalSolverType> local_solver, const GlobalFilterType& filter) :
        BaseClass(section_name, section),
        _local_solver(local_solver),
        _filter(filter),
        _ignore_status(false)
      {
        auto ignore_status_p = section->query("ignore_status");
        if(ignore_status_p.second)
        {
          if(ignore_status_p.first.is_one_of("yes true"))
            this->_ignore_status = true;
          else if(ignore_status_p.first.is_one_of("no false"))
            this->_ignore_status = false;
          else
            throw ParseError(section_name + ".ignore_status", ignore_status_p.first, "yes/no or true/false");
        }
      }
 
      virtual String name() const override
      {
        return "Schwarz";
      }
 
      void set_ignore_status(bool ignore_status)
      {
        this->_ignore_status = ignore_status;
      }
 
      virtual void init_symbolic() override
      {
        _local_solver->init_symbolic();
      }
 
      virtual void init_numeric() override
      {
        _local_solver->init_numeric();
      }
 
      virtual void done_numeric() override
      {
        _local_solver->done_numeric();
      }
 
      virtual void done_symbolic() override
      {
        _local_solver->done_symbolic();
      }
 
      virtual Status apply(GlobalVectorType& vec_cor, const GlobalVectorType& vec_def) override
      {
        Statistics::add_solver_expression(std::make_shared<ExpressionStartSolve>(this->name()));
 
        // apply local solver
        Statistics::add_solver_expression(std::make_shared<ExpressionCallPrecond>(this->name(), this->_local_solver->name()));
        Status status = _local_solver->apply(vec_cor.local(), vec_def.local());
 
        // ignore or synchronize status codes?
        if(this->_ignore_status)
          status = Status::success;
        else
        {
          // synchronize local status over communicator to obtain
          // a consistent status code on all processes
          const Dist::Comm* comm = vec_cor.get_comm();
          if(comm != nullptr)
          {
            // local status: 0: success, 1: failure
            int lstat = (status_success(status) ? 0 : 1);
            int gstat = -1;
 
            // sum up over all processes:
            comm->allreduce(&lstat, &gstat, std::size_t(1), Dist::op_sum);
            XASSERT(gstat >= 0);
 
            // if all processes succeeded, the sum will also be 0
            // if the sum is > 0, then at least one process failed,
            // so the global status will also be failure
            status = (gstat == 0 ? Status::success : Status::aborted);
          }
        }
 
        // did all processes succeed?
        if(status == Status::success)
        {
          // synchronize correction vector
          vec_cor.sync_1();
 
          // apply filter
          _filter.filter_cor(vec_cor);
        }
 
        // okay
        Statistics::add_solver_expression(std::make_shared<ExpressionEndSolve>(this->name(), status, 0));
        return status;
      }
    }; // class SchwarzPrecond<...>
 
    template<typename LocalFilter_, typename Mirror_>
    inline std::shared_ptr<SchwarzPrecond<Global::Vector<typename LocalFilter_::VectorType, Mirror_>, Global::Filter<LocalFilter_, Mirror_>>> new_schwarz_precond(
      std::shared_ptr<SolverBase<typename LocalFilter_::VectorType>> local_solver,
      const Global::Filter<LocalFilter_, Mirror_>& filter)
    {
      return std::make_shared<SchwarzPrecond<Global::Vector<typename LocalFilter_::VectorType, Mirror_>, Global::Filter<LocalFilter_, Mirror_>>>
        (local_solver, filter);
    }
 
    template<typename LocalFilter_, typename Mirror_>
    inline std::shared_ptr
    <
      SchwarzPrecond
      <
        Global::Vector<typename LocalFilter_::VectorType, Mirror_>,
        Global::Filter<LocalFilter_, Mirror_>
      >
    >
    new_schwarz_precond(
      const String& section_name, const PropertyMap* section,
      std::shared_ptr<SolverBase<typename LocalFilter_::VectorType>> local_solver,
      const Global::Filter<LocalFilter_, Mirror_>& filter)
    {
      return std::make_shared<SchwarzPrecond<Global::Vector<typename LocalFilter_::VectorType, Mirror_>, Global::Filter<LocalFilter_, Mirror_>>>
        (section_name, section, local_solver, filter);
    }
  } // namespace Solver
} // namespace FEAT