From 73c6e2399b2379e91c4eed666a8da616833a037a Mon Sep 17 00:00:00 2001 From: kanschat Date: Sun, 28 Nov 2010 23:21:29 +0000 Subject: [PATCH] cleanup git-svn-id: https://svn.dealii.org/trunk@22871 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/deal.II/lac/precondition_block.h | 183 ++++++------------ .../lac/precondition_block.templates.h | 87 +++++++-- 2 files changed, 120 insertions(+), 150 deletions(-) diff --git a/deal.II/include/deal.II/lac/precondition_block.h b/deal.II/include/deal.II/lac/precondition_block.h index 5d75ce1ae3..dfb125d5d2 100644 --- a/deal.II/include/deal.II/lac/precondition_block.h +++ b/deal.II/include/deal.II/lac/precondition_block.h @@ -126,6 +126,22 @@ class PreconditionBlock * are equal to save memory. */ bool same_diagonal; + /** + * Choose the inversion + * method for the blocks. + */ + typename PreconditionBlockBase::Inversion inversion; + + /** + * The if #inversion is SVD, + * the threshold below which + * a singular value will be + * considered zero and thus + * not inverted. This + * parameter is used in the + * call to LAPACKFullMatrix::compute_inverse_svd(). + */ + double threshold; }; @@ -309,6 +325,9 @@ class PreconditionBlock unsigned int block_size () const; /** + * @deprecated Use size() + * instead. + * * The number of blocks of the * matrix. */ @@ -371,11 +390,6 @@ class PreconditionBlock */ double relaxation; - /** - * Number of blocks. - */ - unsigned int nblocks; - /** * The permutation vector */ @@ -547,53 +561,21 @@ class PreconditionBlockJacobi : public virtual Subscriptor, Accessor accessor; }; - /** - * Make initialization function - * publicly available. - */ PreconditionBlock::initialize; - - /** - * Make function of base class public again. - */ PreconditionBlock::clear; - - /** - * Make function of base class public again. - */ PreconditionBlock::empty; - - /** - * Make function of base class public again. - */ PreconditionBlock::el; - - /** - * Make function of base class public again. - */ PreconditionBlock::set_same_diagonal; - - /** - * Make function public. - */ PreconditionBlock::invert_diagblocks; - - /** - * Make function public. - */ PreconditionBlock::block_size; - + PreconditionBlockBase::size; + PreconditionBlockBase::inverse; + PreconditionBlockBase::inverse_householder; + PreconditionBlockBase::inverse_svd; /** - * Make function public. + * @deprecated Use size() instead */ PreconditionBlock::n_blocks; - /** - * Make function accessible to iterator. - */ - PreconditionBlock::inverse; - /** - * Make function accessible. - */ PreconditionBlock::set_permutation; /** @@ -749,42 +731,18 @@ class PreconditionBlockSOR : public virtual Subscriptor, /** * Make type publicly available. */ - using PreconditionBlock::AdditionalData; - - /** - * Make initialization function - * publicly available. - */ - using PreconditionBlock::initialize; - - /** - * Make function of base class public again. - */ - using PreconditionBlock::clear; - - /** - * Make function of base class public again. - */ - using PreconditionBlock::empty; - - /** - * Make function of base class public again. - */ - using PreconditionBlock::el; - - /** - * Make function of base class public again. - */ - using PreconditionBlock::set_same_diagonal; - - /** - * Make function of base class public again. - */ - using PreconditionBlock::invert_diagblocks; - /** - * Make function accessible. - */ - using PreconditionBlock::set_permutation; + PreconditionBlock::AdditionalData; + PreconditionBlock::initialize; + PreconditionBlock::clear; + PreconditionBlock::empty; + PreconditionBlockBase::size; + PreconditionBlockBase::inverse; + PreconditionBlockBase::inverse_householder; + PreconditionBlockBase::inverse_svd; + PreconditionBlock::el; + PreconditionBlock::set_same_diagonal; + PreconditionBlock::invert_diagblocks; + PreconditionBlock::set_permutation; /** * Execute block SOR @@ -957,46 +915,28 @@ class PreconditionBlockSSOR : public virtual Subscriptor, */ PreconditionBlockSSOR (); - /** - * Make type publicly available. - */ - using PreconditionBlockSOR::AdditionalData; + // Keep AdditionalData accessible + PreconditionBlockSOR::AdditionalData; + // The following are the + // functions of the base classes + // which we want to keep + // accessible. /** * Make initialization function * publicly available. */ - using PreconditionBlockSOR::initialize; - - /** - * Make function of base class public again. - */ - using PreconditionBlockSOR::clear; - - /** - * Make function of base class public again. - */ - using PreconditionBlockSOR::set_permutation; - - /** - * Make function of base class public again. - */ - using PreconditionBlockSOR::empty; - - /** - * Make function of base class public again. - */ - using PreconditionBlockSOR::el; - - /** - * Make function of base class public again. - */ - using PreconditionBlockSOR::set_same_diagonal; - - /** - * Make function of base class public again. - */ - using PreconditionBlockSOR::invert_diagblocks; + PreconditionBlockSOR::initialize; + PreconditionBlockSOR::clear; + PreconditionBlockBase::size; + PreconditionBlockBase::inverse; + PreconditionBlockBase::inverse_householder; + PreconditionBlockBase::inverse_svd; + PreconditionBlockSOR::set_permutation; + PreconditionBlockSOR::empty; + PreconditionBlockSOR::el; + PreconditionBlockSOR::set_same_diagonal; + PreconditionBlockSOR::invert_diagblocks; /** * Execute block SSOR @@ -1039,21 +979,6 @@ class PreconditionBlockSSOR : public virtual Subscriptor, #ifndef DOXYGEN -template -inline -PreconditionBlock::AdditionalData:: -AdditionalData (const unsigned int block_size, - const double relaxation, - const bool invert_diagonal, - const bool same_diagonal) - : - relaxation (relaxation), - block_size(block_size), - invert_diagonal(invert_diagonal), - same_diagonal(same_diagonal) -{} - - template inline bool PreconditionBlock::empty () const @@ -1068,7 +993,7 @@ template inline unsigned int PreconditionBlock::n_blocks () const { - return nblocks; + return this->size(); } diff --git a/deal.II/include/deal.II/lac/precondition_block.templates.h b/deal.II/include/deal.II/lac/precondition_block.templates.h index 4e28fa9d23..67b0ab6087 100644 --- a/deal.II/include/deal.II/lac/precondition_block.templates.h +++ b/deal.II/include/deal.II/lac/precondition_block.templates.h @@ -26,6 +26,22 @@ DEAL_II_NAMESPACE_OPEN +template +PreconditionBlock::AdditionalData:: +AdditionalData (const unsigned int block_size, + const double relaxation, + const bool invert_diagonal, + const bool same_diagonal) + : + relaxation (relaxation), + block_size(block_size), + invert_diagonal(invert_diagonal), + same_diagonal(same_diagonal), + inversion(PreconditionBlockBase::gauss_jordan), + threshold(0.) +{} + + template PreconditionBlockBase::~PreconditionBlockBase () {} @@ -67,8 +83,9 @@ void PreconditionBlock::initialize ( Assert (A->m()%bsize==0, ExcWrongBlockSize(bsize, A->m())); blocksize=bsize; relaxation = parameters.relaxation; - nblocks = A->m()/bsize; - this->reinit(nblocks, blocksize, parameters.same_diagonal); + const unsigned int nblocks = A->m()/bsize; + this->reinit(nblocks, blocksize, parameters.same_diagonal, + parameters.inversion); if (parameters.invert_diagonal) invert_diagblocks(); @@ -92,7 +109,7 @@ void PreconditionBlock::initialize ( Assert (A->m()%bsize==0, ExcWrongBlockSize(bsize, A->m())); blocksize=bsize; relaxation = parameters.relaxation; - nblocks = A->m()/bsize; + const unsigned int nblocks = A->m()/bsize; this->reinit(nblocks, blocksize, parameters.same_diagonal); if (parameters.invert_diagonal) @@ -158,7 +175,7 @@ void PreconditionBlock::invert_permuted_diagblocks( // of the unkowns. M_cell = 0; - for (unsigned int cell=0; cellsize(); ++cell) { const unsigned int cell_start = cell*blocksize; for (unsigned int row_cell=0; row_cell::invert_permuted_diagblocks( case PreconditionBlockBase::householder: this->inverse_householder(cell).initialize(M_cell); break; - case PreconditionBlockBase::svd: - this->inverse_svd(cell) = M_cell; - this->inverse_svd(cell).compute_inverse_svd(0.); - break; + case PreconditionBlockBase::svd: + this->inverse_svd(cell) = M_cell; + this->inverse_svd(cell).compute_inverse_svd(0.); + break; default: Assert(false, ExcNotImplemented()); @@ -241,7 +258,7 @@ void PreconditionBlock::forward_step ( FullMatrix M_cell(this->blocksize); // Loop over all blocks - for (unsigned int cell=0; cell < this->nblocks; ++cell) + for (unsigned int cell=0; cell < this->size(); ++cell) { const unsigned int permuted_block_start = permuted ? permutation[block_start] @@ -328,11 +345,11 @@ void PreconditionBlock::backward_step ( // row, column are the global numbering // of the unkowns. unsigned int row, row_cell; - unsigned int block_end=this->blocksize * this->nblocks; + unsigned int block_end=this->blocksize * this->size(); number2 b_cell_row; FullMatrix M_cell(this->blocksize); - for (unsigned int cell=this->nblocks; cell!=0 ;) + for (unsigned int cell=this->size(); cell!=0 ;) { --cell; const unsigned int block_start = block_end - this->blocksize; @@ -433,13 +450,27 @@ void PreconditionBlock::invert_diagblocks() } if (this->store_diagonals()) this->diagonal(0) = M_cell; - this->inverse(0).invert(M_cell); + switch(this->inversion) + { + case PreconditionBlockBase::gauss_jordan: + this->inverse(0).invert(M_cell); + break; + case PreconditionBlockBase::householder: + this->inverse_householder(0).initialize(M_cell); + break; + case PreconditionBlockBase::svd: + this->inverse_svd(0) = M_cell; + this->inverse_svd(0).compute_inverse_svd(1.e-12); + break; + default: + Assert(false, ExcNotImplemented()); + } } else { M_cell = 0; - for (unsigned int cell=0; cellsize(); ++cell) { const unsigned int cell_start = cell*blocksize; for (unsigned int row_cell=0; row_cell::invert_diagblocks() if (this->store_diagonals()) this->diagonal(cell) = M_cell; - this->inverse(cell).invert(M_cell); + switch(this->inversion) + { + case PreconditionBlockBase::gauss_jordan: + this->inverse(cell).invert(M_cell); + break; + case PreconditionBlockBase::householder: + this->inverse_householder(cell).initialize(M_cell); + break; + case PreconditionBlockBase::svd: + this->inverse_svd(cell) = M_cell; + this->inverse_svd(cell).compute_inverse_svd(1.e-12); + break; + default: + Assert(false, ExcNotImplemented()); + } } } this->inverses_computed(true); @@ -543,7 +588,7 @@ void PreconditionBlockJacobi if (!this->inverses_ready()) { FullMatrix M_cell(this->blocksize); - for (unsigned int cell=0; cell < this->nblocks; ++cell) + for (unsigned int cell=0; cell < this->size(); ++cell) { for (row=cell*this->blocksize, row_cell=0; row_cellblocksize; @@ -569,7 +614,7 @@ void PreconditionBlockJacobi } } else - for (unsigned int cell=0; cell < this->nblocks; ++cell) + for (unsigned int cell=0; cell < this->size(); ++cell) { for (row=cell*this->blocksize, row_cell=0; row_cellblocksize; @@ -730,7 +775,7 @@ void PreconditionBlockSOR::forward ( // inverses were not precomputed FullMatrix M_cell(this->blocksize); - for (unsigned int cell=0; cell < this->nblocks; ++cell) + for (unsigned int cell=0; cell < this->size(); ++cell) { const unsigned int permuted_block_start = permuted ? this->permutation[block_start] @@ -833,11 +878,11 @@ void PreconditionBlockSOR::backward ( // row, column are the global numbering // of the unkowns. unsigned int row, row_cell; - unsigned int block_end=this->blocksize * this->nblocks; + unsigned int block_end=this->blocksize * this->size(); number2 b_cell_row; FullMatrix M_cell(this->blocksize); - for (unsigned int cell=this->nblocks; cell!=0 ;) + for (unsigned int cell=this->size(); cell!=0 ;) { --cell; const unsigned int block_start = block_end - this->blocksize; @@ -991,7 +1036,7 @@ void PreconditionBlockSSOR::vmult (Vector &d const double scaling = (2.-this->relaxation)/this->relaxation; // Multiply with diagonal blocks - for (unsigned int cell=0; cell < this->nblocks; ++cell) + for (unsigned int cell=0; cell < this->size(); ++cell) { unsigned int row = cell*this->blocksize; @@ -1022,7 +1067,7 @@ void PreconditionBlockSSOR::Tvmult (Vector & const double scaling = (2.-this->relaxation)/this->relaxation; // Multiply with diagonal blocks - for (unsigned int cell=0; cell < this->nblocks; ++cell) + for (unsigned int cell=0; cell < this->size(); ++cell) { unsigned int row = cell*this->blocksize; -- 2.39.5