codac2_linear_ctc.h

Go to the documentation of this file.

 
#pragma once
 
#include <string>
#include "codac2_Collection.h"
#include "codac2_template_tools.h"
#include "codac2_Matrix.h"
#include "codac2_IntervalVector.h"
#include "codac2_IntervalMatrix.h"
 
namespace codac2
{
  class CtcLinearBase
  {
    public:
 
      virtual void contract(IntervalMatrix& A, IntervalVector& x, IntervalVector& b) const = 0;
      virtual std::shared_ptr<CtcLinearBase> copy() const = 0;
      virtual ~CtcLinearBase() = default;
  };

  class CtcGaussElim : public CtcLinearBase
  {
    public:

      CtcGaussElim()
      { }

      void contract(IntervalMatrix& A, IntervalVector& x, IntervalVector& b) const;
 
      virtual std::shared_ptr<CtcLinearBase> copy() const
      {
        return std::make_shared<CtcGaussElim>(*this);
      }
  };
 

  class CtcGaussSeidel : public CtcLinearBase
  {
    public:

      CtcGaussSeidel()
      { }

      void contract(IntervalMatrix& A, IntervalVector& x, IntervalVector& b) const;
 
      virtual std::shared_ptr<CtcLinearBase> copy() const
      {
        return std::make_shared<CtcGaussSeidel>(*this);
      }
  };
 

  class CtcLinearPrecond : public CtcLinearBase
  {
    public:

      template<typename C>
        requires (std::is_base_of_v<CtcLinearBase,C> || std::is_same_v<std::shared_ptr<CtcLinearBase>,C>)
      CtcLinearPrecond(const C& ctc_no_precond)
        : _ctc_no_precond(ctc_no_precond)
      { }

      void contract(IntervalMatrix& A, IntervalVector& x, IntervalVector& b) const
      {
        assert_release(A.is_squared() && A.rows() == x.size() && A.rows() == b.size());
 
        IntervalMatrix A0 = A.mid();
        IntervalMatrix A0_inv = A.mid().inverse().template cast<Interval>();
        IntervalMatrix Ap = A0_inv*A;
        IntervalVector bp = A0_inv*b;
 
        _ctc_no_precond.front()->contract(Ap,x,bp);
 
        b &= A0*bp;
        A &= A0*Ap;
      }
 
      virtual std::shared_ptr<CtcLinearBase> copy() const
      {
        return std::make_shared<CtcLinearPrecond>(*this);
      }
 
 
    protected:
 
      const Collection<CtcLinearBase> _ctc_no_precond;
  };
}