syrk.c

/*********************************************************************/
/* Copyright 2009, 2010 The University of Texas at Austin. */
/* All rights reserved. */
/* */
/* Redistribution and use in source and binary forms, with or */
/* without modification, are permitted provided that the following */
/* conditions are met: */
/* */
/* 1. Redistributions of source code must retain the above */
/* copyright notice, this list of conditions and the following */
/* disclaimer. */
/* */
/* 2. Redistributions in binary form must reproduce the above */
/* copyright notice, this list of conditions and the following */
/* disclaimer in the documentation and/or other materials */
/* provided with the distribution. */
/* */
/* THIS SOFTWARE IS PROVIDED BY THE UNIVERSITY OF TEXAS AT */
/* AUSTIN ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, */
/* INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF */
/* MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE */
/* DISCLAIMED. IN NO EVENT SHALL THE UNIVERSITY OF TEXAS AT */
/* AUSTIN OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, */
/* INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES */
/* (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE */
/* GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR */
/* BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF */
/* LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT */
/* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT */
/* OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE */
/* POSSIBILITY OF SUCH DAMAGE. */
/* */
/* The views and conclusions contained in the software and */
/* documentation are those of the authors and should not be */
/* interpreted as representing official policies, either expressed */
/* or implied, of The University of Texas at Austin. */
/*********************************************************************/

#include<stdio.h>
#include<ctype.h>
#include"common.h"
#ifdefFUNCTION_PROFILE
#include"functable.h"
#endif

#ifndefCOMPLEX
#ifdefXDOUBLE
#defineERROR_NAME "QSYRK "
#elif defined(DOUBLE)
#defineERROR_NAME "DSYRK "
#else
#defineERROR_NAME "SSYRK "
#endif
#else
#ifndefHEMM
#ifdefXDOUBLE
#defineERROR_NAME "XSYRK "
#elif defined(DOUBLE)
#defineERROR_NAME "ZSYRK "
#else
#defineERROR_NAME "CSYRK "
#endif
#else
#ifdefXDOUBLE
#defineERROR_NAME "XHERK "
#elif defined(DOUBLE)
#defineERROR_NAME "ZHERK "
#else
#defineERROR_NAME "CHERK "
#endif
#endif
#endif

staticint (*syrk[])(blas_arg_t*, BLASLONG*, BLASLONG*, FLOAT*, FLOAT*, BLASLONG) = {
#ifndefHEMM
SYRK_UN, SYRK_UC, SYRK_LN, SYRK_LC,
#if defined(SMP) && !defined(USE_SIMPLE_THREADED_LEVEL3)
SYRK_THREAD_UN, SYRK_THREAD_UC, SYRK_THREAD_LN, SYRK_THREAD_LC,
#endif
#else
HERK_UN, HERK_UC, HERK_LN, HERK_LC,
#if defined(SMP) && !defined(USE_SIMPLE_THREADED_LEVEL3)
HERK_THREAD_UN, HERK_THREAD_UC, HERK_THREAD_LN, HERK_THREAD_LC,
#endif
#endif
};

#ifndefCBLAS

voidNAME(char*UPLO, char*TRANS,
blasint*N, blasint*K,
FLOAT*alpha, FLOAT*a, blasint*ldA,
FLOAT*beta, FLOAT*c, blasint*ldC){

charuplo_arg=*UPLO;
chartrans_arg=*TRANS;

blas_arg_targs;

FLOAT*buffer;
FLOAT*sa, *sb;

#ifdefSMP
#ifdefUSE_SIMPLE_THREADED_LEVEL3
#ifndefCOMPLEX
#ifdefXDOUBLE
intmode=BLAS_XDOUBLE | BLAS_REAL;
#elif defined(DOUBLE)
intmode=BLAS_DOUBLE | BLAS_REAL;
#else
intmode=BLAS_SINGLE | BLAS_REAL;
#endif
#else
#ifdefXDOUBLE
intmode=BLAS_XDOUBLE | BLAS_COMPLEX;
#elif defined(DOUBLE)
intmode=BLAS_DOUBLE | BLAS_COMPLEX;
#else
intmode=BLAS_SINGLE | BLAS_COMPLEX;
#endif
#endif
#endif
#endif

blasintinfo;
intuplo;
inttrans;
intnrowa;

PRINT_DEBUG_NAME;

args.n=*N;
args.k=*K;

args.a= (void*)a;
args.c= (void*)c;

args.lda=*ldA;
args.ldc=*ldC;

args.alpha= (void*)alpha;
args.beta= (void*)beta;

TOUPPER(uplo_arg);
TOUPPER(trans_arg);

uplo=-1;
trans=-1;

if (uplo_arg=='U') uplo=0;
if (uplo_arg=='L') uplo=1;


#ifndefCOMPLEX
if (trans_arg=='N') trans=0;
if (trans_arg=='T') trans=1;
if (trans_arg=='C') trans=1;
#else
#ifdefHEMM
if (trans_arg=='N') trans=0;
if (trans_arg=='C') trans=1;
#else
if (trans_arg=='N') trans=0;
if (trans_arg=='T') trans=1;
#endif

#endif

nrowa=args.n;
if (trans&1) nrowa=args.k;

info=0;

if (args.ldc<MAX(1,args.n)) info=10;
if (args.lda<MAX(1,nrowa)) info=7;
if (args.k<0) info=4;
if (args.n<0) info=3;
if (trans<0) info=2;
if (uplo<0) info=1;

if (info!=0) {
BLASFUNC(xerbla)(ERROR_NAME, &info, sizeof(ERROR_NAME));
return;
 }

#else

voidCNAME(enumCBLAS_ORDERorder, enumCBLAS_UPLOUplo, enumCBLAS_TRANSPOSETrans,
blasintn, blasintk,
#if !defined(COMPLEX) || defined(HEMM)
FLOATalpha,
#else
void*valpha,
#endif
#if !defined(COMPLEX)
FLOAT*a, blasintlda,
#else
void*va, blasintlda,
#endif
#if !defined(COMPLEX) || defined(HEMM)
FLOATbeta,
#else
void*vbeta,
#endif
#if !defined(COMPLEX)
FLOAT*c, blasintldc) {
#else
void*vc, blasintldc) {
#endif

#ifdefCOMPLEX
#if !defined(HEMM)
FLOAT*alpha= (FLOAT*) valpha;
FLOAT*beta= (FLOAT*) vbeta;
#endif
FLOAT*a= (FLOAT*) va;
FLOAT*c= (FLOAT*) vc;
#endif

blas_arg_targs;
intuplo, trans;
blasintinfo, nrowa;

FLOAT*buffer;
FLOAT*sa, *sb;

#ifdefSMP
#ifdefUSE_SIMPLE_THREADED_LEVEL3
#ifndefCOMPLEX
#ifdefXDOUBLE
intmode=BLAS_XDOUBLE | BLAS_REAL;
#elif defined(DOUBLE)
intmode=BLAS_DOUBLE | BLAS_REAL;
#else
intmode=BLAS_SINGLE | BLAS_REAL;
#endif
#else
#ifdefXDOUBLE
intmode=BLAS_XDOUBLE | BLAS_COMPLEX;
#elif defined(DOUBLE)
intmode=BLAS_DOUBLE | BLAS_COMPLEX;
#else
intmode=BLAS_SINGLE | BLAS_COMPLEX;
#endif
#endif
#endif
#endif

PRINT_DEBUG_CNAME;

args.n=n;
args.k=k;

args.a= (void*)a;
args.c= (void*)c;

args.lda=lda;
args.ldc=ldc;

#if !defined(COMPLEX) || defined(HEMM)
args.alpha= (void*)&alpha;
args.beta= (void*)&beta;
#else
args.alpha= (void*)alpha;
args.beta= (void*)beta;
#endif

trans=-1;
uplo=-1;
info=0;

if (order==CblasColMajor) {
if (Uplo==CblasUpper) uplo=0;
if (Uplo==CblasLower) uplo=1;

if (Trans==CblasNoTrans) trans=0;
#ifndefCOMPLEX
if (Trans==CblasTrans) trans=1;
if (Trans==CblasConjNoTrans) trans=0;
if (Trans==CblasConjTrans) trans=1;
#elif !defined(HEMM)
if (Trans==CblasTrans) trans=1;
#else
if (Trans==CblasConjTrans) trans=1;
#endif

info=-1;

nrowa=args.n;
if (trans&1) nrowa=args.k;

if (args.ldc<MAX(1,args.n)) info=10;
if (args.lda<MAX(1,nrowa)) info=7;
if (args.k<0) info=4;
if (args.n<0) info=3;
if (trans<0) info=2;
if (uplo<0) info=1;
 }

if (order==CblasRowMajor) {
if (Uplo==CblasUpper) uplo=1;
if (Uplo==CblasLower) uplo=0;

if (Trans==CblasNoTrans) trans=1;
#ifndefCOMPLEX
if (Trans==CblasTrans) trans=0;
if (Trans==CblasConjNoTrans) trans=1;
if (Trans==CblasConjTrans) trans=0;
#elif !defined(HEMM)
if (Trans==CblasTrans) trans=0;
#else
if (Trans==CblasConjTrans) trans=0;
#endif

info=-1;

nrowa=args.n;
if (trans&1) nrowa=args.k;

if (args.ldc<MAX(1,args.n)) info=10;
if (args.lda<MAX(1,nrowa)) info=7;
if (args.k<0) info=4;
if (args.n<0) info=3;
if (trans<0) info=2;
if (uplo<0) info=1;
 }

if (info >= 0) {
BLASFUNC(xerbla)(ERROR_NAME, &info, sizeof(ERROR_NAME));
return;
 }

#endif

if (args.n==0) return;

IDEBUG_START;

FUNCTION_PROFILE_START();

buffer= (FLOAT*)blas_memory_alloc(0);

sa= (FLOAT*)((BLASLONG)buffer+GEMM_OFFSET_A);
sb= (FLOAT*)(((BLASLONG)sa+ ((GEMM_P*GEMM_Q*COMPSIZE*SIZE+GEMM_ALIGN) & ~GEMM_ALIGN)) +GEMM_OFFSET_B);

#ifdefSMP
#ifdefUSE_SIMPLE_THREADED_LEVEL3
if (!trans){
mode |= (BLAS_TRANSA_N | BLAS_TRANSB_T);
 } else {
mode |= (BLAS_TRANSA_T | BLAS_TRANSB_N);
 }
mode |= (uplo << BLAS_UPLO_SHIFT);
#endif

args.common=NULL;
#ifndefCOMPLEX
#ifdefDOUBLE
if (args.n<100)
#else
if (args.n<200)
#endif
#else
if (args.n<65)
#endif
args.nthreads=1;
else
args.nthreads=num_cpu_avail(3);

if (args.nthreads==1) {
#endif

 (syrk[(uplo << 1) | trans ])(&args, NULL, NULL, sa, sb, 0);

#ifdefSMP

 } else {

#ifndefUSE_SIMPLE_THREADED_LEVEL3

 (syrk[4 | (uplo << 1) | trans ])(&args, NULL, NULL, sa, sb, 0);

#else

syrk_thread(mode, &args, NULL, NULL, syrk[(uplo << 1) | trans ], sa, sb, args.nthreads);

#endif

 }
#endif

blas_memory_free(buffer);

FUNCTION_PROFILE_END(COMPSIZE*COMPSIZE, args.n*args.k+args.n*args.n / 2, args.n*args.n*args.k);

IDEBUG_END;

return;
}