MADNESS 0.10.1
tensor_lapack.h
Go to the documentation of this file.
1/*
2 This file is part of MADNESS.
3
4 Copyright (C) 2007,2010 Oak Ridge National Laboratory
5
6 This program is free software; you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation; either version 2 of the License, or
9 (at your option) any later version.
10
11 This program is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with this program; if not, write to the Free Software
18 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19
20 For more information please contact:
21
22 Robert J. Harrison
23 Oak Ridge National Laboratory
24 One Bethel Valley Road
25 P.O. Box 2008, MS-6367
26
27 email: harrisonrj@ornl.gov
28 tel: 865-241-3937
29 fax: 865-572-0680
30
31
32 $Id$
33*/
34
35
36#ifndef MADNESS_LINALG_TENSOR_LAPACK_H__INCLUDED
37#define MADNESS_LINALG_TENSOR_LAPACK_H__INCLUDED
38
41
42/*!
43 \file tensor_lapack.h
44 \brief Prototypes for a partial interface from Tensor to LAPACK
45 \ingroup linalg
46@{
47*/
48
49namespace madness {
50
51 /// Computes singular value decomposition of matrix
52
53 /// \ingroup linalg
54 template <typename T>
55 void svd(const Tensor<T>& a, Tensor<T>& U,
56 Tensor< typename Tensor<T>::scalar_type >& s, Tensor<T>& VT);
57
58 /// \ingroup linalg
59 template <typename T>
60 void svd_result(Tensor<T>& a, Tensor<T>& U,
61 Tensor< typename Tensor<T>::scalar_type >& s, Tensor<T>& VT, Tensor<T>& work);
62
63 /// SVD - MATLAB syntax
64
65 /// call as
66 /// auto [U,s,VT] = svd(A);
67 /// with A=U*S*VT
68 template <typename T>
69 std::tuple<Tensor<T>, Tensor< typename Tensor<T>::scalar_type >, Tensor<T>>
70 svd(const Tensor<T>& A) {
71 Tensor<T> U,VT;
73 svd(A,U,s,VT);
74 return std::make_tuple(U,s,VT);
75 }
76
77 /// Solves linear equations
78
79 /// \ingroup linalg
80 template <typename T>
81 void gesv(const Tensor<T>& a, const Tensor<T>& b, Tensor<T>& x);
82
83 /// Solves linear equations using least squares
84
85 /// \ingroup linalg
86 template <typename T>
87 void gelss(const Tensor<T>& a, const Tensor<T>& b, double rcond,
88 Tensor<T>& x, Tensor< typename Tensor<T>::scalar_type >& s,
89 long &rank, Tensor<typename Tensor<T>::scalar_type>& sumsq);
90
91 /// Solves symmetric or Hermitian eigenvalue problem
92
93 /// \ingroup linalg
94 template <typename T>
95 void syev(const Tensor<T>& A,
96 Tensor<T>& V, Tensor< typename Tensor<T>::scalar_type >& e);
97
98 /// Solves symmetric or Hermitian eigenvalue problem - MATLAB syntax
99
100 /// call as
101 /// auto [eval, evec] = syev(A);
102 template <typename T>
103 std::tuple<Tensor< typename Tensor<T>::scalar_type >, Tensor<T>>
104 syev(const Tensor<T>& A) {
105 Tensor<T> V;
107 syev(A,V,e);
108 return std::make_tuple(e,V);
109 }
110 // START BRYAN ADDITION
111 /// Solves non-symmetric or non-Hermitian eigenvalue problem
112
113 template <typename T>
114 void geev(const Tensor<T>& A, Tensor<T>& VR, Tensor<std::complex<T>>& e);
115
116 /// Solves non-symmetric or non-Hermitian generalized eigenvalue problem
117
118 template <typename T>
119 void ggev(const Tensor<T>& A, Tensor<T>& B, Tensor<T>& VR,
120 Tensor<std::complex<T>>& e);
121 // END BRYAN ADDITIONS
122 /// Solves symmetric or Hermitian eigenvalue problem - MATLAB syntax
123
124 /// Solves linear equations
125
126 /// \ingroup linalg
127 template <typename T>
128 void gesv(const Tensor<T>& a, const Tensor<T>& b, Tensor<T>& x);
129
130
131 /// Solves symmetric or Hermitian generalized eigenvalue problem
132
133 /// \ingroup linalg
134 template <typename T>
135 void sygv(const Tensor<T>& A, const Tensor<T>& B, int itype,
136 Tensor<T>& V, Tensor< typename Tensor<T>::scalar_type >& e);
137
138 class World; // UGH!
139 /// Solves symmetric or Hermitian generalized eigenvalue problem
140
141
142 // !!!!!!!!!! sygvp and gesvp are now in the ELEMENTAL inteface
143 // /// \ingroup linalg
144 // template <typename T>
145 // void sygvp(World& world, const Tensor<T>& A, const Tensor<T>& B, int itype,
146 // Tensor<T>& V, Tensor< typename Tensor<T>::scalar_type >& e);
147
148 // /// Solves linear equations
149
150 // /// \ingroup linalg
151 // template <typename T>
152 // void gesvp(World& world, const Tensor<T>& a, const Tensor<T>& b, Tensor<T>& x);
153
154 /// Cholesky factorization
155
156 /// \ingroup linalg
157 template <typename T>
158 void cholesky(Tensor<T>& A);
159
160 /// rank-revealing Cholesky factorization
161
162 /// \ingroup linalg
163 template <typename T>
164 void rr_cholesky(Tensor<T>& A, typename Tensor<T>::scalar_type tol, Tensor<integer>& piv, int& rank);
165
166 /// \ingroup linalg
167 template <typename T>
168 Tensor<T> inverse(const Tensor<T>& A);
169
170 /// QR decomposition
171 template<typename T>
172 void qr(Tensor<T>& A, Tensor<T>& R);
173
174 /// LQ decomposition
175 template<typename T>
176 void lq(Tensor<T>& A, Tensor<T>& L);
177 /// LQ decomposition
178 template<typename T>
179 void lq_result(Tensor<T>& A, Tensor<T>& R, Tensor<T>& tau, Tensor<T>& work,bool do_qr);
180
181 template <typename T>
182 void geqp3(Tensor<T>& A, Tensor<T>& tau, Tensor<integer>& jpvt);
183
184 /// orgqr generates an M-by-N complex matrix Q with orthonormal columns
185
186 /// which is defined as the first N columns of a product of K elementary
187 /// reflectors of order M
188 /// Q = H(1) H(2) . . . H(k)
189 /// as returned by ZGEQRF.
190 template <typename T>
191 void orgqr(Tensor<T>& A, const Tensor<T>& tau);
192
193
194 /// Dunno
195
196// /// \ingroup linalg
197// template <typename T>
198// void triangular_solve(const Tensor<T>& L, Tensor<T>& B,
199// const char* side, const char* transa);
200
201 /// Runs the tensor test code, returns true on success
202
203 /// \ingroup linalg
204 bool test_tensor_lapack();
205
206 /// World/MRA initialization calls this before going multithreaded due to static data in \c dlamch
207
208 /// \ingroup linalg
209 void init_tensor_lapack();
210}
211
212#include <madness/tensor/elem.h>
213
214#endif // MADNESS_LINALG_TENSOR_LAPACK_H__INCLUDED
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.
Namespace for all elements and tools of MADNESS.
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 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
Tensor< T > inverse(const Tensor< T > &a_in)
invert general square matrix A
Definition lapack.cc:832
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
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 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