From: janssen Date: Mon, 27 Sep 2010 10:13:50 +0000 (+0000) Subject: some more modifications to use dg and cg in multigrid at the same time X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bc61d7ca8655f30207e0d2ae8a1c616718ad16fa;p=dealii-svn.git some more modifications to use dg and cg in multigrid at the same time git-svn-id: https://svn.dealii.org/trunk@22169 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/multigrid/mg_constrained_dofs.h b/deal.II/deal.II/include/multigrid/mg_constrained_dofs.h index 42ec2f510f..b9d82d14b8 100644 --- a/deal.II/deal.II/include/multigrid/mg_constrained_dofs.h +++ b/deal.II/deal.II/include/multigrid/mg_constrained_dofs.h @@ -125,6 +125,7 @@ class MGConstrainedDoFs : public Subscriptor */ bool set_boundary_values () const; + bool continuity_across_edges () const; private: /** @@ -278,6 +279,20 @@ MGConstrainedDoFs::set_boundary_values () const return boundary_values_need_to_be_set; } +inline +bool +MGConstrainedDoFs::continuity_across_edges () const +{ + bool is_continuous = false; + for(unsigned int l=0; l& M, const std::vector& i1, - const std::vector& i2); + const std::vector& i2, + const unsigned int level = -1); + /** + * Assemble a single matrix + * into a global matrix. + */ + + void assemble_down(MATRIX& G, + const FullMatrix& M, + const std::vector& i1, + const std::vector& i2, + const unsigned int level = -1); /** * Assemble a single matrix * into a global matrix. */ + void assemble_in(MATRIX& G, const FullMatrix& M, const std::vector& i1, @@ -663,6 +676,7 @@ namespace MeshWorker * Assemble a single matrix * into a global matrix. */ + void assemble_out(MATRIX& G, const FullMatrix& M, const std::vector& i1, @@ -1002,14 +1016,14 @@ namespace MeshWorker * the fine and the coarse * level at refinement edges. */ - std::vector interface_down; + std::vector interface_out; /** * The interface matrix between * the coarse and the fine * level at refinement edges. */ - std::vector interface_up; + std::vector interface_in; /** * A pointer to the object containing the block structure. @@ -1562,19 +1576,72 @@ namespace MeshWorker template inline void - MGMatrixSimple::assemble_transpose( + MGMatrixSimple::assemble_up( MATRIX& G, const FullMatrix& M, const std::vector& i1, - const std::vector& i2) + const std::vector& i2, + const unsigned int level) { AssertDimension(M.n(), i1.size()); AssertDimension(M.m(), i2.size()); - for (unsigned int j=0; j= threshold) - G.add(i1[j], i2[k], M(k,j)); + if(mg_constrained_dofs == 0) + { + for (unsigned int j=0; j= threshold) + G.add(i1[j], i2[k], M(k,j)); + } + else + { + for (unsigned int j=0; j= threshold) + { + if(!mg_constrained_dofs->at_refinement_edge(level, i1[j]) || + mg_constrained_dofs->at_refinement_edge(level, i2[k])) + { + if(!mg_constrained_dofs->continuity_across_edges()) + G.add(i1[j], i2[k], M(k,j)); + } + } + } + } + + template + inline void + MGMatrixSimple::assemble_down( + MATRIX& G, + const FullMatrix& M, + const std::vector& i1, + const std::vector& i2, + const unsigned int level) + { + AssertDimension(M.m(), i1.size()); + AssertDimension(M.n(), i2.size()); + + if(mg_constrained_dofs == 0) + { + for (unsigned int j=0; j= threshold) + G.add(i1[j], i2[k], M(j,k)); + } + else + { + for (unsigned int j=0; j= threshold) + { + if(!mg_constrained_dofs->at_refinement_edge(level, i1[j]) || + mg_constrained_dofs->at_refinement_edge(level, i2[k])) + { + if(!mg_constrained_dofs->continuity_across_edges()) + G.add(i1[j], i2[k], M(j,k)); + } + } + } } template @@ -1615,7 +1682,7 @@ namespace MeshWorker G.add(i1[j], i2[k], M(j,k)); } else - G.add(i1[j], i2[k], M(j,k)); + G.add(i1[j], i2[k], M(j,k)); } } } @@ -1671,10 +1738,12 @@ namespace MeshWorker { const unsigned int level = info.cell->level(); assemble((*matrix)[level], info.matrix(0,false).matrix, info.indices, info.indices, level); + if(mg_constrained_dofs != 0) + //if(interface_in != 0 && interface_out != 0 && level>0) { - assemble_out((*interface_out)[level],info.matrix(0,false).matrix, info.indices, info.indices, level); assemble_in((*interface_in)[level], info.matrix(0,false).matrix, info.indices, info.indices, level); + assemble_out((*interface_out)[level],info.matrix(0,false).matrix, info.indices, info.indices, level); } } @@ -1712,8 +1781,13 @@ namespace MeshWorker // which is done by // the coarser cell assemble((*matrix)[level1], info1.matrix(0,false).matrix, info1.indices, info1.indices); - assemble_transpose((*flux_up)[level1],info1.matrix(0,true).matrix, info2.indices, info1.indices); - assemble((*flux_down)[level1], info2.matrix(0,true).matrix, info2.indices, info1.indices); + //assemble_transpose((*flux_up)[level1],info1.matrix(0,true).matrix, info2.indices, info1.indices); + //assemble((*flux_down)[level1], info2.matrix(0,true).matrix, info2.indices, info1.indices); + if(level1>0) + { + assemble_up((*flux_up)[level1],info1.matrix(0,true).matrix, info2.indices, info1.indices, level1); + assemble_down((*flux_down)[level1], info2.matrix(0,true).matrix, info2.indices, info1.indices, level1); + } } } diff --git a/deal.II/deal.II/source/multigrid/mg_tools.cc b/deal.II/deal.II/source/multigrid/mg_tools.cc index 5a7ea92d23..a0706b1a71 100644 --- a/deal.II/deal.II/source/multigrid/mg_tools.cc +++ b/deal.II/deal.II/source/multigrid/mg_tools.cc @@ -1495,6 +1495,9 @@ extract_inner_interface_dofs (const MGDoFHandler &mg_dof_handler, std::fill (cell_dofs.begin(), cell_dofs.end(), false); std::fill (boundary_cell_dofs.begin(), boundary_cell_dofs.end(), false); + if(fe.conforms(FiniteElementData::H1)) + { + deallog << "conforms H1" << std::endl; for (unsigned int face_nr=0; face_nr::faces_per_cell; ++face_nr) { const typename DoFHandler::face_iterator face = cell->face(face_nr); @@ -1509,7 +1512,7 @@ extract_inner_interface_dofs (const MGDoFHandler &mg_dof_handler, if (neighbor->level() < cell->level()) { for (unsigned int j=0; j &mg_dof_handler, boundary_interface_dofs[level][local_dof_indices[i]] = true; } } + } } #endif