#:include "common.fypp" submodule(stdlib_lapack_base) stdlib_lapack_blas_like_l1 implicit none contains #:for ik,it,ii in LINALG_INT_KINDS_TYPES pure module subroutine stdlib${ii}$_clacgv( n, x, incx ) !! CLACGV conjugates a complex vector of length N. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n ! Array Arguments complex(sp), intent(inout) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ioff ! Intrinsic Functions ! Executable Statements if( incx==1_${ik}$ ) then do i = 1, n x( i ) = conjg( x( i ) ) end do else ioff = 1_${ik}$ if( incx<0_${ik}$ )ioff = 1_${ik}$ - ( n-1 )*incx do i = 1, n x( ioff ) = conjg( x( ioff ) ) ioff = ioff + incx end do end if return end subroutine stdlib${ii}$_clacgv pure module subroutine stdlib${ii}$_zlacgv( n, x, incx ) !! ZLACGV conjugates a complex vector of length N. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n ! Array Arguments complex(dp), intent(inout) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ioff ! Intrinsic Functions ! Executable Statements if( incx==1_${ik}$ ) then do i = 1, n x( i ) = conjg( x( i ) ) end do else ioff = 1_${ik}$ if( incx<0_${ik}$ )ioff = 1_${ik}$ - ( n-1 )*incx do i = 1, n x( ioff ) = conjg( x( ioff ) ) ioff = ioff + incx end do end if return end subroutine stdlib${ii}$_zlacgv #:for ck,ct,ci in CMPLX_KINDS_TYPES #:if not ck in ["sp","dp"] pure module subroutine stdlib${ii}$_${ci}$lacgv( n, x, incx ) !! ZLACGV: conjugates a complex vector of length N. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_${ck}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n ! Array Arguments complex(${ck}$), intent(inout) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ioff ! Intrinsic Functions ! Executable Statements if( incx==1_${ik}$ ) then do i = 1, n x( i ) = conjg( x( i ) ) end do else ioff = 1_${ik}$ if( incx<0_${ik}$ )ioff = 1_${ik}$ - ( n-1 )*incx do i = 1, n x( ioff ) = conjg( x( ioff ) ) ioff = ioff + incx end do end if return end subroutine stdlib${ii}$_${ci}$lacgv #:endif #:endfor pure module subroutine stdlib${ii}$_slasrt( id, n, d, info ) !! Sort the numbers in D in increasing order (if ID = 'I') or !! in decreasing order (if ID = 'D' ). !! Use Quick Sort, reverting to Insertion sort on arrays of !! size <= 20. Dimension of STACK limits N to about 2**32. ! -- lapack computational routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments character, intent(in) :: id integer(${ik}$), intent(out) :: info integer(${ik}$), intent(in) :: n ! Array Arguments real(sp), intent(inout) :: d(*) ! ===================================================================== ! Parameters integer(${ik}$), parameter :: select = 20_${ik}$ ! Local Scalars integer(${ik}$) :: dir, endd, i, j, start, stkpnt real(sp) :: d1, d2, d3, dmnmx, tmp ! Local Arrays integer(${ik}$) :: stack(2_${ik}$,32_${ik}$) ! Executable Statements ! test the input parameters. info = 0_${ik}$ dir = -1_${ik}$ if( stdlib_lsame( id, 'D' ) ) then dir = 0_${ik}$ else if( stdlib_lsame( id, 'I' ) ) then dir = 1_${ik}$ end if if( dir==-1_${ik}$ ) then info = -1_${ik}$ else if( n<0_${ik}$ ) then info = -2_${ik}$ end if if( info/=0_${ik}$ ) then call stdlib${ii}$_xerbla( 'SLASRT', -info ) return end if ! quick return if possible if( n<=1 )return stkpnt = 1_${ik}$ stack( 1_${ik}$, 1_${ik}$ ) = 1_${ik}$ stack( 2_${ik}$, 1_${ik}$ ) = n 10 continue start = stack( 1_${ik}$, stkpnt ) endd = stack( 2_${ik}$, stkpnt ) stkpnt = stkpnt - 1_${ik}$ if( endd-start<=select .and. endd-start>0_${ik}$ ) then ! do insertion sort on d( start:endd ) if( dir==0_${ik}$ ) then ! sort into decreasing order loop_30: do i = start + 1, endd do j = i, start + 1, -1 if( d( j )>d( j-1 ) ) then dmnmx = d( j ) d( j ) = d( j-1 ) d( j-1 ) = dmnmx else cycle loop_30 end if end do end do loop_30 else ! sort into increasing order loop_50: do i = start + 1, endd do j = i, start + 1, -1 if( d( j )<d( j-1 ) ) then dmnmx = d( j ) d( j ) = d( j-1 ) d( j-1 ) = dmnmx else cycle loop_50 end if end do end do loop_50 end if else if( endd-start>select ) then ! partition d( start:endd ) and stack parts, largest one first ! choose partition entry as median of 3 d1 = d( start ) d2 = d( endd ) i = ( start+endd ) / 2_${ik}$ d3 = d( i ) if( d1<d2 ) then if( d3<d1 ) then dmnmx = d1 else if( d3<d2 ) then dmnmx = d3 else dmnmx = d2 end if else if( d3<d2 ) then dmnmx = d2 else if( d3<d1 ) then dmnmx = d3 else dmnmx = d1 end if end if if( dir==0_${ik}$ ) then ! sort into decreasing order i = start - 1_${ik}$ j = endd + 1_${ik}$ 60 continue 70 continue j = j - 1_${ik}$ if( d( j )<dmnmx )go to 70 80 continue i = i + 1_${ik}$ if( d( i )>dmnmx )go to 80 if( i<j ) then tmp = d( i ) d( i ) = d( j ) d( j ) = tmp go to 60 end if if( j-start>endd-j-1 ) then stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd else stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j end if else ! sort into increasing order i = start - 1_${ik}$ j = endd + 1_${ik}$ 90 continue 100 continue j = j - 1_${ik}$ if( d( j )>dmnmx )go to 100 110 continue i = i + 1_${ik}$ if( d( i )<dmnmx )go to 110 if( i<j ) then tmp = d( i ) d( i ) = d( j ) d( j ) = tmp go to 90 end if if( j-start>endd-j-1 ) then stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd else stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j end if end if end if if( stkpnt>0 )go to 10 return end subroutine stdlib${ii}$_slasrt pure module subroutine stdlib${ii}$_dlasrt( id, n, d, info ) !! Sort the numbers in D in increasing order (if ID = 'I') or !! in decreasing order (if ID = 'D' ). !! Use Quick Sort, reverting to Insertion sort on arrays of !! size <= 20. Dimension of STACK limits N to about 2**32. ! -- lapack computational routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments character, intent(in) :: id integer(${ik}$), intent(out) :: info integer(${ik}$), intent(in) :: n ! Array Arguments real(dp), intent(inout) :: d(*) ! ===================================================================== ! Parameters integer(${ik}$), parameter :: select = 20_${ik}$ ! Local Scalars integer(${ik}$) :: dir, endd, i, j, start, stkpnt real(dp) :: d1, d2, d3, dmnmx, tmp ! Local Arrays integer(${ik}$) :: stack(2_${ik}$,32_${ik}$) ! Executable Statements ! test the input parameters. info = 0_${ik}$ dir = -1_${ik}$ if( stdlib_lsame( id, 'D' ) ) then dir = 0_${ik}$ else if( stdlib_lsame( id, 'I' ) ) then dir = 1_${ik}$ end if if( dir==-1_${ik}$ ) then info = -1_${ik}$ else if( n<0_${ik}$ ) then info = -2_${ik}$ end if if( info/=0_${ik}$ ) then call stdlib${ii}$_xerbla( 'DLASRT', -info ) return end if ! quick return if possible if( n<=1 )return stkpnt = 1_${ik}$ stack( 1_${ik}$, 1_${ik}$ ) = 1_${ik}$ stack( 2_${ik}$, 1_${ik}$ ) = n 10 continue start = stack( 1_${ik}$, stkpnt ) endd = stack( 2_${ik}$, stkpnt ) stkpnt = stkpnt - 1_${ik}$ if( endd-start<=select .and. endd-start>0_${ik}$ ) then ! do insertion sort on d( start:endd ) if( dir==0_${ik}$ ) then ! sort into decreasing order loop_30: do i = start + 1, endd do j = i, start + 1, -1 if( d( j )>d( j-1 ) ) then dmnmx = d( j ) d( j ) = d( j-1 ) d( j-1 ) = dmnmx else cycle loop_30 end if end do end do loop_30 else ! sort into increasing order loop_50: do i = start + 1, endd do j = i, start + 1, -1 if( d( j )<d( j-1 ) ) then dmnmx = d( j ) d( j ) = d( j-1 ) d( j-1 ) = dmnmx else cycle loop_50 end if end do end do loop_50 end if else if( endd-start>select ) then ! partition d( start:endd ) and stack parts, largest one first ! choose partition entry as median of 3 d1 = d( start ) d2 = d( endd ) i = ( start+endd ) / 2_${ik}$ d3 = d( i ) if( d1<d2 ) then if( d3<d1 ) then dmnmx = d1 else if( d3<d2 ) then dmnmx = d3 else dmnmx = d2 end if else if( d3<d2 ) then dmnmx = d2 else if( d3<d1 ) then dmnmx = d3 else dmnmx = d1 end if end if if( dir==0_${ik}$ ) then ! sort into decreasing order i = start - 1_${ik}$ j = endd + 1_${ik}$ 60 continue 70 continue j = j - 1_${ik}$ if( d( j )<dmnmx )go to 70 80 continue i = i + 1_${ik}$ if( d( i )>dmnmx )go to 80 if( i<j ) then tmp = d( i ) d( i ) = d( j ) d( j ) = tmp go to 60 end if if( j-start>endd-j-1 ) then stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd else stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j end if else ! sort into increasing order i = start - 1_${ik}$ j = endd + 1_${ik}$ 90 continue 100 continue j = j - 1_${ik}$ if( d( j )>dmnmx )go to 100 110 continue i = i + 1_${ik}$ if( d( i )<dmnmx )go to 110 if( i<j ) then tmp = d( i ) d( i ) = d( j ) d( j ) = tmp go to 90 end if if( j-start>endd-j-1 ) then stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd else stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j end if end if end if if( stkpnt>0 )go to 10 return end subroutine stdlib${ii}$_dlasrt #:for rk,rt,ri in REAL_KINDS_TYPES #:if not rk in ["sp","dp"] pure module subroutine stdlib${ii}$_${ri}$lasrt( id, n, d, info ) !! Sort the numbers in D in increasing order (if ID = 'I') or !! in decreasing order (if ID = 'D' ). !! Use Quick Sort, reverting to Insertion sort on arrays of !! size <= 20. Dimension of STACK limits N to about 2**32. ! -- lapack computational routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_${rk}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments character, intent(in) :: id integer(${ik}$), intent(out) :: info integer(${ik}$), intent(in) :: n ! Array Arguments real(${rk}$), intent(inout) :: d(*) ! ===================================================================== ! Parameters integer(${ik}$), parameter :: select = 20_${ik}$ ! Local Scalars integer(${ik}$) :: dir, endd, i, j, start, stkpnt real(${rk}$) :: d1, d2, d3, dmnmx, tmp ! Local Arrays integer(${ik}$) :: stack(2_${ik}$,32_${ik}$) ! Executable Statements ! test the input parameters. info = 0_${ik}$ dir = -1_${ik}$ if( stdlib_lsame( id, 'D' ) ) then dir = 0_${ik}$ else if( stdlib_lsame( id, 'I' ) ) then dir = 1_${ik}$ end if if( dir==-1_${ik}$ ) then info = -1_${ik}$ else if( n<0_${ik}$ ) then info = -2_${ik}$ end if if( info/=0_${ik}$ ) then call stdlib${ii}$_xerbla( 'DLASRT', -info ) return end if ! quick return if possible if( n<=1 )return stkpnt = 1_${ik}$ stack( 1_${ik}$, 1_${ik}$ ) = 1_${ik}$ stack( 2_${ik}$, 1_${ik}$ ) = n 10 continue start = stack( 1_${ik}$, stkpnt ) endd = stack( 2_${ik}$, stkpnt ) stkpnt = stkpnt - 1_${ik}$ if( endd-start<=select .and. endd-start>0_${ik}$ ) then ! do insertion sort on d( start:endd ) if( dir==0_${ik}$ ) then ! sort into decreasing order loop_30: do i = start + 1, endd do j = i, start + 1, -1 if( d( j )>d( j-1 ) ) then dmnmx = d( j ) d( j ) = d( j-1 ) d( j-1 ) = dmnmx else cycle loop_30 end if end do end do loop_30 else ! sort into increasing order loop_50: do i = start + 1, endd do j = i, start + 1, -1 if( d( j )<d( j-1 ) ) then dmnmx = d( j ) d( j ) = d( j-1 ) d( j-1 ) = dmnmx else cycle loop_50 end if end do end do loop_50 end if else if( endd-start>select ) then ! partition d( start:endd ) and stack parts, largest one first ! choose partition entry as median of 3 d1 = d( start ) d2 = d( endd ) i = ( start+endd ) / 2_${ik}$ d3 = d( i ) if( d1<d2 ) then if( d3<d1 ) then dmnmx = d1 else if( d3<d2 ) then dmnmx = d3 else dmnmx = d2 end if else if( d3<d2 ) then dmnmx = d2 else if( d3<d1 ) then dmnmx = d3 else dmnmx = d1 end if end if if( dir==0_${ik}$ ) then ! sort into decreasing order i = start - 1_${ik}$ j = endd + 1_${ik}$ 60 continue 70 continue j = j - 1_${ik}$ if( d( j )<dmnmx )go to 70 80 continue i = i + 1_${ik}$ if( d( i )>dmnmx )go to 80 if( i<j ) then tmp = d( i ) d( i ) = d( j ) d( j ) = tmp go to 60 end if if( j-start>endd-j-1 ) then stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd else stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j end if else ! sort into increasing order i = start - 1_${ik}$ j = endd + 1_${ik}$ 90 continue 100 continue j = j - 1_${ik}$ if( d( j )>dmnmx )go to 100 110 continue i = i + 1_${ik}$ if( d( i )<dmnmx )go to 110 if( i<j ) then tmp = d( i ) d( i ) = d( j ) d( j ) = tmp go to 90 end if if( j-start>endd-j-1 ) then stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd else stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = j + 1_${ik}$ stack( 2_${ik}$, stkpnt ) = endd stkpnt = stkpnt + 1_${ik}$ stack( 1_${ik}$, stkpnt ) = start stack( 2_${ik}$, stkpnt ) = j end if end if end if if( stkpnt>0 )go to 10 return end subroutine stdlib${ii}$_${ri}$lasrt #:endif #:endfor pure module subroutine stdlib${ii}$_slassq( n, x, incx, scl, sumsq ) !! SLASSQ returns the values scl and smsq such that !! ( scl**2 )*smsq = x( 1 )**2 +...+ x( n )**2 + ( scale**2 )*sumsq, !! where x( i ) = X( 1 + ( i - 1 )*INCX ). The value of sumsq is !! assumed to be non-negative. !! scale and sumsq must be supplied in SCALE and SUMSQ and !! scl and smsq are overwritten on SCALE and SUMSQ respectively. !! If scale * sqrt( sumsq ) > tbig then !! we require: scale >= sqrt( TINY*EPS ) / sbig on entry, !! and if 0 < scale * sqrt( sumsq ) < tsml then !! we require: scale <= sqrt( HUGE ) / ssml on entry, !! where !! tbig -- upper threshold for values whose square is representable; !! sbig -- scaling constant for big numbers; \see la_constants.f90 !! tsml -- lower threshold for values whose square is representable; !! ssml -- scaling constant for small numbers; \see la_constants.f90 !! and !! TINY*EPS -- tiniest representable number; !! HUGE -- biggest representable number. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_sp, only: zero, half, one, two, tbig, tsml, ssml, sbig ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(sp), intent(inout) :: scl, sumsq ! Array Arguments real(sp), intent(in) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ix logical(lk) :: notbig real(sp) :: abig, amed, asml, ax, ymax, ymin ! quick return if possible if( ieee_is_nan(scl) .or. ieee_is_nan(sumsq) ) return if( sumsq == zero ) scl = one if( scl == zero ) then scl = one sumsq = zero end if if (n <= 0_${ik}$) then return end if ! compute the sum of squares in 3 accumulators: ! abig -- sums of squares scaled down to avoid overflow ! asml -- sums of squares scaled up to avoid underflow ! amed -- sums of squares that do not require scaling ! the thresholds and multipliers are ! tbig -- values bigger than this are scaled down by sbig ! tsml -- values smaller than this are scaled up by ssml notbig = .true. asml = zero amed = zero abig = zero ix = 1_${ik}$ if( incx < 0_${ik}$ ) ix = 1_${ik}$ - (n-1)*incx do i = 1, n ax = abs(x(ix)) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ix = ix + incx end do ! put the existing sum of squares into one of the accumulators if( sumsq > zero ) then ax = scl*sqrt( sumsq ) if (ax > tbig) then ! we assume scl >= sqrt( tiny*eps ) / sbig abig = abig + (scl*sbig)**2_${ik}$ * sumsq else if (ax < tsml) then ! we assume scl <= sqrt( huge ) / ssml if (notbig) asml = asml + (scl*ssml)**2_${ik}$ * sumsq else amed = amed + scl**2_${ik}$ * sumsq end if end if ! combine abig and amed or amed and asml if more than one ! accumulator was used. if (abig > zero) then ! combine abig and amed if abig > 0. if (amed > zero .or. ieee_is_nan(amed)) then abig = abig + (amed*sbig)*sbig end if scl = one / sbig sumsq = abig else if (asml > zero) then ! combine amed and asml if asml > 0. if (amed > zero .or. ieee_is_nan(amed)) then amed = sqrt(amed) asml = sqrt(asml) / ssml if (asml > amed) then ymin = amed ymax = asml else ymin = asml ymax = amed end if scl = one sumsq = ymax**2_${ik}$*( one + (ymin/ymax)**2_${ik}$ ) else scl = one / ssml sumsq = asml end if else ! otherwise all values are mid-range or zero scl = one sumsq = amed end if return end subroutine stdlib${ii}$_slassq pure module subroutine stdlib${ii}$_dlassq( n, x, incx, scl, sumsq ) !! DLASSQ returns the values scl and smsq such that !! ( scl**2 )*smsq = x( 1 )**2 +...+ x( n )**2 + ( scale**2 )*sumsq, !! where x( i ) = X( 1 + ( i - 1 )*INCX ). The value of sumsq is !! assumed to be non-negative. !! scale and sumsq must be supplied in SCALE and SUMSQ and !! scl and smsq are overwritten on SCALE and SUMSQ respectively. !! If scale * sqrt( sumsq ) > tbig then !! we require: scale >= sqrt( TINY*EPS ) / sbig on entry, !! and if 0 < scale * sqrt( sumsq ) < tsml then !! we require: scale <= sqrt( HUGE ) / ssml on entry, !! where !! tbig -- upper threshold for values whose square is representable; !! sbig -- scaling constant for big numbers; \see la_constants.f90 !! tsml -- lower threshold for values whose square is representable; !! ssml -- scaling constant for small numbers; \see la_constants.f90 !! and !! TINY*EPS -- tiniest representable number; !! HUGE -- biggest representable number. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_dp, only: zero, half, one, two, tbig, tsml, ssml, sbig ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(dp), intent(inout) :: scl, sumsq ! Array Arguments real(dp), intent(in) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ix logical(lk) :: notbig real(dp) :: abig, amed, asml, ax, ymax, ymin ! quick return if possible if( ieee_is_nan(scl) .or. ieee_is_nan(sumsq) ) return if( sumsq == zero ) scl = one if( scl == zero ) then scl = one sumsq = zero end if if (n <= 0_${ik}$) then return end if ! compute the sum of squares in 3 accumulators: ! abig -- sums of squares scaled down to avoid overflow ! asml -- sums of squares scaled up to avoid underflow ! amed -- sums of squares that do not require scaling ! the thresholds and multipliers are ! tbig -- values bigger than this are scaled down by sbig ! tsml -- values smaller than this are scaled up by ssml notbig = .true. asml = zero amed = zero abig = zero ix = 1_${ik}$ if( incx < 0_${ik}$ ) ix = 1_${ik}$ - (n-1)*incx do i = 1, n ax = abs(x(ix)) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ix = ix + incx end do ! put the existing sum of squares into one of the accumulators if( sumsq > zero ) then ax = scl*sqrt( sumsq ) if (ax > tbig) then ! we assume scl >= sqrt( tiny*eps ) / sbig abig = abig + (scl*sbig)**2_${ik}$ * sumsq else if (ax < tsml) then ! we assume scl <= sqrt( huge ) / ssml if (notbig) asml = asml + (scl*ssml)**2_${ik}$ * sumsq else amed = amed + scl**2_${ik}$ * sumsq end if end if ! combine abig and amed or amed and asml if more than one ! accumulator was used. if (abig > zero) then ! combine abig and amed if abig > 0. if (amed > zero .or. ieee_is_nan(amed)) then abig = abig + (amed*sbig)*sbig end if scl = one / sbig sumsq = abig else if (asml > zero) then ! combine amed and asml if asml > 0. if (amed > zero .or. ieee_is_nan(amed)) then amed = sqrt(amed) asml = sqrt(asml) / ssml if (asml > amed) then ymin = amed ymax = asml else ymin = asml ymax = amed end if scl = one sumsq = ymax**2_${ik}$*( one + (ymin/ymax)**2_${ik}$ ) else scl = one / ssml sumsq = asml end if else ! otherwise all values are mid-range or zero scl = one sumsq = amed end if return end subroutine stdlib${ii}$_dlassq #:for rk,rt,ri in REAL_KINDS_TYPES #:if not rk in ["sp","dp"] pure module subroutine stdlib${ii}$_${ri}$lassq( n, x, incx, scl, sumsq ) !! DLASSQ: returns the values scl and smsq such that !! ( scl**2 )*smsq = x( 1 )**2 +...+ x( n )**2 + ( scale**2 )*sumsq, !! where x( i ) = X( 1 + ( i - 1 )*INCX ). The value of sumsq is !! assumed to be non-negative. !! scale and sumsq must be supplied in SCALE and SUMSQ and !! scl and smsq are overwritten on SCALE and SUMSQ respectively. !! If scale * sqrt( sumsq ) > tbig then !! we require: scale >= sqrt( TINY*EPS ) / sbig on entry, !! and if 0 < scale * sqrt( sumsq ) < tsml then !! we require: scale <= sqrt( HUGE ) / ssml on entry, !! where !! tbig -- upper threshold for values whose square is representable; !! sbig -- scaling constant for big numbers; \see la_constants.f90 !! tsml -- lower threshold for values whose square is representable; !! ssml -- scaling constant for small numbers; \see la_constants.f90 !! and !! TINY*EPS -- tiniest representable number; !! HUGE -- biggest representable number. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_${rk}$, only: zero, half, one, two, tbig, tsml, ssml, sbig ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(${rk}$), intent(inout) :: scl, sumsq ! Array Arguments real(${rk}$), intent(in) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ix logical(lk) :: notbig real(${rk}$) :: abig, amed, asml, ax, ymax, ymin ! quick return if possible if( ieee_is_nan(scl) .or. ieee_is_nan(sumsq) ) return if( sumsq == zero ) scl = one if( scl == zero ) then scl = one sumsq = zero end if if (n <= 0_${ik}$) then return end if ! compute the sum of squares in 3 accumulators: ! abig -- sums of squares scaled down to avoid overflow ! asml -- sums of squares scaled up to avoid underflow ! amed -- sums of squares that do not require scaling ! the thresholds and multipliers are ! tbig -- values bigger than this are scaled down by sbig ! tsml -- values smaller than this are scaled up by ssml notbig = .true. asml = zero amed = zero abig = zero ix = 1_${ik}$ if( incx < 0_${ik}$ ) ix = 1_${ik}$ - (n-1)*incx do i = 1, n ax = abs(x(ix)) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ix = ix + incx end do ! put the existing sum of squares into one of the accumulators if( sumsq > zero ) then ax = scl*sqrt( sumsq ) if (ax > tbig) then ! we assume scl >= sqrt( tiny*eps ) / sbig abig = abig + (scl*sbig)**2_${ik}$ * sumsq else if (ax < tsml) then ! we assume scl <= sqrt( huge ) / ssml if (notbig) asml = asml + (scl*ssml)**2_${ik}$ * sumsq else amed = amed + scl**2_${ik}$ * sumsq end if end if ! combine abig and amed or amed and asml if more than one ! accumulator was used. if (abig > zero) then ! combine abig and amed if abig > 0. if (amed > zero .or. ieee_is_nan(amed)) then abig = abig + (amed*sbig)*sbig end if scl = one / sbig sumsq = abig else if (asml > zero) then ! combine amed and asml if asml > 0. if (amed > zero .or. ieee_is_nan(amed)) then amed = sqrt(amed) asml = sqrt(asml) / ssml if (asml > amed) then ymin = amed ymax = asml else ymin = asml ymax = amed end if scl = one sumsq = ymax**2_${ik}$*( one + (ymin/ymax)**2_${ik}$ ) else scl = one / ssml sumsq = asml end if else ! otherwise all values are mid-range or zero scl = one sumsq = amed end if return end subroutine stdlib${ii}$_${ri}$lassq #:endif #:endfor pure module subroutine stdlib${ii}$_classq( n, x, incx, scl, sumsq ) !! CLASSQ returns the values scl and smsq such that !! ( scl**2 )*smsq = x( 1 )**2 +...+ x( n )**2 + ( scale**2 )*sumsq, !! where x( i ) = X( 1 + ( i - 1 )*INCX ). The value of sumsq is !! assumed to be non-negative. !! scale and sumsq must be supplied in SCALE and SUMSQ and !! scl and smsq are overwritten on SCALE and SUMSQ respectively. !! If scale * sqrt( sumsq ) > tbig then !! we require: scale >= sqrt( TINY*EPS ) / sbig on entry, !! and if 0 < scale * sqrt( sumsq ) < tsml then !! we require: scale <= sqrt( HUGE ) / ssml on entry, !! where !! tbig -- upper threshold for values whose square is representable; !! sbig -- scaling constant for big numbers; \see la_constants.f90 !! tsml -- lower threshold for values whose square is representable; !! ssml -- scaling constant for small numbers; \see la_constants.f90 !! and !! TINY*EPS -- tiniest representable number; !! HUGE -- biggest representable number. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_sp, only: zero, half, one, two, tbig, tsml, ssml, sbig ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(sp), intent(inout) :: scl, sumsq ! Array Arguments complex(sp), intent(in) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ix logical(lk) :: notbig real(sp) :: abig, amed, asml, ax, ymax, ymin ! quick return if possible if( ieee_is_nan(scl) .or. ieee_is_nan(sumsq) ) return if( sumsq == zero ) scl = one if( scl == zero ) then scl = one sumsq = zero end if if (n <= 0_${ik}$) then return end if ! compute the sum of squares in 3 accumulators: ! abig -- sums of squares scaled down to avoid overflow ! asml -- sums of squares scaled up to avoid underflow ! amed -- sums of squares that do not require scaling ! the thresholds and multipliers are ! tbig -- values bigger than this are scaled down by sbig ! tsml -- values smaller than this are scaled up by ssml notbig = .true. asml = zero amed = zero abig = zero ix = 1_${ik}$ if( incx < 0_${ik}$ ) ix = 1_${ik}$ - (n-1)*incx do i = 1, n ax = abs(real(x(ix),KIND=sp)) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ax = abs(aimag(x(ix))) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ix = ix + incx end do ! put the existing sum of squares into one of the accumulators if( sumsq > zero ) then ax = scl*sqrt( sumsq ) if (ax > tbig) then ! we assume scl >= sqrt( tiny*eps ) / sbig abig = abig + (scl*sbig)**2_${ik}$ * sumsq else if (ax < tsml) then ! we assume scl <= sqrt( huge ) / ssml if (notbig) asml = asml + (scl*ssml)**2_${ik}$ * sumsq else amed = amed + scl**2_${ik}$ * sumsq end if end if ! combine abig and amed or amed and asml if more than one ! accumulator was used. if (abig > zero) then ! combine abig and amed if abig > 0. if (amed > zero .or. ieee_is_nan(amed)) then abig = abig + (amed*sbig)*sbig end if scl = one / sbig sumsq = abig else if (asml > zero) then ! combine amed and asml if asml > 0. if (amed > zero .or. ieee_is_nan(amed)) then amed = sqrt(amed) asml = sqrt(asml) / ssml if (asml > amed) then ymin = amed ymax = asml else ymin = asml ymax = amed end if scl = one sumsq = ymax**2_${ik}$*( one + (ymin/ymax)**2_${ik}$ ) else scl = one / ssml sumsq = asml end if else ! otherwise all values are mid-range or zero scl = one sumsq = amed end if return end subroutine stdlib${ii}$_classq pure module subroutine stdlib${ii}$_zlassq( n, x, incx, scl, sumsq ) !! ZLASSQ returns the values scl and smsq such that !! ( scl**2 )*smsq = x( 1 )**2 +...+ x( n )**2 + ( scale**2 )*sumsq, !! where x( i ) = X( 1 + ( i - 1 )*INCX ). The value of sumsq is !! assumed to be non-negative. !! scale and sumsq must be supplied in SCALE and SUMSQ and !! scl and smsq are overwritten on SCALE and SUMSQ respectively. !! If scale * sqrt( sumsq ) > tbig then !! we require: scale >= sqrt( TINY*EPS ) / sbig on entry, !! and if 0 < scale * sqrt( sumsq ) < tsml then !! we require: scale <= sqrt( HUGE ) / ssml on entry, !! where !! tbig -- upper threshold for values whose square is representable; !! sbig -- scaling constant for big numbers; \see la_constants.f90 !! tsml -- lower threshold for values whose square is representable; !! ssml -- scaling constant for small numbers; \see la_constants.f90 !! and !! TINY*EPS -- tiniest representable number; !! HUGE -- biggest representable number. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_dp, only: zero, half, one, two, tbig, tsml, ssml, sbig ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(dp), intent(inout) :: scl, sumsq ! Array Arguments complex(dp), intent(in) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ix logical(lk) :: notbig real(dp) :: abig, amed, asml, ax, ymax, ymin ! quick return if possible if( ieee_is_nan(scl) .or. ieee_is_nan(sumsq) ) return if( sumsq == zero ) scl = one if( scl == zero ) then scl = one sumsq = zero end if if (n <= 0_${ik}$) then return end if ! compute the sum of squares in 3 accumulators: ! abig -- sums of squares scaled down to avoid overflow ! asml -- sums of squares scaled up to avoid underflow ! amed -- sums of squares that do not require scaling ! the thresholds and multipliers are ! tbig -- values bigger than this are scaled down by sbig ! tsml -- values smaller than this are scaled up by ssml notbig = .true. asml = zero amed = zero abig = zero ix = 1_${ik}$ if( incx < 0_${ik}$ ) ix = 1_${ik}$ - (n-1)*incx do i = 1, n ax = abs(real(x(ix),KIND=dp)) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ax = abs(aimag(x(ix))) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ix = ix + incx end do ! put the existing sum of squares into one of the accumulators if( sumsq > zero ) then ax = scl*sqrt( sumsq ) if (ax > tbig) then ! we assume scl >= sqrt( tiny*eps ) / sbig abig = abig + (scl*sbig)**2_${ik}$ * sumsq else if (ax < tsml) then ! we assume scl <= sqrt( huge ) / ssml if (notbig) asml = asml + (scl*ssml)**2_${ik}$ * sumsq else amed = amed + scl**2_${ik}$ * sumsq end if end if ! combine abig and amed or amed and asml if more than one ! accumulator was used. if (abig > zero) then ! combine abig and amed if abig > 0. if (amed > zero .or. ieee_is_nan(amed)) then abig = abig + (amed*sbig)*sbig end if scl = one / sbig sumsq = abig else if (asml > zero) then ! combine amed and asml if asml > 0. if (amed > zero .or. ieee_is_nan(amed)) then amed = sqrt(amed) asml = sqrt(asml) / ssml if (asml > amed) then ymin = amed ymax = asml else ymin = asml ymax = amed end if scl = one sumsq = ymax**2_${ik}$*( one + (ymin/ymax)**2_${ik}$ ) else scl = one / ssml sumsq = asml end if else ! otherwise all values are mid-range or zero scl = one sumsq = amed end if return end subroutine stdlib${ii}$_zlassq #:for ck,ct,ci in CMPLX_KINDS_TYPES #:if not ck in ["sp","dp"] pure module subroutine stdlib${ii}$_${ci}$lassq( n, x, incx, scl, sumsq ) !! ZLASSQ: returns the values scl and smsq such that !! ( scl**2 )*smsq = x( 1 )**2 +...+ x( n )**2 + ( scale**2 )*sumsq, !! where x( i ) = X( 1 + ( i - 1 )*INCX ). The value of sumsq is !! assumed to be non-negative. !! scale and sumsq must be supplied in SCALE and SUMSQ and !! scl and smsq are overwritten on SCALE and SUMSQ respectively. !! If scale * sqrt( sumsq ) > tbig then !! we require: scale >= sqrt( TINY*EPS ) / sbig on entry, !! and if 0 < scale * sqrt( sumsq ) < tsml then !! we require: scale <= sqrt( HUGE ) / ssml on entry, !! where !! tbig -- upper threshold for values whose square is representable; !! sbig -- scaling constant for big numbers; \see la_constants.f90 !! tsml -- lower threshold for values whose square is representable; !! ssml -- scaling constant for small numbers; \see la_constants.f90 !! and !! TINY*EPS -- tiniest representable number; !! HUGE -- biggest representable number. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_${ck}$, only: zero, half, one, two, tbig, tsml, ssml, sbig ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(${ck}$), intent(inout) :: scl, sumsq ! Array Arguments complex(${ck}$), intent(in) :: x(*) ! ===================================================================== ! Local Scalars integer(${ik}$) :: i, ix logical(lk) :: notbig real(${ck}$) :: abig, amed, asml, ax, ymax, ymin ! quick return if possible if( ieee_is_nan(scl) .or. ieee_is_nan(sumsq) ) return if( sumsq == zero ) scl = one if( scl == zero ) then scl = one sumsq = zero end if if (n <= 0_${ik}$) then return end if ! compute the sum of squares in 3 accumulators: ! abig -- sums of squares scaled down to avoid overflow ! asml -- sums of squares scaled up to avoid underflow ! amed -- sums of squares that do not require scaling ! the thresholds and multipliers are ! tbig -- values bigger than this are scaled down by sbig ! tsml -- values smaller than this are scaled up by ssml notbig = .true. asml = zero amed = zero abig = zero ix = 1_${ik}$ if( incx < 0_${ik}$ ) ix = 1_${ik}$ - (n-1)*incx do i = 1, n ax = abs(real(x(ix),KIND=${ck}$)) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ax = abs(aimag(x(ix))) if (ax > tbig) then abig = abig + (ax*sbig)**2_${ik}$ notbig = .false. else if (ax < tsml) then if (notbig) asml = asml + (ax*ssml)**2_${ik}$ else amed = amed + ax**2_${ik}$ end if ix = ix + incx end do ! put the existing sum of squares into one of the accumulators if( sumsq > zero ) then ax = scl*sqrt( sumsq ) if (ax > tbig) then ! we assume scl >= sqrt( tiny*eps ) / sbig abig = abig + (scl*sbig)**2_${ik}$ * sumsq else if (ax < tsml) then ! we assume scl <= sqrt( huge ) / ssml if (notbig) asml = asml + (scl*ssml)**2_${ik}$ * sumsq else amed = amed + scl**2_${ik}$ * sumsq end if end if ! combine abig and amed or amed and asml if more than one ! accumulator was used. if (abig > zero) then ! combine abig and amed if abig > 0. if (amed > zero .or. ieee_is_nan(amed)) then abig = abig + (amed*sbig)*sbig end if scl = one / sbig sumsq = abig else if (asml > zero) then ! combine amed and asml if asml > 0. if (amed > zero .or. ieee_is_nan(amed)) then amed = sqrt(amed) asml = sqrt(asml) / ssml if (asml > amed) then ymin = amed ymax = asml else ymin = asml ymax = amed end if scl = one sumsq = ymax**2_${ik}$*( one + (ymin/ymax)**2_${ik}$ ) else scl = one / ssml sumsq = asml end if else ! otherwise all values are mid-range or zero scl = one sumsq = amed end if return end subroutine stdlib${ii}$_${ci}$lassq #:endif #:endfor pure module subroutine stdlib${ii}$_srscl( n, sa, sx, incx ) !! SRSCL multiplies an n-element real vector x by the real scalar 1/a. !! This is done without overflow or underflow as long as !! the final result x/a does not overflow or underflow. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(sp), intent(in) :: sa ! Array Arguments real(sp), intent(inout) :: sx(*) ! ===================================================================== ! Local Scalars logical(lk) :: done real(sp) :: bignum, cden, cden1, cnum, cnum1, mul, smlnum ! Intrinsic Functions ! Executable Statements ! quick return if possible if( n<=0 )return ! get machine parameters smlnum = stdlib${ii}$_slamch( 'S' ) bignum = one / smlnum call stdlib${ii}$_slabad( smlnum, bignum ) ! initialize the denominator to sa and the numerator to 1. cden = sa cnum = one 10 continue cden1 = cden*smlnum cnum1 = cnum / bignum if( abs( cden1 )>abs( cnum ) .and. cnum/=zero ) then ! pre-multiply x by smlnum if cden is large compared to cnum. mul = smlnum done = .false. cden = cden1 else if( abs( cnum1 )>abs( cden ) ) then ! pre-multiply x by bignum if cden is small compared to cnum. mul = bignum done = .false. cnum = cnum1 else ! multiply x by cnum / cden and return. mul = cnum / cden done = .true. end if ! scale the vector x by mul call stdlib${ii}$_sscal( n, mul, sx, incx ) if( .not.done )go to 10 return end subroutine stdlib${ii}$_srscl pure module subroutine stdlib${ii}$_drscl( n, sa, sx, incx ) !! DRSCL multiplies an n-element real vector x by the real scalar 1/a. !! This is done without overflow or underflow as long as !! the final result x/a does not overflow or underflow. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(dp), intent(in) :: sa ! Array Arguments real(dp), intent(inout) :: sx(*) ! ===================================================================== ! Local Scalars logical(lk) :: done real(dp) :: bignum, cden, cden1, cnum, cnum1, mul, smlnum ! Intrinsic Functions ! Executable Statements ! quick return if possible if( n<=0 )return ! get machine parameters smlnum = stdlib${ii}$_dlamch( 'S' ) bignum = one / smlnum call stdlib${ii}$_dlabad( smlnum, bignum ) ! initialize the denominator to sa and the numerator to 1. cden = sa cnum = one 10 continue cden1 = cden*smlnum cnum1 = cnum / bignum if( abs( cden1 )>abs( cnum ) .and. cnum/=zero ) then ! pre-multiply x by smlnum if cden is large compared to cnum. mul = smlnum done = .false. cden = cden1 else if( abs( cnum1 )>abs( cden ) ) then ! pre-multiply x by bignum if cden is small compared to cnum. mul = bignum done = .false. cnum = cnum1 else ! multiply x by cnum / cden and return. mul = cnum / cden done = .true. end if ! scale the vector x by mul call stdlib${ii}$_dscal( n, mul, sx, incx ) if( .not.done )go to 10 return end subroutine stdlib${ii}$_drscl #:for rk,rt,ri in REAL_KINDS_TYPES #:if not rk in ["sp","dp"] pure module subroutine stdlib${ii}$_${ri}$rscl( n, sa, sx, incx ) !! DRSCL: multiplies an n-element real vector x by the real scalar 1/a. !! This is done without overflow or underflow as long as !! the final result x/a does not overflow or underflow. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_${rk}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(${rk}$), intent(in) :: sa ! Array Arguments real(${rk}$), intent(inout) :: sx(*) ! ===================================================================== ! Local Scalars logical(lk) :: done real(${rk}$) :: bignum, cden, cden1, cnum, cnum1, mul, smlnum ! Intrinsic Functions ! Executable Statements ! quick return if possible if( n<=0 )return ! get machine parameters smlnum = stdlib${ii}$_${ri}$lamch( 'S' ) bignum = one / smlnum call stdlib${ii}$_${ri}$labad( smlnum, bignum ) ! initialize the denominator to sa and the numerator to 1. cden = sa cnum = one 10 continue cden1 = cden*smlnum cnum1 = cnum / bignum if( abs( cden1 )>abs( cnum ) .and. cnum/=zero ) then ! pre-multiply x by smlnum if cden is large compared to cnum. mul = smlnum done = .false. cden = cden1 else if( abs( cnum1 )>abs( cden ) ) then ! pre-multiply x by bignum if cden is small compared to cnum. mul = bignum done = .false. cnum = cnum1 else ! multiply x by cnum / cden and return. mul = cnum / cden done = .true. end if ! scale the vector x by mul call stdlib${ii}$_${ri}$scal( n, mul, sx, incx ) if( .not.done )go to 10 return end subroutine stdlib${ii}$_${ri}$rscl #:endif #:endfor pure module subroutine stdlib${ii}$_csrscl( n, sa, sx, incx ) !! CSRSCL multiplies an n-element complex vector x by the real scalar !! 1/a. This is done without overflow or underflow as long as !! the final result x/a does not overflow or underflow. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(sp), intent(in) :: sa ! Array Arguments complex(sp), intent(inout) :: sx(*) ! ===================================================================== ! Local Scalars logical(lk) :: done real(sp) :: bignum, cden, cden1, cnum, cnum1, mul, smlnum ! Intrinsic Functions ! Executable Statements ! quick return if possible if( n<=0 )return ! get machine parameters smlnum = stdlib${ii}$_slamch( 'S' ) bignum = one / smlnum call stdlib${ii}$_slabad( smlnum, bignum ) ! initialize the denominator to sa and the numerator to 1. cden = sa cnum = one 10 continue cden1 = cden*smlnum cnum1 = cnum / bignum if( abs( cden1 )>abs( cnum ) .and. cnum/=zero ) then ! pre-multiply x by smlnum if cden is large compared to cnum. mul = smlnum done = .false. cden = cden1 else if( abs( cnum1 )>abs( cden ) ) then ! pre-multiply x by bignum if cden is small compared to cnum. mul = bignum done = .false. cnum = cnum1 else ! multiply x by cnum / cden and return. mul = cnum / cden done = .true. end if ! scale the vector x by mul call stdlib${ii}$_csscal( n, mul, sx, incx ) if( .not.done )go to 10 return end subroutine stdlib${ii}$_csrscl pure module subroutine stdlib${ii}$_zdrscl( n, sa, sx, incx ) !! ZDRSCL multiplies an n-element complex vector x by the real scalar !! 1/a. This is done without overflow or underflow as long as !! the final result x/a does not overflow or underflow. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(dp), intent(in) :: sa ! Array Arguments complex(dp), intent(inout) :: sx(*) ! ===================================================================== ! Local Scalars logical(lk) :: done real(dp) :: bignum, cden, cden1, cnum, cnum1, mul, smlnum ! Intrinsic Functions ! Executable Statements ! quick return if possible if( n<=0 )return ! get machine parameters smlnum = stdlib${ii}$_dlamch( 'S' ) bignum = one / smlnum call stdlib${ii}$_dlabad( smlnum, bignum ) ! initialize the denominator to sa and the numerator to 1. cden = sa cnum = one 10 continue cden1 = cden*smlnum cnum1 = cnum / bignum if( abs( cden1 )>abs( cnum ) .and. cnum/=zero ) then ! pre-multiply x by smlnum if cden is large compared to cnum. mul = smlnum done = .false. cden = cden1 else if( abs( cnum1 )>abs( cden ) ) then ! pre-multiply x by bignum if cden is small compared to cnum. mul = bignum done = .false. cnum = cnum1 else ! multiply x by cnum / cden and return. mul = cnum / cden done = .true. end if ! scale the vector x by mul call stdlib${ii}$_zdscal( n, mul, sx, incx ) if( .not.done )go to 10 return end subroutine stdlib${ii}$_zdrscl #:for ck,ct,ci in CMPLX_KINDS_TYPES #:if not ck in ["sp","dp"] pure module subroutine stdlib${ii}$_${ci}$drscl( n, sa, sx, incx ) !! ZDRSCL: multiplies an n-element complex vector x by the real scalar !! 1/a. This is done without overflow or underflow as long as !! the final result x/a does not overflow or underflow. ! -- lapack auxiliary routine -- ! -- lapack is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- use stdlib_blas_constants_${ck}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone ! Scalar Arguments integer(${ik}$), intent(in) :: incx, n real(${ck}$), intent(in) :: sa ! Array Arguments complex(${ck}$), intent(inout) :: sx(*) ! ===================================================================== ! Local Scalars logical(lk) :: done real(${ck}$) :: bignum, cden, cden1, cnum, cnum1, mul, smlnum ! Intrinsic Functions ! Executable Statements ! quick return if possible if( n<=0 )return ! get machine parameters smlnum = stdlib${ii}$_${c2ri(ci)}$lamch( 'S' ) bignum = one / smlnum call stdlib${ii}$_${c2ri(ci)}$labad( smlnum, bignum ) ! initialize the denominator to sa and the numerator to 1. cden = sa cnum = one 10 continue cden1 = cden*smlnum cnum1 = cnum / bignum if( abs( cden1 )>abs( cnum ) .and. cnum/=zero ) then ! pre-multiply x by smlnum if cden is large compared to cnum. mul = smlnum done = .false. cden = cden1 else if( abs( cnum1 )>abs( cden ) ) then ! pre-multiply x by bignum if cden is small compared to cnum. mul = bignum done = .false. cnum = cnum1 else ! multiply x by cnum / cden and return. mul = cnum / cden done = .true. end if ! scale the vector x by mul call stdlib${ii}$_${ci}$dscal( n, mul, sx, incx ) if( .not.done )go to 10 return end subroutine stdlib${ii}$_${ci}$drscl #:endif #:endfor #:endfor end submodule stdlib_lapack_blas_like_l1