From 5b75a6ed586722ca05b0fce62978cb4953fe1829 Mon Sep 17 00:00:00 2001 From: hartmann Date: Mon, 3 Feb 2003 10:34:30 +0000 Subject: [PATCH] Make the print_vectors function usable when MG_DEBUG is defined. Implementation of print_vectors for Vector. git-svn-id: https://svn.dealii.org/trunk@7003 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/multigrid/multigrid.h | 22 +++++- .../include/multigrid/multigrid.templates.h | 76 +++++-------------- .../multigrid/multigrid.all_dimensions.cc | 70 +++++++++++++++++ 3 files changed, 110 insertions(+), 58 deletions(-) diff --git a/deal.II/deal.II/include/multigrid/multigrid.h b/deal.II/deal.II/include/multigrid/multigrid.h index 06368d3a47..fc924440ca 100644 --- a/deal.II/deal.II/include/multigrid/multigrid.h +++ b/deal.II/deal.II/include/multigrid/multigrid.h @@ -59,7 +59,7 @@ * returning the solution to the system @p{Ax=b} on the coarsest level * in @p{x}. * - * @author Guido Kanschat, 1999, 2001, 2002, Ralf Hartmann 2002 + * @author Guido Kanschat, 1999, 2001, 2002, Ralf Hartmann 2002, 2003. */ template class Multigrid : public Subscriptor @@ -137,6 +137,13 @@ class Multigrid : public Subscriptor */ void set_edge_matrices (const MGMatrixBase& edge_down, const MGMatrixBase& edge_up); + +#ifdef MG_DEBUG + void print_vector (const unsigned int level, + const VECTOR &v, + const char *name) const; +#endif + private: /** @@ -154,6 +161,7 @@ class Multigrid : public Subscriptor * this level. */ void level_mgstep (const unsigned int level); + /** * Level for coarse grid solution. */ @@ -220,6 +228,15 @@ class Multigrid : public Subscriptor */ SmartPointer > edge_up; + /** + * Pointer to the MGDoFHandler + * given to the constructor. Only + * needed for @p{MG_DEBUG} + * defined. + */ +#ifdef MG_DEBUG + SmartPointer > mg_dof_handler; +#endif /** * Exception. */ @@ -345,6 +362,9 @@ Multigrid::Multigrid (const MGDoFHandler& mg_dof_handler, post_smooth(&post_smooth), edge_down(0), edge_up(0) +#ifdef MG_DEBUG + , mg_dof_handler(&mg_dof_handler) +#endif {} diff --git a/deal.II/deal.II/include/multigrid/multigrid.templates.h b/deal.II/deal.II/include/multigrid/multigrid.templates.h index 60e00fddcf..e235a3a135 100644 --- a/deal.II/deal.II/include/multigrid/multigrid.templates.h +++ b/deal.II/deal.II/include/multigrid/multigrid.templates.h @@ -12,9 +12,19 @@ //---------------------------- multigrid.templates.h --------------------------- #ifndef __deal2__multigrid_templates_h #define __deal2__multigrid_templates_h +#include + +#ifdef MG_DEBUG +#include +#include +#include +#include +#include +#include +#endif + -#include template void @@ -50,7 +60,8 @@ Multigrid::level_mgstep(const unsigned int level) return; } - // smoothing of the residual by modifying s + // smoothing of the residual by + // modifying s pre_smooth->smooth(level, solution[level], defect[level]); #ifdef MG_DEBUG @@ -61,11 +72,11 @@ Multigrid::level_mgstep(const unsigned int level) // t = A*solution[level] matrix->vmult(level, t[level], solution[level]); - // 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. + // 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.; @@ -92,7 +103,7 @@ Multigrid::level_mgstep(const unsigned int level) #ifdef MG_DEBUG sprintf(name, "MG%d-cgc",level); - print_vector(level, t, name); + print_vector(level, t[level], name); #endif solution[level] += t[level]; @@ -137,54 +148,5 @@ Multigrid::vcycle() -// template -// void -// Multigrid::print_vector (const unsigned int /*level*/, -// const Vector& /*v*/, -// const char* /*name*/) const -// { -// if (level!=maxlevel) -// return; - -// const DoFHandler* dof = mg_dof_handler; - -// Vector out_vector(dof->n_dofs()); - -// out_vector = -10000; - -// const unsigned int dofs_per_cell = mg_dof_handler->get_fe().dofs_per_cell; - -// std::vector global_dof_indices (dofs_per_cell); -// std::vector level_dof_indices (dofs_per_cell); - -// typename DoFHandler::active_cell_iterator -// global_cell = dof->begin_active(level); -// typename MGDoFHandler::active_cell_iterator -// level_cell = mg_dof_handler->begin_active(level); -// const typename MGDoFHandler::cell_iterator -// endc = mg_dof_handler->end(level); - -// // traverse all cells and copy the -// // data appropriately to the output -// // vector -// for (; level_cell != endc; ++level_cell, ++global_cell) -// { -// global_cell->get_dof_indices (global_dof_indices); -// level_cell->get_mg_dof_indices(level_dof_indices); - -// // copy level-wise data to -// // global vector -// for (unsigned int i=0; i out; -// out.attach_dof_handler(*dof); -// 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/multigrid.all_dimensions.cc b/deal.II/deal.II/source/multigrid/multigrid.all_dimensions.cc index f734dfc6d5..b634c62f5e 100644 --- a/deal.II/deal.II/source/multigrid/multigrid.all_dimensions.cc +++ b/deal.II/deal.II/source/multigrid/multigrid.all_dimensions.cc @@ -20,6 +20,76 @@ #include + + + +#ifdef MG_DEBUG +template <> +void +Multigrid >::print_vector (const unsigned int level, + const Vector &v, + const char *name) const +{ + if (level!=maxlevel) + return; + const unsigned int dim=deal_II_dimension; + + const DoFHandler *dof = mg_dof_handler; + + Vector out_vector(dof->n_dofs()); + + out_vector = -10000; + + const unsigned int dofs_per_cell = mg_dof_handler->get_fe().dofs_per_cell; + + std::vector global_dof_indices (dofs_per_cell); + std::vector level_dof_indices (dofs_per_cell); + + DoFHandler::active_cell_iterator + global_cell = dof->begin_active(level); + MGDoFHandler::active_cell_iterator + level_cell = mg_dof_handler->begin_active(level); + const MGDoFHandler::cell_iterator + endc = mg_dof_handler->end(level); + + // traverse all cells and copy the + // data appropriately to the output + // vector + for (; level_cell != endc; ++level_cell, ++global_cell) + { + global_cell->get_dof_indices (global_dof_indices); + level_cell->get_mg_dof_indices(level_dof_indices); + + // copy level-wise data to + // global vector + for (unsigned int i=0; i out; + out.attach_dof_handler(*dof); + out.add_data_vector(out_vector, "v"); + out.build_patches(5); + out.write_gnuplot(out_file); +} + + + +template +void +Multigrid::print_vector (const unsigned int level, + const VECTOR &v, + const char *name) const +{ + Assert(false, ExcNotImplemented()); +} +#endif + + + + template MGTransferPrebuilt::~MGTransferPrebuilt () {} -- 2.39.5