]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add tests
authorDaniel Arndt <d.arndt@math.uni-goettingen.de>
Fri, 2 Sep 2016 21:18:39 +0000 (23:18 +0200)
committerDaniel Arndt <d.arndt@math.uni-goettingen.de>
Mon, 3 Oct 2016 16:20:58 +0000 (18:20 +0200)
tests/mpi/fe_tools_extrapolate_01.cc [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_01.with_trilinos=true.mpirun=1.output [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_01.with_trilinos=true.mpirun=3.output [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_02.cc [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_02.with_petsc=true.mpirun=1.output [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_02.with_petsc=true.mpirun=3.output [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_03.cc [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_03.mpirun=1.output [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_03.mpirun=3.output [new file with mode: 0644]
tests/mpi/fe_tools_extrapolate_common.h [new file with mode: 0644]

diff --git a/tests/mpi/fe_tools_extrapolate_01.cc b/tests/mpi/fe_tools_extrapolate_01.cc
new file mode 100644 (file)
index 0000000..a77c7de
--- /dev/null
@@ -0,0 +1,67 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/lac/trilinos_vector.h>
+#include "fe_tools_extrapolate_common.h"
+
+// check FETools::extrapolate on distributed triangulations
+// for TrilinosWrappers::MPI::Vector
+
+template <int dim>
+void
+check (const FiniteElement<dim> &fe1,
+       const FiniteElement<dim> &fe2,
+       const std::string        &name)
+{
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+    deallog << "Checking " << name
+            << " in " << dim << "d:"
+            << std::endl;
+
+  // call main function in .cc files
+  check_this<dim, TrilinosWrappers::MPI::Vector> (fe1, fe2);
+}
+
+#define CHECK(EL1,deg1,EL2,deg2,dim)\
+  { FE_ ## EL1<dim> fe1(deg1);   \
+    FE_ ## EL2<dim> fe2(deg2);   \
+    check(fe1, fe2, #EL1 #deg1 " against " #EL2 #deg2); \
+  }
+
+#define CHECK_ALL(EL1,deg1,EL2,deg2)\
+  { CHECK(EL1,deg1,EL2,deg2,2); \
+    CHECK(EL1,deg1,EL2,deg2,3); \
+  }
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+
+  MPILogInitAll log;
+  deallog.threshold_double(1.e-10);
+
+  CHECK_ALL(Q,1,Q,1);
+  CHECK_ALL(Q,1,Q,2);
+  CHECK_ALL(Q,2,Q,1);
+  CHECK_ALL(Q,2,Q,2);
+
+  CHECK_ALL(DGQ,0,DGQ,0);
+  CHECK_ALL(DGQ,0,DGQ,1);
+  CHECK_ALL(DGQ,1,DGQ,0);
+  CHECK_ALL(DGQ,1,DGQ,1);
+}
+
diff --git a/tests/mpi/fe_tools_extrapolate_01.with_trilinos=true.mpirun=1.output b/tests/mpi/fe_tools_extrapolate_01.with_trilinos=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..39f85f6
--- /dev/null
@@ -0,0 +1,97 @@
+
+DEAL:0::Checking Q1 against Q1 in 2d:
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q1 in 3d:
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 2d:
+DEAL:0::9.125000000 1.682283604 0.5000000000
+DEAL:0::26.88281250 2.700093063 0.5000000000
+DEAL:0::global_error_before: 0.01934481330
+DEAL:0::global_error_after: 0.001208379635
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 3d:
+DEAL:0::66.28125000 5.415995551 0.7500000000
+DEAL:0::342.1987305 11.49092805 0.7500000000
+DEAL:0::global_error_before: 0.03152756518
+DEAL:0::global_error_after: 0.001607487444
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::10.25000000 1.831793861 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.07284222390
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::75.09375000 5.992915740 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1190491594
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 2d:
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 3d:
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 2d:
+DEAL:0::14.00000000 2.828427125 0.8750000000
+DEAL:0::56.00000000 5.678908346 0.9375000000
+DEAL:0::global_error_before: 0.09199750902
+DEAL:0::global_error_after: 0.04599875451
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 3d:
+DEAL:0::72.50000000 7.854139036 1.312500000
+DEAL:0::580.6250000 22.31268054 1.406250000
+DEAL:0::global_error_before: 0.1220351512
+DEAL:0::global_error_after: 0.06040494959
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::14.00000000 2.738612788 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1839950180
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::73.00000000 7.697402159 1.125000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.2440703024
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
diff --git a/tests/mpi/fe_tools_extrapolate_01.with_trilinos=true.mpirun=3.output b/tests/mpi/fe_tools_extrapolate_01.with_trilinos=true.mpirun=3.output
new file mode 100644 (file)
index 0000000..8f5cfcc
--- /dev/null
@@ -0,0 +1,101 @@
+
+DEAL:0::Checking Q1 against Q1 in 2d:
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q1 in 3d:
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 2d:
+DEAL:0::9.125000000 1.682283604 0.5000000000
+DEAL:0::26.88281250 2.700093063 0.5000000000
+DEAL:0::global_error_before: 0.01934481330
+DEAL:0::global_error_after: 0.001208379635
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 3d:
+DEAL:0::66.28125000 5.415995551 0.7500000000
+DEAL:0::342.1987305 11.49092805 0.7500000000
+DEAL:0::global_error_before: 0.03152756518
+DEAL:0::global_error_after: 0.001607487444
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::10.25000000 1.831793861 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.07284222390
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::75.09375000 5.992915740 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1190491594
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 2d:
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 3d:
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 2d:
+DEAL:0::14.00000000 2.828427125 0.8750000000
+DEAL:0::56.00000000 5.678908346 0.9375000000
+DEAL:0::global_error_before: 0.09199750902
+DEAL:0::global_error_after: 0.04599875451
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 3d:
+DEAL:0::72.50000000 7.854139036 1.312500000
+DEAL:0::580.6250000 22.31268054 1.406250000
+DEAL:0::global_error_before: 0.1220351512
+DEAL:0::global_error_after: 0.06040494959
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::14.00000000 2.738612788 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1839950180
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::73.00000000 7.697402159 1.125000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.2440703024
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+
+
+
+
diff --git a/tests/mpi/fe_tools_extrapolate_02.cc b/tests/mpi/fe_tools_extrapolate_02.cc
new file mode 100644 (file)
index 0000000..12b64b6
--- /dev/null
@@ -0,0 +1,67 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/lac/petsc_vector.h>
+#include "fe_tools_extrapolate_common.h"
+
+// check FETools::extrapolate on distributed triangulations
+// for PETScWrappers::MPI::Vector
+
+template <int dim>
+void
+check (const FiniteElement<dim> &fe1,
+       const FiniteElement<dim> &fe2,
+       const std::string        &name)
+{
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+    deallog << "Checking " << name
+            << " in " << dim << "d:"
+            << std::endl;
+
+  // call main function in .cc files
+  check_this<dim, PETScWrappers::MPI::Vector> (fe1, fe2);
+}
+
+#define CHECK(EL1,deg1,EL2,deg2,dim)\
+  { FE_ ## EL1<dim> fe1(deg1);   \
+    FE_ ## EL2<dim> fe2(deg2);   \
+    check(fe1, fe2, #EL1 #deg1 " against " #EL2 #deg2); \
+  }
+
+#define CHECK_ALL(EL1,deg1,EL2,deg2)\
+  { CHECK(EL1,deg1,EL2,deg2,2); \
+    CHECK(EL1,deg1,EL2,deg2,3); \
+  }
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+
+  MPILogInitAll log;
+  deallog.threshold_double(1.e-10);
+
+  CHECK_ALL(Q,1,Q,1);
+  CHECK_ALL(Q,1,Q,2);
+  CHECK_ALL(Q,2,Q,1);
+  CHECK_ALL(Q,2,Q,2);
+
+  CHECK_ALL(DGQ,0,DGQ,0);
+  CHECK_ALL(DGQ,0,DGQ,1);
+  CHECK_ALL(DGQ,1,DGQ,0);
+  CHECK_ALL(DGQ,1,DGQ,1);
+}
+
diff --git a/tests/mpi/fe_tools_extrapolate_02.with_petsc=true.mpirun=1.output b/tests/mpi/fe_tools_extrapolate_02.with_petsc=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..39f85f6
--- /dev/null
@@ -0,0 +1,97 @@
+
+DEAL:0::Checking Q1 against Q1 in 2d:
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q1 in 3d:
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 2d:
+DEAL:0::9.125000000 1.682283604 0.5000000000
+DEAL:0::26.88281250 2.700093063 0.5000000000
+DEAL:0::global_error_before: 0.01934481330
+DEAL:0::global_error_after: 0.001208379635
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 3d:
+DEAL:0::66.28125000 5.415995551 0.7500000000
+DEAL:0::342.1987305 11.49092805 0.7500000000
+DEAL:0::global_error_before: 0.03152756518
+DEAL:0::global_error_after: 0.001607487444
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::10.25000000 1.831793861 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.07284222390
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::75.09375000 5.992915740 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1190491594
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 2d:
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 3d:
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 2d:
+DEAL:0::14.00000000 2.828427125 0.8750000000
+DEAL:0::56.00000000 5.678908346 0.9375000000
+DEAL:0::global_error_before: 0.09199750902
+DEAL:0::global_error_after: 0.04599875451
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 3d:
+DEAL:0::72.50000000 7.854139036 1.312500000
+DEAL:0::580.6250000 22.31268054 1.406250000
+DEAL:0::global_error_before: 0.1220351512
+DEAL:0::global_error_after: 0.06040494959
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::14.00000000 2.738612788 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1839950180
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::73.00000000 7.697402159 1.125000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.2440703024
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
diff --git a/tests/mpi/fe_tools_extrapolate_02.with_petsc=true.mpirun=3.output b/tests/mpi/fe_tools_extrapolate_02.with_petsc=true.mpirun=3.output
new file mode 100644 (file)
index 0000000..8f5cfcc
--- /dev/null
@@ -0,0 +1,101 @@
+
+DEAL:0::Checking Q1 against Q1 in 2d:
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q1 in 3d:
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 2d:
+DEAL:0::9.125000000 1.682283604 0.5000000000
+DEAL:0::26.88281250 2.700093063 0.5000000000
+DEAL:0::global_error_before: 0.01934481330
+DEAL:0::global_error_after: 0.001208379635
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 3d:
+DEAL:0::66.28125000 5.415995551 0.7500000000
+DEAL:0::342.1987305 11.49092805 0.7500000000
+DEAL:0::global_error_before: 0.03152756518
+DEAL:0::global_error_after: 0.001607487444
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::10.25000000 1.831793861 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.07284222390
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::75.09375000 5.992915740 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1190491594
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 2d:
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 3d:
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 2d:
+DEAL:0::14.00000000 2.828427125 0.8750000000
+DEAL:0::56.00000000 5.678908346 0.9375000000
+DEAL:0::global_error_before: 0.09199750902
+DEAL:0::global_error_after: 0.04599875451
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 3d:
+DEAL:0::72.50000000 7.854139036 1.312500000
+DEAL:0::580.6250000 22.31268054 1.406250000
+DEAL:0::global_error_before: 0.1220351512
+DEAL:0::global_error_after: 0.06040494959
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::14.00000000 2.738612788 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1839950180
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::73.00000000 7.697402159 1.125000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.2440703024
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+
+
+
+
diff --git a/tests/mpi/fe_tools_extrapolate_03.cc b/tests/mpi/fe_tools_extrapolate_03.cc
new file mode 100644 (file)
index 0000000..3c805b9
--- /dev/null
@@ -0,0 +1,67 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/lac/la_parallel_vector.h>
+#include "fe_tools_extrapolate_common.h"
+
+// check FETools::extrapolate on distributed triangulations
+// for LinearAlgebra::distributed::Vector<double>
+
+template <int dim>
+void
+check (const FiniteElement<dim> &fe1,
+       const FiniteElement<dim> &fe2,
+       const std::string        &name)
+{
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+    deallog << "Checking " << name
+            << " in " << dim << "d:"
+            << std::endl;
+
+  // call main function in .cc files
+  check_this_dealii<dim, LinearAlgebra::distributed::Vector<double> > (fe1, fe2);
+}
+
+#define CHECK(EL1,deg1,EL2,deg2,dim)\
+  { FE_ ## EL1<dim> fe1(deg1);   \
+    FE_ ## EL2<dim> fe2(deg2);   \
+    check(fe1, fe2, #EL1 #deg1 " against " #EL2 #deg2); \
+  }
+
+#define CHECK_ALL(EL1,deg1,EL2,deg2)\
+  { CHECK(EL1,deg1,EL2,deg2,2); \
+    CHECK(EL1,deg1,EL2,deg2,3); \
+  }
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+
+  MPILogInitAll log;
+  deallog.threshold_double(1.e-10);
+
+  CHECK_ALL(Q,1,Q,1);
+  CHECK_ALL(Q,1,Q,2);
+  CHECK_ALL(Q,2,Q,1);
+  CHECK_ALL(Q,2,Q,2);
+
+  CHECK_ALL(DGQ,0,DGQ,0);
+  CHECK_ALL(DGQ,0,DGQ,1);
+  CHECK_ALL(DGQ,1,DGQ,0);
+  CHECK_ALL(DGQ,1,DGQ,1);
+}
+
diff --git a/tests/mpi/fe_tools_extrapolate_03.mpirun=1.output b/tests/mpi/fe_tools_extrapolate_03.mpirun=1.output
new file mode 100644 (file)
index 0000000..39f85f6
--- /dev/null
@@ -0,0 +1,97 @@
+
+DEAL:0::Checking Q1 against Q1 in 2d:
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q1 in 3d:
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 2d:
+DEAL:0::9.125000000 1.682283604 0.5000000000
+DEAL:0::26.88281250 2.700093063 0.5000000000
+DEAL:0::global_error_before: 0.01934481330
+DEAL:0::global_error_after: 0.001208379635
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 3d:
+DEAL:0::66.28125000 5.415995551 0.7500000000
+DEAL:0::342.1987305 11.49092805 0.7500000000
+DEAL:0::global_error_before: 0.03152756518
+DEAL:0::global_error_after: 0.001607487444
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::10.25000000 1.831793861 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.07284222390
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::75.09375000 5.992915740 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1190491594
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 2d:
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 3d:
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 2d:
+DEAL:0::14.00000000 2.828427125 0.8750000000
+DEAL:0::56.00000000 5.678908346 0.9375000000
+DEAL:0::global_error_before: 0.09199750902
+DEAL:0::global_error_after: 0.04599875451
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 3d:
+DEAL:0::72.50000000 7.854139036 1.312500000
+DEAL:0::580.6250000 22.31268054 1.406250000
+DEAL:0::global_error_before: 0.1220351512
+DEAL:0::global_error_after: 0.06040494959
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::14.00000000 2.738612788 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1839950180
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::73.00000000 7.697402159 1.125000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.2440703024
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
diff --git a/tests/mpi/fe_tools_extrapolate_03.mpirun=3.output b/tests/mpi/fe_tools_extrapolate_03.mpirun=3.output
new file mode 100644 (file)
index 0000000..8f5cfcc
--- /dev/null
@@ -0,0 +1,101 @@
+
+DEAL:0::Checking Q1 against Q1 in 2d:
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::23.00000000 3.944933460 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q1 in 3d:
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::162.3750000 12.71625436 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 2d:
+DEAL:0::9.125000000 1.682283604 0.5000000000
+DEAL:0::26.88281250 2.700093063 0.5000000000
+DEAL:0::global_error_before: 0.01934481330
+DEAL:0::global_error_after: 0.001208379635
+DEAL:0::OK
+DEAL:0::Checking Q1 against Q2 in 3d:
+DEAL:0::66.28125000 5.415995551 0.7500000000
+DEAL:0::342.1987305 11.49092805 0.7500000000
+DEAL:0::global_error_before: 0.03152756518
+DEAL:0::global_error_after: 0.001607487444
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::10.25000000 1.831793861 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.07284222390
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q1 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::75.09375000 5.992915740 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1190491594
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 2d:
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::26.84375000 2.699035020 0.5000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking Q2 against Q2 in 3d:
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::341.5859375 11.47797042 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 2d:
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::56.00000000 10.58300524 2.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ0 in 3d:
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::276.0000000 28.77498914 3.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 2d:
+DEAL:0::14.00000000 2.828427125 0.8750000000
+DEAL:0::56.00000000 5.678908346 0.9375000000
+DEAL:0::global_error_before: 0.09199750902
+DEAL:0::global_error_after: 0.04599875451
+DEAL:0::OK
+DEAL:0::Checking DGQ0 against DGQ1 in 3d:
+DEAL:0::72.50000000 7.854139036 1.312500000
+DEAL:0::580.6250000 22.31268054 1.406250000
+DEAL:0::global_error_before: 0.1220351512
+DEAL:0::global_error_after: 0.06040494959
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::14.00000000 2.738612788 0.7500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.1839950180
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ0 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::73.00000000 7.697402159 1.125000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0.2440703024
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 2d:
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::56.00000000 5.830951895 1.000000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+DEAL:0::Checking DGQ1 against DGQ1 in 3d:
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::580.0000000 22.78157150 1.500000000
+DEAL:0::global_error_before: 0
+DEAL:0::global_error_after: 0
+DEAL:0::OK
+
+
+
+
diff --git a/tests/mpi/fe_tools_extrapolate_common.h b/tests/mpi/fe_tools_extrapolate_common.h
new file mode 100644 (file)
index 0000000..efc3c83
--- /dev/null
@@ -0,0 +1,363 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/base/std_cxx11/unique_ptr.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/fe_tools.h>
+
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <string>
+
+
+template <int dim>
+class  TestFunction : public Function<dim>
+{
+public:
+  TestFunction (unsigned int deg)
+    : Function<dim>(), degree(deg)
+  {}
+
+  virtual double value (const Point<dim>   &p,
+                        const unsigned int  component = 0) const
+  {
+    double return_value = 0.;
+    for (unsigned int d=0; d<dim; ++d)
+      return_value += std::pow(std::abs(.5-p(d)), degree);
+
+    return return_value;
+  }
+private:
+  const unsigned int degree;
+};
+
+
+
+template <int dim>
+parallel::distributed::Triangulation<dim> *make_tria ()
+{
+  parallel::distributed::Triangulation<dim> *tria = new parallel::distributed::Triangulation<dim>(MPI_COMM_WORLD);
+  typename parallel::distributed::Triangulation<dim>::active_cell_iterator cell;
+  GridGenerator::hyper_cube(*tria, 0., 1.);
+  tria->refine_global (2);
+  for (int i=0; i<2; ++i)
+    {
+      cell = tria->begin_active();
+      cell->set_refine_flag();
+      ++cell;
+      cell->set_refine_flag();
+
+      tria->execute_coarsening_and_refinement ();
+    }
+  return tria;
+}
+
+
+
+template <int dim>
+DoFHandler<dim> *make_dof_handler (const parallel::distributed::Triangulation<dim> &tria,
+                                   const FiniteElement<dim> &fe)
+{
+  DoFHandler<dim> *dof_handler = new DoFHandler<dim>(tria);
+  dof_handler->distribute_dofs (fe);
+  return dof_handler;
+}
+
+
+
+// output some indicators for a given vector
+template <unsigned int dim, typename VectorType>
+void
+output_vector (const VectorType &v, const std::string &output_name, const DoFHandler<dim> &dof_handler)
+{
+  //v.print(deallog.get_file_stream());
+
+  DataOut<dim> data_out;
+  data_out.attach_dof_handler (dof_handler);
+  data_out.add_data_vector (v, output_name, DataOut<dim>::type_dof_data);
+
+  const unsigned int degree = dof_handler.get_fe().degree;
+  data_out.build_patches (degree);
+
+  const std::string filename = (output_name +
+                                "." +
+                                Utilities::int_to_string
+                                (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD), 1));
+  std::ofstream output ((filename + ".vtu").c_str());
+  data_out.write_vtu (output);
+
+}
+
+
+
+template <int dim, typename VectorType>
+void
+check_this (const FiniteElement<dim> &fe1,
+            const FiniteElement<dim> &fe2)
+{
+  deallog << std::setprecision (10);
+
+  // only check if both elements have support points. 
+  // otherwise, interpolation doesn't really work
+  if ((fe1.get_unit_support_points().size() == 0) ||
+      (fe2.get_unit_support_points().size() == 0))
+    return;
+  //  likewise for non-primitive elements
+  if (!fe1.is_primitive() || !fe2.is_primitive())
+    return;
+  // we need to have dof_constraints for this test
+  if (!fe2.constraints_are_implemented())
+    return;
+  // we need prolongation matrices in fe2
+  if (!fe2.isotropic_restriction_is_implemented())
+    return;
+
+  std_cxx11::unique_ptr<parallel::distributed::Triangulation<dim> > tria(make_tria<dim>());
+
+  std_cxx11::unique_ptr<DoFHandler<dim> >    dof1(make_dof_handler (*tria, fe1));
+  std_cxx11::unique_ptr<DoFHandler<dim> >    dof2(make_dof_handler (*tria, fe2));
+  ConstraintMatrix cm1;
+  DoFTools::make_hanging_node_constraints (*dof1, cm1);
+  cm1.close ();
+  ConstraintMatrix cm2;
+  DoFTools::make_hanging_node_constraints (*dof2, cm2);
+  cm2.close ();
+
+  IndexSet locally_owned_dofs1 = dof1->locally_owned_dofs();
+  IndexSet locally_relevant_dofs1;
+  DoFTools::extract_locally_relevant_dofs (*dof1, locally_relevant_dofs1);
+  IndexSet locally_owned_dofs2 = dof2->locally_owned_dofs();
+  IndexSet locally_relevant_dofs2;
+  DoFTools::extract_locally_relevant_dofs (*dof2, locally_relevant_dofs2);
+
+  VectorType in_ghosted (locally_owned_dofs1, locally_relevant_dofs1, MPI_COMM_WORLD);
+  VectorType in_distributed (locally_owned_dofs1, MPI_COMM_WORLD);
+  VectorType out_distributed (locally_owned_dofs2, MPI_COMM_WORLD);
+  VectorType out_ghosted (locally_owned_dofs2, locally_relevant_dofs2, MPI_COMM_WORLD);
+  VectorType out_reference (locally_owned_dofs2, MPI_COMM_WORLD);
+
+  // Choose some reference function of sufficiently high polynomial degree.
+  TestFunction<dim> function (std::max(fe1.degree,fe2.degree));
+  VectorTools::interpolate (*dof1, function, in_distributed);
+  cm1.distribute(in_distributed);
+  in_ghosted = in_distributed;
+
+  Vector<double> difference_before;
+  VectorTools::integrate_difference (*dof1, in_ghosted, function, difference_before,
+                                     QGauss<dim>(fe1.degree+2), VectorTools::L2_norm);
+  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<dim, VectorType> (in_ghosted,
+                                  Utilities::int_to_string(fe1.degree,1)+Utilities::int_to_string(dim,1)+std::string("in"),
+                                  *dof1);*/
+  {
+    const double l1_norm = in_distributed.l1_norm();
+    const double l2_norm = in_distributed.l2_norm();
+    const double linfty_norm = in_distributed.linfty_norm();
+    if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+      deallog << l1_norm << ' '
+              << l2_norm << ' '
+              << linfty_norm << std::endl;
+  }
+
+  FETools::extrapolate_parallel (*dof1, in_ghosted, *dof2, cm2, out_distributed);
+  out_distributed.compress(VectorOperation::insert);
+  out_ghosted = out_distributed;
+
+  /*output_vector<dim, VectorType> (out_ghosted,
+                                  Utilities::int_to_string(fe2.degree,1)+Utilities::int_to_string(dim,1)+std::string("out"),
+                                  *dof2);*/
+
+  {
+    const double l1_norm = out_distributed.l1_norm();
+    const double l2_norm = out_distributed.l2_norm();
+    const double linfty_norm = out_distributed.linfty_norm();
+    if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+      deallog << l1_norm << ' '
+              << l2_norm << ' '
+              << linfty_norm << std::endl;
+  }
+
+  Vector<double> difference_after;
+  VectorTools::integrate_difference (*dof2, out_ghosted, function, difference_after,
+                                     QGauss<dim>(fe2.degree+2), VectorTools::L2_norm);
+  const double local_error_after = difference_after.l2_norm();
+  const double global_error_after
+    = std::sqrt(Utilities::MPI::sum(std::pow(local_error_after,2.), MPI_COMM_WORLD));
+
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+    {
+      deallog << "global_error_before: "
+              << global_error_before << std::endl;
+      deallog << "global_error_after: "
+              << global_error_after << std::endl;
+    }
+  if (fe2.degree == fe1.degree)
+    {
+      out_distributed -= in_distributed;
+      AssertThrow(out_distributed.l2_norm()<1.e-8,
+                  ExcInternalError());
+
+    }
+  if (fe2.degree > fe1.degree)
+    {
+      AssertThrow(global_error_after < global_error_before,
+                  ExcInternalError());
+    }
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+    deallog << "OK" << std::endl;
+}
+
+
+template <int dim, typename VectorType>
+void
+check_this_dealii (const FiniteElement<dim> &fe1,
+                   const FiniteElement<dim> &fe2)
+{
+  deallog << std::setprecision (10);
+
+  // only check if both elements have support points. 
+  // otherwise, interpolation doesn't really work
+  if ((fe1.get_unit_support_points().size() == 0) ||
+      (fe2.get_unit_support_points().size() == 0))
+    return;
+  //  likewise for non-primitive elements
+  if (!fe1.is_primitive() || !fe2.is_primitive())
+    return;
+  // we need to have dof_constraints for this test
+  if (!fe2.constraints_are_implemented())
+    return;
+  // we need prolongation matrices in fe2
+  if (!fe2.isotropic_restriction_is_implemented())
+    return;
+
+  std_cxx11::unique_ptr<parallel::distributed::Triangulation<dim> > tria(make_tria<dim>());
+
+  std_cxx11::unique_ptr<DoFHandler<dim> >    dof1(make_dof_handler (*tria, fe1));
+  std_cxx11::unique_ptr<DoFHandler<dim> >    dof2(make_dof_handler (*tria, fe2));
+  ConstraintMatrix cm1;
+  DoFTools::make_hanging_node_constraints (*dof1, cm1);
+  cm1.close ();
+  ConstraintMatrix cm2;
+  DoFTools::make_hanging_node_constraints (*dof2, cm2);
+  cm2.close ();
+
+  IndexSet locally_owned_dofs1 = dof1->locally_owned_dofs();
+  IndexSet locally_relevant_dofs1;
+  DoFTools::extract_locally_relevant_dofs (*dof1, locally_relevant_dofs1);
+  IndexSet locally_owned_dofs2 = dof2->locally_owned_dofs();
+  IndexSet locally_relevant_dofs2;
+  DoFTools::extract_locally_relevant_dofs (*dof2, locally_relevant_dofs2);
+
+  VectorType in_ghosted (locally_owned_dofs1, locally_relevant_dofs1, MPI_COMM_WORLD);
+  VectorType out_ghosted (locally_owned_dofs2, locally_relevant_dofs2, MPI_COMM_WORLD);
+  VectorType out_reference (locally_owned_dofs2, locally_relevant_dofs2, MPI_COMM_WORLD);
+
+  // Choose some reference function of sufficiently high polynomial degree.
+  TestFunction<dim> function (std::max(fe1.degree,fe2.degree));
+  VectorTools::interpolate (*dof1, function, in_ghosted);
+  cm1.distribute(in_ghosted);
+  in_ghosted.update_ghost_values();
+  VectorTools::interpolate (*dof2, function, out_reference);
+  cm2.distribute(out_reference);
+  out_reference.update_ghost_values();
+
+  Vector<double> difference_before;
+  VectorTools::integrate_difference (*dof1, in_ghosted, function, difference_before,
+                                     QGauss<dim>(fe1.degree+2), VectorTools::L2_norm);
+  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<dim, VectorType> (in_ghosted,
+                                  Utilities::int_to_string(fe1.degree,1)+Utilities::int_to_string(dim,1)+std::string("in"),
+                                  *dof1);*/
+  {
+    const double l1_norm = in_ghosted.l1_norm();
+    const double l2_norm = in_ghosted.l2_norm();
+    const double linfty_norm = in_ghosted.linfty_norm();
+    if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+      deallog << l1_norm << ' '
+              << l2_norm << ' '
+              << 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);
+
+  out_ghosted.update_ghost_values();
+
+  /*output_vector<dim, VectorType> (out_ghosted,
+                                  Utilities::int_to_string(fe2.degree,1)+Utilities::int_to_string(dim,1)+std::string("out"),
+                                  *dof2);*/
+  {
+    const double l1_norm = out_ghosted.l1_norm();
+    const double l2_norm = out_ghosted.l2_norm();
+    const double linfty_norm = out_ghosted.linfty_norm();
+    if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+      deallog << l1_norm << ' '
+              << l2_norm << ' '
+              << linfty_norm << std::endl;
+  }
+
+  Vector<double> difference_after;
+  VectorTools::integrate_difference (*dof2, out_ghosted, function, difference_after,
+                                     QGauss<dim>(fe2.degree+2), VectorTools::L2_norm);
+  const double local_error_after = difference_after.l2_norm();
+  const double global_error_after
+    = std::sqrt(Utilities::MPI::sum(std::pow(local_error_after,2.), MPI_COMM_WORLD));
+
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+    {
+      deallog << "global_error_before: "
+              << global_error_before << std::endl;
+      deallog << "global_error_after: "
+              << global_error_after << std::endl;
+    }
+  if (fe2.degree == fe1.degree)
+    {
+      out_ghosted -= in_ghosted;
+      AssertThrow(out_ghosted.l2_norm()<1.e-8,
+                  ExcInternalError());
+
+    }
+  if (fe2.degree > fe1.degree)
+    {
+      AssertThrow(global_error_after < global_error_before,
+                  ExcInternalError());
+    }
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD)==0)
+    deallog << "OK" << std::endl;
+}

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.