partialsolvelt.h

Go to the documentation of this file.

// ========================================================================= //
// This file is part of MyraMath, copyright (c) 2014-2019 by Ryan A Chilton  //
// and distributed by MyraCore, LLC. See LICENSE.txt for license terms.      //
// ========================================================================= //

#ifndef MYRAMATH_MULTIFRONTAL_LU_PARTIALSOLVELT_H
#define MYRAMATH_MULTIFRONTAL_LU_PARTIALSOLVELT_H

#include <myramath/utility/Number.h>

#include <myramath/dense/MatrixRange.h>
#include <myramath/dense/gemm.h>
#include <myramath/dense/detail/nwork.h>

#include <myramath/multifrontal/detail/partialsolveu.h>
#include <myramath/multifrontal/detail/PContainer.h>
#include <myramath/multifrontal/detail/lu/LUContainer.h>

namespace myra {
namespace multifrontal {
namespace lu {
namespace partialsolvelt {

template<class Kernel> class JobGraphBase : public ::myra::multifrontal::detail::partialsolveu::JobGraphBase2<typename ReflectNumber<Kernel>::type>
  {
  public:

    // Structural typedefs for base class.
    typedef typename ReflectNumber<Kernel>::type Number;
    typedef ::myra::multifrontal::detail::partialsolveu::JobGraphBase2<Number> Base;

    // Typedefs related to numeric storage.
    typedef ::myra::multifrontal::detail::lu::LUContainer<Kernel> LUContainer;
    typedef ::myra::multifrontal::detail::PContainer<Number> PContainer;
    typedef CMatrixRange<Number> CRange;
    typedef  MatrixRange<Number>  Range;
    
    // Constructor - requires const reference to LContainer, mutable reference to PContainer, and op (either 'T' or 'H')
    JobGraphBase(const LUContainer* in_lucontainer, const SchurTree* in_stree, PContainer* in_pcontainer, char in_op)
      : Base(in_stree,in_pcontainer), lucontainer(in_lucontainer), op(in_op) { }

    // Virtual copy constructor.
    virtual ::myra::JobGraphBase* clone() const
      { return new JobGraphBase(*this); }

  protected:

    // Accessors for L Kernel's / Block's
    const Kernel& l(int n, int ij) const
      { return lucontainer->lu(this->s2a(n),ij); }
    const CRange l(int n, int i, int j) const
      { return lucontainer->lu(this->s2a(n),i,j); }
  
    // Solves Ln(k,k)'*Xn(k,j) = Bn(k,j)
    virtual uint64_t backsolve(int n, int k, int j)
      {
      const Kernel& Ln_kk = this->l(n,k);
      Range Bn_kj = this->b(n,k,j);
      return Ln_kk.solveL(Bn_kj, 'L', op);
      }

    // Downdates Bn(i,j) -= Ln(k,i)'*Xn(k,j)
    virtual uint64_t downdate(int n, int k, int i, int j) 
      {
      Range Bn_ij = this->b(n,i,j);
      CRange Ln_ki = this->l(n,k,i);
      CRange Xn_kj = this->x(n,k,j);
      Number one(1);
      return gemm_nwork(Bn_ij, Ln_ki, op, Xn_kj, 'N', -one, one);
      }

  private:
  
    // Reference to LUContainer (const).
    const LUContainer* lucontainer;

    // Op modifier, either 'T' or 'H'
    char op;

  }; // class JobGraph

} } } } // namespace

#endif