dhels Subroutine

pure subroutine dhels(a, lda, n, q, b)
!! This routine solves the least squares problem:
!!
!!  $$ \min{\| b - A x \|^2} $$
!!
!! using the factors computed by [[dheqr]]. This is similar to the LINPACK routine [[dgesl]]
!! except that \(A\) is an upper Hessenberg matrix.

   use daskr_kinds, only: rk
   use blas_interfaces, only: daxpy
   implicit none

   real(rk), intent(in) :: a(lda, n)
      !! Output from [[dheqr]] which contains the upper triangular factor R in the QR
      !! decomposition of `a`.
   integer, intent(in) :: lda
      !! Leading dimension of `a`.
   integer, intent(in) :: n
      !! Number of columns in `a` (originally a matrix with shape `(n+1, n)`).
   real(rk), intent(in) :: q(2*n)
      !! Coefficients of the Givens rotations used in decomposing `a`.
   real(rk), intent(inout) :: b(n + 1)
      !! On entry, the right hand side vector. On return, the solution vector \(x\).

   integer :: iq, k, kb, kp1
   real(rk) :: c, s, t, t1, t2

   ! Minimize (B-A*X,B-A*X).
   ! First form Q*B.
   do k = 1, n
      kp1 = k + 1
      iq = 2*(k - 1) + 1
      c = q(iq)
      s = q(iq + 1)
      t1 = b(k)
      t2 = b(kp1)
      b(k) = c*t1 - s*t2
      b(kp1) = s*t1 + c*t2
   end do

   ! Now solve R*X = Q*B.
   do kb = 1, n
      k = n + 1 - kb
      b(k) = b(k)/a(k, k)
      t = -b(k)
      call daxpy(k - 1, t, a(1, k), 1, b(1), 1)
   end do

end subroutine dhels