From: kanschat Date: Fri, 5 Nov 2010 15:30:06 +0000 (+0000) Subject: start some work on SVD X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f161d68f34e1f0d8605c3d9d1ecfc0895b46119d;p=dealii-svn.git start some work on SVD git-svn-id: https://svn.dealii.org/trunk@22611 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/configure.in b/deal.II/configure.in index 7dab4cc12b..f1a5b5c1ab 100644 --- a/deal.II/configure.in +++ b/deal.II/configure.in @@ -621,7 +621,7 @@ dnl exception. AC_CHECK_FUNCS([daxpy_ saxpy_ dgemv_ sgemv_ dgeev_ sgeev_ dgeevx_ sgeevx_]) AC_CHECK_FUNCS([dgesvd_ sgesvd_ dgemm_ sgemm_ dgetrf_ sgetrf_ dgetri_ sgetri_]) AC_CHECK_FUNCS([dgeqrf_ sgeqrf_ dormqr_ sormqr_ dorgqr_ sorgqr_ dtrtrs_ strtrs_]) -AC_CHECK_FUNCS([dgetrs_ sgetrs_ dstev_ sstev_ dsygv_ ssygv_]) +AC_CHECK_FUNCS([dgetrs_ sgetrs_ dstev_ sstev_ dsygv_ ssygv_ dgelsd_ sgelsd_]) dnl ------------------------------------------------------------- dnl set include paths of several libraries diff --git a/deal.II/include/deal.II/lac/lapack_full_matrix.h b/deal.II/include/deal.II/lac/lapack_full_matrix.h index 1d88428a00..52506732b0 100644 --- a/deal.II/include/deal.II/lac/lapack_full_matrix.h +++ b/deal.II/include/deal.II/lac/lapack_full_matrix.h @@ -20,6 +20,7 @@ #include #include +#include #include #include @@ -436,9 +437,9 @@ class LAPACKFullMatrix : public TransposeTable std::vector inv_work; /** - * Real parts of - * eigenvalues. Filled by - * compute_eigenvalues. + * Real parts of eigenvalues or + * the singular values. Filled by + * compute_eigenvalues() or compute_svd(). */ std::vector wr; @@ -460,6 +461,17 @@ class LAPACKFullMatrix : public TransposeTable * can be stored. */ std::vector vr; + + /** + * The matrix U in the singular value decomposition + * USVT. + */ + boost::shared_ptr > svd_u; + /** + * The matrix VT in the singular value decomposition + * USVT. + */ + boost::shared_ptr > svd_vt; }; diff --git a/deal.II/include/deal.II/lac/lapack_support.h b/deal.II/include/deal.II/lac/lapack_support.h index b23eaa14f5..6aa4b9ee05 100644 --- a/deal.II/include/deal.II/lac/lapack_support.h +++ b/deal.II/include/deal.II/lac/lapack_support.h @@ -40,6 +40,8 @@ namespace LAPACKSupport lu, /// Eigenvalue vector is filled eigenvalues, + /// Matrix contains singular value decomposition, + svd, /// Contents is something useless. unusable = 0x8000 }; diff --git a/deal.II/include/deal.II/lac/lapack_templates.h.in b/deal.II/include/deal.II/lac/lapack_templates.h.in index 8ddee8562c..631796654d 100644 --- a/deal.II/include/deal.II/lac/lapack_templates.h.in +++ b/deal.II/include/deal.II/lac/lapack_templates.h.in @@ -77,6 +77,7 @@ void dsygv_ (const int* itype, const char* jobz, const char* uplo, const int* n, double* A, const int* lda, double* B, const int* ldb, double* w, double* work, const int* lwork, int* info); + // Compute singular value decomposition void dgesvd_ (int* jobu, int* jobvt, const int* n, const int* m, double* A, const int* lda, @@ -86,6 +87,14 @@ void dgesvd_ (int* jobu, int* jobvt, double* work, const int* lwork, int* info); +// Solve a least squares problem using SVD +void dgelsd_ (const int* m, const int* n, const int* nrhs, + const double* A, const int* lda, + double* B, const int* ldb, + double* s, const double* rcond, + int* rank, + double* work, const int* lwork, int* iwork, + int* info); // Symmetric tridiagonal matrix void dstev_ (const char* jobz, const int* n,