36 #ifndef MADNESS_LINALG_TENSOR_LAPACK_H__INCLUDED
37 #define MADNESS_LINALG_TENSOR_LAPACK_H__INCLUDED
55 void svd(
const Tensor<T>&
a, Tensor<T>& U,
74 return std::make_tuple(U,s,VT);
81 void gesv(
const Tensor<T>&
a,
const Tensor<T>&
b, Tensor<T>& x);
87 void gelss(
const Tensor<T>&
a,
const Tensor<T>&
b,
double rcond,
95 void syev(
const Tensor<T>&
A,
102 template <
typename T>
103 std::tuple<Tensor< typename Tensor<T>::scalar_type >, Tensor<T>>
108 return std::make_tuple(
e,
V);
113 template <
typename T>
114 void geev(
const Tensor<T>&
A, Tensor<T>& VR, Tensor<std::complex<T>>&
e);
118 template <
typename T>
119 void ggev(
const Tensor<T>&
A, Tensor<T>&
B, Tensor<T>& VR,
120 Tensor<std::complex<T>>&
e);
127 template <
typename T>
128 void gesv(
const Tensor<T>&
a,
const Tensor<T>&
b, Tensor<T>& x);
134 template <
typename T>
135 void sygv(
const Tensor<T>&
A,
const Tensor<T>&
B,
int itype,
157 template <
typename T>
163 template <
typename T>
167 template <
typename T>
168 Tensor<T>
inverse(
const Tensor<T>&
A);
172 void qr(Tensor<T>&
A, Tensor<T>&
R);
176 void lq(Tensor<T>&
A, Tensor<T>&
L);
179 void lq_result(Tensor<T>&
A, Tensor<T>&
R, Tensor<T>& tau, Tensor<T>& work,
bool do_qr);
181 template <
typename T>
182 void geqp3(Tensor<T>&
A, Tensor<T>& tau, Tensor<integer>& jpvt);
190 template <
typename T>
191 void orgqr(Tensor<T>&
A,
const Tensor<T>& tau);
Definition: test_ar.cc:118
Definition: test_ar.cc:141
A tensor is a multidimension array.
Definition: tensor.h:317
TensorTypeData< T >::scalar_type scalar_type
C++ typename of the real type associated with a complex type.
Definition: tensor.h:409
static const double R
Definition: csqrt.cc:46
Correspondence between C++ and Fortran types.
File holds all helper structures necessary for the CC_Operator and CC2 class.
Definition: DFParameters.h:10
void rr_cholesky(Tensor< T > &A, typename Tensor< T >::scalar_type tol, Tensor< integer > &piv, int &rank)
Compute the rank-revealing Cholesky factorization.
Definition: lapack.cc:1203
void svd_result(Tensor< T > &a, Tensor< T > &U, Tensor< typename Tensor< T >::scalar_type > &s, Tensor< T > &VT, Tensor< T > &work)
same as svd, but it optimizes away the tensor construction: a = U * diag(s) * VT
Definition: lapack.cc:773
void sygv(const Tensor< T > &A, const Tensor< T > &B, int itype, Tensor< T > &V, Tensor< typename Tensor< T >::scalar_type > &e)
Generalized real-symmetric or complex-Hermitian eigenproblem.
Definition: lapack.cc:1052
void lq_result(Tensor< T > &A, Tensor< T > &R, Tensor< T > &tau, Tensor< T > &work, bool do_qr)
compute the LQ decomposition of the matrix A = L Q
Definition: lapack.cc:1320
void ggev(const Tensor< T > &A, Tensor< T > &B, Tensor< T > &VR, Tensor< std::complex< T >> &e)
Generalized real-non-symmetric or complex-non-Hermitian eigenproblem.
Definition: lapack.cc:1115
void qr(Tensor< T > &A, Tensor< T > &R)
compute the QR decomposition of the matrix A
Definition: lapack.cc:1273
void lq(Tensor< T > &A, Tensor< T > &R)
compute the LQ decomposition of the matrix A = L Q
Definition: lapack.cc:1297
void cholesky(Tensor< T > &A)
Compute the Cholesky factorization.
Definition: lapack.cc:1174
void init_tensor_lapack()
World/MRA initialization calls this before going multithreaded due to static data in dlamch.
Definition: lapack.cc:1622
void gelss(const Tensor< T > &a, const Tensor< T > &b, double rcond, Tensor< T > &x, Tensor< typename Tensor< T >::scalar_type > &s, long &rank, Tensor< typename Tensor< T >::scalar_type > &sumsq)
Solve Ax = b for general A using the LAPACK *gelss routines.
Definition: lapack.cc:884
void gesv(const Tensor< T > &a, const Tensor< T > &b, Tensor< T > &x)
Solve Ax = b for general A using the LAPACK *gesv routines.
Definition: lapack.cc:804
void orgqr(Tensor< T > &A, const Tensor< T > &tau)
reconstruct the orthogonal matrix Q (e.g. from QR factorization)
Definition: lapack.cc:1382
void svd(const Tensor< T > &a, Tensor< T > &U, Tensor< typename Tensor< T >::scalar_type > &s, Tensor< T > &VT)
Compute the singluar value decomposition of an n-by-m matrix using *gesvd.
Definition: lapack.cc:739
void geev(const Tensor< T > &A, Tensor< T > &VR, Tensor< std::complex< T >> &e)
Real non-symmetric or complex non-Hermitian eigenproblem.
Definition: lapack.cc:1008
Tensor< T > inverse(const Tensor< T > &a_in)
invert general square matrix A
Definition: lapack.cc:832
void geqp3(Tensor< T > &A, Tensor< T > &tau, Tensor< integer > &jpvt)
Compute the QR factorization.
Definition: lapack.cc:1235
bool test_tensor_lapack()
Test the Tensor-LAPACK interface ... currently always returns true!
Definition: lapack.cc:1640
void syev(const Tensor< T > &A, Tensor< T > &V, Tensor< typename Tensor< T >::scalar_type > &e)
Real-symmetric or complex-Hermitian eigenproblem.
Definition: lapack.cc:969
static const double b
Definition: nonlinschro.cc:119
static const double a
Definition: nonlinschro.cc:118
static const double L
Definition: rk.cc:46
static double V(const coordT &r)
Definition: tdse.cc:288
Defines and implements most of Tensor.
void e()
Definition: test_sig.cc:75