From 1dbbb442e267a10a18f6d43ef128ad309919ac8f Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Wed, 18 Dec 2013 17:34:49 +0000 Subject: [PATCH] take over 32048,49 git-svn-id: https://svn.dealii.org/branches/releases/Branch-8-1@32050 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/8.0.0-vs-8.1.0.h | 6 + .../numerics/derivative_approximation.h | 2 +- .../numerics/derivative_approximation.cc | 24 ++-- tests/mpi/derivative_approximation_01.cc | 9 +- ...1.with_trilinos=true.mpirun=3.debug.output | 3 - tests/mpi/derivative_approximation_02.cc | 110 ++++++++++++++++++ ...tion_02.with_trilinos=true.mpirun=1.output | 3 + ...tion_02.with_trilinos=true.mpirun=3.output | 7 ++ ...tion_02.with_trilinos=true.mpirun=7.output | 15 +++ 9 files changed, 160 insertions(+), 19 deletions(-) create mode 100644 tests/mpi/derivative_approximation_02.cc create mode 100644 tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=1.output create mode 100644 tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=3.output create mode 100644 tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=7.output diff --git a/deal.II/doc/news/8.0.0-vs-8.1.0.h b/deal.II/doc/news/8.0.0-vs-8.1.0.h index bf47fc3f32..7be343c267 100644 --- a/deal.II/doc/news/8.0.0-vs-8.1.0.h +++ b/deal.II/doc/news/8.0.0-vs-8.1.0.h @@ -259,6 +259,12 @@ inconvenience this causes.

Specific improvements

    +
  1. Fixed: The DerivativeApproximation class did not work for + parallel programs. This is now fixed. +
    + (Wolfgang Bangerth, 2013/12/18) +
  2. +
  3. Fixed: Move the implementation of Subscriptor::(un)subscribe() to the .cc file so that it is possible to link against the debug library without specifying -DDEBUG diff --git a/deal.II/include/deal.II/numerics/derivative_approximation.h b/deal.II/include/deal.II/numerics/derivative_approximation.h index fb7062f195..002ec01137 100644 --- a/deal.II/include/deal.II/numerics/derivative_approximation.h +++ b/deal.II/include/deal.II/numerics/derivative_approximation.h @@ -568,7 +568,7 @@ private: template class DH, class InputVector, int spacedim> static void - approximate (SynchronousIterators::active_cell_iterator>, + approximate (SynchronousIterators::active_cell_iterator, Vector::iterator> > const &cell, const Mapping &mapping, const DH &dof, diff --git a/deal.II/source/numerics/derivative_approximation.cc b/deal.II/source/numerics/derivative_approximation.cc index dd2a9ca41d..44b5362740 100644 --- a/deal.II/source/numerics/derivative_approximation.cc +++ b/deal.II/source/numerics/derivative_approximation.cc @@ -646,19 +646,17 @@ approximate_derivative (const Mapping &mapping, Assert (component < dof_handler.get_fe().n_components(), ExcIndexRange (component, 0, dof_handler.get_fe().n_components())); - // Only act on the locally owned cells - typedef FilteredIterator::active_cell_iterator> CellFilter; - - typedef std_cxx1x::tuple::iterator> + typedef std_cxx1x::tuple::active_cell_iterator,Vector::iterator> Iterators; - SynchronousIterators begin(Iterators (CellFilter(IteratorFilters::LocallyOwnedCell(), - dof_handler.begin_active()),derivative_norm.begin())), - end(Iterators (CellFilter(IteratorFilters::LocallyOwnedCell(),dof_handler.end()), - derivative_norm.end())); + SynchronousIterators begin(Iterators(dof_handler.begin_active(), + derivative_norm.begin())), + end(Iterators(dof_handler.end(), + derivative_norm.end())); // There is no need for a copier because there is no conflict between threads // to write in derivative_norm. Scratch and CopyData are also useless. - WorkStream::run(begin,end, + WorkStream::run(begin, + end, static_cast const &, internal::Assembler::Scratch const &,internal::Assembler::CopyData &)> > (std_cxx1x::bind(&DerivativeApproximation::template approximate &mapping, template class DH, class InputVector, int spacedim> void -DerivativeApproximation::approximate (SynchronousIterators::active_cell_iterator>,Vector::iterator> > const &cell, +DerivativeApproximation::approximate (SynchronousIterators::active_cell_iterator,Vector::iterator> > const &cell, const Mapping &mapping, const DH &dof_handler, const InputVector &solution, const unsigned int component) { + // if the cell is not locally owned, then there is nothing to do + if (std_cxx1x::get<0>(cell.iterators)->is_locally_owned() == false) + *std_cxx1x::get<1>(cell.iterators) = 0; + else + { typename DerivativeDescription::Derivative derivative; // call the function doing the actual // work on this cell @@ -690,6 +693,7 @@ DerivativeApproximation::approximate (SynchronousIterators(cell.iterators) = DerivativeDescription::derivative_norm (derivative); + } } diff --git a/tests/mpi/derivative_approximation_01.cc b/tests/mpi/derivative_approximation_01.cc index 11d0f1d958..955655f3ac 100644 --- a/tests/mpi/derivative_approximation_01.cc +++ b/tests/mpi/derivative_approximation_01.cc @@ -21,14 +21,14 @@ -------------------------------------------------------- An error occurred in line <3349> of file in function void dealii::DoFCellAccessor::get_dof_values(const InputVector&, ForwardIterator, ForwardIterator) const [with InputVector = dealii::TrilinosWrappers::MPI::Vector, ForwardIterator = double*, DH = dealii::DoFHandler<2>, bool level_dof_access = false] - The violated condition was: + The violated condition was: this->is_artificial() == false The name and call sequence of the exception was: ExcMessage ("Can't ask for DoF indices on artificial cells.") - Additional Information: + Additional Information: Can't ask for DoF indices on artificial cells. -------------------------------------------------------- - + */ #include "../tests.h" @@ -85,8 +85,7 @@ void test() vec_rel, indicators); - deallog <<"output:" << indicators.l2_norm() << std::endl; - + // we got here, so no exception. if (myid == 0) deallog << "OK" << std::endl; } diff --git a/tests/mpi/derivative_approximation_01.with_trilinos=true.mpirun=3.debug.output b/tests/mpi/derivative_approximation_01.with_trilinos=true.mpirun=3.debug.output index 0deadca86b..03cc72804a 100644 --- a/tests/mpi/derivative_approximation_01.with_trilinos=true.mpirun=3.debug.output +++ b/tests/mpi/derivative_approximation_01.with_trilinos=true.mpirun=3.debug.output @@ -1,9 +1,6 @@ -DEAL:0::output:34.4565 DEAL:0::OK -DEAL:1::output:61.6543 -DEAL:2::output:47.0585 diff --git a/tests/mpi/derivative_approximation_02.cc b/tests/mpi/derivative_approximation_02.cc new file mode 100644 index 0000000000..af7e1b11e1 --- /dev/null +++ b/tests/mpi/derivative_approximation_02.cc @@ -0,0 +1,110 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2009 - 2013 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + + + +// DerivativeApproximation didn't work in parallel at all. This test verifies +// that it now does. + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + + +template +class Quadratic : + public Function +{ +public: + double value (const Point &p, const unsigned int) const + { + return p*p; + } +}; + + +template +void test() +{ + const unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + const unsigned int n_processes = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); + + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + GridGenerator::hyper_cube(tr); + tr.refine_global(2); + + const FE_Q fe(1); + DoFHandler dofh(tr); + dofh.distribute_dofs (fe); + + IndexSet locally_relevant_set; + DoFTools::extract_locally_relevant_dofs (dofh, + locally_relevant_set); + + // create a vector representing a function that is independent of the number + // of processors involved + TrilinosWrappers::MPI::Vector vec (dofh.locally_owned_dofs(), MPI_COMM_WORLD); + VectorTools::interpolate (dofh, Quadratic(), vec); + vec.compress(VectorOperation::insert); + + // create such a vector with ghost elements and use it to compute + // derivative information + TrilinosWrappers::MPI::Vector vec_rel ( locally_relevant_set); + vec_rel = vec; + + Vector indicators(tr.n_active_cells()); + DerivativeApproximation::approximate_gradient (dofh, + vec_rel, + indicators); + + // what we get must be a set of derivative indicators, one for each + // cell of the distributed mesh. they need to be the same, no matter + // how many processors we use. So, the sum of absolute values must + // be the same for any processor number + const double sum = Utilities::MPI::sum (indicators.l1_norm(), MPI_COMM_WORLD); + if (myid == 0) + deallog << sum << std::endl; +} + + + +int main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + MPILogInitAll log; + + test<2> (); + test<3> (); +} + diff --git a/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=1.output b/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=1.output new file mode 100644 index 0000000000..380b2ea54c --- /dev/null +++ b/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=1.output @@ -0,0 +1,3 @@ + +DEAL:0::23.4892 +DEAL:0::116.543 diff --git a/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=3.output b/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=3.output new file mode 100644 index 0000000000..ece4add40a --- /dev/null +++ b/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=3.output @@ -0,0 +1,7 @@ + +DEAL:0::23.4892 +DEAL:0::116.543 + + + + diff --git a/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=7.output b/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=7.output new file mode 100644 index 0000000000..ce82c719f1 --- /dev/null +++ b/tests/mpi/derivative_approximation_02.with_trilinos=true.mpirun=7.output @@ -0,0 +1,15 @@ + +DEAL:0::23.4892 +DEAL:0::116.543 + + + + + + + + + + + + -- 2.39.5