stdlib_lapack_blas_like_base.fypp Source File


Source Code

#:include "common.fypp" 
submodule(stdlib_lapack_base) stdlib_lapack_blas_like_base
  implicit none


  contains
#:for ik,it,ii in LINALG_INT_KINDS_TYPES

     pure module subroutine stdlib${ii}$_slaset( uplo, m, n, alpha, beta, a, lda )
     !! SLASET initializes an m-by-n matrix A to BETA on the diagonal and
     !! ALPHA on the offdiagonals.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, m, n
           real(sp), intent(in) :: alpha, beta
           ! Array Arguments 
           real(sp), intent(out) :: a(lda,*)
       ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              ! set the strictly upper triangular or trapezoidal part of the
              ! array to alpha.
              do j = 2, n
                 do i = 1, min( j-1, m )
                    a( i, j ) = alpha
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              ! set the strictly lower triangular or trapezoidal part of the
              ! array to alpha.
              do j = 1, min( m, n )
                 do i = j + 1, m
                    a( i, j ) = alpha
                 end do
              end do
           else
              ! set the leading m-by-n submatrix to alpha.
              do j = 1, n
                 do i = 1, m
                    a( i, j ) = alpha
                 end do
              end do
           end if
           ! set the first min(m,n) diagonal elements to beta.
           do i = 1, min( m, n )
              a( i, i ) = beta
           end do
           return
     end subroutine stdlib${ii}$_slaset

     pure module subroutine stdlib${ii}$_dlaset( uplo, m, n, alpha, beta, a, lda )
     !! DLASET initializes an m-by-n matrix A to BETA on the diagonal and
     !! ALPHA on the offdiagonals.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, m, n
           real(dp), intent(in) :: alpha, beta
           ! Array Arguments 
           real(dp), intent(out) :: a(lda,*)
       ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              ! set the strictly upper triangular or trapezoidal part of the
              ! array to alpha.
              do j = 2, n
                 do i = 1, min( j-1, m )
                    a( i, j ) = alpha
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              ! set the strictly lower triangular or trapezoidal part of the
              ! array to alpha.
              do j = 1, min( m, n )
                 do i = j + 1, m
                    a( i, j ) = alpha
                 end do
              end do
           else
              ! set the leading m-by-n submatrix to alpha.
              do j = 1, n
                 do i = 1, m
                    a( i, j ) = alpha
                 end do
              end do
           end if
           ! set the first min(m,n) diagonal elements to beta.
           do i = 1, min( m, n )
              a( i, i ) = beta
           end do
           return
     end subroutine stdlib${ii}$_dlaset

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ri}$laset( uplo, m, n, alpha, beta, a, lda )
     !! DLASET: initializes an m-by-n matrix A to BETA on the diagonal and
     !! ALPHA on the offdiagonals.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, m, n
           real(${rk}$), intent(in) :: alpha, beta
           ! Array Arguments 
           real(${rk}$), intent(out) :: a(lda,*)
       ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              ! set the strictly upper triangular or trapezoidal part of the
              ! array to alpha.
              do j = 2, n
                 do i = 1, min( j-1, m )
                    a( i, j ) = alpha
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              ! set the strictly lower triangular or trapezoidal part of the
              ! array to alpha.
              do j = 1, min( m, n )
                 do i = j + 1, m
                    a( i, j ) = alpha
                 end do
              end do
           else
              ! set the leading m-by-n submatrix to alpha.
              do j = 1, n
                 do i = 1, m
                    a( i, j ) = alpha
                 end do
              end do
           end if
           ! set the first min(m,n) diagonal elements to beta.
           do i = 1, min( m, n )
              a( i, i ) = beta
           end do
           return
     end subroutine stdlib${ii}$_${ri}$laset

#:endif
#:endfor

     pure module subroutine stdlib${ii}$_claset( uplo, m, n, alpha, beta, a, lda )
     !! CLASET initializes a 2-D array A to BETA on the diagonal and
     !! ALPHA on the offdiagonals.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, m, n
           complex(sp), intent(in) :: alpha, beta
           ! Array Arguments 
           complex(sp), intent(out) :: a(lda,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              ! set the diagonal to beta and the strictly upper triangular
              ! part of the array to alpha.
              do j = 2, n
                 do i = 1, min( j-1, m )
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( n, m )
                 a( i, i ) = beta
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              ! set the diagonal to beta and the strictly lower triangular
              ! part of the array to alpha.
              do j = 1, min( m, n )
                 do i = j + 1, m
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( n, m )
                 a( i, i ) = beta
              end do
           else
              ! set the array to beta on the diagonal and alpha on the
              ! offdiagonal.
              do j = 1, n
                 do i = 1, m
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( m, n )
                 a( i, i ) = beta
              end do
           end if
           return
     end subroutine stdlib${ii}$_claset

     pure module subroutine stdlib${ii}$_zlaset( uplo, m, n, alpha, beta, a, lda )
     !! ZLASET initializes a 2-D array A to BETA on the diagonal and
     !! ALPHA on the offdiagonals.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, m, n
           complex(dp), intent(in) :: alpha, beta
           ! Array Arguments 
           complex(dp), intent(out) :: a(lda,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              ! set the diagonal to beta and the strictly upper triangular
              ! part of the array to alpha.
              do j = 2, n
                 do i = 1, min( j-1, m )
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( n, m )
                 a( i, i ) = beta
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              ! set the diagonal to beta and the strictly lower triangular
              ! part of the array to alpha.
              do j = 1, min( m, n )
                 do i = j + 1, m
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( n, m )
                 a( i, i ) = beta
              end do
           else
              ! set the array to beta on the diagonal and alpha on the
              ! offdiagonal.
              do j = 1, n
                 do i = 1, m
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( m, n )
                 a( i, i ) = beta
              end do
           end if
           return
     end subroutine stdlib${ii}$_zlaset

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ci}$laset( uplo, m, n, alpha, beta, a, lda )
     !! ZLASET: initializes a 2-D array A to BETA on the diagonal and
     !! ALPHA on the offdiagonals.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, m, n
           complex(${ck}$), intent(in) :: alpha, beta
           ! Array Arguments 
           complex(${ck}$), intent(out) :: a(lda,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              ! set the diagonal to beta and the strictly upper triangular
              ! part of the array to alpha.
              do j = 2, n
                 do i = 1, min( j-1, m )
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( n, m )
                 a( i, i ) = beta
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              ! set the diagonal to beta and the strictly lower triangular
              ! part of the array to alpha.
              do j = 1, min( m, n )
                 do i = j + 1, m
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( n, m )
                 a( i, i ) = beta
              end do
           else
              ! set the array to beta on the diagonal and alpha on the
              ! offdiagonal.
              do j = 1, n
                 do i = 1, m
                    a( i, j ) = alpha
                 end do
              end do
              do i = 1, min( m, n )
                 a( i, i ) = beta
              end do
           end if
           return
     end subroutine stdlib${ii}$_${ci}$laset

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_slarnv( idist, iseed, n, x )
     !! SLARNV returns a vector of n random real numbers from a uniform or
     !! normal distribution.
        ! -- 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: one, two
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: idist, n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           real(sp), intent(out) :: x(*)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           real(sp), parameter :: twopi = 6.28318530717958647692528676655900576839e+0_sp
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, il, il2, iv
           ! Local Arrays 
           real(sp) :: u(lv)
           ! Intrinsic Functions 
           ! Executable Statements 
           do 40 iv = 1, n, lv / 2
              il = min( lv / 2_${ik}$, n-iv+1 )
              if( idist==3_${ik}$ ) then
                 il2 = 2_${ik}$*il
              else
                 il2 = il
              end if
              ! call stdlib${ii}$_slaruv to generate il2 numbers from a uniform (0,1)
              ! distribution (il2 <= lv)
              call stdlib${ii}$_slaruv( iseed, il2, u )
              if( idist==1_${ik}$ ) then
                 ! copy generated numbers
                 do i = 1, il
                    x( iv+i-1 ) = u( i )
                 end do
              else if( idist==2_${ik}$ ) then
                 ! convert generated numbers to uniform (-1,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = two*u( i ) - one
                 end do
              else if( idist==3_${ik}$ ) then
                 ! convert generated numbers to normal (0,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( -two*log( u( 2_${ik}$*i-1 ) ) )*cos( twopi*u( 2_${ik}$*i ) )
                 end do
              end if
              40 continue
           return
     end subroutine stdlib${ii}$_slarnv

     pure module subroutine stdlib${ii}$_dlarnv( idist, iseed, n, x )
     !! DLARNV returns a vector of n random real numbers from a uniform or
     !! normal distribution.
        ! -- 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: one, two
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: idist, n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           real(dp), intent(out) :: x(*)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           real(dp), parameter :: twopi = 6.28318530717958647692528676655900576839e+0_dp
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, il, il2, iv
           ! Local Arrays 
           real(dp) :: u(lv)
           ! Intrinsic Functions 
           ! Executable Statements 
           do 40 iv = 1, n, lv / 2
              il = min( lv / 2_${ik}$, n-iv+1 )
              if( idist==3_${ik}$ ) then
                 il2 = 2_${ik}$*il
              else
                 il2 = il
              end if
              ! call stdlib${ii}$_dlaruv to generate il2 numbers from a uniform (0,1)
              ! distribution (il2 <= lv)
              call stdlib${ii}$_dlaruv( iseed, il2, u )
              if( idist==1_${ik}$ ) then
                 ! copy generated numbers
                 do i = 1, il
                    x( iv+i-1 ) = u( i )
                 end do
              else if( idist==2_${ik}$ ) then
                 ! convert generated numbers to uniform (-1,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = two*u( i ) - one
                 end do
              else if( idist==3_${ik}$ ) then
                 ! convert generated numbers to normal (0,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( -two*log( u( 2_${ik}$*i-1 ) ) )*cos( twopi*u( 2_${ik}$*i ) )
                 end do
              end if
              40 continue
           return
     end subroutine stdlib${ii}$_dlarnv

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ri}$larnv( idist, iseed, n, x )
     !! DLARNV: returns a vector of n random real numbers from a uniform or
     !! normal distribution.
        ! -- 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: one, two
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: idist, n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           real(${rk}$), intent(out) :: x(*)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           real(${rk}$), parameter :: twopi = 6.28318530717958647692528676655900576839e+0_${rk}$
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, il, il2, iv
           ! Local Arrays 
           real(${rk}$) :: u(lv)
           ! Intrinsic Functions 
           ! Executable Statements 
           do 40 iv = 1, n, lv / 2
              il = min( lv / 2_${ik}$, n-iv+1 )
              if( idist==3_${ik}$ ) then
                 il2 = 2_${ik}$*il
              else
                 il2 = il
              end if
              ! call stdlib${ii}$_${ri}$laruv to generate il2 numbers from a uniform (0,1)
              ! distribution (il2 <= lv)
              call stdlib${ii}$_${ri}$laruv( iseed, il2, u )
              if( idist==1_${ik}$ ) then
                 ! copy generated numbers
                 do i = 1, il
                    x( iv+i-1 ) = u( i )
                 end do
              else if( idist==2_${ik}$ ) then
                 ! convert generated numbers to uniform (-1,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = two*u( i ) - one
                 end do
              else if( idist==3_${ik}$ ) then
                 ! convert generated numbers to normal (0,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( -two*log( u( 2_${ik}$*i-1 ) ) )*cos( twopi*u( 2_${ik}$*i ) )
                 end do
              end if
              40 continue
           return
     end subroutine stdlib${ii}$_${ri}$larnv

#:endif
#:endfor

     pure module subroutine stdlib${ii}$_clarnv( idist, iseed, n, x )
     !! CLARNV returns a vector of n random complex numbers from a uniform or
     !! normal distribution.
        ! -- 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, one, two
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: idist, n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           complex(sp), intent(out) :: x(*)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           real(sp), parameter :: twopi = 6.28318530717958647692528676655900576839e+0_sp
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, il, iv
           ! Local Arrays 
           real(sp) :: u(lv)
           ! Intrinsic Functions 
           ! Executable Statements 
           do 60 iv = 1, n, lv / 2
              il = min( lv / 2_${ik}$, n-iv+1 )
              ! call stdlib${ii}$_slaruv to generate 2*il realnumbers from a uniform (0,1,KIND=sp)
              ! distribution (2*il <= lv)
              call stdlib${ii}$_slaruv( iseed, 2_${ik}$*il, u )
              if( idist==1_${ik}$ ) then
                 ! copy generated numbers
                 do i = 1, il
                    x( iv+i-1 ) = cmplx( u( 2_${ik}$*i-1 ), u( 2_${ik}$*i ),KIND=sp)
                 end do
              else if( idist==2_${ik}$ ) then
                 ! convert generated numbers to uniform (-1,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = cmplx( two*u( 2_${ik}$*i-1 )-one,two*u( 2_${ik}$*i )-one,KIND=sp)
                 end do
              else if( idist==3_${ik}$ ) then
                 ! convert generated numbers to normal (0,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( -two*log( u( 2_${ik}$*i-1 ) ) )*exp( cmplx( zero, twopi*u( 2_${ik}$*i ),&
                              KIND=sp) )
                 end do
              else if( idist==4_${ik}$ ) then
                 ! convert generated numbers to complex numbers uniformly
                 ! distributed on the unit disk
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( u( 2_${ik}$*i-1 ) )*exp( cmplx( zero, twopi*u( 2_${ik}$*i ),KIND=sp) )
                              
                 end do
              else if( idist==5_${ik}$ ) then
                 ! convert generated numbers to complex numbers uniformly
                 ! distributed on the unit circle
                 do i = 1, il
                    x( iv+i-1 ) = exp( cmplx( zero, twopi*u( 2_${ik}$*i ),KIND=sp) )
                 end do
              end if
              60 continue
           return
     end subroutine stdlib${ii}$_clarnv

     pure module subroutine stdlib${ii}$_zlarnv( idist, iseed, n, x )
     !! ZLARNV returns a vector of n random complex numbers from a uniform or
     !! normal distribution.
        ! -- 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, one, two
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: idist, n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           complex(dp), intent(out) :: x(*)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           real(dp), parameter :: twopi = 6.28318530717958647692528676655900576839e+0_dp
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, il, iv
           ! Local Arrays 
           real(dp) :: u(lv)
           ! Intrinsic Functions 
           ! Executable Statements 
           do 60 iv = 1, n, lv / 2
              il = min( lv / 2_${ik}$, n-iv+1 )
              ! call stdlib${ii}$_dlaruv to generate 2*il realnumbers from a uniform (0,1,KIND=dp)
              ! distribution (2*il <= lv)
              call stdlib${ii}$_dlaruv( iseed, 2_${ik}$*il, u )
              if( idist==1_${ik}$ ) then
                 ! copy generated numbers
                 do i = 1, il
                    x( iv+i-1 ) = cmplx( u( 2_${ik}$*i-1 ), u( 2_${ik}$*i ),KIND=dp)
                 end do
              else if( idist==2_${ik}$ ) then
                 ! convert generated numbers to uniform (-1,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = cmplx( two*u( 2_${ik}$*i-1 )-one,two*u( 2_${ik}$*i )-one,KIND=dp)
                 end do
              else if( idist==3_${ik}$ ) then
                 ! convert generated numbers to normal (0,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( -two*log( u( 2_${ik}$*i-1 ) ) )*exp( cmplx( zero, twopi*u( 2_${ik}$*i ),&
                              KIND=dp) )
                 end do
              else if( idist==4_${ik}$ ) then
                 ! convert generated numbers to complex numbers uniformly
                 ! distributed on the unit disk
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( u( 2_${ik}$*i-1 ) )*exp( cmplx( zero, twopi*u( 2_${ik}$*i ),KIND=dp) )
                              
                 end do
              else if( idist==5_${ik}$ ) then
                 ! convert generated numbers to complex numbers uniformly
                 ! distributed on the unit circle
                 do i = 1, il
                    x( iv+i-1 ) = exp( cmplx( zero, twopi*u( 2_${ik}$*i ),KIND=dp) )
                 end do
              end if
              60 continue
           return
     end subroutine stdlib${ii}$_zlarnv

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ci}$larnv( idist, iseed, n, x )
     !! ZLARNV: returns a vector of n random complex numbers from a uniform or
     !! normal distribution.
        ! -- 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, one, two
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: idist, n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           complex(${ck}$), intent(out) :: x(*)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           real(${ck}$), parameter :: twopi = 6.28318530717958647692528676655900576839e+0_${ck}$
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, il, iv
           ! Local Arrays 
           real(${ck}$) :: u(lv)
           ! Intrinsic Functions 
           ! Executable Statements 
           do 60 iv = 1, n, lv / 2
              il = min( lv / 2_${ik}$, n-iv+1 )
              ! call stdlib${ii}$_${c2ri(ci)}$laruv to generate 2*il realnumbers from a uniform (0,1,KIND=${ck}$)
              ! distribution (2*il <= lv)
              call stdlib${ii}$_${c2ri(ci)}$laruv( iseed, 2_${ik}$*il, u )
              if( idist==1_${ik}$ ) then
                 ! copy generated numbers
                 do i = 1, il
                    x( iv+i-1 ) = cmplx( u( 2_${ik}$*i-1 ), u( 2_${ik}$*i ),KIND=${ck}$)
                 end do
              else if( idist==2_${ik}$ ) then
                 ! convert generated numbers to uniform (-1,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = cmplx( two*u( 2_${ik}$*i-1 )-one,two*u( 2_${ik}$*i )-one,KIND=${ck}$)
                 end do
              else if( idist==3_${ik}$ ) then
                 ! convert generated numbers to normal (0,1) distribution
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( -two*log( u( 2_${ik}$*i-1 ) ) )*exp( cmplx( zero, twopi*u( 2_${ik}$*i ),&
                              KIND=${ck}$) )
                 end do
              else if( idist==4_${ik}$ ) then
                 ! convert generated numbers to complex numbers uniformly
                 ! distributed on the unit disk
                 do i = 1, il
                    x( iv+i-1 ) = sqrt( u( 2_${ik}$*i-1 ) )*exp( cmplx( zero, twopi*u( 2_${ik}$*i ),KIND=${ck}$) )
                              
                 end do
              else if( idist==5_${ik}$ ) then
                 ! convert generated numbers to complex numbers uniformly
                 ! distributed on the unit circle
                 do i = 1, il
                    x( iv+i-1 ) = exp( cmplx( zero, twopi*u( 2_${ik}$*i ),KIND=${ck}$) )
                 end do
              end if
              60 continue
           return
     end subroutine stdlib${ii}$_${ci}$larnv

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_slaruv( iseed, n, x )
     !! SLARUV returns a vector of n random real numbers from a uniform (0,1)
     !! distribution (n <= 128).
     !! This is an auxiliary routine called by SLARNV and CLARNV.
        ! -- 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: one
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           real(sp), intent(out) :: x(n)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           integer(${ik}$), parameter :: ipw2 = 4096_${ik}$
           real(sp), parameter :: r = one/ipw2
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, i1, i2, i3, i4, it1, it2, it3, it4
           ! Local Arrays 
           integer(${ik}$) :: mm(lv,4_${ik}$)
           ! Intrinsic Functions 
           ! Data Statements 
           mm(1_${ik}$,1_${ik}$:4_${ik}$)=[494_${ik}$,322_${ik}$,2508_${ik}$,2549_${ik}$]
           mm(2_${ik}$,1_${ik}$:4_${ik}$)=[2637_${ik}$,789_${ik}$,3754_${ik}$,1145_${ik}$]
           mm(3_${ik}$,1_${ik}$:4_${ik}$)=[255_${ik}$,1440_${ik}$,1766_${ik}$,2253_${ik}$]
           mm(4_${ik}$,1_${ik}$:4_${ik}$)=[2008_${ik}$,752_${ik}$,3572_${ik}$,305_${ik}$]
           mm(5_${ik}$,1_${ik}$:4_${ik}$)=[1253_${ik}$,2859_${ik}$,2893_${ik}$,3301_${ik}$]
           mm(6_${ik}$,1_${ik}$:4_${ik}$)=[3344_${ik}$,123_${ik}$,307_${ik}$,1065_${ik}$]
           mm(7_${ik}$,1_${ik}$:4_${ik}$)=[4084_${ik}$,1848_${ik}$,1297_${ik}$,3133_${ik}$]
           mm(8_${ik}$,1_${ik}$:4_${ik}$)=[1739_${ik}$,643_${ik}$,3966_${ik}$,2913_${ik}$]
           mm(9_${ik}$,1_${ik}$:4_${ik}$)=[3143_${ik}$,2405_${ik}$,758_${ik}$,3285_${ik}$]
           mm(10_${ik}$,1_${ik}$:4_${ik}$)=[3468_${ik}$,2638_${ik}$,2598_${ik}$,1241_${ik}$]
           mm(11_${ik}$,1_${ik}$:4_${ik}$)=[688_${ik}$,2344_${ik}$,3406_${ik}$,1197_${ik}$]
           mm(12_${ik}$,1_${ik}$:4_${ik}$)=[1657_${ik}$,46_${ik}$,2922_${ik}$,3729_${ik}$]
           mm(13_${ik}$,1_${ik}$:4_${ik}$)=[1238_${ik}$,3814_${ik}$,1038_${ik}$,2501_${ik}$]
           mm(14_${ik}$,1_${ik}$:4_${ik}$)=[3166_${ik}$,913_${ik}$,2934_${ik}$,1673_${ik}$]
           mm(15_${ik}$,1_${ik}$:4_${ik}$)=[1292_${ik}$,3649_${ik}$,2091_${ik}$,541_${ik}$]
           mm(16_${ik}$,1_${ik}$:4_${ik}$)=[3422_${ik}$,339_${ik}$,2451_${ik}$,2753_${ik}$]
           mm(17_${ik}$,1_${ik}$:4_${ik}$)=[1270_${ik}$,3808_${ik}$,1580_${ik}$,949_${ik}$]
           mm(18_${ik}$,1_${ik}$:4_${ik}$)=[2016_${ik}$,822_${ik}$,1958_${ik}$,2361_${ik}$]
           mm(19_${ik}$,1_${ik}$:4_${ik}$)=[154_${ik}$,2832_${ik}$,2055_${ik}$,1165_${ik}$]
           mm(20_${ik}$,1_${ik}$:4_${ik}$)=[2862_${ik}$,3078_${ik}$,1507_${ik}$,4081_${ik}$]
           mm(21_${ik}$,1_${ik}$:4_${ik}$)=[697_${ik}$,3633_${ik}$,1078_${ik}$,2725_${ik}$]
           mm(22_${ik}$,1_${ik}$:4_${ik}$)=[1706_${ik}$,2970_${ik}$,3273_${ik}$,3305_${ik}$]
           mm(23_${ik}$,1_${ik}$:4_${ik}$)=[491_${ik}$,637_${ik}$,17_${ik}$,3069_${ik}$]
           mm(24_${ik}$,1_${ik}$:4_${ik}$)=[931_${ik}$,2249_${ik}$,854_${ik}$,3617_${ik}$]
           mm(25_${ik}$,1_${ik}$:4_${ik}$)=[1444_${ik}$,2081_${ik}$,2916_${ik}$,3733_${ik}$]
           mm(26_${ik}$,1_${ik}$:4_${ik}$)=[444_${ik}$,4019_${ik}$,3971_${ik}$,409_${ik}$]
           mm(27_${ik}$,1_${ik}$:4_${ik}$)=[3577_${ik}$,1478_${ik}$,2889_${ik}$,2157_${ik}$]
           mm(28_${ik}$,1_${ik}$:4_${ik}$)=[3944_${ik}$,242_${ik}$,3831_${ik}$,1361_${ik}$]
           mm(29_${ik}$,1_${ik}$:4_${ik}$)=[2184_${ik}$,481_${ik}$,2621_${ik}$,3973_${ik}$]
           mm(30_${ik}$,1_${ik}$:4_${ik}$)=[1661_${ik}$,2075_${ik}$,1541_${ik}$,1865_${ik}$]
           mm(31_${ik}$,1_${ik}$:4_${ik}$)=[3482_${ik}$,4058_${ik}$,893_${ik}$,2525_${ik}$]
           mm(32_${ik}$,1_${ik}$:4_${ik}$)=[657_${ik}$,622_${ik}$,736_${ik}$,1409_${ik}$]
           mm(33_${ik}$,1_${ik}$:4_${ik}$)=[3023_${ik}$,3376_${ik}$,3992_${ik}$,3445_${ik}$]
           mm(34_${ik}$,1_${ik}$:4_${ik}$)=[3618_${ik}$,812_${ik}$,787_${ik}$,3577_${ik}$]
           mm(35_${ik}$,1_${ik}$:4_${ik}$)=[1267_${ik}$,234_${ik}$,2125_${ik}$,77_${ik}$]
           mm(36_${ik}$,1_${ik}$:4_${ik}$)=[1828_${ik}$,641_${ik}$,2364_${ik}$,3761_${ik}$]
           mm(37_${ik}$,1_${ik}$:4_${ik}$)=[164_${ik}$,4005_${ik}$,2460_${ik}$,2149_${ik}$]
           mm(38_${ik}$,1_${ik}$:4_${ik}$)=[3798_${ik}$,1122_${ik}$,257_${ik}$,1449_${ik}$]
           mm(39_${ik}$,1_${ik}$:4_${ik}$)=[3087_${ik}$,3135_${ik}$,1574_${ik}$,3005_${ik}$]
           mm(40_${ik}$,1_${ik}$:4_${ik}$)=[2400_${ik}$,2640_${ik}$,3912_${ik}$,225_${ik}$]
           mm(41_${ik}$,1_${ik}$:4_${ik}$)=[2870_${ik}$,2302_${ik}$,1216_${ik}$,85_${ik}$]
           mm(42_${ik}$,1_${ik}$:4_${ik}$)=[3876_${ik}$,40_${ik}$,3248_${ik}$,3673_${ik}$]
           mm(43_${ik}$,1_${ik}$:4_${ik}$)=[1905_${ik}$,1832_${ik}$,3401_${ik}$,3117_${ik}$]
           mm(44_${ik}$,1_${ik}$:4_${ik}$)=[1593_${ik}$,2247_${ik}$,2124_${ik}$,3089_${ik}$]
           mm(45_${ik}$,1_${ik}$:4_${ik}$)=[1797_${ik}$,2034_${ik}$,2762_${ik}$,1349_${ik}$]
           mm(46_${ik}$,1_${ik}$:4_${ik}$)=[1234_${ik}$,2637_${ik}$,149_${ik}$,2057_${ik}$]
           mm(47_${ik}$,1_${ik}$:4_${ik}$)=[3460_${ik}$,1287_${ik}$,2245_${ik}$,413_${ik}$]
           mm(48_${ik}$,1_${ik}$:4_${ik}$)=[328_${ik}$,1691_${ik}$,166_${ik}$,65_${ik}$]
           mm(49_${ik}$,1_${ik}$:4_${ik}$)=[2861_${ik}$,496_${ik}$,466_${ik}$,1845_${ik}$]
           mm(50_${ik}$,1_${ik}$:4_${ik}$)=[1950_${ik}$,1597_${ik}$,4018_${ik}$,697_${ik}$]
           mm(51_${ik}$,1_${ik}$:4_${ik}$)=[617_${ik}$,2394_${ik}$,1399_${ik}$,3085_${ik}$]
           mm(52_${ik}$,1_${ik}$:4_${ik}$)=[2070_${ik}$,2584_${ik}$,190_${ik}$,3441_${ik}$]
           mm(53_${ik}$,1_${ik}$:4_${ik}$)=[3331_${ik}$,1843_${ik}$,2879_${ik}$,1573_${ik}$]
           mm(54_${ik}$,1_${ik}$:4_${ik}$)=[769_${ik}$,336_${ik}$,153_${ik}$,3689_${ik}$]
           mm(55_${ik}$,1_${ik}$:4_${ik}$)=[1558_${ik}$,1472_${ik}$,2320_${ik}$,2941_${ik}$]
           mm(56_${ik}$,1_${ik}$:4_${ik}$)=[2412_${ik}$,2407_${ik}$,18_${ik}$,929_${ik}$]
           mm(57_${ik}$,1_${ik}$:4_${ik}$)=[2800_${ik}$,433_${ik}$,712_${ik}$,533_${ik}$]
           mm(58_${ik}$,1_${ik}$:4_${ik}$)=[189_${ik}$,2096_${ik}$,2159_${ik}$,2841_${ik}$]
           mm(59_${ik}$,1_${ik}$:4_${ik}$)=[287_${ik}$,1761_${ik}$,2318_${ik}$,4077_${ik}$]
           mm(60_${ik}$,1_${ik}$:4_${ik}$)=[2045_${ik}$,2810_${ik}$,2091_${ik}$,721_${ik}$]
           mm(61_${ik}$,1_${ik}$:4_${ik}$)=[1227_${ik}$,566_${ik}$,3443_${ik}$,2821_${ik}$]
           mm(62_${ik}$,1_${ik}$:4_${ik}$)=[2838_${ik}$,442_${ik}$,1510_${ik}$,2249_${ik}$]
           mm(63_${ik}$,1_${ik}$:4_${ik}$)=[209_${ik}$,41_${ik}$,449_${ik}$,2397_${ik}$]
           mm(64_${ik}$,1_${ik}$:4_${ik}$)=[2770_${ik}$,1238_${ik}$,1956_${ik}$,2817_${ik}$]
           mm(65_${ik}$,1_${ik}$:4_${ik}$)=[3654_${ik}$,1086_${ik}$,2201_${ik}$,245_${ik}$]
           mm(66_${ik}$,1_${ik}$:4_${ik}$)=[3993_${ik}$,603_${ik}$,3137_${ik}$,1913_${ik}$]
           mm(67_${ik}$,1_${ik}$:4_${ik}$)=[192_${ik}$,840_${ik}$,3399_${ik}$,1997_${ik}$]
           mm(68_${ik}$,1_${ik}$:4_${ik}$)=[2253_${ik}$,3168_${ik}$,1321_${ik}$,3121_${ik}$]
           mm(69_${ik}$,1_${ik}$:4_${ik}$)=[3491_${ik}$,1499_${ik}$,2271_${ik}$,997_${ik}$]
           mm(70_${ik}$,1_${ik}$:4_${ik}$)=[2889_${ik}$,1084_${ik}$,3667_${ik}$,1833_${ik}$]
           mm(71_${ik}$,1_${ik}$:4_${ik}$)=[2857_${ik}$,3438_${ik}$,2703_${ik}$,2877_${ik}$]
           mm(72_${ik}$,1_${ik}$:4_${ik}$)=[2094_${ik}$,2408_${ik}$,629_${ik}$,1633_${ik}$]
           mm(73_${ik}$,1_${ik}$:4_${ik}$)=[1818_${ik}$,1589_${ik}$,2365_${ik}$,981_${ik}$]
           mm(74_${ik}$,1_${ik}$:4_${ik}$)=[688_${ik}$,2391_${ik}$,2431_${ik}$,2009_${ik}$]
           mm(75_${ik}$,1_${ik}$:4_${ik}$)=[1407_${ik}$,288_${ik}$,1113_${ik}$,941_${ik}$]
           mm(76_${ik}$,1_${ik}$:4_${ik}$)=[634_${ik}$,26_${ik}$,3922_${ik}$,2449_${ik}$]
           mm(77_${ik}$,1_${ik}$:4_${ik}$)=[3231_${ik}$,512_${ik}$,2554_${ik}$,197_${ik}$]
           mm(78_${ik}$,1_${ik}$:4_${ik}$)=[815_${ik}$,1456_${ik}$,184_${ik}$,2441_${ik}$]
           mm(79_${ik}$,1_${ik}$:4_${ik}$)=[3524_${ik}$,171_${ik}$,2099_${ik}$,285_${ik}$]
           mm(80_${ik}$,1_${ik}$:4_${ik}$)=[1914_${ik}$,1677_${ik}$,3228_${ik}$,1473_${ik}$]
           mm(81_${ik}$,1_${ik}$:4_${ik}$)=[516_${ik}$,2657_${ik}$,4012_${ik}$,2741_${ik}$]
           mm(82_${ik}$,1_${ik}$:4_${ik}$)=[164_${ik}$,2270_${ik}$,1921_${ik}$,3129_${ik}$]
           mm(83_${ik}$,1_${ik}$:4_${ik}$)=[303_${ik}$,2587_${ik}$,3452_${ik}$,909_${ik}$]
           mm(84_${ik}$,1_${ik}$:4_${ik}$)=[2144_${ik}$,2961_${ik}$,3901_${ik}$,2801_${ik}$]
           mm(85_${ik}$,1_${ik}$:4_${ik}$)=[3480_${ik}$,1970_${ik}$,572_${ik}$,421_${ik}$]
           mm(86_${ik}$,1_${ik}$:4_${ik}$)=[119_${ik}$,1817_${ik}$,3309_${ik}$,4073_${ik}$]
           mm(87_${ik}$,1_${ik}$:4_${ik}$)=[3357_${ik}$,676_${ik}$,3171_${ik}$,2813_${ik}$]
           mm(88_${ik}$,1_${ik}$:4_${ik}$)=[837_${ik}$,1410_${ik}$,817_${ik}$,2337_${ik}$]
           mm(89_${ik}$,1_${ik}$:4_${ik}$)=[2826_${ik}$,3723_${ik}$,3039_${ik}$,1429_${ik}$]
           mm(90_${ik}$,1_${ik}$:4_${ik}$)=[2332_${ik}$,2803_${ik}$,1696_${ik}$,1177_${ik}$]
           mm(91_${ik}$,1_${ik}$:4_${ik}$)=[2089_${ik}$,3185_${ik}$,1256_${ik}$,1901_${ik}$]
           mm(92_${ik}$,1_${ik}$:4_${ik}$)=[3780_${ik}$,184_${ik}$,3715_${ik}$,81_${ik}$]
           mm(93_${ik}$,1_${ik}$:4_${ik}$)=[1700_${ik}$,663_${ik}$,2077_${ik}$,1669_${ik}$]
           mm(94_${ik}$,1_${ik}$:4_${ik}$)=[3712_${ik}$,499_${ik}$,3019_${ik}$,2633_${ik}$]
           mm(95_${ik}$,1_${ik}$:4_${ik}$)=[150_${ik}$,3784_${ik}$,1497_${ik}$,2269_${ik}$]
           mm(96_${ik}$,1_${ik}$:4_${ik}$)=[2000_${ik}$,1631_${ik}$,1101_${ik}$,129_${ik}$]
           mm(97_${ik}$,1_${ik}$:4_${ik}$)=[3375_${ik}$,1925_${ik}$,717_${ik}$,1141_${ik}$]
           mm(98_${ik}$,1_${ik}$:4_${ik}$)=[1621_${ik}$,3912_${ik}$,51_${ik}$,249_${ik}$]
           mm(99_${ik}$,1_${ik}$:4_${ik}$)=[3090_${ik}$,1398_${ik}$,981_${ik}$,3917_${ik}$]
           mm(100_${ik}$,1_${ik}$:4_${ik}$)=[3765_${ik}$,1349_${ik}$,1978_${ik}$,2481_${ik}$]
           mm(101_${ik}$,1_${ik}$:4_${ik}$)=[1149_${ik}$,1441_${ik}$,1813_${ik}$,3941_${ik}$]
           mm(102_${ik}$,1_${ik}$:4_${ik}$)=[3146_${ik}$,2224_${ik}$,3881_${ik}$,2217_${ik}$]
           mm(103_${ik}$,1_${ik}$:4_${ik}$)=[33_${ik}$,2411_${ik}$,76_${ik}$,2749_${ik}$]
           mm(104_${ik}$,1_${ik}$:4_${ik}$)=[3082_${ik}$,1907_${ik}$,3846_${ik}$,3041_${ik}$]
           mm(105_${ik}$,1_${ik}$:4_${ik}$)=[2741_${ik}$,3192_${ik}$,3694_${ik}$,1877_${ik}$]
           mm(106_${ik}$,1_${ik}$:4_${ik}$)=[359_${ik}$,2786_${ik}$,1682_${ik}$,345_${ik}$]
           mm(107_${ik}$,1_${ik}$:4_${ik}$)=[3316_${ik}$,382_${ik}$,124_${ik}$,2861_${ik}$]
           mm(108_${ik}$,1_${ik}$:4_${ik}$)=[1749_${ik}$,37_${ik}$,1660_${ik}$,1809_${ik}$]
           mm(109_${ik}$,1_${ik}$:4_${ik}$)=[185_${ik}$,759_${ik}$,3997_${ik}$,3141_${ik}$]
           mm(110_${ik}$,1_${ik}$:4_${ik}$)=[2784_${ik}$,2948_${ik}$,479_${ik}$,2825_${ik}$]
           mm(111_${ik}$,1_${ik}$:4_${ik}$)=[2202_${ik}$,1862_${ik}$,1141_${ik}$,157_${ik}$]
           mm(112_${ik}$,1_${ik}$:4_${ik}$)=[2199_${ik}$,3802_${ik}$,886_${ik}$,2881_${ik}$]
           mm(113_${ik}$,1_${ik}$:4_${ik}$)=[1364_${ik}$,2423_${ik}$,3514_${ik}$,3637_${ik}$]
           mm(114_${ik}$,1_${ik}$:4_${ik}$)=[1244_${ik}$,2051_${ik}$,1301_${ik}$,1465_${ik}$]
           mm(115_${ik}$,1_${ik}$:4_${ik}$)=[2020_${ik}$,2295_${ik}$,3604_${ik}$,2829_${ik}$]
           mm(116_${ik}$,1_${ik}$:4_${ik}$)=[3160_${ik}$,1332_${ik}$,1888_${ik}$,2161_${ik}$]
           mm(117_${ik}$,1_${ik}$:4_${ik}$)=[2785_${ik}$,1832_${ik}$,1836_${ik}$,3365_${ik}$]
           mm(118_${ik}$,1_${ik}$:4_${ik}$)=[2772_${ik}$,2405_${ik}$,1990_${ik}$,361_${ik}$]
           mm(119_${ik}$,1_${ik}$:4_${ik}$)=[1217_${ik}$,3638_${ik}$,2058_${ik}$,2685_${ik}$]
           mm(120_${ik}$,1_${ik}$:4_${ik}$)=[1822_${ik}$,3661_${ik}$,692_${ik}$,3745_${ik}$]
           mm(121_${ik}$,1_${ik}$:4_${ik}$)=[1245_${ik}$,327_${ik}$,1194_${ik}$,2325_${ik}$]
           mm(122_${ik}$,1_${ik}$:4_${ik}$)=[2252_${ik}$,3660_${ik}$,20_${ik}$,3609_${ik}$]
           mm(123_${ik}$,1_${ik}$:4_${ik}$)=[3904_${ik}$,716_${ik}$,3285_${ik}$,3821_${ik}$]
           mm(124_${ik}$,1_${ik}$:4_${ik}$)=[2774_${ik}$,1842_${ik}$,2046_${ik}$,3537_${ik}$]
           mm(125_${ik}$,1_${ik}$:4_${ik}$)=[997_${ik}$,3987_${ik}$,2107_${ik}$,517_${ik}$]
           mm(126_${ik}$,1_${ik}$:4_${ik}$)=[2573_${ik}$,1368_${ik}$,3508_${ik}$,3017_${ik}$]
           mm(127_${ik}$,1_${ik}$:4_${ik}$)=[1148_${ik}$,1848_${ik}$,3525_${ik}$,2141_${ik}$]
           mm(128_${ik}$,1_${ik}$:4_${ik}$)=[545_${ik}$,2366_${ik}$,3801_${ik}$,1537_${ik}$]
           ! Executable Statements 
           i1 = iseed( 1_${ik}$ )
           i2 = iseed( 2_${ik}$ )
           i3 = iseed( 3_${ik}$ )
           i4 = iseed( 4_${ik}$ )
           loop_10: do i = 1, min( n, lv )
           20 continue
              ! multiply the seed by i-th power of the multiplier modulo 2**48
              it4 = i4*mm( i, 4_${ik}$ )
              it3 = it4 / ipw2
              it4 = it4 - ipw2*it3
              it3 = it3 + i3*mm( i, 4_${ik}$ ) + i4*mm( i, 3_${ik}$ )
              it2 = it3 / ipw2
              it3 = it3 - ipw2*it2
              it2 = it2 + i2*mm( i, 4_${ik}$ ) + i3*mm( i, 3_${ik}$ ) + i4*mm( i, 2_${ik}$ )
              it1 = it2 / ipw2
              it2 = it2 - ipw2*it1
              it1 = it1 + i1*mm( i, 4_${ik}$ ) + i2*mm( i, 3_${ik}$ ) + i3*mm( i, 2_${ik}$ ) +i4*mm( i, 1_${ik}$ )
              it1 = mod( it1, ipw2 )
              ! convert 48-bit integer to a realnumber in the interval (0,1,KIND=sp)
              x( i ) = r*( real( it1,KIND=sp)+r*( real( it2,KIND=sp)+r*( real( it3,KIND=sp)+&
                        r*real( it4,KIND=sp) ) ) )
              if (x( i )==one) then
                 ! if a real number has n bits of precision, and the first
                 ! n bits of the 48-bit integer above happen to be all 1 (which
                 ! will occur about once every 2**n calls), then x( i ) will
                 ! be rounded to exactly one. in ieee single precision arithmetic,
                 ! this will happen relatively often since n = 24.
                 ! since x( i ) is not supposed to return exactly 0.0_sp or one,
                 ! the statistically correct thing to do in this situation is
                 ! simply to iterate again.
                 ! n.b. the case x( i ) = 0.0_sp should not be possible.
                 i1 = i1 + 2_${ik}$
                 i2 = i2 + 2_${ik}$
                 i3 = i3 + 2_${ik}$
                 i4 = i4 + 2_${ik}$
                 goto 20
              end if
           end do loop_10
           ! return final value of seed
           iseed( 1_${ik}$ ) = it1
           iseed( 2_${ik}$ ) = it2
           iseed( 3_${ik}$ ) = it3
           iseed( 4_${ik}$ ) = it4
           return
     end subroutine stdlib${ii}$_slaruv

     pure module subroutine stdlib${ii}$_dlaruv( iseed, n, x )
     !! DLARUV returns a vector of n random real numbers from a uniform (0,1)
     !! distribution (n <= 128).
     !! This is an auxiliary routine called by DLARNV and ZLARNV.
        ! -- 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: one
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           real(dp), intent(out) :: x(n)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           integer(${ik}$), parameter :: ipw2 = 4096_${ik}$
           real(dp), parameter :: r = one/ipw2
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, i1, i2, i3, i4, it1, it2, it3, it4
           ! Local Arrays 
           integer(${ik}$) :: mm(lv,4_${ik}$)
           ! Intrinsic Functions 
           ! Data Statements 
           mm(1_${ik}$,1_${ik}$:4_${ik}$)=[494_${ik}$,322_${ik}$,2508_${ik}$,2549_${ik}$]
           mm(2_${ik}$,1_${ik}$:4_${ik}$)=[2637_${ik}$,789_${ik}$,3754_${ik}$,1145_${ik}$]
           mm(3_${ik}$,1_${ik}$:4_${ik}$)=[255_${ik}$,1440_${ik}$,1766_${ik}$,2253_${ik}$]
           mm(4_${ik}$,1_${ik}$:4_${ik}$)=[2008_${ik}$,752_${ik}$,3572_${ik}$,305_${ik}$]
           mm(5_${ik}$,1_${ik}$:4_${ik}$)=[1253_${ik}$,2859_${ik}$,2893_${ik}$,3301_${ik}$]
           mm(6_${ik}$,1_${ik}$:4_${ik}$)=[3344_${ik}$,123_${ik}$,307_${ik}$,1065_${ik}$]
           mm(7_${ik}$,1_${ik}$:4_${ik}$)=[4084_${ik}$,1848_${ik}$,1297_${ik}$,3133_${ik}$]
           mm(8_${ik}$,1_${ik}$:4_${ik}$)=[1739_${ik}$,643_${ik}$,3966_${ik}$,2913_${ik}$]
           mm(9_${ik}$,1_${ik}$:4_${ik}$)=[3143_${ik}$,2405_${ik}$,758_${ik}$,3285_${ik}$]
           mm(10_${ik}$,1_${ik}$:4_${ik}$)=[3468_${ik}$,2638_${ik}$,2598_${ik}$,1241_${ik}$]
           mm(11_${ik}$,1_${ik}$:4_${ik}$)=[688_${ik}$,2344_${ik}$,3406_${ik}$,1197_${ik}$]
           mm(12_${ik}$,1_${ik}$:4_${ik}$)=[1657_${ik}$,46_${ik}$,2922_${ik}$,3729_${ik}$]
           mm(13_${ik}$,1_${ik}$:4_${ik}$)=[1238_${ik}$,3814_${ik}$,1038_${ik}$,2501_${ik}$]
           mm(14_${ik}$,1_${ik}$:4_${ik}$)=[3166_${ik}$,913_${ik}$,2934_${ik}$,1673_${ik}$]
           mm(15_${ik}$,1_${ik}$:4_${ik}$)=[1292_${ik}$,3649_${ik}$,2091_${ik}$,541_${ik}$]
           mm(16_${ik}$,1_${ik}$:4_${ik}$)=[3422_${ik}$,339_${ik}$,2451_${ik}$,2753_${ik}$]
           mm(17_${ik}$,1_${ik}$:4_${ik}$)=[1270_${ik}$,3808_${ik}$,1580_${ik}$,949_${ik}$]
           mm(18_${ik}$,1_${ik}$:4_${ik}$)=[2016_${ik}$,822_${ik}$,1958_${ik}$,2361_${ik}$]
           mm(19_${ik}$,1_${ik}$:4_${ik}$)=[154_${ik}$,2832_${ik}$,2055_${ik}$,1165_${ik}$]
           mm(20_${ik}$,1_${ik}$:4_${ik}$)=[2862_${ik}$,3078_${ik}$,1507_${ik}$,4081_${ik}$]
           mm(21_${ik}$,1_${ik}$:4_${ik}$)=[697_${ik}$,3633_${ik}$,1078_${ik}$,2725_${ik}$]
           mm(22_${ik}$,1_${ik}$:4_${ik}$)=[1706_${ik}$,2970_${ik}$,3273_${ik}$,3305_${ik}$]
           mm(23_${ik}$,1_${ik}$:4_${ik}$)=[491_${ik}$,637_${ik}$,17_${ik}$,3069_${ik}$]
           mm(24_${ik}$,1_${ik}$:4_${ik}$)=[931_${ik}$,2249_${ik}$,854_${ik}$,3617_${ik}$]
           mm(25_${ik}$,1_${ik}$:4_${ik}$)=[1444_${ik}$,2081_${ik}$,2916_${ik}$,3733_${ik}$]
           mm(26_${ik}$,1_${ik}$:4_${ik}$)=[444_${ik}$,4019_${ik}$,3971_${ik}$,409_${ik}$]
           mm(27_${ik}$,1_${ik}$:4_${ik}$)=[3577_${ik}$,1478_${ik}$,2889_${ik}$,2157_${ik}$]
           mm(28_${ik}$,1_${ik}$:4_${ik}$)=[3944_${ik}$,242_${ik}$,3831_${ik}$,1361_${ik}$]
           mm(29_${ik}$,1_${ik}$:4_${ik}$)=[2184_${ik}$,481_${ik}$,2621_${ik}$,3973_${ik}$]
           mm(30_${ik}$,1_${ik}$:4_${ik}$)=[1661_${ik}$,2075_${ik}$,1541_${ik}$,1865_${ik}$]
           mm(31_${ik}$,1_${ik}$:4_${ik}$)=[3482_${ik}$,4058_${ik}$,893_${ik}$,2525_${ik}$]
           mm(32_${ik}$,1_${ik}$:4_${ik}$)=[657_${ik}$,622_${ik}$,736_${ik}$,1409_${ik}$]
           mm(33_${ik}$,1_${ik}$:4_${ik}$)=[3023_${ik}$,3376_${ik}$,3992_${ik}$,3445_${ik}$]
           mm(34_${ik}$,1_${ik}$:4_${ik}$)=[3618_${ik}$,812_${ik}$,787_${ik}$,3577_${ik}$]
           mm(35_${ik}$,1_${ik}$:4_${ik}$)=[1267_${ik}$,234_${ik}$,2125_${ik}$,77_${ik}$]
           mm(36_${ik}$,1_${ik}$:4_${ik}$)=[1828_${ik}$,641_${ik}$,2364_${ik}$,3761_${ik}$]
           mm(37_${ik}$,1_${ik}$:4_${ik}$)=[164_${ik}$,4005_${ik}$,2460_${ik}$,2149_${ik}$]
           mm(38_${ik}$,1_${ik}$:4_${ik}$)=[3798_${ik}$,1122_${ik}$,257_${ik}$,1449_${ik}$]
           mm(39_${ik}$,1_${ik}$:4_${ik}$)=[3087_${ik}$,3135_${ik}$,1574_${ik}$,3005_${ik}$]
           mm(40_${ik}$,1_${ik}$:4_${ik}$)=[2400_${ik}$,2640_${ik}$,3912_${ik}$,225_${ik}$]
           mm(41_${ik}$,1_${ik}$:4_${ik}$)=[2870_${ik}$,2302_${ik}$,1216_${ik}$,85_${ik}$]
           mm(42_${ik}$,1_${ik}$:4_${ik}$)=[3876_${ik}$,40_${ik}$,3248_${ik}$,3673_${ik}$]
           mm(43_${ik}$,1_${ik}$:4_${ik}$)=[1905_${ik}$,1832_${ik}$,3401_${ik}$,3117_${ik}$]
           mm(44_${ik}$,1_${ik}$:4_${ik}$)=[1593_${ik}$,2247_${ik}$,2124_${ik}$,3089_${ik}$]
           mm(45_${ik}$,1_${ik}$:4_${ik}$)=[1797_${ik}$,2034_${ik}$,2762_${ik}$,1349_${ik}$]
           mm(46_${ik}$,1_${ik}$:4_${ik}$)=[1234_${ik}$,2637_${ik}$,149_${ik}$,2057_${ik}$]
           mm(47_${ik}$,1_${ik}$:4_${ik}$)=[3460_${ik}$,1287_${ik}$,2245_${ik}$,413_${ik}$]
           mm(48_${ik}$,1_${ik}$:4_${ik}$)=[328_${ik}$,1691_${ik}$,166_${ik}$,65_${ik}$]
           mm(49_${ik}$,1_${ik}$:4_${ik}$)=[2861_${ik}$,496_${ik}$,466_${ik}$,1845_${ik}$]
           mm(50_${ik}$,1_${ik}$:4_${ik}$)=[1950_${ik}$,1597_${ik}$,4018_${ik}$,697_${ik}$]
           mm(51_${ik}$,1_${ik}$:4_${ik}$)=[617_${ik}$,2394_${ik}$,1399_${ik}$,3085_${ik}$]
           mm(52_${ik}$,1_${ik}$:4_${ik}$)=[2070_${ik}$,2584_${ik}$,190_${ik}$,3441_${ik}$]
           mm(53_${ik}$,1_${ik}$:4_${ik}$)=[3331_${ik}$,1843_${ik}$,2879_${ik}$,1573_${ik}$]
           mm(54_${ik}$,1_${ik}$:4_${ik}$)=[769_${ik}$,336_${ik}$,153_${ik}$,3689_${ik}$]
           mm(55_${ik}$,1_${ik}$:4_${ik}$)=[1558_${ik}$,1472_${ik}$,2320_${ik}$,2941_${ik}$]
           mm(56_${ik}$,1_${ik}$:4_${ik}$)=[2412_${ik}$,2407_${ik}$,18_${ik}$,929_${ik}$]
           mm(57_${ik}$,1_${ik}$:4_${ik}$)=[2800_${ik}$,433_${ik}$,712_${ik}$,533_${ik}$]
           mm(58_${ik}$,1_${ik}$:4_${ik}$)=[189_${ik}$,2096_${ik}$,2159_${ik}$,2841_${ik}$]
           mm(59_${ik}$,1_${ik}$:4_${ik}$)=[287_${ik}$,1761_${ik}$,2318_${ik}$,4077_${ik}$]
           mm(60_${ik}$,1_${ik}$:4_${ik}$)=[2045_${ik}$,2810_${ik}$,2091_${ik}$,721_${ik}$]
           mm(61_${ik}$,1_${ik}$:4_${ik}$)=[1227_${ik}$,566_${ik}$,3443_${ik}$,2821_${ik}$]
           mm(62_${ik}$,1_${ik}$:4_${ik}$)=[2838_${ik}$,442_${ik}$,1510_${ik}$,2249_${ik}$]
           mm(63_${ik}$,1_${ik}$:4_${ik}$)=[209_${ik}$,41_${ik}$,449_${ik}$,2397_${ik}$]
           mm(64_${ik}$,1_${ik}$:4_${ik}$)=[2770_${ik}$,1238_${ik}$,1956_${ik}$,2817_${ik}$]
           mm(65_${ik}$,1_${ik}$:4_${ik}$)=[3654_${ik}$,1086_${ik}$,2201_${ik}$,245_${ik}$]
           mm(66_${ik}$,1_${ik}$:4_${ik}$)=[3993_${ik}$,603_${ik}$,3137_${ik}$,1913_${ik}$]
           mm(67_${ik}$,1_${ik}$:4_${ik}$)=[192_${ik}$,840_${ik}$,3399_${ik}$,1997_${ik}$]
           mm(68_${ik}$,1_${ik}$:4_${ik}$)=[2253_${ik}$,3168_${ik}$,1321_${ik}$,3121_${ik}$]
           mm(69_${ik}$,1_${ik}$:4_${ik}$)=[3491_${ik}$,1499_${ik}$,2271_${ik}$,997_${ik}$]
           mm(70_${ik}$,1_${ik}$:4_${ik}$)=[2889_${ik}$,1084_${ik}$,3667_${ik}$,1833_${ik}$]
           mm(71_${ik}$,1_${ik}$:4_${ik}$)=[2857_${ik}$,3438_${ik}$,2703_${ik}$,2877_${ik}$]
           mm(72_${ik}$,1_${ik}$:4_${ik}$)=[2094_${ik}$,2408_${ik}$,629_${ik}$,1633_${ik}$]
           mm(73_${ik}$,1_${ik}$:4_${ik}$)=[1818_${ik}$,1589_${ik}$,2365_${ik}$,981_${ik}$]
           mm(74_${ik}$,1_${ik}$:4_${ik}$)=[688_${ik}$,2391_${ik}$,2431_${ik}$,2009_${ik}$]
           mm(75_${ik}$,1_${ik}$:4_${ik}$)=[1407_${ik}$,288_${ik}$,1113_${ik}$,941_${ik}$]
           mm(76_${ik}$,1_${ik}$:4_${ik}$)=[634_${ik}$,26_${ik}$,3922_${ik}$,2449_${ik}$]
           mm(77_${ik}$,1_${ik}$:4_${ik}$)=[3231_${ik}$,512_${ik}$,2554_${ik}$,197_${ik}$]
           mm(78_${ik}$,1_${ik}$:4_${ik}$)=[815_${ik}$,1456_${ik}$,184_${ik}$,2441_${ik}$]
           mm(79_${ik}$,1_${ik}$:4_${ik}$)=[3524_${ik}$,171_${ik}$,2099_${ik}$,285_${ik}$]
           mm(80_${ik}$,1_${ik}$:4_${ik}$)=[1914_${ik}$,1677_${ik}$,3228_${ik}$,1473_${ik}$]
           mm(81_${ik}$,1_${ik}$:4_${ik}$)=[516_${ik}$,2657_${ik}$,4012_${ik}$,2741_${ik}$]
           mm(82_${ik}$,1_${ik}$:4_${ik}$)=[164_${ik}$,2270_${ik}$,1921_${ik}$,3129_${ik}$]
           mm(83_${ik}$,1_${ik}$:4_${ik}$)=[303_${ik}$,2587_${ik}$,3452_${ik}$,909_${ik}$]
           mm(84_${ik}$,1_${ik}$:4_${ik}$)=[2144_${ik}$,2961_${ik}$,3901_${ik}$,2801_${ik}$]
           mm(85_${ik}$,1_${ik}$:4_${ik}$)=[3480_${ik}$,1970_${ik}$,572_${ik}$,421_${ik}$]
           mm(86_${ik}$,1_${ik}$:4_${ik}$)=[119_${ik}$,1817_${ik}$,3309_${ik}$,4073_${ik}$]
           mm(87_${ik}$,1_${ik}$:4_${ik}$)=[3357_${ik}$,676_${ik}$,3171_${ik}$,2813_${ik}$]
           mm(88_${ik}$,1_${ik}$:4_${ik}$)=[837_${ik}$,1410_${ik}$,817_${ik}$,2337_${ik}$]
           mm(89_${ik}$,1_${ik}$:4_${ik}$)=[2826_${ik}$,3723_${ik}$,3039_${ik}$,1429_${ik}$]
           mm(90_${ik}$,1_${ik}$:4_${ik}$)=[2332_${ik}$,2803_${ik}$,1696_${ik}$,1177_${ik}$]
           mm(91_${ik}$,1_${ik}$:4_${ik}$)=[2089_${ik}$,3185_${ik}$,1256_${ik}$,1901_${ik}$]
           mm(92_${ik}$,1_${ik}$:4_${ik}$)=[3780_${ik}$,184_${ik}$,3715_${ik}$,81_${ik}$]
           mm(93_${ik}$,1_${ik}$:4_${ik}$)=[1700_${ik}$,663_${ik}$,2077_${ik}$,1669_${ik}$]
           mm(94_${ik}$,1_${ik}$:4_${ik}$)=[3712_${ik}$,499_${ik}$,3019_${ik}$,2633_${ik}$]
           mm(95_${ik}$,1_${ik}$:4_${ik}$)=[150_${ik}$,3784_${ik}$,1497_${ik}$,2269_${ik}$]
           mm(96_${ik}$,1_${ik}$:4_${ik}$)=[2000_${ik}$,1631_${ik}$,1101_${ik}$,129_${ik}$]
           mm(97_${ik}$,1_${ik}$:4_${ik}$)=[3375_${ik}$,1925_${ik}$,717_${ik}$,1141_${ik}$]
           mm(98_${ik}$,1_${ik}$:4_${ik}$)=[1621_${ik}$,3912_${ik}$,51_${ik}$,249_${ik}$]
           mm(99_${ik}$,1_${ik}$:4_${ik}$)=[3090_${ik}$,1398_${ik}$,981_${ik}$,3917_${ik}$]
           mm(100_${ik}$,1_${ik}$:4_${ik}$)=[3765_${ik}$,1349_${ik}$,1978_${ik}$,2481_${ik}$]
           mm(101_${ik}$,1_${ik}$:4_${ik}$)=[1149_${ik}$,1441_${ik}$,1813_${ik}$,3941_${ik}$]
           mm(102_${ik}$,1_${ik}$:4_${ik}$)=[3146_${ik}$,2224_${ik}$,3881_${ik}$,2217_${ik}$]
           mm(103_${ik}$,1_${ik}$:4_${ik}$)=[33_${ik}$,2411_${ik}$,76_${ik}$,2749_${ik}$]
           mm(104_${ik}$,1_${ik}$:4_${ik}$)=[3082_${ik}$,1907_${ik}$,3846_${ik}$,3041_${ik}$]
           mm(105_${ik}$,1_${ik}$:4_${ik}$)=[2741_${ik}$,3192_${ik}$,3694_${ik}$,1877_${ik}$]
           mm(106_${ik}$,1_${ik}$:4_${ik}$)=[359_${ik}$,2786_${ik}$,1682_${ik}$,345_${ik}$]
           mm(107_${ik}$,1_${ik}$:4_${ik}$)=[3316_${ik}$,382_${ik}$,124_${ik}$,2861_${ik}$]
           mm(108_${ik}$,1_${ik}$:4_${ik}$)=[1749_${ik}$,37_${ik}$,1660_${ik}$,1809_${ik}$]
           mm(109_${ik}$,1_${ik}$:4_${ik}$)=[185_${ik}$,759_${ik}$,3997_${ik}$,3141_${ik}$]
           mm(110_${ik}$,1_${ik}$:4_${ik}$)=[2784_${ik}$,2948_${ik}$,479_${ik}$,2825_${ik}$]
           mm(111_${ik}$,1_${ik}$:4_${ik}$)=[2202_${ik}$,1862_${ik}$,1141_${ik}$,157_${ik}$]
           mm(112_${ik}$,1_${ik}$:4_${ik}$)=[2199_${ik}$,3802_${ik}$,886_${ik}$,2881_${ik}$]
           mm(113_${ik}$,1_${ik}$:4_${ik}$)=[1364_${ik}$,2423_${ik}$,3514_${ik}$,3637_${ik}$]
           mm(114_${ik}$,1_${ik}$:4_${ik}$)=[1244_${ik}$,2051_${ik}$,1301_${ik}$,1465_${ik}$]
           mm(115_${ik}$,1_${ik}$:4_${ik}$)=[2020_${ik}$,2295_${ik}$,3604_${ik}$,2829_${ik}$]
           mm(116_${ik}$,1_${ik}$:4_${ik}$)=[3160_${ik}$,1332_${ik}$,1888_${ik}$,2161_${ik}$]
           mm(117_${ik}$,1_${ik}$:4_${ik}$)=[2785_${ik}$,1832_${ik}$,1836_${ik}$,3365_${ik}$]
           mm(118_${ik}$,1_${ik}$:4_${ik}$)=[2772_${ik}$,2405_${ik}$,1990_${ik}$,361_${ik}$]
           mm(119_${ik}$,1_${ik}$:4_${ik}$)=[1217_${ik}$,3638_${ik}$,2058_${ik}$,2685_${ik}$]
           mm(120_${ik}$,1_${ik}$:4_${ik}$)=[1822_${ik}$,3661_${ik}$,692_${ik}$,3745_${ik}$]
           mm(121_${ik}$,1_${ik}$:4_${ik}$)=[1245_${ik}$,327_${ik}$,1194_${ik}$,2325_${ik}$]
           mm(122_${ik}$,1_${ik}$:4_${ik}$)=[2252_${ik}$,3660_${ik}$,20_${ik}$,3609_${ik}$]
           mm(123_${ik}$,1_${ik}$:4_${ik}$)=[3904_${ik}$,716_${ik}$,3285_${ik}$,3821_${ik}$]
           mm(124_${ik}$,1_${ik}$:4_${ik}$)=[2774_${ik}$,1842_${ik}$,2046_${ik}$,3537_${ik}$]
           mm(125_${ik}$,1_${ik}$:4_${ik}$)=[997_${ik}$,3987_${ik}$,2107_${ik}$,517_${ik}$]
           mm(126_${ik}$,1_${ik}$:4_${ik}$)=[2573_${ik}$,1368_${ik}$,3508_${ik}$,3017_${ik}$]
           mm(127_${ik}$,1_${ik}$:4_${ik}$)=[1148_${ik}$,1848_${ik}$,3525_${ik}$,2141_${ik}$]
           mm(128_${ik}$,1_${ik}$:4_${ik}$)=[545_${ik}$,2366_${ik}$,3801_${ik}$,1537_${ik}$]
           ! Executable Statements 
           i1 = iseed( 1_${ik}$ )
           i2 = iseed( 2_${ik}$ )
           i3 = iseed( 3_${ik}$ )
           i4 = iseed( 4_${ik}$ )
           loop_10: do i = 1, min( n, lv )
           20 continue
              ! multiply the seed by i-th power of the multiplier modulo 2**48
              it4 = i4*mm( i, 4_${ik}$ )
              it3 = it4 / ipw2
              it4 = it4 - ipw2*it3
              it3 = it3 + i3*mm( i, 4_${ik}$ ) + i4*mm( i, 3_${ik}$ )
              it2 = it3 / ipw2
              it3 = it3 - ipw2*it2
              it2 = it2 + i2*mm( i, 4_${ik}$ ) + i3*mm( i, 3_${ik}$ ) + i4*mm( i, 2_${ik}$ )
              it1 = it2 / ipw2
              it2 = it2 - ipw2*it1
              it1 = it1 + i1*mm( i, 4_${ik}$ ) + i2*mm( i, 3_${ik}$ ) + i3*mm( i, 2_${ik}$ ) +i4*mm( i, 1_${ik}$ )
              it1 = mod( it1, ipw2 )
              ! convert 48-bit integer to a realnumber in the interval (0,1,KIND=dp)
              x( i ) = r*( real( it1,KIND=dp)+r*( real( it2,KIND=dp)+r*( real( it3,KIND=dp)+&
                        r*real( it4,KIND=dp) ) ) )
              if (x( i )==one) then
                 ! if a real number has n bits of precision, and the first
                 ! n bits of the 48-bit integer above happen to be all 1 (which
                 ! will occur about once every 2**n calls), then x( i ) will
                 ! be rounded to exactly one.
                 ! since x( i ) is not supposed to return exactly 0.0_dp or one,
                 ! the statistically correct thing to do in this situation is
                 ! simply to iterate again.
                 ! n.b. the case x( i ) = 0.0_dp should not be possible.
                 i1 = i1 + 2_${ik}$
                 i2 = i2 + 2_${ik}$
                 i3 = i3 + 2_${ik}$
                 i4 = i4 + 2_${ik}$
                 goto 20
              end if
           end do loop_10
           ! return final value of seed
           iseed( 1_${ik}$ ) = it1
           iseed( 2_${ik}$ ) = it2
           iseed( 3_${ik}$ ) = it3
           iseed( 4_${ik}$ ) = it4
           return
     end subroutine stdlib${ii}$_dlaruv

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ri}$laruv( iseed, n, x )
     !! DLARUV: returns a vector of n random real numbers from a uniform (0,1)
     !! distribution (n <= 128).
     !! This is an auxiliary routine called by DLARNV and ZLARNV.
        ! -- 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: one
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: n
           ! Array Arguments 
           integer(${ik}$), intent(inout) :: iseed(4_${ik}$)
           real(${rk}$), intent(out) :: x(n)
        ! =====================================================================
           ! Parameters 
           integer(${ik}$), parameter :: lv = 128_${ik}$
           integer(${ik}$), parameter :: ipw2 = 4096_${ik}$
           real(${rk}$), parameter :: r = one/ipw2
           
           
           
           ! Local Scalars 
           integer(${ik}$) :: i, i1, i2, i3, i4, it1, it2, it3, it4, j
           ! Local Arrays 
           integer(${ik}$) :: mm(lv,4_${ik}$)
           ! Intrinsic Functions 
           ! Data Statements 
           mm(1_${ik}$,1_${ik}$:4_${ik}$)=[494_${ik}$,322_${ik}$,2508_${ik}$,2549_${ik}$]
           mm(2_${ik}$,1_${ik}$:4_${ik}$)=[2637_${ik}$,789_${ik}$,3754_${ik}$,1145_${ik}$]
           mm(3_${ik}$,1_${ik}$:4_${ik}$)=[255_${ik}$,1440_${ik}$,1766_${ik}$,2253_${ik}$]
           mm(4_${ik}$,1_${ik}$:4_${ik}$)=[2008_${ik}$,752_${ik}$,3572_${ik}$,305_${ik}$]
           mm(5_${ik}$,1_${ik}$:4_${ik}$)=[1253_${ik}$,2859_${ik}$,2893_${ik}$,3301_${ik}$]
           mm(6_${ik}$,1_${ik}$:4_${ik}$)=[3344_${ik}$,123_${ik}$,307_${ik}$,1065_${ik}$]
           mm(7_${ik}$,1_${ik}$:4_${ik}$)=[4084_${ik}$,1848_${ik}$,1297_${ik}$,3133_${ik}$]
           mm(8_${ik}$,1_${ik}$:4_${ik}$)=[1739_${ik}$,643_${ik}$,3966_${ik}$,2913_${ik}$]
           mm(9_${ik}$,1_${ik}$:4_${ik}$)=[3143_${ik}$,2405_${ik}$,758_${ik}$,3285_${ik}$]
           mm(10_${ik}$,1_${ik}$:4_${ik}$)=[3468_${ik}$,2638_${ik}$,2598_${ik}$,1241_${ik}$]
           mm(11_${ik}$,1_${ik}$:4_${ik}$)=[688_${ik}$,2344_${ik}$,3406_${ik}$,1197_${ik}$]
           mm(12_${ik}$,1_${ik}$:4_${ik}$)=[1657_${ik}$,46_${ik}$,2922_${ik}$,3729_${ik}$]
           mm(13_${ik}$,1_${ik}$:4_${ik}$)=[1238_${ik}$,3814_${ik}$,1038_${ik}$,2501_${ik}$]
           mm(14_${ik}$,1_${ik}$:4_${ik}$)=[3166_${ik}$,913_${ik}$,2934_${ik}$,1673_${ik}$]
           mm(15_${ik}$,1_${ik}$:4_${ik}$)=[1292_${ik}$,3649_${ik}$,2091_${ik}$,541_${ik}$]
           mm(16_${ik}$,1_${ik}$:4_${ik}$)=[3422_${ik}$,339_${ik}$,2451_${ik}$,2753_${ik}$]
           mm(17_${ik}$,1_${ik}$:4_${ik}$)=[1270_${ik}$,3808_${ik}$,1580_${ik}$,949_${ik}$]
           mm(18_${ik}$,1_${ik}$:4_${ik}$)=[2016_${ik}$,822_${ik}$,1958_${ik}$,2361_${ik}$]
           mm(19_${ik}$,1_${ik}$:4_${ik}$)=[154_${ik}$,2832_${ik}$,2055_${ik}$,1165_${ik}$]
           mm(20_${ik}$,1_${ik}$:4_${ik}$)=[2862_${ik}$,3078_${ik}$,1507_${ik}$,4081_${ik}$]
           mm(21_${ik}$,1_${ik}$:4_${ik}$)=[697_${ik}$,3633_${ik}$,1078_${ik}$,2725_${ik}$]
           mm(22_${ik}$,1_${ik}$:4_${ik}$)=[1706_${ik}$,2970_${ik}$,3273_${ik}$,3305_${ik}$]
           mm(23_${ik}$,1_${ik}$:4_${ik}$)=[491_${ik}$,637_${ik}$,17_${ik}$,3069_${ik}$]
           mm(24_${ik}$,1_${ik}$:4_${ik}$)=[931_${ik}$,2249_${ik}$,854_${ik}$,3617_${ik}$]
           mm(25_${ik}$,1_${ik}$:4_${ik}$)=[1444_${ik}$,2081_${ik}$,2916_${ik}$,3733_${ik}$]
           mm(26_${ik}$,1_${ik}$:4_${ik}$)=[444_${ik}$,4019_${ik}$,3971_${ik}$,409_${ik}$]
           mm(27_${ik}$,1_${ik}$:4_${ik}$)=[3577_${ik}$,1478_${ik}$,2889_${ik}$,2157_${ik}$]
           mm(28_${ik}$,1_${ik}$:4_${ik}$)=[3944_${ik}$,242_${ik}$,3831_${ik}$,1361_${ik}$]
           mm(29_${ik}$,1_${ik}$:4_${ik}$)=[2184_${ik}$,481_${ik}$,2621_${ik}$,3973_${ik}$]
           mm(30_${ik}$,1_${ik}$:4_${ik}$)=[1661_${ik}$,2075_${ik}$,1541_${ik}$,1865_${ik}$]
           mm(31_${ik}$,1_${ik}$:4_${ik}$)=[3482_${ik}$,4058_${ik}$,893_${ik}$,2525_${ik}$]
           mm(32_${ik}$,1_${ik}$:4_${ik}$)=[657_${ik}$,622_${ik}$,736_${ik}$,1409_${ik}$]
           mm(33_${ik}$,1_${ik}$:4_${ik}$)=[3023_${ik}$,3376_${ik}$,3992_${ik}$,3445_${ik}$]
           mm(34_${ik}$,1_${ik}$:4_${ik}$)=[3618_${ik}$,812_${ik}$,787_${ik}$,3577_${ik}$]
           mm(35_${ik}$,1_${ik}$:4_${ik}$)=[1267_${ik}$,234_${ik}$,2125_${ik}$,77_${ik}$]
           mm(36_${ik}$,1_${ik}$:4_${ik}$)=[1828_${ik}$,641_${ik}$,2364_${ik}$,3761_${ik}$]
           mm(37_${ik}$,1_${ik}$:4_${ik}$)=[164_${ik}$,4005_${ik}$,2460_${ik}$,2149_${ik}$]
           mm(38_${ik}$,1_${ik}$:4_${ik}$)=[3798_${ik}$,1122_${ik}$,257_${ik}$,1449_${ik}$]
           mm(39_${ik}$,1_${ik}$:4_${ik}$)=[3087_${ik}$,3135_${ik}$,1574_${ik}$,3005_${ik}$]
           mm(40_${ik}$,1_${ik}$:4_${ik}$)=[2400_${ik}$,2640_${ik}$,3912_${ik}$,225_${ik}$]
           mm(41_${ik}$,1_${ik}$:4_${ik}$)=[2870_${ik}$,2302_${ik}$,1216_${ik}$,85_${ik}$]
           mm(42_${ik}$,1_${ik}$:4_${ik}$)=[3876_${ik}$,40_${ik}$,3248_${ik}$,3673_${ik}$]
           mm(43_${ik}$,1_${ik}$:4_${ik}$)=[1905_${ik}$,1832_${ik}$,3401_${ik}$,3117_${ik}$]
           mm(44_${ik}$,1_${ik}$:4_${ik}$)=[1593_${ik}$,2247_${ik}$,2124_${ik}$,3089_${ik}$]
           mm(45_${ik}$,1_${ik}$:4_${ik}$)=[1797_${ik}$,2034_${ik}$,2762_${ik}$,1349_${ik}$]
           mm(46_${ik}$,1_${ik}$:4_${ik}$)=[1234_${ik}$,2637_${ik}$,149_${ik}$,2057_${ik}$]
           mm(47_${ik}$,1_${ik}$:4_${ik}$)=[3460_${ik}$,1287_${ik}$,2245_${ik}$,413_${ik}$]
           mm(48_${ik}$,1_${ik}$:4_${ik}$)=[328_${ik}$,1691_${ik}$,166_${ik}$,65_${ik}$]
           mm(49_${ik}$,1_${ik}$:4_${ik}$)=[2861_${ik}$,496_${ik}$,466_${ik}$,1845_${ik}$]
           mm(50_${ik}$,1_${ik}$:4_${ik}$)=[1950_${ik}$,1597_${ik}$,4018_${ik}$,697_${ik}$]
           mm(51_${ik}$,1_${ik}$:4_${ik}$)=[617_${ik}$,2394_${ik}$,1399_${ik}$,3085_${ik}$]
           mm(52_${ik}$,1_${ik}$:4_${ik}$)=[2070_${ik}$,2584_${ik}$,190_${ik}$,3441_${ik}$]
           mm(53_${ik}$,1_${ik}$:4_${ik}$)=[3331_${ik}$,1843_${ik}$,2879_${ik}$,1573_${ik}$]
           mm(54_${ik}$,1_${ik}$:4_${ik}$)=[769_${ik}$,336_${ik}$,153_${ik}$,3689_${ik}$]
           mm(55_${ik}$,1_${ik}$:4_${ik}$)=[1558_${ik}$,1472_${ik}$,2320_${ik}$,2941_${ik}$]
           mm(56_${ik}$,1_${ik}$:4_${ik}$)=[2412_${ik}$,2407_${ik}$,18_${ik}$,929_${ik}$]
           mm(57_${ik}$,1_${ik}$:4_${ik}$)=[2800_${ik}$,433_${ik}$,712_${ik}$,533_${ik}$]
           mm(58_${ik}$,1_${ik}$:4_${ik}$)=[189_${ik}$,2096_${ik}$,2159_${ik}$,2841_${ik}$]
           mm(59_${ik}$,1_${ik}$:4_${ik}$)=[287_${ik}$,1761_${ik}$,2318_${ik}$,4077_${ik}$]
           mm(60_${ik}$,1_${ik}$:4_${ik}$)=[2045_${ik}$,2810_${ik}$,2091_${ik}$,721_${ik}$]
           mm(61_${ik}$,1_${ik}$:4_${ik}$)=[1227_${ik}$,566_${ik}$,3443_${ik}$,2821_${ik}$]
           mm(62_${ik}$,1_${ik}$:4_${ik}$)=[2838_${ik}$,442_${ik}$,1510_${ik}$,2249_${ik}$]
           mm(63_${ik}$,1_${ik}$:4_${ik}$)=[209_${ik}$,41_${ik}$,449_${ik}$,2397_${ik}$]
           mm(64_${ik}$,1_${ik}$:4_${ik}$)=[2770_${ik}$,1238_${ik}$,1956_${ik}$,2817_${ik}$]
           mm(65_${ik}$,1_${ik}$:4_${ik}$)=[3654_${ik}$,1086_${ik}$,2201_${ik}$,245_${ik}$]
           mm(66_${ik}$,1_${ik}$:4_${ik}$)=[3993_${ik}$,603_${ik}$,3137_${ik}$,1913_${ik}$]
           mm(67_${ik}$,1_${ik}$:4_${ik}$)=[192_${ik}$,840_${ik}$,3399_${ik}$,1997_${ik}$]
           mm(68_${ik}$,1_${ik}$:4_${ik}$)=[2253_${ik}$,3168_${ik}$,1321_${ik}$,3121_${ik}$]
           mm(69_${ik}$,1_${ik}$:4_${ik}$)=[3491_${ik}$,1499_${ik}$,2271_${ik}$,997_${ik}$]
           mm(70_${ik}$,1_${ik}$:4_${ik}$)=[2889_${ik}$,1084_${ik}$,3667_${ik}$,1833_${ik}$]
           mm(71_${ik}$,1_${ik}$:4_${ik}$)=[2857_${ik}$,3438_${ik}$,2703_${ik}$,2877_${ik}$]
           mm(72_${ik}$,1_${ik}$:4_${ik}$)=[2094_${ik}$,2408_${ik}$,629_${ik}$,1633_${ik}$]
           mm(73_${ik}$,1_${ik}$:4_${ik}$)=[1818_${ik}$,1589_${ik}$,2365_${ik}$,981_${ik}$]
           mm(74_${ik}$,1_${ik}$:4_${ik}$)=[688_${ik}$,2391_${ik}$,2431_${ik}$,2009_${ik}$]
           mm(75_${ik}$,1_${ik}$:4_${ik}$)=[1407_${ik}$,288_${ik}$,1113_${ik}$,941_${ik}$]
           mm(76_${ik}$,1_${ik}$:4_${ik}$)=[634_${ik}$,26_${ik}$,3922_${ik}$,2449_${ik}$]
           mm(77_${ik}$,1_${ik}$:4_${ik}$)=[3231_${ik}$,512_${ik}$,2554_${ik}$,197_${ik}$]
           mm(78_${ik}$,1_${ik}$:4_${ik}$)=[815_${ik}$,1456_${ik}$,184_${ik}$,2441_${ik}$]
           mm(79_${ik}$,1_${ik}$:4_${ik}$)=[3524_${ik}$,171_${ik}$,2099_${ik}$,285_${ik}$]
           mm(80_${ik}$,1_${ik}$:4_${ik}$)=[1914_${ik}$,1677_${ik}$,3228_${ik}$,1473_${ik}$]
           mm(81_${ik}$,1_${ik}$:4_${ik}$)=[516_${ik}$,2657_${ik}$,4012_${ik}$,2741_${ik}$]
           mm(82_${ik}$,1_${ik}$:4_${ik}$)=[164_${ik}$,2270_${ik}$,1921_${ik}$,3129_${ik}$]
           mm(83_${ik}$,1_${ik}$:4_${ik}$)=[303_${ik}$,2587_${ik}$,3452_${ik}$,909_${ik}$]
           mm(84_${ik}$,1_${ik}$:4_${ik}$)=[2144_${ik}$,2961_${ik}$,3901_${ik}$,2801_${ik}$]
           mm(85_${ik}$,1_${ik}$:4_${ik}$)=[3480_${ik}$,1970_${ik}$,572_${ik}$,421_${ik}$]
           mm(86_${ik}$,1_${ik}$:4_${ik}$)=[119_${ik}$,1817_${ik}$,3309_${ik}$,4073_${ik}$]
           mm(87_${ik}$,1_${ik}$:4_${ik}$)=[3357_${ik}$,676_${ik}$,3171_${ik}$,2813_${ik}$]
           mm(88_${ik}$,1_${ik}$:4_${ik}$)=[837_${ik}$,1410_${ik}$,817_${ik}$,2337_${ik}$]
           mm(89_${ik}$,1_${ik}$:4_${ik}$)=[2826_${ik}$,3723_${ik}$,3039_${ik}$,1429_${ik}$]
           mm(90_${ik}$,1_${ik}$:4_${ik}$)=[2332_${ik}$,2803_${ik}$,1696_${ik}$,1177_${ik}$]
           mm(91_${ik}$,1_${ik}$:4_${ik}$)=[2089_${ik}$,3185_${ik}$,1256_${ik}$,1901_${ik}$]
           mm(92_${ik}$,1_${ik}$:4_${ik}$)=[3780_${ik}$,184_${ik}$,3715_${ik}$,81_${ik}$]
           mm(93_${ik}$,1_${ik}$:4_${ik}$)=[1700_${ik}$,663_${ik}$,2077_${ik}$,1669_${ik}$]
           mm(94_${ik}$,1_${ik}$:4_${ik}$)=[3712_${ik}$,499_${ik}$,3019_${ik}$,2633_${ik}$]
           mm(95_${ik}$,1_${ik}$:4_${ik}$)=[150_${ik}$,3784_${ik}$,1497_${ik}$,2269_${ik}$]
           mm(96_${ik}$,1_${ik}$:4_${ik}$)=[2000_${ik}$,1631_${ik}$,1101_${ik}$,129_${ik}$]
           mm(97_${ik}$,1_${ik}$:4_${ik}$)=[3375_${ik}$,1925_${ik}$,717_${ik}$,1141_${ik}$]
           mm(98_${ik}$,1_${ik}$:4_${ik}$)=[1621_${ik}$,3912_${ik}$,51_${ik}$,249_${ik}$]
           mm(99_${ik}$,1_${ik}$:4_${ik}$)=[3090_${ik}$,1398_${ik}$,981_${ik}$,3917_${ik}$]
           mm(100_${ik}$,1_${ik}$:4_${ik}$)=[3765_${ik}$,1349_${ik}$,1978_${ik}$,2481_${ik}$]
           mm(101_${ik}$,1_${ik}$:4_${ik}$)=[1149_${ik}$,1441_${ik}$,1813_${ik}$,3941_${ik}$]
           mm(102_${ik}$,1_${ik}$:4_${ik}$)=[3146_${ik}$,2224_${ik}$,3881_${ik}$,2217_${ik}$]
           mm(103_${ik}$,1_${ik}$:4_${ik}$)=[33_${ik}$,2411_${ik}$,76_${ik}$,2749_${ik}$]
           mm(104_${ik}$,1_${ik}$:4_${ik}$)=[3082_${ik}$,1907_${ik}$,3846_${ik}$,3041_${ik}$]
           mm(105_${ik}$,1_${ik}$:4_${ik}$)=[2741_${ik}$,3192_${ik}$,3694_${ik}$,1877_${ik}$]
           mm(106_${ik}$,1_${ik}$:4_${ik}$)=[359_${ik}$,2786_${ik}$,1682_${ik}$,345_${ik}$]
           mm(107_${ik}$,1_${ik}$:4_${ik}$)=[3316_${ik}$,382_${ik}$,124_${ik}$,2861_${ik}$]
           mm(108_${ik}$,1_${ik}$:4_${ik}$)=[1749_${ik}$,37_${ik}$,1660_${ik}$,1809_${ik}$]
           mm(109_${ik}$,1_${ik}$:4_${ik}$)=[185_${ik}$,759_${ik}$,3997_${ik}$,3141_${ik}$]
           mm(110_${ik}$,1_${ik}$:4_${ik}$)=[2784_${ik}$,2948_${ik}$,479_${ik}$,2825_${ik}$]
           mm(111_${ik}$,1_${ik}$:4_${ik}$)=[2202_${ik}$,1862_${ik}$,1141_${ik}$,157_${ik}$]
           mm(112_${ik}$,1_${ik}$:4_${ik}$)=[2199_${ik}$,3802_${ik}$,886_${ik}$,2881_${ik}$]
           mm(113_${ik}$,1_${ik}$:4_${ik}$)=[1364_${ik}$,2423_${ik}$,3514_${ik}$,3637_${ik}$]
           mm(114_${ik}$,1_${ik}$:4_${ik}$)=[1244_${ik}$,2051_${ik}$,1301_${ik}$,1465_${ik}$]
           mm(115_${ik}$,1_${ik}$:4_${ik}$)=[2020_${ik}$,2295_${ik}$,3604_${ik}$,2829_${ik}$]
           mm(116_${ik}$,1_${ik}$:4_${ik}$)=[3160_${ik}$,1332_${ik}$,1888_${ik}$,2161_${ik}$]
           mm(117_${ik}$,1_${ik}$:4_${ik}$)=[2785_${ik}$,1832_${ik}$,1836_${ik}$,3365_${ik}$]
           mm(118_${ik}$,1_${ik}$:4_${ik}$)=[2772_${ik}$,2405_${ik}$,1990_${ik}$,361_${ik}$]
           mm(119_${ik}$,1_${ik}$:4_${ik}$)=[1217_${ik}$,3638_${ik}$,2058_${ik}$,2685_${ik}$]
           mm(120_${ik}$,1_${ik}$:4_${ik}$)=[1822_${ik}$,3661_${ik}$,692_${ik}$,3745_${ik}$]
           mm(121_${ik}$,1_${ik}$:4_${ik}$)=[1245_${ik}$,327_${ik}$,1194_${ik}$,2325_${ik}$]
           mm(122_${ik}$,1_${ik}$:4_${ik}$)=[2252_${ik}$,3660_${ik}$,20_${ik}$,3609_${ik}$]
           mm(123_${ik}$,1_${ik}$:4_${ik}$)=[3904_${ik}$,716_${ik}$,3285_${ik}$,3821_${ik}$]
           mm(124_${ik}$,1_${ik}$:4_${ik}$)=[2774_${ik}$,1842_${ik}$,2046_${ik}$,3537_${ik}$]
           mm(125_${ik}$,1_${ik}$:4_${ik}$)=[997_${ik}$,3987_${ik}$,2107_${ik}$,517_${ik}$]
           mm(126_${ik}$,1_${ik}$:4_${ik}$)=[2573_${ik}$,1368_${ik}$,3508_${ik}$,3017_${ik}$]
           mm(127_${ik}$,1_${ik}$:4_${ik}$)=[1148_${ik}$,1848_${ik}$,3525_${ik}$,2141_${ik}$]
           mm(128_${ik}$,1_${ik}$:4_${ik}$)=[545_${ik}$,2366_${ik}$,3801_${ik}$,1537_${ik}$]
           ! Executable Statements 
           i1 = iseed( 1_${ik}$ )
           i2 = iseed( 2_${ik}$ )
           i3 = iseed( 3_${ik}$ )
           i4 = iseed( 4_${ik}$ )
           loop_10: do i = 1, min( n, lv )
           20 continue
              ! multiply the seed by i-th power of the multiplier modulo 2**48
              it4 = i4*mm( i, 4_${ik}$ )
              it3 = it4 / ipw2
              it4 = it4 - ipw2*it3
              it3 = it3 + i3*mm( i, 4_${ik}$ ) + i4*mm( i, 3_${ik}$ )
              it2 = it3 / ipw2
              it3 = it3 - ipw2*it2
              it2 = it2 + i2*mm( i, 4_${ik}$ ) + i3*mm( i, 3_${ik}$ ) + i4*mm( i, 2_${ik}$ )
              it1 = it2 / ipw2
              it2 = it2 - ipw2*it1
              it1 = it1 + i1*mm( i, 4_${ik}$ ) + i2*mm( i, 3_${ik}$ ) + i3*mm( i, 2_${ik}$ ) +i4*mm( i, 1_${ik}$ )
              it1 = mod( it1, ipw2 )
              ! convert 48-bit integer to a realnumber in the interval (0,1,KIND=${rk}$)
              x( i ) = r*( real( it1,KIND=${rk}$)+r*( real( it2,KIND=${rk}$)+r*( real( it3,KIND=${rk}$)+&
                        r*real( it4,KIND=${rk}$) ) ) )
              if (x( i )==one) then
                 ! if a real number has n bits of precision, and the first
                 ! n bits of the 48-bit integer above happen to be all 1 (which
                 ! will occur about once every 2**n calls), then x( i ) will
                 ! be rounded to exactly one.
                 ! since x( i ) is not supposed to return exactly 0.0_${rk}$ or one,
                 ! the statistically correct thing to do in this situation is
                 ! simply to iterate again.
                 ! n.b. the case x( i ) = 0.0_${rk}$ should not be possible.
                 i1 = i1 + 2_${ik}$
                 i2 = i2 + 2_${ik}$
                 i3 = i3 + 2_${ik}$
                 i4 = i4 + 2_${ik}$
                 goto 20
              end if
           end do loop_10
           ! return final value of seed
           iseed( 1_${ik}$ ) = it1
           iseed( 2_${ik}$ ) = it2
           iseed( 3_${ik}$ ) = it3
           iseed( 4_${ik}$ ) = it4
           return
     end subroutine stdlib${ii}$_${ri}$laruv

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_slacpy( uplo, m, n, a, lda, b, ldb )
     !! SLACPY copies all or part of a two-dimensional matrix A to another
     !! matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           real(sp), intent(in) :: a(lda,*)
           real(sp), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_slacpy

     pure module subroutine stdlib${ii}$_dlacpy( uplo, m, n, a, lda, b, ldb )
     !! DLACPY copies all or part of a two-dimensional matrix A to another
     !! matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           real(dp), intent(in) :: a(lda,*)
           real(dp), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_dlacpy

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ri}$lacpy( uplo, m, n, a, lda, b, ldb )
     !! DLACPY: copies all or part of a two-dimensional matrix A to another
     !! matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           real(${rk}$), intent(in) :: a(lda,*)
           real(${rk}$), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_${ri}$lacpy

#:endif
#:endfor

     pure module subroutine stdlib${ii}$_clacpy( uplo, m, n, a, lda, b, ldb )
     !! CLACPY copies all or part of a two-dimensional matrix A to another
     !! matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           complex(sp), intent(in) :: a(lda,*)
           complex(sp), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_clacpy

     pure module subroutine stdlib${ii}$_zlacpy( uplo, m, n, a, lda, b, ldb )
     !! ZLACPY copies all or part of a two-dimensional matrix A to another
     !! matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           complex(dp), intent(in) :: a(lda,*)
           complex(dp), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_zlacpy

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ci}$lacpy( uplo, m, n, a, lda, b, ldb )
     !! ZLACPY: copies all or part of a two-dimensional matrix A to another
     !! matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           complex(${ck}$), intent(in) :: a(lda,*)
           complex(${ck}$), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_${ci}$lacpy

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_clacp2( uplo, m, n, a, lda, b, ldb )
     !! CLACP2 copies all or part of a real two-dimensional matrix A to a
     !! complex matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           real(sp), intent(in) :: a(lda,*)
           complex(sp), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_clacp2

     pure module subroutine stdlib${ii}$_zlacp2( uplo, m, n, a, lda, b, ldb )
     !! ZLACP2 copies all or part of a real two-dimensional matrix A to a
     !! complex matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           real(dp), intent(in) :: a(lda,*)
           complex(dp), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_zlacp2

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ci}$lacp2( uplo, m, n, a, lda, b, ldb )
     !! ZLACP2: copies all or part of a real two-dimensional matrix A to a
     !! complex matrix B.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: uplo
           integer(${ik}$), intent(in) :: lda, ldb, m, n
           ! Array Arguments 
           real(${ck}$), intent(in) :: a(lda,*)
           complex(${ck}$), intent(out) :: b(ldb,*)
        ! =====================================================================
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( stdlib_lsame( uplo, 'U' ) ) then
              do j = 1, n
                 do i = 1, min( j, m )
                    b( i, j ) = a( i, j )
                 end do
              end do
           else if( stdlib_lsame( uplo, 'L' ) ) then
              do j = 1, n
                 do i = j, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           else
              do j = 1, n
                 do i = 1, m
                    b( i, j ) = a( i, j )
                 end do
              end do
           end if
           return
     end subroutine stdlib${ii}$_${ci}$lacp2

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_stfttp( transr, uplo, n, arf, ap, info )
     !! STFTTP copies a triangular matrix A from rectangular full packed
     !! format (TF) to standard packed format (TP).
        ! -- lapack computational routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: transr, uplo
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: n
           ! Array Arguments 
           real(sp), intent(out) :: ap(0_${ik}$:*)
           real(sp), intent(in) :: arf(0_${ik}$:*)
        ! =====================================================================
           ! Parameters 
           ! Local Scalars 
           logical(lk) :: lower, nisodd, normaltransr
           integer(${ik}$) :: n1, n2, k, nt
           integer(${ik}$) :: i, j, ij
           integer(${ik}$) :: ijp, jp, lda, js
           ! Executable Statements 
           ! test the input parameters.
           info = 0_${ik}$
           normaltransr = stdlib_lsame( transr, 'N' )
           lower = stdlib_lsame( uplo, 'L' )
           if( .not.normaltransr .and. .not.stdlib_lsame( transr, 'T' ) ) then
              info = -1_${ik}$
           else if( .not.lower .and. .not.stdlib_lsame( uplo, 'U' ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'STFTTP', -info )
              return
           end if
           ! quick return if possible
           if( n==0 )return
           if( n==1_${ik}$ ) then
              if( normaltransr ) then
                 ap( 0_${ik}$ ) = arf( 0_${ik}$ )
              else
                 ap( 0_${ik}$ ) = arf( 0_${ik}$ )
              end if
              return
           end if
           ! size of array arf(0:nt-1)
           nt = n*( n+1 ) / 2_${ik}$
           ! set n1 and n2 depending on lower
           if( lower ) then
              n2 = n / 2_${ik}$
              n1 = n - n2
           else
              n1 = n / 2_${ik}$
              n2 = n - n1
           end if
           ! if n is odd, set nisodd = .true.
           ! if n is even, set k = n/2 and nisodd = .false.
           ! set lda of arf^c; arf^c is (0:(n+1)/2-1,0:n-noe)
           ! where noe = 0 if n is even, noe = 1 if n is odd
           if( mod( n, 2_${ik}$ )==0_${ik}$ ) then
              k = n / 2_${ik}$
              nisodd = .false.
              lda = n + 1_${ik}$
           else
              nisodd = .true.
              lda = n
           end if
           ! arf^c has lda rows and n+1-noe cols
           if( .not.normaltransr )lda = ( n+1 ) / 2_${ik}$
           ! start execution: there are eight cases
           if( nisodd ) then
              ! n is odd
              if( normaltransr ) then
                 ! n is odd and transr = 'n'
                 if( lower ) then
                   ! srpa for lower, normal and n is odd ( a(0:n-1,0:n1-1) )
                   ! t1 -> a(0,0), t2 -> a(0,1), s -> a(n1,0)
                   ! t1 -> a(0), t2 -> a(n), s -> a(n1); lda = n
                    ijp = 0_${ik}$
                    jp = 0_${ik}$
                    do j = 0, n2
                       do i = j, n - 1
                          ij = i + jp
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       jp = jp + lda
                    end do
                    do i = 0, n2 - 1
                       do j = 1 + i, n2
                          ij = i + j*lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 else
                   ! srpa for upper, normal and n is odd ( a(0:n-1,0:n2-1)
                   ! t1 -> a(n1+1,0), t2 -> a(n1,0), s -> a(0,0)
                   ! t1 -> a(n2), t2 -> a(n1), s -> a(0)
                    ijp = 0_${ik}$
                    do j = 0, n1 - 1
                       ij = n2 + j
                       do i = 0, j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                          ij = ij + lda
                       end do
                    end do
                    js = 0_${ik}$
                    do j = n1, n - 1
                       ij = js
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                 end if
              else
                 ! n is odd and transr = 't'
                 if( lower ) then
                    ! srpa for lower, transpose and n is odd
                    ! t1 -> a(0,0) , t2 -> a(1,0) , s -> a(0,n1)
                    ! t1 -> a(0+0) , t2 -> a(1+0) , s -> a(0+n1*n1); lda=n1
                    ijp = 0_${ik}$
                    do i = 0, n2
                       do ij = i*( lda+1 ), n*lda - 1, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                    js = 1_${ik}$
                    do j = 0, n2 - 1
                       do ij = js, js + n2 - j - 1
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda + 1_${ik}$
                    end do
                 else
                    ! srpa for upper, transpose and n is odd
                    ! t1 -> a(0,n1+1), t2 -> a(0,n1), s -> a(0,0)
                    ! t1 -> a(n2*n2), t2 -> a(n1*n2), s -> a(0); lda = n2
                    ijp = 0_${ik}$
                    js = n2*lda
                    do j = 0, n1 - 1
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                    do i = 0, n1
                       do ij = i, i + ( n1+i )*lda, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 end if
              end if
           else
              ! n is even
              if( normaltransr ) then
                 ! n is even and transr = 'n'
                 if( lower ) then
                    ! srpa for lower, normal, and n is even ( a(0:n,0:k-1) )
                    ! t1 -> a(1,0), t2 -> a(0,0), s -> a(k+1,0)
                    ! t1 -> a(1), t2 -> a(0), s -> a(k+1)
                    ijp = 0_${ik}$
                    jp = 0_${ik}$
                    do j = 0, k - 1
                       do i = j, n - 1
                          ij = 1_${ik}$ + i + jp
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       jp = jp + lda
                    end do
                    do i = 0, k - 1
                       do j = i, k - 1
                          ij = i + j*lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 else
                    ! srpa for upper, normal, and n is even ( a(0:n,0:k-1) )
                    ! t1 -> a(k+1,0) ,  t2 -> a(k,0),   s -> a(0,0)
                    ! t1 -> a(k+1), t2 -> a(k), s -> a(0)
                    ijp = 0_${ik}$
                    do j = 0, k - 1
                       ij = k + 1_${ik}$ + j
                       do i = 0, j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                          ij = ij + lda
                       end do
                    end do
                    js = 0_${ik}$
                    do j = k, n - 1
                       ij = js
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                 end if
              else
                 ! n is even and transr = 't'
                 if( lower ) then
                    ! srpa for lower, transpose and n is even (see paper)
                    ! t1 -> b(0,1), t2 -> b(0,0), s -> b(0,k+1)
                    ! t1 -> a(0+k), t2 -> a(0+0), s -> a(0+k*(k+1)); lda=k
                    ijp = 0_${ik}$
                    do i = 0, k - 1
                       do ij = i + ( i+1 )*lda, ( n+1 )*lda - 1, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                    js = 0_${ik}$
                    do j = 0, k - 1
                       do ij = js, js + k - j - 1
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda + 1_${ik}$
                    end do
                 else
                    ! srpa for upper, transpose and n is even (see paper)
                    ! t1 -> b(0,k+1),     t2 -> b(0,k),   s -> b(0,0)
                    ! t1 -> a(0+k*(k+1)), t2 -> a(0+k*k), s -> a(0+0)); lda=k
                    ijp = 0_${ik}$
                    js = ( k+1 )*lda
                    do j = 0, k - 1
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                    do i = 0, k - 1
                       do ij = i, i + ( k+i )*lda, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 end if
              end if
           end if
           return
     end subroutine stdlib${ii}$_stfttp

     pure module subroutine stdlib${ii}$_dtfttp( transr, uplo, n, arf, ap, info )
     !! DTFTTP copies a triangular matrix A from rectangular full packed
     !! format (TF) to standard packed format (TP).
        ! -- lapack computational routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: transr, uplo
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: n
           ! Array Arguments 
           real(dp), intent(out) :: ap(0_${ik}$:*)
           real(dp), intent(in) :: arf(0_${ik}$:*)
        ! =====================================================================
           ! Parameters 
           ! Local Scalars 
           logical(lk) :: lower, nisodd, normaltransr
           integer(${ik}$) :: n1, n2, k, nt
           integer(${ik}$) :: i, j, ij
           integer(${ik}$) :: ijp, jp, lda, js
           ! Executable Statements 
           ! test the input parameters.
           info = 0_${ik}$
           normaltransr = stdlib_lsame( transr, 'N' )
           lower = stdlib_lsame( uplo, 'L' )
           if( .not.normaltransr .and. .not.stdlib_lsame( transr, 'T' ) ) then
              info = -1_${ik}$
           else if( .not.lower .and. .not.stdlib_lsame( uplo, 'U' ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'DTFTTP', -info )
              return
           end if
           ! quick return if possible
           if( n==0 )return
           if( n==1_${ik}$ ) then
              if( normaltransr ) then
                 ap( 0_${ik}$ ) = arf( 0_${ik}$ )
              else
                 ap( 0_${ik}$ ) = arf( 0_${ik}$ )
              end if
              return
           end if
           ! size of array arf(0:nt-1)
           nt = n*( n+1 ) / 2_${ik}$
           ! set n1 and n2 depending on lower
           if( lower ) then
              n2 = n / 2_${ik}$
              n1 = n - n2
           else
              n1 = n / 2_${ik}$
              n2 = n - n1
           end if
           ! if n is odd, set nisodd = .true.
           ! if n is even, set k = n/2 and nisodd = .false.
           ! set lda of arf^c; arf^c is (0:(n+1)/2-1,0:n-noe)
           ! where noe = 0 if n is even, noe = 1 if n is odd
           if( mod( n, 2_${ik}$ )==0_${ik}$ ) then
              k = n / 2_${ik}$
              nisodd = .false.
              lda = n + 1_${ik}$
           else
              nisodd = .true.
              lda = n
           end if
           ! arf^c has lda rows and n+1-noe cols
           if( .not.normaltransr )lda = ( n+1 ) / 2_${ik}$
           ! start execution: there are eight cases
           if( nisodd ) then
              ! n is odd
              if( normaltransr ) then
                 ! n is odd and transr = 'n'
                 if( lower ) then
                   ! srpa for lower, normal and n is odd ( a(0:n-1,0:n1-1) )
                   ! t1 -> a(0,0), t2 -> a(0,1), s -> a(n1,0)
                   ! t1 -> a(0), t2 -> a(n), s -> a(n1); lda = n
                    ijp = 0_${ik}$
                    jp = 0_${ik}$
                    do j = 0, n2
                       do i = j, n - 1
                          ij = i + jp
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       jp = jp + lda
                    end do
                    do i = 0, n2 - 1
                       do j = 1 + i, n2
                          ij = i + j*lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 else
                   ! srpa for upper, normal and n is odd ( a(0:n-1,0:n2-1)
                   ! t1 -> a(n1+1,0), t2 -> a(n1,0), s -> a(0,0)
                   ! t1 -> a(n2), t2 -> a(n1), s -> a(0)
                    ijp = 0_${ik}$
                    do j = 0, n1 - 1
                       ij = n2 + j
                       do i = 0, j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                          ij = ij + lda
                       end do
                    end do
                    js = 0_${ik}$
                    do j = n1, n - 1
                       ij = js
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                 end if
              else
                 ! n is odd and transr = 't'
                 if( lower ) then
                    ! srpa for lower, transpose and n is odd
                    ! t1 -> a(0,0) , t2 -> a(1,0) , s -> a(0,n1)
                    ! t1 -> a(0+0) , t2 -> a(1+0) , s -> a(0+n1*n1); lda=n1
                    ijp = 0_${ik}$
                    do i = 0, n2
                       do ij = i*( lda+1 ), n*lda - 1, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                    js = 1_${ik}$
                    do j = 0, n2 - 1
                       do ij = js, js + n2 - j - 1
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda + 1_${ik}$
                    end do
                 else
                    ! srpa for upper, transpose and n is odd
                    ! t1 -> a(0,n1+1), t2 -> a(0,n1), s -> a(0,0)
                    ! t1 -> a(n2*n2), t2 -> a(n1*n2), s -> a(0); lda = n2
                    ijp = 0_${ik}$
                    js = n2*lda
                    do j = 0, n1 - 1
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                    do i = 0, n1
                       do ij = i, i + ( n1+i )*lda, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 end if
              end if
           else
              ! n is even
              if( normaltransr ) then
                 ! n is even and transr = 'n'
                 if( lower ) then
                    ! srpa for lower, normal, and n is even ( a(0:n,0:k-1) )
                    ! t1 -> a(1,0), t2 -> a(0,0), s -> a(k+1,0)
                    ! t1 -> a(1), t2 -> a(0), s -> a(k+1)
                    ijp = 0_${ik}$
                    jp = 0_${ik}$
                    do j = 0, k - 1
                       do i = j, n - 1
                          ij = 1_${ik}$ + i + jp
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       jp = jp + lda
                    end do
                    do i = 0, k - 1
                       do j = i, k - 1
                          ij = i + j*lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 else
                    ! srpa for upper, normal, and n is even ( a(0:n,0:k-1) )
                    ! t1 -> a(k+1,0) ,  t2 -> a(k,0),   s -> a(0,0)
                    ! t1 -> a(k+1), t2 -> a(k), s -> a(0)
                    ijp = 0_${ik}$
                    do j = 0, k - 1
                       ij = k + 1_${ik}$ + j
                       do i = 0, j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                          ij = ij + lda
                       end do
                    end do
                    js = 0_${ik}$
                    do j = k, n - 1
                       ij = js
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                 end if
              else
                 ! n is even and transr = 't'
                 if( lower ) then
                    ! srpa for lower, transpose and n is even (see paper)
                    ! t1 -> b(0,1), t2 -> b(0,0), s -> b(0,k+1)
                    ! t1 -> a(0+k), t2 -> a(0+0), s -> a(0+k*(k+1)); lda=k
                    ijp = 0_${ik}$
                    do i = 0, k - 1
                       do ij = i + ( i+1 )*lda, ( n+1 )*lda - 1, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                    js = 0_${ik}$
                    do j = 0, k - 1
                       do ij = js, js + k - j - 1
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda + 1_${ik}$
                    end do
                 else
                    ! srpa for upper, transpose and n is even (see paper)
                    ! t1 -> b(0,k+1),     t2 -> b(0,k),   s -> b(0,0)
                    ! t1 -> a(0+k*(k+1)), t2 -> a(0+k*k), s -> a(0+0)); lda=k
                    ijp = 0_${ik}$
                    js = ( k+1 )*lda
                    do j = 0, k - 1
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                    do i = 0, k - 1
                       do ij = i, i + ( k+i )*lda, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 end if
              end if
           end if
           return
     end subroutine stdlib${ii}$_dtfttp

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ri}$tfttp( transr, uplo, n, arf, ap, info )
     !! DTFTTP: copies a triangular matrix A from rectangular full packed
     !! format (TF) to standard packed format (TP).
        ! -- lapack computational routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: transr, uplo
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: n
           ! Array Arguments 
           real(${rk}$), intent(out) :: ap(0_${ik}$:*)
           real(${rk}$), intent(in) :: arf(0_${ik}$:*)
        ! =====================================================================
           ! Parameters 
           ! Local Scalars 
           logical(lk) :: lower, nisodd, normaltransr
           integer(${ik}$) :: n1, n2, k, nt
           integer(${ik}$) :: i, j, ij
           integer(${ik}$) :: ijp, jp, lda, js
           ! Executable Statements 
           ! test the input parameters.
           info = 0_${ik}$
           normaltransr = stdlib_lsame( transr, 'N' )
           lower = stdlib_lsame( uplo, 'L' )
           if( .not.normaltransr .and. .not.stdlib_lsame( transr, 'T' ) ) then
              info = -1_${ik}$
           else if( .not.lower .and. .not.stdlib_lsame( uplo, 'U' ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'DTFTTP', -info )
              return
           end if
           ! quick return if possible
           if( n==0 )return
           if( n==1_${ik}$ ) then
              if( normaltransr ) then
                 ap( 0_${ik}$ ) = arf( 0_${ik}$ )
              else
                 ap( 0_${ik}$ ) = arf( 0_${ik}$ )
              end if
              return
           end if
           ! size of array arf(0:nt-1)
           nt = n*( n+1 ) / 2_${ik}$
           ! set n1 and n2 depending on lower
           if( lower ) then
              n2 = n / 2_${ik}$
              n1 = n - n2
           else
              n1 = n / 2_${ik}$
              n2 = n - n1
           end if
           ! if n is odd, set nisodd = .true.
           ! if n is even, set k = n/2 and nisodd = .false.
           ! set lda of arf^c; arf^c is (0:(n+1)/2-1,0:n-noe)
           ! where noe = 0 if n is even, noe = 1 if n is odd
           if( mod( n, 2_${ik}$ )==0_${ik}$ ) then
              k = n / 2_${ik}$
              nisodd = .false.
              lda = n + 1_${ik}$
           else
              nisodd = .true.
              lda = n
           end if
           ! arf^c has lda rows and n+1-noe cols
           if( .not.normaltransr )lda = ( n+1 ) / 2_${ik}$
           ! start execution: there are eight cases
           if( nisodd ) then
              ! n is odd
              if( normaltransr ) then
                 ! n is odd and transr = 'n'
                 if( lower ) then
                   ! srpa for lower, normal and n is odd ( a(0:n-1,0:n1-1) )
                   ! t1 -> a(0,0), t2 -> a(0,1), s -> a(n1,0)
                   ! t1 -> a(0), t2 -> a(n), s -> a(n1); lda = n
                    ijp = 0_${ik}$
                    jp = 0_${ik}$
                    do j = 0, n2
                       do i = j, n - 1
                          ij = i + jp
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       jp = jp + lda
                    end do
                    do i = 0, n2 - 1
                       do j = 1 + i, n2
                          ij = i + j*lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 else
                   ! srpa for upper, normal and n is odd ( a(0:n-1,0:n2-1)
                   ! t1 -> a(n1+1,0), t2 -> a(n1,0), s -> a(0,0)
                   ! t1 -> a(n2), t2 -> a(n1), s -> a(0)
                    ijp = 0_${ik}$
                    do j = 0, n1 - 1
                       ij = n2 + j
                       do i = 0, j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                          ij = ij + lda
                       end do
                    end do
                    js = 0_${ik}$
                    do j = n1, n - 1
                       ij = js
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                 end if
              else
                 ! n is odd and transr = 't'
                 if( lower ) then
                    ! srpa for lower, transpose and n is odd
                    ! t1 -> a(0,0) , t2 -> a(1,0) , s -> a(0,n1)
                    ! t1 -> a(0+0) , t2 -> a(1+0) , s -> a(0+n1*n1); lda=n1
                    ijp = 0_${ik}$
                    do i = 0, n2
                       do ij = i*( lda+1 ), n*lda - 1, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                    js = 1_${ik}$
                    do j = 0, n2 - 1
                       do ij = js, js + n2 - j - 1
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda + 1_${ik}$
                    end do
                 else
                    ! srpa for upper, transpose and n is odd
                    ! t1 -> a(0,n1+1), t2 -> a(0,n1), s -> a(0,0)
                    ! t1 -> a(n2*n2), t2 -> a(n1*n2), s -> a(0); lda = n2
                    ijp = 0_${ik}$
                    js = n2*lda
                    do j = 0, n1 - 1
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                    do i = 0, n1
                       do ij = i, i + ( n1+i )*lda, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 end if
              end if
           else
              ! n is even
              if( normaltransr ) then
                 ! n is even and transr = 'n'
                 if( lower ) then
                    ! srpa for lower, normal, and n is even ( a(0:n,0:k-1) )
                    ! t1 -> a(1,0), t2 -> a(0,0), s -> a(k+1,0)
                    ! t1 -> a(1), t2 -> a(0), s -> a(k+1)
                    ijp = 0_${ik}$
                    jp = 0_${ik}$
                    do j = 0, k - 1
                       do i = j, n - 1
                          ij = 1_${ik}$ + i + jp
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       jp = jp + lda
                    end do
                    do i = 0, k - 1
                       do j = i, k - 1
                          ij = i + j*lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 else
                    ! srpa for upper, normal, and n is even ( a(0:n,0:k-1) )
                    ! t1 -> a(k+1,0) ,  t2 -> a(k,0),   s -> a(0,0)
                    ! t1 -> a(k+1), t2 -> a(k), s -> a(0)
                    ijp = 0_${ik}$
                    do j = 0, k - 1
                       ij = k + 1_${ik}$ + j
                       do i = 0, j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                          ij = ij + lda
                       end do
                    end do
                    js = 0_${ik}$
                    do j = k, n - 1
                       ij = js
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                 end if
              else
                 ! n is even and transr = 't'
                 if( lower ) then
                    ! srpa for lower, transpose and n is even (see paper)
                    ! t1 -> b(0,1), t2 -> b(0,0), s -> b(0,k+1)
                    ! t1 -> a(0+k), t2 -> a(0+0), s -> a(0+k*(k+1)); lda=k
                    ijp = 0_${ik}$
                    do i = 0, k - 1
                       do ij = i + ( i+1 )*lda, ( n+1 )*lda - 1, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                    js = 0_${ik}$
                    do j = 0, k - 1
                       do ij = js, js + k - j - 1
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda + 1_${ik}$
                    end do
                 else
                    ! srpa for upper, transpose and n is even (see paper)
                    ! t1 -> b(0,k+1),     t2 -> b(0,k),   s -> b(0,0)
                    ! t1 -> a(0+k*(k+1)), t2 -> a(0+k*k), s -> a(0+0)); lda=k
                    ijp = 0_${ik}$
                    js = ( k+1 )*lda
                    do j = 0, k - 1
                       do ij = js, js + j
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                       js = js + lda
                    end do
                    do i = 0, k - 1
                       do ij = i, i + ( k+i )*lda, lda
                          ap( ijp ) = arf( ij )
                          ijp = ijp + 1_${ik}$
                       end do
                    end do
                 end if
              end if
           end if
           return
     end subroutine stdlib${ii}$_${ri}$tfttp

#:endif
#:endfor

     pure module subroutine stdlib${ii}$_ctfttp( transr, uplo, n, arf, ap, info )
     !! CTFTTP copies a triangular matrix A from rectangular full packed
     !! format (TF) to standard packed format (TP).
        ! -- lapack computational routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           ! Scalar Arguments 
           character, intent(in) :: transr, uplo
           integer(${ik}$), intent(out) :: info
           integer(${ik}$), intent(in) :: n
           ! Array Arguments 
           complex(sp), intent(out) :: ap(0_${ik}$:*)
           complex(sp), intent(in) :: arf(0_${ik}$:*)
        ! =====================================================================
           ! Parameters 
           ! Local Scalars 
           logical(lk) :: lower, nisodd, normaltransr
           integer(${ik}$) :: n1, n2, k, nt
           integer(${ik}$) :: i, j, ij
           integer(${ik}$) :: ijp, jp, lda, js
           ! Intrinsic Functions 
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input parameters.
           info = 0_${ik}$
           normaltransr = stdlib_lsame( transr, 'N' )
           lower = stdlib_lsame( uplo, 'L' )
           if( .not.normaltransr .and. .not.stdlib_lsame( transr, 'C' ) ) then
              info = -1_${ik}$
           else if( .not.lower .and. .not.stdlib_lsame( uplo, 'U' ) ) then
              info = -2_${ik}$
           else if( n<0_${ik}$ ) then
              info = -3_${ik}$
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'CTFTTP', -info )
              return
           end if
           ! quick return if possible
           if( n==0 )return
           if( n==1_${ik}$ ) then
              if( normaltransr ) then
                 ap( 0_${ik}$ ) = arf( 0_${ik}$ )
              else
                 ap( 0_${ik}$ ) = conjg( arf( 0_${ik}$ ) )
              end if
              return
           end if
           ! size of array arf(0:nt-1)
           nt = n*( n+1 ) / 2_${ik}$
           ! set n1 and n2 depending on lower
           if( lower ) then
              n2 = n / 2_${ik}$
              n1 = n - n2
           else
              n1 = n / 2_${ik}$
              n2 = n - n1
           end if
           ! if n is odd, set nisodd = .true.
           ! if n is even, set k = n/2 and nisodd = .false.
           ! set lda of arf^c; arf^c is (0:(n+1)/2-1,0:n-noe)