From: Daniel Arndt Date: Wed, 28 Sep 2016 16:29:44 +0000 (+0200) Subject: Hide parallel implementation and use guards for vtk test output X-Git-Tag: v8.5.0-rc1~615^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=651e24ac9f714438cd8b292c2bdcc807db106be9;p=dealii.git Hide parallel implementation and use guards for vtk test output --- diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 0e7dab6dee..203013761a 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -783,28 +783,6 @@ namespace FETools const ConstraintMatrix &constraints, OutVector &z2); - /** - * Same as above for external parallel VectorTypes - * on a parallel::distributed::Triangulation - */ - template - void extrapolate_parallel(const DoFHandler &dof1, - const VectorType &u1, - const DoFHandler &dof2, - const ConstraintMatrix &constraints2, - VectorType &u2); - - /** - * Same as above for deal.II parallel Vectors - * on a parallel::distributed::Triangulation - */ - template - void extrapolate_parallel(const DoFHandler &dof1, - const LinearAlgebra::distributed::Vector &u1, - const DoFHandler &dof2, - const ConstraintMatrix &constraints2, - LinearAlgebra::distributed::Vector &u2); - //@} /** * The numbering of the degrees of freedom in continuous finite elements is diff --git a/source/fe/fe_tools_extrapolate.cc b/source/fe/fe_tools_extrapolate.cc index 9dae497d7c..13e85c8ca2 100755 --- a/source/fe/fe_tools_extrapolate.cc +++ b/source/fe/fe_tools_extrapolate.cc @@ -39,7 +39,27 @@ namespace FETools namespace internal { -#ifdef DEAL_II_WITH_P4EST +#ifndef DEAL_II_WITH_P4EST + // Dummy implementation in case p4est is not available. + template + class ExtrapolateImplementation + { + public: + + ExtrapolateImplementation () + { + Assert(false, ExcNotImplemented()); + }; + + template + void extrapolate_parallel (const InVector &/*u2_relevant*/, + const DoFHandler &/*dof2*/, + OutVector &/*u2*/) + { + Assert(false, ExcNotImplemented()) + }; + }; +#else // Implementation of the @p extrapolate function // on parallel distributed grids. template @@ -1395,6 +1415,110 @@ namespace FETools namespace { + template + void reinit_distributed(const DH &/*dh*/, + VectorType &/*vector*/) + { + Assert(false, ExcNotImplemented()); + } + +#ifdef DEAL_II_WITH_PETSC + template + void reinit_distributed(const DoFHandler &dh, + PETScWrappers::MPI::Vector &vector) + { + const parallel::distributed::Triangulation *parallel_tria = + dynamic_cast*>(&dh.get_tria()); + Assert (parallel_tria !=0, ExcNotImplemented()); + + const IndexSet locally_owned_dofs = dh.locally_owned_dofs(); + vector.reinit(locally_owned_dofs, parallel_tria->get_communicator()); + } +#endif //DEAL_II_WITH_PETSC + +#ifdef DEAL_II_WITH_TRILINOS + template + void reinit_distributed(const DoFHandler &dh, + TrilinosWrappers::MPI::Vector &vector) + { + const parallel::distributed::Triangulation *parallel_tria = + dynamic_cast*>(&dh.get_tria()); + Assert (parallel_tria !=0, ExcNotImplemented()); + + const IndexSet locally_owned_dofs = dh.locally_owned_dofs(); + vector.reinit(locally_owned_dofs, parallel_tria->get_communicator()); + } +#endif //DEAL_II_WITH_TRILINOS + + template + void reinit_distributed(const DoFHandler &dh, + LinearAlgebra::distributed::Vector &vector) + { + const parallel::distributed::Triangulation *parallel_tria = + dynamic_cast*>(&dh.get_tria()); + Assert (parallel_tria !=0, ExcNotImplemented()); + + const IndexSet locally_owned_dofs = dh.locally_owned_dofs(); + vector.reinit(locally_owned_dofs, parallel_tria->get_communicator()); + } + + + + template + void reinit_ghosted(const DH &/*dh*/, + VectorType &/*vector*/) + { + Assert(false, ExcNotImplemented()); + } + +#ifdef DEAL_II_WITH_PETSC + template + void reinit_ghosted(const DoFHandler &dh, + PETScWrappers::MPI::Vector &vector) + { + const parallel::distributed::Triangulation *parallel_tria = + dynamic_cast*>(&dh.get_tria()); + Assert(parallel_tria !=0, ExcNotImplemented()); + const IndexSet locally_owned_dofs = dh.locally_owned_dofs(); + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs (dh, locally_relevant_dofs); + vector.reinit (locally_owned_dofs, locally_relevant_dofs, + parallel_tria->get_communicator()); + } +#endif //DEAL_II_WITH_PETSC + +#ifdef DEAL_II_WITH_TRILINOS + template + void reinit_ghosted(const DoFHandler &dh, + TrilinosWrappers::MPI::Vector &vector) + { + const parallel::distributed::Triangulation *parallel_tria = + dynamic_cast*>(&dh.get_tria()); + Assert (parallel_tria !=0, ExcNotImplemented()); + const IndexSet locally_owned_dofs = dh.locally_owned_dofs(); + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs (dh, locally_relevant_dofs); + vector.reinit (locally_owned_dofs, locally_relevant_dofs, + parallel_tria->get_communicator()); + } +#endif //DEAL_II_WITH_TRILINOS + + template + void reinit_ghosted(const DoFHandler &dh, + LinearAlgebra::distributed::Vector &vector) + { + const parallel::distributed::Triangulation *parallel_tria = + dynamic_cast*>(&dh.get_tria()); + Assert (parallel_tria !=0, ExcNotImplemented()); + const IndexSet locally_owned_dofs = dh.locally_owned_dofs(); + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs (dh, locally_relevant_dofs); + vector.reinit (locally_owned_dofs, locally_relevant_dofs, + parallel_tria->get_communicator()); + } + + + template void extrapolate_serial(const InVector &u3, const DoFHandler &dof2, @@ -1438,26 +1562,6 @@ namespace FETools } } } - - - -#ifdef DEAL_II_WITH_P4EST - template - void extrapolate_parallel(const InVector &u2_relevant, - const DoFHandler &dof2, - OutVector &u2) - { - if (dynamic_cast*>(&dof2.get_tria()) != 0) - { - ExtrapolateImplementation implementation; - implementation.extrapolate_parallel (u2_relevant, dof2, u2); - } - else - { - extrapolate_serial (u2_relevant, dof2, u2); - } - } -#endif //DEAL_II_WITH_P4EST } } @@ -1487,79 +1591,38 @@ namespace FETools Assert(u1.size()==dof1.n_dofs(), ExcDimensionMismatch(u1.size(), dof1.n_dofs())); Assert(u2.size()==dof2.n_dofs(), ExcDimensionMismatch(u2.size(), dof2.n_dofs())); - // make sure that each cell on the - // coarsest level is at least once - // refined, otherwise, these cells - // can't be treated and would - // generate a bogus result + // make sure that each cell on the coarsest level is at least once refined, otherwise, these cells + // can't be treated and would generate a bogus result { typename DoFHandler::cell_iterator cell = dof2.begin(0), endc = dof2.end(0); for (; cell!=endc; ++cell) - Assert (cell->has_children(), ExcGridNotRefinedAtLeastOnce()); + Assert (cell->has_children() || cell->is_artificial(), ExcGridNotRefinedAtLeastOnce()); } - OutVector u3; - u3.reinit(u2); - interpolate(dof1, u1, dof2, constraints, u3); - - internal::extrapolate_serial (u3, dof2, u2); - - constraints.distribute(u2); - } - - -#ifdef DEAL_II_WITH_P4EST - template - void extrapolate_parallel(const DoFHandler &dof1, - const VectorType &u1, - const DoFHandler &dof2, - const ConstraintMatrix &constraints2, - VectorType &u2) - { - IndexSet dof2_locally_owned_dofs = dof2.locally_owned_dofs(); - IndexSet dof2_locally_relevant_dofs; - DoFTools::extract_locally_relevant_dofs (dof2, - dof2_locally_relevant_dofs); - - VectorType u3 (dof2_locally_owned_dofs, u1.get_mpi_communicator ()); - interpolate (dof1, u1, dof2, constraints2, u3); - VectorType u3_relevant (dof2_locally_owned_dofs, - dof2_locally_relevant_dofs, - u1.get_mpi_communicator ()); - u3_relevant = u3; - - internal::extrapolate_parallel (u3_relevant, dof2, u2); - - constraints2.distribute(u2); - } - + OutVector u3; + internal::reinit_distributed(dof2, u3); + if (dynamic_cast*>(&dof2.get_tria()) != 0) + { + interpolate(dof1, u1, dof2, constraints, u3); - template - void extrapolate_parallel(const DoFHandler &dof1, - const LinearAlgebra::distributed::Vector &u1, - const DoFHandler &dof2, - const ConstraintMatrix &constraints2, - LinearAlgebra::distributed::Vector &u2) - { - IndexSet dof2_locally_owned_dofs = dof2.locally_owned_dofs(); - IndexSet dof2_locally_relevant_dofs; - DoFTools::extract_locally_relevant_dofs (dof2, - dof2_locally_relevant_dofs); + OutVector u3_relevant; + internal::reinit_ghosted(dof2, u3_relevant); + u3_relevant = u3; - LinearAlgebra::distributed::Vector u3(dof2_locally_owned_dofs, - dof2_locally_relevant_dofs, - u1.get_mpi_communicator ()); - interpolate (dof1, u1, dof2, constraints2, u3); - u3.update_ghost_values(); + internal::ExtrapolateImplementation implementation; + implementation.extrapolate_parallel (u3_relevant, dof2, u2); + } + else + { + interpolate(dof1, u1, dof2, constraints, u3); - internal::extrapolate_parallel (u3, dof2, u2); + internal::extrapolate_serial (u3, dof2, u2); + } - constraints2.distribute(u2); + constraints.distribute(u2); } -#endif //DEAL_II_WITH_P4EST - } // end of namespace FETools diff --git a/source/fe/fe_tools_extrapolate.inst.in b/source/fe/fe_tools_extrapolate.inst.in index dbcb9375d3..30dd33def6 100755 --- a/source/fe/fe_tools_extrapolate.inst.in +++ b/source/fe/fe_tools_extrapolate.inst.in @@ -33,58 +33,3 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS #endif \} } - - - -for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) -{ -#ifdef DEAL_II_WITH_P4EST - namespace FETools - \{ -#if deal_II_dimension == deal_II_space_dimension - -#ifdef DEAL_II_WITH_PETSC - template - void extrapolate_parallel - (const DoFHandler &, - const PETScWrappers::MPI::Vector &, - const DoFHandler &, - const ConstraintMatrix &, - PETScWrappers::MPI::Vector &); -#endif - -#ifdef DEAL_II_WITH_TRILINOS - template - void extrapolate_parallel - (const DoFHandler &, - const TrilinosWrappers::MPI::Vector &, - const DoFHandler &, - const ConstraintMatrix &, - TrilinosWrappers::MPI::Vector &); -#endif - -#endif - \} -#endif //DEAL_II_WITH_P4EST -} - - - -for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; scalar : REAL_SCALARS) -{ -#ifdef DEAL_II_WITH_P4EST - namespace FETools - \{ -#if deal_II_dimension == deal_II_space_dimension - template - void extrapolate_parallel - (const DoFHandler &, - const LinearAlgebra::distributed::Vector &, - const DoFHandler &, - const ConstraintMatrix &, - LinearAlgebra::distributed::Vector &); -#endif - \} -#endif //DEAL_II_WITH_P4EST -} - diff --git a/tests/mpi/fe_tools_extrapolate_common.h b/tests/mpi/fe_tools_extrapolate_common.h index 5cd255fbf3..3c9f7bf51d 100755 --- a/tests/mpi/fe_tools_extrapolate_common.h +++ b/tests/mpi/fe_tools_extrapolate_common.h @@ -37,6 +37,7 @@ #include #include +//#define DEBUG_OUTPUT_VTK template class TestFunction : public Function @@ -176,9 +177,11 @@ check_this (const FiniteElement &fe1, const double global_error_before = std::sqrt(Utilities::MPI::sum(std::pow(local_error_before,2.), MPI_COMM_WORLD)); - /*output_vector (in_ghosted, +#ifdef DEBUG_OUTPUT_VTK + output_vector (in_ghosted, Utilities::int_to_string(fe1.degree,1)+Utilities::int_to_string(dim,1)+std::string("in"), - *dof1);*/ + *dof1); +#endif { const double l1_norm = in_distributed.l1_norm(); const double l2_norm = in_distributed.l2_norm(); @@ -189,13 +192,14 @@ check_this (const FiniteElement &fe1, << linfty_norm << std::endl; } - FETools::extrapolate_parallel (*dof1, in_ghosted, *dof2, cm2, out_distributed); + FETools::extrapolate (*dof1, in_ghosted, *dof2, cm2, out_distributed); out_distributed.compress(VectorOperation::insert); out_ghosted = out_distributed; - - /*output_vector (out_ghosted, +#ifdef DEBUG_OUTPUT_VTK + output_vector (out_ghosted, Utilities::int_to_string(fe2.degree,1)+Utilities::int_to_string(dim,1)+std::string("out"), - *dof2);*/ + *dof2); +#endif { const double l1_norm = out_distributed.l1_norm(); @@ -297,10 +301,11 @@ check_this_dealii (const FiniteElement &fe1, const double local_error_before = difference_before.l2_norm(); const double global_error_before = std::sqrt(Utilities::MPI::sum(std::pow(local_error_before,2.), MPI_COMM_WORLD)); - - /*output_vector (in_ghosted, +#ifdef DEBUG_OUTPUT_VTK + output_vector (in_ghosted, Utilities::int_to_string(fe1.degree,1)+Utilities::int_to_string(dim,1)+std::string("in"), - *dof1);*/ + *dof1); +#endif { const double l1_norm = in_ghosted.l1_norm(); const double l2_norm = in_ghosted.l2_norm(); @@ -311,16 +316,14 @@ check_this_dealii (const FiniteElement &fe1, << linfty_norm << std::endl; } - if (Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD)==1) - FETools::extrapolate (*dof1, in_ghosted, *dof2, cm2, out_ghosted); - else - FETools::extrapolate_parallel (*dof1, in_ghosted, *dof2, cm2, out_ghosted); + FETools::extrapolate (*dof1, in_ghosted, *dof2, cm2, out_ghosted); out_ghosted.update_ghost_values(); - - /*output_vector (out_ghosted, +#ifdef DEBUG_OUTPUT_VTK + output_vector (out_ghosted, Utilities::int_to_string(fe2.degree,1)+Utilities::int_to_string(dim,1)+std::string("out"), - *dof2);*/ + *dof2); +#endif { const double l1_norm = out_ghosted.l1_norm(); const double l2_norm = out_ghosted.l2_norm();