From c19cd2c1485fa24811ac72e708678f36a8ff3491 Mon Sep 17 00:00:00 2001 From: "denis.davydov" Date: Thu, 20 Mar 2014 16:57:24 +0000 Subject: [PATCH] introduced shared parallel triangulation git-svn-id: https://svn.dealii.org/branches/branch_sharedtria@32674 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/distributed/tria.h | 41 +++++++ .../deal.II/dofs/dof_accessor.templates.h | 6 +- .../deal.II/grid/tria_accessor.templates.h | 35 ++++++ deal.II/source/distributed/tria.cc | 103 ++++++++++++++++++ deal.II/source/distributed/tria.inst.in | 11 ++ deal.II/source/dofs/dof_tools.cc | 7 +- 6 files changed, 198 insertions(+), 5 deletions(-) diff --git a/deal.II/include/deal.II/distributed/tria.h b/deal.II/include/deal.II/distributed/tria.h index b00a01b224..2e3bd67ed6 100644 --- a/deal.II/include/deal.II/distributed/tria.h +++ b/deal.II/include/deal.II/distributed/tria.h @@ -931,6 +931,47 @@ namespace parallel #endif +namespace parallel +{ + namespace shared + { + template + class Triangulation : public dealii::Triangulation + { + public: + /** + * Constructor. + */ + Triangulation (MPI_Comm mpi_communicator); + + /** + * Destructor. + */ + virtual ~Triangulation (); + + types::subdomain_id locally_owned_subdomain () const; + + /** + * Return the MPI communicator used by this triangulation. + */ + MPI_Comm get_communicator () const; + + virtual void execute_coarsening_and_refinement (); + + virtual void create_triangulation (const std::vector< Point< spacedim > > &vertices, + const std::vector< CellData< dim > > &cells, + const SubCellData &subcelldata); + + + private: + MPI_Comm mpi_communicator; + types::subdomain_id my_subdomain; + types::subdomain_id num_subdomains; + + void mark_artificial(); + }; + } +} DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/include/deal.II/dofs/dof_accessor.templates.h b/deal.II/include/deal.II/dofs/dof_accessor.templates.h index 42446a05c0..cd33e62451 100644 --- a/deal.II/include/deal.II/dofs/dof_accessor.templates.h +++ b/deal.II/include/deal.II/dofs/dof_accessor.templates.h @@ -3261,8 +3261,10 @@ void DoFCellAccessor::get_dof_indices (std::vector &dof_indices) const { Assert (this->active(), ExcMessage ("get_dof_indices() only works on active cells.")); - Assert (this->is_artificial() == false, - ExcMessage ("Can't ask for DoF indices on artificial cells.")); + + Assert ( (dynamic_cast*>(this->tria)==0) + ||(this->is_artificial() == false), + ExcMessage ("Can't ask for DoF indices on artificial cells of distributed triangulation.")); AssertDimension (dof_indices.size(), this->get_fe().dofs_per_cell); const types::global_dof_index *cache diff --git a/deal.II/include/deal.II/grid/tria_accessor.templates.h b/deal.II/include/deal.II/grid/tria_accessor.templates.h index 19bb7ad61a..5f51543613 100644 --- a/deal.II/include/deal.II/grid/tria_accessor.templates.h +++ b/deal.II/include/deal.II/grid/tria_accessor.templates.h @@ -39,6 +39,11 @@ namespace parallel { template class Triangulation; } + + namespace shared + { + template class Triangulation; + } } @@ -2877,6 +2882,13 @@ bool CellAccessor::is_locally_owned () const { Assert (this->active(), ExcMessage("is_locally_owned() only works on active cells!")); + + const parallel::shared::Triangulation *pst + = dynamic_cast *>(this->tria); + + if (pst != 0) + return (this->subdomain_id() == pst->locally_owned_subdomain()); + #ifndef DEAL_II_WITH_P4EST return true; #else @@ -2899,6 +2911,12 @@ inline bool CellAccessor::is_locally_owned_on_level () const { + const parallel::shared::Triangulation *pst + = dynamic_cast *>(this->tria); + + if (pst != 0) + return (this->level_subdomain_id() == pst->locally_owned_subdomain()); + #ifndef DEAL_II_WITH_P4EST return true; #else @@ -2919,6 +2937,16 @@ bool CellAccessor::is_ghost () const { Assert (this->active(), ExcMessage("is_ghost() only works on active cells!")); + + if (is_artificial() || this->has_children()) + return false; + + const parallel::shared::Triangulation *pst + = dynamic_cast *>(this->tria); + + if (pst != 0) + return (this->subdomain_id() != pst->locally_owned_subdomain()); + #ifndef DEAL_II_WITH_P4EST return false; #else @@ -2943,6 +2971,13 @@ bool CellAccessor::is_artificial () const { Assert (this->active(), ExcMessage("is_artificial() only works on active cells!")); + + const parallel::shared::Triangulation *pst + = dynamic_cast *>(this->tria); + + if (pst != 0) + return this->subdomain_id() == numbers::artificial_subdomain_id;//could probably keep this line only? + #ifndef DEAL_II_WITH_P4EST return false; #else diff --git a/deal.II/source/distributed/tria.cc b/deal.II/source/distributed/tria.cc index 4085903ce0..ef29d175da 100644 --- a/deal.II/source/distributed/tria.cc +++ b/deal.II/source/distributed/tria.cc @@ -4115,6 +4115,109 @@ namespace parallel #endif // DEAL_II_WITH_P4EST +namespace parallel +{ + namespace shared + { + template + Triangulation::Triangulation (MPI_Comm mpi_communicator): + dealii::Triangulation(), + mpi_communicator (Utilities::MPI:: + duplicate_communicator(mpi_communicator)), + my_subdomain (Utilities::MPI::this_mpi_process (this->mpi_communicator)), + num_subdomains(Utilities::MPI::n_mpi_processes(mpi_communicator)) + { + + } + + + template + Triangulation::~Triangulation () + { + + } + + + + template + types::subdomain_id + Triangulation::locally_owned_subdomain () const + { + return my_subdomain; + } + + template + void + Triangulation::execute_coarsening_and_refinement () { + dealii::Triangulation::execute_coarsening_and_refinement (); + dealii::GridTools::partition_triangulation (num_subdomains, *this); + mark_artificial(); + } + + template + void + Triangulation::create_triangulation (const std::vector< Point< spacedim > > &vertices, + const std::vector< CellData< dim > > &cells, + const SubCellData &subcelldata) { + try + { + dealii::Triangulation:: + create_triangulation (vertices, cells, subcelldata); + } + catch (const typename dealii::Triangulation::DistortedCellList &) + { + // the underlying triangulation should not be checking for distorted + // cells + AssertThrow (false, ExcInternalError()); + } + dealii::GridTools::partition_triangulation (num_subdomains, *this); + mark_artificial(); + } + + + template + MPI_Comm + Triangulation::get_communicator () const + { + return mpi_communicator; + } + + template + void + Triangulation::mark_artificial() + { + + std::vector marked_vertices(this->n_vertices(),false); + for (typename dealii::Triangulation::active_cell_iterator + cell = this->begin_active(); + cell != this->end(); ++cell) + if (cell->subdomain_id() == this->locally_owned_subdomain()) + for (unsigned int v=0; v::vertices_per_cell; ++v) + marked_vertices[cell->vertex_index(v)] = true; + + + + for (typename dealii::Triangulation::active_cell_iterator + cell = this->begin_active(); + cell != this->end(); ++cell) + if (cell->subdomain_id() != this->locally_owned_subdomain()) { + bool is_ghost = false; + for (unsigned int v=0; v::vertices_per_cell; ++v) + if (marked_vertices[cell->vertex_index(v)] == true) + { + is_ghost = true; + break; + } + if (is_ghost) + continue; + + cell->set_subdomain_id(numbers::artificial_subdomain_id); + } + } + + + } +} /*-------------- Explicit Instantiations -------------------------------*/ diff --git a/deal.II/source/distributed/tria.inst.in b/deal.II/source/distributed/tria.inst.in index b5199699f3..3b80440cb5 100644 --- a/deal.II/source/distributed/tria.inst.in +++ b/deal.II/source/distributed/tria.inst.in @@ -61,6 +61,17 @@ for (deal_II_dimension : DIMENSIONS) # if deal_II_dimension < 3 template class Triangulation; # endif +# if deal_II_dimension < 2 + template class Triangulation; +# endif + \} + + namespace shared + \{ + template class Triangulation; +# if deal_II_dimension < 3 + template class Triangulation; +# endif # if deal_II_dimension < 2 template class Triangulation; # endif diff --git a/deal.II/source/dofs/dof_tools.cc b/deal.II/source/dofs/dof_tools.cc index 3ea5499334..1402055720 100644 --- a/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/source/dofs/dof_tools.cc @@ -1096,9 +1096,10 @@ namespace DoFTools endc = dof_handler.end(); for (; cell!=endc; ++cell) { - Assert (cell->is_artificial() == false, - ExcMessage ("You can't call this function for meshes that " - "have artificial cells.")); + //can comment it since the case of distributed triangulations is handled above. + //Assert (cell->is_artificial() == false, + // ExcMessage ("You can't call this function for meshes that " + // "have artificial cells.")); const types::subdomain_id subdomain_id = cell->subdomain_id(); const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell; -- 2.39.5