#:include "common.fypp"
submodule(stdlib_blas) stdlib_blas_level2_tri
implicit none
contains
#:for ik,it,ii in LINALG_INT_KINDS_TYPES
pure module subroutine stdlib${ii}$_strmv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_sp
!! STRMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(sp), intent(in) :: a(lda,*)
real(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('STRMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
do i = 1,j - 1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do i = 1,j - 1
x(ix) = x(ix) + temp*a(i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
do i = n,j + 1,-1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do i = n,j + 1,-1
x(ix) = x(ix) + temp*a(i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx - incx
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
temp = temp + a(i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
ix = ix - incx
temp = temp + a(i,j)*x(ix)
end do
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
temp = temp + a(i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
ix = ix + incx
temp = temp + a(i,j)*x(ix)
end do
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_strmv
pure module subroutine stdlib${ii}$_dtrmv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_dp
!! DTRMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(dp), intent(in) :: a(lda,*)
real(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTRMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
do i = 1,j - 1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do i = 1,j - 1
x(ix) = x(ix) + temp*a(i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
do i = n,j + 1,-1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do i = n,j + 1,-1
x(ix) = x(ix) + temp*a(i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx - incx
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
temp = temp + a(i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
ix = ix - incx
temp = temp + a(i,j)*x(ix)
end do
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
temp = temp + a(i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
ix = ix + incx
temp = temp + a(i,j)*x(ix)
end do
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_dtrmv
#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ri}$trmv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_${rk}$
!! DTRMV: performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(${rk}$), intent(in) :: a(lda,*)
real(${rk}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(${rk}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTRMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
do i = 1,j - 1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do i = 1,j - 1
x(ix) = x(ix) + temp*a(i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
do i = n,j + 1,-1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do i = n,j + 1,-1
x(ix) = x(ix) + temp*a(i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx - incx
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
temp = temp + a(i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
ix = ix - incx
temp = temp + a(i,j)*x(ix)
end do
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
temp = temp + a(i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
ix = ix + incx
temp = temp + a(i,j)*x(ix)
end do
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ri}$trmv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_ctrmv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_sp
!! CTRMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(sp), intent(in) :: a(lda,*)
complex(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('CTRMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
do i = 1,j - 1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do i = 1,j - 1
x(ix) = x(ix) + temp*a(i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
do i = n,j + 1,-1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do i = n,j + 1,-1
x(ix) = x(ix) + temp*a(i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx - incx
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
temp = temp + a(i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j - 1,1,-1
temp = temp + conjg(a(i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
ix = ix - incx
temp = temp + a(i,j)*x(ix)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j - 1,1,-1
ix = ix - incx
temp = temp + conjg(a(i,j))*x(ix)
end do
end if
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
temp = temp + a(i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j + 1,n
temp = temp + conjg(a(i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
ix = ix + incx
temp = temp + a(i,j)*x(ix)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j + 1,n
ix = ix + incx
temp = temp + conjg(a(i,j))*x(ix)
end do
end if
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ctrmv
pure module subroutine stdlib${ii}$_ztrmv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_dp
!! ZTRMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(dp), intent(in) :: a(lda,*)
complex(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTRMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
do i = 1,j - 1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do i = 1,j - 1
x(ix) = x(ix) + temp*a(i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
do i = n,j + 1,-1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do i = n,j + 1,-1
x(ix) = x(ix) + temp*a(i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx - incx
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
temp = temp + a(i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j - 1,1,-1
temp = temp + conjg(a(i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
ix = ix - incx
temp = temp + a(i,j)*x(ix)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j - 1,1,-1
ix = ix - incx
temp = temp + conjg(a(i,j))*x(ix)
end do
end if
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
temp = temp + a(i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j + 1,n
temp = temp + conjg(a(i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
ix = ix + incx
temp = temp + a(i,j)*x(ix)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j + 1,n
ix = ix + incx
temp = temp + conjg(a(i,j))*x(ix)
end do
end if
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ztrmv
#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ci}$trmv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_${ck}$
!! ZTRMV: performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(${ck}$), intent(in) :: a(lda,*)
complex(${ck}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(${ck}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTRMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
do i = 1,j - 1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do i = 1,j - 1
x(ix) = x(ix) + temp*a(i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
do i = n,j + 1,-1
x(i) = x(i) + temp*a(i,j)
end do
if (nounit) x(j) = x(j)*a(j,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do i = n,j + 1,-1
x(ix) = x(ix) + temp*a(i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(j,j)
end if
jx = jx - incx
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
temp = temp + a(i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j - 1,1,-1
temp = temp + conjg(a(i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j - 1,1,-1
ix = ix - incx
temp = temp + a(i,j)*x(ix)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j - 1,1,-1
ix = ix - incx
temp = temp + conjg(a(i,j))*x(ix)
end do
end if
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
temp = temp + a(i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j + 1,n
temp = temp + conjg(a(i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*a(j,j)
do i = j + 1,n
ix = ix + incx
temp = temp + a(i,j)*x(ix)
end do
else
if (nounit) temp = temp*conjg(a(j,j))
do i = j + 1,n
ix = ix + incx
temp = temp + conjg(a(i,j))*x(ix)
end do
end if
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ci}$trmv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_stbmv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_sp
!! STBMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular band matrix, with ( k + 1 ) diagonals.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(sp), intent(in) :: a(lda,*)
real(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('STBMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(kplus1,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(kplus1,j)
end if
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(1,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(1,j)
end if
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = kplus1 - j
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(i)
end do
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
kx = kx - incx
ix = kx
l = kplus1 - j
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(ix)
ix = ix - incx
end do
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
l = 1 - j
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
kx = kx + incx
ix = kx
l = 1 - j
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(ix)
ix = ix + incx
end do
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_stbmv
pure module subroutine stdlib${ii}$_dtbmv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_dp
!! DTBMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular band matrix, with ( k + 1 ) diagonals.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(dp), intent(in) :: a(lda,*)
real(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTBMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(kplus1,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(kplus1,j)
end if
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(1,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(1,j)
end if
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = kplus1 - j
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(i)
end do
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
kx = kx - incx
ix = kx
l = kplus1 - j
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(ix)
ix = ix - incx
end do
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
l = 1 - j
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
kx = kx + incx
ix = kx
l = 1 - j
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(ix)
ix = ix + incx
end do
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_dtbmv
#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ri}$tbmv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_${rk}$
!! DTBMV: performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular band matrix, with ( k + 1 ) diagonals.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(${rk}$), intent(in) :: a(lda,*)
real(${rk}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(${rk}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTBMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(kplus1,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(kplus1,j)
end if
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(1,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(1,j)
end if
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = kplus1 - j
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(i)
end do
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
kx = kx - incx
ix = kx
l = kplus1 - j
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(ix)
ix = ix - incx
end do
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
l = 1 - j
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(i)
end do
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
kx = kx + incx
ix = kx
l = 1 - j
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(ix)
ix = ix + incx
end do
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ri}$tbmv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_ctbmv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_sp
!! CTBMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular band matrix, with ( k + 1 ) diagonals.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(sp), intent(in) :: a(lda,*)
complex(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('CTBMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(kplus1,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(kplus1,j)
end if
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(1,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(1,j)
end if
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = kplus1 - j
if (noconj) then
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(kplus1,j))
do i = j - 1,max(1,j-k),-1
temp = temp + conjg(a(l+i,j))*x(i)
end do
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
kx = kx - incx
ix = kx
l = kplus1 - j
if (noconj) then
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(ix)
ix = ix - incx
end do
else
if (nounit) temp = temp*conjg(a(kplus1,j))
do i = j - 1,max(1,j-k),-1
temp = temp + conjg(a(l+i,j))*x(ix)
ix = ix - incx
end do
end if
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
l = 1 - j
if (noconj) then
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(1,j))
do i = j + 1,min(n,j+k)
temp = temp + conjg(a(l+i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
kx = kx + incx
ix = kx
l = 1 - j
if (noconj) then
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(ix)
ix = ix + incx
end do
else
if (nounit) temp = temp*conjg(a(1,j))
do i = j + 1,min(n,j+k)
temp = temp + conjg(a(l+i,j))*x(ix)
ix = ix + incx
end do
end if
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ctbmv
pure module subroutine stdlib${ii}$_ztbmv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_dp
!! ZTBMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular band matrix, with ( k + 1 ) diagonals.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(dp), intent(in) :: a(lda,*)
complex(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTBMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(kplus1,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(kplus1,j)
end if
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(1,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(1,j)
end if
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = kplus1 - j
if (noconj) then
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(kplus1,j))
do i = j - 1,max(1,j-k),-1
temp = temp + conjg(a(l+i,j))*x(i)
end do
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
kx = kx - incx
ix = kx
l = kplus1 - j
if (noconj) then
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(ix)
ix = ix - incx
end do
else
if (nounit) temp = temp*conjg(a(kplus1,j))
do i = j - 1,max(1,j-k),-1
temp = temp + conjg(a(l+i,j))*x(ix)
ix = ix - incx
end do
end if
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
l = 1 - j
if (noconj) then
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(1,j))
do i = j + 1,min(n,j+k)
temp = temp + conjg(a(l+i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
kx = kx + incx
ix = kx
l = 1 - j
if (noconj) then
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(ix)
ix = ix + incx
end do
else
if (nounit) temp = temp*conjg(a(1,j))
do i = j + 1,min(n,j+k)
temp = temp + conjg(a(l+i,j))*x(ix)
ix = ix + incx
end do
end if
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ztbmv
#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ci}$tbmv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_${ck}$
!! ZTBMV: performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular band matrix, with ( k + 1 ) diagonals.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(${ck}$), intent(in) :: a(lda,*)
complex(${ck}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(${ck}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTBMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := a*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(kplus1,j)
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*a(kplus1,j)
end if
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(i) = x(i) + temp*a(l+i,j)
end do
if (nounit) x(j) = x(j)*a(1,j)
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
x(ix) = x(ix) + temp*a(l+i,j)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*a(1,j)
end if
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = kplus1 - j
if (noconj) then
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(kplus1,j))
do i = j - 1,max(1,j-k),-1
temp = temp + conjg(a(l+i,j))*x(i)
end do
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
kx = kx - incx
ix = kx
l = kplus1 - j
if (noconj) then
if (nounit) temp = temp*a(kplus1,j)
do i = j - 1,max(1,j-k),-1
temp = temp + a(l+i,j)*x(ix)
ix = ix - incx
end do
else
if (nounit) temp = temp*conjg(a(kplus1,j))
do i = j - 1,max(1,j-k),-1
temp = temp + conjg(a(l+i,j))*x(ix)
ix = ix - incx
end do
end if
x(jx) = temp
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
temp = x(j)
l = 1 - j
if (noconj) then
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(i)
end do
else
if (nounit) temp = temp*conjg(a(1,j))
do i = j + 1,min(n,j+k)
temp = temp + conjg(a(l+i,j))*x(i)
end do
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
kx = kx + incx
ix = kx
l = 1 - j
if (noconj) then
if (nounit) temp = temp*a(1,j)
do i = j + 1,min(n,j+k)
temp = temp + a(l+i,j)*x(ix)
ix = ix + incx
end do
else
if (nounit) temp = temp*conjg(a(1,j))
do i = j + 1,min(n,j+k)
temp = temp + conjg(a(l+i,j))*x(ix)
ix = ix + incx
end do
end if
x(jx) = temp
jx = jx + incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ci}$tbmv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_stpmv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_sp
!! STPMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix, supplied in packed form.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(sp), intent(in) :: ap(*)
real(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: nounit
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('STPMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with one pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x:= a*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
k = kk
do i = 1,j - 1
x(i) = x(i) + temp*ap(k)
k = k + 1
end do
if (nounit) x(j) = x(j)*ap(kk+j-1)
end if
kk = kk + j
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
x(ix) = x(ix) + temp*ap(k)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*ap(kk+j-1)
end if
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
k = kk
do i = n,j + 1,-1
x(i) = x(i) + temp*ap(k)
k = k - 1
end do
if (nounit) x(j) = x(j)*ap(kk-n+j)
end if
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
x(ix) = x(ix) + temp*ap(k)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*ap(kk-n+j)
end if
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (nounit) temp = temp*ap(kk)
k = kk - 1
do i = j - 1,1,-1
temp = temp + ap(k)*x(i)
k = k - 1
end do
x(j) = temp
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (nounit) temp = temp*ap(kk)
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + ap(k)*x(ix)
end do
x(jx) = temp
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
if (nounit) temp = temp*ap(kk)
k = kk + 1
do i = j + 1,n
temp = temp + ap(k)*x(i)
k = k + 1
end do
x(j) = temp
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (nounit) temp = temp*ap(kk)
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + ap(k)*x(ix)
end do
x(jx) = temp
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_stpmv
pure module subroutine stdlib${ii}$_dtpmv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_dp
!! DTPMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix, supplied in packed form.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(dp), intent(in) :: ap(*)
real(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: nounit
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTPMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with one pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x:= a*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
k = kk
do i = 1,j - 1
x(i) = x(i) + temp*ap(k)
k = k + 1
end do
if (nounit) x(j) = x(j)*ap(kk+j-1)
end if
kk = kk + j
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
x(ix) = x(ix) + temp*ap(k)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*ap(kk+j-1)
end if
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
k = kk
do i = n,j + 1,-1
x(i) = x(i) + temp*ap(k)
k = k - 1
end do
if (nounit) x(j) = x(j)*ap(kk-n+j)
end if
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
x(ix) = x(ix) + temp*ap(k)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*ap(kk-n+j)
end if
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (nounit) temp = temp*ap(kk)
k = kk - 1
do i = j - 1,1,-1
temp = temp + ap(k)*x(i)
k = k - 1
end do
x(j) = temp
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (nounit) temp = temp*ap(kk)
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + ap(k)*x(ix)
end do
x(jx) = temp
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
if (nounit) temp = temp*ap(kk)
k = kk + 1
do i = j + 1,n
temp = temp + ap(k)*x(i)
k = k + 1
end do
x(j) = temp
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (nounit) temp = temp*ap(kk)
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + ap(k)*x(ix)
end do
x(jx) = temp
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_dtpmv
#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ri}$tpmv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_${rk}$
!! DTPMV: performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix, supplied in packed form.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(${rk}$), intent(in) :: ap(*)
real(${rk}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(${rk}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: nounit
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTPMV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with one pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x:= a*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
temp = x(j)
k = kk
do i = 1,j - 1
x(i) = x(i) + temp*ap(k)
k = k + 1
end do
if (nounit) x(j) = x(j)*ap(kk+j-1)
end if
kk = kk + j
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
x(ix) = x(ix) + temp*ap(k)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*ap(kk+j-1)
end if
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
temp = x(j)
k = kk
do i = n,j + 1,-1
x(i) = x(i) + temp*ap(k)
k = k - 1
end do
if (nounit) x(j) = x(j)*ap(kk-n+j)
end if
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=zero) then
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
x(ix) = x(ix) + temp*ap(k)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*ap(kk-n+j)
end if
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
else
! form x := a**t*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (nounit) temp = temp*ap(kk)
k = kk - 1
do i = j - 1,1,-1
temp = temp + ap(k)*x(i)
k = k - 1
end do
x(j) = temp
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (nounit) temp = temp*ap(kk)
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + ap(k)*x(ix)
end do
x(jx) = temp
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
if (nounit) temp = temp*ap(kk)
k = kk + 1
do i = j + 1,n
temp = temp + ap(k)*x(i)
k = k + 1
end do
x(j) = temp
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (nounit) temp = temp*ap(kk)
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + ap(k)*x(ix)
end do
x(jx) = temp
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ri}$tpmv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_ctpmv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_sp
!! CTPMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix, supplied in packed form.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(sp), intent(in) :: ap(*)
complex(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('CTPMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with cone pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x:= a*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
k = kk
do i = 1,j - 1
x(i) = x(i) + temp*ap(k)
k = k + 1
end do
if (nounit) x(j) = x(j)*ap(kk+j-1)
end if
kk = kk + j
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
x(ix) = x(ix) + temp*ap(k)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*ap(kk+j-1)
end if
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
k = kk
do i = n,j + 1,-1
x(i) = x(i) + temp*ap(k)
k = k - 1
end do
if (nounit) x(j) = x(j)*ap(kk-n+j)
end if
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
x(ix) = x(ix) + temp*ap(k)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*ap(kk-n+j)
end if
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk - 1
if (noconj) then
if (nounit) temp = temp*ap(kk)
do i = j - 1,1,-1
temp = temp + ap(k)*x(i)
k = k - 1
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do i = j - 1,1,-1
temp = temp + conjg(ap(k))*x(i)
k = k - 1
end do
end if
x(j) = temp
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*ap(kk)
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + ap(k)*x(ix)
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + conjg(ap(k))*x(ix)
end do
end if
x(jx) = temp
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk + 1
if (noconj) then
if (nounit) temp = temp*ap(kk)
do i = j + 1,n
temp = temp + ap(k)*x(i)
k = k + 1
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do i = j + 1,n
temp = temp + conjg(ap(k))*x(i)
k = k + 1
end do
end if
x(j) = temp
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*ap(kk)
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + ap(k)*x(ix)
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + conjg(ap(k))*x(ix)
end do
end if
x(jx) = temp
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ctpmv
pure module subroutine stdlib${ii}$_ztpmv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_dp
!! ZTPMV performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix, supplied in packed form.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(dp), intent(in) :: ap(*)
complex(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTPMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with cone pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x:= a*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
k = kk
do i = 1,j - 1
x(i) = x(i) + temp*ap(k)
k = k + 1
end do
if (nounit) x(j) = x(j)*ap(kk+j-1)
end if
kk = kk + j
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
x(ix) = x(ix) + temp*ap(k)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*ap(kk+j-1)
end if
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
k = kk
do i = n,j + 1,-1
x(i) = x(i) + temp*ap(k)
k = k - 1
end do
if (nounit) x(j) = x(j)*ap(kk-n+j)
end if
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
x(ix) = x(ix) + temp*ap(k)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*ap(kk-n+j)
end if
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk - 1
if (noconj) then
if (nounit) temp = temp*ap(kk)
do i = j - 1,1,-1
temp = temp + ap(k)*x(i)
k = k - 1
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do i = j - 1,1,-1
temp = temp + conjg(ap(k))*x(i)
k = k - 1
end do
end if
x(j) = temp
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*ap(kk)
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + ap(k)*x(ix)
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + conjg(ap(k))*x(ix)
end do
end if
x(jx) = temp
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk + 1
if (noconj) then
if (nounit) temp = temp*ap(kk)
do i = j + 1,n
temp = temp + ap(k)*x(i)
k = k + 1
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do i = j + 1,n
temp = temp + conjg(ap(k))*x(i)
k = k + 1
end do
end if
x(j) = temp
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*ap(kk)
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + ap(k)*x(ix)
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + conjg(ap(k))*x(ix)
end do
end if
x(jx) = temp
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ztpmv
#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ci}$tpmv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_${ck}$
!! ZTPMV: performs one of the matrix-vector operations
!! x := A*x, or x := A**T*x, or x := A**H*x,
!! where x is an n element vector and A is an n by n unit, or non-unit,
!! upper or lower triangular matrix, supplied in packed form.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(${ck}$), intent(in) :: ap(*)
complex(${ck}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(${ck}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTPMV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with cone pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x:= a*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
temp = x(j)
k = kk
do i = 1,j - 1
x(i) = x(i) + temp*ap(k)
k = k + 1
end do
if (nounit) x(j) = x(j)*ap(kk+j-1)
end if
kk = kk + j
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
x(ix) = x(ix) + temp*ap(k)
ix = ix + incx
end do
if (nounit) x(jx) = x(jx)*ap(kk+j-1)
end if
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
temp = x(j)
k = kk
do i = n,j + 1,-1
x(i) = x(i) + temp*ap(k)
k = k - 1
end do
if (nounit) x(j) = x(j)*ap(kk-n+j)
end if
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
if (x(jx)/=czero) then
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
x(ix) = x(ix) + temp*ap(k)
ix = ix - incx
end do
if (nounit) x(jx) = x(jx)*ap(kk-n+j)
end if
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
else
! form x := a**t*x or x := a**h*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk - 1
if (noconj) then
if (nounit) temp = temp*ap(kk)
do i = j - 1,1,-1
temp = temp + ap(k)*x(i)
k = k - 1
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do i = j - 1,1,-1
temp = temp + conjg(ap(k))*x(i)
k = k - 1
end do
end if
x(j) = temp
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*ap(kk)
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + ap(k)*x(ix)
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
temp = temp + conjg(ap(k))*x(ix)
end do
end if
x(jx) = temp
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk + 1
if (noconj) then
if (nounit) temp = temp*ap(kk)
do i = j + 1,n
temp = temp + ap(k)*x(i)
k = k + 1
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do i = j + 1,n
temp = temp + conjg(ap(k))*x(i)
k = k + 1
end do
end if
x(j) = temp
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = jx
if (noconj) then
if (nounit) temp = temp*ap(kk)
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + ap(k)*x(ix)
end do
else
if (nounit) temp = temp*conjg(ap(kk))
do k = kk + 1,kk + n - j
ix = ix + incx
temp = temp + conjg(ap(k))*x(ix)
end do
end if
x(jx) = temp
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ci}$tpmv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_strsv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_sp
!! STRSV solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(sp), intent(in) :: a(lda,*)
real(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('STRSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j - 1,1,-1
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j - 1,1,-1
ix = ix - incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j + 1,n
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j + 1,n
ix = ix + incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
temp = x(j)
do i = 1,j - 1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
do i = 1,j - 1
temp = temp - a(i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(j,j)
x(jx) = temp
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
do i = n,j + 1,-1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
do i = n,j + 1,-1
temp = temp - a(i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(j,j)
x(jx) = temp
jx = jx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_strsv
pure module subroutine stdlib${ii}$_dtrsv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_dp
!! DTRSV solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level1 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(dp), intent(in) :: a(lda,*)
real(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTRSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j - 1,1,-1
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j - 1,1,-1
ix = ix - incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j + 1,n
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j + 1,n
ix = ix + incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
temp = x(j)
do i = 1,j - 1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
do i = 1,j - 1
temp = temp - a(i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(j,j)
x(jx) = temp
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
do i = n,j + 1,-1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
do i = n,j + 1,-1
temp = temp - a(i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(j,j)
x(jx) = temp
jx = jx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_dtrsv
#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ri}$trsv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_${rk}$
!! DTRSV: solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level1 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(${rk}$), intent(in) :: a(lda,*)
real(${rk}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(${rk}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTRSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j - 1,1,-1
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j - 1,1,-1
ix = ix - incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j + 1,n
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j + 1,n
ix = ix + incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
temp = x(j)
do i = 1,j - 1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
do i = 1,j - 1
temp = temp - a(i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(j,j)
x(jx) = temp
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
do i = n,j + 1,-1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
do i = n,j + 1,-1
temp = temp - a(i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(j,j)
x(jx) = temp
jx = jx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ri}$trsv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_ctrsv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_sp
!! CTRSV solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(sp), intent(in) :: a(lda,*)
complex(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('CTRSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j - 1,1,-1
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j - 1,1,-1
ix = ix - incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j + 1,n
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j + 1,n
ix = ix + incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
temp = x(j)
if (noconj) then
do i = 1,j - 1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
else
do i = 1,j - 1
temp = temp - conjg(a(i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
ix = kx
temp = x(jx)
if (noconj) then
do i = 1,j - 1
temp = temp - a(i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(j,j)
else
do i = 1,j - 1
temp = temp - conjg(a(i,j))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(jx) = temp
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (noconj) then
do i = n,j + 1,-1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
else
do i = n,j + 1,-1
temp = temp - conjg(a(i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
ix = kx
temp = x(jx)
if (noconj) then
do i = n,j + 1,-1
temp = temp - a(i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(j,j)
else
do i = n,j + 1,-1
temp = temp - conjg(a(i,j))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(jx) = temp
jx = jx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ctrsv
pure module subroutine stdlib${ii}$_ztrsv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_dp
!! ZTRSV solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(dp), intent(in) :: a(lda,*)
complex(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTRSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j - 1,1,-1
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j - 1,1,-1
ix = ix - incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j + 1,n
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j + 1,n
ix = ix + incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
temp = x(j)
if (noconj) then
do i = 1,j - 1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
else
do i = 1,j - 1
temp = temp - conjg(a(i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
ix = kx
temp = x(jx)
if (noconj) then
do i = 1,j - 1
temp = temp - a(i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(j,j)
else
do i = 1,j - 1
temp = temp - conjg(a(i,j))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(jx) = temp
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (noconj) then
do i = n,j + 1,-1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
else
do i = n,j + 1,-1
temp = temp - conjg(a(i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
ix = kx
temp = x(jx)
if (noconj) then
do i = n,j + 1,-1
temp = temp - a(i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(j,j)
else
do i = n,j + 1,-1
temp = temp - conjg(a(i,j))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(jx) = temp
jx = jx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ztrsv
#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ci}$trsv(uplo,trans,diag,n,a,lda,x,incx)
use stdlib_blas_constants_${ck}$
!! ZTRSV: solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(${ck}$), intent(in) :: a(lda,*)
complex(${ck}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(${ck}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (lda<max(1,n)) then
info = 6
else if (incx==0) then
info = 8
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTRSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j - 1,1,-1
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j - 1,1,-1
ix = ix - incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/a(j,j)
temp = x(j)
do i = j + 1,n
x(i) = x(i) - temp*a(i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/a(j,j)
temp = x(jx)
ix = jx
do i = j + 1,n
ix = ix + incx
x(ix) = x(ix) - temp*a(i,j)
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
if (incx==1) then
do j = 1,n
temp = x(j)
if (noconj) then
do i = 1,j - 1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
else
do i = 1,j - 1
temp = temp - conjg(a(i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
ix = kx
temp = x(jx)
if (noconj) then
do i = 1,j - 1
temp = temp - a(i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(j,j)
else
do i = 1,j - 1
temp = temp - conjg(a(i,j))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(jx) = temp
jx = jx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
if (noconj) then
do i = n,j + 1,-1
temp = temp - a(i,j)*x(i)
end do
if (nounit) temp = temp/a(j,j)
else
do i = n,j + 1,-1
temp = temp - conjg(a(i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
ix = kx
temp = x(jx)
if (noconj) then
do i = n,j + 1,-1
temp = temp - a(i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(j,j)
else
do i = n,j + 1,-1
temp = temp - conjg(a(i,j))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(a(j,j))
end if
x(jx) = temp
jx = jx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ci}$trsv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_stbsv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_sp
!! STBSV solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular band matrix, with ( k + 1 )
!! diagonals.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(sp), intent(in) :: a(lda,*)
real(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('STBSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed by sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
l = kplus1 - j
if (nounit) x(j) = x(j)/a(kplus1,j)
temp = x(j)
do i = j - 1,max(1,j-k),-1
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
kx = kx - incx
if (x(jx)/=zero) then
ix = kx
l = kplus1 - j
if (nounit) x(jx) = x(jx)/a(kplus1,j)
temp = x(jx)
do i = j - 1,max(1,j-k),-1
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix - incx
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
l = 1 - j
if (nounit) x(j) = x(j)/a(1,j)
temp = x(j)
do i = j + 1,min(n,j+k)
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
kx = kx + incx
if (x(jx)/=zero) then
ix = kx
l = 1 - j
if (nounit) x(jx) = x(jx)/a(1,j)
temp = x(jx)
do i = j + 1,min(n,j+k)
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix + incx
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t)*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(kplus1,j)
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(kplus1,j)
x(jx) = temp
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(1,j)
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(1,j)
x(jx) = temp
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_stbsv
pure module subroutine stdlib${ii}$_dtbsv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_dp
!! DTBSV solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular band matrix, with ( k + 1 )
!! diagonals.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(dp), intent(in) :: a(lda,*)
real(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTBSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed by sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
l = kplus1 - j
if (nounit) x(j) = x(j)/a(kplus1,j)
temp = x(j)
do i = j - 1,max(1,j-k),-1
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
kx = kx - incx
if (x(jx)/=zero) then
ix = kx
l = kplus1 - j
if (nounit) x(jx) = x(jx)/a(kplus1,j)
temp = x(jx)
do i = j - 1,max(1,j-k),-1
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix - incx
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
l = 1 - j
if (nounit) x(j) = x(j)/a(1,j)
temp = x(j)
do i = j + 1,min(n,j+k)
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
kx = kx + incx
if (x(jx)/=zero) then
ix = kx
l = 1 - j
if (nounit) x(jx) = x(jx)/a(1,j)
temp = x(jx)
do i = j + 1,min(n,j+k)
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix + incx
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t)*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(kplus1,j)
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(kplus1,j)
x(jx) = temp
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(1,j)
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(1,j)
x(jx) = temp
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_dtbsv
#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ri}$tbsv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_${rk}$
!! DTBSV: solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular band matrix, with ( k + 1 )
!! diagonals.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(${rk}$), intent(in) :: a(lda,*)
real(${rk}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(${rk}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: nounit
! Intrinsic Functions
intrinsic :: max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTBSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed by sequentially with one pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
l = kplus1 - j
if (nounit) x(j) = x(j)/a(kplus1,j)
temp = x(j)
do i = j - 1,max(1,j-k),-1
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
kx = kx - incx
if (x(jx)/=zero) then
ix = kx
l = kplus1 - j
if (nounit) x(jx) = x(jx)/a(kplus1,j)
temp = x(jx)
do i = j - 1,max(1,j-k),-1
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix - incx
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
l = 1 - j
if (nounit) x(j) = x(j)/a(1,j)
temp = x(j)
do i = j + 1,min(n,j+k)
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
kx = kx + incx
if (x(jx)/=zero) then
ix = kx
l = 1 - j
if (nounit) x(jx) = x(jx)/a(1,j)
temp = x(jx)
do i = j + 1,min(n,j+k)
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix + incx
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t)*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
temp = x(j)
l = kplus1 - j
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(kplus1,j)
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
l = kplus1 - j
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(kplus1,j)
x(jx) = temp
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = 1 - j
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(1,j)
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
l = 1 - j
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(1,j)
x(jx) = temp
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ri}$tbsv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_ctbsv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_sp
!! CTBSV solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular band matrix, with ( k + 1 )
!! diagonals.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(sp), intent(in) :: a(lda,*)
complex(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('CTBSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed by sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
l = kplus1 - j
if (nounit) x(j) = x(j)/a(kplus1,j)
temp = x(j)
do i = j - 1,max(1,j-k),-1
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
kx = kx - incx
if (x(jx)/=czero) then
ix = kx
l = kplus1 - j
if (nounit) x(jx) = x(jx)/a(kplus1,j)
temp = x(jx)
do i = j - 1,max(1,j-k),-1
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix - incx
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
l = 1 - j
if (nounit) x(j) = x(j)/a(1,j)
temp = x(j)
do i = j + 1,min(n,j+k)
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
kx = kx + incx
if (x(jx)/=czero) then
ix = kx
l = 1 - j
if (nounit) x(jx) = x(jx)/a(1,j)
temp = x(jx)
do i = j + 1,min(n,j+k)
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix + incx
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
temp = x(j)
l = kplus1 - j
if (noconj) then
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(kplus1,j)
else
do i = max(1,j-k),j - 1
temp = temp - conjg(a(l+i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(kplus1,j))
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
l = kplus1 - j
if (noconj) then
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(kplus1,j)
else
do i = max(1,j-k),j - 1
temp = temp - conjg(a(l+i,j))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(a(kplus1,j))
end if
x(jx) = temp
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = 1 - j
if (noconj) then
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(1,j)
else
do i = min(n,j+k),j + 1,-1
temp = temp - conjg(a(l+i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(1,j))
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
l = 1 - j
if (noconj) then
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(1,j)
else
do i = min(n,j+k),j + 1,-1
temp = temp - conjg(a(l+i,j))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(a(1,j))
end if
x(jx) = temp
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ctbsv
pure module subroutine stdlib${ii}$_ztbsv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_dp
!! ZTBSV solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular band matrix, with ( k + 1 )
!! diagonals.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(dp), intent(in) :: a(lda,*)
complex(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTBSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed by sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
l = kplus1 - j
if (nounit) x(j) = x(j)/a(kplus1,j)
temp = x(j)
do i = j - 1,max(1,j-k),-1
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
kx = kx - incx
if (x(jx)/=czero) then
ix = kx
l = kplus1 - j
if (nounit) x(jx) = x(jx)/a(kplus1,j)
temp = x(jx)
do i = j - 1,max(1,j-k),-1
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix - incx
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
l = 1 - j
if (nounit) x(j) = x(j)/a(1,j)
temp = x(j)
do i = j + 1,min(n,j+k)
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
kx = kx + incx
if (x(jx)/=czero) then
ix = kx
l = 1 - j
if (nounit) x(jx) = x(jx)/a(1,j)
temp = x(jx)
do i = j + 1,min(n,j+k)
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix + incx
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
temp = x(j)
l = kplus1 - j
if (noconj) then
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(kplus1,j)
else
do i = max(1,j-k),j - 1
temp = temp - conjg(a(l+i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(kplus1,j))
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
l = kplus1 - j
if (noconj) then
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(kplus1,j)
else
do i = max(1,j-k),j - 1
temp = temp - conjg(a(l+i,j))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(a(kplus1,j))
end if
x(jx) = temp
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = 1 - j
if (noconj) then
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(1,j)
else
do i = min(n,j+k),j + 1,-1
temp = temp - conjg(a(l+i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(1,j))
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
l = 1 - j
if (noconj) then
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(1,j)
else
do i = min(n,j+k),j + 1,-1
temp = temp - conjg(a(l+i,j))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(a(1,j))
end if
x(jx) = temp
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ztbsv
#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ci}$tbsv(uplo,trans,diag,n,k,a,lda,x,incx)
use stdlib_blas_constants_${ck}$
!! ZTBSV: solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular band matrix, with ( k + 1 )
!! diagonals.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, k, lda, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(${ck}$), intent(in) :: a(lda,*)
complex(${ck}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(${ck}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, kplus1, kx, l
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg,max,min
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (k<0) then
info = 5
else if (lda< (k+1)) then
info = 7
else if (incx==0) then
info = 9
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTBSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of a are
! accessed by sequentially with cone pass through a.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
l = kplus1 - j
if (nounit) x(j) = x(j)/a(kplus1,j)
temp = x(j)
do i = j - 1,max(1,j-k),-1
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
kx = kx - incx
if (x(jx)/=czero) then
ix = kx
l = kplus1 - j
if (nounit) x(jx) = x(jx)/a(kplus1,j)
temp = x(jx)
do i = j - 1,max(1,j-k),-1
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix - incx
end do
end if
jx = jx - incx
end do
end if
else
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
l = 1 - j
if (nounit) x(j) = x(j)/a(1,j)
temp = x(j)
do i = j + 1,min(n,j+k)
x(i) = x(i) - temp*a(l+i,j)
end do
end if
end do
else
jx = kx
do j = 1,n
kx = kx + incx
if (x(jx)/=czero) then
ix = kx
l = 1 - j
if (nounit) x(jx) = x(jx)/a(1,j)
temp = x(jx)
do i = j + 1,min(n,j+k)
x(ix) = x(ix) - temp*a(l+i,j)
ix = ix + incx
end do
end if
jx = jx + incx
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
kplus1 = k + 1
if (incx==1) then
do j = 1,n
temp = x(j)
l = kplus1 - j
if (noconj) then
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(kplus1,j)
else
do i = max(1,j-k),j - 1
temp = temp - conjg(a(l+i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(kplus1,j))
end if
x(j) = temp
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
l = kplus1 - j
if (noconj) then
do i = max(1,j-k),j - 1
temp = temp - a(l+i,j)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/a(kplus1,j)
else
do i = max(1,j-k),j - 1
temp = temp - conjg(a(l+i,j))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(a(kplus1,j))
end if
x(jx) = temp
jx = jx + incx
if (j>k) kx = kx + incx
end do
end if
else
if (incx==1) then
do j = n,1,-1
temp = x(j)
l = 1 - j
if (noconj) then
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(i)
end do
if (nounit) temp = temp/a(1,j)
else
do i = min(n,j+k),j + 1,-1
temp = temp - conjg(a(l+i,j))*x(i)
end do
if (nounit) temp = temp/conjg(a(1,j))
end if
x(j) = temp
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
l = 1 - j
if (noconj) then
do i = min(n,j+k),j + 1,-1
temp = temp - a(l+i,j)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/a(1,j)
else
do i = min(n,j+k),j + 1,-1
temp = temp - conjg(a(l+i,j))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(a(1,j))
end if
x(jx) = temp
jx = jx - incx
if ((n-j)>=k) kx = kx - incx
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ci}$tbsv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_stpsv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_sp
!! STPSV solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix, supplied in packed form.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(sp), intent(in) :: ap(*)
real(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: nounit
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('STPSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with one pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk - 1
do i = j - 1,1,-1
x(i) = x(i) - temp*ap(k)
k = k - 1
end do
end if
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk + 1
do i = j + 1,n
x(i) = x(i) - temp*ap(k)
k = k + 1
end do
end if
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk + 1,kk + n - j
ix = ix + incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
else
! form x := inv( a**t )*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk
do i = 1,j - 1
temp = temp - ap(k)*x(i)
k = k + 1
end do
if (nounit) temp = temp/ap(kk+j-1)
x(j) = temp
kk = kk + j
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
temp = temp - ap(k)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/ap(kk+j-1)
x(jx) = temp
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk
do i = n,j + 1,-1
temp = temp - ap(k)*x(i)
k = k - 1
end do
if (nounit) temp = temp/ap(kk-n+j)
x(j) = temp
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
temp = temp - ap(k)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/ap(kk-n+j)
x(jx) = temp
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_stpsv
pure module subroutine stdlib${ii}$_dtpsv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_dp
!! DTPSV solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix, supplied in packed form.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(dp), intent(in) :: ap(*)
real(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: nounit
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTPSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with one pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk - 1
do i = j - 1,1,-1
x(i) = x(i) - temp*ap(k)
k = k - 1
end do
end if
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk + 1
do i = j + 1,n
x(i) = x(i) - temp*ap(k)
k = k + 1
end do
end if
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk + 1,kk + n - j
ix = ix + incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
else
! form x := inv( a**t )*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk
do i = 1,j - 1
temp = temp - ap(k)*x(i)
k = k + 1
end do
if (nounit) temp = temp/ap(kk+j-1)
x(j) = temp
kk = kk + j
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
temp = temp - ap(k)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/ap(kk+j-1)
x(jx) = temp
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk
do i = n,j + 1,-1
temp = temp - ap(k)*x(i)
k = k - 1
end do
if (nounit) temp = temp/ap(kk-n+j)
x(j) = temp
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
temp = temp - ap(k)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/ap(kk-n+j)
x(jx) = temp
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_dtpsv
#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ri}$tpsv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_${rk}$
!! DTPSV: solves one of the systems of equations
!! A*x = b, or A**T*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix, supplied in packed form.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
real(${rk}$), intent(in) :: ap(*)
real(${rk}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
real(${rk}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: nounit
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('DTPSV ',info)
return
end if
! quick return if possible.
if (n==0) return
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with one pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk - 1
do i = j - 1,1,-1
x(i) = x(i) - temp*ap(k)
k = k - 1
end do
end if
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=zero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk + 1
do i = j + 1,n
x(i) = x(i) - temp*ap(k)
k = k + 1
end do
end if
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
if (x(jx)/=zero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk + 1,kk + n - j
ix = ix + incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
else
! form x := inv( a**t )*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk
do i = 1,j - 1
temp = temp - ap(k)*x(i)
k = k + 1
end do
if (nounit) temp = temp/ap(kk+j-1)
x(j) = temp
kk = kk + j
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
do k = kk,kk + j - 2
temp = temp - ap(k)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/ap(kk+j-1)
x(jx) = temp
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk
do i = n,j + 1,-1
temp = temp - ap(k)*x(i)
k = k - 1
end do
if (nounit) temp = temp/ap(kk-n+j)
x(j) = temp
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
do k = kk,kk - (n- (j+1)),-1
temp = temp - ap(k)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/ap(kk-n+j)
x(jx) = temp
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ri}$tpsv
#:endif
#:endfor
pure module subroutine stdlib${ii}$_ctpsv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_sp
!! CTPSV solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix, supplied in packed form.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(sp), intent(in) :: ap(*)
complex(sp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(sp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('CTPSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with cone pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk - 1
do i = j - 1,1,-1
x(i) = x(i) - temp*ap(k)
k = k - 1
end do
end if
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk + 1
do i = j + 1,n
x(i) = x(i) - temp*ap(k)
k = k + 1
end do
end if
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk + 1,kk + n - j
ix = ix + incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk
if (noconj) then
do i = 1,j - 1
temp = temp - ap(k)*x(i)
k = k + 1
end do
if (nounit) temp = temp/ap(kk+j-1)
else
do i = 1,j - 1
temp = temp - conjg(ap(k))*x(i)
k = k + 1
end do
if (nounit) temp = temp/conjg(ap(kk+j-1))
end if
x(j) = temp
kk = kk + j
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
if (noconj) then
do k = kk,kk + j - 2
temp = temp - ap(k)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/ap(kk+j-1)
else
do k = kk,kk + j - 2
temp = temp - conjg(ap(k))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(ap(kk+j-1))
end if
x(jx) = temp
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk
if (noconj) then
do i = n,j + 1,-1
temp = temp - ap(k)*x(i)
k = k - 1
end do
if (nounit) temp = temp/ap(kk-n+j)
else
do i = n,j + 1,-1
temp = temp - conjg(ap(k))*x(i)
k = k - 1
end do
if (nounit) temp = temp/conjg(ap(kk-n+j))
end if
x(j) = temp
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
if (noconj) then
do k = kk,kk - (n- (j+1)),-1
temp = temp - ap(k)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/ap(kk-n+j)
else
do k = kk,kk - (n- (j+1)),-1
temp = temp - conjg(ap(k))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(ap(kk-n+j))
end if
x(jx) = temp
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ctpsv
pure module subroutine stdlib${ii}$_ztpsv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_dp
!! ZTPSV solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix, supplied in packed form.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(dp), intent(in) :: ap(*)
complex(dp), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(dp) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTPSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with cone pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk - 1
do i = j - 1,1,-1
x(i) = x(i) - temp*ap(k)
k = k - 1
end do
end if
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk + 1
do i = j + 1,n
x(i) = x(i) - temp*ap(k)
k = k + 1
end do
end if
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk + 1,kk + n - j
ix = ix + incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk
if (noconj) then
do i = 1,j - 1
temp = temp - ap(k)*x(i)
k = k + 1
end do
if (nounit) temp = temp/ap(kk+j-1)
else
do i = 1,j - 1
temp = temp - conjg(ap(k))*x(i)
k = k + 1
end do
if (nounit) temp = temp/conjg(ap(kk+j-1))
end if
x(j) = temp
kk = kk + j
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
if (noconj) then
do k = kk,kk + j - 2
temp = temp - ap(k)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/ap(kk+j-1)
else
do k = kk,kk + j - 2
temp = temp - conjg(ap(k))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(ap(kk+j-1))
end if
x(jx) = temp
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk
if (noconj) then
do i = n,j + 1,-1
temp = temp - ap(k)*x(i)
k = k - 1
end do
if (nounit) temp = temp/ap(kk-n+j)
else
do i = n,j + 1,-1
temp = temp - conjg(ap(k))*x(i)
k = k - 1
end do
if (nounit) temp = temp/conjg(ap(kk-n+j))
end if
x(j) = temp
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
if (noconj) then
do k = kk,kk - (n- (j+1)),-1
temp = temp - ap(k)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/ap(kk-n+j)
else
do k = kk,kk - (n- (j+1)),-1
temp = temp - conjg(ap(k))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(ap(kk-n+j))
end if
x(jx) = temp
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_ztpsv
#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
pure module subroutine stdlib${ii}$_${ci}$tpsv(uplo,trans,diag,n,ap,x,incx)
use stdlib_blas_constants_${ck}$
!! ZTPSV: solves one of the systems of equations
!! A*x = b, or A**T*x = b, or A**H*x = b,
!! where b and x are n element vectors and A is an n by n unit, or
!! non-unit, upper or lower triangular matrix, supplied in packed form.
!! No test for singularity or near-singularity is included in this
!! routine. Such tests must be performed before calling this routine.
! -- reference blas level2 routine --
! -- reference blas is a software package provided by univ. of tennessee, --
! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
! Scalar Arguments
integer(${ik}$), intent(in) :: incx, n
character, intent(in) :: diag, trans, uplo
! Array Arguments
complex(${ck}$), intent(in) :: ap(*)
complex(${ck}$), intent(inout) :: x(*)
! =====================================================================
! Local Scalars
complex(${ck}$) :: temp
integer(${ik}$) :: i, info, ix, j, jx, k, kk, kx
logical(lk) :: noconj, nounit
! Intrinsic Functions
intrinsic :: conjg
! test the input parameters.
info = 0
if (.not.stdlib_lsame(uplo,'U') .and. .not.stdlib_lsame(uplo,'L')) then
info = 1
else if (.not.stdlib_lsame(trans,'N') .and. .not.stdlib_lsame(trans,'T') &
.and..not.stdlib_lsame(trans,'C')) then
info = 2
else if (.not.stdlib_lsame(diag,'U') .and. .not.stdlib_lsame(diag,'N')) then
info = 3
else if (n<0) then
info = 4
else if (incx==0) then
info = 7
end if
if (info/=0) then
call stdlib${ii}$_xerbla('ZTPSV ',info)
return
end if
! quick return if possible.
if (n==0) return
noconj = stdlib_lsame(trans,'T')
nounit = stdlib_lsame(diag,'N')
! set up the start point in x if the increment is not unity. this
! will be ( n - 1 )*incx too small for descending loops.
if (incx<=0) then
kx = 1 - (n-1)*incx
else if (incx/=1) then
kx = 1
end if
! start the operations. in this version the elements of ap are
! accessed sequentially with cone pass through ap.
if (stdlib_lsame(trans,'N')) then
! form x := inv( a )*x.
if (stdlib_lsame(uplo,'U')) then
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk - 1
do i = j - 1,1,-1
x(i) = x(i) - temp*ap(k)
k = k - 1
end do
end if
kk = kk - j
end do
else
jx = kx + (n-1)*incx
do j = n,1,-1
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk - 1,kk - j + 1,-1
ix = ix - incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx - incx
kk = kk - j
end do
end if
else
kk = 1
if (incx==1) then
do j = 1,n
if (x(j)/=czero) then
if (nounit) x(j) = x(j)/ap(kk)
temp = x(j)
k = kk + 1
do i = j + 1,n
x(i) = x(i) - temp*ap(k)
k = k + 1
end do
end if
kk = kk + (n-j+1)
end do
else
jx = kx
do j = 1,n
if (x(jx)/=czero) then
if (nounit) x(jx) = x(jx)/ap(kk)
temp = x(jx)
ix = jx
do k = kk + 1,kk + n - j
ix = ix + incx
x(ix) = x(ix) - temp*ap(k)
end do
end if
jx = jx + incx
kk = kk + (n-j+1)
end do
end if
end if
else
! form x := inv( a**t )*x or x := inv( a**h )*x.
if (stdlib_lsame(uplo,'U')) then
kk = 1
if (incx==1) then
do j = 1,n
temp = x(j)
k = kk
if (noconj) then
do i = 1,j - 1
temp = temp - ap(k)*x(i)
k = k + 1
end do
if (nounit) temp = temp/ap(kk+j-1)
else
do i = 1,j - 1
temp = temp - conjg(ap(k))*x(i)
k = k + 1
end do
if (nounit) temp = temp/conjg(ap(kk+j-1))
end if
x(j) = temp
kk = kk + j
end do
else
jx = kx
do j = 1,n
temp = x(jx)
ix = kx
if (noconj) then
do k = kk,kk + j - 2
temp = temp - ap(k)*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/ap(kk+j-1)
else
do k = kk,kk + j - 2
temp = temp - conjg(ap(k))*x(ix)
ix = ix + incx
end do
if (nounit) temp = temp/conjg(ap(kk+j-1))
end if
x(jx) = temp
jx = jx + incx
kk = kk + j
end do
end if
else
kk = (n* (n+1))/2
if (incx==1) then
do j = n,1,-1
temp = x(j)
k = kk
if (noconj) then
do i = n,j + 1,-1
temp = temp - ap(k)*x(i)
k = k - 1
end do
if (nounit) temp = temp/ap(kk-n+j)
else
do i = n,j + 1,-1
temp = temp - conjg(ap(k))*x(i)
k = k - 1
end do
if (nounit) temp = temp/conjg(ap(kk-n+j))
end if
x(j) = temp
kk = kk - (n-j+1)
end do
else
kx = kx + (n-1)*incx
jx = kx
do j = n,1,-1
temp = x(jx)
ix = kx
if (noconj) then
do k = kk,kk - (n- (j+1)),-1
temp = temp - ap(k)*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/ap(kk-n+j)
else
do k = kk,kk - (n- (j+1)),-1
temp = temp - conjg(ap(k))*x(ix)
ix = ix - incx
end do
if (nounit) temp = temp/conjg(ap(kk-n+j))
end if
x(jx) = temp
jx = jx - incx
kk = kk - (n-j+1)
end do
end if
end if
end if
return
end subroutine stdlib${ii}$_${ci}$tpsv
#:endif
#:endfor
#:endfor
end submodule stdlib_blas_level2_tri
