From: kanschat Date: Tue, 21 Aug 2007 21:28:55 +0000 (+0000) Subject: remove functions in FullMatrix deprecated for about 2 years X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=2be63c0cb5c1e51fe5344f6716649b3472db2646;p=dealii-svn.git remove functions in FullMatrix deprecated for about 2 years git-svn-id: https://svn.dealii.org/trunk@15008 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/multigrid/mg_coarse.h b/deal.II/deal.II/include/multigrid/mg_coarse.h index bf0c217ec7..3cc764de65 100644 --- a/deal.II/deal.II/include/multigrid/mg_coarse.h +++ b/deal.II/deal.II/include/multigrid/mg_coarse.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -16,6 +16,7 @@ #include #include +#include #include DEAL_II_NAMESPACE_OPEN @@ -156,21 +157,10 @@ class MGCoarseGridHouseholder : public MGCoarseGridBase const VECTOR &src) const; private: - /** - * Pointer to the coarse grid - * matrix. - */ - SmartPointer > matrix; - /** * Matrix for QR-factorization. */ - mutable FullMatrix work; - - /** - * Auxiliary vector. - */ - mutable VECTOR aux; + Householder householder; }; /*@}*/ @@ -274,9 +264,9 @@ MGCoarseGridLACIteration template MGCoarseGridHouseholder::MGCoarseGridHouseholder( const FullMatrix* A) - : - matrix(A, typeid(*this).name()) -{} +{ + if (A != 0) householder.initialize(*A); +} @@ -285,7 +275,7 @@ void MGCoarseGridHouseholder::initialize( const FullMatrix& A) { - matrix = &A; + householder.initialize(A); } @@ -293,9 +283,7 @@ MGCoarseGridHouseholder::initialize( template void MGCoarseGridHouseholder::clear() -{ - matrix = 0; -} +{} @@ -306,9 +294,7 @@ MGCoarseGridHouseholder::operator() ( VECTOR &dst, const VECTOR &src) const { - work = *matrix; - aux = src; - work.least_squares(dst, aux); + householder.least_squares(dst, src); } #endif // DOXYGEN diff --git a/deal.II/lac/include/lac/full_matrix.h b/deal.II/lac/include/lac/full_matrix.h index 7d6e1cda8f..0848a37e6c 100644 --- a/deal.II/lac/include/lac/full_matrix.h +++ b/deal.II/lac/include/lac/full_matrix.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -485,13 +485,7 @@ class FullMatrix : public Table<2,number> * the vector space. */ number frobenius_norm () const; - - /** - * @deprecated Old name for - * FullMatrix::frobenius_norm(). - */ - number norm2 () const; - + /** * Compute the relative norm of * the skew-symmetric part. The @@ -691,23 +685,7 @@ class FullMatrix : public Table<2,number> const FullMatrix &B, const number c, const FullMatrix &C); - - /** - * @deprecated Simple addition of - * a scaled matrix, - * i.e. *this += a*A. - * - * This function is - * deprecated. Use add - * instead, since this has the - * same interface as the other - * matrix and vector classes in - * the library. - */ - template - void add_scaled (const number a, - const FullMatrix &A); - + /** * Add rectangular block. * @@ -917,26 +895,6 @@ class FullMatrix : public Table<2,number> void invert (const FullMatrix &M); - /** - * @deprecated Use the class Householder - * to compute a QR-factorization which - * can be applied to several vectors. - * - * QR-factorization of a matrix. - * The orthogonal transformation - * Q is applied to the vector y - * and this matrix. - * - * After execution of - * householder, the upper - * triangle contains the - * resulting matrix R, the lower - * the incomplete factorization - * matrices. - */ - template - void householder (Vector &y); - //@} ///@name Multiplications //@{ @@ -1090,7 +1048,7 @@ class FullMatrix : public Table<2,number> /** * Forward elimination of lower * triangle. Inverts the lower - * triangle of a quadratic matrix + * triangle of a rectangular matrix * for a given right hand side. * * If the matrix has more columns @@ -1100,28 +1058,9 @@ class FullMatrix : public Table<2,number> * rows, the upper quadratic part * of the matrix is considered. * - * Note that this function does - * not fit into this class at - * all, since it assumes that the - * elements of this object do not - * represent a matrix, but rather - * a decomposition into two - * factors. Therefore, if this - * assumption holds, all - * functions like multiplication - * by matrices or vectors, norms, - * etc, have no meaning any - * more. Conversely, if these - * functions have a meaning on - * this object, then the - * forward() function has no - * meaning. This bifacial - * property of this class is - * probably a design mistake and - * may once go away by separating - * the forward() and backward() - * functions into a class of - * their own. + * @note It is safe to use the + * same object for @p dst and @p + * src. */ template void forward (Vector &dst, @@ -1132,24 +1071,15 @@ class FullMatrix : public Table<2,number> * triangle. * * See forward() + * + * @note It is safe to use the + * same object for @p dst and @p + * src. */ template void backward (Vector &dst, const Vector &src) const; - /** - * @deprecated Use the class - * Householder to solve least - * squares problems. - * - * Least-Squares-Approximation by - * QR-factorization. The return - * value is the Euclidean norm of - * the approximation error. - */ - template - double least_squares (Vector &dst, - Vector &src); //@} /** @addtogroup Exceptions diff --git a/deal.II/lac/include/lac/full_matrix.templates.h b/deal.II/lac/include/lac/full_matrix.templates.h index 893aeb3401..abaead2eca 100644 --- a/deal.II/lac/include/lac/full_matrix.templates.h +++ b/deal.II/lac/include/lac/full_matrix.templates.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -1609,20 +1609,6 @@ void FullMatrix::Tadd (const FullMatrix &src, -template -template -void -FullMatrix::add_scaled (const number s, - const FullMatrix &src) -{ - // this function is - // deprecated. forward to the other - // one - add (s, src); -} - - - template template void @@ -1956,15 +1942,6 @@ FullMatrix::frobenius_norm () const -template -number -FullMatrix::norm2 () const -{ - return frobenius_norm(); -} - - - template number FullMatrix::relative_symmetry_norm2 () const @@ -2321,71 +2298,6 @@ FullMatrix::gauss_jordan () } } -// QR-transformation cf. Stoer 1 4.8.2 (p. 191) - -template -template -void -FullMatrix::householder(Vector& src) -{ - Assert (!this->empty(), ExcEmptyMatrix()); - - // m > n, src.n() = m - Assert (this->n_cols() <= this->n_rows(), - ExcDimensionMismatch(this->n_cols(), this->n_rows())); - Assert (src.size() == this->n_rows(), - ExcDimensionMismatch(src.size(), this->n_rows())); - - for (unsigned int j=0 ; jel(i,j)*this->el(i,j); - if (std::fabs(sigma) < 1.e-15) return; - number2 s = this->el(j,j); - s = (s<0) ? std::sqrt(sigma) : -std::sqrt(sigma); - number2 dj = s; - - number2 beta = 1./(s*this->el(j,j)-sigma); - this->el(j,j) -= s; - - for (unsigned int k=j+1 ; kel(i,j)*this->el(i,k); - sum *= beta; - - for (i=j ; iel(i,k) += sum*this->el(i,j); - } - - number2 sum = 0.; - for (i=j ; iel(i,j)*src(i); - sum *= beta; - - for (i=j ; iel(i,j); - this->el(j,j) = dj; - } -} - - -template -template -double -FullMatrix::least_squares (Vector& dst, - Vector& src) -{ - Assert (!this->empty(), ExcEmptyMatrix()); - - // m > n, m = src.n, n = dst.n - - householder(src); - backward(dst, src); - - number2 sum = 0.; - for (unsigned int i=n() ; i diff --git a/deal.II/lac/include/lac/householder.h b/deal.II/lac/include/lac/householder.h index f44a4c4426..943548816d 100644 --- a/deal.II/lac/include/lac/householder.h +++ b/deal.II/lac/include/lac/householder.h @@ -51,6 +51,11 @@ template class Householder : private FullMatrix { public: + /** + * Create an empty object. + */ + Householder (); + /** * Create an object holding the * QR-decomposition of a matrix. @@ -58,6 +63,14 @@ class Householder : private FullMatrix template Householder (const FullMatrix&); + /** + * Compute the QR-decomposition + * of another matrix. + */ + template + void + initialize (const FullMatrix&); + /** * Solve the least-squares * problem for the right hand @@ -77,7 +90,7 @@ class Householder : private FullMatrix */ template double least_squares (Vector &dst, - Vector &src) const; + const Vector &src) const; private: /** @@ -95,15 +108,20 @@ class Householder : private FullMatrix // QR-transformation cf. Stoer 1 4.8.2 (p. 191) +template +Householder::Householder() +{} + + template template -Householder::Householder(const FullMatrix& M) - : - FullMatrix(M), - diagonal(M.n_rows()) +void +Householder::initialize(const FullMatrix& M) { -// Assert (!this->empty(), ExcEmptyMatrix()); - + this->reinit(M.n_rows(), M.n_cols()); + this->fill(M); + diagonal.resize(M.n_rows()); + Assert (!this->empty(), typename FullMatrix::ExcEmptyMatrix()); // m > n, src.n() = m Assert (this->n_cols() <= this->n_rows(), ExcDimensionMismatch(this->n_cols(), this->n_rows())); @@ -149,35 +167,44 @@ Householder::Householder(const FullMatrix& M) } +template +template +Householder::Householder(const FullMatrix& M) +{ + initialize(M); +} + + template template double Householder::least_squares (Vector& dst, - Vector& src) const + const Vector& src) const { -// Assert (!this->empty(), ExcEmptyMatrix()); - + Assert (!this->empty(), typename FullMatrix::ExcEmptyMatrix()); + dst = src; // m > n, m = src.n, n = dst.n // Multiply Q_n ... Q_2 Q_1 src // Where Q_i = I-v_i v_i^T for (unsigned int j=0;jn();++j) { - // sum = v_i^T src - number2 sum = diagonal[j]*src(j); + // sum = v_i^T dst + number2 sum = diagonal[j]*dst(j); for (unsigned int i=j+1 ; im() ; ++i) - sum += this->el(i,j)*src(i); - // src -= v * sum - src(j) -= sum*diagonal[j]; + sum += this->el(i,j)*dst(i); + // dst -= v * sum + dst(j) -= sum*diagonal[j]; for (unsigned int i=j+1 ; im() ; ++i) - src(i) -= sum*this->el(i,j); + dst(i) -= sum*this->el(i,j); } - - this->backward(dst, src); - + // Compute norm of residual number2 sum = 0.; for (unsigned int i=this->n() ; im() ; ++i) - sum += src(i) * src(i); + sum += dst(i) * dst(i); + // Compute solution + this->backward(dst, dst); + return std::sqrt(sum); } diff --git a/deal.II/lac/include/lac/precondition_block.templates.h b/deal.II/lac/include/lac/precondition_block.templates.h index fd56919f6a..255be19c3e 100644 --- a/deal.II/lac/include/lac/precondition_block.templates.h +++ b/deal.II/lac/include/lac/precondition_block.templates.h @@ -18,6 +18,7 @@ #include #include #include +#include #include #include #include @@ -446,8 +447,8 @@ void PreconditionBlockJacobi column_cellblocksize; ++column_cell, ++column) M_cell(row_cell,column_cell)=M(row,column); } - M_cell.householder(b_cell); - M_cell.backward(x_cell,b_cell); + Householder house(M_cell); + house.least_squares(x_cell,b_cell); // distribute x_cell to dst for (row=cell*this->blocksize, row_cell=0; row_cellblocksize; @@ -637,8 +638,8 @@ void PreconditionBlockSOR::forward ( } else { - M_cell.householder(b_cell); - M_cell.backward(x_cell,b_cell); + Householder house(M_cell); + house.least_squares(x_cell,b_cell); } // distribute x_cell to dst @@ -751,8 +752,8 @@ void PreconditionBlockSOR::backward ( } else { - M_cell.householder(b_cell); - M_cell.backward(x_cell,b_cell); + Householder house(M_cell); + house.least_squares(x_cell,b_cell); } diff --git a/deal.II/lac/include/lac/solver_gmres.h b/deal.II/lac/include/lac/solver_gmres.h index ce8ad80c63..9256bb7f08 100644 --- a/deal.II/lac/include/lac/solver_gmres.h +++ b/deal.II/lac/include/lac/solver_gmres.h @@ -18,6 +18,7 @@ #include #include #include +#include #include #include #include @@ -930,8 +931,8 @@ SolverFGMRES::solve ( y.reinit(j); projected_rhs(0) = beta; H1.fill(H); - - double res = H1.least_squares(y, projected_rhs); + Householder house(H1); + double res = house.least_squares(y, projected_rhs); iteration_state = this->control().check(++accumulated_iterations, res); if (iteration_state != SolverControl::iterate) break; diff --git a/deal.II/lac/source/full_matrix.double.cc b/deal.II/lac/source/full_matrix.double.cc index 6a0ca614c1..29c5e79095 100644 --- a/deal.II/lac/source/full_matrix.double.cc +++ b/deal.II/lac/source/full_matrix.double.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -39,7 +39,6 @@ template void FullMatrix::add (const TYPEMAT, const FullMatri template void FullMatrix::add (const TYPEMAT, const FullMatrix&, const TYPEMAT, const FullMatrix&, const TYPEMAT, const FullMatrix&); -template void FullMatrix::add_scaled (const TYPEMAT, const FullMatrix&); template void FullMatrix::add ( const FullMatrix&, double, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::Tadd (const TYPEMAT, const FullMatrix&); @@ -77,9 +76,6 @@ template void FullMatrix::forward( Vector&, const Vector&) const; template void FullMatrix::backward( Vector&, const Vector&) const; -template void FullMatrix::householder(Vector&); -template double FullMatrix::least_squares( - Vector&, Vector&); template void FullMatrix::precondition_Jacobi ( @@ -106,9 +102,6 @@ template void FullMatrix::forward( Vector&, const Vector&) const; template void FullMatrix::backward( Vector&, const Vector&) const; -template void FullMatrix::householder(Vector&); -template double FullMatrix::least_squares( - Vector&, Vector&); template void FullMatrix::precondition_Jacobi ( Vector &, const Vector &, const TYPEMAT) const; diff --git a/deal.II/lac/source/full_matrix.float.cc b/deal.II/lac/source/full_matrix.float.cc index 14754a8b08..12030d4481 100644 --- a/deal.II/lac/source/full_matrix.float.cc +++ b/deal.II/lac/source/full_matrix.float.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -40,7 +40,6 @@ template void FullMatrix::add (const TYPEMAT, const FullMatri template void FullMatrix::add (const TYPEMAT, const FullMatrix&, const TYPEMAT, const FullMatrix&, const TYPEMAT, const FullMatrix&); -template void FullMatrix::add_scaled (const TYPEMAT, const FullMatrix&); template void FullMatrix::add ( const FullMatrix&, double, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::Tadd (const TYPEMAT, const FullMatrix&); @@ -78,9 +77,6 @@ template void FullMatrix::forward( Vector&, const Vector&) const; template void FullMatrix::backward( Vector&, const Vector&) const; -template void FullMatrix::householder(Vector&); -template double FullMatrix::least_squares( - Vector&, Vector&); template void FullMatrix::precondition_Jacobi ( @@ -107,9 +103,6 @@ template void FullMatrix::forward( Vector&, const Vector&) const; template void FullMatrix::backward( Vector&, const Vector&) const; -template void FullMatrix::householder(Vector&); -template double FullMatrix::least_squares( - Vector&, Vector&); template void FullMatrix::precondition_Jacobi ( Vector &, const Vector &, const TYPEMAT) const;