From: janssen Date: Sun, 6 Nov 2011 15:48:08 +0000 (+0000) Subject: multigrid for continuous elements of systems X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=cc834032146c06acb72609324d86fc9f6af2ecb9;p=dealii-svn.git multigrid for continuous elements of systems git-svn-id: https://svn.dealii.org/trunk@24732 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/multigrid/mg_transfer_block.h b/deal.II/include/deal.II/multigrid/mg_transfer_block.h index f2906d6cc5..6715471ee8 100644 --- a/deal.II/include/deal.II/multigrid/mg_transfer_block.h +++ b/deal.II/include/deal.II/multigrid/mg_transfer_block.h @@ -26,6 +26,7 @@ #include #include +#include #include @@ -60,6 +61,19 @@ template class MGDoFHandler; class MGTransferBlockBase { public: + /** + * Constructor without constraint + * matrices. Use this constructor + * only with discontinuous finite + * elements or with no local + * refinement. + */ + MGTransferBlockBase (); + /** + * Constructor with constraint matrices as well as mg_constrained_dofs. + */ + MGTransferBlockBase (const ConstraintMatrix& constraints, + const MGConstrainedDoFs& mg_constrained_dofs); /** * Memory used by this object. */ @@ -159,6 +173,17 @@ class MGTransferBlockBase */ std::vector > > > copy_indices; + /** + * The constraints of the global + * system. + */ + SmartPointer constraints; + /** + * The mg_constrained_dofs of the level + * systems. + */ + + SmartPointer mg_constrained_dofs; }; /** @@ -315,6 +340,19 @@ class MGTransferBlockSelect : public MGTransferBase >, private MGTransferBlockBase { public: + /** + * Constructor without constraint + * matrices. Use this constructor + * only with discontinuous finite + * elements or with no local + * refinement. + */ + MGTransferBlockSelect (); + /** + * Constructor with constraint matrices as well as mg_constrained_dofs. + */ + MGTransferBlockSelect (const ConstraintMatrix& constraints, + const MGConstrainedDoFs& mg_constrained_dofs); /** * Destructor. */ diff --git a/deal.II/source/multigrid/mg_transfer_block.cc b/deal.II/source/multigrid/mg_transfer_block.cc index 1a6a3cdfdb..533cab4ab8 100644 --- a/deal.II/source/multigrid/mg_transfer_block.cc +++ b/deal.II/source/multigrid/mg_transfer_block.cc @@ -305,6 +305,7 @@ void MGTransferBlockBase::build_matrices ( // need to take care of cells on // the coarser level which have // children + for (unsigned int level=0; levelset_boundary_values()) + { + std::vector constrain_indices; + std::vector > constraints_per_block (n_blocks); + for (int level=n_levels-2; level>=0; --level) + { + if (mg_constrained_dofs->get_boundary_indices()[level].size() == 0) + continue; + + // need to delete all the columns in the + // matrix that are on the boundary. to achive + // this, create an array as long as there are + // matrix columns, and find which columns we + // need to filter away. + constrain_indices.resize (0); + constrain_indices.resize (prolongation_matrices[level]->n(), 0); + std::set::const_iterator dof + = mg_constrained_dofs->get_boundary_indices()[level].begin(), + endd = mg_constrained_dofs->get_boundary_indices()[level].end(); + for (; dof != endd; ++dof) + constrain_indices[*dof] = 1; + + unsigned int index = 0; + for(unsigned int block=0; blockblock(block, block).m(); + constraints_per_block[block].resize(0); + constraints_per_block[block].resize(n_dofs, 0); + for (unsigned int i=0; i::iterator + start_row = prolongation_matrices[level]->block(block, block).begin(i), + end_row = prolongation_matrices[level]->block(block, block).end(i); + for(; start_row != end_row; ++start_row) + { + if (constraints_per_block[block][start_row->column()]) + start_row->value() = 0.; + } + } + } + } } } @@ -427,6 +477,7 @@ void MGTransferBlockSelect::build_matrices ( std::vector temp_copy_indices; std::vector global_dof_indices (fe.dofs_per_cell); std::vector level_dof_indices (fe.dofs_per_cell); + for (int level=dof.get_tria().n_levels()-1; level>=0; --level) { typename MGDoFHandler::active_cell_iterator @@ -454,8 +505,17 @@ void MGTransferBlockSelect::build_matrices ( { const unsigned int block = fe.system_to_block_index(i).first; if (selected[block]) + { + if(mg_constrained_dofs != 0) + { + if(!mg_constrained_dofs->at_refinement_edge(level,level_dof_indices[i])) + temp_copy_indices[level_dof_indices[i] - mg_block_start[level][block]] + = global_dof_indices[i] - block_start[block]; + } + else temp_copy_indices[level_dof_indices[i] - mg_block_start[level][block]] = global_dof_indices[i] - block_start[block]; + } } } diff --git a/deal.II/source/multigrid/mg_transfer_prebuilt.cc b/deal.II/source/multigrid/mg_transfer_prebuilt.cc index 56de732bc7..2e658172cb 100644 --- a/deal.II/source/multigrid/mg_transfer_prebuilt.cc +++ b/deal.II/source/multigrid/mg_transfer_prebuilt.cc @@ -137,7 +137,8 @@ MGTransferPrebuilt::copy_to_mg ( // have fine level basis // functions if (!first) - restrict_and_add (level+1, dst[level], dst[level+1]); + restrict_and_add (level+1, dst[level], dst[level+1]); + first = false; } } diff --git a/deal.II/source/multigrid/multigrid.cc b/deal.II/source/multigrid/multigrid.cc index a7029bbaca..ed4803c321 100644 --- a/deal.II/source/multigrid/multigrid.cc +++ b/deal.II/source/multigrid/multigrid.cc @@ -145,6 +145,19 @@ void MGTransferPrebuilt::restrict_and_add ( prolongation_matrices[from_level-1]->Tvmult_add (dst, src); } + +MGTransferBlockBase::MGTransferBlockBase () +{} + + +MGTransferBlockBase::MGTransferBlockBase ( + const ConstraintMatrix &c, const MGConstrainedDoFs& mg_c) + : + constraints(&c), + mg_constrained_dofs(&mg_c) +{} + + template MGTransferBlock::MGTransferBlock () : @@ -333,6 +346,17 @@ void MGTransferSelect::restrict_and_add ( //----------------------------------------------------------------------// +template +MGTransferBlockSelect::MGTransferBlockSelect () +{} + + +template +MGTransferBlockSelect::MGTransferBlockSelect ( + const ConstraintMatrix &c, const MGConstrainedDoFs& mg_c) +: MGTransferBlockBase(c, mg_c) +{} + template MGTransferBlockSelect::~MGTransferBlockSelect () {}