/* QBLAS.H
*
* a header file for QBLAS library
*
* Last Modified 8/2009
* Original Version 12/2008
*
* (c) Elegant Mathematics, Ltd. 2006-2008
* Hanauer Muehle 2
* 66564 Ottweiler-Fuerth
* Germany
* e-mail: info@elegant-mathematics.com
* Tel: +49-163-7414473
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see .
*
* PURPOSE:
* This is BLAS (www.netlib.org/blas) implementation using quad precision
*
* Quad presicion is implemented as two double words (a,b) that a>=b*eps,
* where eps is machine presicion
*
* USAGE:
* This is complete C version of BLAS as it is implemented in ATLAS package
* (www.netlib.org/atlas), with support of complex valuers.
*
* The calling convention is similar to ATLAS, so, that
*
* cblas_d* and cblas_z* subroutins correspond to dq* and zq* subroutines.
*
*
* LICENSE:
* This package is distributed under GPL 1.0 license. In case if you want
* You should have received a copy of the GNU General Public License
* along with this program. If not, see .
*
* In case if you want to incorporate in commercial package, please,
* contact Elegant Mathematics Ltd. for commercial version of this
* package or derivatives, for example, iterative linear solvers.
*/
#ifndef CBLAS_ENUM_DEFINED_H
#define CBLAS_ENUM_DEFINED_H
// we need a CBlas/Clapack definitions for compatibility with
// ATLAS (www.netlib.org/atlas) package
enum CBLAS_ORDER {CblasRowMajor=101, CblasColMajor=102 };
enum CBLAS_TRANSPOSE {CblasNoTrans=111, CblasTrans=112, CblasConjTrans=113,
AtlasConj=114};
enum CBLAS_UPLO {CblasUpper=121, CblasLower=122};
enum CBLAS_DIAG {CblasNonUnit=131, CblasUnit=132};
enum CBLAS_SIDE {CblasLeft=141, CblasRight=142};
#endif
#ifndef _QBLAS_H_
#define _QBLAS_H_
// in case if you need to compile this package for any other precision,
// for examle, for floats, you should change the following two lines
// so, that FLOAT represents to your actual floating point type and
// _SPLITTER_ should be equal to 2^{k/2}+1, where k is the total abount
// of bits in your floating point representation.
//
#define _SPLITTER_ 134217729.0
#define FLOAT double
//
// If you are going to use floats, please, uncomment the next two lines.
//
// #define _SPLITTER_ 4097.0
// #define FLOAT float
typedef struct
{ FLOAT Hi, Lo;
} FLOATDD;
typedef struct
{ FLOATDD Re;
FLOATDD Im;
} zZAHL;
extern const zZAHL zZERO;
extern const zZAHL zONE;
extern const FLOATDD dZERO;
extern const FLOATDD dONE;
int iszero_ZZ (zZAHL A);
int isone_ZZ (zZAHL A);
zZAHL split_ZZ (zZAHL A);
zZAHL cmplx_ZZ (FLOATDD A, FLOATDD B);
FLOATDD abs_ZZ (zZAHL A);
zZAHL sq_ZZ (zZAHL A);
zZAHL sq_ZZSS (zZAHL A, zZAHL ASS);
FLOATDD aca_ZZ (zZAHL A);
FLOATDD aca_ZZSS (zZAHL A, zZAHL ASS);
zZAHL conj_ZZ (zZAHL A);
void aapc_ZZ_ZZ (zZAHL A, zZAHL *C);
void aapc_ZZSS_ZZ (zZAHL A, zZAHL ASS, zZAHL *C);
void acapc_ZZ_DD (zZAHL A, FLOATDD *C);
void acapc_ZZSS_DD (zZAHL A, zZAHL ASS, FLOATDD *C);
void abpc_ZZ_ZZ_ZZ (zZAHL A, zZAHL B, zZAHL *Y);
void abpc_ZZSS_ZZ_ZZ (zZAHL A, zZAHL ASS, zZAHL B, zZAHL *Y);
void abpc_ZZSS_ZZSS_ZZ (zZAHL A, zZAHL ASS, zZAHL B, zZAHL BSS, zZAHL *Y);
zZAHL sum_ZZ_ZZ (zZAHL A, zZAHL B);
zZAHL prod_ZZSS_ZZSS (zZAHL A, zZAHL ASS, zZAHL B, zZAHL BSS);
zZAHL prod_ZZSS_ZZ (zZAHL A, zZAHL ASS, zZAHL B);
zZAHL prod_ZZ_ZZ (zZAHL A, zZAHL B);
zZAHL div_ZZSS_ZZSS (zZAHL A, zZAHL ASS, zZAHL B, zZAHL BSS);
zZAHL div_ZZSS_ZZ (zZAHL A, zZAHL ASS, zZAHL B);
zZAHL div_ZZ_ZZSS (zZAHL A, zZAHL B, zZAHL BSS);
zZAHL div_ZZ_ZZ (zZAHL A, zZAHL B);
void qrscal (const int N, const FLOATDD Alpha, zZAHL *X, const int IncX);
// a, d --- float numbers, a>=d
// DD --- float[2] numbers, DD[0]>DD[1]*SPLITTER^2
// SS --- float[2] numbers, SS[0]>SS[1]*SPLITTER
// dSS --- d and SS arguments
//
FLOATDD split_d (FLOAT a);
FLOATDD split_DD (FLOATDD A);
FLOATDD sign_DD (FLOATDD A);
FLOATDD sq_dSS (FLOAT a, FLOATDD ASS);
FLOATDD sq_d (FLOAT a);
FLOATDD sq_DD (FLOATDD A);
FLOATDD sq_DDSS (FLOATDD A, FLOATDD ASS);
FLOATDD sqrt_DD (FLOATDD A);
FLOATDD sqrt_d (FLOAT a);
int iszero_DD (FLOATDD A);
int isone_DD (FLOATDD A);
FLOATDD abs_DD (FLOATDD A);
FLOATDD aca_DD (FLOATDD A);
FLOATDD aca_DDSS (FLOATDD A, FLOATDD ASS);
FLOATDD conj_DD (FLOATDD A);
FLOATDD sum_a_d (FLOAT a, FLOAT b);
FLOATDD sum_d_d (FLOAT a, FLOAT b);
FLOATDD sum_DD_d (FLOATDD A, FLOAT d);
FLOATDD sum_DD_DD (FLOATDD A, FLOATDD B);
FLOATDD prod_dSS_dSS (FLOAT a, FLOATDD ASS, FLOAT b, FLOATDD BSS);
FLOATDD prod_d_dSS (FLOAT a, FLOAT b, FLOATDD BSS);
FLOATDD prod_d_d (FLOAT a, FLOAT b);
FLOATDD prod_DD_d (FLOATDD A, FLOAT b);
FLOATDD prod_DD_dSS (FLOATDD A, FLOAT b, FLOATDD BSS);
FLOATDD prod_DDSS_d (FLOATDD A, FLOATDD ASS, FLOAT b);
FLOATDD prod_DDSS_dSS (FLOATDD A, FLOATDD ASS, FLOAT b, FLOATDD BSS);
FLOATDD prod_DD_DD (FLOATDD A, FLOATDD B);
FLOATDD prod_DDSS_DD (FLOATDD A, FLOATDD ASS, FLOATDD B);
FLOATDD prod_DD_DDSS (FLOATDD A, FLOATDD B, FLOATDD BSS);
FLOATDD prod_DDSS_DDSS (FLOATDD A, FLOATDD ASS, FLOATDD B, FLOATDD BSS);
FLOATDD div_DDSS_DDSS (FLOATDD A, FLOATDD ASS, FLOATDD B, FLOATDD BSS);
FLOATDD div_DDSS_DD (FLOATDD A, FLOATDD ASS, FLOATDD B);
FLOATDD div_DD_DDSS (FLOATDD A, FLOATDD B, FLOATDD BSS);
FLOATDD div_DD_DD (FLOATDD A, FLOATDD B);
FLOATDD div_DD_d (FLOATDD A, FLOAT b);
int cmp_DD_DD (FLOATDD A, FLOATDD B);
int cmp_DD_d (FLOATDD A, FLOAT b);
void abpc_DD_DD_DDD (FLOATDD A, FLOATDD B, FLOAT *C); // 88 ops
void aapc_DD_DD (FLOATDD A, FLOATDD *C);
void abpc_DD_DD_DD (FLOATDD A, FLOATDD B, FLOATDD *Y);
void abpc_DDSS_DD_DD (FLOATDD A, FLOATDD ASS, FLOATDD B, FLOATDD *Y);
void abpc_DDSS_DDSS_DD (FLOATDD A, FLOATDD ASS, FLOATDD B, FLOATDD BSS, FLOATDD *C); // 22 ops
void abpc_CCSS_CCSS_CC (FLOATDD AR, FLOATDD ARSS, FLOATDD AI, FLOATDD AISS,
FLOATDD BR, FLOATDD BRSS, FLOATDD BI, FLOATDD BISS, FLOATDD *C); // 88 ops
void aapc_DDSS_DD (FLOATDD A, FLOATDD ASS, FLOATDD *C);
void acapc_DD_DD (FLOATDD A, FLOATDD *C);
void acapc_DDSS_DD (FLOATDD A, FLOATDD ASS, FLOATDD *C);
// Real functions
void dqcopy (const int N, const FLOATDD *X, const int IncX,
FLOATDD *Y, const int IncY);
void dqswap (const int N, FLOATDD *X, const int IncX,
FLOATDD *Y, const int IncY);
FLOATDD dqdotu (const int N, const FLOATDD *X, const int IncX,
const FLOATDD *Y, const int IncY);
FLOATDD dqdot (const int N, const FLOATDD *X, const int IncX,
const FLOATDD *Y, const int IncY);
FLOATDD dqnrm2 (const int N, const FLOATDD *X, const int IncX);
FLOATDD dqasum (const int N, const FLOATDD *X, const int IncX);
int dqamax (const int N, const FLOATDD *X, const int IncX);
void dqaxpy (const int N, const FLOATDD Alpha, const FLOATDD *X, const int IncX,
FLOATDD *Y, const int IncY);
void dqset (const int N, const FLOATDD Alpha, FLOATDD *X, const int IncX);
void dqaxpby (const int N, const FLOATDD Alpha, const FLOATDD *X, const int IncX,
const FLOATDD Beta, FLOATDD *Y, const int IncY);
void dqscal (const int N, const FLOATDD Alpha, FLOATDD *X, const int IncX);
void dqgemv (const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
const FLOATDD Alpha, const FLOATDD *A, const int LdA, const FLOATDD *X,
const int IncX, const FLOATDD Beta, FLOATDD *Y, const int IncY);
void dqgbmv (const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
const int KL, const int KU, const FLOATDD Alpha, const FLOATDD *A,
const int LdA, const FLOATDD *X, const int IncX, const FLOATDD Beta,
FLOATDD *Y, const int IncY);
void dqtrmv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const FLOATDD *A,
const int LdA, FLOATDD *X, const int IncX);
void dqtbmv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const int K,
const FLOATDD *A, const int LdA, FLOATDD *X, const int IncX);
void dqtpmv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const FLOATDD *AP,
FLOATDD *X, const int IncX);
void dqtrsv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const FLOATDD *A,
const int LdA, FLOATDD *X, const int IncX);
void dqtbsv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const int K,
const FLOATDD *A, const int LdA, FLOATDD *X, const int IncX);
void dqtpsv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const FLOATDD *AP,
FLOATDD *X, const int IncX);
void dqhemv (const enum CBLAS_UPLO Uplo, const int N, FLOATDD Alpha,
const FLOATDD *A, const int LdA, const FLOATDD *X, const int IncX,
FLOATDD Beta, FLOATDD *Y, const int IncY);
void dqhbmv (const enum CBLAS_UPLO Uplo, const int N, const int K,
FLOATDD Alpha, const FLOATDD *A, const int LdA, const FLOATDD *X,
const int IncX, FLOATDD Beta, FLOATDD *Y, const int IncY);
void dqhpmv (const enum CBLAS_UPLO Uplo, const int N, FLOATDD Alpha,
const FLOATDD *AP, const FLOATDD *X, const int IncX, FLOATDD Beta,
FLOATDD *Y, const int IncY);
void dqher (const enum CBLAS_UPLO Uplo, const int N, const FLOATDD Alpha,
const FLOATDD *X, const int IncX, FLOATDD *A, const int LdA);
void dqhpr (const enum CBLAS_UPLO Uplo, const int N, const FLOATDD Alpha,
const FLOATDD *X, const int IncX, FLOATDD *A);
void dqher2 (const enum CBLAS_UPLO Uplo, const int N, FLOATDD Alpha,
const FLOATDD *X, const int IncX, const FLOATDD *Y,
const int IncY, FLOATDD *A, const int LdA);
void dqhpr2 (const enum CBLAS_UPLO Uplo, const int N, FLOATDD Alpha,
const FLOATDD *X, const int IncX, const FLOATDD *Y,
const int IncY, FLOATDD *AP);
void dqgeru (const int M, const int N, FLOATDD Alpha, const FLOATDD *X,
const int IncX, const FLOATDD *Y, const int IncY,
FLOATDD *A, const int LdA);
void dqgerc (const int M, const int N, FLOATDD Alpha, const FLOATDD *X,
const int IncX, const FLOATDD *Y, const int IncY,
FLOATDD *A, const int LdA);
void dqgemm (const enum CBLAS_TRANSPOSE TransA, const enum
CBLAS_TRANSPOSE TransB, const int M, const int N, const int K,
const FLOATDD Alpha, const FLOATDD *A, const int LdA, const FLOATDD *B,
const int LdB, const FLOATDD Beta, FLOATDD *C, const int LdC);
void dqsymm (const enum CBLAS_TRANSPOSE ConjA, const enum CBLAS_SIDE Side,
const enum CBLAS_UPLO Uplo, const int M, const int N,
const FLOATDD Alpha, const FLOATDD *A, const int LdA, const FLOATDD *B,
const int LdB, const FLOATDD Beta, FLOATDD *C, const int LdC);
void dqsyrk (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE Trans,
const int N, const int K, const FLOATDD Alpha, const FLOATDD *A,
const int LdA, const FLOATDD Beta, FLOATDD *C, const int LdC);
void dqsyr2k (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE Trans,
const int N, const int K, const FLOATDD Alpha, const FLOATDD *A,
const int LdA, const FLOATDD *B, const int LdB, const FLOATDD Beta,
FLOATDD *C, const int LdC);
void dqtrmm (const enum CBLAS_SIDE Side, const enum CBLAS_UPLO Uplo,
const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
const int M, const int N, const FLOATDD Alpha, const FLOATDD *A,
const int LdA, FLOATDD *B, const int LdB);
// Complex functions
void zqcopy (const int N, const zZAHL *X, const int IncX,
zZAHL *Y, const int IncY);
void zqswap (const int N, zZAHL *X, const int IncX,
zZAHL *Y, const int IncY);
zZAHL zqdotu (const int N, const zZAHL *X, const int IncX,
const zZAHL *Y, const int IncY);
zZAHL zqdot (const int N, const zZAHL *X, const int IncX,
const zZAHL *Y, const int IncY);
FLOATDD zqnrm2 (const int N, const zZAHL *X, const int IncX);
FLOATDD zqasum (const int N, const zZAHL *X, const int IncX);
int zqamax (const int N, const zZAHL *X, const int IncX);
void zqaxpy (const int N, const zZAHL Alpha, const zZAHL *X, const int IncX,
zZAHL *Y, const int IncY);
void zqset (const int N, const zZAHL Alpha, zZAHL *X, const int IncX);
void zqaxpby (const int N, const zZAHL Alpha, const zZAHL *X, const int IncX,
const zZAHL Beta, zZAHL *Y, const int IncY);
void zqscal (const int N, const zZAHL Alpha, zZAHL *X, const int IncX);
void zqgemv (const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
const zZAHL Alpha, const zZAHL *A, const int LdA, const zZAHL *X,
const int IncX, const zZAHL Beta, zZAHL *Y, const int IncY);
void zqgbmv (const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
const int KL, const int KU, const zZAHL Alpha, const zZAHL *A,
const int LdA, const zZAHL *X, const int IncX, const zZAHL Beta,
zZAHL *Y, const int IncY);
void zqtrmv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const zZAHL *A,
const int LdA, zZAHL *X, const int IncX);
void zqtbmv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const int K,
const zZAHL *A, const int LdA, zZAHL *X, const int IncX);
void zqtpmv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const zZAHL *AP,
zZAHL *X, const int IncX);
void zqtrsv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const zZAHL *A,
const int LdA, zZAHL *X, const int IncX);
void zqtbsv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const int K,
const zZAHL *A, const int LdA, zZAHL *X, const int IncX);
void zqtpsv (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
const enum CBLAS_DIAG Diag, const int N, const zZAHL *AP,
zZAHL *X, const int IncX);
void zqhemv (const enum CBLAS_UPLO Uplo, const int N, zZAHL Alpha,
const zZAHL *A, const int LdA, const zZAHL *X, const int IncX,
zZAHL Beta, zZAHL *Y, const int IncY);
void zqhbmv (const enum CBLAS_UPLO Uplo, const int N, const int K,
zZAHL Alpha, const zZAHL *A, const int LdA, const zZAHL *X,
const int IncX, zZAHL Beta, zZAHL *Y, const int IncY);
void zqhpmv (const enum CBLAS_UPLO Uplo, const int N, zZAHL Alpha,
const zZAHL *AP, const zZAHL *X, const int IncX, zZAHL Beta,
zZAHL *Y, const int IncY);
void zqher (const enum CBLAS_UPLO Uplo, const int N, const zZAHL Alpha,
const zZAHL *X, const int IncX, zZAHL *A, const int LdA);
void zqhpr (const enum CBLAS_UPLO Uplo, const int N, const zZAHL Alpha,
const zZAHL *X, const int IncX, zZAHL *A);
void zqher2 (const enum CBLAS_UPLO Uplo, const int N, zZAHL Alpha,
const zZAHL *X, const int IncX, const zZAHL *Y,
const int IncY, zZAHL *A, const int LdA);
void zqhpr2 (const enum CBLAS_UPLO Uplo, const int N, zZAHL Alpha,
const zZAHL *X, const int IncX, const zZAHL *Y,
const int IncY, zZAHL *AP);
void zqgeru (const int M, const int N, zZAHL Alpha, const zZAHL *X,
const int IncX, const zZAHL *Y, const int IncY,
zZAHL *A, const int LdA);
void zqgerc (const int M, const int N, zZAHL Alpha, const zZAHL *X,
const int IncX, const zZAHL *Y, const int IncY,
zZAHL *A, const int LdA);
void zqgemm (const enum CBLAS_TRANSPOSE TransA, const enum
CBLAS_TRANSPOSE TransB, const int M, const int N, const int K,
const zZAHL Alpha, const zZAHL *A, const int LdA, const zZAHL *B,
const int LdB, const zZAHL Beta, zZAHL *C, const int LdC);
void zqsymm (const enum CBLAS_TRANSPOSE ConjA, const enum CBLAS_SIDE Side,
const enum CBLAS_UPLO Uplo, const int M, const int N,
const zZAHL Alpha, const zZAHL *A, const int LdA, const zZAHL *B,
const int LdB, const zZAHL Beta, zZAHL *C, const int LdC);
void zqsyrk (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE Trans,
const int N, const int K, const zZAHL Alpha, const zZAHL *A,
const int LdA, const zZAHL Beta, zZAHL *C, const int LdC);
void zqsyr2k (const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE Trans,
const int N, const int K, const zZAHL Alpha, const zZAHL *A,
const int LdA, const zZAHL *B, const int LdB, const zZAHL Beta,
zZAHL *C, const int LdC);
void zqtrmm (const enum CBLAS_SIDE Side, const enum CBLAS_UPLO Uplo,
const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
const int M, const int N, const zZAHL Alpha, const zZAHL *A,
const int LdA, zZAHL *B, const int LdB);
//
// The following functions are not yet implemented and will be supported in the nearest future
//
// void qtrsm (const enum CBLAS_SIDE Side, const enum CBLAS_UPLO Uplo,
// const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
// const int M, const int N, const FLOATDD Alpha, const FLOATDD *A,
// const int LdA, FLOATDD *B, const int LdB);
//
// void qrotmg (FLOATDD *d1, FLOATDD *d2, FLOATDD *b1, const FLOATDD b2, FLOATDD *P);
// void qrotm (const int N, FLOATDD *X, const int IncX, FLOATDD *Y,
// const int IncY, const FLOATDD *P);
// void qrotg (FLOATDD *a, FLOATDD *b, FLOATDD *c, FLOATDD *s);
// void qrot (const int N, FLOATDD *X, const int IncX, FLOATDD *Y,
// const int IncY, const FLOATDD c, const FLOATDD s);
#endif