]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add a test that truly tests hp::DoFHandler on shared triangulations. 4593/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 10 Jul 2017 22:45:29 +0000 (16:45 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Mon, 10 Jul 2017 22:45:29 +0000 (16:45 -0600)
In particular, extend a previously added test to really do use different
finite elements on different cells.

Comparing the output for these two tests shows that indeed
the test with different (and higher order) elements is more
accurate, as one might have surmised.

tests/vector_tools/integrate_difference_04_hp_02.cc [new file with mode: 0644]
tests/vector_tools/integrate_difference_04_hp_02.with_trilinos=true.with_metis=true.mpirun=2.output [new file with mode: 0644]
tests/vector_tools/integrate_difference_04_hp_02.with_trilinos=true.with_metis=true.mpirun=3.output [new file with mode: 0644]

diff --git a/tests/vector_tools/integrate_difference_04_hp_02.cc b/tests/vector_tools/integrate_difference_04_hp_02.cc
new file mode 100644 (file)
index 0000000..2c69e4e
--- /dev/null
@@ -0,0 +1,177 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2003 - 2016 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Test integrate_difference and compute_global_error in parallel
+// see integrate_difference_02.cc for the serial version
+//
+// like the integrate_difference_04_hp test, but this time really do
+// use multiple different finite elements on different
+// cells. comparing the output for these two tests shows that indeed
+// the test with different (and higher order) elements is more
+// accurate, as one might have surmised.
+
+#include "../tests.h"
+#include <deal.II/base/function.h>
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/numerics/vector_tools.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/trilinos_vector.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/hp/q_collection.h>
+#include <deal.II/distributed/shared_tria.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_refinement.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+
+using namespace dealii;
+
+
+// x+y+z, x^2+y^2, z+xy
+// div = 1+2y+1
+template <int dim>
+class Ref : public Function<dim>
+{
+public:
+  Ref()
+    :Function<dim>(dim)
+  {}
+
+  double value (const Point<dim> &p, const unsigned int c) const
+  {
+    if (c==0)
+      return p[0]+p[1]+((dim==3)?p[2]:0.0);
+    if (c==1)
+      return p[0]*p[0]+p[1]*p[1];
+    if (c==2)
+      return p[2]+p[0]*p[1];
+    return 0.0;
+  }
+};
+
+
+
+template <int dim>
+void test(VectorTools::NormType norm, double value, double exp = 2.0)
+{
+  parallel::shared::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(3);
+
+  hp::FECollection<dim> fe;
+  fe.push_back (FESystem<dim> (FE_Q<dim>(4),dim));
+  fe.push_back (FESystem<dim> (FE_Q<dim>(5),dim));
+  fe.push_back (FESystem<dim> (FE_Q<dim>(6),dim));
+  hp::DoFHandler<dim> dofh(tria);
+
+  // assign FEs mostly randomly to each cell
+  for (auto cell : dofh.active_cell_iterators())
+    cell->set_active_fe_index (cell->active_cell_index() %
+                               fe.size());
+  dofh.distribute_dofs(fe);
+
+  TrilinosWrappers::MPI::Vector interpolated(dofh.locally_owned_dofs(),
+                                             MPI_COMM_WORLD);
+
+  VectorTools::interpolate(dofh, Ref<dim>(), interpolated);
+
+  IndexSet relevant_set;
+  DoFTools::extract_locally_relevant_dofs (dofh, relevant_set);
+  TrilinosWrappers::MPI::Vector solution(relevant_set, MPI_COMM_WORLD);
+  solution = interpolated;
+
+  Vector<double> cellwise_errors (tria.n_active_cells());
+  hp::QCollection<dim> quadrature;
+  quadrature.push_back (QIterated<dim> (QTrapez<1>(), 5));
+  quadrature.push_back (QIterated<dim> (QTrapez<1>(), 6));
+  quadrature.push_back (QIterated<dim> (QTrapez<1>(), 7));
+
+  const dealii::Function<dim,double> *w = nullptr;
+  VectorTools::integrate_difference (dofh,
+                                     solution,
+                                     ZeroFunction<dim>(dim),
+                                     cellwise_errors,
+                                     quadrature,
+                                     norm,
+                                     w,
+                                     exp);
+
+  const double error
+    = VectorTools::compute_global_error(tria, cellwise_errors, norm, exp);
+
+  const double difference = std::abs(error-value);
+  deallog << "computed: " << error
+          << " expected: " << value
+          << " difference: " << difference
+          << std::endl;
+  Assert(difference<2e-3, ExcMessage("Error in integrate_difference"));
+}
+
+template <int dim>
+void test()
+{
+  deallog << "Hdiv_seminorm:" << std::endl;
+  // sqrt(\int (div f)^2 = sqrt(\int (1+2y+1)^2)
+  test<dim>(VectorTools::Hdiv_seminorm, 2.0*std::sqrt(7.0/3.0));
+
+  deallog << "L2_norm:" << std::endl;
+  // sqrt(\int_\Omega f^2) = sqrt(\int (x+y+z)^2+(x^2+y^2)^2+(z+xy)^2)
+  test<dim>(VectorTools::L2_norm, std::sqrt(229.0/60.0));
+
+  deallog << "H1_seminorm:" << std::endl;
+  // sqrt( \int sum_k | d/dxi f_k |_0^2 )
+  // = sqrt( \int 3+ (2x)^2+(2y)^2 + y^2+x^2+1
+  // = sqrt( 22/3  )
+  test<dim>(VectorTools::H1_seminorm, std::sqrt(22.0/3.0));
+
+  deallog << "H1_norm:" << std::endl;
+  test<dim>(VectorTools::H1_norm, std::sqrt(229.0/60.0+22.0/3.0));
+
+  deallog << "L1_norm:" << std::endl;
+  test<dim>(VectorTools::L1_norm, 35.0/12.0);
+
+  deallog << "Linfty_norm:" << std::endl;
+  test<dim>(VectorTools::Linfty_norm, 3.0);
+
+  deallog << "mean:" << std::endl;
+  // int -(x+y+z + x^2+y^2 + z+xy)
+  test<dim>(VectorTools::mean, -35.0/12.0);
+
+  deallog << "Lp_norm:" << std::endl;
+  // (int (x+y+z)^p+(x^2+y^2)^p+(z+xy)^p) ) ^ 1./p
+  // = std::pow(9937.0/1680.0, 1.0/3.0)
+  test<dim>(VectorTools::Lp_norm, std::pow(9937.0/1680.0, 1.0/3.0), 3.0);
+
+  deallog << "W1p_seminorm:" << std::endl;
+  // ( \int_K sum_k | d/dxi f_k |_0^2^(p/2) )^1/p
+  // ( integrate 3^(3/2) + ((2x)^2+(2y)^2)^(3/2) + (y^2+x^2+1)^(3/2) ) ^1/p
+  // = (12.4164) ^1/3 = 2.31560
+  test<dim>(VectorTools::W1p_seminorm, 2.31560, 3.0);
+
+  deallog << "OK" << std::endl;
+}
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+  MPILogInitAll log;
+  test<3>();
+}
diff --git a/tests/vector_tools/integrate_difference_04_hp_02.with_trilinos=true.with_metis=true.mpirun=2.output b/tests/vector_tools/integrate_difference_04_hp_02.with_trilinos=true.with_metis=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..5c766ac
--- /dev/null
@@ -0,0 +1,41 @@
+
+DEAL:0::Hdiv_seminorm:
+DEAL:0::computed: 3.05510 expected: 3.05505 difference: 5.01434e-05
+DEAL:0::L2_norm:
+DEAL:0::computed: 1.95383 expected: 1.95363 difference: 0.000196056
+DEAL:0::H1_seminorm:
+DEAL:0::computed: 2.70815 expected: 2.70801 difference: 0.000141421
+DEAL:0::H1_norm:
+DEAL:0::computed: 3.33939 expected: 3.33916 difference: 0.000229397
+DEAL:0::L1_norm:
+DEAL:0::computed: 2.91682 expected: 2.91667 difference: 0.000153192
+DEAL:0::Linfty_norm:
+DEAL:0::computed: 3.00000 expected: 3.00000 difference: 0.00000
+DEAL:0::mean:
+DEAL:0::computed: -2.91682 expected: -2.91667 difference: 0.000153192
+DEAL:0::Lp_norm:
+DEAL:0::computed: 1.80871 expected: 1.80849 difference: 0.000224088
+DEAL:0::W1p_seminorm:
+DEAL:0::computed: 2.31576 expected: 2.31560 difference: 0.000161338
+DEAL:0::OK
+
+DEAL:1::Hdiv_seminorm:
+DEAL:1::computed: 3.05510 expected: 3.05505 difference: 5.01434e-05
+DEAL:1::L2_norm:
+DEAL:1::computed: 1.95383 expected: 1.95363 difference: 0.000196056
+DEAL:1::H1_seminorm:
+DEAL:1::computed: 2.70815 expected: 2.70801 difference: 0.000141421
+DEAL:1::H1_norm:
+DEAL:1::computed: 3.33939 expected: 3.33916 difference: 0.000229397
+DEAL:1::L1_norm:
+DEAL:1::computed: 2.91682 expected: 2.91667 difference: 0.000153192
+DEAL:1::Linfty_norm:
+DEAL:1::computed: 3.00000 expected: 3.00000 difference: 0.00000
+DEAL:1::mean:
+DEAL:1::computed: -2.91682 expected: -2.91667 difference: 0.000153192
+DEAL:1::Lp_norm:
+DEAL:1::computed: 1.80871 expected: 1.80849 difference: 0.000224088
+DEAL:1::W1p_seminorm:
+DEAL:1::computed: 2.31576 expected: 2.31560 difference: 0.000161338
+DEAL:1::OK
+
diff --git a/tests/vector_tools/integrate_difference_04_hp_02.with_trilinos=true.with_metis=true.mpirun=3.output b/tests/vector_tools/integrate_difference_04_hp_02.with_trilinos=true.with_metis=true.mpirun=3.output
new file mode 100644 (file)
index 0000000..e96aee4
--- /dev/null
@@ -0,0 +1,62 @@
+
+DEAL:0::Hdiv_seminorm:
+DEAL:0::computed: 3.05510 expected: 3.05505 difference: 5.01434e-05
+DEAL:0::L2_norm:
+DEAL:0::computed: 1.95383 expected: 1.95363 difference: 0.000196056
+DEAL:0::H1_seminorm:
+DEAL:0::computed: 2.70815 expected: 2.70801 difference: 0.000141421
+DEAL:0::H1_norm:
+DEAL:0::computed: 3.33939 expected: 3.33916 difference: 0.000229397
+DEAL:0::L1_norm:
+DEAL:0::computed: 2.91682 expected: 2.91667 difference: 0.000153192
+DEAL:0::Linfty_norm:
+DEAL:0::computed: 3.00000 expected: 3.00000 difference: 0.00000
+DEAL:0::mean:
+DEAL:0::computed: -2.91682 expected: -2.91667 difference: 0.000153192
+DEAL:0::Lp_norm:
+DEAL:0::computed: 1.80871 expected: 1.80849 difference: 0.000224088
+DEAL:0::W1p_seminorm:
+DEAL:0::computed: 2.31576 expected: 2.31560 difference: 0.000161338
+DEAL:0::OK
+
+DEAL:1::Hdiv_seminorm:
+DEAL:1::computed: 3.05510 expected: 3.05505 difference: 5.01434e-05
+DEAL:1::L2_norm:
+DEAL:1::computed: 1.95383 expected: 1.95363 difference: 0.000196056
+DEAL:1::H1_seminorm:
+DEAL:1::computed: 2.70815 expected: 2.70801 difference: 0.000141421
+DEAL:1::H1_norm:
+DEAL:1::computed: 3.33939 expected: 3.33916 difference: 0.000229397
+DEAL:1::L1_norm:
+DEAL:1::computed: 2.91682 expected: 2.91667 difference: 0.000153192
+DEAL:1::Linfty_norm:
+DEAL:1::computed: 3.00000 expected: 3.00000 difference: 0.00000
+DEAL:1::mean:
+DEAL:1::computed: -2.91682 expected: -2.91667 difference: 0.000153192
+DEAL:1::Lp_norm:
+DEAL:1::computed: 1.80871 expected: 1.80849 difference: 0.000224088
+DEAL:1::W1p_seminorm:
+DEAL:1::computed: 2.31576 expected: 2.31560 difference: 0.000161338
+DEAL:1::OK
+
+
+DEAL:2::Hdiv_seminorm:
+DEAL:2::computed: 3.05510 expected: 3.05505 difference: 5.01434e-05
+DEAL:2::L2_norm:
+DEAL:2::computed: 1.95383 expected: 1.95363 difference: 0.000196056
+DEAL:2::H1_seminorm:
+DEAL:2::computed: 2.70815 expected: 2.70801 difference: 0.000141421
+DEAL:2::H1_norm:
+DEAL:2::computed: 3.33939 expected: 3.33916 difference: 0.000229397
+DEAL:2::L1_norm:
+DEAL:2::computed: 2.91682 expected: 2.91667 difference: 0.000153192
+DEAL:2::Linfty_norm:
+DEAL:2::computed: 3.00000 expected: 3.00000 difference: 0.00000
+DEAL:2::mean:
+DEAL:2::computed: -2.91682 expected: -2.91667 difference: 0.000153192
+DEAL:2::Lp_norm:
+DEAL:2::computed: 1.80871 expected: 1.80849 difference: 0.000224088
+DEAL:2::W1p_seminorm:
+DEAL:2::computed: 2.31576 expected: 2.31560 difference: 0.000161338
+DEAL:2::OK
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.