[136] | 1 | // This file is part of Eigen, a lightweight C++ template library
|
---|
| 2 | // for linear algebra.
|
---|
| 3 | //
|
---|
| 4 | // Copyright (C) 2009-2010 Gael Guennebaud <gael.guennebaud@inria.fr>
|
---|
| 5 | //
|
---|
| 6 | // This Source Code Form is subject to the terms of the Mozilla
|
---|
| 7 | // Public License v. 2.0. If a copy of the MPL was not distributed
|
---|
| 8 | // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
|
---|
| 9 |
|
---|
| 10 | #ifndef EIGEN_BLAS_COMMON_H
|
---|
| 11 | #define EIGEN_BLAS_COMMON_H
|
---|
| 12 |
|
---|
| 13 | #include <Eigen/Core>
|
---|
| 14 | #include <Eigen/Jacobi>
|
---|
| 15 |
|
---|
| 16 | #include <iostream>
|
---|
| 17 | #include <complex>
|
---|
| 18 |
|
---|
| 19 | #ifndef SCALAR
|
---|
| 20 | #error the token SCALAR must be defined to compile this file
|
---|
| 21 | #endif
|
---|
| 22 |
|
---|
| 23 | #include <Eigen/src/misc/blas.h>
|
---|
| 24 |
|
---|
| 25 |
|
---|
| 26 | #define NOTR 0
|
---|
| 27 | #define TR 1
|
---|
| 28 | #define ADJ 2
|
---|
| 29 |
|
---|
| 30 | #define LEFT 0
|
---|
| 31 | #define RIGHT 1
|
---|
| 32 |
|
---|
| 33 | #define UP 0
|
---|
| 34 | #define LO 1
|
---|
| 35 |
|
---|
| 36 | #define NUNIT 0
|
---|
| 37 | #define UNIT 1
|
---|
| 38 |
|
---|
| 39 | #define INVALID 0xff
|
---|
| 40 |
|
---|
| 41 | #define OP(X) ( ((X)=='N' || (X)=='n') ? NOTR \
|
---|
| 42 | : ((X)=='T' || (X)=='t') ? TR \
|
---|
| 43 | : ((X)=='C' || (X)=='c') ? ADJ \
|
---|
| 44 | : INVALID)
|
---|
| 45 |
|
---|
| 46 | #define SIDE(X) ( ((X)=='L' || (X)=='l') ? LEFT \
|
---|
| 47 | : ((X)=='R' || (X)=='r') ? RIGHT \
|
---|
| 48 | : INVALID)
|
---|
| 49 |
|
---|
| 50 | #define UPLO(X) ( ((X)=='U' || (X)=='u') ? UP \
|
---|
| 51 | : ((X)=='L' || (X)=='l') ? LO \
|
---|
| 52 | : INVALID)
|
---|
| 53 |
|
---|
| 54 | #define DIAG(X) ( ((X)=='N' || (X)=='n') ? NUNIT \
|
---|
| 55 | : ((X)=='U' || (X)=='u') ? UNIT \
|
---|
| 56 | : INVALID)
|
---|
| 57 |
|
---|
| 58 |
|
---|
| 59 | inline bool check_op(const char* op)
|
---|
| 60 | {
|
---|
| 61 | return OP(*op)!=0xff;
|
---|
| 62 | }
|
---|
| 63 |
|
---|
| 64 | inline bool check_side(const char* side)
|
---|
| 65 | {
|
---|
| 66 | return SIDE(*side)!=0xff;
|
---|
| 67 | }
|
---|
| 68 |
|
---|
| 69 | inline bool check_uplo(const char* uplo)
|
---|
| 70 | {
|
---|
| 71 | return UPLO(*uplo)!=0xff;
|
---|
| 72 | }
|
---|
| 73 |
|
---|
| 74 |
|
---|
| 75 | namespace Eigen {
|
---|
| 76 | #include "BandTriangularSolver.h"
|
---|
| 77 | #include "GeneralRank1Update.h"
|
---|
| 78 | #include "PackedSelfadjointProduct.h"
|
---|
| 79 | #include "PackedTriangularMatrixVector.h"
|
---|
| 80 | #include "PackedTriangularSolverVector.h"
|
---|
| 81 | #include "Rank2Update.h"
|
---|
| 82 | }
|
---|
| 83 |
|
---|
| 84 | using namespace Eigen;
|
---|
| 85 |
|
---|
| 86 | typedef SCALAR Scalar;
|
---|
| 87 | typedef NumTraits<Scalar>::Real RealScalar;
|
---|
| 88 | typedef std::complex<RealScalar> Complex;
|
---|
| 89 |
|
---|
| 90 | enum
|
---|
| 91 | {
|
---|
| 92 | IsComplex = Eigen::NumTraits<SCALAR>::IsComplex,
|
---|
| 93 | Conj = IsComplex
|
---|
| 94 | };
|
---|
| 95 |
|
---|
| 96 | typedef Matrix<Scalar,Dynamic,Dynamic,ColMajor> PlainMatrixType;
|
---|
| 97 | typedef Map<Matrix<Scalar,Dynamic,Dynamic,ColMajor>, 0, OuterStride<> > MatrixType;
|
---|
| 98 | typedef Map<Matrix<Scalar,Dynamic,1>, 0, InnerStride<Dynamic> > StridedVectorType;
|
---|
| 99 | typedef Map<Matrix<Scalar,Dynamic,1> > CompactVectorType;
|
---|
| 100 |
|
---|
| 101 | template<typename T>
|
---|
| 102 | Map<Matrix<T,Dynamic,Dynamic,ColMajor>, 0, OuterStride<> >
|
---|
| 103 | matrix(T* data, int rows, int cols, int stride)
|
---|
| 104 | {
|
---|
| 105 | return Map<Matrix<T,Dynamic,Dynamic,ColMajor>, 0, OuterStride<> >(data, rows, cols, OuterStride<>(stride));
|
---|
| 106 | }
|
---|
| 107 |
|
---|
| 108 | template<typename T>
|
---|
| 109 | Map<Matrix<T,Dynamic,1>, 0, InnerStride<Dynamic> > vector(T* data, int size, int incr)
|
---|
| 110 | {
|
---|
| 111 | return Map<Matrix<T,Dynamic,1>, 0, InnerStride<Dynamic> >(data, size, InnerStride<Dynamic>(incr));
|
---|
| 112 | }
|
---|
| 113 |
|
---|
| 114 | template<typename T>
|
---|
| 115 | Map<Matrix<T,Dynamic,1> > vector(T* data, int size)
|
---|
| 116 | {
|
---|
| 117 | return Map<Matrix<T,Dynamic,1> >(data, size);
|
---|
| 118 | }
|
---|
| 119 |
|
---|
| 120 | template<typename T>
|
---|
| 121 | T* get_compact_vector(T* x, int n, int incx)
|
---|
| 122 | {
|
---|
| 123 | if(incx==1)
|
---|
| 124 | return x;
|
---|
| 125 |
|
---|
| 126 | T* ret = new Scalar[n];
|
---|
| 127 | if(incx<0) vector(ret,n) = vector(x,n,-incx).reverse();
|
---|
| 128 | else vector(ret,n) = vector(x,n, incx);
|
---|
| 129 | return ret;
|
---|
| 130 | }
|
---|
| 131 |
|
---|
| 132 | template<typename T>
|
---|
| 133 | T* copy_back(T* x_cpy, T* x, int n, int incx)
|
---|
| 134 | {
|
---|
| 135 | if(x_cpy==x)
|
---|
| 136 | return 0;
|
---|
| 137 |
|
---|
| 138 | if(incx<0) vector(x,n,-incx).reverse() = vector(x_cpy,n);
|
---|
| 139 | else vector(x,n, incx) = vector(x_cpy,n);
|
---|
| 140 | return x_cpy;
|
---|
| 141 | }
|
---|
| 142 |
|
---|
| 143 | #define EIGEN_BLAS_FUNC(X) EIGEN_CAT(SCALAR_SUFFIX,X##_)
|
---|
| 144 |
|
---|
| 145 | #endif // EIGEN_BLAS_COMMON_H
|
---|