lapack.h

Go to the documentation of this file.

#ifndef LAPACK_H
#define LAPACK_H

#include <complex>
#include <vector>

#include "matpack.h"

namespace lapack_help {
  
  template <typename T> struct Inverse {
    mutable int n;
    mutable int lwork;
    mutable std::vector<T> work;
    mutable std::vector<int> ipiv;
    
    explicit Inverse(Index N) noexcept : n(int(N)), lwork(int(N)), work(N), ipiv(N) {}
    
    bool resize_if_smaller(Index N) const {
      if (N>n) {
        n=int(N);
        lwork=int(N);
        work.resize(N);
        ipiv.resize(N);
        return true;
      } else {
        return false;
      }
    }
    
    int* size() const noexcept {return &n;}
    int* lsize() const noexcept {return &lwork;}
    T* workdata() const noexcept {return work.data();}
    int* ipivdata() const noexcept {return ipiv.data();}
  };
};

namespace lapack {


extern "C" void dgetrf_(
    int *m, int *n, double *A, int *lda, int *ipiv, int *info);

extern "C" void zgetrf_(
    int *m, int *n, std::complex<double> *A, int *lda, int *ipiv, int *info);


extern "C" void dgetrs_(char *trans,
                        int *n,
                        int *nrhs,
                        double *A,
                        int *lda,
                        int *ipiv,
                        double *b,
                        int *ldb,
                        int *info);

//

extern "C" void dgetri_(int *n,
                        double *A,
                        int *lda,
                        int *ipiv,
                        double *work,
                        int *lwork,
                        int *info);

extern "C" void zgetri_(int *n,
                        std::complex<double> *A,
                        int *lda,
                        int *ipiv,
                        std::complex<double> *work,
                        int *lwork,
                        int *info);


extern "C" void ilaenv_(
    int *ispec, char *name, char *opts, int *n1, int *n2, int *n3, int *n4);


extern "C" void dgesvx_(char *fact,
                        char *trans,
                        int *n,
                        int *nrhs,
                        double *A,
                        int *lda,
                        double *AF,
                        int *ldaf,
                        int *ipiv,
                        char *equed,
                        double *R,
                        double *C,
                        double *B,
                        int *ldb,
                        double *X,
                        int *ldx,
                        double *RCOND,
                        double *FERR,
                        double *BERR,
                        double *WORK,
                        int *IWORK,
                        int *info);

/* Computes eigenvalues and eigenvectors for the Complex n-by-n Matrix A
 * Use-case example diag(VR*A*inv(VR)) = W.

    \param[in] jobvl calculate left eigenvectors if 'V', otherwise 'N'
    \param[in] jobvl calculate right eigen vectors if 'V', otherwise 'N'
    \param[in] n Dimensionality of the system.
    \param[in,out] A matrix to find eigenvalues and eigenvectors.  Changed in execution.
    \param[in] lda The leading dimension of A.
    \param[out] W/WR/WI The eigenvalues in a vector (double matrix can have complex eigenvalues)
    \param[out] VL The left eigenvectors of A.
    \param[in] ldvl The leading dimension of VL.
    \param[out] VR The right eigenvectors of A.
    \param[in] ldvr The leading dimension of VR.
    \param[out] WORK
    \param[out] IWORK
    \param[out] RWORK
    \param[out] info
 */

extern "C" void dgeev_(char *jobvl,
                       char *jobvr,
                       int *n,
                       double *A,
                       int *lda,
                       double *WR,
                       double *WI,
                       double *VL,
                       int *ldvl,
                       double *VR,
                       int *ldvr,
                       double *work,
                       int *lwork,
                       double *rwork,
                       int *info);

extern "C" void zgeev_(char *jobvl,
                       char *jobvr,
                       int *n,
                       std::complex<double> *A,
                       int *lda,
                       std::complex<double> *W,
                       std::complex<double> *VL,
                       int *ldvl,
                       std::complex<double> *VR,
                       int *ldvr,
                       std::complex<double> *work,
                       int *lwork,
                       double *rwork,
                       int *info);
}  // namespace lapack

#endif  // LAPACK_H