From 1d7a13d715c5c8cef8a1f1263768f9ce65e05994 Mon Sep 17 00:00:00 2001 From: young Date: Mon, 31 Mar 2014 07:14:35 +0000 Subject: [PATCH] Add an assignment operator from a regular SparseMatrix. git-svn-id: https://svn.dealii.org/trunk@32706 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/deal.II/lac/lapack_full_matrix.h | 21 +++++- deal.II/source/lac/lapack_full_matrix.cc | 67 +++++++++---------- deal.II/source/lac/lapack_full_matrix.inst.in | 3 + 3 files changed, 52 insertions(+), 39 deletions(-) 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 92e0dd09e8..ee8339ba60 100644 --- a/deal.II/include/deal.II/lac/lapack_full_matrix.h +++ b/deal.II/include/deal.II/lac/lapack_full_matrix.h @@ -34,6 +34,7 @@ DEAL_II_NAMESPACE_OPEN template class Vector; template class BlockVector; template class FullMatrix; +template class SparseMatrix; /** @@ -70,6 +71,7 @@ public: */ explicit LAPACKFullMatrix (const size_type n = 0); + /** * Constructor. Initialize the matrix as a rectangular matrix. */ @@ -95,14 +97,23 @@ public: operator = (const LAPACKFullMatrix &); /** - * Assignment operator for a regular FullMatrix. Note that since LAPACK - * expects matrices in transposed order, this transposition is included - * here. + * Assignment operator from a regular FullMatrix. @note Since LAPACK + * expects matrices in transposed order, this transposition is + * included here. */ template LAPACKFullMatrix & operator = (const FullMatrix &); + /** + * Assignment operator from a regular SparseMatrix. @note Since + * LAPACK expects matrices in transposed order, this transposition + * is included here. + */ + template + LAPACKFullMatrix & + operator = (const SparseMatrix &); + /** * This operator assigns a scalar to a matrix. To avoid confusion with * constructors, zero is the only value allowed for d @@ -544,6 +555,10 @@ public: const double threshold = 0.) const; private: + /** + * n_rows + */ + /** * Since LAPACK operations notoriously change the meaning of the matrix * entries, we record the current state after the last operation here. diff --git a/deal.II/source/lac/lapack_full_matrix.cc b/deal.II/source/lac/lapack_full_matrix.cc index 2cf6cbee9a..98aa0a1d52 100644 --- a/deal.II/source/lac/lapack_full_matrix.cc +++ b/deal.II/source/lac/lapack_full_matrix.cc @@ -19,6 +19,7 @@ #include #include #include +#include #include #include @@ -30,34 +31,30 @@ DEAL_II_NAMESPACE_OPEN using namespace LAPACKSupport; template -LAPACKFullMatrix::LAPACKFullMatrix(const size_type n) +LAPACKFullMatrix::LAPACKFullMatrix (const size_type n) : TransposeTable (n,n), - state(matrix) + state (matrix) {} - template -LAPACKFullMatrix::LAPACKFullMatrix( - const size_type m, - const size_type n) +LAPACKFullMatrix::LAPACKFullMatrix (const size_type m, + const size_type n) : TransposeTable (m,n), - state(matrix) + state (matrix) {} - template -LAPACKFullMatrix::LAPACKFullMatrix(const LAPACKFullMatrix &M) +LAPACKFullMatrix::LAPACKFullMatrix (const LAPACKFullMatrix &M) : TransposeTable (M), - state(matrix) + state (matrix) {} - template LAPACKFullMatrix & LAPACKFullMatrix::operator = (const LAPACKFullMatrix &M) @@ -68,7 +65,6 @@ LAPACKFullMatrix::operator = (const LAPACKFullMatrix &M) } - template template LAPACKFullMatrix & @@ -85,6 +81,21 @@ LAPACKFullMatrix::operator = (const FullMatrix &M) } +template +template +LAPACKFullMatrix & +LAPACKFullMatrix::operator = (const SparseMatrix &M) +{ + Assert (this->n_rows() == M.n(), ExcDimensionMismatch(this->n_rows(), M.n())); + Assert (this->n_cols() == M.m(), ExcDimensionMismatch(this->n_cols(), M.m())); + for (size_type i=0; in_rows(); ++i) + for (size_type j=0; jn_cols(); ++j) + (*this)(i,j) = M.el(i,j); + + state = LAPACKSupport::matrix; + return *this; +} + template LAPACKFullMatrix & @@ -100,7 +111,6 @@ LAPACKFullMatrix::operator = (const double d) } - template void LAPACKFullMatrix::vmult ( @@ -159,7 +169,6 @@ LAPACKFullMatrix::vmult ( } - template void LAPACKFullMatrix::Tvmult ( @@ -220,7 +229,6 @@ LAPACKFullMatrix::Tvmult ( } - template void LAPACKFullMatrix::mmult(LAPACKFullMatrix &C, @@ -244,7 +252,6 @@ LAPACKFullMatrix::mmult(LAPACKFullMatrix &C, } - template void LAPACKFullMatrix::mmult(FullMatrix &C, @@ -293,7 +300,6 @@ LAPACKFullMatrix::Tmmult(LAPACKFullMatrix &C, } - template void LAPACKFullMatrix::Tmmult(FullMatrix &C, @@ -318,7 +324,6 @@ LAPACKFullMatrix::Tmmult(FullMatrix &C, } - template void LAPACKFullMatrix::mTmult(LAPACKFullMatrix &C, @@ -367,7 +372,6 @@ LAPACKFullMatrix::mTmult(FullMatrix &C, } - template void LAPACKFullMatrix::TmTmult(LAPACKFullMatrix &C, @@ -391,7 +395,6 @@ LAPACKFullMatrix::TmTmult(LAPACKFullMatrix &C, } - template void LAPACKFullMatrix::TmTmult(FullMatrix &C, @@ -416,7 +419,6 @@ LAPACKFullMatrix::TmTmult(FullMatrix &C, } - template void LAPACKFullMatrix::compute_lu_factorization() @@ -436,7 +438,6 @@ LAPACKFullMatrix::compute_lu_factorization() } - template void LAPACKFullMatrix::compute_svd() @@ -487,7 +488,6 @@ LAPACKFullMatrix::compute_svd() } - template void LAPACKFullMatrix::compute_inverse_svd(const double threshold) @@ -509,7 +509,6 @@ LAPACKFullMatrix::compute_inverse_svd(const double threshold) } - template void LAPACKFullMatrix::invert() @@ -542,7 +541,6 @@ LAPACKFullMatrix::invert() } - template void LAPACKFullMatrix::apply_lu_factorization(Vector &v, @@ -565,7 +563,6 @@ LAPACKFullMatrix::apply_lu_factorization(Vector &v, } - template void LAPACKFullMatrix::apply_lu_factorization(LAPACKFullMatrix &B, @@ -588,12 +585,10 @@ LAPACKFullMatrix::apply_lu_factorization(LAPACKFullMatrix &B, } - template void -LAPACKFullMatrix::compute_eigenvalues( - const bool right, - const bool left) +LAPACKFullMatrix::compute_eigenvalues(const bool right, + const bool left) { Assert(state == matrix, ExcState(state)); const int nn = this->n_cols(); @@ -662,12 +657,11 @@ LAPACKFullMatrix::compute_eigenvalues( template void -LAPACKFullMatrix::compute_eigenvalues_symmetric( - const number lower_bound, - const number upper_bound, - const number abs_accuracy, - Vector &eigenvalues, - FullMatrix &eigenvectors) +LAPACKFullMatrix::compute_eigenvalues_symmetric(const number lower_bound, + const number upper_bound, + const number abs_accuracy, + Vector &eigenvalues, + FullMatrix &eigenvectors) { Assert(state == matrix, ExcState(state)); const int nn = (this->n_cols() > 0 ? this->n_cols() : 1); @@ -968,6 +962,7 @@ LAPACKFullMatrix::Tvmult_add ( Tvmult(w, v, true); } + template void LAPACKFullMatrix::print_formatted ( diff --git a/deal.II/source/lac/lapack_full_matrix.inst.in b/deal.II/source/lac/lapack_full_matrix.inst.in index 9db46db4a0..23d91651bd 100644 --- a/deal.II/source/lac/lapack_full_matrix.inst.in +++ b/deal.II/source/lac/lapack_full_matrix.inst.in @@ -26,4 +26,7 @@ for (S1, S2 : REAL_SCALARS) { template LAPACKFullMatrix & LAPACKFullMatrix::operator = (const FullMatrix &M); + + template LAPACKFullMatrix & + LAPACKFullMatrix::operator = (const SparseMatrix &M); } -- 2.39.5