]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add test.
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 30 Jan 2018 10:19:01 +0000 (11:19 +0100)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 30 Jan 2018 10:19:01 +0000 (11:19 +0100)
tests/multigrid/transfer_matrix_free_08.cc
tests/multigrid/transfer_matrix_free_12.cc [new file with mode: 0644]
tests/multigrid/transfer_matrix_free_12.output [new file with mode: 0644]

index 472a0e6e8fffa3616e244aa35df02b3b2ad8be2b..bf079cf6848ce2c618e1b954e6dfe6b7fa80c432 100644 (file)
@@ -81,8 +81,9 @@ void check(const FiniteElement<dim> &fe)
   for (unsigned int level=1; level<tr.n_global_levels(); ++level)
     transfer.prolongate(level, vec[level], vec[level-1]);
   vec.back() -= vref;
+  const Number tolerance = 1000. * std::numeric_limits<Number>::epsilon();
   deallog << "Error after prolongation: "
-          << (double)vec.back().linfty_norm() << std::endl;
+          << filter_out_small_numbers(vec.back().linfty_norm(), tolerance) << std::endl;
 
   // unfortunately, no completely trivial expression of what should happen
   // during restriction
diff --git a/tests/multigrid/transfer_matrix_free_12.cc b/tests/multigrid/transfer_matrix_free_12.cc
new file mode 100644 (file)
index 0000000..dd6f9fa
--- /dev/null
@@ -0,0 +1,116 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2018 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Similar to transfer_matrix_free_08 but using a system of equations: high
+// polynomial degrees beyond 10 that are compared against a linear function
+// that gets prolongated between the mesh levels. Since the function is
+// linear, it should be exactly represented on the finest level. Then also
+// look into the result of restriction
+
+#include "../tests.h"
+#include <deal.II/base/function_lib.h>
+#include <deal.II/lac/la_parallel_vector.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/multigrid/mg_transfer.h>
+#include <deal.II/multigrid/mg_transfer_matrix_free.h>
+#include <deal.II/numerics/vector_tools.h>
+
+
+template <int dim, typename Number>
+void check(const FiniteElement<dim> &fe_scalar)
+{
+  FESystem<dim> fe(fe_scalar, dim);
+  deallog << "FE: " << fe.get_name() << std::endl;
+
+  // run a few different sizes...
+
+  Triangulation<dim> tr(Triangulation<dim>::limit_level_difference_at_vertices);
+  GridGenerator::hyper_cube(tr);
+  tr.refine_global(2);
+
+  Triangulation<dim> trcoarse(Triangulation<dim>::limit_level_difference_at_vertices);
+  GridGenerator::hyper_cube(trcoarse);
+  DoFHandler<dim> dofcoarse(trcoarse);
+  dofcoarse.distribute_dofs(fe);
+  Tensor<1,dim> exponents_monomial;
+  for (unsigned int d=0; d<dim; ++d)
+    exponents_monomial[d] = 1;
+  LinearAlgebra::distributed::Vector<double> vrefcoarse;
+  vrefcoarse.reinit(dofcoarse.n_dofs());
+  VectorTools::interpolate(dofcoarse, Functions::Monomial<dim>(exponents_monomial, dim), vrefcoarse);
+
+  deallog << "no. cells: " << tr.n_global_active_cells() << std::endl;
+
+  DoFHandler<dim> mgdof(tr);
+  mgdof.distribute_dofs(fe);
+  mgdof.distribute_mg_dofs(fe);
+
+  // build matrix-free transfer
+  MGTransferMatrixFree<dim, Number> transfer;
+  transfer.build(mgdof);
+
+  LinearAlgebra::distributed::Vector<double> vrefdouble;
+  LinearAlgebra::distributed::Vector<Number> vref;
+  AssertDimension(mgdof.n_dofs(tr.n_global_levels()-1), mgdof.n_dofs());
+  vrefdouble.reinit(mgdof.n_dofs());
+  VectorTools::interpolate(mgdof, Functions::Monomial<dim>(exponents_monomial, dim), vrefdouble);
+
+  vref.reinit(mgdof.n_dofs());
+  vref = vrefdouble;
+  std::vector<LinearAlgebra::distributed::Vector<Number> > vec(tr.n_global_levels());
+  for (unsigned int level=0; level<tr.n_global_levels(); ++level)
+    vec[level].reinit(mgdof.n_dofs(level));
+  vec.back() = vref;
+
+  // prolongate monomial from coarse to fine, should be exact for monomial
+  vec[0] = vrefcoarse;
+  for (unsigned int level=1; level<tr.n_global_levels(); ++level)
+    transfer.prolongate(level, vec[level], vec[level-1]);
+  vec.back() -= vref;
+  const Number tolerance = 1000. * std::numeric_limits<Number>::epsilon();
+  deallog << "Error after prolongation: "
+          << filter_out_small_numbers(vec.back().linfty_norm(), tolerance) << std::endl;
+
+  // unfortunately, no completely trivial expression of what should happen
+  // during restriction
+  vec.back() = vref;
+  for (unsigned int level=tr.n_global_levels()-1; level > 0 ; --level)
+    transfer.restrict_and_add(level, vec[level-1], vec[level]);
+  deallog << "Norm after restriction: " << vec[0].l2_norm() << std::endl;
+}
+
+
+int main(int argc, char **argv)
+{
+  // no threading in this test...
+  Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+  initlog();
+
+  check<2,double>(FE_DGQ<2>(3));
+  check<2,double>(FE_DGQ<2>(8));
+  check<2,double>(FE_DGQ<2>(16));
+  check<2,double>(FE_Q<2>(16));
+  check<3,double>(FE_DGQ<3>(11));
+  check<2,float> (FE_DGQ<2>(3));
+  check<2,float> (FE_DGQ<2>(8));
+  check<2,float> (FE_DGQ<2>(16));
+  check<2,float> (FE_Q<2>(16));
+  check<3,float> (FE_DGQ<3>(11));
+}
diff --git a/tests/multigrid/transfer_matrix_free_12.output b/tests/multigrid/transfer_matrix_free_12.output
new file mode 100644 (file)
index 0000000..b3cf0ec
--- /dev/null
@@ -0,0 +1,41 @@
+
+DEAL::FE: FESystem<2>[FE_DGQ<2>(3)^2]
+DEAL::no. cells: 16
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 48.2543
+DEAL::FE: FESystem<2>[FE_DGQ<2>(8)^2]
+DEAL::no. cells: 16
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 101.354
+DEAL::FE: FESystem<2>[FE_DGQ<2>(16)^2]
+DEAL::no. cells: 16
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 214.405
+DEAL::FE: FESystem<2>[FE_Q<2>(16)^2]
+DEAL::no. cells: 16
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 189.129
+DEAL::FE: FESystem<3>[FE_DGQ<3>(11)^3]
+DEAL::no. cells: 64
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 1340.82
+DEAL::FE: FESystem<2>[FE_DGQ<2>(3)^2]
+DEAL::no. cells: 16
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 48.2544
+DEAL::FE: FESystem<2>[FE_DGQ<2>(8)^2]
+DEAL::no. cells: 16
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 101.354
+DEAL::FE: FESystem<2>[FE_DGQ<2>(16)^2]
+DEAL::no. cells: 16
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 214.405
+DEAL::FE: FESystem<2>[FE_Q<2>(16)^2]
+DEAL::no. cells: 16
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 189.129
+DEAL::FE: FESystem<3>[FE_DGQ<3>(11)^3]
+DEAL::no. cells: 64
+DEAL::Error after prolongation: 0.00000
+DEAL::Norm after restriction: 1340.82

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.