Merge pull request #26 from xianyi/develop

rebase

common.h

@@ -416,6 +416,15 @@ please https://github.com/xianyi/OpenBLAS/issues/246
 #include "common_alpha.h"
 #endif
 
+#if (defined(ARCH_X86) || defined(ARCH_X86_64)) && defined(__CET__) && defined(__has_include)
+#if __has_include(<cet.h>)
+#include <cet.h>
+#endif
+#endif
+#ifndef _CET_ENDBR
+#define _CET_ENDBR
+#endif
+
 #ifdef ARCH_X86
 #include "common_x86.h"
 #endif

common_x86.h

@@ -340,7 +340,8 @@ REALNAME:
 	.align 16; \
 	.globl REALNAME ;\
        .type REALNAME, @function; \
-REALNAME:
+REALNAME: \
+	_CET_ENDBR
 
 #ifdef PROFILE
 #define PROFCODE call mcount

common_x86_64.h

@@ -451,7 +451,8 @@ REALNAME:
 	.align 512; \
 	.globl REALNAME ;\
        .type REALNAME, @function; \
-REALNAME:
+REALNAME: \
+	_CET_ENDBR
 
 #ifdef PROFILE
 #define PROFCODE call *mcount@GOTPCREL(%rip)

ctest/c_cblas2.c

@@ -20,7 +20,7 @@ void F77_cgemv(int *order, char *transp, int *m, int *n,
   get_transpose_type(transp, &trans);
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A  = (CBLAS_TEST_COMPLEX *)malloc( (*m)*LDA*sizeof( CBLAS_TEST_COMPLEX) );
+     A  = (CBLAS_TEST_COMPLEX *)malloc( (*m)*(size_t)LDA*sizeof( CBLAS_TEST_COMPLEX) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ){
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -50,7 +50,7 @@ void F77_cgbmv(int *order, char *transp, int *m, int *n, int *kl, int *ku,
   get_transpose_type(transp, &trans);
   if (*order == TEST_ROW_MJR) {
      LDA = *ku+*kl+2;
-     A=( CBLAS_TEST_COMPLEX* )malloc((*n+*kl)*LDA*sizeof(CBLAS_TEST_COMPLEX));
+     A=( CBLAS_TEST_COMPLEX* )malloc((*n+*kl)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX));
      for( i=0; i<*ku; i++ ){
         irow=*ku+*kl-i;
         jcol=(*ku)-i;
@@ -94,7 +94,7 @@ void F77_cgeru(int *order, int *m, int *n, CBLAS_TEST_COMPLEX *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A=(CBLAS_TEST_COMPLEX*)malloc((*m)*LDA*sizeof(CBLAS_TEST_COMPLEX));
+     A=(CBLAS_TEST_COMPLEX*)malloc((*m)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX));
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ){
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -122,7 +122,7 @@ void F77_cgerc(int *order, int *m, int *n, CBLAS_TEST_COMPLEX *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A=(CBLAS_TEST_COMPLEX* )malloc((*m)*LDA*sizeof(CBLAS_TEST_COMPLEX ) );
+     A=(CBLAS_TEST_COMPLEX* )malloc((*m)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ){
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -154,7 +154,7 @@ void F77_chemv(int *order, char *uplow, int *n, CBLAS_TEST_COMPLEX *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A = (CBLAS_TEST_COMPLEX *)malloc((*n)*LDA*sizeof(CBLAS_TEST_COMPLEX));
+     A = (CBLAS_TEST_COMPLEX *)malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX));
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ ){
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -190,7 +190,7 @@ int i,irow,j,jcol,LDA;
 		 *incx, beta, y, *incy );
      else {
         LDA = *k+2;
-        A =(CBLAS_TEST_COMPLEX*)malloc((*n+*k)*LDA*sizeof(CBLAS_TEST_COMPLEX));
+        A =(CBLAS_TEST_COMPLEX*)malloc((*n+*k)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX));
         if (uplo == CblasUpper) {
            for( i=0; i<*k; i++ ){
               irow=*k-i;
@@ -251,8 +251,8 @@ void F77_chpmv(int *order, char *uplow, int *n, CBLAS_TEST_COMPLEX *alpha,
 	         beta, y, *incy);
      else {
         LDA = *n;
-        A = (CBLAS_TEST_COMPLEX* )malloc(LDA*LDA*sizeof(CBLAS_TEST_COMPLEX ));
+        A = (CBLAS_TEST_COMPLEX* )malloc((size_t)LDA*LDA*sizeof(CBLAS_TEST_COMPLEX ));
-        AP = (CBLAS_TEST_COMPLEX* )malloc( (((LDA+1)*LDA)/2)*
+        AP = (CBLAS_TEST_COMPLEX* )malloc( ((((size_t)LDA+1)*LDA)/2)*
 	        sizeof( CBLAS_TEST_COMPLEX ));
         if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -311,7 +311,7 @@ void F77_ctbmv(int *order, char *uplow, char *transp, char *diagn,
 	x, *incx);
      else {
         LDA = *k+2;
-        A=(CBLAS_TEST_COMPLEX *)malloc((*n+*k)*LDA*sizeof(CBLAS_TEST_COMPLEX));
+        A=(CBLAS_TEST_COMPLEX *)malloc((*n+*k)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX));
         if (uplo == CblasUpper) {
            for( i=0; i<*k; i++ ){
               irow=*k-i;
@@ -375,7 +375,7 @@ void F77_ctbsv(int *order, char *uplow, char *transp, char *diagn,
 	         *incx);
      else {
         LDA = *k+2;
-        A=(CBLAS_TEST_COMPLEX*)malloc((*n+*k)*LDA*sizeof(CBLAS_TEST_COMPLEX ));
+        A=(CBLAS_TEST_COMPLEX*)malloc((*n+*k)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX ));
         if (uplo == CblasUpper) {
            for( i=0; i<*k; i++ ){
               irow=*k-i;
@@ -436,8 +436,8 @@ void F77_ctpmv(int *order, char *uplow, char *transp, char *diagn,
         cblas_ctpmv( CblasRowMajor, UNDEFINED, trans, diag, *n, ap, x, *incx );
      else {
         LDA = *n;
-        A=(CBLAS_TEST_COMPLEX*)malloc(LDA*LDA*sizeof(CBLAS_TEST_COMPLEX));
+        A=(CBLAS_TEST_COMPLEX*)malloc((size_t)LDA*LDA*sizeof(CBLAS_TEST_COMPLEX));
-        AP=(CBLAS_TEST_COMPLEX*)malloc((((LDA+1)*LDA)/2)*
+        AP=(CBLAS_TEST_COMPLEX*)malloc(((((size_t)LDA+1)*LDA)/2)*
 	 	sizeof(CBLAS_TEST_COMPLEX));
         if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -491,8 +491,8 @@ void F77_ctpsv(int *order, char *uplow, char *transp, char *diagn,
         cblas_ctpsv( CblasRowMajor, UNDEFINED, trans, diag, *n, ap, x, *incx );
      else {
         LDA = *n;
-        A=(CBLAS_TEST_COMPLEX*)malloc(LDA*LDA*sizeof(CBLAS_TEST_COMPLEX));
+        A=(CBLAS_TEST_COMPLEX*)malloc((size_t)LDA*LDA*sizeof(CBLAS_TEST_COMPLEX));
-        AP=(CBLAS_TEST_COMPLEX*)malloc((((LDA+1)*LDA)/2)*
+        AP=(CBLAS_TEST_COMPLEX*)malloc(((((size_t)LDA+1)*LDA)/2)*
 		sizeof(CBLAS_TEST_COMPLEX));
      	if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -544,7 +544,7 @@ void F77_ctrmv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA=*n+1;
-     A=(CBLAS_TEST_COMPLEX*)malloc((*n)*LDA*sizeof(CBLAS_TEST_COMPLEX));
+     A=(CBLAS_TEST_COMPLEX*)malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX));
      for( i=0; i<*n; i++ )
        for( j=0; j<*n; j++ ) {
 	  A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -573,7 +573,7 @@ void F77_ctrsv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A =(CBLAS_TEST_COMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_COMPLEX ) );
+     A =(CBLAS_TEST_COMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ ) {
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -601,8 +601,8 @@ void F77_chpr(int *order, char *uplow, int *n, float *alpha,
         cblas_chpr(CblasRowMajor, UNDEFINED, *n, *alpha, x, *incx, ap );
      else {
         LDA = *n;
-        A = (CBLAS_TEST_COMPLEX* )malloc(LDA*LDA*sizeof(CBLAS_TEST_COMPLEX ) );
+        A = (CBLAS_TEST_COMPLEX* )malloc((size_t)LDA*LDA*sizeof(CBLAS_TEST_COMPLEX ) );
-        AP = ( CBLAS_TEST_COMPLEX* )malloc( (((LDA+1)*LDA)/2)*
+        AP = ( CBLAS_TEST_COMPLEX* )malloc( ((((size_t)LDA+1)*LDA)/2)*
 		sizeof( CBLAS_TEST_COMPLEX ));
         if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -678,8 +678,8 @@ void F77_chpr2(int *order, char *uplow, int *n, CBLAS_TEST_COMPLEX *alpha,
 		     *incy, ap );
      else {
         LDA = *n;
-        A=(CBLAS_TEST_COMPLEX*)malloc( LDA*LDA*sizeof(CBLAS_TEST_COMPLEX ) );
+        A=(CBLAS_TEST_COMPLEX*)malloc( (size_t)LDA*LDA*sizeof(CBLAS_TEST_COMPLEX ) );
-        AP=(CBLAS_TEST_COMPLEX*)malloc( (((LDA+1)*LDA)/2)*
+        AP=(CBLAS_TEST_COMPLEX*)malloc( ((((size_t)LDA+1)*LDA)/2)*
 	sizeof( CBLAS_TEST_COMPLEX ));
         if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -750,7 +750,7 @@ void F77_cher(int *order, char *uplow, int *n, float *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A=(CBLAS_TEST_COMPLEX*)malloc((*n)*LDA*sizeof( CBLAS_TEST_COMPLEX ));
+     A=(CBLAS_TEST_COMPLEX*)malloc((*n)*(size_t)LDA*sizeof( CBLAS_TEST_COMPLEX ));
 
      for( i=0; i<*n; i++ )
        for( j=0; j<*n; j++ ) {
@@ -784,7 +784,7 @@ void F77_cher2(int *order, char *uplow, int *n, CBLAS_TEST_COMPLEX *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A= ( CBLAS_TEST_COMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_COMPLEX ) );
+     A= ( CBLAS_TEST_COMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_COMPLEX ) );
 
      for( i=0; i<*n; i++ )
        for( j=0; j<*n; j++ ) {

ctest/c_dblas2.c

@@ -19,7 +19,7 @@ void F77_dgemv(int *order, char *transp, int *m, int *n, double *alpha,
   get_transpose_type(transp, &trans);
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( double* )malloc( (*m)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*m)*(size_t)LDA*sizeof( double ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -43,7 +43,7 @@ void F77_dger(int *order, int *m, int *n, double *alpha, double *x, int *incx,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( double* )malloc( (*m)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*m)*(size_t)LDA*sizeof( double ) );
 
      for( i=0; i<*m; i++ ) {
        for( j=0; j<*n; j++ )
@@ -74,7 +74,7 @@ void F77_dtrmv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
      for( i=0; i<*n; i++ )
        for( j=0; j<*n; j++ )
          A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -102,7 +102,7 @@ void F77_dtrsv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -123,7 +123,7 @@ void F77_dsymv(int *order, char *uplow, int *n, double *alpha, double *a,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -146,7 +146,7 @@ void F77_dsyr(int *order, char *uplow, int *n, double *alpha, double *x,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -170,7 +170,7 @@ void F77_dsyr2(int *order, char *uplow, int *n, double *alpha, double *x,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -196,7 +196,7 @@ void F77_dgbmv(int *order, char *transp, int *m, int *n, int *kl, int *ku,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *ku+*kl+2;
-     A   = ( double* )malloc( (*n+*kl)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*n+*kl)*(size_t)LDA*sizeof( double ) );
      for( i=0; i<*ku; i++ ){
         irow=*ku+*kl-i;
         jcol=(*ku)-i;
@@ -236,7 +236,7 @@ void F77_dtbmv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *k+1;
-     A = ( double* )malloc( (*n+*k)*LDA*sizeof( double ) );
+     A = ( double* )malloc( (*n+*k)*(size_t)LDA*sizeof( double ) );
      if (uplo == CblasUpper) {
         for( i=0; i<*k; i++ ){
            irow=*k-i;
@@ -282,7 +282,7 @@ void F77_dtbsv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *k+1;
-     A = ( double* )malloc( (*n+*k)*LDA*sizeof( double ) );
+     A = ( double* )malloc( (*n+*k)*(size_t)LDA*sizeof( double ) );
      if (uplo == CblasUpper) {
         for( i=0; i<*k; i++ ){
         irow=*k-i;
@@ -325,7 +325,7 @@ void F77_dsbmv(int *order, char *uplow, int *n, int *k, double *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *k+1;
-     A   = ( double* )malloc( (*n+*k)*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (*n+*k)*(size_t)LDA*sizeof( double ) );
      if (uplo == CblasUpper) {
         for( i=0; i<*k; i++ ){
            irow=*k-i;
@@ -369,8 +369,8 @@ void F77_dspmv(int *order, char *uplow, int *n, double *alpha, double *ap,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( double* )malloc( LDA*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (size_t)LDA*LDA*sizeof( double ) );
-     AP  = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) );
+     AP  = ( double* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( double ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )
@@ -411,8 +411,8 @@ void F77_dtpmv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( double* )malloc( LDA*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (size_t)LDA*LDA*sizeof( double ) );
-     AP  = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) );
+     AP  = ( double* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( double ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )
@@ -451,8 +451,8 @@ void F77_dtpsv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( double* )malloc( LDA*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (size_t)LDA*LDA*sizeof( double ) );
-     AP  = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) );
+     AP  = ( double* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( double ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )
@@ -488,8 +488,8 @@ void F77_dspr(int *order, char *uplow, int *n, double *alpha, double *x,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( double* )malloc( LDA*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (size_t)LDA*LDA*sizeof( double ) );
-     AP  = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) );
+     AP  = ( double* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( double ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )
@@ -540,8 +540,8 @@ void F77_dspr2(int *order, char *uplow, int *n, double *alpha, double *x,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( double* )malloc( LDA*LDA*sizeof( double ) );
+     A   = ( double* )malloc( (size_t)LDA*LDA*sizeof( double ) );
-     AP  = ( double* )malloc( (((LDA+1)*LDA)/2)*sizeof( double ) );
+     AP  = ( double* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( double ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )

ctest/c_dblas3.c

@@ -26,34 +26,34 @@ void F77_dgemm(int *order, char *transpa, char *transpb, int *m, int *n,
   if (*order == TEST_ROW_MJR) {
      if (transa == CblasNoTrans) {
         LDA = *k+1;
-        A = (double *)malloc( (*m)*LDA*sizeof( double ) );
+        A = (double *)malloc( (*m)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*k; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else {
         LDA = *m+1;
-        A   = ( double* )malloc( LDA*(*k)*sizeof( double ) );
+        A   = ( double* )malloc( (size_t)LDA*(*k)*sizeof( double ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      if (transb == CblasNoTrans) {
         LDB = *n+1;
-        B   = ( double* )malloc( (*k)*LDB*sizeof( double ) );
+        B   = ( double* )malloc( (*k)*(size_t)LDB*sizeof( double ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ )
               B[i*LDB+j]=b[j*(*ldb)+i];
      }
      else {
         LDB = *k+1;
-        B   = ( double* )malloc( LDB*(*n)*sizeof( double ) );
+        B   = ( double* )malloc( (size_t)LDB*(*n)*sizeof( double ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ )
               B[i*LDB+j]=b[j*(*ldb)+i];
      }
      LDC = *n+1;
-     C   = ( double* )malloc( (*m)*LDC*sizeof( double ) );
+     C   = ( double* )malloc( (*m)*(size_t)LDC*sizeof( double ) );
      for( j=0; j<*n; j++ )
         for( i=0; i<*m; i++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -89,25 +89,25 @@ void F77_dsymm(int *order, char *rtlf, char *uplow, int *m, int *n,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A   = ( double* )malloc( (*m)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*m)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDB = *n+1;
-     B   = ( double* )malloc( (*m)*LDB*sizeof( double ) );
+     B   = ( double* )malloc( (*m)*(size_t)LDB*sizeof( double ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            B[i*LDB+j]=b[j*(*ldb)+i];
      LDC = *n+1;
-     C   = ( double* )malloc( (*m)*LDC*sizeof( double ) );
+     C   = ( double* )malloc( (*m)*(size_t)LDC*sizeof( double ) );
      for( j=0; j<*n; j++ )
         for( i=0; i<*m; i++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -143,20 +143,20 @@ void F77_dsyrk(int *order, char *uplow, char *transp, int *n, int *k,
   if (*order == TEST_ROW_MJR) {
      if (trans == CblasNoTrans) {
         LDA = *k+1;
-        A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A   = ( double* )malloc( (*k)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*k)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDC = *n+1;
-     C   = ( double* )malloc( (*n)*LDC*sizeof( double ) );
+     C   = ( double* )malloc( (*n)*(size_t)LDC*sizeof( double ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -191,8 +191,8 @@ void F77_dsyr2k(int *order, char *uplow, char *transp, int *n, int *k,
      if (trans == CblasNoTrans) {
         LDA = *k+1;
         LDB = *k+1;
-        A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
-        B   = ( double* )malloc( (*n)*LDB*sizeof( double ) );
+        B   = ( double* )malloc( (*n)*(size_t)LDB*sizeof( double ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ ) {
               A[i*LDA+j]=a[j*(*lda)+i];
@@ -202,8 +202,8 @@ void F77_dsyr2k(int *order, char *uplow, char *transp, int *n, int *k,
      else {
         LDA = *n+1;
         LDB = *n+1;
-        A   = ( double* )malloc( LDA*(*k)*sizeof( double ) );
+        A   = ( double* )malloc( (size_t)LDA*(*k)*sizeof( double ) );
-        B   = ( double* )malloc( LDB*(*k)*sizeof( double ) );
+        B   = ( double* )malloc( (size_t)LDB*(*k)*sizeof( double ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ ){
               A[i*LDA+j]=a[j*(*lda)+i];
@@ -211,7 +211,7 @@ void F77_dsyr2k(int *order, char *uplow, char *transp, int *n, int *k,
            }
      }
      LDC = *n+1;
-     C   = ( double* )malloc( (*n)*LDC*sizeof( double ) );
+     C   = ( double* )malloc( (*n)*(size_t)LDC*sizeof( double ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -249,20 +249,20 @@ void F77_dtrmm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A   = ( double* )malloc( (*m)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*m)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDB = *n+1;
-     B   = ( double* )malloc( (*m)*LDB*sizeof( double ) );
+     B   = ( double* )malloc( (*m)*(size_t)LDB*sizeof( double ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            B[i*LDB+j]=b[j*(*ldb)+i];
@@ -300,20 +300,20 @@ void F77_dtrsm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A   = ( double* )malloc( (*m)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*m)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A   = ( double* )malloc( (*n)*LDA*sizeof( double ) );
+        A   = ( double* )malloc( (*n)*(size_t)LDA*sizeof( double ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDB = *n+1;
-     B   = ( double* )malloc( (*m)*LDB*sizeof( double ) );
+     B   = ( double* )malloc( (*m)*(size_t)LDB*sizeof( double ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            B[i*LDB+j]=b[j*(*ldb)+i];

ctest/c_sblas2.c

@@ -19,7 +19,7 @@ void F77_sgemv(int *order, char *transp, int *m, int *n, float *alpha,
   get_transpose_type(transp, &trans);
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*m)*(size_t)LDA*sizeof( float ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -43,7 +43,7 @@ void F77_sger(int *order, int *m, int *n, float *alpha, float *x, int *incx,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*m)*(size_t)LDA*sizeof( float ) );
 
      for( i=0; i<*m; i++ ) {
        for( j=0; j<*n; j++ )
@@ -74,7 +74,7 @@ void F77_strmv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
      for( i=0; i<*n; i++ )
        for( j=0; j<*n; j++ )
          A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -102,7 +102,7 @@ void F77_strsv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -123,7 +123,7 @@ void F77_ssymv(int *order, char *uplow, int *n, float *alpha, float *a,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -146,7 +146,7 @@ void F77_ssyr(int *order, char *uplow, int *n, float *alpha, float *x,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -170,7 +170,7 @@ void F77_ssyr2(int *order, char *uplow, int *n, float *alpha, float *x,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            A[ LDA*i+j ]=a[ (*lda)*j+i ];
@@ -196,7 +196,7 @@ void F77_sgbmv(int *order, char *transp, int *m, int *n, int *kl, int *ku,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *ku+*kl+2;
-     A   = ( float* )malloc( (*n+*kl)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*n+*kl)*(size_t)LDA*sizeof( float ) );
      for( i=0; i<*ku; i++ ){
         irow=*ku+*kl-i;
         jcol=(*ku)-i;
@@ -236,7 +236,7 @@ void F77_stbmv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *k+1;
-     A = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
+     A = ( float* )malloc( (*n+*k)*(size_t)LDA*sizeof( float ) );
      if (uplo == CblasUpper) {
         for( i=0; i<*k; i++ ){
            irow=*k-i;
@@ -282,7 +282,7 @@ void F77_stbsv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *k+1;
-     A = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
+     A = ( float* )malloc( (*n+*k)*(size_t)LDA*sizeof( float ) );
      if (uplo == CblasUpper) {
         for( i=0; i<*k; i++ ){
         irow=*k-i;
@@ -325,7 +325,7 @@ void F77_ssbmv(int *order, char *uplow, int *n, int *k, float *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *k+1;
-     A   = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (*n+*k)*(size_t)LDA*sizeof( float ) );
      if (uplo == CblasUpper) {
         for( i=0; i<*k; i++ ){
            irow=*k-i;
@@ -369,8 +369,8 @@ void F77_sspmv(int *order, char *uplow, int *n, float *alpha, float *ap,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (size_t)LDA*LDA*sizeof( float ) );
-     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
+     AP  = ( float* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( float ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )
@@ -410,8 +410,8 @@ void F77_stpmv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (size_t)LDA*LDA*sizeof( float ) );
-     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
+     AP  = ( float* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( float ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )
@@ -449,8 +449,8 @@ void F77_stpsv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (size_t)LDA*LDA*sizeof( float ) );
-     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
+     AP  = ( float* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( float ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )
@@ -485,8 +485,8 @@ void F77_sspr(int *order, char *uplow, int *n, float *alpha, float *x,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (size_t)LDA*LDA*sizeof( float ) );
-     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
+     AP  = ( float* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( float ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )
@@ -536,8 +536,8 @@ void F77_sspr2(int *order, char *uplow, int *n, float *alpha, float *x,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n;
-     A   = ( float* )malloc( LDA*LDA*sizeof( float ) );
+     A   = ( float* )malloc( (size_t)LDA*LDA*sizeof( float ) );
-     AP  = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
+     AP  = ( float* )malloc( ((((size_t)LDA+1)*LDA)/2)*sizeof( float ) );
      if (uplo == CblasUpper) {
         for( j=0, k=0; j<*n; j++ )
            for( i=0; i<j+1; i++, k++ )

ctest/c_sblas3.c

@@ -23,34 +23,34 @@ void F77_sgemm(int *order, char *transpa, char *transpb, int *m, int *n,
   if (*order == TEST_ROW_MJR) {
      if (transa == CblasNoTrans) {
         LDA = *k+1;
-        A = (float *)malloc( (*m)*LDA*sizeof( float ) );
+        A = (float *)malloc( (*m)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*k; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else {
         LDA = *m+1;
-        A   = ( float* )malloc( LDA*(*k)*sizeof( float ) );
+        A   = ( float* )malloc( (size_t)LDA*(*k)*sizeof( float ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      if (transb == CblasNoTrans) {
         LDB = *n+1;
-        B   = ( float* )malloc( (*k)*LDB*sizeof( float ) );
+        B   = ( float* )malloc( (*k)*(size_t)LDB*sizeof( float ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ )
               B[i*LDB+j]=b[j*(*ldb)+i];
      }
      else {
         LDB = *k+1;
-        B   = ( float* )malloc( LDB*(*n)*sizeof( float ) );
+        B   = ( float* )malloc( (size_t)LDB*(*n)*sizeof( float ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ )
               B[i*LDB+j]=b[j*(*ldb)+i];
      }
      LDC = *n+1;
-     C   = ( float* )malloc( (*m)*LDC*sizeof( float ) );
+     C   = ( float* )malloc( (*m)*(size_t)LDC*sizeof( float ) );
      for( j=0; j<*n; j++ )
         for( i=0; i<*m; i++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -85,25 +85,25 @@ void F77_ssymm(int *order, char *rtlf, char *uplow, int *m, int *n,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*m)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDB = *n+1;
-     B   = ( float* )malloc( (*m)*LDB*sizeof( float ) );
+     B   = ( float* )malloc( (*m)*(size_t)LDB*sizeof( float ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            B[i*LDB+j]=b[j*(*ldb)+i];
      LDC = *n+1;
-     C   = ( float* )malloc( (*m)*LDC*sizeof( float ) );
+     C   = ( float* )malloc( (*m)*(size_t)LDC*sizeof( float ) );
      for( j=0; j<*n; j++ )
         for( i=0; i<*m; i++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -139,20 +139,20 @@ void F77_ssyrk(int *order, char *uplow, char *transp, int *n, int *k,
   if (*order == TEST_ROW_MJR) {
      if (trans == CblasNoTrans) {
         LDA = *k+1;
-        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A   = ( float* )malloc( (*k)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*k)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDC = *n+1;
-     C   = ( float* )malloc( (*n)*LDC*sizeof( float ) );
+     C   = ( float* )malloc( (*n)*(size_t)LDC*sizeof( float ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -187,8 +187,8 @@ void F77_ssyr2k(int *order, char *uplow, char *transp, int *n, int *k,
      if (trans == CblasNoTrans) {
         LDA = *k+1;
         LDB = *k+1;
-        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
-        B   = ( float* )malloc( (*n)*LDB*sizeof( float ) );
+        B   = ( float* )malloc( (*n)*(size_t)LDB*sizeof( float ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ ) {
               A[i*LDA+j]=a[j*(*lda)+i];
@@ -198,8 +198,8 @@ void F77_ssyr2k(int *order, char *uplow, char *transp, int *n, int *k,
      else {
         LDA = *n+1;
         LDB = *n+1;
-        A   = ( float* )malloc( LDA*(*k)*sizeof( float ) );
+        A   = ( float* )malloc( (size_t)LDA*(*k)*sizeof( float ) );
-        B   = ( float* )malloc( LDB*(*k)*sizeof( float ) );
+        B   = ( float* )malloc( (size_t)LDB*(*k)*sizeof( float ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ ){
               A[i*LDA+j]=a[j*(*lda)+i];
@@ -207,7 +207,7 @@ void F77_ssyr2k(int *order, char *uplow, char *transp, int *n, int *k,
            }
      }
      LDC = *n+1;
-     C   = ( float* )malloc( (*n)*LDC*sizeof( float ) );
+     C   = ( float* )malloc( (*n)*(size_t)LDC*sizeof( float ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -245,20 +245,20 @@ void F77_strmm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*m)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDB = *n+1;
-     B   = ( float* )malloc( (*m)*LDB*sizeof( float ) );
+     B   = ( float* )malloc( (*m)*(size_t)LDB*sizeof( float ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            B[i*LDB+j]=b[j*(*ldb)+i];
@@ -296,20 +296,20 @@ void F77_strsm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A   = ( float* )malloc( (*m)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*m)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A   = ( float* )malloc( (*n)*LDA*sizeof( float ) );
+        A   = ( float* )malloc( (*n)*(size_t)LDA*sizeof( float ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDB = *n+1;
-     B   = ( float* )malloc( (*m)*LDB*sizeof( float ) );
+     B   = ( float* )malloc( (*m)*(size_t)LDB*sizeof( float ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            B[i*LDB+j]=b[j*(*ldb)+i];

ctest/c_zblas2.c

@@ -20,7 +20,7 @@ void F77_zgemv(int *order, char *transp, int *m, int *n,
   get_transpose_type(transp, &trans);
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A  = (CBLAS_TEST_ZOMPLEX *)malloc( (*m)*LDA*sizeof( CBLAS_TEST_ZOMPLEX) );
+     A  = (CBLAS_TEST_ZOMPLEX *)malloc( (*m)*(size_t)LDA*sizeof( CBLAS_TEST_ZOMPLEX) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ){
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -50,7 +50,7 @@ void F77_zgbmv(int *order, char *transp, int *m, int *n, int *kl, int *ku,
   get_transpose_type(transp, &trans);
   if (*order == TEST_ROW_MJR) {
      LDA = *ku+*kl+2;
-     A=( CBLAS_TEST_ZOMPLEX* )malloc((*n+*kl)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+     A=( CBLAS_TEST_ZOMPLEX* )malloc((*n+*kl)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
      for( i=0; i<*ku; i++ ){
         irow=*ku+*kl-i;
         jcol=(*ku)-i;
@@ -94,7 +94,7 @@ void F77_zgeru(int *order, int *m, int *n, CBLAS_TEST_ZOMPLEX *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A=(CBLAS_TEST_ZOMPLEX*)malloc((*m)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+     A=(CBLAS_TEST_ZOMPLEX*)malloc((*m)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ){
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -122,7 +122,7 @@ void F77_zgerc(int *order, int *m, int *n, CBLAS_TEST_ZOMPLEX *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+     A=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ){
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -154,7 +154,7 @@ void F77_zhemv(int *order, char *uplow, int *n, CBLAS_TEST_ZOMPLEX *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A = (CBLAS_TEST_ZOMPLEX *)malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+     A = (CBLAS_TEST_ZOMPLEX *)malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ ){
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -190,7 +190,7 @@ int i,irow,j,jcol,LDA;
 		 *incx, beta, y, *incy );
      else {
         LDA = *k+2;
-        A =(CBLAS_TEST_ZOMPLEX*)malloc((*n+*k)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A =(CBLAS_TEST_ZOMPLEX*)malloc((*n+*k)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         if (uplo == CblasUpper) {
            for( i=0; i<*k; i++ ){
               irow=*k-i;
@@ -251,8 +251,8 @@ void F77_zhpmv(int *order, char *uplow, int *n, CBLAS_TEST_ZOMPLEX *alpha,
 	         beta, y, *incy);
      else {
         LDA = *n;
-        A = (CBLAS_TEST_ZOMPLEX* )malloc(LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX ));
+        A = (CBLAS_TEST_ZOMPLEX* )malloc((size_t)LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX ));
-        AP = (CBLAS_TEST_ZOMPLEX* )malloc( (((LDA+1)*LDA)/2)*
+        AP = (CBLAS_TEST_ZOMPLEX* )malloc( ((((size_t)LDA+1)*LDA)/2)*
 	        sizeof( CBLAS_TEST_ZOMPLEX ));
         if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -311,7 +311,7 @@ void F77_ztbmv(int *order, char *uplow, char *transp, char *diagn,
 	x, *incx);
      else {
         LDA = *k+2;
-        A=(CBLAS_TEST_ZOMPLEX *)malloc((*n+*k)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX *)malloc((*n+*k)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         if (uplo == CblasUpper) {
            for( i=0; i<*k; i++ ){
               irow=*k-i;
@@ -375,7 +375,7 @@ void F77_ztbsv(int *order, char *uplow, char *transp, char *diagn,
 	         *incx);
      else {
         LDA = *k+2;
-        A=(CBLAS_TEST_ZOMPLEX*)malloc((*n+*k)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ));
+        A=(CBLAS_TEST_ZOMPLEX*)malloc((*n+*k)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ));
         if (uplo == CblasUpper) {
            for( i=0; i<*k; i++ ){
               irow=*k-i;
@@ -436,8 +436,8 @@ void F77_ztpmv(int *order, char *uplow, char *transp, char *diagn,
         cblas_ztpmv( CblasRowMajor, UNDEFINED, trans, diag, *n, ap, x, *incx );
      else {
         LDA = *n;
-        A=(CBLAS_TEST_ZOMPLEX*)malloc(LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX*)malloc((size_t)LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
-        AP=(CBLAS_TEST_ZOMPLEX*)malloc((((LDA+1)*LDA)/2)*
+        AP=(CBLAS_TEST_ZOMPLEX*)malloc(((((size_t)LDA+1)*LDA)/2)*
 	 	sizeof(CBLAS_TEST_ZOMPLEX));
         if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -491,8 +491,8 @@ void F77_ztpsv(int *order, char *uplow, char *transp, char *diagn,
         cblas_ztpsv( CblasRowMajor, UNDEFINED, trans, diag, *n, ap, x, *incx );
      else {
         LDA = *n;
-        A=(CBLAS_TEST_ZOMPLEX*)malloc(LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX*)malloc((size_t)LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
-        AP=(CBLAS_TEST_ZOMPLEX*)malloc((((LDA+1)*LDA)/2)*
+        AP=(CBLAS_TEST_ZOMPLEX*)malloc(((((size_t)LDA+1)*LDA)/2)*
 		sizeof(CBLAS_TEST_ZOMPLEX));
      	if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -544,7 +544,7 @@ void F77_ztrmv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA=*n+1;
-     A=(CBLAS_TEST_ZOMPLEX*)malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+     A=(CBLAS_TEST_ZOMPLEX*)malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
      for( i=0; i<*n; i++ )
        for( j=0; j<*n; j++ ) {
 	  A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -573,7 +573,7 @@ void F77_ztrsv(int *order, char *uplow, char *transp, char *diagn,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A =(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+     A =(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ ) {
            A[ LDA*i+j ].real=a[ (*lda)*j+i ].real;
@@ -601,8 +601,8 @@ void F77_zhpr(int *order, char *uplow, int *n, double *alpha,
         cblas_zhpr(CblasRowMajor, UNDEFINED, *n, *alpha, x, *incx, ap );
      else {
         LDA = *n;
-        A = (CBLAS_TEST_ZOMPLEX* )malloc(LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A = (CBLAS_TEST_ZOMPLEX* )malloc((size_t)LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
-        AP = ( CBLAS_TEST_ZOMPLEX* )malloc( (((LDA+1)*LDA)/2)*
+        AP = ( CBLAS_TEST_ZOMPLEX* )malloc( ((((size_t)LDA+1)*LDA)/2)*
 		sizeof( CBLAS_TEST_ZOMPLEX ));
         if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -678,8 +678,8 @@ void F77_zhpr2(int *order, char *uplow, int *n, CBLAS_TEST_ZOMPLEX *alpha,
 		     *incy, ap );
      else {
         LDA = *n;
-        A=(CBLAS_TEST_ZOMPLEX*)malloc( LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A=(CBLAS_TEST_ZOMPLEX*)malloc( (size_t)LDA*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
-        AP=(CBLAS_TEST_ZOMPLEX*)malloc( (((LDA+1)*LDA)/2)*
+        AP=(CBLAS_TEST_ZOMPLEX*)malloc( ((((size_t)LDA+1)*LDA)/2)*
 	sizeof( CBLAS_TEST_ZOMPLEX ));
         if (uplo == CblasUpper) {
            for( j=0, k=0; j<*n; j++ )
@@ -750,7 +750,7 @@ void F77_zher(int *order, char *uplow, int *n, double *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A=(CBLAS_TEST_ZOMPLEX*)malloc((*n)*LDA*sizeof( CBLAS_TEST_ZOMPLEX ));
+     A=(CBLAS_TEST_ZOMPLEX*)malloc((*n)*(size_t)LDA*sizeof( CBLAS_TEST_ZOMPLEX ));
 
      for( i=0; i<*n; i++ )
        for( j=0; j<*n; j++ ) {
@@ -784,7 +784,7 @@ void F77_zher2(int *order, char *uplow, int *n, CBLAS_TEST_ZOMPLEX *alpha,
 
   if (*order == TEST_ROW_MJR) {
      LDA = *n+1;
-     A= ( CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+     A= ( CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
 
      for( i=0; i<*n; i++ )
        for( j=0; j<*n; j++ ) {

ctest/c_zblas3.c

@@ -26,7 +26,7 @@ void F77_zgemm(int *order, char *transpa, char *transpb, int *m, int *n,
   if (*order == TEST_ROW_MJR) {
      if (transa == CblasNoTrans) {
         LDA = *k+1;
-        A=(CBLAS_TEST_ZOMPLEX*)malloc((*m)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX*)malloc((*m)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*m; i++ )
            for( j=0; j<*k; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -35,7 +35,7 @@ void F77_zgemm(int *order, char *transpa, char *transpb, int *m, int *n,
      }
      else {
         LDA = *m+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc(LDA*(*k)*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((size_t)LDA*(*k)*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*k; i++ )
            for( j=0; j<*m; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -45,7 +45,7 @@ void F77_zgemm(int *order, char *transpa, char *transpb, int *m, int *n,
 
      if (transb == CblasNoTrans) {
         LDB = *n+1;
-        B=(CBLAS_TEST_ZOMPLEX* )malloc((*k)*LDB*sizeof(CBLAS_TEST_ZOMPLEX) );
+        B=(CBLAS_TEST_ZOMPLEX* )malloc((*k)*(size_t)LDB*sizeof(CBLAS_TEST_ZOMPLEX) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ ) {
               B[i*LDB+j].real=b[j*(*ldb)+i].real;
@@ -54,7 +54,7 @@ void F77_zgemm(int *order, char *transpa, char *transpb, int *m, int *n,
      }
      else {
         LDB = *k+1;
-        B=(CBLAS_TEST_ZOMPLEX* )malloc(LDB*(*n)*sizeof(CBLAS_TEST_ZOMPLEX));
+        B=(CBLAS_TEST_ZOMPLEX* )malloc((size_t)LDB*(*n)*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ ) {
               B[i*LDB+j].real=b[j*(*ldb)+i].real;
@@ -63,7 +63,7 @@ void F77_zgemm(int *order, char *transpa, char *transpb, int *m, int *n,
      }
 
      LDC = *n+1;
-     C=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDC*sizeof(CBLAS_TEST_ZOMPLEX));
+     C=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDC*sizeof(CBLAS_TEST_ZOMPLEX));
      for( j=0; j<*n; j++ )
         for( i=0; i<*m; i++ ) {
            C[i*LDC+j].real=c[j*(*ldc)+i].real;
@@ -103,7 +103,7 @@ void F77_zhemm(int *order, char *rtlf, char *uplow, int *m, int *n,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A= (CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A= (CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -112,7 +112,7 @@ void F77_zhemm(int *order, char *rtlf, char *uplow, int *m, int *n,
      }
      else{
         LDA = *n+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -120,14 +120,14 @@ void F77_zhemm(int *order, char *rtlf, char *uplow, int *m, int *n,
            }
      }
      LDB = *n+1;
-     B=(CBLAS_TEST_ZOMPLEX* )malloc( (*m)*LDB*sizeof(CBLAS_TEST_ZOMPLEX ) );
+     B=(CBLAS_TEST_ZOMPLEX* )malloc( (*m)*(size_t)LDB*sizeof(CBLAS_TEST_ZOMPLEX ) );
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ) {
            B[i*LDB+j].real=b[j*(*ldb)+i].real;
            B[i*LDB+j].imag=b[j*(*ldb)+i].imag;
         }
      LDC = *n+1;
-     C=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDC*sizeof(CBLAS_TEST_ZOMPLEX ) );
+     C=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDC*sizeof(CBLAS_TEST_ZOMPLEX ) );
      for( j=0; j<*n; j++ )
         for( i=0; i<*m; i++ ) {
            C[i*LDC+j].real=c[j*(*ldc)+i].real;
@@ -167,25 +167,25 @@ void F77_zsymm(int *order, char *rtlf, char *uplow, int *m, int *n,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      else{
         LDA = *n+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ )
               A[i*LDA+j]=a[j*(*lda)+i];
      }
      LDB = *n+1;
-     B=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDB*sizeof(CBLAS_TEST_ZOMPLEX ));
+     B=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDB*sizeof(CBLAS_TEST_ZOMPLEX ));
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ )
            B[i*LDB+j]=b[j*(*ldb)+i];
      LDC = *n+1;
-     C=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDC*sizeof(CBLAS_TEST_ZOMPLEX));
+     C=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDC*sizeof(CBLAS_TEST_ZOMPLEX));
      for( j=0; j<*n; j++ )
         for( i=0; i<*m; i++ )
            C[i*LDC+j]=c[j*(*ldc)+i];
@@ -221,7 +221,7 @@ void F77_zherk(int *order, char *uplow, char *transp, int *n, int *k,
   if (*order == TEST_ROW_MJR) {
      if (trans == CblasNoTrans) {
         LDA = *k+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -230,7 +230,7 @@ void F77_zherk(int *order, char *uplow, char *transp, int *n, int *k,
      }
      else{
         LDA = *n+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*k)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*k)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -238,7 +238,7 @@ void F77_zherk(int *order, char *uplow, char *transp, int *n, int *k,
            }
      }
      LDC = *n+1;
-     C=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDC*sizeof(CBLAS_TEST_ZOMPLEX ) );
+     C=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDC*sizeof(CBLAS_TEST_ZOMPLEX ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ ) {
            C[i*LDC+j].real=c[j*(*ldc)+i].real;
@@ -277,7 +277,7 @@ void F77_zsyrk(int *order, char *uplow, char *transp, int *n, int *k,
   if (*order == TEST_ROW_MJR) {
      if (trans == CblasNoTrans) {
         LDA = *k+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -286,7 +286,7 @@ void F77_zsyrk(int *order, char *uplow, char *transp, int *n, int *k,
      }
      else{
         LDA = *n+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*k)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*k)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -294,7 +294,7 @@ void F77_zsyrk(int *order, char *uplow, char *transp, int *n, int *k,
            }
      }
      LDC = *n+1;
-     C=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDC*sizeof(CBLAS_TEST_ZOMPLEX ) );
+     C=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDC*sizeof(CBLAS_TEST_ZOMPLEX ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ ) {
            C[i*LDC+j].real=c[j*(*ldc)+i].real;
@@ -333,8 +333,8 @@ void F77_zher2k(int *order, char *uplow, char *transp, int *n, int *k,
      if (trans == CblasNoTrans) {
         LDA = *k+1;
         LDB = *k+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ));
-        B=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDB*sizeof(CBLAS_TEST_ZOMPLEX ));
+        B=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDB*sizeof(CBLAS_TEST_ZOMPLEX ));
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -346,8 +346,8 @@ void F77_zher2k(int *order, char *uplow, char *transp, int *n, int *k,
      else {
         LDA = *n+1;
         LDB = *n+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc( LDA*(*k)*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A=(CBLAS_TEST_ZOMPLEX* )malloc( (size_t)LDA*(*k)*sizeof(CBLAS_TEST_ZOMPLEX ) );
-        B=(CBLAS_TEST_ZOMPLEX* )malloc( LDB*(*k)*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        B=(CBLAS_TEST_ZOMPLEX* )malloc( (size_t)LDB*(*k)*sizeof(CBLAS_TEST_ZOMPLEX ) );
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ ){
 	      A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -357,7 +357,7 @@ void F77_zher2k(int *order, char *uplow, char *transp, int *n, int *k,
            }
      }
      LDC = *n+1;
-     C=(CBLAS_TEST_ZOMPLEX* )malloc( (*n)*LDC*sizeof(CBLAS_TEST_ZOMPLEX ) );
+     C=(CBLAS_TEST_ZOMPLEX* )malloc( (*n)*(size_t)LDC*sizeof(CBLAS_TEST_ZOMPLEX ) );
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ ) {
            C[i*LDC+j].real=c[j*(*ldc)+i].real;
@@ -397,8 +397,8 @@ void F77_zsyr2k(int *order, char *uplow, char *transp, int *n, int *k,
      if (trans == CblasNoTrans) {
         LDA = *k+1;
         LDB = *k+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
-        B=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDB*sizeof(CBLAS_TEST_ZOMPLEX));
+        B=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDB*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*n; i++ )
            for( j=0; j<*k; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -410,8 +410,8 @@ void F77_zsyr2k(int *order, char *uplow, char *transp, int *n, int *k,
      else {
         LDA = *n+1;
         LDB = *n+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc(LDA*(*k)*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((size_t)LDA*(*k)*sizeof(CBLAS_TEST_ZOMPLEX));
-        B=(CBLAS_TEST_ZOMPLEX* )malloc(LDB*(*k)*sizeof(CBLAS_TEST_ZOMPLEX));
+        B=(CBLAS_TEST_ZOMPLEX* )malloc((size_t)LDB*(*k)*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*k; i++ )
            for( j=0; j<*n; j++ ){
 	      A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -421,7 +421,7 @@ void F77_zsyr2k(int *order, char *uplow, char *transp, int *n, int *k,
            }
      }
      LDC = *n+1;
-     C=(CBLAS_TEST_ZOMPLEX* )malloc( (*n)*LDC*sizeof(CBLAS_TEST_ZOMPLEX));
+     C=(CBLAS_TEST_ZOMPLEX* )malloc( (*n)*(size_t)LDC*sizeof(CBLAS_TEST_ZOMPLEX));
      for( i=0; i<*n; i++ )
         for( j=0; j<*n; j++ ) {
            C[i*LDC+j].real=c[j*(*ldc)+i].real;
@@ -463,7 +463,7 @@ void F77_ztrmm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -472,7 +472,7 @@ void F77_ztrmm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
      }
      else{
         LDA = *n+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -480,7 +480,7 @@ void F77_ztrmm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
            }
      }
      LDB = *n+1;
-     B=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDB*sizeof(CBLAS_TEST_ZOMPLEX));
+     B=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDB*sizeof(CBLAS_TEST_ZOMPLEX));
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ) {
            B[i*LDB+j].real=b[j*(*ldb)+i].real;
@@ -522,7 +522,7 @@ void F77_ztrsm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
   if (*order == TEST_ROW_MJR) {
      if (side == CblasLeft) {
         LDA = *m+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc( (*m)*LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
+        A=(CBLAS_TEST_ZOMPLEX* )malloc( (*m)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX ) );
         for( i=0; i<*m; i++ )
            for( j=0; j<*m; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -531,7 +531,7 @@ void F77_ztrsm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
      }
      else{
         LDA = *n+1;
-        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*LDA*sizeof(CBLAS_TEST_ZOMPLEX));
+        A=(CBLAS_TEST_ZOMPLEX* )malloc((*n)*(size_t)LDA*sizeof(CBLAS_TEST_ZOMPLEX));
         for( i=0; i<*n; i++ )
            for( j=0; j<*n; j++ ) {
               A[i*LDA+j].real=a[j*(*lda)+i].real;
@@ -539,7 +539,7 @@ void F77_ztrsm(int *order, char *rtlf, char *uplow, char *transp, char *diagn,
            }
      }
      LDB = *n+1;
-     B=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*LDB*sizeof(CBLAS_TEST_ZOMPLEX));
+     B=(CBLAS_TEST_ZOMPLEX* )malloc((*m)*(size_t)LDB*sizeof(CBLAS_TEST_ZOMPLEX));
      for( i=0; i<*m; i++ )
         for( j=0; j<*n; j++ ) {
            B[i*LDB+j].real=b[j*(*ldb)+i].real;

interface/imatcopy.c

@@ -150,9 +150,9 @@ void CNAME( enum CBLAS_ORDER CORDER, enum CBLAS_TRANSPOSE CTRANS, blasint crows,
 #endif
 
 	if ( *lda >  *ldb )
-		msize = (*lda) * (*ldb)  * sizeof(FLOAT);
+		msize = (size_t)(*lda) * (*ldb)  * sizeof(FLOAT);
 	else
-		msize = (*ldb) * (*ldb)  * sizeof(FLOAT);
+		msize = (size_t)(*ldb) * (*ldb)  * sizeof(FLOAT);
 
 	b = malloc(msize);
 	if ( b == NULL )

interface/zimatcopy.c

@@ -172,9 +172,9 @@ void CNAME( enum CBLAS_ORDER CORDER, enum CBLAS_TRANSPOSE CTRANS, blasint crows,
 #endif
 
 	if ( *lda >  *ldb )
-                msize = (*lda) * (*ldb)  * sizeof(FLOAT) * 2;
+                msize = (size_t)(*lda) * (*ldb)  * sizeof(FLOAT) * 2;
         else
-                msize = (*ldb) * (*ldb)  * sizeof(FLOAT) * 2;
+                msize = (size_t)(*ldb) * (*ldb)  * sizeof(FLOAT) * 2;
 
         b = malloc(msize);
         if ( b == NULL )

lapack-netlib/LAPACKE/include/lapack.h

@@ -4768,7 +4768,7 @@ void LAPACK_chegst(
     lapack_int const* itype, char const* uplo,
     lapack_int const* n,
     lapack_complex_float* A, lapack_int const* lda,
-    lapack_complex_float* B, lapack_int const* ldb,
+    const lapack_complex_float* B, lapack_int const* ldb,
     lapack_int* info );
 
 #define LAPACK_zhegst LAPACK_GLOBAL(zhegst,ZHEGST)
@@ -4776,7 +4776,7 @@ void LAPACK_zhegst(
     lapack_int const* itype, char const* uplo,
     lapack_int const* n,
     lapack_complex_double* A, lapack_int const* lda,
-    lapack_complex_double* B, lapack_int const* ldb,
+    const lapack_complex_double* B, lapack_int const* ldb,
     lapack_int* info );
 
 #define LAPACK_chegv LAPACK_GLOBAL(chegv,CHEGV)
@@ -11556,7 +11556,7 @@ void LAPACK_zsytrs(
 void LAPACK_csytrs2(
     char const* uplo,
     lapack_int const* n, lapack_int const* nrhs,
-    lapack_complex_float* A, lapack_int const* lda, lapack_int const* ipiv,
+    const lapack_complex_float* A, lapack_int const* lda, lapack_int const* ipiv,
     lapack_complex_float* B, lapack_int const* ldb,
     lapack_complex_float* work,
     lapack_int* info );
@@ -11565,7 +11565,7 @@ void LAPACK_csytrs2(
 void LAPACK_dsytrs2(
     char const* uplo,
     lapack_int const* n, lapack_int const* nrhs,
-    double* A, lapack_int const* lda, lapack_int const* ipiv,
+    const double* A, lapack_int const* lda, lapack_int const* ipiv,
     double* B, lapack_int const* ldb,
     double* work,
     lapack_int* info );
@@ -11574,7 +11574,7 @@ void LAPACK_dsytrs2(
 void LAPACK_ssytrs2(
     char const* uplo,
     lapack_int const* n, lapack_int const* nrhs,
-    float* A, lapack_int const* lda, lapack_int const* ipiv,
+    const float* A, lapack_int const* lda, lapack_int const* ipiv,
     float* B, lapack_int const* ldb,
     float* work,
     lapack_int* info );
@@ -11583,7 +11583,7 @@ void LAPACK_ssytrs2(
 void LAPACK_zsytrs2(
     char const* uplo,
     lapack_int const* n, lapack_int const* nrhs,
-    lapack_complex_double* A, lapack_int const* lda, lapack_int const* ipiv,
+    const lapack_complex_double* A, lapack_int const* lda, lapack_int const* ipiv,
     lapack_complex_double* B, lapack_int const* ldb,
     lapack_complex_double* work,
     lapack_int* info );

lapack-netlib/LAPACKE/include/lapacke.h

@@ -1867,11 +1867,11 @@ lapack_int LAPACKE_zheevx( int matrix_layout, char jobz, char range, char uplo,
 
 lapack_int LAPACKE_chegst( int matrix_layout, lapack_int itype, char uplo,
                            lapack_int n, lapack_complex_float* a,
-                           lapack_int lda, lapack_complex_float* b,
+                           lapack_int lda, const lapack_complex_float* b,
                            lapack_int ldb );
 lapack_int LAPACKE_zhegst( int matrix_layout, lapack_int itype, char uplo,
                            lapack_int n, lapack_complex_double* a,
-                           lapack_int lda, lapack_complex_double* b,
+                           lapack_int lda, const lapack_complex_double* b,
                            lapack_int ldb );
 
 lapack_int LAPACKE_chegv( int matrix_layout, lapack_int itype, char jobz,
@@ -6932,11 +6932,11 @@ lapack_int LAPACKE_zheevx_work( int matrix_layout, char jobz, char range,
 
 lapack_int LAPACKE_chegst_work( int matrix_layout, lapack_int itype, char uplo,
                                 lapack_int n, lapack_complex_float* a,
-                                lapack_int lda, lapack_complex_float* b,
+                                lapack_int lda, const lapack_complex_float* b,
                                 lapack_int ldb );
 lapack_int LAPACKE_zhegst_work( int matrix_layout, lapack_int itype, char uplo,
                                 lapack_int n, lapack_complex_double* a,
-                                lapack_int lda, lapack_complex_double* b,
+                                lapack_int lda, const lapack_complex_double* b,
                                 lapack_int ldb );
 
 lapack_int LAPACKE_chegv_work( int matrix_layout, lapack_int itype, char jobz,
@@ -10553,11 +10553,11 @@ lapack_int LAPACKE_csytri2x_work( int matrix_layout, char uplo, lapack_int n,
                                   const lapack_int* ipiv,
                                   lapack_complex_float* work, lapack_int nb );
 lapack_int LAPACKE_csytrs2( int matrix_layout, char uplo, lapack_int n,
-                            lapack_int nrhs, lapack_complex_float* a,
+                            lapack_int nrhs, const lapack_complex_float* a,
                             lapack_int lda, const lapack_int* ipiv,
                             lapack_complex_float* b, lapack_int ldb );
 lapack_int LAPACKE_csytrs2_work( int matrix_layout, char uplo, lapack_int n,
-                                 lapack_int nrhs, lapack_complex_float* a,
+                                 lapack_int nrhs, const lapack_complex_float* a,
                                  lapack_int lda, const lapack_int* ipiv,
                                  lapack_complex_float* b, lapack_int ldb,
                                  lapack_complex_float* work );
@@ -10718,10 +10718,10 @@ lapack_int LAPACKE_dsytri2x_work( int matrix_layout, char uplo, lapack_int n,
                                   const lapack_int* ipiv, double* work,
                                   lapack_int nb );
 lapack_int LAPACKE_dsytrs2( int matrix_layout, char uplo, lapack_int n,
-                            lapack_int nrhs, double* a, lapack_int lda,
+                            lapack_int nrhs, const double* a, lapack_int lda,
                             const lapack_int* ipiv, double* b, lapack_int ldb );
 lapack_int LAPACKE_dsytrs2_work( int matrix_layout, char uplo, lapack_int n,
-                                 lapack_int nrhs, double* a,
+                                 lapack_int nrhs, const double* a,
                                  lapack_int lda, const lapack_int* ipiv,
                                  double* b, lapack_int ldb, double* work );
 lapack_int LAPACKE_sbbcsd( int matrix_layout, char jobu1, char jobu2,
@@ -10813,10 +10813,10 @@ lapack_int LAPACKE_ssytri2x_work( int matrix_layout, char uplo, lapack_int n,
                                   const lapack_int* ipiv, float* work,
                                   lapack_int nb );
 lapack_int LAPACKE_ssytrs2( int matrix_layout, char uplo, lapack_int n,
-                            lapack_int nrhs, float* a, lapack_int lda,
+                            lapack_int nrhs, const float* a, lapack_int lda,
                             const lapack_int* ipiv, float* b, lapack_int ldb );
 lapack_int LAPACKE_ssytrs2_work( int matrix_layout, char uplo, lapack_int n,
-                                 lapack_int nrhs, float* a,
+                                 lapack_int nrhs, const float* a,
                                  lapack_int lda, const lapack_int* ipiv,
                                  float* b, lapack_int ldb, float* work );
 lapack_int LAPACKE_zbbcsd( int matrix_layout, char jobu1, char jobu2,
@@ -10898,11 +10898,11 @@ lapack_int LAPACKE_zsytri2x_work( int matrix_layout, char uplo, lapack_int n,
                                   const lapack_int* ipiv,
                                   lapack_complex_double* work, lapack_int nb );
 lapack_int LAPACKE_zsytrs2( int matrix_layout, char uplo, lapack_int n,
-                            lapack_int nrhs, lapack_complex_double* a,
+                            lapack_int nrhs, const lapack_complex_double* a,
                             lapack_int lda, const lapack_int* ipiv,
                             lapack_complex_double* b, lapack_int ldb );
 lapack_int LAPACKE_zsytrs2_work( int matrix_layout, char uplo, lapack_int n,
-                                 lapack_int nrhs, lapack_complex_double* a,
+                                 lapack_int nrhs, const lapack_complex_double* a,
                                  lapack_int lda, const lapack_int* ipiv,
                                  lapack_complex_double* b, lapack_int ldb,
                                  lapack_complex_double* work );

lapack-netlib/LAPACKE/src/lapacke_cgesvd_work.c

@@ -56,6 +56,8 @@ lapack_int LAPACKE_cgesvd_work( int matrix_layout, char jobu, char jobvt,
                              ( LAPACKE_lsame( jobu, 's' ) ? MIN(m,n) : 1);
         lapack_int nrows_vt = LAPACKE_lsame( jobvt, 'a' ) ? n :
                               ( LAPACKE_lsame( jobvt, 's' ) ? MIN(m,n) : 1);
+        lapack_int ncols_vt = ( LAPACKE_lsame( jobvt, 'a' ) ||
+                               LAPACKE_lsame( jobvt, 's' ) ) ? n : 1;
         lapack_int lda_t = MAX(1,m);
         lapack_int ldu_t = MAX(1,nrows_u);
         lapack_int ldvt_t = MAX(1,nrows_vt);
@@ -73,7 +75,7 @@ lapack_int LAPACKE_cgesvd_work( int matrix_layout, char jobu, char jobvt,
             LAPACKE_xerbla( "LAPACKE_cgesvd_work", info );
             return info;
         }
-        if( ldvt < n ) {
+        if( ldvt < ncols_vt ) {
             info = -12;
             LAPACKE_xerbla( "LAPACKE_cgesvd_work", info );
             return info;

lapack-netlib/LAPACKE/src/lapacke_cheev_work.c

@@ -78,7 +78,7 @@ lapack_int LAPACKE_cheev_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_cge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_che_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_cheevd_2stage_work.c

@@ -79,7 +79,7 @@ lapack_int LAPACKE_cheevd_2stage_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_cge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_che_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda ); 

lapack-netlib/LAPACKE/src/lapacke_cheevd_work.c

@@ -79,7 +79,7 @@ lapack_int LAPACKE_cheevd_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_cge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else { 
             LAPACKE_che_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_chegst.c

@@ -35,7 +35,7 @@
 
 lapack_int LAPACKE_chegst( int matrix_layout, lapack_int itype, char uplo,
                            lapack_int n, lapack_complex_float* a,
-                           lapack_int lda, lapack_complex_float* b,
+                           lapack_int lda, const lapack_complex_float* b,
                            lapack_int ldb )
 {
     if( matrix_layout != LAPACK_COL_MAJOR && matrix_layout != LAPACK_ROW_MAJOR ) {

lapack-netlib/LAPACKE/src/lapacke_chegst_work.c

@@ -35,7 +35,7 @@
 
 lapack_int LAPACKE_chegst_work( int matrix_layout, lapack_int itype, char uplo,
                                 lapack_int n, lapack_complex_float* a,
-                                lapack_int lda, lapack_complex_float* b,
+                                lapack_int lda, const lapack_complex_float* b,
                                 lapack_int ldb )
 {
     lapack_int info = 0;

lapack-netlib/LAPACKE/src/lapacke_chegv.c

@@ -50,10 +50,10 @@ lapack_int LAPACKE_chegv( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_chegv_2stage.c

@@ -50,10 +50,10 @@ lapack_int LAPACKE_chegv_2stage( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_chegvd.c

@@ -55,10 +55,10 @@ lapack_int LAPACKE_chegvd( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_chegvx.c

@@ -60,7 +60,7 @@ lapack_int LAPACKE_chegvx( int matrix_layout, lapack_int itype, char jobz,
         if( LAPACKE_s_nancheck( 1, &abstol, 1 ) ) {
             return -15;
         }
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -9;
         }
         if( LAPACKE_lsame( range, 'v' ) ) {

lapack-netlib/LAPACKE/src/lapacke_chetri2x.c

@@ -46,7 +46,7 @@ lapack_int LAPACKE_chetri2x( int matrix_layout, char uplo, lapack_int n,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -4;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_clacpy_work.c

@@ -42,9 +42,6 @@ lapack_int LAPACKE_clacpy_work( int matrix_layout, char uplo, lapack_int m,
     if( matrix_layout == LAPACK_COL_MAJOR ) {
         /* Call LAPACK function and adjust info */
         LAPACK_clacpy( &uplo, &m, &n, a, &lda, b, &ldb );
-        if( info < 0 ) {
-            info = info - 1;
-        }
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
         lapack_int lda_t = MAX(1,m);
         lapack_int ldb_t = MAX(1,m);

lapack-netlib/LAPACKE/src/lapacke_clantr_work.c

@@ -41,45 +41,46 @@ float LAPACKE_clantr_work( int matrix_layout, char norm, char uplo,
     lapack_int info = 0;
     float res = 0.;
     if( matrix_layout == LAPACK_COL_MAJOR ) {
-        /* Call LAPACK function and adjust info */
+        /* Call LAPACK function */
         res = LAPACK_clantr( &norm, &uplo, &diag, &m, &n, a, &lda, work );
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
-        lapack_int lda_t = MAX(1,m);
-        lapack_complex_float* a_t = NULL;
         float* work_lapack = NULL;
+        char norm_lapack;
+        char uplo_lapack;
         /* Check leading dimension(s) */
         if( lda < n ) {
             info = -8;
             LAPACKE_xerbla( "LAPACKE_clantr_work", info );
             return info;
         }
-        /* Allocate memory for temporary array(s) */
+        if( LAPACKE_lsame( norm, '1' ) || LAPACKE_lsame( norm, 'o' ) ) {
-        a_t = (lapack_complex_float*)
+            norm_lapack = 'i';
-            LAPACKE_malloc( sizeof(lapack_complex_float) * lda_t * MAX(1,MAX(m,n)) );
+        } else if( LAPACKE_lsame( norm, 'i' ) ) {
-        if( a_t == NULL ) {
+            norm_lapack = '1';
-            info = LAPACK_TRANSPOSE_MEMORY_ERROR;
+        } else {
-            goto exit_level_0;
+            norm_lapack = norm;
+        }
+        if( LAPACKE_lsame( uplo, 'u' ) ) {
+            uplo_lapack = 'l';
+        } else {
+            uplo_lapack = 'u';
         }
         /* Allocate memory for work array(s) */
-        if( LAPACKE_lsame( norm, 'i' ) ) {
+        if( LAPACKE_lsame( norm_lapack, 'i' ) ) {
-            work_lapack = (float*)LAPACKE_malloc( sizeof(float) * MAX(1,m) );
+            work_lapack = (float*)LAPACKE_malloc( sizeof(float) * MAX(1,n) );
             if( work_lapack == NULL ) {
                 info = LAPACK_WORK_MEMORY_ERROR;
-                goto exit_level_1;
+                goto exit_level_0;
             }
         }
-        /* Transpose input matrices */
+        /* Call LAPACK function */
-        LAPACKE_ctr_trans( matrix_layout, uplo, diag, MAX(m,n), a, lda, a_t, lda_t );
+        res = LAPACK_clantr( &norm_lapack, &uplo_lapack, &diag, &n, &m, a, &lda, work_lapack );
-        /* Call LAPACK function and adjust info */
-        res = LAPACK_clantr( &norm, &uplo, &diag, &m, &n, a_t, &lda_t, work_lapack );
         /* Release memory and exit */
         if( work_lapack ) {
             LAPACKE_free( work_lapack );
         }
-exit_level_1:
-        LAPACKE_free( a_t );
 exit_level_0:
-        if( info == LAPACK_TRANSPOSE_MEMORY_ERROR ) {
+        if( info == LAPACK_WORK_MEMORY_ERROR ) {
             LAPACKE_xerbla( "LAPACKE_clantr_work", info );
         }
     } else {

lapack-netlib/LAPACKE/src/lapacke_claset_work.c

@@ -42,9 +42,6 @@ lapack_int LAPACKE_claset_work( int matrix_layout, char uplo, lapack_int m,
     if( matrix_layout == LAPACK_COL_MAJOR ) {
         /* Call LAPACK function and adjust info */
         LAPACK_claset( &uplo, &m, &n, &alpha, &beta, a, &lda );
-        if( info < 0 ) {
-            info = info - 1;
-        }
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
         lapack_int lda_t = MAX(1,m);
         lapack_complex_float* a_t = NULL;

lapack-netlib/LAPACKE/src/lapacke_csyconv.c

@@ -45,7 +45,7 @@ lapack_int LAPACKE_csyconv( int matrix_layout, char uplo, char way, lapack_int n
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_csy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -5;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_csytrs2.c

@@ -34,7 +34,7 @@
 #include "lapacke_utils.h"
 
 lapack_int LAPACKE_csytrs2( int matrix_layout, char uplo, lapack_int n,
-                            lapack_int nrhs, lapack_complex_float* a,
+                            lapack_int nrhs, const lapack_complex_float* a,
                             lapack_int lda, const lapack_int* ipiv,
                             lapack_complex_float* b, lapack_int ldb )
 {

lapack-netlib/LAPACKE/src/lapacke_csytrs2_work.c

@@ -34,7 +34,7 @@
 #include "lapacke_utils.h"
 
 lapack_int LAPACKE_csytrs2_work( int matrix_layout, char uplo, lapack_int n,
-                                 lapack_int nrhs, lapack_complex_float* a,
+                                 lapack_int nrhs, const lapack_complex_float* a,
                                  lapack_int lda, const lapack_int* ipiv,
                                  lapack_complex_float* b, lapack_int ldb,
                                  lapack_complex_float* work )

lapack-netlib/LAPACKE/src/lapacke_ctrttf.c

@@ -44,7 +44,7 @@ lapack_int LAPACKE_ctrttf( int matrix_layout, char transr, char uplo,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ctr_nancheck( matrix_layout, uplo, 'n', n, a, lda ) ) {
             return -5;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_ctrttp.c

@@ -44,7 +44,7 @@ lapack_int LAPACKE_ctrttp( int matrix_layout, char uplo, lapack_int n,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ctr_nancheck( matrix_layout, uplo, 'n', n, a, lda ) ) {
             return -4;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_cungtr.c

@@ -48,7 +48,7 @@ lapack_int LAPACKE_cungtr( int matrix_layout, char uplo, lapack_int n,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_cge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -4;
         }
         if( LAPACKE_c_nancheck( n-1, tau, 1 ) ) {

lapack-netlib/LAPACKE/src/lapacke_cunmtr.c

@@ -52,7 +52,7 @@ lapack_int LAPACKE_cunmtr( int matrix_layout, char side, char uplo, char trans,
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
         r = LAPACKE_lsame( side, 'l' ) ? m : n;
-        if( LAPACKE_cge_nancheck( matrix_layout, r, r, a, lda ) ) {
+        if( LAPACKE_che_nancheck( matrix_layout, uplo, r, a, lda ) ) {
             return -7;
         }
         if( LAPACKE_cge_nancheck( matrix_layout, m, n, c, ldc ) ) {

lapack-netlib/LAPACKE/src/lapacke_dgesvd_work.c

@@ -54,6 +54,8 @@ lapack_int LAPACKE_dgesvd_work( int matrix_layout, char jobu, char jobvt,
                              ( LAPACKE_lsame( jobu, 's' ) ? MIN(m,n) : 1);
         lapack_int nrows_vt = LAPACKE_lsame( jobvt, 'a' ) ? n :
                               ( LAPACKE_lsame( jobvt, 's' ) ? MIN(m,n) : 1);
+        lapack_int ncols_vt = ( LAPACKE_lsame( jobvt, 'a' ) ||
+                               LAPACKE_lsame( jobvt, 's' ) ) ? n : 1;
         lapack_int lda_t = MAX(1,m);
         lapack_int ldu_t = MAX(1,nrows_u);
         lapack_int ldvt_t = MAX(1,nrows_vt);
@@ -71,7 +73,7 @@ lapack_int LAPACKE_dgesvd_work( int matrix_layout, char jobu, char jobvt,
             LAPACKE_xerbla( "LAPACKE_dgesvd_work", info );
             return info;
         }
-        if( ldvt < n ) {
+        if( ldvt < ncols_vt ) {
             info = -12;
             LAPACKE_xerbla( "LAPACKE_dgesvd_work", info );
             return info;

lapack-netlib/LAPACKE/src/lapacke_dlacpy_work.c

@@ -41,9 +41,6 @@ lapack_int LAPACKE_dlacpy_work( int matrix_layout, char uplo, lapack_int m,
     if( matrix_layout == LAPACK_COL_MAJOR ) {
         /* Call LAPACK function and adjust info */
         LAPACK_dlacpy( &uplo, &m, &n, a, &lda, b, &ldb );
-        if( info < 0 ) {
-            info = info - 1;
-        }
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
         lapack_int lda_t = MAX(1,m);
         lapack_int ldb_t = MAX(1,m);

lapack-netlib/LAPACKE/src/lapacke_dlantr_work.c

@@ -40,44 +40,46 @@ double LAPACKE_dlantr_work( int matrix_layout, char norm, char uplo,
     lapack_int info = 0;
     double res = 0.;
     if( matrix_layout == LAPACK_COL_MAJOR ) {
-        /* Call LAPACK function and adjust info */
+        /* Call LAPACK function */
         res = LAPACK_dlantr( &norm, &uplo, &diag, &m, &n, a, &lda, work );
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
-        lapack_int lda_t = MAX(1,m);
-        double* a_t = NULL;
         double* work_lapack = NULL;
+        char norm_lapack;
+        char uplo_lapack;
         /* Check leading dimension(s) */
         if( lda < n ) {
             info = -8;
             LAPACKE_xerbla( "LAPACKE_dlantr_work", info );
             return info;
         }
-        /* Allocate memory for temporary array(s) */
+        if( LAPACKE_lsame( norm, '1' ) || LAPACKE_lsame( norm, 'o' ) ) {
-        a_t = (double*)LAPACKE_malloc( sizeof(double) * lda_t * MAX(1,MAX(m,n)) );
+            norm_lapack = 'i';
-        if( a_t == NULL ) {
+        } else if( LAPACKE_lsame( norm, 'i' ) ) {
-            info = LAPACK_TRANSPOSE_MEMORY_ERROR;
+            norm_lapack = '1';
-            goto exit_level_0;
+        } else {
+            norm_lapack = norm;
+        }
+        if( LAPACKE_lsame( uplo, 'u' ) ) {
+            uplo_lapack = 'l';
+        } else {
+            uplo_lapack = 'u';
         }
         /* Allocate memory for work array(s) */
-        if( LAPACKE_lsame( norm, 'i' ) ) {
+        if( LAPACKE_lsame( norm_lapack, 'i' ) ) {
-            work_lapack = (double*)LAPACKE_malloc( sizeof(double) * MAX(1,m) );
+            work_lapack = (double*)LAPACKE_malloc( sizeof(double) * MAX(1,n) );
             if( work_lapack == NULL ) {
                 info = LAPACK_WORK_MEMORY_ERROR;
-                goto exit_level_1;
+                goto exit_level_0;
             }
         }
-        /* Transpose input matrices */
+        /* Call LAPACK function */
-        LAPACKE_dtr_trans( matrix_layout, uplo, diag, MAX(m,n), a, lda, a_t, lda_t );
+        res = LAPACK_dlantr( &norm_lapack, &uplo_lapack, &diag, &n, &m, a, &lda, work_lapack );
-        /* Call LAPACK function and adjust info */
-        res = LAPACK_dlantr( &norm, &uplo, &diag, &m, &n, a_t, &lda_t, work_lapack );
         /* Release memory and exit */
         if( work_lapack ) {
             LAPACKE_free( work_lapack );
         }
-exit_level_1:
-        LAPACKE_free( a_t );
 exit_level_0:
-        if( info == LAPACK_TRANSPOSE_MEMORY_ERROR ) {
+        if( info == LAPACK_WORK_MEMORY_ERROR ) {
             LAPACKE_xerbla( "LAPACKE_dlantr_work", info );
         }
     } else {

lapack-netlib/LAPACKE/src/lapacke_dlaset_work.c

@@ -41,9 +41,6 @@ lapack_int LAPACKE_dlaset_work( int matrix_layout, char uplo, lapack_int m,
     if( matrix_layout == LAPACK_COL_MAJOR ) {
         /* Call LAPACK function and adjust info */
         LAPACK_dlaset( &uplo, &m, &n, &alpha, &beta, a, &lda );
-        if( info < 0 ) {
-            info = info - 1;
-        }
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
         lapack_int lda_t = MAX(1,m);
         double* a_t = NULL;

lapack-netlib/LAPACKE/src/lapacke_dorgtr.c

@@ -47,7 +47,7 @@ lapack_int LAPACKE_dorgtr( int matrix_layout, char uplo, lapack_int n, double* a
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -4;
         }
         if( LAPACKE_d_nancheck( n-1, tau, 1 ) ) {

lapack-netlib/LAPACKE/src/lapacke_dormtr.c

@@ -51,7 +51,7 @@ lapack_int LAPACKE_dormtr( int matrix_layout, char side, char uplo, char trans,
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
         r = LAPACKE_lsame( side, 'l' ) ? m : n;
-        if( LAPACKE_dge_nancheck( matrix_layout, r, r, a, lda ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, r, a, lda ) ) {
             return -7;
         }
         if( LAPACKE_dge_nancheck( matrix_layout, m, n, c, ldc ) ) {

lapack-netlib/LAPACKE/src/lapacke_dsyconv.c

@@ -43,7 +43,7 @@ lapack_int LAPACKE_dsyconv( int matrix_layout, char uplo, char way, lapack_int n
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -5;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_dsyev_work.c

@@ -72,7 +72,7 @@ lapack_int LAPACKE_dsyev_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_dsy_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_dsyevd_2stage_work.c

@@ -76,7 +76,7 @@ lapack_int LAPACKE_dsyevd_2stage_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_dsy_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_dsyevd_work.c

@@ -76,7 +76,7 @@ lapack_int LAPACKE_dsyevd_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_dge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_dsy_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_dsygst.c

@@ -47,7 +47,7 @@ lapack_int LAPACKE_dsygst( int matrix_layout, lapack_int itype, char uplo,
         if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -5;
         }
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -7;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_dsygv.c

@@ -48,10 +48,10 @@ lapack_int LAPACKE_dsygv( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_dsygv_2stage.c

@@ -48,10 +48,10 @@ lapack_int LAPACKE_dsygv_2stage( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_dsygvd.c

@@ -51,10 +51,10 @@ lapack_int LAPACKE_dsygvd( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_dsygvx.c

@@ -58,7 +58,7 @@ lapack_int LAPACKE_dsygvx( int matrix_layout, lapack_int itype, char jobz,
         if( LAPACKE_d_nancheck( 1, &abstol, 1 ) ) {
             return -15;
         }
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_dsy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -9;
         }
         if( LAPACKE_lsame( range, 'v' ) ) {

lapack-netlib/LAPACKE/src/lapacke_dsytrs2.c

@@ -34,7 +34,7 @@
 #include "lapacke_utils.h"
 
 lapack_int LAPACKE_dsytrs2( int matrix_layout, char uplo, lapack_int n,
-                            lapack_int nrhs, double* a, lapack_int lda,
+                            lapack_int nrhs, const double* a, lapack_int lda,
                             const lapack_int* ipiv, double* b, lapack_int ldb )
 {
     lapack_int info = 0;

lapack-netlib/LAPACKE/src/lapacke_dsytrs2_work.c

@@ -34,7 +34,7 @@
 #include "lapacke_utils.h"
 
 lapack_int LAPACKE_dsytrs2_work( int matrix_layout, char uplo, lapack_int n,
-                                 lapack_int nrhs, double* a,
+                                 lapack_int nrhs, const double* a,
                                  lapack_int lda, const lapack_int* ipiv,
                                  double* b, lapack_int ldb, double* work )
 {

lapack-netlib/LAPACKE/src/lapacke_dtrttf.c

@@ -44,7 +44,7 @@ lapack_int LAPACKE_dtrttf( int matrix_layout, char transr, char uplo,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_dtr_nancheck( matrix_layout, uplo, 'n', n, a, lda ) ) {
             return -5;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_dtrttp.c

@@ -43,7 +43,7 @@ lapack_int LAPACKE_dtrttp( int matrix_layout, char uplo, lapack_int n,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_dge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_dtr_nancheck( matrix_layout, uplo, 'n', n, a, lda ) ) {
             return -4;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_sgesvd_work.c

@@ -54,6 +54,8 @@ lapack_int LAPACKE_sgesvd_work( int matrix_layout, char jobu, char jobvt,
                              ( LAPACKE_lsame( jobu, 's' ) ? MIN(m,n) : 1);
         lapack_int nrows_vt = LAPACKE_lsame( jobvt, 'a' ) ? n :
                               ( LAPACKE_lsame( jobvt, 's' ) ? MIN(m,n) : 1);
+        lapack_int ncols_vt = ( LAPACKE_lsame( jobvt, 'a' ) ||
+                               LAPACKE_lsame( jobvt, 's' ) ) ? n : 1;
         lapack_int lda_t = MAX(1,m);
         lapack_int ldu_t = MAX(1,nrows_u);
         lapack_int ldvt_t = MAX(1,nrows_vt);
@@ -71,7 +73,7 @@ lapack_int LAPACKE_sgesvd_work( int matrix_layout, char jobu, char jobvt,
             LAPACKE_xerbla( "LAPACKE_sgesvd_work", info );
             return info;
         }
-        if( ldvt < n ) {
+        if( ldvt < ncols_vt ) {
             info = -12;
             LAPACKE_xerbla( "LAPACKE_sgesvd_work", info );
             return info;

lapack-netlib/LAPACKE/src/lapacke_slacpy_work.c

@@ -41,9 +41,6 @@ lapack_int LAPACKE_slacpy_work( int matrix_layout, char uplo, lapack_int m,
     if( matrix_layout == LAPACK_COL_MAJOR ) {
         /* Call LAPACK function and adjust info */
         LAPACK_slacpy( &uplo, &m, &n, a, &lda, b, &ldb );
-        if( info < 0 ) {
-            info = info - 1;
-        }
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
         lapack_int lda_t = MAX(1,m);
         lapack_int ldb_t = MAX(1,m);

lapack-netlib/LAPACKE/src/lapacke_slantr_work.c

@@ -40,44 +40,46 @@ float LAPACKE_slantr_work( int matrix_layout, char norm, char uplo,
     lapack_int info = 0;
     float res = 0.;
     if( matrix_layout == LAPACK_COL_MAJOR ) {
-        /* Call LAPACK function and adjust info */
+        /* Call LAPACK function */
         res = LAPACK_slantr( &norm, &uplo, &diag, &m, &n, a, &lda, work );
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
-        lapack_int lda_t = MAX(1,m);
-        float* a_t = NULL;
         float* work_lapack = NULL;
+        char norm_lapack;
+        char uplo_lapack;
         /* Check leading dimension(s) */
         if( lda < n ) {
             info = -8;
             LAPACKE_xerbla( "LAPACKE_slantr_work", info );
             return info;
         }
-        /* Allocate memory for temporary array(s) */
+        if( LAPACKE_lsame( norm, '1' ) || LAPACKE_lsame( norm, 'o' ) ) {
-        a_t = (float*)LAPACKE_malloc( sizeof(float) * lda_t * MAX(1,MAX(m,n)) );
+            norm_lapack = 'i';
-        if( a_t == NULL ) {
+        } else if( LAPACKE_lsame( norm, 'i' ) ) {
-            info = LAPACK_TRANSPOSE_MEMORY_ERROR;
+            norm_lapack = '1';
-            goto exit_level_0;
+        } else {
+            norm_lapack = norm;
+        }
+        if( LAPACKE_lsame( uplo, 'u' ) ) {
+            uplo_lapack = 'l';
+        } else {
+            uplo_lapack = 'u';
         }
         /* Allocate memory for work array(s) */
-        if( LAPACKE_lsame( norm, 'i' ) ) {
+        if( LAPACKE_lsame( norm_lapack, 'i' ) ) {
-            work_lapack = (float*)LAPACKE_malloc( sizeof(float) * MAX(1,m) );
+            work_lapack = (float*)LAPACKE_malloc( sizeof(float) * MAX(1,n) );
             if( work_lapack == NULL ) {
                 info = LAPACK_WORK_MEMORY_ERROR;
-                goto exit_level_1;
+                goto exit_level_0;
             }
         }
-        /* Transpose input matrices */
+        /* Call LAPACK function */
-        LAPACKE_str_trans( matrix_layout, uplo, diag, MAX(m,n), a, lda, a_t, lda_t );
+        res = LAPACK_slantr( &norm_lapack, &uplo_lapack, &diag, &n, &m, a, &lda, work_lapack );
-        /* Call LAPACK function and adjust info */
-        res = LAPACK_slantr( &norm, &uplo, &diag, &m, &n, a_t, &lda_t, work_lapack );
         /* Release memory and exit */
         if( work_lapack ) {
             LAPACKE_free( work_lapack );
         }
-exit_level_1:
-        LAPACKE_free( a_t );
 exit_level_0:
-        if( info == LAPACK_TRANSPOSE_MEMORY_ERROR ) {
+        if( info == LAPACK_WORK_MEMORY_ERROR ) {
             LAPACKE_xerbla( "LAPACKE_slantr_work", info );
         }
     } else {

lapack-netlib/LAPACKE/src/lapacke_slaset_work.c

@@ -41,9 +41,6 @@ lapack_int LAPACKE_slaset_work( int matrix_layout, char uplo, lapack_int m,
     if( matrix_layout == LAPACK_COL_MAJOR ) {
         /* Call LAPACK function and adjust info */
         LAPACK_slaset( &uplo, &m, &n, &alpha, &beta, a, &lda );
-        if( info < 0 ) {
-            info = info - 1;
-        }
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
         lapack_int lda_t = MAX(1,m);
         float* a_t = NULL;

lapack-netlib/LAPACKE/src/lapacke_sorgtr.c

@@ -47,7 +47,7 @@ lapack_int LAPACKE_sorgtr( int matrix_layout, char uplo, lapack_int n, float* a,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -4;
         }
         if( LAPACKE_s_nancheck( n-1, tau, 1 ) ) {

lapack-netlib/LAPACKE/src/lapacke_sormtr.c

@@ -51,7 +51,7 @@ lapack_int LAPACKE_sormtr( int matrix_layout, char side, char uplo, char trans,
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
         r = LAPACKE_lsame( side, 'l' ) ? m : n;
-        if( LAPACKE_sge_nancheck( matrix_layout, r, r, a, lda ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, r, a, lda ) ) {
             return -7;
         }
         if( LAPACKE_sge_nancheck( matrix_layout, m, n, c, ldc ) ) {

lapack-netlib/LAPACKE/src/lapacke_ssyconv.c

@@ -43,7 +43,7 @@ lapack_int LAPACKE_ssyconv( int matrix_layout, char uplo, char way, lapack_int n
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -5;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_ssyev_work.c

@@ -72,7 +72,7 @@ lapack_int LAPACKE_ssyev_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_ssy_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_ssyevd_2stage_work.c

@@ -76,7 +76,7 @@ lapack_int LAPACKE_ssyevd_2stage_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_ssy_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda ); 

lapack-netlib/LAPACKE/src/lapacke_ssyevd_work.c

@@ -76,7 +76,7 @@ lapack_int LAPACKE_ssyevd_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_sge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_ssy_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_ssygst.c

@@ -47,7 +47,7 @@ lapack_int LAPACKE_ssygst( int matrix_layout, lapack_int itype, char uplo,
         if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -5;
         }
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -7;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_ssygv.c

@@ -48,10 +48,10 @@ lapack_int LAPACKE_ssygv( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_ssygv_2stage.c

@@ -48,10 +48,10 @@ lapack_int LAPACKE_ssygv_2stage( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_ssygvd.c

@@ -51,10 +51,10 @@ lapack_int LAPACKE_ssygvd( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_ssygvx.c

@@ -58,7 +58,7 @@ lapack_int LAPACKE_ssygvx( int matrix_layout, lapack_int itype, char jobz,
         if( LAPACKE_s_nancheck( 1, &abstol, 1 ) ) {
             return -15;
         }
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_ssy_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -9;
         }
         if( LAPACKE_lsame( range, 'v' ) ) {

lapack-netlib/LAPACKE/src/lapacke_ssytrs2.c

@@ -34,7 +34,7 @@
 #include "lapacke_utils.h"
 
 lapack_int LAPACKE_ssytrs2( int matrix_layout, char uplo, lapack_int n,
-                            lapack_int nrhs, float* a, lapack_int lda,
+                            lapack_int nrhs, const float* a, lapack_int lda,
                             const lapack_int* ipiv, float* b, lapack_int ldb )
 {
     lapack_int info = 0;

lapack-netlib/LAPACKE/src/lapacke_ssytrs2_work.c

@@ -34,7 +34,7 @@
 #include "lapacke_utils.h"
 
 lapack_int LAPACKE_ssytrs2_work( int matrix_layout, char uplo, lapack_int n,
-                                 lapack_int nrhs, float* a,
+                                 lapack_int nrhs, const float* a,
                                  lapack_int lda, const lapack_int* ipiv,
                                  float* b, lapack_int ldb, float* work )
 {

lapack-netlib/LAPACKE/src/lapacke_strttf.c

@@ -44,7 +44,7 @@ lapack_int LAPACKE_strttf( int matrix_layout, char transr, char uplo,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_str_nancheck( matrix_layout, uplo, 'n', n, a, lda ) ) {
             return -5;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_strttp.c

@@ -43,7 +43,7 @@ lapack_int LAPACKE_strttp( int matrix_layout, char uplo, lapack_int n,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_sge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_str_nancheck( matrix_layout, uplo, 'n', n, a, lda ) ) {
             return -4;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_zgesvd_work.c

@@ -56,6 +56,8 @@ lapack_int LAPACKE_zgesvd_work( int matrix_layout, char jobu, char jobvt,
                              ( LAPACKE_lsame( jobu, 's' ) ? MIN(m,n) : 1);
         lapack_int nrows_vt = LAPACKE_lsame( jobvt, 'a' ) ? n :
                               ( LAPACKE_lsame( jobvt, 's' ) ? MIN(m,n) : 1);
+        lapack_int ncols_vt = ( LAPACKE_lsame( jobvt, 'a' ) ||
+                               LAPACKE_lsame( jobvt, 's' ) ) ? n : 1;
         lapack_int lda_t = MAX(1,m);
         lapack_int ldu_t = MAX(1,nrows_u);
         lapack_int ldvt_t = MAX(1,nrows_vt);
@@ -73,7 +75,7 @@ lapack_int LAPACKE_zgesvd_work( int matrix_layout, char jobu, char jobvt,
             LAPACKE_xerbla( "LAPACKE_zgesvd_work", info );
             return info;
         }
-        if( ldvt < n ) {
+        if( ldvt < ncols_vt ) {
             info = -12;
             LAPACKE_xerbla( "LAPACKE_zgesvd_work", info );
             return info;

lapack-netlib/LAPACKE/src/lapacke_zheev_work.c

@@ -78,7 +78,7 @@ lapack_int LAPACKE_zheev_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_zhe_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_zheevd_2stage_work.c

@@ -79,7 +79,7 @@ lapack_int LAPACKE_zheevd_2stage_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else {
             LAPACKE_zhe_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_zheevd_work.c

@@ -79,7 +79,7 @@ lapack_int LAPACKE_zheevd_work( int matrix_layout, char jobz, char uplo,
             info = info - 1;
         }
         /* Transpose output matrices */
-        if ( jobz == 'V') {
+        if ( jobz == 'V' || jobz == 'v' ) {
             LAPACKE_zge_trans( LAPACK_COL_MAJOR, n, n, a_t, lda_t, a, lda );
         } else { 
             LAPACKE_zhe_trans( LAPACK_COL_MAJOR, uplo, n, a_t, lda_t, a, lda );

lapack-netlib/LAPACKE/src/lapacke_zhegst.c

@@ -35,7 +35,7 @@
 
 lapack_int LAPACKE_zhegst( int matrix_layout, lapack_int itype, char uplo,
                            lapack_int n, lapack_complex_double* a,
-                           lapack_int lda, lapack_complex_double* b,
+                           lapack_int lda, const lapack_complex_double* b,
                            lapack_int ldb )
 {
     if( matrix_layout != LAPACK_COL_MAJOR && matrix_layout != LAPACK_ROW_MAJOR ) {

lapack-netlib/LAPACKE/src/lapacke_zhegst_work.c

@@ -35,7 +35,7 @@
 
 lapack_int LAPACKE_zhegst_work( int matrix_layout, lapack_int itype, char uplo,
                                 lapack_int n, lapack_complex_double* a,
-                                lapack_int lda, lapack_complex_double* b,
+                                lapack_int lda, const lapack_complex_double* b,
                                 lapack_int ldb )
 {
     lapack_int info = 0;

lapack-netlib/LAPACKE/src/lapacke_zhegv.c

@@ -50,10 +50,10 @@ lapack_int LAPACKE_zhegv( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_zhegv_2stage.c

@@ -50,10 +50,10 @@ lapack_int LAPACKE_zhegv_2stage( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_zhegvd.c

@@ -55,10 +55,10 @@ lapack_int LAPACKE_zhegvd( int matrix_layout, lapack_int itype, char jobz,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -6;
         }
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -8;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_zhegvx.c

@@ -61,7 +61,7 @@ lapack_int LAPACKE_zhegvx( int matrix_layout, lapack_int itype, char jobz,
         if( LAPACKE_d_nancheck( 1, &abstol, 1 ) ) {
             return -15;
         }
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, b, ldb ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, b, ldb ) ) {
             return -9;
         }
         if( LAPACKE_lsame( range, 'v' ) ) {

lapack-netlib/LAPACKE/src/lapacke_zhetri2x.c

@@ -46,7 +46,7 @@ lapack_int LAPACKE_zhetri2x( int matrix_layout, char uplo, lapack_int n,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -4;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_zlacpy_work.c

@@ -42,9 +42,6 @@ lapack_int LAPACKE_zlacpy_work( int matrix_layout, char uplo, lapack_int m,
     if( matrix_layout == LAPACK_COL_MAJOR ) {
         /* Call LAPACK function and adjust info */
         LAPACK_zlacpy( &uplo, &m, &n, a, &lda, b, &ldb );
-        if( info < 0 ) {
-            info = info - 1;
-        }
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
         lapack_int lda_t = MAX(1,m);
         lapack_int ldb_t = MAX(1,m);

lapack-netlib/LAPACKE/src/lapacke_zlantr_work.c

@@ -41,45 +41,46 @@ double LAPACKE_zlantr_work( int matrix_layout, char norm, char uplo,
     lapack_int info = 0;
     double res = 0.;
     if( matrix_layout == LAPACK_COL_MAJOR ) {
-        /* Call LAPACK function and adjust info */
+        /* Call LAPACK function */
         res = LAPACK_zlantr( &norm, &uplo, &diag, &m, &n, a, &lda, work );
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
-        lapack_int lda_t = MAX(1,m);
-        lapack_complex_double* a_t = NULL;
         double* work_lapack = NULL;
+        char norm_lapack;
+        char uplo_lapack;
         /* Check leading dimension(s) */
         if( lda < n ) {
             info = -8;
             LAPACKE_xerbla( "LAPACKE_zlantr_work", info );
             return info;
         }
-        /* Allocate memory for temporary array(s) */
+        if( LAPACKE_lsame( norm, '1' ) || LAPACKE_lsame( norm, 'o' ) ) {
-        a_t = (lapack_complex_double*)
+            norm_lapack = 'i';
-            LAPACKE_malloc( sizeof(lapack_complex_double) * lda_t * MAX(1,MAX(m,n)) );
+        } else if( LAPACKE_lsame( norm, 'i' ) ) {
-        if( a_t == NULL ) {
+            norm_lapack = '1';
-            info = LAPACK_TRANSPOSE_MEMORY_ERROR;
+        } else {
-            goto exit_level_0;
+            norm_lapack = norm;
+        }
+        if( LAPACKE_lsame( uplo, 'u' ) ) {
+            uplo_lapack = 'l';
+        } else {
+            uplo_lapack = 'u';
         }
         /* Allocate memory for work array(s) */
-        if( LAPACKE_lsame( norm, 'i' ) ) {
+        if( LAPACKE_lsame( norm_lapack, 'i' ) ) {
-            work_lapack = (double*)LAPACKE_malloc( sizeof(double) * MAX(1,m) );
+            work_lapack = (double*)LAPACKE_malloc( sizeof(double) * MAX(1,n) );
             if( work_lapack == NULL ) {
                 info = LAPACK_WORK_MEMORY_ERROR;
-                goto exit_level_1;
+                goto exit_level_0;
             }
         }
-        /* Transpose input matrices */
+        /* Call LAPACK function */
-        LAPACKE_ztr_trans( matrix_layout, uplo, diag, MAX(m,n), a, lda, a_t, lda_t );
+        res = LAPACK_zlantr( &norm_lapack, &uplo_lapack, &diag, &n, &m, a, &lda, work_lapack );
-        /* Call LAPACK function and adjust info */
-        res = LAPACK_zlantr( &norm, &uplo, &diag, &m, &n, a_t, &lda_t, work_lapack );
         /* Release memory and exit */
         if( work_lapack ) {
             LAPACKE_free( work_lapack );
         }
-exit_level_1:
-        LAPACKE_free( a_t );
 exit_level_0:
-        if( info == LAPACK_TRANSPOSE_MEMORY_ERROR ) {
+        if( info == LAPACK_WORK_MEMORY_ERROR ) {
             LAPACKE_xerbla( "LAPACKE_zlantr_work", info );
         }
     } else {

lapack-netlib/LAPACKE/src/lapacke_zlaset_work.c

@@ -42,9 +42,6 @@ lapack_int LAPACKE_zlaset_work( int matrix_layout, char uplo, lapack_int m,
     if( matrix_layout == LAPACK_COL_MAJOR ) {
         /* Call LAPACK function and adjust info */
         LAPACK_zlaset( &uplo, &m, &n, &alpha, &beta, a, &lda );
-        if( info < 0 ) {
-            info = info - 1;
-        }
     } else if( matrix_layout == LAPACK_ROW_MAJOR ) {
         lapack_int lda_t = MAX(1,m);
         lapack_complex_double* a_t = NULL;

lapack-netlib/LAPACKE/src/lapacke_zsyconv.c

@@ -45,7 +45,7 @@ lapack_int LAPACKE_zsyconv( int matrix_layout, char uplo, char way, lapack_int n
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_zsy_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -5;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_zsytrs2.c

@@ -34,7 +34,7 @@
 #include "lapacke_utils.h"
 
 lapack_int LAPACKE_zsytrs2( int matrix_layout, char uplo, lapack_int n,
-                            lapack_int nrhs, lapack_complex_double* a,
+                            lapack_int nrhs, const lapack_complex_double* a,
                             lapack_int lda, const lapack_int* ipiv,
                             lapack_complex_double* b, lapack_int ldb )
 {

lapack-netlib/LAPACKE/src/lapacke_zsytrs2_work.c

@@ -35,7 +35,7 @@
 
 lapack_int LAPACKE_zsytrs2_work( int matrix_layout, char uplo, lapack_int n,
                                  lapack_int nrhs,
-                                 lapack_complex_double* a, lapack_int lda,
+                                 const lapack_complex_double* a, lapack_int lda,
                                  const lapack_int* ipiv,
                                  lapack_complex_double* b, lapack_int ldb,
                                  lapack_complex_double* work )

lapack-netlib/LAPACKE/src/lapacke_ztrttf.c

@@ -44,7 +44,7 @@ lapack_int LAPACKE_ztrttf( int matrix_layout, char transr, char uplo,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ztr_nancheck( matrix_layout, uplo, 'n', n, a, lda ) ) {
             return -5;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_ztrttp.c

@@ -44,7 +44,7 @@ lapack_int LAPACKE_ztrttp( int matrix_layout, char uplo, lapack_int n,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_ztr_nancheck( matrix_layout, uplo, 'n', n, a, lda ) ) {
             return -4;
         }
     }

lapack-netlib/LAPACKE/src/lapacke_zungtr.c

@@ -48,7 +48,7 @@ lapack_int LAPACKE_zungtr( int matrix_layout, char uplo, lapack_int n,
 #ifndef LAPACK_DISABLE_NAN_CHECK
     if( LAPACKE_get_nancheck() ) {
         /* Optionally check input matrices for NaNs */
-        if( LAPACKE_zge_nancheck( matrix_layout, n, n, a, lda ) ) {
+        if( LAPACKE_zhe_nancheck( matrix_layout, uplo, n, a, lda ) ) {
             return -4;
         }
         if( LAPACKE_z_nancheck( n-1, tau, 1 ) ) {

lapack-netlib/SRC/cgesdd.f

@@ -281,9 +281,9 @@
      $                   CUNGQR, CUNMBR, SBDSDC, SLASCL, XERBLA
 *     ..
 *     .. External Functions ..
-      LOGICAL            LSAME
+      LOGICAL            LSAME, SISNAN
       REAL               SLAMCH, CLANGE
-      EXTERNAL           LSAME, SLAMCH, CLANGE
+      EXTERNAL           LSAME, SLAMCH, CLANGE, SISNAN
 *     ..
 *     .. Intrinsic Functions ..
       INTRINSIC          INT, MAX, MIN, SQRT
@@ -647,6 +647,10 @@
 *     Scale A if max element outside range [SMLNUM,BIGNUM]
 *
       ANRM = CLANGE( 'M', M, N, A, LDA, DUM )
+      IF( SISNAN ( ANRM ) ) THEN
+          INFO = -4
+          RETURN
+      END IF
       ISCL = 0
       IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
          ISCL = 1

lapack-netlib/SRC/cggglm.f

@@ -271,8 +271,15 @@
 *
 *     Quick return if possible
 *
-      IF( N.EQ.0 )
+      IF( N.EQ.0 ) THEN
-     $   RETURN
+         DO I = 1, M
+            X(I) = CZERO
+         END DO
+         DO I = 1, P
+            Y(I) = CZERO
+         END DO
+         RETURN
+      END IF
 *
 *     Compute the GQR factorization of matrices A and B:
 *

lapack-netlib/SRC/chseqr.f

@@ -320,10 +320,10 @@
 *     .    CLAHQR because of insufficient subdiagonal scratch space.
 *     .    (This is a hard limit.) ====
       INTEGER            NTINY
-      PARAMETER          ( NTINY = 11 )
+      PARAMETER          ( NTINY = 15 )
 *
 *     ==== NL allocates some local workspace to help small matrices
-*     .    through a rare CLAHQR failure.  NL > NTINY = 11 is
+*     .    through a rare CLAHQR failure.  NL > NTINY = 15 is
 *     .    required and NL <= NMIN = ILAENV(ISPEC=12,...) is recom-
 *     .    mended.  (The default value of NMIN is 75.)  Using NL = 49
 *     .    allows up to six simultaneous shifts and a 16-by-16

lapack-netlib/SRC/claqr0.f

@@ -260,7 +260,7 @@
 *     .    CLAHQR because of insufficient subdiagonal scratch space.
 *     .    (This is a hard limit.) ====
       INTEGER            NTINY
-      PARAMETER          ( NTINY = 11 )
+      PARAMETER          ( NTINY = 15 )
 *
 *     ==== Exceptional deflation windows:  try to cure rare
 *     .    slow convergence by varying the size of the
@@ -355,22 +355,22 @@
          END IF
 *
 *        ==== NWR = recommended deflation window size.  At this
-*        .    point,  N .GT. NTINY = 11, so there is enough
+*        .    point,  N .GT. NTINY = 15, so there is enough
 *        .    subdiagonal workspace for NWR.GE.2 as required.
 *        .    (In fact, there is enough subdiagonal space for
-*        .    NWR.GE.3.) ====
+*        .    NWR.GE.4.) ====
 *
          NWR = ILAENV( 13, 'CLAQR0', JBCMPZ, N, ILO, IHI, LWORK )
          NWR = MAX( 2, NWR )
          NWR = MIN( IHI-ILO+1, ( N-1 ) / 3, NWR )
 *
 *        ==== NSR = recommended number of simultaneous shifts.
-*        .    At this point N .GT. NTINY = 11, so there is at
+*        .    At this point N .GT. NTINY = 15, so there is at
 *        .    enough subdiagonal workspace for NSR to be even
 *        .    and greater than or equal to two as required. ====
 *
          NSR = ILAENV( 15, 'CLAQR0', JBCMPZ, N, ILO, IHI, LWORK )
-         NSR = MIN( NSR, ( N+6 ) / 9, IHI-ILO )
+         NSR = MIN( NSR, ( N-3 ) / 6, IHI-ILO )
          NSR = MAX( 2, NSR-MOD( NSR, 2 ) )
 *
 *        ==== Estimate optimal workspace ====
@@ -418,7 +418,7 @@
 *        ==== NSMAX = the Largest number of simultaneous shifts
 *        .    for which there is sufficient workspace. ====
 *
-         NSMAX = MIN( ( N+6 ) / 9, 2*LWORK / 3 )
+         NSMAX = MIN( ( N-3 ) / 6, 2*LWORK / 3 )
          NSMAX = NSMAX - MOD( NSMAX, 2 )
 *
 *        ==== NDFL: an iteration count restarted at deflation. ====
@@ -558,7 +558,7 @@
 *
 *                 ==== Got NS/2 or fewer shifts? Use CLAQR4 or
 *                 .    CLAHQR on a trailing principal submatrix to
-*                 .    get more. (Since NS.LE.NSMAX.LE.(N+6)/9,
+*                 .    get more. (Since NS.LE.NSMAX.LE.(N-3)/6,
 *                 .    there is enough space below the subdiagonal
 *                 .    to fit an NS-by-NS scratch array.) ====
 *
@@ -659,7 +659,7 @@
 *              .      (NVE-by-KDU) vertical work WV arrow along
 *              .      the left-hand-edge. ====
 *
-               KDU = 3*NS - 3
+               KDU = 2*NS
                KU = N - KDU + 1
                KWH = KDU + 1
                NHO = ( N-KDU+1-4 ) - ( KDU+1 ) + 1

lapack-netlib/SRC/claqr4.f

@@ -270,7 +270,7 @@
 *     .    CLAHQR because of insufficient subdiagonal scratch space.
 *     .    (This is a hard limit.) ====
       INTEGER            NTINY
-      PARAMETER          ( NTINY = 11 )
+      PARAMETER          ( NTINY = 15 )
 *
 *     ==== Exceptional deflation windows:  try to cure rare
 *     .    slow convergence by varying the size of the
@@ -365,22 +365,22 @@
          END IF
 *
 *        ==== NWR = recommended deflation window size.  At this
-*        .    point,  N .GT. NTINY = 11, so there is enough
+*        .    point,  N .GT. NTINY = 15, so there is enough
 *        .    subdiagonal workspace for NWR.GE.2 as required.
 *        .    (In fact, there is enough subdiagonal space for
-*        .    NWR.GE.3.) ====
+*        .    NWR.GE.4.) ====
 *
          NWR = ILAENV( 13, 'CLAQR4', JBCMPZ, N, ILO, IHI, LWORK )
          NWR = MAX( 2, NWR )
          NWR = MIN( IHI-ILO+1, ( N-1 ) / 3, NWR )
 *
 *        ==== NSR = recommended number of simultaneous shifts.
-*        .    At this point N .GT. NTINY = 11, so there is at
+*        .    At this point N .GT. NTINY = 15, so there is at
 *        .    enough subdiagonal workspace for NSR to be even
 *        .    and greater than or equal to two as required. ====
 *
          NSR = ILAENV( 15, 'CLAQR4', JBCMPZ, N, ILO, IHI, LWORK )
-         NSR = MIN( NSR, ( N+6 ) / 9, IHI-ILO )
+         NSR = MIN( NSR, ( N-3 ) / 6, IHI-ILO )
          NSR = MAX( 2, NSR-MOD( NSR, 2 ) )
 *
 *        ==== Estimate optimal workspace ====
@@ -428,7 +428,7 @@
 *        ==== NSMAX = the Largest number of simultaneous shifts
 *        .    for which there is sufficient workspace. ====
 *
-         NSMAX = MIN( ( N+6 ) / 9, 2*LWORK / 3 )
+         NSMAX = MIN( ( N-3 ) / 6, 2*LWORK / 3 )
          NSMAX = NSMAX - MOD( NSMAX, 2 )
 *
 *        ==== NDFL: an iteration count restarted at deflation. ====
@@ -568,7 +568,7 @@
 *
 *                 ==== Got NS/2 or fewer shifts? Use CLAHQR
 *                 .    on a trailing principal submatrix to
-*                 .    get more. (Since NS.LE.NSMAX.LE.(N+6)/9,
+*                 .    get more. (Since NS.LE.NSMAX.LE.(N-3)/6,
 *                 .    there is enough space below the subdiagonal
 *                 .    to fit an NS-by-NS scratch array.) ====
 *
@@ -663,7 +663,7 @@
 *              .      (NVE-by-KDU) vertical work WV arrow along
 *              .      the left-hand-edge. ====
 *
-               KDU = 3*NS - 3
+               KDU = 2*NS
                KU = N - KDU + 1
                KWH = KDU + 1
                NHO = ( N-KDU+1-4 ) - ( KDU+1 ) + 1

lapack-netlib/SRC/claqr5.f

@@ -69,10 +69,9 @@
 *>             matrix entries.
 *>        = 1: CLAQR5 accumulates reflections and uses matrix-matrix
 *>             multiply to update the far-from-diagonal matrix entries.
-*>        = 2: CLAQR5 accumulates reflections, uses matrix-matrix
+*>        = 2: Same as KACC22 = 1. This option used to enable exploiting
-*>             multiply to update the far-from-diagonal matrix entries,
+*>             the 2-by-2 structure during matrix multiplications, but
-*>             and takes advantage of 2-by-2 block structure during
+*>             this is no longer supported.
-*>             matrix multiplies.
 *> \endverbatim
 *>
 *> \param[in] N
@@ -170,14 +169,14 @@
 *>
 *> \param[out] U
 *> \verbatim
-*>          U is COMPLEX array, dimension (LDU,3*NSHFTS-3)
+*>          U is COMPLEX array, dimension (LDU,2*NSHFTS)
 *> \endverbatim
 *>
 *> \param[in] LDU
 *> \verbatim
 *>          LDU is INTEGER
 *>             LDU is the leading dimension of U just as declared in the
-*>             in the calling subroutine.  LDU >= 3*NSHFTS-3.
+*>             in the calling subroutine.  LDU >= 2*NSHFTS.
 *> \endverbatim
 *>
 *> \param[in] NV
@@ -189,7 +188,7 @@
 *>
 *> \param[out] WV
 *> \verbatim
-*>          WV is COMPLEX array, dimension (LDWV,3*NSHFTS-3)
+*>          WV is COMPLEX array, dimension (LDWV,2*NSHFTS)
 *> \endverbatim
 *>
 *> \param[in] LDWV
@@ -215,7 +214,7 @@
 *> \verbatim
 *>          LDWH is INTEGER
 *>             Leading dimension of WH just as declared in the
-*>             calling procedure.  LDWH >= 3*NSHFTS-3.
+*>             calling procedure.  LDWH >= 2*NSHFTS.
 *> \endverbatim
 *>
 *  Authors:
@@ -226,7 +225,7 @@
 *> \author Univ. of Colorado Denver
 *> \author NAG Ltd.
 *
-*> \date June 2016
+*> \date January 2021
 *
 *> \ingroup complexOTHERauxiliary
 *
@@ -235,6 +234,11 @@
 *>
 *>       Karen Braman and Ralph Byers, Department of Mathematics,
 *>       University of Kansas, USA
+*>
+*>       Lars Karlsson, Daniel Kressner, and Bruno Lang
+*>
+*>       Thijs Steel, Department of Computer science,
+*>       KU Leuven, Belgium
 *
 *> \par References:
 *  ================
@@ -244,10 +248,15 @@
 *>       Performance, SIAM Journal of Matrix Analysis, volume 23, pages
 *>       929--947, 2002.
 *>
+*>       Lars Karlsson, Daniel Kressner, and Bruno Lang, Optimally packed
+*>       chains of bulges in multishift QR algorithms.
+*>       ACM Trans. Math. Softw. 40, 2, Article 12 (February 2014).
+*>
 *  =====================================================================
       SUBROUTINE CLAQR5( WANTT, WANTZ, KACC22, N, KTOP, KBOT, NSHFTS, S,
      $                   H, LDH, ILOZ, IHIZ, Z, LDZ, V, LDV, U, LDU, NV,
      $                   WV, LDWV, NH, WH, LDWH )
+      IMPLICIT NONE
 *
 *  -- LAPACK auxiliary routine (version 3.7.1) --
 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
@@ -276,11 +285,11 @@
       COMPLEX            ALPHA, BETA, CDUM, REFSUM
       REAL               H11, H12, H21, H22, SAFMAX, SAFMIN, SCL,
      $                   SMLNUM, TST1, TST2, ULP
-      INTEGER            I2, I4, INCOL, J, J2, J4, JBOT, JCOL, JLEN,
+      INTEGER            I2, I4, INCOL, J, JBOT, JCOL, JLEN,
-     $                   JROW, JTOP, K, K1, KDU, KMS, KNZ, KRCOL, KZS,
+     $                   JROW, JTOP, K, K1, KDU, KMS, KRCOL,
-     $                   M, M22, MBOT, MEND, MSTART, MTOP, NBMPS, NDCOL,
+     $                   M, M22, MBOT, MTOP, NBMPS, NDCOL,
      $                   NS, NU
-      LOGICAL            ACCUM, BLK22, BMP22
+      LOGICAL            ACCUM, BMP22
 *     ..
 *     .. External Functions ..
       REAL               SLAMCH
@@ -334,10 +343,6 @@
 *
       ACCUM = ( KACC22.EQ.1 ) .OR. ( KACC22.EQ.2 )
 *
-*     ==== If so, exploit the 2-by-2 block structure? ====
-*
-      BLK22 = ( NS.GT.2 ) .AND. ( KACC22.EQ.2 )
-*
 *     ==== clear trash ====
 *
       IF( KTOP+2.LE.KBOT )
@@ -349,28 +354,39 @@
 *
 *     ==== KDU = width of slab ====
 *
-      KDU = 6*NBMPS - 3
+      KDU = 4*NBMPS
 *
 *     ==== Create and chase chains of NBMPS bulges ====
 *
-      DO 210 INCOL = 3*( 1-NBMPS ) + KTOP - 1, KBOT - 2, 3*NBMPS - 2
+      DO 180 INCOL = KTOP - 2*NBMPS + 1, KBOT - 2, 2*NBMPS
+*
+*        JTOP = Index from which updates from the right start.
+*
+         IF( ACCUM ) THEN
+            JTOP = MAX( KTOP, INCOL )
+         ELSE IF( WANTT ) THEN
+            JTOP = 1
+         ELSE
+            JTOP = KTOP
+         END IF
+*
          NDCOL = INCOL + KDU
          IF( ACCUM )
      $      CALL CLASET( 'ALL', KDU, KDU, ZERO, ONE, U, LDU )
 *
 *        ==== Near-the-diagonal bulge chase.  The following loop
 *        .    performs the near-the-diagonal part of a small bulge
-*        .    multi-shift QR sweep.  Each 6*NBMPS-2 column diagonal
+*        .    multi-shift QR sweep.  Each 4*NBMPS column diagonal
 *        .    chunk extends from column INCOL to column NDCOL
 *        .    (including both column INCOL and column NDCOL). The
-*        .    following loop chases a 3*NBMPS column long chain of
+*        .    following loop chases a 2*NBMPS+1 column long chain of
-*        .    NBMPS bulges 3*NBMPS-2 columns to the right.  (INCOL
+*        .    NBMPS bulges 2*NBMPS columns to the right.  (INCOL
 *        .    may be less than KTOP and and NDCOL may be greater than
 *        .    KBOT indicating phantom columns from which to chase
 *        .    bulges before they are actually introduced or to which
 *        .    to chase bulges beyond column KBOT.)  ====
 *
-         DO 140 KRCOL = INCOL, MIN( INCOL+3*NBMPS-3, KBOT-2 )
+         DO 145 KRCOL = INCOL, MIN( INCOL+2*NBMPS-1, KBOT-2 )
 *
 *           ==== Bulges number MTOP to MBOT are active double implicit
 *           .    shift bulges.  There may or may not also be small
@@ -379,24 +395,156 @@
 *           .    down the diagonal to make room.  The phantom matrix
 *           .    paradigm described above helps keep track.  ====
 *
-            MTOP = MAX( 1, ( ( KTOP-1 )-KRCOL+2 ) / 3+1 )
+            MTOP = MAX( 1, ( KTOP-KRCOL ) / 2+1 )
-            MBOT = MIN( NBMPS, ( KBOT-KRCOL ) / 3 )
+            MBOT = MIN( NBMPS, ( KBOT-KRCOL-1 ) / 2 )
             M22 = MBOT + 1
-            BMP22 = ( MBOT.LT.NBMPS ) .AND. ( KRCOL+3*( M22-1 ) ).EQ.
+            BMP22 = ( MBOT.LT.NBMPS ) .AND. ( KRCOL+2*( M22-1 ) ).EQ.
      $              ( KBOT-2 )
 *
 *           ==== Generate reflections to chase the chain right
 *           .    one column.  (The minimum value of K is KTOP-1.) ====
 *
-            DO 10 M = MTOP, MBOT
+            IF ( BMP22 ) THEN
-               K = KRCOL + 3*( M-1 )
+*
+*              ==== Special case: 2-by-2 reflection at bottom treated
+*              .    separately ====
+*
+               K = KRCOL + 2*( M22-1 )
+               IF( K.EQ.KTOP-1 ) THEN
+                  CALL CLAQR1( 2, H( K+1, K+1 ), LDH, S( 2*M22-1 ),
+     $                         S( 2*M22 ), V( 1, M22 ) )
+                  BETA = V( 1, M22 )
+                  CALL CLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) )
+               ELSE
+                  BETA = H( K+1, K )
+                  V( 2, M22 ) = H( K+2, K )
+                  CALL CLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) )
+                  H( K+1, K ) = BETA
+                  H( K+2, K ) = ZERO
+               END IF
+
+*
+*              ==== Perform update from right within 
+*              .    computational window. ====
+*
+               DO 30 J = JTOP, MIN( KBOT, K+3 )
+                  REFSUM = V( 1, M22 )*( H( J, K+1 )+V( 2, M22 )*
+     $                     H( J, K+2 ) )
+                  H( J, K+1 ) = H( J, K+1 ) - REFSUM
+                  H( J, K+2 ) = H( J, K+2 ) -
+     $                          REFSUM*CONJG( V( 2, M22 ) )
+   30          CONTINUE
+*
+*              ==== Perform update from left within 
+*              .    computational window. ====
+*
+               IF( ACCUM ) THEN
+                  JBOT = MIN( NDCOL, KBOT )
+               ELSE IF( WANTT ) THEN
+                  JBOT = N
+               ELSE
+                  JBOT = KBOT
+               END IF
+               DO 40 J = K+1, JBOT
+                  REFSUM = CONJG( V( 1, M22 ) )*
+     $                     ( H( K+1, J )+CONJG( V( 2, M22 ) )*
+     $                     H( K+2, J ) )
+                  H( K+1, J ) = H( K+1, J ) - REFSUM
+                  H( K+2, J ) = H( K+2, J ) - REFSUM*V( 2, M22 )
+   40          CONTINUE
+*
+*              ==== The following convergence test requires that
+*              .    the tradition small-compared-to-nearby-diagonals
+*              .    criterion and the Ahues & Tisseur (LAWN 122, 1997)
+*              .    criteria both be satisfied.  The latter improves
+*              .    accuracy in some examples. Falling back on an
+*              .    alternate convergence criterion when TST1 or TST2
+*              .    is zero (as done here) is traditional but probably
+*              .    unnecessary. ====
+*
+               IF( K.GE.KTOP) THEN
+                  IF( H( K+1, K ).NE.ZERO ) THEN
+                     TST1 = CABS1( H( K, K ) ) + CABS1( H( K+1, K+1 ) )
+                     IF( TST1.EQ.RZERO ) THEN
+                        IF( K.GE.KTOP+1 )
+     $                     TST1 = TST1 + CABS1( H( K, K-1 ) )
+                        IF( K.GE.KTOP+2 )
+     $                     TST1 = TST1 + CABS1( H( K, K-2 ) )
+                        IF( K.GE.KTOP+3 )
+     $                     TST1 = TST1 + CABS1( H( K, K-3 ) )
+                        IF( K.LE.KBOT-2 )
+     $                     TST1 = TST1 + CABS1( H( K+2, K+1 ) )
+                        IF( K.LE.KBOT-3 )
+     $                     TST1 = TST1 + CABS1( H( K+3, K+1 ) )
+                        IF( K.LE.KBOT-4 )
+     $                     TST1 = TST1 + CABS1( H( K+4, K+1 ) )
+                     END IF
+                     IF( CABS1( H( K+1, K ) )
+     $                   .LE.MAX( SMLNUM, ULP*TST1 ) ) THEN
+                        H12 = MAX( CABS1( H( K+1, K ) ),
+     $                        CABS1( H( K, K+1 ) ) )
+                        H21 = MIN( CABS1( H( K+1, K ) ),
+     $                        CABS1( H( K, K+1 ) ) )
+                        H11 = MAX( CABS1( H( K+1, K+1 ) ),
+     $                        CABS1( H( K, K )-H( K+1, K+1 ) ) )
+                        H22 = MIN( CABS1( H( K+1, K+1 ) ),
+     $                        CABS1( H( K, K )-H( K+1, K+1 ) ) )
+                        SCL = H11 + H12
+                        TST2 = H22*( H11 / SCL )
+*
+                        IF( TST2.EQ.RZERO .OR. H21*( H12 / SCL ).LE.
+     $                      MAX( SMLNUM, ULP*TST2 ) )H( K+1, K ) = ZERO
+                     END IF
+                  END IF
+               END IF
+*
+*              ==== Accumulate orthogonal transformations. ====
+*
+               IF( ACCUM ) THEN
+                  KMS = K - INCOL
+                  DO 50 J = MAX( 1, KTOP-INCOL ), KDU
+                     REFSUM = V( 1, M22 )*( U( J, KMS+1 )+
+     $                        V( 2, M22 )*U( J, KMS+2 ) )
+                     U( J, KMS+1 ) = U( J, KMS+1 ) - REFSUM
+                     U( J, KMS+2 ) = U( J, KMS+2 ) -
+     $                               REFSUM*CONJG( V( 2, M22 ) )
+  50                 CONTINUE
+               ELSE IF( WANTZ ) THEN
+                  DO 60 J = ILOZ, IHIZ
+                     REFSUM = V( 1, M22 )*( Z( J, K+1 )+V( 2, M22 )*
+     $                        Z( J, K+2 ) )
+                     Z( J, K+1 ) = Z( J, K+1 ) - REFSUM
+                     Z( J, K+2 ) = Z( J, K+2 ) -
+     $                             REFSUM*CONJG( V( 2, M22 ) )
+  60              CONTINUE
+               END IF
+            END IF
+*
+*           ==== Normal case: Chain of 3-by-3 reflections ====
+*
+            DO 80 M = MBOT, MTOP, -1
+               K = KRCOL + 2*( M-1 )
                IF( K.EQ.KTOP-1 ) THEN
                   CALL CLAQR1( 3, H( KTOP, KTOP ), LDH, S( 2*M-1 ),
      $                         S( 2*M ), V( 1, M ) )
                   ALPHA = V( 1, M )
                   CALL CLARFG( 3, ALPHA, V( 2, M ), 1, V( 1, M ) )
                ELSE
-                  BETA = H( K+1, K )
+*
+*                 ==== Perform delayed transformation of row below
+*                 .    Mth bulge. Exploit fact that first two elements
+*                 .    of row are actually zero. ====
+*
+                  REFSUM = V( 1, M )*V( 3, M )*H( K+3, K+2 )
+                  H( K+3, K   ) = -REFSUM
+                  H( K+3, K+1 ) = -REFSUM*CONJG( V( 2, M ) )
+                  H( K+3, K+2 ) = H( K+3, K+2 ) -
+     $                            REFSUM*CONJG( V( 3, M ) )
+*
+*                 ==== Calculate reflection to move
+*                 .    Mth bulge one step. ====
+*
+                  BETA      = H( K+1, K )
                   V( 2, M ) = H( K+2, K )
                   V( 3, M ) = H( K+3, K )
                   CALL CLARFG( 3, BETA, V( 2, M ), 1, V( 1, M ) )
@@ -444,7 +592,7 @@
                         H( K+3, K ) = ZERO
                      ELSE
 *
-*                       ==== Stating a new bulge here would
+*                       ==== Starting a new bulge here would
 *                       .    create only negligible fill.
 *                       .    Replace the old reflector with
 *                       .    the new one. ====
@@ -458,163 +606,32 @@
                      END IF
                   END IF
                END IF
-   10       CONTINUE
 *
-*           ==== Generate a 2-by-2 reflection, if needed. ====
+*              ====  Apply reflection from the right and
+*              .     the first column of update from the left.
+*              .     These updates are required for the vigilant
+*              .     deflation check. We still delay most of the
+*              .     updates from the left for efficiency. ====
 *
-            K = KRCOL + 3*( M22-1 )
+               DO 70 J = JTOP, MIN( KBOT, K+3 )
-            IF( BMP22 ) THEN
+                  REFSUM = V( 1, M )*( H( J, K+1 )+V( 2, M )*
-               IF( K.EQ.KTOP-1 ) THEN
+     $                     H( J, K+2 )+V( 3, M )*H( J, K+3 ) )
-                  CALL CLAQR1( 2, H( K+1, K+1 ), LDH, S( 2*M22-1 ),
+                  H( J, K+1 ) = H( J, K+1 ) - REFSUM
-     $                         S( 2*M22 ), V( 1, M22 ) )
+                  H( J, K+2 ) = H( J, K+2 ) -
-                  BETA = V( 1, M22 )
+     $                          REFSUM*CONJG( V( 2, M ) )
-                  CALL CLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) )
+                  H( J, K+3 ) = H( J, K+3 ) -
-               ELSE
+     $                          REFSUM*CONJG( V( 3, M ) )
-                  BETA = H( K+1, K )
+   70          CONTINUE
-                  V( 2, M22 ) = H( K+2, K )
-                  CALL CLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) )
-                  H( K+1, K ) = BETA
-                  H( K+2, K ) = ZERO
-               END IF
-            END IF
 *
-*           ==== Multiply H by reflections from the left ====
+*              ==== Perform update from left for subsequent
+*              .    column. ====
 *
-            IF( ACCUM ) THEN
+               REFSUM =  CONJG( V( 1, M ) )*( H( K+1, K+1 )
-               JBOT = MIN( NDCOL, KBOT )
+     $                  +CONJG( V( 2, M ) )*H( K+2, K+1 )
-            ELSE IF( WANTT ) THEN
+     $                  +CONJG( V( 3, M ) )*H( K+3, K+1 ) )
-               JBOT = N
+               H( K+1, K+1 ) = H( K+1, K+1 ) - REFSUM
-            ELSE
+               H( K+2, K+1 ) = H( K+2, K+1 ) - REFSUM*V( 2, M )
-               JBOT = KBOT
+               H( K+3, K+1 ) = H( K+3, K+1 ) - REFSUM*V( 3, M )
-            END IF
-            DO 30 J = MAX( KTOP, KRCOL ), JBOT
-               MEND = MIN( MBOT, ( J-KRCOL+2 ) / 3 )
-               DO 20 M = MTOP, MEND
-                  K = KRCOL + 3*( M-1 )
-                  REFSUM = CONJG( V( 1, M ) )*
-     $                     ( H( K+1, J )+CONJG( V( 2, M ) )*H( K+2, J )+
-     $                     CONJG( V( 3, M ) )*H( K+3, J ) )
-                  H( K+1, J ) = H( K+1, J ) - REFSUM
-                  H( K+2, J ) = H( K+2, J ) - REFSUM*V( 2, M )
-                  H( K+3, J ) = H( K+3, J ) - REFSUM*V( 3, M )
-   20          CONTINUE
-   30       CONTINUE
-            IF( BMP22 ) THEN
-               K = KRCOL + 3*( M22-1 )
-               DO 40 J = MAX( K+1, KTOP ), JBOT
-                  REFSUM = CONJG( V( 1, M22 ) )*
-     $                     ( H( K+1, J )+CONJG( V( 2, M22 ) )*
-     $                     H( K+2, J ) )
-                  H( K+1, J ) = H( K+1, J ) - REFSUM
-                  H( K+2, J ) = H( K+2, J ) - REFSUM*V( 2, M22 )
-   40          CONTINUE
-            END IF
-*
-*           ==== Multiply H by reflections from the right.
-*           .    Delay filling in the last row until the
-*           .    vigilant deflation check is complete. ====
-*
-            IF( ACCUM ) THEN
-               JTOP = MAX( KTOP, INCOL )
-            ELSE IF( WANTT ) THEN
-               JTOP = 1
-            ELSE
-               JTOP = KTOP
-            END IF
-            DO 80 M = MTOP, MBOT
-               IF( V( 1, M ).NE.ZERO ) THEN
-                  K = KRCOL + 3*( M-1 )
-                  DO 50 J = JTOP, MIN( KBOT, K+3 )
-                     REFSUM = V( 1, M )*( H( J, K+1 )+V( 2, M )*
-     $                        H( J, K+2 )+V( 3, M )*H( J, K+3 ) )
-                     H( J, K+1 ) = H( J, K+1 ) - REFSUM
-                     H( J, K+2 ) = H( J, K+2 ) -
-     $                             REFSUM*CONJG( V( 2, M ) )
-                     H( J, K+3 ) = H( J, K+3 ) -
-     $                             REFSUM*CONJG( V( 3, M ) )
-   50             CONTINUE
-*
-                  IF( ACCUM ) THEN
-*
-*                    ==== Accumulate U. (If necessary, update Z later
-*                    .    with with an efficient matrix-matrix
-*                    .    multiply.) ====
-*
-                     KMS = K - INCOL
-                     DO 60 J = MAX( 1, KTOP-INCOL ), KDU
-                        REFSUM = V( 1, M )*( U( J, KMS+1 )+V( 2, M )*
-     $                           U( J, KMS+2 )+V( 3, M )*U( J, KMS+3 ) )
-                        U( J, KMS+1 ) = U( J, KMS+1 ) - REFSUM
-                        U( J, KMS+2 ) = U( J, KMS+2 ) -
-     $                                  REFSUM*CONJG( V( 2, M ) )
-                        U( J, KMS+3 ) = U( J, KMS+3 ) -
-     $                                  REFSUM*CONJG( V( 3, M ) )
-   60                CONTINUE
-                  ELSE IF( WANTZ ) THEN
-*
-*                    ==== U is not accumulated, so update Z
-*                    .    now by multiplying by reflections
-*                    .    from the right. ====
-*
-                     DO 70 J = ILOZ, IHIZ
-                        REFSUM = V( 1, M )*( Z( J, K+1 )+V( 2, M )*
-     $                           Z( J, K+2 )+V( 3, M )*Z( J, K+3 ) )
-                        Z( J, K+1 ) = Z( J, K+1 ) - REFSUM
-                        Z( J, K+2 ) = Z( J, K+2 ) -
-     $                                REFSUM*CONJG( V( 2, M ) )
-                        Z( J, K+3 ) = Z( J, K+3 ) -
-     $                                REFSUM*CONJG( V( 3, M ) )
-   70                CONTINUE
-                  END IF
-               END IF
-   80       CONTINUE
-*
-*           ==== Special case: 2-by-2 reflection (if needed) ====
-*
-            K = KRCOL + 3*( M22-1 )
-            IF( BMP22 ) THEN
-               IF ( V( 1, M22 ).NE.ZERO ) THEN
-                  DO 90 J = JTOP, MIN( KBOT, K+3 )
-                     REFSUM = V( 1, M22 )*( H( J, K+1 )+V( 2, M22 )*
-     $                        H( J, K+2 ) )
-                     H( J, K+1 ) = H( J, K+1 ) - REFSUM
-                     H( J, K+2 ) = H( J, K+2 ) -
-     $                             REFSUM*CONJG( V( 2, M22 ) )
-   90             CONTINUE
-*
-                  IF( ACCUM ) THEN
-                     KMS = K - INCOL
-                     DO 100 J = MAX( 1, KTOP-INCOL ), KDU
-                        REFSUM = V( 1, M22 )*( U( J, KMS+1 )+
-     $                           V( 2, M22 )*U( J, KMS+2 ) )
-                        U( J, KMS+1 ) = U( J, KMS+1 ) - REFSUM
-                        U( J, KMS+2 ) = U( J, KMS+2 ) -
-     $                                  REFSUM*CONJG( V( 2, M22 ) )
-  100                CONTINUE
-                  ELSE IF( WANTZ ) THEN
-                     DO 110 J = ILOZ, IHIZ
-                        REFSUM = V( 1, M22 )*( Z( J, K+1 )+V( 2, M22 )*
-     $                           Z( J, K+2 ) )
-                        Z( J, K+1 ) = Z( J, K+1 ) - REFSUM
-                        Z( J, K+2 ) = Z( J, K+2 ) -
-     $                                REFSUM*CONJG( V( 2, M22 ) )
-  110                CONTINUE
-                  END IF
-               END IF
-            END IF
-*
-*           ==== Vigilant deflation check ====
-*
-            MSTART = MTOP
-            IF( KRCOL+3*( MSTART-1 ).LT.KTOP )
-     $         MSTART = MSTART + 1
-            MEND = MBOT
-            IF( BMP22 )
-     $         MEND = MEND + 1
-            IF( KRCOL.EQ.KBOT-2 )
-     $         MEND = MEND + 1
-            DO 120 M = MSTART, MEND
-               K = MIN( KBOT-1, KRCOL+3*( M-1 ) )
 *
 *              ==== The following convergence test requires that
 *              .    the tradition small-compared-to-nearby-diagonals
@@ -625,6 +642,8 @@
 *              .    is zero (as done here) is traditional but probably
 *              .    unnecessary. ====
 *
+               IF( K.LT.KTOP)
+     $              CYCLE
                IF( H( K+1, K ).NE.ZERO ) THEN
                   TST1 = CABS1( H( K, K ) ) + CABS1( H( K+1, K+1 ) )
                   IF( TST1.EQ.RZERO ) THEN
@@ -658,22 +677,77 @@
      $                   MAX( SMLNUM, ULP*TST2 ) )H( K+1, K ) = ZERO
                   END IF
                END IF
-  120       CONTINUE
+   80       CONTINUE
 *
-*           ==== Fill in the last row of each bulge. ====
+*           ==== Multiply H by reflections from the left ====
 *
-            MEND = MIN( NBMPS, ( KBOT-KRCOL-1 ) / 3 )
+            IF( ACCUM ) THEN
-            DO 130 M = MTOP, MEND
+               JBOT = MIN( NDCOL, KBOT )
-               K = KRCOL + 3*( M-1 )
+            ELSE IF( WANTT ) THEN
-               REFSUM = V( 1, M )*V( 3, M )*H( K+4, K+3 )
+               JBOT = N
-               H( K+4, K+1 ) = -REFSUM
+            ELSE
-               H( K+4, K+2 ) = -REFSUM*CONJG( V( 2, M ) )
+               JBOT = KBOT
-               H( K+4, K+3 ) = H( K+4, K+3 ) - REFSUM*CONJG( V( 3, M ) )
+            END IF
-  130       CONTINUE
+*
+            DO 100 M = MBOT, MTOP, -1
+               K = KRCOL + 2*( M-1 )
+               DO 90 J = MAX( KTOP, KRCOL + 2*M ), JBOT
+                  REFSUM = CONJG( V( 1, M ) )*
+     $                     ( H( K+1, J )+CONJG( V( 2, M ) )*
+     $                     H( K+2, J )+CONJG( V( 3, M ) )*H( K+3, J ) )
+                  H( K+1, J ) = H( K+1, J ) - REFSUM
+                  H( K+2, J ) = H( K+2, J ) - REFSUM*V( 2, M )
+                  H( K+3, J ) = H( K+3, J ) - REFSUM*V( 3, M )
+   90          CONTINUE
+  100       CONTINUE
+*
+*           ==== Accumulate orthogonal transformations. ====
+*
+            IF( ACCUM ) THEN
+*
+*              ==== Accumulate U. (If needed, update Z later
+*              .    with an efficient matrix-matrix
+*              .    multiply.) ====
+*
+               DO 120 M = MBOT, MTOP, -1
+                  K = KRCOL + 2*( M-1 )
+                  KMS = K - INCOL
+                  I2 = MAX( 1, KTOP-INCOL )
+                  I2 = MAX( I2, KMS-(KRCOL-INCOL)+1 )
+                  I4 = MIN( KDU, KRCOL + 2*( MBOT-1 ) - INCOL + 5 )
+                  DO 110 J = I2, I4
+                     REFSUM = V( 1, M )*( U( J, KMS+1 )+V( 2, M )*
+     $                        U( J, KMS+2 )+V( 3, M )*U( J, KMS+3 ) )
+                     U( J, KMS+1 ) = U( J, KMS+1 ) - REFSUM
+                     U( J, KMS+2 ) = U( J, KMS+2 ) -
+     $                               REFSUM*CONJG( V( 2, M ) )
+                     U( J, KMS+3 ) = U( J, KMS+3 ) -
+     $                               REFSUM*CONJG( V( 3, M ) )
+  110             CONTINUE
+  120          CONTINUE
+            ELSE IF( WANTZ ) THEN
+*
+*              ==== U is not accumulated, so update Z
+*              .    now by multiplying by reflections
+*              .    from the right. ====
+*
+               DO 140 M = MBOT, MTOP, -1
+                  K = KRCOL + 2*( M-1 )
+                  DO 130 J = ILOZ, IHIZ
+                     REFSUM = V( 1, M )*( Z( J, K+1 )+V( 2, M )*
+     $                        Z( J, K+2 )+V( 3, M )*Z( J, K+3 ) )
+                     Z( J, K+1 ) = Z( J, K+1 ) - REFSUM
+                     Z( J, K+2 ) = Z( J, K+2 ) -
+     $                             REFSUM*CONJG( V( 2, M ) )
+                     Z( J, K+3 ) = Z( J, K+3 ) -
+     $                             REFSUM*CONJG( V( 3, M ) )
+  130             CONTINUE
+  140          CONTINUE
+            END IF
 *
 *           ==== End of near-the-diagonal bulge chase. ====
 *
-  140    CONTINUE
+  145    CONTINUE
 *
 *        ==== Use U (if accumulated) to update far-from-diagonal
 *        .    entries in H.  If required, use U to update Z as
@@ -687,220 +761,45 @@
                JTOP = KTOP
                JBOT = KBOT
             END IF
-            IF( ( .NOT.BLK22 ) .OR. ( INCOL.LT.KTOP ) .OR.
+            K1 = MAX( 1, KTOP-INCOL )
-     $          ( NDCOL.GT.KBOT ) .OR. ( NS.LE.2 ) ) THEN
+            NU = ( KDU-MAX( 0, NDCOL-KBOT ) ) - K1 + 1
 *
-*              ==== Updates not exploiting the 2-by-2 block
+*           ==== Horizontal Multiply ====
-*              .    structure of U.  K1 and NU keep track of
-*              .    the location and size of U in the special
-*              .    cases of introducing bulges and chasing
-*              .    bulges off the bottom.  In these special
-*              .    cases and in case the number of shifts
-*              .    is NS = 2, there is no 2-by-2 block
-*              .    structure to exploit.  ====
 *
-               K1 = MAX( 1, KTOP-INCOL )
+            DO 150 JCOL = MIN( NDCOL, KBOT ) + 1, JBOT, NH
-               NU = ( KDU-MAX( 0, NDCOL-KBOT ) ) - K1 + 1
+               JLEN = MIN( NH, JBOT-JCOL+1 )
+               CALL CGEMM( 'C', 'N', NU, JLEN, NU, ONE, U( K1, K1 ),
+     $                     LDU, H( INCOL+K1, JCOL ), LDH, ZERO, WH,
+     $                     LDWH )
+               CALL CLACPY( 'ALL', NU, JLEN, WH, LDWH,
+     $                      H( INCOL+K1, JCOL ), LDH )
+  150       CONTINUE
 *
-*              ==== Horizontal Multiply ====
+*           ==== Vertical multiply ====
 *
-               DO 150 JCOL = MIN( NDCOL, KBOT ) + 1, JBOT, NH
+            DO 160 JROW = JTOP, MAX( KTOP, INCOL ) - 1, NV
-                  JLEN = MIN( NH, JBOT-JCOL+1 )
+               JLEN = MIN( NV, MAX( KTOP, INCOL )-JROW )
-                  CALL CGEMM( 'C', 'N', NU, JLEN, NU, ONE, U( K1, K1 ),
+               CALL CGEMM( 'N', 'N', JLEN, NU, NU, ONE,
-     $                        LDU, H( INCOL+K1, JCOL ), LDH, ZERO, WH,
+     $                     H( JROW, INCOL+K1 ), LDH, U( K1, K1 ),
-     $                        LDWH )
+     $                     LDU, ZERO, WV, LDWV )
-                  CALL CLACPY( 'ALL', NU, JLEN, WH, LDWH,
+               CALL CLACPY( 'ALL', JLEN, NU, WV, LDWV,
-     $                         H( INCOL+K1, JCOL ), LDH )
+     $                      H( JROW, INCOL+K1 ), LDH )
-  150          CONTINUE
+  160       CONTINUE
 *
-*              ==== Vertical multiply ====
+*           ==== Z multiply (also vertical) ====
 *
-               DO 160 JROW = JTOP, MAX( KTOP, INCOL ) - 1, NV
+            IF( WANTZ ) THEN
-                  JLEN = MIN( NV, MAX( KTOP, INCOL )-JROW )
+               DO 170 JROW = ILOZ, IHIZ, NV
+                  JLEN = MIN( NV, IHIZ-JROW+1 )
                   CALL CGEMM( 'N', 'N', JLEN, NU, NU, ONE,
-     $                        H( JROW, INCOL+K1 ), LDH, U( K1, K1 ),
+     $                        Z( JROW, INCOL+K1 ), LDZ, U( K1, K1 ),
      $                        LDU, ZERO, WV, LDWV )
                   CALL CLACPY( 'ALL', JLEN, NU, WV, LDWV,
-     $                         H( JROW, INCOL+K1 ), LDH )
+     $                         Z( JROW, INCOL+K1 ), LDZ )
-  160          CONTINUE
+  170          CONTINUE
-*
-*              ==== Z multiply (also vertical) ====
-*
-               IF( WANTZ ) THEN
-                  DO 170 JROW = ILOZ, IHIZ, NV
-                     JLEN = MIN( NV, IHIZ-JROW+1 )
-                     CALL CGEMM( 'N', 'N', JLEN, NU, NU, ONE,
-     $                           Z( JROW, INCOL+K1 ), LDZ, U( K1, K1 ),
-     $                           LDU, ZERO, WV, LDWV )
-                     CALL CLACPY( 'ALL', JLEN, NU, WV, LDWV,
-     $                            Z( JROW, INCOL+K1 ), LDZ )
-  170             CONTINUE
-               END IF
-            ELSE
-*
-*              ==== Updates exploiting U's 2-by-2 block structure.
-*              .    (I2, I4, J2, J4 are the last rows and columns
-*              .    of the blocks.) ====
-*
-               I2 = ( KDU+1 ) / 2
-               I4 = KDU
-               J2 = I4 - I2
-               J4 = KDU
-*
-*              ==== KZS and KNZ deal with the band of zeros
-*              .    along the diagonal of one of the triangular
-*              .    blocks. ====
-*
-               KZS = ( J4-J2 ) - ( NS+1 )
-               KNZ = NS + 1
-*
-*              ==== Horizontal multiply ====
-*
-               DO 180 JCOL = MIN( NDCOL, KBOT ) + 1, JBOT, NH
-                  JLEN = MIN( NH, JBOT-JCOL+1 )
-*
-*                 ==== Copy bottom of H to top+KZS of scratch ====
-*                  (The first KZS rows get multiplied by zero.) ====
-*
-                  CALL CLACPY( 'ALL', KNZ, JLEN, H( INCOL+1+J2, JCOL ),
-     $                         LDH, WH( KZS+1, 1 ), LDWH )
-*
-*                 ==== Multiply by U21**H ====
-*
-                  CALL CLASET( 'ALL', KZS, JLEN, ZERO, ZERO, WH, LDWH )
-                  CALL CTRMM( 'L', 'U', 'C', 'N', KNZ, JLEN, ONE,
-     $                        U( J2+1, 1+KZS ), LDU, WH( KZS+1, 1 ),
-     $                        LDWH )
-*
-*                 ==== Multiply top of H by U11**H ====
-*
-                  CALL CGEMM( 'C', 'N', I2, JLEN, J2, ONE, U, LDU,
-     $                        H( INCOL+1, JCOL ), LDH, ONE, WH, LDWH )
-*
-*                 ==== Copy top of H to bottom of WH ====
-*
-                  CALL CLACPY( 'ALL', J2, JLEN, H( INCOL+1, JCOL ), LDH,
-     $                         WH( I2+1, 1 ), LDWH )
-*
-*                 ==== Multiply by U21**H ====
-*
-                  CALL CTRMM( 'L', 'L', 'C', 'N', J2, JLEN, ONE,
-     $                        U( 1, I2+1 ), LDU, WH( I2+1, 1 ), LDWH )
-*
-*                 ==== Multiply by U22 ====
-*
-                  CALL CGEMM( 'C', 'N', I4-I2, JLEN, J4-J2, ONE,
-     $                        U( J2+1, I2+1 ), LDU,
-     $                        H( INCOL+1+J2, JCOL ), LDH, ONE,
-     $                        WH( I2+1, 1 ), LDWH )
-*
-*                 ==== Copy it back ====
-*
-                  CALL CLACPY( 'ALL', KDU, JLEN, WH, LDWH,
-     $                         H( INCOL+1, JCOL ), LDH )
-  180          CONTINUE
-*
-*              ==== Vertical multiply ====
-*
-               DO 190 JROW = JTOP, MAX( INCOL, KTOP ) - 1, NV
-                  JLEN = MIN( NV, MAX( INCOL, KTOP )-JROW )
-*
-*                 ==== Copy right of H to scratch (the first KZS
-*                 .    columns get multiplied by zero) ====
-*
-                  CALL CLACPY( 'ALL', JLEN, KNZ, H( JROW, INCOL+1+J2 ),
-     $                         LDH, WV( 1, 1+KZS ), LDWV )
-*
-*                 ==== Multiply by U21 ====
-*
-                  CALL CLASET( 'ALL', JLEN, KZS, ZERO, ZERO, WV, LDWV )
-                  CALL CTRMM( 'R', 'U', 'N', 'N', JLEN, KNZ, ONE,
-     $                        U( J2+1, 1+KZS ), LDU, WV( 1, 1+KZS ),
-     $                        LDWV )
-*
-*                 ==== Multiply by U11 ====
-*
-                  CALL CGEMM( 'N', 'N', JLEN, I2, J2, ONE,
-     $                        H( JROW, INCOL+1 ), LDH, U, LDU, ONE, WV,
-     $                        LDWV )
-*
-*                 ==== Copy left of H to right of scratch ====
-*
-                  CALL CLACPY( 'ALL', JLEN, J2, H( JROW, INCOL+1 ), LDH,
-     $                         WV( 1, 1+I2 ), LDWV )
-*
-*                 ==== Multiply by U21 ====
-*
-                  CALL CTRMM( 'R', 'L', 'N', 'N', JLEN, I4-I2, ONE,
-     $                        U( 1, I2+1 ), LDU, WV( 1, 1+I2 ), LDWV )
-*
-*                 ==== Multiply by U22 ====
-*
-                  CALL CGEMM( 'N', 'N', JLEN, I4-I2, J4-J2, ONE,
-     $                        H( JROW, INCOL+1+J2 ), LDH,
-     $                        U( J2+1, I2+1 ), LDU, ONE, WV( 1, 1+I2 ),
-     $                        LDWV )
-*
-*                 ==== Copy it back ====
-*
-                  CALL CLACPY( 'ALL', JLEN, KDU, WV, LDWV,
-     $                         H( JROW, INCOL+1 ), LDH )
-  190          CONTINUE
-*
-*              ==== Multiply Z (also vertical) ====
-*
-               IF( WANTZ ) THEN
-                  DO 200 JROW = ILOZ, IHIZ, NV
-                     JLEN = MIN( NV, IHIZ-JROW+1 )
-*
-*                    ==== Copy right of Z to left of scratch (first
-*                    .     KZS columns get multiplied by zero) ====
-*
-                     CALL CLACPY( 'ALL', JLEN, KNZ,
-     $                            Z( JROW, INCOL+1+J2 ), LDZ,
-     $                            WV( 1, 1+KZS ), LDWV )
-*
-*                    ==== Multiply by U12 ====
-*
-                     CALL CLASET( 'ALL', JLEN, KZS, ZERO, ZERO, WV,
-     $                            LDWV )
-                     CALL CTRMM( 'R', 'U', 'N', 'N', JLEN, KNZ, ONE,
-     $                           U( J2+1, 1+KZS ), LDU, WV( 1, 1+KZS ),
-     $                           LDWV )
-*
-*                    ==== Multiply by U11 ====
-*
-                     CALL CGEMM( 'N', 'N', JLEN, I2, J2, ONE,
-     $                           Z( JROW, INCOL+1 ), LDZ, U, LDU, ONE,
-     $                           WV, LDWV )
-*
-*                    ==== Copy left of Z to right of scratch ====
-*
-                     CALL CLACPY( 'ALL', JLEN, J2, Z( JROW, INCOL+1 ),
-     $                            LDZ, WV( 1, 1+I2 ), LDWV )
-*
-*                    ==== Multiply by U21 ====
-*
-                     CALL CTRMM( 'R', 'L', 'N', 'N', JLEN, I4-I2, ONE,
-     $                           U( 1, I2+1 ), LDU, WV( 1, 1+I2 ),
-     $                           LDWV )
-*
-*                    ==== Multiply by U22 ====
-*
-                     CALL CGEMM( 'N', 'N', JLEN, I4-I2, J4-J2, ONE,
-     $                           Z( JROW, INCOL+1+J2 ), LDZ,
-     $                           U( J2+1, I2+1 ), LDU, ONE,
-     $                           WV( 1, 1+I2 ), LDWV )
-*
-*                    ==== Copy the result back to Z ====
-*
-                     CALL CLACPY( 'ALL', JLEN, KDU, WV, LDWV,
-     $                            Z( JROW, INCOL+1 ), LDZ )
-  200             CONTINUE
-               END IF
             END IF
          END IF
-  210 CONTINUE
+  180 CONTINUE
 *
 *     ==== End of CLAQR5 ====
 *

lapack-netlib/SRC/ctgsja.f

@@ -401,7 +401,7 @@
 *     .. Parameters ..
       INTEGER            MAXIT
       PARAMETER          ( MAXIT = 40 )
-      REAL               ZERO, ONE
+      REAL               ZERO, ONE, HUGENUM
       PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0 )
       COMPLEX            CZERO, CONE
       PARAMETER          ( CZERO = ( 0.0E+0, 0.0E+0 ),
@@ -424,7 +424,8 @@
      $                   SLARTG, XERBLA
 *     ..
 *     .. Intrinsic Functions ..
-      INTRINSIC          ABS, CONJG, MAX, MIN, REAL
+      INTRINSIC          ABS, CONJG, MAX, MIN, REAL, HUGE
+      PARAMETER          ( HUGENUM = HUGE(ZERO) )
 *     ..
 *     .. Executable Statements ..
 *
@@ -610,9 +611,9 @@
 *
          A1 = REAL( A( K+I, N-L+I ) )
          B1 = REAL( B( I, N-L+I ) )
+         GAMMA = B1 / A1
 *
-         IF( A1.NE.ZERO ) THEN
+         IF( (GAMMA.LE.HUGENUM).AND.(GAMMA.GE.-HUGENUM) ) THEN
-            GAMMA = B1 / A1
 *
             IF( GAMMA.LT.ZERO ) THEN
                CALL CSSCAL( L-I+1, -ONE, B( I, N-L+I ), LDB )