unmhr.hpp

Go to the documentation of this file.

 
//
// Copyright (c) 2021-2023, University of Colorado Denver. All rights reserved.
//
// This file is part of <T>LAPACK.
// <T>LAPACK is free software: you can redistribute it and/or modify it under
// the terms of the BSD 3-Clause license. See the accompanying LICENSE file.
 
#ifndef TLAPACK_UNMHR_HH
#define TLAPACK_UNMHR_HH
 
#include "tlapack/base/utils.hpp"
#include "tlapack/lapack/ung2r.hpp"
#include "tlapack/lapack/unm2r.hpp"
 
namespace tlapack {
 
template <class T, TLAPACK_SMATRIX matrix_t, TLAPACK_SVECTOR vector_t>
constexpr WorkInfo unmhr_worksize(Side side,
                                  Op trans,
                                  size_type<matrix_t> ilo,
                                  size_type<matrix_t> ihi,
                                  const matrix_t& A,
                                  const vector_t& tau,
                                  const matrix_t& C)
{
    using idx_t = size_type<matrix_t>;
    using range = pair<idx_t, idx_t>;
 
    auto&& A_s = slice(A, range{ilo + 1, ihi}, range{ilo, ihi - 1});
    auto&& tau_s = slice(tau, range{ilo, ihi - 1});
    auto&& C_s = (side == Side::Left)
                     ? slice(C, range{ilo + 1, ihi}, range{0, ncols(C)})
                     : slice(C, range{0, nrows(C)}, range{ilo + 1, ihi});
 
    return unmq_level2_worksize<T>(side, trans, FORWARD, COLUMNWISE_STORAGE,
                                   A_s, tau_s, C_s);
}
 
template <TLAPACK_SMATRIX matrix_t,
          TLAPACK_SVECTOR vector_t,
          TLAPACK_WORKSPACE work_t>
int unmhr_work(Side side,
               Op trans,
               size_type<matrix_t> ilo,
               size_type<matrix_t> ihi,
               const matrix_t& A,
               const vector_t& tau,
               matrix_t& C,
               work_t& work)
{
    using idx_t = size_type<matrix_t>;
    using range = pair<idx_t, idx_t>;
 
    auto A_s = slice(A, range{ilo + 1, ihi}, range{ilo, ihi - 1});
    auto tau_s = slice(tau, range{ilo, ihi - 1});
    auto C_s = (side == Side::Left)
                   ? slice(C, range{ilo + 1, ihi}, range{0, ncols(C)})
                   : slice(C, range{0, nrows(C)}, range{ilo + 1, ihi});
 
    return unmq_level2_work(side, trans, FORWARD, COLUMNWISE_STORAGE, A_s,
                            tau_s, C_s, work);
}
 
template <TLAPACK_SMATRIX matrix_t, TLAPACK_SVECTOR vector_t>
int unmhr(Side side,
          Op trans,
          size_type<matrix_t> ilo,
          size_type<matrix_t> ihi,
          const matrix_t& A,
          const vector_t& tau,
          matrix_t& C)
{
    using T = type_t<matrix_t>;
 
    // functor
    Create<matrix_t> new_matrix;
 
    // Allocates workspace
    WorkInfo workinfo = unmhr_worksize<T>(side, trans, ilo, ihi, A, tau, C);
    std::vector<T> work_;
    auto work = new_matrix(work_, workinfo.m, workinfo.n);
 
    return unmhr_work(side, trans, ilo, ihi, A, tau, C, work);
}
}  // namespace tlapack
 
#endif  // TLAPACK_UNMHR_HH