Fork me on GitHub
 All Classes Files Functions Variables Groups Pages
lapack.h
Go to the documentation of this file.
1 /*
2  * Copyright 2008-2009 NVIDIA Corporation
3  *
4  * Licensed under the Apache License, Version 2.0 (the "License");
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
8  * http://www.apache.org/licenses/LICENSE-2.0
9  *
10  * Unless required by applicable law or agreed to in writing, software
11  * distributed under the License is distributed on an "AS IS" BASIS,
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16 
21 #pragma once
22 
23 #include <cusp/detail/config.h>
24 
25 #include <cusp/lapack/detail/defs.h>
26 
27 namespace cusp
28 {
29 namespace lapack
30 {
31 
41 template<DerivedPolicy, typename Array2d, typename Array1d>
42 void getrf(const thrust::detail::execution_policy_base<DerivedPolicy> &exec,
43  Array2d& A, Array1d& piv );
90 template<typename Array2d, typename Array1d>
91 void getrf( Array2d& A, Array1d& piv );
92 
94 template<DerivedPolicy, typename Array2d>
95 void potrf(const thrust::detail::execution_policy_base<DerivedPolicy> &exec,
96  Array2d& A, char uplo = 'U' );
141 template<typename Array2d>
142 void potrf( Array2d& A, char uplo = 'U' );
143 
145 template<DerivedPolicy, typename Array2d, typename Array1d>
146 void sytrf(const thrust::detail::execution_policy_base<DerivedPolicy> &exec,
147  Array2d& A, Array1d& piv, char uplo = 'U' );
200 template<typename Array2d, typename Array1d>
201 void sytrf( Array2d& A, Array1d& piv, char uplo = 'U' );
202 
204 template<DerivedPolicy, typename Array2d, typename Array1d>
205 void getrs(const thrust::detail::execution_policy_base<DerivedPolicy> &exec,
206  Array2d& A, Array1d& piv, Array2d& B, char trans = 'N' );
263 template<typename Array2d, typename Array1d>
264 void getrs( const Array2d& A, const Array1d& piv, Array2d& B, char trans = 'N' );
265 
325 template<typename Array2d>
326 void potrs( const Array2d& A, Array2d& B, char uplo = 'U');
327 
390 template<typename Array2d, typename Array1d>
391 void sytrs( const Array2d& A, const Array1d& piv, Array2d& B, char uplo = 'U' );
392 
443 template<typename Array2d>
444 void trtrs( const Array2d& A, Array2d& B, char uplo = 'U', char trans = 'N', char diag = 'N' );
445 
479 template<typename Array2d>
480 void trtri( Array2d& A, char uplo = 'U', char diag = 'N' );
481 
526 template<typename Array2d, typename Array1d>
527 void syev( const Array2d& A, Array1d& eigvals, Array2d& eigvecs, char uplo = 'U' );
528 
575 template<typename Array1d1, typename Array1d2, typename Array1d3, typename Array2d>
576 void stev( const Array1d1& alphas, const Array1d2& betas,
577  Array1d3& eigvals, Array2d& eigvecs, char job = 'V' );
578 
629 template<typename Array2d1, typename Array2d2, typename Array1d, typename Array2d3>
630 void sygv( const Array2d1& A, const Array2d2& B, Array1d& eigvals, Array2d3& eigvecs );
631 
682 template<typename Array2d, typename Array1d>
683 void gesv( const Array2d& A, Array2d& B, Array1d& piv );
684 
688 } // end namespace lapack
689 } // end namespace cusp
690 
691 #include <cusp/lapack/detail/lapack.inl>
void sygv(const Array2d1 &A, const Array2d2 &B, Array1d &eigvals, Array2d3 &eigvecs)
Computes all eigenvalues and, optionally, eigenvectors of a real generalized symmetric definite eigen...
void syev(const Array2d &A, Array1d &eigvals, Array2d &eigvecs, char uplo= 'U')
Computes all eigenvalues and eigenvectors of a real symmetric matrix.
void gesv(const Array2d &A, Array2d &B, Array1d &piv)
Computes the solution of a system of linear equations with a square matrix A and multiple right-hand ...
void potrs(const Array2d &A, Array2d &B, char uplo= 'U')
Solves a system of linear equations with a Cholesky-factored symmetric (Hermitian) positive-definite ...
void sytrs(const Array2d &A, const Array1d &piv, Array2d &B, char uplo= 'U')
Solves a system of linear equations with a UDU- or LDL-factored symmetric matrix. ...
void getrf(Array2d &A, Array1d &piv)
Compute LU factorization of matrix.
void stev(const Array1d1 &alphas, const Array1d2 &betas, Array1d3 &eigvals, Array2d &eigvecs, char job= 'V')
Computes all eigenvalues and, optionally, eigenvectors of a real symmetric tridiagonal matrix...
void trtri(Array2d &A, char uplo= 'U', char diag= 'N')
This routine computes the inverse of a triangular matrix.
void sytrf(Array2d &A, Array1d &piv, char uplo= 'U')
Computes the Bunch-Kaufman factorization of a symmetric matrix.
void potrf(Array2d &A, char uplo= 'U')
Computes the Cholesky factorization of a symmetric (Hermitian) positive-definite matrix.
void trtrs(const Array2d &A, Array2d &B, char uplo= 'U', char trans= 'N', char diag= 'N')
Solves a system of linear equations with a triangular matrix, with multiple right-hand sides...
void getrs(const Array2d &A, const Array1d &piv, Array2d &B, char trans= 'N')
Solves a system of linear equations with an LU-factored square matrix, with multiple right-hand sides...