From: Timo Heister Date: Tue, 2 Aug 2016 19:04:01 +0000 (+0200) Subject: implement parallel RelaxationBlock X-Git-Tag: v8.5.0-rc1~792^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3f3d4612bf90c4bfa6efeae3462acd3f7e3754f8;p=dealii.git implement parallel RelaxationBlock - use a ghosted vector in RelaxationBlock - introduce VectorType in base class RelaxationBlock - update test for parallel RelaxationBlock - add documentation --- diff --git a/doc/news/changes.h b/doc/news/changes.h index 2e59d52a6b..6444cdfac0 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -428,6 +428,12 @@ SparsityPattern.

Specific improvements

    +
  1. New: RelaxationBlock classes for geometric multigrid now support parallel + computations using Trilinos. +
    + (Timo Heister, Guido Kanschat, 2016/08/08) +
  2. +
  3. Improved: The regular and hp versions of DoFTools::make_flux_sparsity_pattern() no longer use the user flags of the underlying triangulation to determine if entries along a certain face have been diff --git a/include/deal.II/lac/relaxation_block.h b/include/deal.II/lac/relaxation_block.h index ce5597255d..696fea9206 100644 --- a/include/deal.II/lac/relaxation_block.h +++ b/include/deal.II/lac/relaxation_block.h @@ -43,13 +43,18 @@ DEAL_II_NAMESPACE_OPEN * implementation relies on a straight forward implementation of the Gauss- * Seidel process. * + * Parallel computations require you to specify an initialized + * ghost vector in AdditionalData::temp_ghost_vector. + * * @ingroup Preconditioners * @author Guido Kanschat * @date 2010 */ -template +template > class RelaxationBlock : - protected PreconditionBlockBase + protected PreconditionBlockBase { private: /** @@ -60,7 +65,7 @@ private: /** * Value type for inverse matrices. */ - typedef inverse_type value_type; + typedef InverseNumberType value_type; public: /** @@ -82,7 +87,11 @@ public: */ AdditionalData (const double relaxation = 1., const bool invert_diagonal = true, - const bool same_diagonal = false); + const bool same_diagonal = false, + const typename PreconditionBlockBase::Inversion inversion + = PreconditionBlockBase::gauss_jordan, + const double threshold = 0., + VectorType *temp_ghost_vector = NULL); /** * The mapping from indices to blocks. Each row of this pattern enumerates @@ -112,10 +121,11 @@ public: * particular if their sizes differ. */ bool same_diagonal; + /** * Choose the inversion method for the blocks. */ - typename PreconditionBlockBase::Inversion inversion; + typename PreconditionBlockBase::Inversion inversion; /** * If #inversion is SVD, we can compute the Penrose-Moore inverse of the @@ -147,6 +157,17 @@ public: *
*/ std::vector > order; + + /** + * Temporary ghost vector that is used in the relaxation method when + * performing parallel MPI computations. The user is required to have this + * point to an initialized vector that contains all indices + * that appear in the @p block_list sa ghost values. Typically, this the + * set of locally active level DoFs. Unused when VectorType is a serial + * vector type like Vector. + */ + mutable VectorType *temp_ghost_vector; + /** * Return the memory allocated in this object. */ @@ -210,24 +231,24 @@ protected: * vectors). For the Jacobi step, the calling function must copy @p dst to * @p pref after this. */ - template void do_step ( - VECTOR &dst, - const VECTOR &prev, - const VECTOR &src, + VectorType &dst, + const VectorType &prev, + const VectorType &src, const bool backward) const; + /** * Pointer to the matrix. Make sure that the matrix exists as long as this * class needs it, i.e. until calling @p invert_diagblocks, or (if the * inverse matrices should not be stored) until the last call of the * preconditioning @p vmult function of the derived classes. */ - SmartPointer > A; + SmartPointer > A; /** * Control information. */ - SmartPointer > additional_data; + SmartPointer > additional_data; }; @@ -246,9 +267,11 @@ protected: * @author Guido Kanschat * @date 2010 */ -template +template > class RelaxationBlockJacobi : public virtual Subscriptor, - protected RelaxationBlock + protected RelaxationBlock { public: /** @@ -264,50 +287,51 @@ public: /** * Make type publicly available. */ - using typename RelaxationBlock::AdditionalData; + using typename RelaxationBlock::AdditionalData; /** * Make initialization function publicly available. */ - using RelaxationBlock::initialize; + using RelaxationBlock::initialize; /** * Make function of base class public again. */ - using RelaxationBlock::clear; + using RelaxationBlock::clear; /** * Make function of base class public again. */ - using RelaxationBlock::empty; + using RelaxationBlock::empty; + /** + * Make function of base class public again. + */ + using RelaxationBlock::size; /** * Make function of base class public again. */ - using RelaxationBlock::size; + using RelaxationBlock::inverse; /** * Make function of base class public again. */ - using RelaxationBlock::inverse; + using RelaxationBlock::inverse_householder; /** * Make function of base class public again. */ - using RelaxationBlock::inverse_householder; + using RelaxationBlock::inverse_svd; /** * Make function of base class public again. */ - using RelaxationBlock::inverse_svd; - using PreconditionBlockBase::log_statistics; + using PreconditionBlockBase::log_statistics; /** * Perform one step of the Jacobi iteration. */ - template - void step (VECTOR &dst, const VECTOR &rhs) const; + void step (VectorType &dst, const VectorType &rhs) const; /** * Perform one step of the Jacobi iteration. */ - template - void Tstep (VECTOR &dst, const VECTOR &rhs) const; + void Tstep (VectorType &dst, const VectorType &rhs) const; /** * Return the memory allocated in this object. @@ -331,9 +355,11 @@ public: * @author Guido Kanschat * @date 2010 */ -template +template > class RelaxationBlockSOR : public virtual Subscriptor, - protected RelaxationBlock + protected RelaxationBlock { public: /** @@ -349,50 +375,51 @@ public: /** * Make type publicly available. */ - using typename RelaxationBlock::AdditionalData; + using typename RelaxationBlock::AdditionalData; /** * Make initialization function publicly available. */ - using RelaxationBlock::initialize; + using RelaxationBlock::initialize; /** * Make function of base class public again. */ - using RelaxationBlock::clear; + using RelaxationBlock::clear; /** * Make function of base class public again. */ - using RelaxationBlock::empty; + using RelaxationBlock::empty; + /** + * Make function of base class public again. + */ + using RelaxationBlock::size; /** * Make function of base class public again. */ - using RelaxationBlock::size; + using RelaxationBlock::inverse; /** * Make function of base class public again. */ - using RelaxationBlock::inverse; + using RelaxationBlock::inverse_householder; /** * Make function of base class public again. */ - using RelaxationBlock::inverse_householder; + using RelaxationBlock::inverse_svd; /** * Make function of base class public again. */ - using RelaxationBlock::inverse_svd; - using PreconditionBlockBase::log_statistics; + using PreconditionBlockBase::log_statistics; /** * Perform one step of the SOR iteration. */ - template - void step (VECTOR &dst, const VECTOR &rhs) const; + void step (VectorType &dst, const VectorType &rhs) const; /** * Perform one step of the transposed SOR iteration. */ - template - void Tstep (VECTOR &dst, const VECTOR &rhs) const; + void Tstep (VectorType &dst, const VectorType &rhs) const; }; @@ -411,9 +438,11 @@ public: * @author Guido Kanschat * @date 2010 */ -template +template > class RelaxationBlockSSOR : public virtual Subscriptor, - protected RelaxationBlock + protected RelaxationBlock { public: /** @@ -424,51 +453,52 @@ public: /** * Make type publicly available. */ - using typename RelaxationBlock::AdditionalData; + using typename RelaxationBlock::AdditionalData; /** * Make initialization function publicly available. */ - using RelaxationBlock::initialize; + using RelaxationBlock::initialize; /** * Make function of base class public again. */ - using RelaxationBlock::clear; + using RelaxationBlock::clear; /** * Make function of base class public again. */ - using RelaxationBlock::empty; + using RelaxationBlock::empty; /** * Make function of base class public again. */ - using RelaxationBlock::size; + using RelaxationBlock::size; + /** + * Make function of base class public again. + */ + using RelaxationBlock::inverse; /** * Make function of base class public again. */ - using RelaxationBlock::inverse; + using RelaxationBlock::inverse_householder; /** * Make function of base class public again. */ - using RelaxationBlock::inverse_householder; + using RelaxationBlock::inverse_svd; /** * Make function of base class public again. */ - using RelaxationBlock::inverse_svd; - using PreconditionBlockBase::log_statistics; + using PreconditionBlockBase::log_statistics; /** * Perform one step of the SSOR iteration. */ - template - void step (VECTOR &dst, const VECTOR &rhs) const; + void step (VectorType &dst, const VectorType &rhs) const; /** * Perform one step of the transposed SSOR iteration. */ - template - void Tstep (VECTOR &dst, const VECTOR &rhs) const; + void Tstep (VectorType &dst, const VectorType &rhs) const; }; diff --git a/include/deal.II/lac/relaxation_block.templates.h b/include/deal.II/lac/relaxation_block.templates.h index 70a99df7a7..84eafb84ae 100644 --- a/include/deal.II/lac/relaxation_block.templates.h +++ b/include/deal.II/lac/relaxation_block.templates.h @@ -19,28 +19,33 @@ #include #include #include +#include DEAL_II_NAMESPACE_OPEN -template +template inline -RelaxationBlock::AdditionalData::AdditionalData +RelaxationBlock::AdditionalData::AdditionalData (const double relaxation, const bool invert_diagonal, - const bool same_diagonal) + const bool same_diagonal, + const typename PreconditionBlockBase::Inversion inversion, + const double threshold, + VectorType *temp_ghost_vector) : relaxation(relaxation), invert_diagonal(invert_diagonal), same_diagonal(same_diagonal), - inversion(PreconditionBlockBase::gauss_jordan), - threshold(0.) + inversion(inversion), + threshold(threshold), + temp_ghost_vector (temp_ghost_vector) {} -template +template inline std::size_t -RelaxationBlock::AdditionalData::memory_consumption() const +RelaxationBlock::AdditionalData::memory_consumption() const { std::size_t result = sizeof(*this) + block_list.memory_consumption() @@ -51,11 +56,11 @@ RelaxationBlock::AdditionalData::memory_consumption() c } -template +template inline void -RelaxationBlock::initialize (const MatrixType &M, - const AdditionalData ¶meters) +RelaxationBlock::initialize (const MatrixType &M, + const AdditionalData ¶meters) { Assert (parameters.invert_diagonal, ExcNotImplemented()); @@ -73,24 +78,24 @@ RelaxationBlock::initialize (const MatrixType &M, } -template +template inline void -RelaxationBlock::clear () +RelaxationBlock::clear () { A = 0; additional_data = 0; - PreconditionBlockBase::clear (); + PreconditionBlockBase::clear (); } -template +template inline void -RelaxationBlock::invert_diagblocks () +RelaxationBlock::invert_diagblocks () { const MatrixType &M=*A; - FullMatrix M_cell; + FullMatrix M_cell; if (this->same_diagonal()) { @@ -131,14 +136,14 @@ RelaxationBlock::invert_diagblocks () } switch (this->inversion) { - case PreconditionBlockBase::gauss_jordan: + case PreconditionBlockBase::gauss_jordan: this->inverse(block).reinit(bs, bs); this->inverse(block).invert(M_cell); break; - case PreconditionBlockBase::householder: + case PreconditionBlockBase::householder: this->inverse_householder(block).initialize(M_cell); break; - case PreconditionBlockBase::svd: + case PreconditionBlockBase::svd: this->inverse_svd(block).reinit(bs, bs); this->inverse_svd(block) = M_cell; this->inverse_svd(block).compute_inverse_svd(additional_data->threshold); @@ -151,20 +156,59 @@ RelaxationBlock::invert_diagblocks () this->inverses_computed(true); } - -template -template +namespace internal +{ + /** + * Default implementation for serial vectors. Here we don't need to make a + * copy into a ghosted vector, so just return a reference to @p prev. + */ + template + const VectorType & + prepare_ghost_vector( + const VectorType &prev, + VectorType *other) + { + // If the following Assertion triggers, you either set temp_ghost_vector + // for a serial computation (don't!), or nobody implemented, instantiated, and + // tested the parallel version for your vector type. + Assert(other==NULL, ExcNotImplemented()); + (void)other; + return prev; + } + + /** + * Specialization for Trilinos. Use the ghosted vector. + */ + template <> + const TrilinosWrappers::MPI::Vector & + prepare_ghost_vector( + const TrilinosWrappers::MPI::Vector &prev, + TrilinosWrappers::MPI::Vector *other) + { + Assert(other!=NULL, + ExcMessage("You need to provide a ghosted vector in RelaxationBlock::AdditionalData::temp_trilinos_ghost_vector.")); + Assert(other->size()==prev.size(), ExcInternalError()); + + // import ghost values: + *other = prev; + return *other; + } +} // end namespace internal + +template inline void -RelaxationBlock::do_step (VECTOR &dst, - const VECTOR &prev, - const VECTOR &src, - const bool backward) const +RelaxationBlock::do_step (VectorType &dst, + const VectorType &prev, + const VectorType &src, + const bool backward) const { Assert (additional_data->invert_diagonal, ExcNotImplemented()); + const VectorType &ghosted_prev = internal::prepare_ghost_vector(prev, additional_data->temp_ghost_vector); + const MatrixType &M=*this->A; - Vector b_cell, x_cell; + Vector b_cell, x_cell; const bool permutation_empty = additional_data->order.size() == 0; const unsigned int n_permutations = (permutation_empty) @@ -197,7 +241,7 @@ RelaxationBlock::do_step (VECTOR &dst, b_cell(row_cell) = src(row->column()); for (typename MatrixType::const_iterator entry = M.begin(row->column()); entry != M.end(row->column()); ++entry) - b_cell(row_cell) -= entry->value() * prev(entry->column()); + b_cell(row_cell) -= entry->value() * ghosted_prev(entry->column()); } // Apply inverse diagonal this->inverse_vmult(block, x_cell, b_cell); @@ -208,7 +252,7 @@ RelaxationBlock::do_step (VECTOR &dst, } #endif // Store in result vector - row=additional_data->block_list.begin(block); + row = additional_data->block_list.begin(block); for (size_type row_cell=0; row_cellcolumn()) += additional_data->relaxation * x_cell(row_cell); } @@ -219,28 +263,26 @@ RelaxationBlock::do_step (VECTOR &dst, //----------------------------------------------------------------------// -template -template -void RelaxationBlockJacobi::step -(VECTOR &dst, - const VECTOR &src) const +template +void RelaxationBlockJacobi::step +(VectorType &dst, + const VectorType &src) const { - GrowingVectorMemory mem; - typename VectorMemory::Pointer aux = mem; + GrowingVectorMemory mem; + typename VectorMemory::Pointer aux = mem; aux->reinit(dst, false); *aux = dst; this->do_step(dst, *aux, src, false); } -template -template -void RelaxationBlockJacobi::Tstep -(VECTOR &dst, - const VECTOR &src) const +template +void RelaxationBlockJacobi::Tstep +(VectorType &dst, + const VectorType &src) const { - GrowingVectorMemory mem; - typename VectorMemory::Pointer aux = mem; + GrowingVectorMemory mem; + typename VectorMemory::Pointer aux = mem; aux->reinit(dst, false); *aux = dst; this->do_step(dst, *aux, src, true); @@ -249,21 +291,19 @@ void RelaxationBlockJacobi::Tstep //----------------------------------------------------------------------// -template -template -void RelaxationBlockSOR::step -(VECTOR &dst, - const VECTOR &src) const +template +void RelaxationBlockSOR::step +(VectorType &dst, + const VectorType &src) const { this->do_step(dst, dst, src, false); } -template -template -void RelaxationBlockSOR::Tstep -(VECTOR &dst, - const VECTOR &src) const +template +void RelaxationBlockSOR::Tstep +(VectorType &dst, + const VectorType &src) const { this->do_step(dst, dst, src, true); } @@ -271,22 +311,20 @@ void RelaxationBlockSOR::Tstep //----------------------------------------------------------------------// -template -template -void RelaxationBlockSSOR::step -(VECTOR &dst, - const VECTOR &src) const +template +void RelaxationBlockSSOR::step +(VectorType &dst, + const VectorType &src) const { this->do_step(dst, dst, src, false); this->do_step(dst, dst, src, true); } -template -template -void RelaxationBlockSSOR::Tstep -(VECTOR &dst, - const VECTOR &src) const +template +void RelaxationBlockSSOR::Tstep +(VectorType &dst, + const VectorType &src) const { this->do_step(dst, dst, src, true); this->do_step(dst, dst, src, false); diff --git a/source/lac/relaxation_block.inst.in b/source/lac/relaxation_block.inst.in index 47cf6b86fb..2cc7accd18 100644 --- a/source/lac/relaxation_block.inst.in +++ b/source/lac/relaxation_block.inst.in @@ -15,66 +15,22 @@ -for (S1, S2 : REAL_SCALARS) +for (S1, S2, S3 : REAL_SCALARS) { - template class RelaxationBlock, S2>; - template class RelaxationBlockJacobi, S2>; - template class RelaxationBlockSOR, S2>; - template class RelaxationBlockSSOR, S2>; + template class RelaxationBlock, S2, Vector >; + template class RelaxationBlockJacobi, S2, Vector >; + template class RelaxationBlockSOR, S2, Vector >; + template class RelaxationBlockSSOR, S2, Vector >; } -for (S1, S2, S3: REAL_SCALARS) -{ -// ------------ RelaxationBlockJacobi ----------------- - template - void RelaxationBlockJacobi, S2>::step - (Vector &, const Vector &) const; - template - void RelaxationBlockJacobi, S2>::Tstep - (Vector &, const Vector &) const; - -// ------------ RelaxationBlockSOR ----------------- - template - void RelaxationBlockSOR, S2>::step - (Vector &, const Vector &) const; - template - void RelaxationBlockSOR, S2>::Tstep - (Vector &, const Vector &) const; - -// ------------ RelaxationBlockSSOR ----------------- - template - void RelaxationBlockSSOR, S2>::step - (Vector &, const Vector &) const; - template - void RelaxationBlockSSOR, S2>::Tstep - (Vector &, const Vector &) const; -} - - -for (S2 : REAL_SCALARS) +for (S1 : REAL_SCALARS) { #ifdef DEAL_II_WITH_TRILINOS - template class RelaxationBlock; - template class RelaxationBlockJacobi; - template class RelaxationBlockSOR; - template class RelaxationBlockSSOR; - - template void RelaxationBlockJacobi::step - (TrilinosWrappers::MPI::Vector &, const TrilinosWrappers::MPI::Vector &) const; - template void RelaxationBlockJacobi::Tstep - (TrilinosWrappers::MPI::Vector &, const TrilinosWrappers::MPI::Vector &) const; - - template void RelaxationBlockSOR::step - (TrilinosWrappers::MPI::Vector &, const TrilinosWrappers::MPI::Vector &) const; - template void RelaxationBlockSOR::Tstep - (TrilinosWrappers::MPI::Vector &, const TrilinosWrappers::MPI::Vector &) const; - -// ------------ RelaxationBlockSSOR ----------------- - template void RelaxationBlockSSOR::step - (TrilinosWrappers::MPI::Vector &, const TrilinosWrappers::MPI::Vector &) const; - template void RelaxationBlockSSOR::Tstep - (TrilinosWrappers::MPI::Vector &, const TrilinosWrappers::MPI::Vector &) const; + template class RelaxationBlock; + template class RelaxationBlockJacobi; + template class RelaxationBlockSOR; + template class RelaxationBlockSSOR; #endif } diff --git a/tests/mpi/step-39-block.cc b/tests/mpi/step-39-block.cc index e24b17e9a3..c82f3ae28b 100644 --- a/tests/mpi/step-39-block.cc +++ b/tests/mpi/step-39-block.cc @@ -14,7 +14,7 @@ // --------------------------------------------------------------------- /* - * Author: Guido Kanschat, Texas A&M University, 2009 + * Test RelaxationBlockJacobi in parallel */ #include "../tests.h" @@ -565,24 +565,35 @@ namespace Step39 MGCoarseGridLACIteration,TrilinosWrappers::MPI::Vector> coarse_grid_solver(coarse_solver, coarse_matrix, identity); - typedef RelaxationBlockJacobi + typedef RelaxationBlockJacobi Smoother; MGLevelObject smoother_data; smoother_data.resize(0, triangulation.n_levels() - 1); mg::SmootherRelaxation mg_smoother; + MGLevelObject temp_vectors(0, triangulation.n_levels() - 1); + + for (unsigned int l = smoother_data.min_level() + 1; l <= smoother_data.max_level(); ++l) { DoFTools::make_cell_patches(smoother_data[l].block_list, this->dof_handler, l); - smoother_data[l].block_list.compress(); + if (smoother_data[l].block_list.n_rows()>0) + smoother_data[l].block_list.compress(); smoother_data[l].relaxation = 0.7; smoother_data[l].inversion = PreconditionBlockBase::svd; + TrilinosWrappers::MPI::Vector *ghost = &(temp_vectors[l]); + IndexSet relevant_dofs; + DoFTools::extract_locally_relevant_level_dofs(dof_handler, l, relevant_dofs); + ghost->reinit(dof_handler.locally_owned_mg_dofs(l), + relevant_dofs, + MPI_COMM_WORLD); + smoother_data[l].temp_ghost_vector = ghost; } mg_smoother.initialize(mg_matrix, smoother_data); - mg_smoother.set_steps(1); - mg_smoother.set_variable(true); + mg_smoother.set_steps(2); + mg_smoother.set_variable(false); mg::Matrix mgmatrix(mg_matrix); mg::Matrix mgdown(mg_matrix_dg_down); @@ -712,35 +723,8 @@ int main(int argc, char *argv[]) Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads()); MPILogInitAll log; - try - { - FE_DGQ<2> fe1(2); - InteriorPenaltyProblem<2> test1(fe1); - test1.run(6); - } - catch (std::exception &exc) - { - std::cerr << std::endl << std::endl - << "----------------------------------------------------" - << std::endl; - std::cerr << "Exception on processing: " << std::endl - << exc.what() << std::endl - << "Aborting!" << std::endl - << "----------------------------------------------------" - << std::endl; - return 1; - } - catch (...) - { - std::cerr << std::endl << std::endl - << "----------------------------------------------------" - << std::endl; - std::cerr << "Unknown exception!" << std::endl - << "Aborting!" << std::endl - << "----------------------------------------------------" - << std::endl; - return 1; - } - + FE_DGQ<2> fe1(2); + InteriorPenaltyProblem<2> test1(fe1); + test1.run(6); return 0; } diff --git a/tests/mpi/step-39-block.with_trilinos=true.mpirun=3.output b/tests/mpi/step-39-block.with_trilinos=true.mpirun=3.output index 0519ecba6e..75cd43ea51 100644 --- a/tests/mpi/step-39-block.with_trilinos=true.mpirun=3.output +++ b/tests/mpi/step-39-block.with_trilinos=true.mpirun=3.output @@ -1 +1,188 @@ - \ No newline at end of file + +DEAL:0::Element: FE_DGQ<2>(2) +DEAL:0::Step 0 +DEAL:0::Triangulation 12 cells, 2 levels +DEAL:0::DoFHandler 108 dofs, level dofs 27 108 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 20.6576 +DEAL:0:cg::Convergence step 15 value 4.82452e-13 +DEAL:0::Estimate 1.10045 +DEAL:0::Step 1 +DEAL:0::Triangulation 18 cells, 3 levels +DEAL:0::DoFHandler 162 dofs, level dofs 27 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 21.1045 +DEAL:0:cg::Convergence step 15 value 9.35135e-13 +DEAL:0::Estimate 1.03046 +DEAL:0::Step 2 +DEAL:0::Triangulation 27 cells, 4 levels +DEAL:0::DoFHandler 243 dofs, level dofs 27 108 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 21.3244 +DEAL:0:cg::Convergence step 16 value 7.03391e-13 +DEAL:0::Estimate 0.836764 +DEAL:0::Step 3 +DEAL:0::Triangulation 39 cells, 5 levels +DEAL:0::DoFHandler 351 dofs, level dofs 27 108 144 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 21.4584 +DEAL:0:cg::Convergence step 17 value 2.88305e-13 +DEAL:0::Estimate 0.644552 +DEAL:0::Step 4 +DEAL:0::Triangulation 54 cells, 6 levels +DEAL:0::DoFHandler 486 dofs, level dofs 27 108 180 144 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 21.5437 +DEAL:0:cg::Convergence step 16 value 6.40351e-13 +DEAL:0::Estimate 0.474987 +DEAL:0::Step 5 +DEAL:0::Triangulation 84 cells, 7 levels +DEAL:0::DoFHandler 756 dofs, level dofs 27 108 324 216 144 108 72 +DEAL:0::Assemble matrix +DEAL:0::Assemble multilevel matrix +DEAL:0::Assemble right hand side +DEAL:0::Solve +DEAL:0:cg::Starting value 25.5435 +DEAL:0:cg::Convergence step 16 value 4.60775e-13 +DEAL:0::Estimate 0.356044 + +DEAL:1::Element: FE_DGQ<2>(2) +DEAL:1::Step 0 +DEAL:1::Triangulation 12 cells, 2 levels +DEAL:1::DoFHandler 108 dofs, level dofs 27 108 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 20.6576 +DEAL:1:cg::Convergence step 15 value 4.82452e-13 +DEAL:1::Estimate 1.10045 +DEAL:1::Step 1 +DEAL:1::Triangulation 18 cells, 3 levels +DEAL:1::DoFHandler 162 dofs, level dofs 27 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 21.1045 +DEAL:1:cg::Convergence step 15 value 9.35135e-13 +DEAL:1::Estimate 1.03046 +DEAL:1::Step 2 +DEAL:1::Triangulation 27 cells, 4 levels +DEAL:1::DoFHandler 243 dofs, level dofs 27 108 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 21.3244 +DEAL:1:cg::Convergence step 16 value 7.03391e-13 +DEAL:1::Estimate 0.836764 +DEAL:1::Step 3 +DEAL:1::Triangulation 39 cells, 5 levels +DEAL:1::DoFHandler 351 dofs, level dofs 27 108 144 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 21.4584 +DEAL:1:cg::Convergence step 17 value 2.88305e-13 +DEAL:1::Estimate 0.644552 +DEAL:1::Step 4 +DEAL:1::Triangulation 54 cells, 6 levels +DEAL:1::DoFHandler 486 dofs, level dofs 27 108 180 144 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 21.5437 +DEAL:1:cg::Convergence step 16 value 6.40351e-13 +DEAL:1::Estimate 0.474987 +DEAL:1::Step 5 +DEAL:1::Triangulation 84 cells, 7 levels +DEAL:1::DoFHandler 756 dofs, level dofs 27 108 324 216 144 108 72 +DEAL:1::Assemble matrix +DEAL:1::Assemble multilevel matrix +DEAL:1::Assemble right hand side +DEAL:1::Solve +DEAL:1:cg::Starting value 25.5435 +DEAL:1:cg::Convergence step 16 value 4.60775e-13 +DEAL:1::Estimate 0.356044 + + +DEAL:2::Element: FE_DGQ<2>(2) +DEAL:2::Step 0 +DEAL:2::Triangulation 12 cells, 2 levels +DEAL:2::DoFHandler 108 dofs, level dofs 27 108 +DEAL:2::Assemble matrix +DEAL:2::Assemble multilevel matrix +DEAL:2::Assemble right hand side +DEAL:2::Solve +DEAL:2:cg::Starting value 20.6576 +DEAL:2:cg::Convergence step 15 value 4.82452e-13 +DEAL:2::Estimate 1.10045 +DEAL:2::Step 1 +DEAL:2::Triangulation 18 cells, 3 levels +DEAL:2::DoFHandler 162 dofs, level dofs 27 108 72 +DEAL:2::Assemble matrix +DEAL:2::Assemble multilevel matrix +DEAL:2::Assemble right hand side +DEAL:2::Solve +DEAL:2:cg::Starting value 21.1045 +DEAL:2:cg::Convergence step 15 value 9.35135e-13 +DEAL:2::Estimate 1.03046 +DEAL:2::Step 2 +DEAL:2::Triangulation 27 cells, 4 levels +DEAL:2::DoFHandler 243 dofs, level dofs 27 108 108 72 +DEAL:2::Assemble matrix +DEAL:2::Assemble multilevel matrix +DEAL:2::Assemble right hand side +DEAL:2::Solve +DEAL:2:cg::Starting value 21.3244 +DEAL:2:cg::Convergence step 16 value 7.03391e-13 +DEAL:2::Estimate 0.836764 +DEAL:2::Step 3 +DEAL:2::Triangulation 39 cells, 5 levels +DEAL:2::DoFHandler 351 dofs, level dofs 27 108 144 108 72 +DEAL:2::Assemble matrix +DEAL:2::Assemble multilevel matrix +DEAL:2::Assemble right hand side +DEAL:2::Solve +DEAL:2:cg::Starting value 21.4584 +DEAL:2:cg::Convergence step 17 value 2.88305e-13 +DEAL:2::Estimate 0.644552 +DEAL:2::Step 4 +DEAL:2::Triangulation 54 cells, 6 levels +DEAL:2::DoFHandler 486 dofs, level dofs 27 108 180 144 108 72 +DEAL:2::Assemble matrix +DEAL:2::Assemble multilevel matrix +DEAL:2::Assemble right hand side +DEAL:2::Solve +DEAL:2:cg::Starting value 21.5437 +DEAL:2:cg::Convergence step 16 value 6.40351e-13 +DEAL:2::Estimate 0.474987 +DEAL:2::Step 5 +DEAL:2::Triangulation 84 cells, 7 levels +DEAL:2::DoFHandler 756 dofs, level dofs 27 108 324 216 144 108 72 +DEAL:2::Assemble matrix +DEAL:2::Assemble multilevel matrix +DEAL:2::Assemble right hand side +DEAL:2::Solve +DEAL:2:cg::Starting value 25.5435 +DEAL:2:cg::Convergence step 16 value 4.60775e-13 +DEAL:2::Estimate 0.356044 +