This routine solves the least squares problem:
using the factors computed by dheqr. This is similar to the LINPACK routine dgesl except that is an upper Hessenberg matrix.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real(kind=rk), | intent(in) | :: | a(lda,n) |
Output from dheqr which contains the upper triangular factor R in the QR
decomposition of |
||
integer, | intent(in) | :: | lda |
Leading dimension of |
||
integer, | intent(in) | :: | n |
Number of columns in |
||
real(kind=rk), | intent(in) | :: | q(2*n) |
Coefficients of the Givens rotations used in decomposing |
||
real(kind=rk), | intent(inout) | :: | b(n+1) |
On entry, the right hand side vector. On return, the solution vector . |
Type | Visibility | Attributes | Name | Initial | |||
---|---|---|---|---|---|---|---|
integer, | public | :: | iq | ||||
integer, | public | :: | k | ||||
integer, | public | :: | kb | ||||
integer, | public | :: | kp1 | ||||
real(kind=rk), | public | :: | c | ||||
real(kind=rk), | public | :: | s | ||||
real(kind=rk), | public | :: | t | ||||
real(kind=rk), | public | :: | t1 | ||||
real(kind=rk), | public | :: | t2 |
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