stdlib_lapack_blas_like_l3.fypp Source File


Source Code

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


  contains
#:for ik,it,ii in LINALG_INT_KINDS_TYPES

     pure module subroutine stdlib${ii}$_slagtm( trans, n, nrhs, alpha, dl, d, du, x, ldx, beta,b, ldb )
     !! SLAGTM performs a matrix-vector product of the form
     !! B := alpha * A * X + beta * B
     !! where A is a tridiagonal matrix of order N, B and X are N by NRHS
     !! matrices, and alpha and beta are real scalars, each of which may be
     !! 0., 1., or -1.
               
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: trans
           integer(${ik}$), intent(in) :: ldb, ldx, n, nrhs
           real(sp), intent(in) :: alpha, beta
           ! Array Arguments 
           real(sp), intent(inout) :: b(ldb,*)
           real(sp), intent(in) :: d(*), dl(*), du(*), x(ldx,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Executable Statements 
           if( n==0 )return
           ! multiply b by beta if beta/=1.
           if( beta==zero ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = zero
                 end do
              end do
           else if( beta==-one ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = -b( i, j )
                 end do
              end do
           end if
           if( alpha==one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b + a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + dl( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + dl( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else
                 ! compute b := b + a**t*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + du( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + du( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           else if( alpha==-one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b - a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - dl( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - dl( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else
                 ! compute b := b - a**t*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - du( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - du( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           end if
           return
     end subroutine stdlib${ii}$_slagtm

     pure module subroutine stdlib${ii}$_dlagtm( trans, n, nrhs, alpha, dl, d, du, x, ldx, beta,b, ldb )
     !! DLAGTM performs a matrix-vector product of the form
     !! B := alpha * A * X + beta * B
     !! where A is a tridiagonal matrix of order N, B and X are N by NRHS
     !! matrices, and alpha and beta are real scalars, each of which may be
     !! 0., 1., or -1.
               
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: trans
           integer(${ik}$), intent(in) :: ldb, ldx, n, nrhs
           real(dp), intent(in) :: alpha, beta
           ! Array Arguments 
           real(dp), intent(inout) :: b(ldb,*)
           real(dp), intent(in) :: d(*), dl(*), du(*), x(ldx,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Executable Statements 
           if( n==0 )return
           ! multiply b by beta if beta/=1.
           if( beta==zero ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = zero
                 end do
              end do
           else if( beta==-one ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = -b( i, j )
                 end do
              end do
           end if
           if( alpha==one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b + a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + dl( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + dl( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else
                 ! compute b := b + a**t*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + du( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + du( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           else if( alpha==-one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b - a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - dl( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - dl( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else
                 ! compute b := b - a**t*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - du( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - du( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           end if
           return
     end subroutine stdlib${ii}$_dlagtm

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if not rk in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ri}$lagtm( trans, n, nrhs, alpha, dl, d, du, x, ldx, beta,b, ldb )
     !! DLAGTM: performs a matrix-vector product of the form
     !! B := alpha * A * X + beta * B
     !! where A is a tridiagonal matrix of order N, B and X are N by NRHS
     !! matrices, and alpha and beta are real scalars, each of which may be
     !! 0., 1., or -1.
               
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_${rk}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: trans
           integer(${ik}$), intent(in) :: ldb, ldx, n, nrhs
           real(${rk}$), intent(in) :: alpha, beta
           ! Array Arguments 
           real(${rk}$), intent(inout) :: b(ldb,*)
           real(${rk}$), intent(in) :: d(*), dl(*), du(*), x(ldx,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Executable Statements 
           if( n==0 )return
           ! multiply b by beta if beta/=1.
           if( beta==zero ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = zero
                 end do
              end do
           else if( beta==-one ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = -b( i, j )
                 end do
              end do
           end if
           if( alpha==one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b + a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + dl( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + dl( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else
                 ! compute b := b + a**t*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + du( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + du( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           else if( alpha==-one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b - a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - dl( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - dl( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else
                 ! compute b := b - a**t*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - du( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - du( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           end if
           return
     end subroutine stdlib${ii}$_${ri}$lagtm

#:endif
#:endfor

     pure module subroutine stdlib${ii}$_clagtm( trans, n, nrhs, alpha, dl, d, du, x, ldx, beta,b, ldb )
     !! CLAGTM performs a matrix-vector product of the form
     !! B := alpha * A * X + beta * B
     !! where A is a tridiagonal matrix of order N, B and X are N by NRHS
     !! matrices, and alpha and beta are real scalars, each of which may be
     !! 0., 1., or -1.
               
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: trans
           integer(${ik}$), intent(in) :: ldb, ldx, n, nrhs
           real(sp), intent(in) :: alpha, beta
           ! Array Arguments 
           complex(sp), intent(inout) :: b(ldb,*)
           complex(sp), intent(in) :: d(*), dl(*), du(*), x(ldx,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( n==0 )return
           ! multiply b by beta if beta/=1.
           if( beta==zero ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = zero
                 end do
              end do
           else if( beta==-one ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = -b( i, j )
                 end do
              end do
           end if
           if( alpha==one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b + a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + dl( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + dl( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'T' ) ) then
                 ! compute b := b + a**t * x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + du( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + du( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'C' ) ) then
                 ! compute b := b + a**h * x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j ) +conjg( dl( 1_${ik}$ ) )*x( 2_${ik}$, &
                                 j )
                       b( n, j ) = b( n, j ) + conjg( du( n-1 ) )*x( n-1, j ) + conjg( d( n ) )*x(&
                                  n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + conjg( du( i-1 ) )*x( i-1, j ) + conjg( d( i ) )&
                                    *x( i, j ) + conjg( dl( i ) )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           else if( alpha==-one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b - a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - dl( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - dl( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'T' ) ) then
                 ! compute b := b - a**t*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - du( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - du( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'C' ) ) then
                 ! compute b := b - a**h*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j ) -conjg( dl( 1_${ik}$ ) )*x( 2_${ik}$, &
                                 j )
                       b( n, j ) = b( n, j ) - conjg( du( n-1 ) )*x( n-1, j ) - conjg( d( n ) )*x(&
                                  n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - conjg( du( i-1 ) )*x( i-1, j ) - conjg( d( i ) )&
                                    *x( i, j ) - conjg( dl( i ) )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           end if
           return
     end subroutine stdlib${ii}$_clagtm

     pure module subroutine stdlib${ii}$_zlagtm( trans, n, nrhs, alpha, dl, d, du, x, ldx, beta,b, ldb )
     !! ZLAGTM performs a matrix-vector product of the form
     !! B := alpha * A * X + beta * B
     !! where A is a tridiagonal matrix of order N, B and X are N by NRHS
     !! matrices, and alpha and beta are real scalars, each of which may be
     !! 0., 1., or -1.
               
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: trans
           integer(${ik}$), intent(in) :: ldb, ldx, n, nrhs
           real(dp), intent(in) :: alpha, beta
           ! Array Arguments 
           complex(dp), intent(inout) :: b(ldb,*)
           complex(dp), intent(in) :: d(*), dl(*), du(*), x(ldx,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( n==0 )return
           ! multiply b by beta if beta/=1.
           if( beta==zero ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = zero
                 end do
              end do
           else if( beta==-one ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = -b( i, j )
                 end do
              end do
           end if
           if( alpha==one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b + a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + dl( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + dl( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'T' ) ) then
                 ! compute b := b + a**t * x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + du( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + du( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'C' ) ) then
                 ! compute b := b + a**h * x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j ) +conjg( dl( 1_${ik}$ ) )*x( 2_${ik}$, &
                                 j )
                       b( n, j ) = b( n, j ) + conjg( du( n-1 ) )*x( n-1, j ) + conjg( d( n ) )*x(&
                                  n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + conjg( du( i-1 ) )*x( i-1, j ) + conjg( d( i ) )&
                                    *x( i, j ) + conjg( dl( i ) )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           else if( alpha==-one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b - a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - dl( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - dl( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'T' ) ) then
                 ! compute b := b - a**t *x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - du( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - du( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'C' ) ) then
                 ! compute b := b - a**h *x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j ) -conjg( dl( 1_${ik}$ ) )*x( 2_${ik}$, &
                                 j )
                       b( n, j ) = b( n, j ) - conjg( du( n-1 ) )*x( n-1, j ) - conjg( d( n ) )*x(&
                                  n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - conjg( du( i-1 ) )*x( i-1, j ) - conjg( d( i ) )&
                                    *x( i, j ) - conjg( dl( i ) )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           end if
           return
     end subroutine stdlib${ii}$_zlagtm

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ci}$lagtm( trans, n, nrhs, alpha, dl, d, du, x, ldx, beta,b, ldb )
     !! ZLAGTM: performs a matrix-vector product of the form
     !! B := alpha * A * X + beta * B
     !! where A is a tridiagonal matrix of order N, B and X are N by NRHS
     !! matrices, and alpha and beta are real scalars, each of which may be
     !! 0., 1., or -1.
               
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_${ck}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: trans
           integer(${ik}$), intent(in) :: ldb, ldx, n, nrhs
           real(${ck}$), intent(in) :: alpha, beta
           ! Array Arguments 
           complex(${ck}$), intent(inout) :: b(ldb,*)
           complex(${ck}$), intent(in) :: d(*), dl(*), du(*), x(ldx,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           if( n==0 )return
           ! multiply b by beta if beta/=1.
           if( beta==zero ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = zero
                 end do
              end do
           else if( beta==-one ) then
              do j = 1, nrhs
                 do i = 1, n
                    b( i, j ) = -b( i, j )
                 end do
              end do
           end if
           if( alpha==one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b + a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + dl( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + dl( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'T' ) ) then
                 ! compute b := b + a**t * x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + d( 1_${ik}$ )*x( 1_${ik}$, j ) +dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) + du( n-1 )*x( n-1, j ) +d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + du( i-1 )*x( i-1, j ) +d( i )*x( i, j ) + dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'C' ) ) then
                 ! compute b := b + a**h * x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) + conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j ) +conjg( dl( 1_${ik}$ ) )*x( 2_${ik}$, &
                                 j )
                       b( n, j ) = b( n, j ) + conjg( du( n-1 ) )*x( n-1, j ) + conjg( d( n ) )*x(&
                                  n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) + conjg( du( i-1 ) )*x( i-1, j ) + conjg( d( i ) )&
                                    *x( i, j ) + conjg( dl( i ) )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           else if( alpha==-one ) then
              if( stdlib_lsame( trans, 'N' ) ) then
                 ! compute b := b - a*x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -du( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - dl( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - dl( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - du( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'T' ) ) then
                 ! compute b := b - a**t *x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - d( 1_${ik}$ )*x( 1_${ik}$, j ) -dl( 1_${ik}$ )*x( 2_${ik}$, j )
                       b( n, j ) = b( n, j ) - du( n-1 )*x( n-1, j ) -d( n )*x( n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - du( i-1 )*x( i-1, j ) -d( i )*x( i, j ) - dl( i &
                                    )*x( i+1, j )
                       end do
                    end if
                 end do
              else if( stdlib_lsame( trans, 'C' ) ) then
                 ! compute b := b - a**h *x
                 do j = 1, nrhs
                    if( n==1_${ik}$ ) then
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j )
                    else
                       b( 1_${ik}$, j ) = b( 1_${ik}$, j ) - conjg( d( 1_${ik}$ ) )*x( 1_${ik}$, j ) -conjg( dl( 1_${ik}$ ) )*x( 2_${ik}$, &
                                 j )
                       b( n, j ) = b( n, j ) - conjg( du( n-1 ) )*x( n-1, j ) - conjg( d( n ) )*x(&
                                  n, j )
                       do i = 2, n - 1
                          b( i, j ) = b( i, j ) - conjg( du( i-1 ) )*x( i-1, j ) - conjg( d( i ) )&
                                    *x( i, j ) - conjg( dl( i ) )*x( i+1, j )
                       end do
                    end if
                 end do
              end if
           end if
           return
     end subroutine stdlib${ii}$_${ci}$lagtm

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_clacrm( m, n, a, lda, b, ldb, c, ldc, rwork )
     !! CLACRM performs a very simple matrix-matrix multiplication:
     !! C := A * B,
     !! where A is M by N and complex; B is N by N and real;
     !! C is M by N and complex.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: lda, ldb, ldc, m, n
           ! Array Arguments 
           real(sp), intent(in) :: b(ldb,*)
           real(sp), intent(out) :: rwork(*)
           complex(sp), intent(in) :: a(lda,*)
           complex(sp), intent(out) :: c(ldc,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j, l
           ! Intrinsic Functions 
           ! Executable Statements 
           ! quick return if possible.
           if( ( m==0 ) .or. ( n==0 ) )return
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = real( a( i, j ),KIND=sp)
              end do
           end do
           l = m*n + 1_${ik}$
           call stdlib${ii}$_sgemm( 'N', 'N', m, n, n, one, rwork, m, b, ldb, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = rwork( l+( j-1 )*m+i-1 )
              end do
           end do
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = aimag( a( i, j ) )
              end do
           end do
           call stdlib${ii}$_sgemm( 'N', 'N', m, n, n, one, rwork, m, b, ldb, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = cmplx( real( c( i, j ),KIND=sp),rwork( l+( j-1 )*m+i-1 ),KIND=sp)
                           
              end do
           end do
           return
     end subroutine stdlib${ii}$_clacrm

     pure module subroutine stdlib${ii}$_zlacrm( m, n, a, lda, b, ldb, c, ldc, rwork )
     !! ZLACRM performs a very simple matrix-matrix multiplication:
     !! C := A * B,
     !! where A is M by N and complex; B is N by N and real;
     !! C is M by N and complex.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: lda, ldb, ldc, m, n
           ! Array Arguments 
           real(dp), intent(in) :: b(ldb,*)
           real(dp), intent(out) :: rwork(*)
           complex(dp), intent(in) :: a(lda,*)
           complex(dp), intent(out) :: c(ldc,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j, l
           ! Intrinsic Functions 
           ! Executable Statements 
           ! quick return if possible.
           if( ( m==0 ) .or. ( n==0 ) )return
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = real( a( i, j ),KIND=dp)
              end do
           end do
           l = m*n + 1_${ik}$
           call stdlib${ii}$_dgemm( 'N', 'N', m, n, n, one, rwork, m, b, ldb, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = rwork( l+( j-1 )*m+i-1 )
              end do
           end do
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = aimag( a( i, j ) )
              end do
           end do
           call stdlib${ii}$_dgemm( 'N', 'N', m, n, n, one, rwork, m, b, ldb, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = cmplx( real( c( i, j ),KIND=dp),rwork( l+( j-1 )*m+i-1 ),KIND=dp)
                           
              end do
           end do
           return
     end subroutine stdlib${ii}$_zlacrm

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ci}$lacrm( m, n, a, lda, b, ldb, c, ldc, rwork )
     !! ZLACRM: performs a very simple matrix-matrix multiplication:
     !! C := A * B,
     !! where A is M by N and complex; B is N by N and real;
     !! C is M by N and complex.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_${ck}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: lda, ldb, ldc, m, n
           ! Array Arguments 
           real(${ck}$), intent(in) :: b(ldb,*)
           real(${ck}$), intent(out) :: rwork(*)
           complex(${ck}$), intent(in) :: a(lda,*)
           complex(${ck}$), intent(out) :: c(ldc,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j, l
           ! Intrinsic Functions 
           ! Executable Statements 
           ! quick return if possible.
           if( ( m==0 ) .or. ( n==0 ) )return
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = real( a( i, j ),KIND=${ck}$)
              end do
           end do
           l = m*n + 1_${ik}$
           call stdlib${ii}$_${c2ri(ci)}$gemm( 'N', 'N', m, n, n, one, rwork, m, b, ldb, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = rwork( l+( j-1 )*m+i-1 )
              end do
           end do
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = aimag( a( i, j ) )
              end do
           end do
           call stdlib${ii}$_${c2ri(ci)}$gemm( 'N', 'N', m, n, n, one, rwork, m, b, ldb, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = cmplx( real( c( i, j ),KIND=${ck}$),rwork( l+( j-1 )*m+i-1 ),KIND=${ck}$)
                           
              end do
           end do
           return
     end subroutine stdlib${ii}$_${ci}$lacrm

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_clarcm( m, n, a, lda, b, ldb, c, ldc, rwork )
     !! CLARCM performs a very simple matrix-matrix multiplication:
     !! C := A * B,
     !! where A is M by M and real; B is M by N and complex;
     !! C is M by N and complex.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: lda, ldb, ldc, m, n
           ! Array Arguments 
           real(sp), intent(in) :: a(lda,*)
           real(sp), intent(out) :: rwork(*)
           complex(sp), intent(in) :: b(ldb,*)
           complex(sp), intent(out) :: c(ldc,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j, l
           ! Intrinsic Functions 
           ! Executable Statements 
           ! quick return if possible.
           if( ( m==0 ) .or. ( n==0 ) )return
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = real( b( i, j ),KIND=sp)
              end do
           end do
           l = m*n + 1_${ik}$
           call stdlib${ii}$_sgemm( 'N', 'N', m, n, m, one, a, lda, rwork, m, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = rwork( l+( j-1 )*m+i-1 )
              end do
           end do
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = aimag( b( i, j ) )
              end do
           end do
           call stdlib${ii}$_sgemm( 'N', 'N', m, n, m, one, a, lda, rwork, m, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = cmplx( real( c( i, j ),KIND=sp),rwork( l+( j-1 )*m+i-1 ),KIND=sp)
                           
              end do
           end do
           return
     end subroutine stdlib${ii}$_clarcm

     pure module subroutine stdlib${ii}$_zlarcm( m, n, a, lda, b, ldb, c, ldc, rwork )
     !! ZLARCM performs a very simple matrix-matrix multiplication:
     !! C := A * B,
     !! where A is M by M and real; B is M by N and complex;
     !! C is M by N and complex.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: lda, ldb, ldc, m, n
           ! Array Arguments 
           real(dp), intent(in) :: a(lda,*)
           real(dp), intent(out) :: rwork(*)
           complex(dp), intent(in) :: b(ldb,*)
           complex(dp), intent(out) :: c(ldc,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j, l
           ! Intrinsic Functions 
           ! Executable Statements 
           ! quick return if possible.
           if( ( m==0 ) .or. ( n==0 ) )return
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = real( b( i, j ),KIND=dp)
              end do
           end do
           l = m*n + 1_${ik}$
           call stdlib${ii}$_dgemm( 'N', 'N', m, n, m, one, a, lda, rwork, m, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = rwork( l+( j-1 )*m+i-1 )
              end do
           end do
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = aimag( b( i, j ) )
              end do
           end do
           call stdlib${ii}$_dgemm( 'N', 'N', m, n, m, one, a, lda, rwork, m, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = cmplx( real( c( i, j ),KIND=dp),rwork( l+( j-1 )*m+i-1 ),KIND=dp)
                           
              end do
           end do
           return
     end subroutine stdlib${ii}$_zlarcm

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ci}$larcm( m, n, a, lda, b, ldb, c, ldc, rwork )
     !! ZLARCM: performs a very simple matrix-matrix multiplication:
     !! C := A * B,
     !! where A is M by M and real; B is M by N and complex;
     !! C is M by N and complex.
        ! -- lapack auxiliary routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_${ck}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           integer(${ik}$), intent(in) :: lda, ldb, ldc, m, n
           ! Array Arguments 
           real(${ck}$), intent(in) :: a(lda,*)
           real(${ck}$), intent(out) :: rwork(*)
           complex(${ck}$), intent(in) :: b(ldb,*)
           complex(${ck}$), intent(out) :: c(ldc,*)
        ! =====================================================================
           
           ! Local Scalars 
           integer(${ik}$) :: i, j, l
           ! Intrinsic Functions 
           ! Executable Statements 
           ! quick return if possible.
           if( ( m==0 ) .or. ( n==0 ) )return
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = real( b( i, j ),KIND=${ck}$)
              end do
           end do
           l = m*n + 1_${ik}$
           call stdlib${ii}$_${c2ri(ci)}$gemm( 'N', 'N', m, n, m, one, a, lda, rwork, m, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = rwork( l+( j-1 )*m+i-1 )
              end do
           end do
           do j = 1, n
              do i = 1, m
                 rwork( ( j-1 )*m+i ) = aimag( b( i, j ) )
              end do
           end do
           call stdlib${ii}$_${c2ri(ci)}$gemm( 'N', 'N', m, n, m, one, a, lda, rwork, m, zero,rwork( l ), m )
                     
           do j = 1, n
              do i = 1, m
                 c( i, j ) = cmplx( real( c( i, j ),KIND=${ck}$),rwork( l+( j-1 )*m+i-1 ),KIND=${ck}$)
                           
              end do
           end do
           return
     end subroutine stdlib${ii}$_${ci}$larcm

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_chfrk( transr, uplo, trans, n, k, alpha, a, lda, beta,c )
     !! Level 3 BLAS like routine for C in RFP Format.
     !! CHFRK performs one of the Hermitian rank--k operations
     !! C := alpha*A*A**H + beta*C,
     !! or
     !! C := alpha*A**H*A + beta*C,
     !! where alpha and beta are real scalars, C is an n--by--n Hermitian
     !! matrix and A is an n--by--k matrix in the first case and a k--by--n
     !! matrix in the second case.
        ! -- lapack computational routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           real(sp), intent(in) :: alpha, beta
           integer(${ik}$), intent(in) :: k, lda, n
           character, intent(in) :: trans, transr, uplo
           ! Array Arguments 
           complex(sp), intent(in) :: a(lda,*)
           complex(sp), intent(inout) :: c(*)
        ! =====================================================================
           
           
           ! Local Scalars 
           logical(lk) :: lower, normaltransr, nisodd, notrans
           integer(${ik}$) :: info, nrowa, j, nk, n1, n2
           complex(sp) :: calpha, cbeta
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input parameters.
           info = 0_${ik}$
           normaltransr = stdlib_lsame( transr, 'N' )
           lower = stdlib_lsame( uplo, 'L' )
           notrans = stdlib_lsame( trans, 'N' )
           if( notrans ) then
              nrowa = n
           else
              nrowa = k
           end if
           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( .not.notrans .and. .not.stdlib_lsame( trans, 'C' ) ) then
              info = -3_${ik}$
           else if( n<0_${ik}$ ) then
              info = -4_${ik}$
           else if( k<0_${ik}$ ) then
              info = -5_${ik}$
           else if( lda<max( 1_${ik}$, nrowa ) ) then
              info = -8_${ik}$
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'CHFRK ', -info )
              return
           end if
           ! quick return if possible.
           ! the quick return case: ((alpha==0).and.(beta/=zero)) is not
           ! done (it is in stdlib${ii}$_cherk for example) and left in the general case.
           if( ( n==0_${ik}$ ) .or. ( ( ( alpha==zero ) .or. ( k==0_${ik}$ ) ) .and.( beta==one ) ) )&
                     return
           if( ( alpha==zero ) .and. ( beta==zero ) ) then
              do j = 1, ( ( n*( n+1 ) ) / 2 )
                 c( j ) = czero
              end do
              return
           end if
           calpha = cmplx( alpha, zero,KIND=sp)
           cbeta = cmplx( beta, zero,KIND=sp)
           ! c is n-by-n.
           ! if n is odd, set nisodd = .true., and n1 and n2.
           ! if n is even, nisodd = .false., and nk.
           if( mod( n, 2_${ik}$ )==0_${ik}$ ) then
              nisodd = .false.
              nk = n / 2_${ik}$
           else
              nisodd = .true.
              if( lower ) then
                 n2 = n / 2_${ik}$
                 n1 = n - n2
              else
                 n1 = n / 2_${ik}$
                 n2 = n - n1
              end if
           end if
           if( nisodd ) then
              ! n is odd
              if( normaltransr ) then
                 ! n is odd and transr = 'n'
                 if( lower ) then
                    ! n is odd, transr = 'n', and uplo = 'l'
                    if( notrans ) then
                       ! n is odd, transr = 'n', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_cherk( 'L', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n )
                                 
                       call stdlib${ii}$_cherk( 'U', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( n+1 )&
                                 , n )
                       call stdlib${ii}$_cgemm( 'N', 'C', n2, n1, k, calpha, a( n1+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( n1+1 ), n )
                    else
                       ! n is odd, transr = 'n', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_cherk( 'L', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n )
                                 
                       call stdlib${ii}$_cherk( 'U', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( n+1 )&
                                 , n )
                       call stdlib${ii}$_cgemm( 'C', 'N', n2, n1, k, calpha, a( 1_${ik}$, n1+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( n1+1 ), n )
                    end if
                 else
                    ! n is odd, transr = 'n', and uplo = 'u'
                    if( notrans ) then
                       ! n is odd, transr = 'n', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_cherk( 'L', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2+1 ), &
                                 n )
                       call stdlib${ii}$_cherk( 'U', 'N', n2, k, alpha, a( n2, 1_${ik}$ ), lda,beta, c( n1+1 ),&
                                  n )
                       call stdlib${ii}$_cgemm( 'N', 'C', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( n2, 1_${ik}$ ), &
                                 lda, cbeta, c( 1_${ik}$ ), n )
                    else
                       ! n is odd, transr = 'n', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_cherk( 'L', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2+1 ), &
                                 n )
                       call stdlib${ii}$_cherk( 'U', 'C', n2, k, alpha, a( 1_${ik}$, n2 ), lda,beta, c( n1+1 ),&
                                  n )
                       call stdlib${ii}$_cgemm( 'C', 'N', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, n2 ), &
                                 lda, cbeta, c( 1_${ik}$ ), n )
                    end if
                 end if
              else
                 ! n is odd, and transr = 'c'
                 if( lower ) then
                    ! n is odd, transr = 'c', and uplo = 'l'
                    if( notrans ) then
                       ! n is odd, transr = 'c', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_cherk( 'U', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n1 &
                                 )
                       call stdlib${ii}$_cherk( 'L', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), &
                                 n1 )
                       call stdlib${ii}$_cgemm( 'N', 'C', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( n1+1, 1_${ik}$ )&
                                 , lda, cbeta,c( n1*n1+1 ), n1 )
                    else
                       ! n is odd, transr = 'c', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_cherk( 'U', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n1 &
                                 )
                       call stdlib${ii}$_cherk( 'L', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( 2_${ik}$ ), &
                                 n1 )
                       call stdlib${ii}$_cgemm( 'C', 'N', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, n1+1 )&
                                 , lda, cbeta,c( n1*n1+1 ), n1 )
                    end if
                 else
                    ! n is odd, transr = 'c', and uplo = 'u'
                    if( notrans ) then
                       ! n is odd, transr = 'c', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_cherk( 'U', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2*n2+1 &
                                 ), n2 )
                       call stdlib${ii}$_cherk( 'L', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( &
                                 n1*n2+1 ), n2 )
                       call stdlib${ii}$_cgemm( 'N', 'C', n2, n1, k, calpha, a( n1+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), n2 )
                    else
                       ! n is odd, transr = 'c', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_cherk( 'U', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2*n2+1 &
                                 ), n2 )
                       call stdlib${ii}$_cherk( 'L', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( &
                                 n1*n2+1 ), n2 )
                       call stdlib${ii}$_cgemm( 'C', 'N', n2, n1, k, calpha, a( 1_${ik}$, n1+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), n2 )
                    end if
                 end if
              end if
           else
              ! n is even
              if( normaltransr ) then
                 ! n is even and transr = 'n'
                 if( lower ) then
                    ! n is even, transr = 'n', and uplo = 'l'
                    if( notrans ) then
                       ! n is even, transr = 'n', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_cherk( 'L', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), n+&
                                 1_${ik}$ )
                       call stdlib${ii}$_cherk( 'U', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), &
                                 n+1 )
                       call stdlib${ii}$_cgemm( 'N', 'C', nk, nk, k, calpha, a( nk+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( nk+2 ),n+1 )
                    else
                       ! n is even, transr = 'n', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_cherk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), n+&
                                 1_${ik}$ )
                       call stdlib${ii}$_cherk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( 1_${ik}$ ), &
                                 n+1 )
                       call stdlib${ii}$_cgemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, nk+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( nk+2 ),n+1 )
                    end if
                 else
                    ! n is even, transr = 'n', and uplo = 'u'
                    if( notrans ) then
                       ! n is even, transr = 'n', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_cherk( 'L', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+2 ), &
                                 n+1 )
                       call stdlib${ii}$_cherk( 'U', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( nk+1 &
                                 ), n+1 )
                       call stdlib${ii}$_cgemm( 'N', 'C', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( nk+1, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ),n+1 )
                    else
                       ! n is even, transr = 'n', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_cherk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+2 ), &
                                 n+1 )
                       call stdlib${ii}$_cherk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( nk+1 &
                                 ), n+1 )
                       call stdlib${ii}$_cgemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, nk+1 )&
                                 , lda, cbeta, c( 1_${ik}$ ),n+1 )
                    end if
                 end if
              else
                 ! n is even, and transr = 'c'
                 if( lower ) then
                    ! n is even, transr = 'c', and uplo = 'l'
                    if( notrans ) then
                       ! n is even, transr = 'c', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_cherk( 'U', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+1 ), &
                                 nk )
                       call stdlib${ii}$_cherk( 'L', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), &
                                 nk )
                       call stdlib${ii}$_cgemm( 'N', 'C', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( nk+1, 1_${ik}$ )&
                                 , lda, cbeta,c( ( ( nk+1 )*nk )+1_${ik}$ ), nk )
                    else
                       ! n is even, transr = 'c', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_cherk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+1 ), &
                                 nk )
                       call stdlib${ii}$_cherk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( 1_${ik}$ ), &
                                 nk )
                       call stdlib${ii}$_cgemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, nk+1 )&
                                 , lda, cbeta,c( ( ( nk+1 )*nk )+1_${ik}$ ), nk )
                    end if
                 else
                    ! n is even, transr = 'c', and uplo = 'u'
                    if( notrans ) then
                       ! n is even, transr = 'c', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_cherk( 'U', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk*( nk+&
                                 1_${ik}$ )+1_${ik}$ ), nk )
                       call stdlib${ii}$_cherk( 'L', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( &
                                 nk*nk+1 ), nk )
                       call stdlib${ii}$_cgemm( 'N', 'C', nk, nk, k, calpha, a( nk+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), nk )
                    else
                       ! n is even, transr = 'c', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_cherk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk*( nk+&
                                 1_${ik}$ )+1_${ik}$ ), nk )
                       call stdlib${ii}$_cherk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( &
                                 nk*nk+1 ), nk )
                       call stdlib${ii}$_cgemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, nk+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), nk )
                    end if
                 end if
              end if
           end if
           return
     end subroutine stdlib${ii}$_chfrk

     pure module subroutine stdlib${ii}$_zhfrk( transr, uplo, trans, n, k, alpha, a, lda, beta,c )
     !! Level 3 BLAS like routine for C in RFP Format.
     !! ZHFRK performs one of the Hermitian rank--k operations
     !! C := alpha*A*A**H + beta*C,
     !! or
     !! C := alpha*A**H*A + beta*C,
     !! where alpha and beta are real scalars, C is an n--by--n Hermitian
     !! matrix and A is an n--by--k matrix in the first case and a k--by--n
     !! matrix in the second case.
        ! -- lapack computational routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_dp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           real(dp), intent(in) :: alpha, beta
           integer(${ik}$), intent(in) :: k, lda, n
           character, intent(in) :: trans, transr, uplo
           ! Array Arguments 
           complex(dp), intent(in) :: a(lda,*)
           complex(dp), intent(inout) :: c(*)
        ! =====================================================================
           
           
           ! Local Scalars 
           logical(lk) :: lower, normaltransr, nisodd, notrans
           integer(${ik}$) :: info, nrowa, j, nk, n1, n2
           complex(dp) :: calpha, cbeta
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input parameters.
           info = 0_${ik}$
           normaltransr = stdlib_lsame( transr, 'N' )
           lower = stdlib_lsame( uplo, 'L' )
           notrans = stdlib_lsame( trans, 'N' )
           if( notrans ) then
              nrowa = n
           else
              nrowa = k
           end if
           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( .not.notrans .and. .not.stdlib_lsame( trans, 'C' ) ) then
              info = -3_${ik}$
           else if( n<0_${ik}$ ) then
              info = -4_${ik}$
           else if( k<0_${ik}$ ) then
              info = -5_${ik}$
           else if( lda<max( 1_${ik}$, nrowa ) ) then
              info = -8_${ik}$
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'ZHFRK ', -info )
              return
           end if
           ! quick return if possible.
           ! the quick return case: ((alpha==0).and.(beta/=zero)) is not
           ! done (it is in stdlib${ii}$_zherk for example) and left in the general case.
           if( ( n==0_${ik}$ ) .or. ( ( ( alpha==zero ) .or. ( k==0_${ik}$ ) ) .and.( beta==one ) ) )&
                     return
           if( ( alpha==zero ) .and. ( beta==zero ) ) then
              do j = 1, ( ( n*( n+1 ) ) / 2 )
                 c( j ) = czero
              end do
              return
           end if
           calpha = cmplx( alpha, zero,KIND=dp)
           cbeta = cmplx( beta, zero,KIND=dp)
           ! c is n-by-n.
           ! if n is odd, set nisodd = .true., and n1 and n2.
           ! if n is even, nisodd = .false., and nk.
           if( mod( n, 2_${ik}$ )==0_${ik}$ ) then
              nisodd = .false.
              nk = n / 2_${ik}$
           else
              nisodd = .true.
              if( lower ) then
                 n2 = n / 2_${ik}$
                 n1 = n - n2
              else
                 n1 = n / 2_${ik}$
                 n2 = n - n1
              end if
           end if
           if( nisodd ) then
              ! n is odd
              if( normaltransr ) then
                 ! n is odd and transr = 'n'
                 if( lower ) then
                    ! n is odd, transr = 'n', and uplo = 'l'
                    if( notrans ) then
                       ! n is odd, transr = 'n', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_zherk( 'L', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n )
                                 
                       call stdlib${ii}$_zherk( 'U', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( n+1 )&
                                 , n )
                       call stdlib${ii}$_zgemm( 'N', 'C', n2, n1, k, calpha, a( n1+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( n1+1 ), n )
                    else
                       ! n is odd, transr = 'n', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_zherk( 'L', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n )
                                 
                       call stdlib${ii}$_zherk( 'U', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( n+1 )&
                                 , n )
                       call stdlib${ii}$_zgemm( 'C', 'N', n2, n1, k, calpha, a( 1_${ik}$, n1+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( n1+1 ), n )
                    end if
                 else
                    ! n is odd, transr = 'n', and uplo = 'u'
                    if( notrans ) then
                       ! n is odd, transr = 'n', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_zherk( 'L', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2+1 ), &
                                 n )
                       call stdlib${ii}$_zherk( 'U', 'N', n2, k, alpha, a( n2, 1_${ik}$ ), lda,beta, c( n1+1 ),&
                                  n )
                       call stdlib${ii}$_zgemm( 'N', 'C', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( n2, 1_${ik}$ ), &
                                 lda, cbeta, c( 1_${ik}$ ), n )
                    else
                       ! n is odd, transr = 'n', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_zherk( 'L', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2+1 ), &
                                 n )
                       call stdlib${ii}$_zherk( 'U', 'C', n2, k, alpha, a( 1_${ik}$, n2 ), lda,beta, c( n1+1 ),&
                                  n )
                       call stdlib${ii}$_zgemm( 'C', 'N', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, n2 ), &
                                 lda, cbeta, c( 1_${ik}$ ), n )
                    end if
                 end if
              else
                 ! n is odd, and transr = 'c'
                 if( lower ) then
                    ! n is odd, transr = 'c', and uplo = 'l'
                    if( notrans ) then
                       ! n is odd, transr = 'c', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_zherk( 'U', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n1 &
                                 )
                       call stdlib${ii}$_zherk( 'L', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), &
                                 n1 )
                       call stdlib${ii}$_zgemm( 'N', 'C', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( n1+1, 1_${ik}$ )&
                                 , lda, cbeta,c( n1*n1+1 ), n1 )
                    else
                       ! n is odd, transr = 'c', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_zherk( 'U', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n1 &
                                 )
                       call stdlib${ii}$_zherk( 'L', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( 2_${ik}$ ), &
                                 n1 )
                       call stdlib${ii}$_zgemm( 'C', 'N', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, n1+1 )&
                                 , lda, cbeta,c( n1*n1+1 ), n1 )
                    end if
                 else
                    ! n is odd, transr = 'c', and uplo = 'u'
                    if( notrans ) then
                       ! n is odd, transr = 'c', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_zherk( 'U', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2*n2+1 &
                                 ), n2 )
                       call stdlib${ii}$_zherk( 'L', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( &
                                 n1*n2+1 ), n2 )
                       call stdlib${ii}$_zgemm( 'N', 'C', n2, n1, k, calpha, a( n1+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), n2 )
                    else
                       ! n is odd, transr = 'c', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_zherk( 'U', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2*n2+1 &
                                 ), n2 )
                       call stdlib${ii}$_zherk( 'L', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( &
                                 n1*n2+1 ), n2 )
                       call stdlib${ii}$_zgemm( 'C', 'N', n2, n1, k, calpha, a( 1_${ik}$, n1+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), n2 )
                    end if
                 end if
              end if
           else
              ! n is even
              if( normaltransr ) then
                 ! n is even and transr = 'n'
                 if( lower ) then
                    ! n is even, transr = 'n', and uplo = 'l'
                    if( notrans ) then
                       ! n is even, transr = 'n', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_zherk( 'L', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), n+&
                                 1_${ik}$ )
                       call stdlib${ii}$_zherk( 'U', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), &
                                 n+1 )
                       call stdlib${ii}$_zgemm( 'N', 'C', nk, nk, k, calpha, a( nk+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( nk+2 ),n+1 )
                    else
                       ! n is even, transr = 'n', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_zherk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), n+&
                                 1_${ik}$ )
                       call stdlib${ii}$_zherk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( 1_${ik}$ ), &
                                 n+1 )
                       call stdlib${ii}$_zgemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, nk+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( nk+2 ),n+1 )
                    end if
                 else
                    ! n is even, transr = 'n', and uplo = 'u'
                    if( notrans ) then
                       ! n is even, transr = 'n', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_zherk( 'L', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+2 ), &
                                 n+1 )
                       call stdlib${ii}$_zherk( 'U', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( nk+1 &
                                 ), n+1 )
                       call stdlib${ii}$_zgemm( 'N', 'C', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( nk+1, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ),n+1 )
                    else
                       ! n is even, transr = 'n', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_zherk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+2 ), &
                                 n+1 )
                       call stdlib${ii}$_zherk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( nk+1 &
                                 ), n+1 )
                       call stdlib${ii}$_zgemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, nk+1 )&
                                 , lda, cbeta, c( 1_${ik}$ ),n+1 )
                    end if
                 end if
              else
                 ! n is even, and transr = 'c'
                 if( lower ) then
                    ! n is even, transr = 'c', and uplo = 'l'
                    if( notrans ) then
                       ! n is even, transr = 'c', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_zherk( 'U', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+1 ), &
                                 nk )
                       call stdlib${ii}$_zherk( 'L', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), &
                                 nk )
                       call stdlib${ii}$_zgemm( 'N', 'C', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( nk+1, 1_${ik}$ )&
                                 , lda, cbeta,c( ( ( nk+1 )*nk )+1_${ik}$ ), nk )
                    else
                       ! n is even, transr = 'c', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_zherk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+1 ), &
                                 nk )
                       call stdlib${ii}$_zherk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( 1_${ik}$ ), &
                                 nk )
                       call stdlib${ii}$_zgemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, nk+1 )&
                                 , lda, cbeta,c( ( ( nk+1 )*nk )+1_${ik}$ ), nk )
                    end if
                 else
                    ! n is even, transr = 'c', and uplo = 'u'
                    if( notrans ) then
                       ! n is even, transr = 'c', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_zherk( 'U', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk*( nk+&
                                 1_${ik}$ )+1_${ik}$ ), nk )
                       call stdlib${ii}$_zherk( 'L', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( &
                                 nk*nk+1 ), nk )
                       call stdlib${ii}$_zgemm( 'N', 'C', nk, nk, k, calpha, a( nk+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), nk )
                    else
                       ! n is even, transr = 'c', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_zherk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk*( nk+&
                                 1_${ik}$ )+1_${ik}$ ), nk )
                       call stdlib${ii}$_zherk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( &
                                 nk*nk+1 ), nk )
                       call stdlib${ii}$_zgemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, nk+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), nk )
                    end if
                 end if
              end if
           end if
           return
     end subroutine stdlib${ii}$_zhfrk

#:for ck,ct,ci in CMPLX_KINDS_TYPES
#:if not ck in ["sp","dp"]
     pure module subroutine stdlib${ii}$_${ci}$hfrk( transr, uplo, trans, n, k, alpha, a, lda, beta,c )
     !! Level 3 BLAS like routine for C in RFP Format.
     !! ZHFRK: performs one of the Hermitian rank--k operations
     !! C := alpha*A*A**H + beta*C,
     !! or
     !! C := alpha*A**H*A + beta*C,
     !! where alpha and beta are real scalars, C is an n--by--n Hermitian
     !! matrix and A is an n--by--k matrix in the first case and a k--by--n
     !! matrix in the second case.
        ! -- lapack computational routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_${ck}$, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           real(${ck}$), intent(in) :: alpha, beta
           integer(${ik}$), intent(in) :: k, lda, n
           character, intent(in) :: trans, transr, uplo
           ! Array Arguments 
           complex(${ck}$), intent(in) :: a(lda,*)
           complex(${ck}$), intent(inout) :: c(*)
        ! =====================================================================
           
           
           ! Local Scalars 
           logical(lk) :: lower, normaltransr, nisodd, notrans
           integer(${ik}$) :: info, nrowa, j, nk, n1, n2
           complex(${ck}$) :: calpha, cbeta
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input parameters.
           info = 0_${ik}$
           normaltransr = stdlib_lsame( transr, 'N' )
           lower = stdlib_lsame( uplo, 'L' )
           notrans = stdlib_lsame( trans, 'N' )
           if( notrans ) then
              nrowa = n
           else
              nrowa = k
           end if
           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( .not.notrans .and. .not.stdlib_lsame( trans, 'C' ) ) then
              info = -3_${ik}$
           else if( n<0_${ik}$ ) then
              info = -4_${ik}$
           else if( k<0_${ik}$ ) then
              info = -5_${ik}$
           else if( lda<max( 1_${ik}$, nrowa ) ) then
              info = -8_${ik}$
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'ZHFRK ', -info )
              return
           end if
           ! quick return if possible.
           ! the quick return case: ((alpha==0).and.(beta/=zero)) is not
           ! done (it is in stdlib${ii}$_${ci}$herk for example) and left in the general case.
           if( ( n==0_${ik}$ ) .or. ( ( ( alpha==zero ) .or. ( k==0_${ik}$ ) ) .and.( beta==one ) ) )&
                     return
           if( ( alpha==zero ) .and. ( beta==zero ) ) then
              do j = 1, ( ( n*( n+1 ) ) / 2 )
                 c( j ) = czero
              end do
              return
           end if
           calpha = cmplx( alpha, zero,KIND=${ck}$)
           cbeta = cmplx( beta, zero,KIND=${ck}$)
           ! c is n-by-n.
           ! if n is odd, set nisodd = .true., and n1 and n2.
           ! if n is even, nisodd = .false., and nk.
           if( mod( n, 2_${ik}$ )==0_${ik}$ ) then
              nisodd = .false.
              nk = n / 2_${ik}$
           else
              nisodd = .true.
              if( lower ) then
                 n2 = n / 2_${ik}$
                 n1 = n - n2
              else
                 n1 = n / 2_${ik}$
                 n2 = n - n1
              end if
           end if
           if( nisodd ) then
              ! n is odd
              if( normaltransr ) then
                 ! n is odd and transr = 'n'
                 if( lower ) then
                    ! n is odd, transr = 'n', and uplo = 'l'
                    if( notrans ) then
                       ! n is odd, transr = 'n', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_${ci}$herk( 'L', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n )
                                 
                       call stdlib${ii}$_${ci}$herk( 'U', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( n+1 )&
                                 , n )
                       call stdlib${ii}$_${ci}$gemm( 'N', 'C', n2, n1, k, calpha, a( n1+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( n1+1 ), n )
                    else
                       ! n is odd, transr = 'n', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_${ci}$herk( 'L', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n )
                                 
                       call stdlib${ii}$_${ci}$herk( 'U', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( n+1 )&
                                 , n )
                       call stdlib${ii}$_${ci}$gemm( 'C', 'N', n2, n1, k, calpha, a( 1_${ik}$, n1+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( n1+1 ), n )
                    end if
                 else
                    ! n is odd, transr = 'n', and uplo = 'u'
                    if( notrans ) then
                       ! n is odd, transr = 'n', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_${ci}$herk( 'L', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2+1 ), &
                                 n )
                       call stdlib${ii}$_${ci}$herk( 'U', 'N', n2, k, alpha, a( n2, 1_${ik}$ ), lda,beta, c( n1+1 ),&
                                  n )
                       call stdlib${ii}$_${ci}$gemm( 'N', 'C', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( n2, 1_${ik}$ ), &
                                 lda, cbeta, c( 1_${ik}$ ), n )
                    else
                       ! n is odd, transr = 'n', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_${ci}$herk( 'L', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2+1 ), &
                                 n )
                       call stdlib${ii}$_${ci}$herk( 'U', 'C', n2, k, alpha, a( 1_${ik}$, n2 ), lda,beta, c( n1+1 ),&
                                  n )
                       call stdlib${ii}$_${ci}$gemm( 'C', 'N', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, n2 ), &
                                 lda, cbeta, c( 1_${ik}$ ), n )
                    end if
                 end if
              else
                 ! n is odd, and transr = 'c'
                 if( lower ) then
                    ! n is odd, transr = 'c', and uplo = 'l'
                    if( notrans ) then
                       ! n is odd, transr = 'c', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_${ci}$herk( 'U', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n1 &
                                 )
                       call stdlib${ii}$_${ci}$herk( 'L', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), &
                                 n1 )
                       call stdlib${ii}$_${ci}$gemm( 'N', 'C', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( n1+1, 1_${ik}$ )&
                                 , lda, cbeta,c( n1*n1+1 ), n1 )
                    else
                       ! n is odd, transr = 'c', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_${ci}$herk( 'U', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), n1 &
                                 )
                       call stdlib${ii}$_${ci}$herk( 'L', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( 2_${ik}$ ), &
                                 n1 )
                       call stdlib${ii}$_${ci}$gemm( 'C', 'N', n1, n2, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, n1+1 )&
                                 , lda, cbeta,c( n1*n1+1 ), n1 )
                    end if
                 else
                    ! n is odd, transr = 'c', and uplo = 'u'
                    if( notrans ) then
                       ! n is odd, transr = 'c', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_${ci}$herk( 'U', 'N', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2*n2+1 &
                                 ), n2 )
                       call stdlib${ii}$_${ci}$herk( 'L', 'N', n2, k, alpha, a( n1+1, 1_${ik}$ ), lda,beta, c( &
                                 n1*n2+1 ), n2 )
                       call stdlib${ii}$_${ci}$gemm( 'N', 'C', n2, n1, k, calpha, a( n1+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), n2 )
                    else
                       ! n is odd, transr = 'c', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_${ci}$herk( 'U', 'C', n1, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( n2*n2+1 &
                                 ), n2 )
                       call stdlib${ii}$_${ci}$herk( 'L', 'C', n2, k, alpha, a( 1_${ik}$, n1+1 ), lda,beta, c( &
                                 n1*n2+1 ), n2 )
                       call stdlib${ii}$_${ci}$gemm( 'C', 'N', n2, n1, k, calpha, a( 1_${ik}$, n1+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), n2 )
                    end if
                 end if
              end if
           else
              ! n is even
              if( normaltransr ) then
                 ! n is even and transr = 'n'
                 if( lower ) then
                    ! n is even, transr = 'n', and uplo = 'l'
                    if( notrans ) then
                       ! n is even, transr = 'n', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_${ci}$herk( 'L', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), n+&
                                 1_${ik}$ )
                       call stdlib${ii}$_${ci}$herk( 'U', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), &
                                 n+1 )
                       call stdlib${ii}$_${ci}$gemm( 'N', 'C', nk, nk, k, calpha, a( nk+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( nk+2 ),n+1 )
                    else
                       ! n is even, transr = 'n', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_${ci}$herk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( 2_${ik}$ ), n+&
                                 1_${ik}$ )
                       call stdlib${ii}$_${ci}$herk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( 1_${ik}$ ), &
                                 n+1 )
                       call stdlib${ii}$_${ci}$gemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, nk+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( nk+2 ),n+1 )
                    end if
                 else
                    ! n is even, transr = 'n', and uplo = 'u'
                    if( notrans ) then
                       ! n is even, transr = 'n', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_${ci}$herk( 'L', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+2 ), &
                                 n+1 )
                       call stdlib${ii}$_${ci}$herk( 'U', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( nk+1 &
                                 ), n+1 )
                       call stdlib${ii}$_${ci}$gemm( 'N', 'C', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( nk+1, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ),n+1 )
                    else
                       ! n is even, transr = 'n', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_${ci}$herk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+2 ), &
                                 n+1 )
                       call stdlib${ii}$_${ci}$herk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( nk+1 &
                                 ), n+1 )
                       call stdlib${ii}$_${ci}$gemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, nk+1 )&
                                 , lda, cbeta, c( 1_${ik}$ ),n+1 )
                    end if
                 end if
              else
                 ! n is even, and transr = 'c'
                 if( lower ) then
                    ! n is even, transr = 'c', and uplo = 'l'
                    if( notrans ) then
                       ! n is even, transr = 'c', uplo = 'l', and trans = 'n'
                       call stdlib${ii}$_${ci}$herk( 'U', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+1 ), &
                                 nk )
                       call stdlib${ii}$_${ci}$herk( 'L', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( 1_${ik}$ ), &
                                 nk )
                       call stdlib${ii}$_${ci}$gemm( 'N', 'C', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( nk+1, 1_${ik}$ )&
                                 , lda, cbeta,c( ( ( nk+1 )*nk )+1_${ik}$ ), nk )
                    else
                       ! n is even, transr = 'c', uplo = 'l', and trans = 'c'
                       call stdlib${ii}$_${ci}$herk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk+1 ), &
                                 nk )
                       call stdlib${ii}$_${ci}$herk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( 1_${ik}$ ), &
                                 nk )
                       call stdlib${ii}$_${ci}$gemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, 1_${ik}$ ),lda, a( 1_${ik}$, nk+1 )&
                                 , lda, cbeta,c( ( ( nk+1 )*nk )+1_${ik}$ ), nk )
                    end if
                 else
                    ! n is even, transr = 'c', and uplo = 'u'
                    if( notrans ) then
                       ! n is even, transr = 'c', uplo = 'u', and trans = 'n'
                       call stdlib${ii}$_${ci}$herk( 'U', 'N', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk*( nk+&
                                 1_${ik}$ )+1_${ik}$ ), nk )
                       call stdlib${ii}$_${ci}$herk( 'L', 'N', nk, k, alpha, a( nk+1, 1_${ik}$ ), lda,beta, c( &
                                 nk*nk+1 ), nk )
                       call stdlib${ii}$_${ci}$gemm( 'N', 'C', nk, nk, k, calpha, a( nk+1, 1_${ik}$ ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), nk )
                    else
                       ! n is even, transr = 'c', uplo = 'u', and trans = 'c'
                       call stdlib${ii}$_${ci}$herk( 'U', 'C', nk, k, alpha, a( 1_${ik}$, 1_${ik}$ ), lda,beta, c( nk*( nk+&
                                 1_${ik}$ )+1_${ik}$ ), nk )
                       call stdlib${ii}$_${ci}$herk( 'L', 'C', nk, k, alpha, a( 1_${ik}$, nk+1 ), lda,beta, c( &
                                 nk*nk+1 ), nk )
                       call stdlib${ii}$_${ci}$gemm( 'C', 'N', nk, nk, k, calpha, a( 1_${ik}$, nk+1 ),lda, a( 1_${ik}$, 1_${ik}$ )&
                                 , lda, cbeta, c( 1_${ik}$ ), nk )
                    end if
                 end if
              end if
           end if
           return
     end subroutine stdlib${ii}$_${ci}$hfrk

#:endif
#:endfor



     pure module subroutine stdlib${ii}$_stfsm( transr, side, uplo, trans, diag, m, n, alpha, a,b, ldb )
     !! Level 3 BLAS like routine for A in RFP Format.
     !! STFSM solves the matrix equation
     !! op( A )*X = alpha*B  or  X*op( A ) = alpha*B
     !! where alpha is a scalar, X and B are m by n matrices, A is a unit, or
     !! non-unit,  upper or lower triangular matrix  and  op( A )  is one  of
     !! op( A ) = A   or   op( A ) = A**T.
     !! A is in Rectangular Full Packed (RFP) Format.
     !! The matrix X is overwritten on B.
               
        ! -- lapack computational routine --
        ! -- lapack is a software package provided by univ. of tennessee,    --
        ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..--
           use stdlib_blas_constants_sp, only: negone, zero, half, one, two, three, four, eight, ten, czero, chalf, cone, cnegone
           ! Scalar Arguments 
           character, intent(in) :: transr, diag, side, trans, uplo
           integer(${ik}$), intent(in) :: ldb, m, n
           real(sp), intent(in) :: alpha
           ! Array Arguments 
           real(sp), intent(in) :: a(0_${ik}$:*)
           real(sp), intent(inout) :: b(0_${ik}$:ldb-1,0_${ik}$:*)
        ! =====================================================================
           
           ! Local Scalars 
           logical(lk) :: lower, lside, misodd, nisodd, normaltransr, notrans
           integer(${ik}$) :: m1, m2, n1, n2, k, info, i, j
           ! Intrinsic Functions 
           ! Executable Statements 
           ! test the input parameters.
           info = 0_${ik}$
           normaltransr = stdlib_lsame( transr, 'N' )
           lside = stdlib_lsame( side, 'L' )
           lower = stdlib_lsame( uplo, 'L' )
           notrans = stdlib_lsame( trans, 'N' )
           if( .not.normaltransr .and. .not.stdlib_lsame( transr, 'T' ) ) then
              info = -1_${ik}$
           else if( .not.lside .and. .not.stdlib_lsame( side, 'R' ) ) then
              info = -2_${ik}$
           else if( .not.lower .and. .not.stdlib_lsame( uplo, 'U' ) ) then
              info = -3_${ik}$
           else if( .not.notrans .and. .not.stdlib_lsame( trans, 'T' ) ) then
              info = -4_${ik}$
           else if( .not.stdlib_lsame( diag, 'N' ) .and. .not.stdlib_lsame( diag, 'U' ) )&
                     then
              info = -5_${ik}$
           else if( m<0_${ik}$ ) then
              info = -6_${ik}$
           else if( n<0_${ik}$ ) then
              info = -7_${ik}$
           else if( ldb<max( 1_${ik}$, m ) ) then
              info = -11_${ik}$
           end if
           if( info/=0_${ik}$ ) then
              call stdlib${ii}$_xerbla( 'STFSM ', -info )
              return
           end if
           ! quick return when ( (n==0).or.(m==0) )
           if( ( m==0 ) .or. ( n==0 ) )return
           ! quick return when alpha==(0e+0_sp)
           if( alpha==zero ) then
              do j = 0, n - 1
                 do i = 0, m - 1
                    b( i, j ) = zero
                 end do
              end do
              return
           end if
           if( lside ) then
              ! side = 'l'
              ! a is m-by-m.
              ! if m is odd, set nisodd = .true., and m1 and m2.
              ! if m is even, nisodd = .false., and m.
              if( mod( m, 2_${ik}$ )==0_${ik}$ ) then
                 misodd = .false.
                 k = m / 2_${ik}$
              else
                 misodd = .true.
                 if( lower ) then
                    m2 = m / 2_${ik}$
                    m1 = m - m2
                 else
                    m1 = m / 2_${ik}$
                    m2 = m - m1
                 end if
              end if
              if( misodd ) then
                 ! side = 'l' and n is odd
                 if( normaltransr ) then
                    ! side = 'l', n is odd, and transr = 'n'
                    if( lower ) then
                       ! side  ='l', n is odd, transr = 'n', and uplo = 'l'
                       if( notrans ) then
                          ! side  ='l', n is odd, transr = 'n', uplo = 'l', and
                          ! trans = 'n'
                          if( m==1_${ik}$ ) then
                             call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, m1, n, alpha,a, m, b, ldb )
                                       
                          else
                             call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, m1, n, alpha,a( 0_${ik}$ ), m, b, &
                                       ldb )
                             call stdlib${ii}$_sgemm( 'N', 'N', m2, n, m1, -one, a( m1 ),m, b, ldb, &
                                       alpha, b( m1, 0_${ik}$ ), ldb )
                             call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, m2, n, one,a( m ), m, b( m1, &
                                       0_${ik}$ ), ldb )
                          end if
                       else
                          ! side  ='l', n is odd, transr = 'n', uplo = 'l', and
                          ! trans = 't'
                          if( m==1_${ik}$ ) then
                             call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, m1, n, alpha,a( 0_${ik}$ ), m, b, &
                                       ldb )
                          else
                             call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, m2, n, alpha,a( m ), m, b( &
                                       m1, 0_${ik}$ ), ldb )
                             call stdlib${ii}$_sgemm( 'T', 'N', m1, n, m2, -one, a( m1 ),m, b( m1, 0_${ik}$ ), &
                                       ldb, alpha, b, ldb )
                             call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, m1, n, one,a( 0_${ik}$ ), m, b, ldb &
                                       )
                          end if
                       end if
                    else
                       ! side  ='l', n is odd, transr = 'n', and uplo = 'u'
                       if( .not.notrans ) then
                          ! side  ='l', n is odd, transr = 'n', uplo = 'u', and
                          ! trans = 'n'
                          call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, m1, n, alpha,a( m2 ), m, b, ldb &
                                    )
                          call stdlib${ii}$_sgemm( 'T', 'N', m2, n, m1, -one, a( 0_${ik}$ ), m,b, ldb, alpha, &
                                    b( m1, 0_${ik}$ ), ldb )
                          call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, m2, n, one,a( m1 ), m, b( m1, 0_${ik}$ &
                                    ), ldb )
                       else
                          ! side  ='l', n is odd, transr = 'n', uplo = 'u', and
                          ! trans = 't'
                          call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, m2, n, alpha,a( m1 ), m, b( m1, &
                                    0_${ik}$ ), ldb )
                          call stdlib${ii}$_sgemm( 'N', 'N', m1, n, m2, -one, a( 0_${ik}$ ), m,b( m1, 0_${ik}$ ), ldb,&
                                     alpha, b, ldb )
                          call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, m1, n, one,a( m2 ), m, b, ldb )
                                    
                       end if
                    end if
                 else
                    ! side = 'l', n is odd, and transr = 't'
                    if( lower ) then
                       ! side  ='l', n is odd, transr = 't', and uplo = 'l'
                       if( notrans ) then
                          ! side  ='l', n is odd, transr = 't', uplo = 'l', and
                          ! trans = 'n'
                          if( m==1_${ik}$ ) then
                             call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, m1, n, alpha,a( 0_${ik}$ ), m1, b, &
                                       ldb )
                          else
                             call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, m1, n, alpha,a( 0_${ik}$ ), m1, b, &
                                       ldb )
                             call stdlib${ii}$_sgemm( 'T', 'N', m2, n, m1, -one,a( m1*m1 ), m1, b, ldb, &
                                       alpha,b( m1, 0_${ik}$ ), ldb )
                             call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, m2, n, one,a( 1_${ik}$ ), m1, b( m1,&
                                        0_${ik}$ ), ldb )
                          end if
                       else
                          ! side  ='l', n is odd, transr = 't', uplo = 'l', and
                          ! trans = 't'
                          if( m==1_${ik}$ ) then
                             call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, m1, n, alpha,a( 0_${ik}$ ), m1, b, &
                                       ldb )
                          else
                             call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, m2, n, alpha,a( 1_${ik}$ ), m1, b( &
                                       m1, 0_${ik}$ ), ldb )
                             call stdlib${ii}$_sgemm( 'N', 'N', m1, n, m2, -one,a( m1*m1 ), m1, b( m1, &
                                       0_${ik}$ ), ldb,alpha, b, ldb )
                             call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, m1, n, one,a( 0_${ik}$ ), m1, b, &
                                       ldb )
                          end if
                       end if
                    else
                       ! side  ='l', n is odd, transr = 't', and uplo = 'u'
                       if( .not.notrans ) then
                          ! side  ='l', n is odd, transr = 't', uplo = 'u', and
                          ! trans = 'n'
                          call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, m1, n, alpha,a( m2*m2 ), m2, b, &
                                    ldb )
                          call stdlib${ii}$_sgemm( 'N', 'N', m2, n, m1, -one, a( 0_${ik}$ ), m2,b, ldb, alpha, &
                                    b( m1, 0_${ik}$ ), ldb )
                          call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, m2, n, one,a( m1*m2 ), m2, b( &
                                    m1, 0_${ik}$ ), ldb )
                       else
                          ! side  ='l', n is odd, transr = 't', uplo = 'u', and
                          ! trans = 't'
                          call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, m2, n, alpha,a( m1*m2 ), m2, b( &
                                    m1, 0_${ik}$ ), ldb )
                          call stdlib${ii}$_sgemm( 'T', 'N', m1, n, m2, -one, a( 0_${ik}$ ), m2,b( m1, 0_${ik}$ ), &
                                    ldb, alpha, b, ldb )
                          call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, m1, n, one,a( m2*m2 ), m2, b, &
                                    ldb )
                       end if
                    end if
                 end if
              else
                 ! side = 'l' and n is even
                 if( normaltransr ) then
                    ! side = 'l', n is even, and transr = 'n'
                    if( lower ) then
                       ! side  ='l', n is even, transr = 'n', and uplo = 'l'
                       if( notrans ) then
                          ! side  ='l', n is even, transr = 'n', uplo = 'l',
                          ! and trans = 'n'
                          call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, k, n, alpha,a( 1_${ik}$ ), m+1, b, ldb &
                                    )
                          call stdlib${ii}$_sgemm( 'N', 'N', k, n, k, -one, a( k+1 ),m+1, b, ldb, alpha,&
                                     b( k, 0_${ik}$ ), ldb )
                          call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, k, n, one,a( 0_${ik}$ ), m+1, b( k, 0_${ik}$ )&
                                    , ldb )
                       else
                          ! side  ='l', n is even, transr = 'n', uplo = 'l',
                          ! and trans = 't'
                          call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, k, n, alpha,a( 0_${ik}$ ), m+1, b( k, &
                                    0_${ik}$ ), ldb )
                          call stdlib${ii}$_sgemm( 'T', 'N', k, n, k, -one, a( k+1 ),m+1, b( k, 0_${ik}$ ), &
                                    ldb, alpha, b, ldb )
                          call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, k, n, one,a( 1_${ik}$ ), m+1, b, ldb )
                                    
                       end if
                    else
                       ! side  ='l', n is even, transr = 'n', and uplo = 'u'
                       if( .not.notrans ) then
                          ! side  ='l', n is even, transr = 'n', uplo = 'u',
                          ! and trans = 'n'
                          call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, k, n, alpha,a( k+1 ), m+1, b, &
                                    ldb )
                          call stdlib${ii}$_sgemm( 'T', 'N', k, n, k, -one, a( 0_${ik}$ ), m+1,b, ldb, alpha, &
                                    b( k, 0_${ik}$ ), ldb )
                          call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, k, n, one,a( k ), m+1, b( k, 0_${ik}$ )&
                                    , ldb )
                       else
                          ! side  ='l', n is even, transr = 'n', uplo = 'u',
                          ! and trans = 't'
                          call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, k, n, alpha,a( k ), m+1, b( k, &
                                    0_${ik}$ ), ldb )
                          call stdlib${ii}$_sgemm( 'N', 'N', k, n, k, -one, a( 0_${ik}$ ), m+1,b( k, 0_${ik}$ ), ldb, &
                                    alpha, b, ldb )
                          call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, k, n, one,a( k+1 ), m+1, b, ldb &
                                    )
                       end if
                    end if
                 else
                    ! side = 'l', n is even, and transr = 't'
                    if( lower ) then
                       ! side  ='l', n is even, transr = 't', and uplo = 'l'
                       if( notrans ) then
                          ! side  ='l', n is even, transr = 't', uplo = 'l',
                          ! and trans = 'n'
                          call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, k, n, alpha,a( k ), k, b, ldb )
                                    
                          call stdlib${ii}$_sgemm( 'T', 'N', k, n, k, -one,a( k*( k+1 ) ), k, b, ldb, &
                                    alpha,b( k, 0_${ik}$ ), ldb )
                          call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, k, n, one,a( 0_${ik}$ ), k, b( k, 0_${ik}$ ), &
                                    ldb )
                       else
                          ! side  ='l', n is even, transr = 't', uplo = 'l',
                          ! and trans = 't'
                          call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, k, n, alpha,a( 0_${ik}$ ), k, b( k, 0_${ik}$ )&
                                    , ldb )
                          call stdlib${ii}$_sgemm( 'N', 'N', k, n, k, -one,a( k*( k+1 ) ), k, b( k, 0_${ik}$ ),&
                                     ldb,alpha, b, ldb )
                          call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, k, n, one,a( k ), k, b, ldb )
                                    
                       end if
                    else
                       ! side  ='l', n is even, transr = 't', and uplo = 'u'
                       if( .not.notrans ) then
                          ! side  ='l', n is even, transr = 't', uplo = 'u',
                          ! and trans = 'n'
                          call stdlib${ii}$_strsm( 'L', 'U', 'T', diag, k, n, alpha,a( k*( k+1 ) ), k, &
                                    b, ldb )
                          call stdlib${ii}$_sgemm( 'N', 'N', k, n, k, -one, a( 0_${ik}$ ), k, b,ldb, alpha, b( &
                                    k, 0_${ik}$ ), ldb )
                          call stdlib${ii}$_strsm( 'L', 'L', 'N', diag, k, n, one,a( k*k ), k, b( k, 0_${ik}$ )&
                                    , ldb )
                       else
                          ! side  ='l', n is even, transr = 't', uplo = 'u',
                          ! and trans = 't'
                          call stdlib${ii}$_strsm( 'L', 'L', 'T', diag, k, n, alpha,a( k*k ), k, b( k, &
                                    0_${ik}$ ), ldb )
                          call stdlib${ii}$_sgemm( 'T', 'N', k, n, k, -one, a( 0_${ik}$ ), k,b( k, 0_${ik}$ ), ldb, &
                                    alpha, b, ldb )
                          call stdlib${ii}$_strsm( 'L', 'U', 'N', diag, k, n, one,a( k*( k+1 ) ), k, b, &
                                    ldb )
                       end if
                    end if
                 end if
              end if
           else
              ! side = 'r'
              ! a is n-by-n.
              ! if n is odd, set nisodd = .true., and n1 and n2.
              ! if n is even, nisodd = .false., and k.
              if( mod( n, 2_${ik}$ )==0_${ik}$ ) then
                 nisodd = .false.
                 k = n / 2_${ik}$
              else
                 nisodd = .true.
                 if( lower ) then
                    n2 = n / 2_${ik}$
                    n1 = n - n2
                 else
                    n1 = n / 2_${ik}$
                    n2 = n - n1
                 end if
              end if
              if( nisodd ) then
                 ! side = 'r' and n is odd
                 if( normaltransr ) then
                    ! side = 'r', n is odd, and transr = 'n'
                    if( lower ) then
                       ! side  ='r', n is odd, transr = 'n', and uplo = 'l'
                       if( notrans ) then
                          ! side  ='r', n is odd, transr = 'n', uplo = 'l', and
                          ! trans = 'n'
                          call stdlib${ii}$_strsm( 'R', 'U', 'T', diag, m, n2, alpha,a( n ), n, b( 0_${ik}$, &
                                    n1 ), ldb )
                          call stdlib${ii}$_sgemm( 'N', 'N', m, n1, n2, -one, b( 0_${ik}$, n1 ),ldb, a( n1 ), &
                                    n, alpha, b( 0_${ik}$, 0_${ik}$ ),ldb )
                          call stdlib${ii}$_strsm( 'R', 'L', 'N', diag, m, n1, one,a( 0_${ik}$ ), n, b( 0_${ik}$, 0_${ik}$ ),&
                                     ldb )
                       else
                          ! side  ='r', n is odd, transr = 'n', uplo = 'l', and
                          ! trans = 't'
                          call stdlib${ii}$_strsm( 'R', 'L', 'T', diag, m, n1, alpha,a( 0_${ik}$ ), n, b( 0_${ik}$, 0_${ik}$ &
                                    ), ldb )
                          call stdlib${ii}$_sgemm( 'N', 'T', m, n2, n1, -one, b( 0_${ik}$, 0_${ik}$ ),ldb, a( n1 ), n,&
                                     alpha, b( 0_${ik}$, n1 ),ldb )
                          call stdlib${ii}$_strsm( 'R', 'U', 'N', diag, m, n2, one,a( n ), n, b( 0_${ik}$, n1 )&
                                    , ldb )
                       end if
                    else
                       ! side  ='r', n is odd, transr = 'n', and uplo = 'u'
                       if( notrans ) then
                          ! side  ='r', n is odd, transr = 'n', uplo = 'u', and
                          ! trans = 'n'
                          call stdlib${ii}$_strsm( 'R', 'L', 'T', diag, m, n1, alpha,a( n2 ), n, b( 0_${ik}$, &
                                    0_${ik}$ ), ldb )
                          call stdlib${ii}$_sgemm( 'N', 'N', m, n2, n1, -one, b( 0_${ik}$, 0_${ik}$ ),ldb, a( 0_${ik}$ ), n, &
                                    alpha, b( 0_${ik}$, n1 ),ldb )
                          call stdlib${ii}$_strsm( 'R', 'U', 'N', diag, m, n2, one,a( n1 ), n, b( 0_${ik}$, n1 &
                                    ), ldb )
                       else
                          ! side  ='r', n is odd, transr = 'n', uplo = 'u', and
                          ! trans = 't'
                          call stdlib${ii}$_strsm( 'R', 'U', 'T', diag, m, n2, alpha,a( n1 ), n, b( 0_${ik}$, &
                                    n1 ), ldb )
                          call stdlib${ii}$_sgemm( 'N', 'T',