df_transforming.h

Go to the documentation of this file.

/* 
 * SPDX-FileCopyrightText: 2025 Gesellschaft fuer Anlagen- und Reaktorsicherheit gGmbH
 * SPDX-License-Identifier: EUPL-1.2
 * SPDX-FileContributor: Arne Naegel
 * SPDX-FileContributor: Goethe Universität Frankfurt
 * SPDX-FileType: SOURCE
 * 
 * This file is part of d3f++.
 * d3f++ is an extension for UG4. Licensing information and citation requirements of UG4 are provided in LICENSES/UG4-LGPL_2.1
 */
 
#ifndef __H__UG__LIB_DISC__OPERATOR__LINEAR_OPERATOR__TRANSFORMING__
#define __H__UG__LIB_DISC__OPERATOR__LINEAR_OPERATOR__TRANSFORMING__
 
#include "lib_algebra/operator/interface/preconditioner.h"
#include "lib_algebra/algebra_common/sparsematrix_util.h"
#include "lib_algebra/operator/preconditioner/preconditioners.h"
#include "lib_algebra/operator/debug_writer.h"
 
#include "df_transform.h"
#ifdef UG_PARALLEL
  #include "lib_algebra/parallelization/parallelization.h"
  #include "lib_algebra/parallelization/parallel_matrix_overlap_impl.h"
#endif
 
namespace ug{
 
 
template <typename TAlgebra, typename TTransform>
class LeftTransformIterator :
//public virtual ILinearIterator<typename TAlgebra::vector_type>,
//public DebugWritingObject <TAlgebra>
public ITransformingIteration<TAlgebra, LeftTransformIterator<TAlgebra, TTransform> >
{
  public:
    typedef TAlgebra algebra_type;
 
    typedef typename TAlgebra::vector_type vector_type;
 
    typedef typename TAlgebra::matrix_type matrix_type;
 
    typedef ILinearIterator<typename TAlgebra::vector_type> base_type;
 
    typedef TTransform transform_type;
 
  public:
  //  Constructor
    LeftTransformIterator() :  m_transform()
 
    {
      MatrixOperator<matrix_type, vector_type> *ptr =new MatrixOperator<matrix_type, vector_type>();
       m_spOpTransformed = SmartPtr<MatrixOperator<matrix_type, vector_type> >(ptr);
    };
 
  //  Clone
    virtual SmartPtr<ILinearIterator<vector_type> > clone()
    {
      SmartPtr<LeftTransformIterator<algebra_type, transform_type> > newInst(new LeftTransformIterator<algebra_type, transform_type>());
      newInst->set_damp(this->damping());
      return newInst;
    }
 
private:
      bool init(SmartPtr<MatrixOperator<matrix_type, vector_type> > pOp)
      {
        // transform
        get_transform().transform_operator(pOp->get_matrix(), transformed_matrix()->get_matrix());
 
        // write debug
        this->write_debug(pOp->get_matrix(), "DebugOriginalMatrix");
        this->write_debug(transformed_matrix()->get_matrix(), "DebugTransformedMatrix");
 
        // forward init
        return get_preconditioner()->init(transformed_matrix());
      }
 
public:
 
 
        bool init(SmartPtr<ILinearOperator<vector_type> > J, const vector_type& u)
        {
          m_spOriginalOp = J;
 
          //  cast to matrix based operator
          SmartPtr<MatrixOperator<matrix_type, vector_type> > pOp =
              J.template cast_dynamic<MatrixOperator<matrix_type, vector_type> >();
 
          //  Check that matrix if of correct type
          if(pOp.invalid())
            UG_THROW(name() << "::init': Passed Operator is "
                "not based on matrix. This Preconditioner can only "
                "handle matrix-based operators.");
 
          //  forward request to matrix based implementation
          get_transform().init(pOp->get_matrix(), &u);
          return init(pOp);
        }
 
 
        bool init(SmartPtr<ILinearOperator<vector_type> > L)
        {
          m_spOriginalOp = L;
          //  cast to matrix based operator
          SmartPtr<MatrixOperator<matrix_type, vector_type> > pOp =
            L.template cast_dynamic<MatrixOperator<matrix_type, vector_type> >();
 
          //  Check that matrix if of correct type
            if(pOp.invalid())
              UG_THROW(name() << "::init': Passed Operator is "
              "not based on matrix. This Preconditioner can only "
              "handle matrix-based operators.");
 
          //  forward request to matrix based implementation
          get_transform().init(pOp->get_matrix(), NULL);
          return init(pOp);
        };
 
 
        bool transform_defect(vector_type &c, const vector_type& d)
        { return get_transform().transform_rhs(d, m_dTransformed); }
 
        bool apply_transformed(vector_type &c, const vector_type &d)
        { return get_preconditioner()->apply(c, m_dTransformed); }
 
        bool untransform_correction(vector_type& c, const vector_type& d)
        { return true; }
 
        SmartPtr<ILinearOperator<vector_type> > original_operator()
        {return m_spOriginalOp;}
 
 
        void set_preconditioner(SmartPtr<base_type> pprecond) {m_pprecond=pprecond;}
 
 
 
  protected:
  // reference to 'real' preconditioner
 
    SmartPtr<base_type> get_preconditioner() {return m_pprecond;}
    transform_type& get_transform() {return m_transform;}
 
    SmartPtr <MatrixOperator<matrix_type, vector_type> > transformed_matrix()
    {return m_spOpTransformed;}
 
    vector_type& transformed_rhs()
    {return m_dTransformed; }
 
    virtual bool supports_parallel() const
    {
      UG_LOG("WARNING: Not sure if LeftTransformIterator supports parallel execution!\n");
      return true;
    }
 
  //  Name of preconditioner
    virtual const char* name() const {return "LeftTransform";}
 
 
 
 
  protected:
    transform_type m_transform;
    SmartPtr<base_type> m_pprecond;
 
    vector_type m_dTransformed;
    SmartPtr<MatrixOperator<matrix_type, vector_type> > m_spOpTransformed;
 
    SmartPtr<ILinearOperator<vector_type> > m_spOriginalOp;
 
};
 
 
 
} // end namespace ug
 
#endif // __H__UG__LIB_DISC__OPERATOR__LINEAR_OPERATOR__TRANSFORMING__