From 328b7ecfc09447d4433875552dd82b2ef606106b Mon Sep 17 00:00:00 2001 From: guido Date: Thu, 22 Aug 2002 15:14:55 +0000 Subject: [PATCH] preparation for local multigrid git-svn-id: https://svn.dealii.org/trunk@6338 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/multigrid/mg_base.h | 78 ++++++++++ .../deal.II/include/multigrid/mg_dof_tools.h | 9 +- .../include/multigrid/mg_tools.templates.h | 2 +- deal.II/deal.II/include/multigrid/multigrid.h | 143 ++---------------- .../include/multigrid/multigrid.templates.h | 143 ++++++++++++++++-- .../deal.II/source/multigrid/mg_dof_tools.cc | 4 +- .../multigrid/multigrid.all_dimensions.cc | 6 + 7 files changed, 236 insertions(+), 149 deletions(-) diff --git a/deal.II/deal.II/include/multigrid/mg_base.h b/deal.II/deal.II/include/multigrid/mg_base.h index c858ae747e..968e55b585 100644 --- a/deal.II/deal.II/include/multigrid/mg_base.h +++ b/deal.II/deal.II/include/multigrid/mg_base.h @@ -130,6 +130,29 @@ class MGMatrixBase : public Subscriptor */ virtual void vmult(unsigned int level, VECTOR& dst, const VECTOR& src) const = 0; + + /** + * Adding matrix-vector-multiplication on + * a certain level. + */ + virtual void vmult_add(unsigned int level, VECTOR& dst, + const VECTOR& src) const = 0; + + /** + * Transpose + * matrix-vector-multiplication on + * a certain level. + */ + virtual void Tvmult(unsigned int level, VECTOR& dst, + const VECTOR& src) const = 0; + + /** + * Adding transpose + * matrix-vector-multiplication on + * a certain level. + */ + virtual void Tvmult_add(unsigned int level, VECTOR& dst, + const VECTOR& src) const = 0; }; @@ -159,6 +182,28 @@ class MGMatrix : public MGMatrixBase, virtual void vmult(unsigned int level, VECTOR& dst, const VECTOR& src) const; + /** + * Adding matrix-vector-multiplication on + * a certain level. + */ + virtual void vmult_add(unsigned int level, VECTOR& dst, + const VECTOR& src) const; + + /** + * Transpose + * matrix-vector-multiplication on + * a certain level. + */ + virtual void Tvmult(unsigned int level, VECTOR& dst, + const VECTOR& src) const; + + /** + * Adding transpose + * matrix-vector-multiplication on + * a certain level. + */ + virtual void Tvmult_add(unsigned int level, VECTOR& dst, + const VECTOR& src) const; }; @@ -373,6 +418,39 @@ MGMatrix::vmult (unsigned int level, } +template +void +MGMatrix::vmult_add (unsigned int level, + VECTOR& dst, + const VECTOR& src) const +{ + const MGLevelObject& m = **this; + m[level].vmult_add(dst, src); +} + + +template +void +MGMatrix::Tvmult (unsigned int level, + VECTOR& dst, + const VECTOR& src) const +{ + const MGLevelObject& m = **this; + m[level].Tvmult(dst, src); +} + + +template +void +MGMatrix::Tvmult_add (unsigned int level, + VECTOR& dst, + const VECTOR& src) const +{ + const MGLevelObject& m = **this; + m[level].Tvmult_add(dst, src); +} + + /*---------------------------- mgbase.h ---------------------------*/ #endif diff --git a/deal.II/deal.II/include/multigrid/mg_dof_tools.h b/deal.II/deal.II/include/multigrid/mg_dof_tools.h index a7d2233fb4..54ea5318aa 100644 --- a/deal.II/deal.II/include/multigrid/mg_dof_tools.h +++ b/deal.II/deal.II/include/multigrid/mg_dof_tools.h @@ -62,7 +62,7 @@ class MGTools * Make a sparsity pattern including fluxes * of discontinuous Galerkin methods. * @ref make_sparsity_pattern - * $see DoFTools + * @ref DoFTools */ template static void @@ -73,9 +73,12 @@ class MGTools /** * Create sparsity pattern for * the fluxes at refinement - * edges. + * edges. The matrix maps a + * function of the fine level + * space @p{level} to the coarser + * space. + * * @ref{make_flux_sparsity_pattern} - * @ref{DoFTools} */ template static void diff --git a/deal.II/deal.II/include/multigrid/mg_tools.templates.h b/deal.II/deal.II/include/multigrid/mg_tools.templates.h index efeb9e2f15..abc271eb36 100644 --- a/deal.II/deal.II/include/multigrid/mg_tools.templates.h +++ b/deal.II/deal.II/include/multigrid/mg_tools.templates.h @@ -16,7 +16,7 @@ #include #include -#include +#include #include #include diff --git a/deal.II/deal.II/include/multigrid/multigrid.h b/deal.II/deal.II/include/multigrid/multigrid.h index 7a787fbc9d..ad3afd3643 100644 --- a/deal.II/deal.II/include/multigrid/multigrid.h +++ b/deal.II/deal.II/include/multigrid/multigrid.h @@ -195,7 +195,15 @@ class Multigrid : public Subscriptor */ SmartPointer > post_smooth; - + /** + * Edge matrix from fine to coarse. + */ + SmartPointer > edge_down; + + /** + * Transpose edge matrix from coarse to fine. + */ + SmartPointer > edge_up; /** * Exception. @@ -290,142 +298,9 @@ class PreconditionMG : public Subscriptor }; - - /* --------------------------- inline functions --------------------- */ -template -template -Multigrid::Multigrid (const MGDoFHandler& mg_dof_handler, - const MGMatrixBase& matrix, - const MGCoarseGrid& coarse, - const MGTransfer& transfer, - const MGSmoother& pre_smooth, - const MGSmoother& post_smooth, - const unsigned int min_level, - const unsigned int max_level) - : - minlevel(min_level), - maxlevel(std::min(mg_dof_handler.get_tria().n_levels()-1, - max_level)), - defect(minlevel,maxlevel), - solution(minlevel,maxlevel), - t(minlevel,maxlevel), - matrix(&matrix), - coarse(&coarse), - transfer(&transfer), - pre_smooth(&pre_smooth), - post_smooth(&post_smooth) -{ -}; - - -template -void -Multigrid::level_mgstep(const unsigned int level) -{ -#ifdef MG_DEBUG - char *name = new char[100]; - sprintf(name, "MG%d-defect",level); - print_vector(level, defect[level], name); -#endif - - solution[level] = 0.; - - if (level == minlevel) - { - (*coarse)(level, solution[level], defect[level]); -#ifdef MG_DEBUG - sprintf(name, "MG%d-solution",level); - print_vector(level, solution[level], name); -#endif - return; - } - - // smoothing of the residual by modifying s - pre_smooth->smooth(level, solution[level], defect[level]); - -#ifdef MG_DEBUG - sprintf(name, "MG%d-pre",level); - print_vector(level, solution[level], name); -#endif - - // t = -A*solution[level] - matrix->vmult(level, t[level], solution[level]); - t[level].scale(-1); - - // make t rhs of lower level - // The non-refined parts of the - // coarse-level defect already contain - // the global defect, the refined parts - // its restriction. - for (unsigned int l = level;l>minlevel;--l) - { - t[l-1] = 0.; - transfer->restrict_and_add (l, t[l-1], t[l]); - defect[l-1] += t[l-1]; - } - - // add additional DG contribution -// edge_vmult(level, defect[level-1], defect[level]); - - // do recursion - level_mgstep(level-1); - - // reset size of the auxiliary - // vector, since it has been - // resized in the recursive call to - // level_mgstep directly above - t[level] = 0.; - - // do coarse grid correction - - transfer->prolongate(level, t[level], solution[level-1]); - -#ifdef MG_DEBUG - sprintf(name, "MG%d-cgc",level); - print_vector(level, t, name); -#endif - - solution[level] += t[level]; - - // post-smoothing - post_smooth->smooth(level, solution[level], defect[level]); - -#ifdef MG_DEBUG - sprintf(name, "MG%d-post",level); - print_vector(level, solution[level], name); - - delete[] name; -#endif -} - - -template -void -Multigrid::vcycle() -{ - // The defect vector has been - // initialized by copy_to_mg. Now - // adjust the other vectors. - solution.resize(minlevel, maxlevel); - t.resize(minlevel, maxlevel); - - for (unsigned int level=minlevel; level<=maxlevel;++level) - { - solution[level].reinit(defect[level]); - t[level].reinit(defect[level]); - } - - level_mgstep (maxlevel); -// abort (); -} - - -/* ------------------------------------------------------------------- */ - - template PreconditionMG ::PreconditionMG(const MGDoFHandler& mg_dof_handler, diff --git a/deal.II/deal.II/include/multigrid/multigrid.templates.h b/deal.II/deal.II/include/multigrid/multigrid.templates.h index c23bf90d5d..dbd9e42384 100644 --- a/deal.II/deal.II/include/multigrid/multigrid.templates.h +++ b/deal.II/deal.II/include/multigrid/multigrid.templates.h @@ -14,21 +14,146 @@ #define __deal2__multigrid_templates_h -#include -//#include #include -#include -#include -#include +/* --------------------------- inline functions --------------------- */ +template template +Multigrid::Multigrid (const MGDoFHandler& mg_dof_handler, + const MGMatrixBase& matrix, + const MGCoarseGrid& coarse, + const MGTransfer& transfer, + const MGSmoother& pre_smooth, + const MGSmoother& post_smooth, + const unsigned int min_level, + const unsigned int max_level) + : + minlevel(min_level), + maxlevel(std::min(mg_dof_handler.get_tria().n_levels()-1, + max_level)), + defect(minlevel,maxlevel), + solution(minlevel,maxlevel), + t(minlevel,maxlevel), + matrix(&matrix), + coarse(&coarse), + transfer(&transfer), + pre_smooth(&pre_smooth), + post_smooth(&post_smooth) +{ +}; + + +template +void +Multigrid::level_mgstep(const unsigned int level) +{ +#ifdef MG_DEBUG + char *name = new char[100]; + sprintf(name, "MG%d-defect",level); + print_vector(level, defect[level], name); +#endif + + solution[level] = 0.; + + if (level == minlevel) + { + (*coarse)(level, solution[level], defect[level]); +#ifdef MG_DEBUG + sprintf(name, "MG%d-solution",level); + print_vector(level, solution[level], name); +#endif + return; + } + + // smoothing of the residual by modifying s + pre_smooth->smooth(level, solution[level], defect[level]); + +#ifdef MG_DEBUG + sprintf(name, "MG%d-pre",level); + print_vector(level, solution[level], name); +#endif + + // t = -A*solution[level] + matrix->vmult(level, t[level], solution[level]); + t[level].scale(-1); + + // make t rhs of lower level + // The non-refined parts of the + // coarse-level defect already contain + // the global defect, the refined parts + // its restriction. + for (unsigned int l = level;l>minlevel;--l) + { + t[l-1] = 0.; + transfer->restrict_and_add (l, t[l-1], t[l]); + defect[l-1] += t[l-1]; + } + + // add additional DG contribution +// edge_vmult(level, defect[level-1], defect[level]); + + // do recursion + level_mgstep(level-1); + + // reset size of the auxiliary + // vector, since it has been + // resized in the recursive call to + // level_mgstep directly above + t[level] = 0.; + + // do coarse grid correction + + transfer->prolongate(level, t[level], solution[level-1]); + +#ifdef MG_DEBUG + sprintf(name, "MG%d-cgc",level); + print_vector(level, t, name); +#endif + + solution[level] += t[level]; + + // post-smoothing + post_smooth->smooth(level, solution[level], defect[level]); + +#ifdef MG_DEBUG + sprintf(name, "MG%d-post",level); + print_vector(level, solution[level], name); + + delete[] name; +#endif +} + + +template void -Multigrid::print_vector (const unsigned int /*level*/, - const Vector& /*v*/, - const char* /*name*/) const +Multigrid::vcycle() { + // The defect vector has been + // initialized by copy_to_mg. Now + // adjust the other vectors. + solution.resize(minlevel, maxlevel); + t.resize(minlevel, maxlevel); + + for (unsigned int level=minlevel; level<=maxlevel;++level) + { + solution[level].reinit(defect[level]); + t[level].reinit(defect[level]); + } + + level_mgstep (maxlevel); +// abort (); +} + + + +// template +// void +// Multigrid::print_vector (const unsigned int /*level*/, +// const Vector& /*v*/, +// const char* /*name*/) const +// { // if (level!=maxlevel) // return; @@ -71,6 +196,6 @@ Multigrid::print_vector (const unsigned int /*level*/, // out.add_data_vector(out_vector, "v"); // out.build_patches(5); // out.write_gnuplot(out_file); -} +//} #endif diff --git a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc index c3bc1ba9da..80843d7cba 100644 --- a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc +++ b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc @@ -121,8 +121,8 @@ MGTools::make_flux_sparsity_pattern (const MGDoFHandler &dof, template void MGTools::make_flux_sparsity_pattern_edge (const MGDoFHandler &dof, - SparsityPattern &sparsity, - const unsigned int level) + SparsityPattern &sparsity, + const unsigned int level) { Assert ((level>=1) && (level #include #include +#include template @@ -135,6 +136,11 @@ void MGTransferSelect::restrict_and_add ( // Explicit instantiations +template class Multigrid >; +template class Multigrid >; +template class Multigrid >; +template class Multigrid >; + template class MGTransferPrebuilt; template class MGTransferPrebuilt; template class MGTransferBlock; -- 2.39.5