MatrixTranslator.h

Go to the documentation of this file.

#pragma once

#include <memory>

#include "TimeDiscretization.h"
#include "Types.h"

#include "MathLib/LinAlg/LinAlg.h"

namespace NumLib
{

template <ODESystemTag ODETag>
class MatrixTranslator;

template <>
class MatrixTranslator<ODESystemTag::FirstOrderImplicitQuasilinear>
{
public:
    virtual void computeA(GlobalMatrix const& M, GlobalMatrix const& K,
                          GlobalMatrix& A) const = 0;

    virtual void computeRhs(const GlobalMatrix& M, const GlobalMatrix& K,
                            const GlobalVector& b, const GlobalVector& x_prev,
                            GlobalVector& rhs) const = 0;

    virtual void computeResidual(GlobalMatrix const& M, GlobalMatrix const& K,
                                 GlobalVector const& b,
                                 GlobalVector const& x_curr,
                                 GlobalVector const& xdot,
                                 GlobalVector& res) const = 0;

    virtual void computeJacobian(GlobalMatrix const& Jac_in,
                                 GlobalMatrix& Jac_out) const = 0;

    virtual ~MatrixTranslator() = default;
};

template <ODESystemTag ODETag>
class MatrixTranslatorGeneral;

template <>
class MatrixTranslatorGeneral<ODESystemTag::FirstOrderImplicitQuasilinear>
    : public MatrixTranslator<ODESystemTag::FirstOrderImplicitQuasilinear>
{
public:
    explicit MatrixTranslatorGeneral(TimeDiscretization const& timeDisc)
        : _time_disc(timeDisc)
    {
    }

    void computeA(GlobalMatrix const& M, GlobalMatrix const& K,
                  GlobalMatrix& A) const override;

    void computeRhs(const GlobalMatrix& M, const GlobalMatrix& /*K*/,
                    const GlobalVector& b, const GlobalVector& x_prev,
                    GlobalVector& rhs) const override;

    void computeResidual(GlobalMatrix const& M, GlobalMatrix const& K,
                         GlobalVector const& b,
                         GlobalVector const& x_curr,
                         GlobalVector const& xdot,
                         GlobalVector& res) const override;

    void computeJacobian(GlobalMatrix const& Jac_in,
                         GlobalMatrix& Jac_out) const override;

private:
    TimeDiscretization const&
        _time_disc;  

    mutable std::size_t _tmp_id = 0u;
};

template <ODESystemTag ODETag>
std::unique_ptr<MatrixTranslator<ODETag>> createMatrixTranslator(
    TimeDiscretization const& timeDisc)
{
    return std::unique_ptr<MatrixTranslator<ODETag>>(
        new MatrixTranslatorGeneral<ODETag>(timeDisc));
}

}  // namespace NumLib