saf_utility_veclib.h

Go to the documentation of this file.

/*
 * Copyright 2016-2018 Leo McCormack
 *
 * Permission to use, copy, modify, and/or distribute this software for any
 * purpose with or without fee is hereby granted, provided that the above
 * copyright notice and this permission notice appear in all copies.
 *
 * THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES WITH
 * REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY
 * AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT,
 * INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
 * LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR
 * OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR
 * PERFORMANCE OF THIS SOFTWARE.
 */
 
#ifndef SAF_VECLIB_H_INCLUDED
#define SAF_VECLIB_H_INCLUDED
 
#ifdef __cplusplus
extern "C" {
#endif /* __cplusplus */
 
#include "saf_utility_complex.h"
 
typedef enum {
  NO_CONJ = 1,  
  CONJ = 2      
}CONJ_FLAG;
 
#ifdef SAF_USE_BUILT_IN_NAIVE_CBLAS
 enum { CblasRowMajor = 101, CblasColMajor = 102 } CBLAS_ORDER;
 enum { CblasNoTrans = 111, CblasTrans = 112, CblasConjTrans = 113 } CBLAS_TRANSPOSE;
#endif
 
/* KEY:
 * ? -> s -> single floating-point precision
 * ? -> d -> double floating-point precision
 * ? -> c -> complex single floating-point precision
 * ? -> z -> complex double floating-point precision
 *
 * s -> scalar
 * v -> vector
 * m -> matrix */
 
/* ========================================================================== */
/*                     Built-in CBLAS Functions (Level 0)                     */
/* ========================================================================== */
 
#ifdef SAF_USE_BUILT_IN_NAIVE_CBLAS
void cblas_scopy(const int N, const float *X, const int incX,
                 float *Y, const int incY);
void cblas_dcopy(const int N, const double *X, const int incX,
                 double *Y, const int incY);
void cblas_ccopy(const int N, const void *X, const int incX,
                 void *Y, const int incY);
void cblas_zcopy(const int N, const void *X, const int incX,
                 void *Y, const int incY);
 
void cblas_saxpy(const int N, const float alpha, const float* X,
                 const int incX, float* Y, const int incY);
void cblas_daxpy(const int N, const double alpha, const double* X,
                 const int incX, double* Y, const int incY);
void cblas_caxpy(const int N, const void* alpha, const void* X,
                 const int incX, void* Y, const int incY);
void cblas_zaxpy(const int N, const void* alpha, const void* X,
                 const int incX, void* Y, const int incY);
#endif
 
 
/* ========================================================================== */
/*                     Built-in CBLAS Functions (Level 3)                     */
/* ========================================================================== */
 
#ifdef SAF_USE_BUILT_IN_NAIVE_CBLAS
void cblas_sgemm(const CBLAS_LAYOUT Layout, const CBLAS_TRANSPOSE TransA,
                     const CBLAS_TRANSPOSE TransB, const int M, const int N,
                     const int K, const float alpha, const float *A,
                     const int lda, const float *B, const int ldb,
                     const float beta, float *C, const int ldc);
#endif
 
 
/* ========================================================================== */
/*                     Find Index of Min-Abs-Value (?iminv)                   */
/* ========================================================================== */
 
void utility_siminv(/* Input Arguments */
                    const float* a,
                    const int len,
                    /* Output Arguments */
                    int* index);
 
void utility_ciminv(/* Input Arguments */
                    const float_complex* a,
                    const int len,
                    /* Output Arguments */
                    int* index);
 
void utility_diminv(/* Input Arguments */
                    const double* a,
                    const int len,
                    /* Output Arguments */
                    int* index);
 
void utility_ziminv(/* Input Arguments */
                    const double_complex* a,
                    const int len,
                    /* Output Arguments */
                    int* index);
 
 
/* ========================================================================== */
/*                     Find Index of Max-Abs-Value (?imaxv)                   */
/* ========================================================================== */
 
void utility_simaxv(/* Input Arguments */
                    const float* a,
                    const int len,
                    /* Output Arguments */
                    int* index);
 
void utility_cimaxv(/* Input Arguments */
                    const float_complex* a,
                    const int len,
                    /* Output Arguments */
                    int* index);
 
void utility_dimaxv(/* Input Arguments */
                    const double* a,
                    const int len,
                    /* Output Arguments */
                    int* index);
 
void utility_zimaxv(/* Input Arguments */
                    const double_complex* a,
                    const int len,
                    /* Output Arguments */
                    int* index);
 
 
/* ========================================================================== */
/*                              Vector-Abs (?vabs)                            */
/* ========================================================================== */
 
void utility_svabs(/* Input Arguments */
                   const float* a,
                   const int len,
                   /* Output Arguments */
                   float* c);
 
void utility_cvabs(/* Input Arguments */
                   const float_complex* a,
                   const int len,
                   /* Output Arguments */
                   float* c);
 
 
/* ========================================================================== */
/*                            Vector-Modulus (?vmod)                          */
/* ========================================================================== */
 
void utility_svmod(/* Input Arguments */
                   const float* a,
                   const float* b,
                   const int len,
                   /* Output Arguments */
                   float* c);
 
 
/* ========================================================================== */
/*                          Vector-Reciprocal (?vrecip)                       */
/* ========================================================================== */
 
void utility_svrecip(const float* a,
                     const int len,
                     float* c);
 
 
/* ========================================================================== */
/*                           Vector-Conjugate (?vconj)                        */
/* ========================================================================== */
 
void utility_cvconj(const float_complex* a,
                    const int len,
                    float_complex* c);
 
void utility_zvconj(const double_complex* a,
                    const int len,
                    double_complex* c);
 
 
/* ========================================================================== */
/*                        Vector-Vector Copy (?vvcopy)                        */
/* ========================================================================== */
 
void utility_svvcopy(/* Input Arguments */
                     const float* a,
                     const int len,
                     /* Output Arguments */
                     float* c);
 
void utility_cvvcopy(/* Input Arguments */
                     const float_complex* a,
                     const int len,
                     /* Output Arguments */
                     float_complex* c);
 
void utility_dvvcopy(/* Input Arguments */
                     const double* a,
                     const int len,
                     /* Output Arguments */
                     double* c);
 
void utility_zvvcopy(/* Input Arguments */
                     const double_complex* a,
                     const int len,
                     /* Output Arguments */
                     double_complex* c);
 
 
/* ========================================================================== */
/*                       Vector-Vector Addition (?vvadd)                      */
/* ========================================================================== */
 
void utility_svvadd(/* Input Arguments */
                    const float* a,
                    const float* b,
                    const int len,
                    /* Output Arguments */
                    float* c);
 
void utility_cvvadd(/* Input Arguments */
                    const float_complex* a,
                    const float_complex* b,
                    const int len,
                    /* Output Arguments */
                    float_complex* c);
 
void utility_dvvadd(/* Input Arguments */
                    const double* a,
                    const double* b,
                    const int len,
                    /* Output Arguments */
                    double* c);
 
void utility_zvvadd(/* Input Arguments */
                    const double_complex* a,
                    const double_complex* b,
                    const int len,
                    /* Output Arguments */
                    double_complex* c);
 
 
/* ========================================================================== */
/*                     Vector-Vector Subtraction (?vvsub)                     */
/* ========================================================================== */
 
void utility_svvsub(/* Input Arguments */
                    const float* a,
                    const float* b,
                    const int len,
                    /* Output Arguments */
                    float* c);
 
void utility_cvvsub(/* Input Arguments */
                    const float_complex* a,
                    const float_complex* b,
                    const int len,
                    /* Output Arguments */
                    float_complex* c);
 
void utility_dvvsub(/* Input Arguments */
                    const double* a,
                    const double* b,
                    const int len,
                    /* Output Arguments */
                    double* c);
 
void utility_zvvsub(/* Input Arguments */
                    const double_complex* a,
                    const double_complex* b,
                    const int len,
                    /* Output Arguments */
                    double_complex* c);
 
 
/* ========================================================================== */
/*                    Vector-Vector Multiplication (?vvmul)                   */
/* ========================================================================== */
 
void utility_svvmul(/* Input Arguments */
                    const float* a,
                    const float* b,
                    const int len,
                    /* Output Arguments */
                    float* c);
 
void utility_cvvmul(/* Input Arguments */
                    const float_complex* a,
                    const float_complex* b,
                    const int len,
                    /* Output Arguments */
                    float_complex* c);
 
 
/* ========================================================================== */
/*                     Vector-Vector Dot Product (?vvdot)                     */
/* ========================================================================== */
 
void utility_svvdot(/* Input Arguments */
                    const float* a,
                    const float* b,
                    const int len,
                    /* Output Arguments */
                    float* c);
 
void utility_cvvdot(/* Input Arguments */
                    const float_complex* a,
                    const float_complex* b,
                    const int len,
                    CONJ_FLAG flag,
                    /* Output Arguments */
                    float_complex* c);
 
 
/* ========================================================================== */
/*                       Vector-Scalar Product (?vsmul)                       */
/* ========================================================================== */
 
void utility_svsmul(/* Input Arguments */
                    float* a,
                    const float* s,
                    const int len,
                    /* Output Arguments */
                    float* c);
 
void utility_cvsmul(/* Input Arguments */
                    float_complex* a,
                    const float_complex* s,
                    const int len,
                    /* Output Arguments */
                    float_complex* c);
 
void utility_dvsmul(/* Input Arguments */
                    double* a,
                    const double* s,
                    const int len,
                    /* Output Arguments */
                    double* c);
 
void utility_zvsmul(/* Input Arguments */
                    double_complex* a,
                    const double_complex* s,
                    const int len,
                    /* Output Arguments */
                    double_complex* c);
 
 
/* ========================================================================== */
/*                       Vector-Scalar Division (?vsdiv)                      */
/* ========================================================================== */
 
void utility_svsdiv(/* Input Arguments */
                    const float* a,
                    const float* s,
                    const int len,
                    /* Output Arguments */
                    float* c);
 
 
/* ========================================================================== */
/*                       Vector-Scalar Addition (?vsadd)                      */
/* ========================================================================== */
 
void utility_svsadd(/* Input Arguments */
                    float* a,
                    const float* s,
                    const int len,
                    /* Output Arguments */
                    float* c);
 
 
/* ========================================================================== */
/*                     Vector-Scalar Subtraction (?vssub)                     */
/* ========================================================================== */
    
void utility_svssub(/* Input Arguments */
                    float* a,
                    const float* s,
                    const int len,
                    /* Output Arguments */
                    float* c);
 
 
/* ========================================================================== */
/*      Sparse-Vector to Compressed-Vector (Known Indices) (?sv2cv_inds)      */
/* ========================================================================== */
 
void utility_ssv2cv_inds(/* Input Arguments */
                         const float* sv,
                         const int* inds,
                         const int len,
                         /* Output Arguments */
                         float* cv);
 
void utility_csv2cv_inds(/* Input Arguments */
                         const float_complex* sv,
                         const int* inds,
                         const int len,
                         /* Output Arguments */
                         float_complex* cv);
 
void utility_dsv2cv_inds(/* Input Arguments */
                         const double* sv,
                         const int* inds,
                         const int len,
                         /* Output Arguments */
                         double* cv);
 
void utility_zsv2cv_inds(/* Input Arguments */
                         const double_complex* sv,
                         const int* inds,
                         const int len,
                         /* Output Arguments */
                         double_complex* cv);
 
 
/* ========================================================================== */
/*                     Singular-Value Decomposition (?svd)                    */
/* ========================================================================== */
 
void utility_ssvd_create(void ** const phWork, int maxDim1, int maxDim2);
 
void utility_ssvd_destroy(void ** const phWork);
 
void utility_ssvd(/* Input Arguments */
                  void* const hWork,
                  const float* A,
                  const int dim1,
                  const int dim2,
                  /* Output Arguments */
                  float* U,
                  float* S,
                  float* V,
                  float* sing);
 
void utility_csvd_create(void ** const phWork, int maxDim1, int maxDim2);
 
void utility_csvd_destroy(void ** const phWork);
 
void utility_csvd(/* Input Arguments */
                  void* const hWork,
                  const float_complex* A,
                  const int dim1,
                  const int dim2,
                  /* Output Arguments */
                  float_complex* U,
                  float_complex* S,
                  float_complex* V,
                  float* sing);
 
 
/* ========================================================================== */
/*                 Symmetric Eigenvalue Decomposition (?seig)                 */
/* ========================================================================== */
 
void utility_sseig_create(void ** const phWork, int maxDim);
 
void utility_sseig_destroy(void ** const phWork);
 
void utility_sseig(/* Input Arguments */
                   void* const hWork,
                   const float* A,
                   const int dim,
                   int sortDecFLAG,
                   /* Output Arguments */
                   float* V,
                   float* D,
                   float* eig);
 
void utility_cseig_create(void ** const phWork,
                          int maxDim);
 
void utility_cseig_destroy(void ** const phWork);
 
void utility_cseig(/* Input Arguments */
                   void* const hWork,
                   const float_complex* A,
                   const int dim,
                   int sortDecFLAG,
                   /* Output Arguments */
                   float_complex* V,
                   float_complex* D,
                   float* eig);
 
 
/* ========================================================================== */
/*                     Eigenvalues of Matrix Pair (?eigmp)                    */
/* ========================================================================== */
 
void utility_ceigmp_create(void** const phWork, int maxDim);
 
void utility_ceigmp_destroy(void ** const phWork);
 
void utility_ceigmp(/* Input Arguments */
                    void* const hWork,
                    const float_complex* A,
                    const float_complex* B,
                    const int dim,
                    /* Output Arguments */
                    float_complex* VL,
                    float_complex* VR,
                    float_complex* D);
 
void utility_zeigmp_create(void** const phWork, int maxDim);
 
void utility_zeigmp_destroy(void ** const phWork);
 
void utility_zeigmp(/* Input Arguments */
                    void* const hWork,
                    const double_complex* A,
                    const double_complex* B,
                    const int dim,
                    /* Output Arguments */
                    double_complex* VL,
                    double_complex* VR,
                    double_complex* D);
 
 
/* ========================================================================== */
/*                       Eigenvalue Decomposition (?eig)                      */
/* ========================================================================== */
 
void utility_ceig_create(void ** const phWork, int maxDim);
 
void utility_ceig_destroy(void ** const phWork);
 
void utility_ceig(/* Input Arguments */
                  void* const hWork,
                  const float_complex* A,
                  const int dim, 
                  /* Output Arguments */
                  float_complex* VL,
                  float_complex* VR,
                  float_complex* D,
                  float_complex* eig);
 
void utility_zeig_create(void ** const phWork, int maxDim);
 
void utility_zeig_destroy(void ** const phWork);
 
void utility_zeig(/* Input Arguments */
                  void* const hWork,
                  const double_complex* A,
                  const int dim,
                  /* Output Arguments */
                  double_complex* VL,
                  double_complex* VR,
                  double_complex* D,
                  double_complex* eig);
 
 
/* ========================================================================== */
/*                       General Linear Solver (?glslv)                       */
/* ========================================================================== */
 
void utility_sglslv_create(void ** const phWork,
                           int maxDim,
                           int maxNCol);
 
void utility_sglslv_destroy(void ** const phWork);
 
void utility_sglslv(/* Input Arguments */
                    void* const hWork,
                    const float* A,
                    const int dim,
                    float* B,
                    int nCol,
                    /* Output Arguments */
                    float* X);
 
void utility_cglslv_create(void ** const phWork, int maxDim, int maxNCol);
 
void utility_cglslv_destroy(void ** const phWork);
 
void utility_cglslv(/* Input Arguments */
                    void* const hWork,
                    const float_complex* A,
                    const int dim,
                    float_complex* B,
                    int nCol,
                    /* Output Arguments */
                    float_complex* X);
 
void utility_dglslv_create(void ** const phWork, int maxDim, int maxNCol);
 
void utility_dglslv_destroy(void ** const phWork);
 
void utility_dglslv(/* Input Arguments */
                    void* const hWork,
                    const double* A,
                    const int dim,
                    double* B,
                    int nCol,
                    /* Output Arguments */
                    double* X);
 
void utility_zglslv_create(void ** const phWork, int maxDim, int maxNCol);
 
void utility_zglslv_destroy(void ** const phWork);
 
void utility_zglslv(/* Input Arguments */
                    void* const hWork,
                    const double_complex* A,
                    const int dim,
                    double_complex* B,
                    int nCol,
                    /* Output Arguments */
                    double_complex* X);
 
 
/* ========================================================================== */
/*                      General Linear Solver (?glslvt)                       */
/* ========================================================================== */
 
void utility_sglslvt_create(void ** const phWork,
                            int maxDim,
                            int maxNCol);
 
void utility_sglslvt_destroy(void ** const phWork);
 
void utility_sglslvt(/* Input Arguments */
                     void* const hWork,
                     const float* A,
                     const int dim,
                     float* B,
                     int nCol,
                     /* Output Arguments */
                     float* X);
 
 
/* ========================================================================== */
/*                      Symmetric Linear Solver (?slslv)                      */
/* ========================================================================== */
 
void utility_sslslv_create(void ** const phWork, int maxDim, int maxNCol);
 
void utility_sslslv_destroy(void ** const phWork);
 
void utility_sslslv(/* Input Arguments */
                    void* const hWork,
                    const float* A,
                    const int dim,
                    float* B,
                    int nCol,
                    /* Output Arguments */
                    float* X);
 
void utility_cslslv_create(void ** const phWork, int maxDim, int maxNCol);
 
void utility_cslslv_destroy(void ** const phWork);
 
void utility_cslslv(/* Input Arguments */
                    void* const hWork,
                    const float_complex* A,
                    const int dim,
                    float_complex* B,
                    int nCol,
                    /* Output Arguments */
                    float_complex* X);
 
 
/* ========================================================================== */
/*                        Matrix Pseudo-Inverse (?pinv)                       */
/* ========================================================================== */
 
void utility_spinv_create(void ** const phWork, int maxDim1, int maxDim2);
 
void utility_spinv_destroy(void ** const phWork);
 
void utility_spinv(/* Input Arguments */
                   void* const hWork,
                   const float* A,
                   const int dim1,
                   const int dim2,
                   /* Output Arguments */
                   float* B);
 
void utility_cpinv_create(void ** const phWork, int maxDim1, int maxDim2);
 
void utility_cpinv_destroy(void ** const phWork);
 
void utility_cpinv(/* Input Arguments */
                   void* const hWork,
                   const float_complex* A,
                   const int dim1,
                   const int dim2,
                   /* Output Arguments */
                   float_complex* B);
 
void utility_dpinv_create(void ** const phWork, int maxDim1, int maxDim2);
 
void utility_dpinv_destroy(void ** const phWork);
 
void utility_dpinv(/* Input Arguments */
                   void* const hWork,
                   const double* A,
                   const int dim1,
                   const int dim2,
                   /* Output Arguments */
                   double* B);
 
void utility_zpinv_create(void ** const phWork, int maxDim1, int maxDim2);
 
void utility_zpinv_destroy(void ** const phWork);
 
void utility_zpinv(/* Input Arguments */
                   void* const hWork,
                   const double_complex* A,
                   const int dim1,
                   const int dim2,
                   /* Output Arguments */
                   double_complex* B);
 
 
/* ========================================================================== */
/*                       Cholesky Factorisation (?chol)                       */
/* ========================================================================== */
 
void utility_schol_create(void ** const phWork, int maxDim);
 
void utility_schol_destroy(void ** const phWork);
 
void utility_schol(/* Input Arguments */
                   void* const hWork,
                   const float* A,
                   const int dim,
                   /* Output Arguments */
                   float* X);
 
void utility_cchol_create(void ** const phWork, int maxDim);
 
void utility_cchol_destroy(void ** const phWork);
 
void utility_cchol(/* Input Arguments */
                   void* const hWork,
                   const float_complex* A,
                   const int dim,
                   /* Output Arguments */
                   float_complex* X);
 
 
/* ========================================================================== */
/*                        Determinant of a Matrix (?det)                      */
/* ========================================================================== */
 
void utility_sdet_create(void ** const phWork, int maxN);
 
void utility_sdet_destroy(void ** const phWork);
 
float utility_sdet(void* const hWork,
                   float* A,
                   int N);
 
void utility_ddet_create(void ** const phWork, int maxN);
 
void utility_ddet_destroy(void ** const phWork);
 
double utility_ddet(void* const hWork,
                    double* A,
                    int N);
 
    
/* ========================================================================== */
/*                           Matrix Inversion (?inv)                          */
/* ========================================================================== */
 
void utility_sinv_create(void ** const phWork, int maxDim);
 
void utility_sinv_destroy(void ** const phWork);
 
void utility_sinv(void* const hWork,
                  float* A,
                  float* B,
                  const int dim);
 
void utility_dinv_create(void ** const phWork, int maxDim);
 
void utility_dinv_destroy(void ** const phWork);
 
void utility_dinv(void* const hWork,
                  double* A,
                  double* B,
                  const int dim);
 
void utility_cinv_create(void ** const phWork, int maxDim);
 
void utility_cinv_destroy(void ** const phWork);
 
void utility_cinv(void* const hWork,
                  float_complex* A,
                  float_complex* B,
                  const int dim);
 
 
#ifdef __cplusplus
}/* extern "C" */
#endif  /* __cplusplus */
 
#endif /* SAF_VECLIB_H_INCLUDED */
 
 /* doxygen addtogroup Utilities */