From 871806b7c5e3372a36fbdd08002b57f1b4372dbe Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Thu, 30 Mar 2017 14:35:03 +0200 Subject: [PATCH] FEEvaluation: Fix collocation evaluation path --- .../deal.II/matrix_free/evaluation_kernels.h | 2 +- .../matrix_free/get_functions_variants_gl.cc | 209 ++++++++++++++++++ .../get_functions_variants_gl.output | 64 ++++++ 3 files changed, 274 insertions(+), 1 deletion(-) create mode 100644 tests/matrix_free/get_functions_variants_gl.cc create mode 100644 tests/matrix_free/get_functions_variants_gl.output diff --git a/include/deal.II/matrix_free/evaluation_kernels.h b/include/deal.II/matrix_free/evaluation_kernels.h index 54b374b7c6..2dc040df16 100644 --- a/include/deal.II/matrix_free/evaluation_kernels.h +++ b/include/deal.II/matrix_free/evaluation_kernels.h @@ -807,7 +807,7 @@ namespace internal } if (evaluate_hessians == true) { - eval.template hessians<0,true,false> (values_quad[c], hessians_quad[c][0]); + eval.template hessians<0,true,false> (values_dofs[c], hessians_quad[c][0]); if (dim > 1) { // re-use grad_x already in gradients diff --git a/tests/matrix_free/get_functions_variants_gl.cc b/tests/matrix_free/get_functions_variants_gl.cc new file mode 100644 index 0000000000..93847b5436 --- /dev/null +++ b/tests/matrix_free/get_functions_variants_gl.cc @@ -0,0 +1,209 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2013 - 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. +// +// --------------------------------------------------------------------- + + + +// same as get_functions_variants for testing values, gradients, and +// Laplacians individually but now on Gauss-Lobatto elements + +#include "../tests.h" + +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +std::ofstream logfile("output"); + + +template +class MatrixFreeTest +{ +public: + typedef Vector VectorType; + + MatrixFreeTest(const MatrixFree &data_in): + data (data_in) + {}; + + void operator () (const MatrixFree &data, + VectorType &dst, + const VectorType &src, + const std::pair &cell_range) const; + + void test_functions (const VectorType &src) const + { + for (unsigned int i=0; i<5; ++i) + errors[i] = 0; + data.cell_loop (&MatrixFreeTest::operator(), this, + const_cast(src), src); + + deallog << "Error val, function values alone: " + << errors[0] << std::endl; + deallog << "Error grad, function gradients alone: " + << errors[1] << std::endl; + deallog << "Error val, function values and gradients alone: " + << errors[2] << std::endl; + deallog << "Error grad, function values and gradients alone: " + << errors[3] << std::endl; + deallog << "Error Lapl, function Laplacians alone: " + << errors[4] << std::endl; + }; + +private: + const MatrixFree &data; + mutable double errors [5]; +}; + + + + +template +void MatrixFreeTest:: +operator () (const MatrixFree &data, + VectorType &, + const VectorType &src, + const std::pair &cell_range) const +{ + FEEvaluation fe_eval (data); + FEEvaluation fe_eval2 (data); + FEEvaluation fe_eval3 (data); + FEEvaluation fe_eval4 (data); + FEEvaluation fe_eval5 (data); + for (unsigned int cell=cell_range.first; cell::n_array_elements; ++j) + { + errors[0] += std::fabs(fe_eval.get_value(q)[j]- + fe_eval2.get_value(q)[j]); + errors[2] += std::fabs(fe_eval.get_value(q)[j]- + fe_eval4.get_value(q)[j]); + for (unsigned int d=0; d +void test () +{ + Triangulation tria; + GridGenerator::hyper_cube (tria); + tria.refine_global(1); + + FE_Q fe (fe_degree); + DoFHandler dof (tria); + dof.distribute_dofs(fe); + + deallog << "Testing " << fe.get_name() << std::endl; + //std::cout << "Number of degrees of freedom: " << dof.n_dofs() << std::endl; + + Vector solution_dist (dof.n_dofs()); + + // create vector with random entries + for (unsigned int i=0; i mf_data; + { + const QGaussLobatto<1> quad (fe_degree+1); + typename MatrixFree::AdditionalData data; + data.tasks_parallel_scheme = MatrixFree::AdditionalData::none; + data.mapping_update_flags = update_gradients | update_hessians; + mf_data.reinit (dof, constraints, quad, data); + } + + MatrixFreeTest mf (mf_data); + mf.test_functions(solution_dist); + deallog << std::endl; +} + + +int main () +{ + deallog.attach(logfile); + + deallog.threshold_double(1.e-14); + deallog << std::setprecision (3); + + { + deallog.push("2d"); + test<2,1>(); + test<2,2>(); + test<2,3>(); + test<2,4>(); + test<2,5>(); + deallog.pop(); + deallog.push("3d"); + test<3,1>(); + test<3,2>(); + test<3,3>(); + test<3,4>(); + deallog.pop(); + } +} diff --git a/tests/matrix_free/get_functions_variants_gl.output b/tests/matrix_free/get_functions_variants_gl.output new file mode 100644 index 0000000000..9ac254fb9f --- /dev/null +++ b/tests/matrix_free/get_functions_variants_gl.output @@ -0,0 +1,64 @@ + +DEAL:2d::Testing FE_Q<2>(1) +DEAL:2d::Error val, function values alone: 0 +DEAL:2d::Error grad, function gradients alone: 0 +DEAL:2d::Error val, function values and gradients alone: 0 +DEAL:2d::Error grad, function values and gradients alone: 0 +DEAL:2d::Error Lapl, function Laplacians alone: 0 +DEAL:2d:: +DEAL:2d::Testing FE_Q<2>(2) +DEAL:2d::Error val, function values alone: 0 +DEAL:2d::Error grad, function gradients alone: 0 +DEAL:2d::Error val, function values and gradients alone: 0 +DEAL:2d::Error grad, function values and gradients alone: 0 +DEAL:2d::Error Lapl, function Laplacians alone: 0 +DEAL:2d:: +DEAL:2d::Testing FE_Q<2>(3) +DEAL:2d::Error val, function values alone: 0 +DEAL:2d::Error grad, function gradients alone: 0 +DEAL:2d::Error val, function values and gradients alone: 0 +DEAL:2d::Error grad, function values and gradients alone: 0 +DEAL:2d::Error Lapl, function Laplacians alone: 0 +DEAL:2d:: +DEAL:2d::Testing FE_Q<2>(4) +DEAL:2d::Error val, function values alone: 0 +DEAL:2d::Error grad, function gradients alone: 0 +DEAL:2d::Error val, function values and gradients alone: 0 +DEAL:2d::Error grad, function values and gradients alone: 0 +DEAL:2d::Error Lapl, function Laplacians alone: 0 +DEAL:2d:: +DEAL:2d::Testing FE_Q<2>(5) +DEAL:2d::Error val, function values alone: 0 +DEAL:2d::Error grad, function gradients alone: 0 +DEAL:2d::Error val, function values and gradients alone: 0 +DEAL:2d::Error grad, function values and gradients alone: 0 +DEAL:2d::Error Lapl, function Laplacians alone: 0 +DEAL:2d:: +DEAL:3d::Testing FE_Q<3>(1) +DEAL:3d::Error val, function values alone: 0 +DEAL:3d::Error grad, function gradients alone: 0 +DEAL:3d::Error val, function values and gradients alone: 0 +DEAL:3d::Error grad, function values and gradients alone: 0 +DEAL:3d::Error Lapl, function Laplacians alone: 0 +DEAL:3d:: +DEAL:3d::Testing FE_Q<3>(2) +DEAL:3d::Error val, function values alone: 0 +DEAL:3d::Error grad, function gradients alone: 0 +DEAL:3d::Error val, function values and gradients alone: 0 +DEAL:3d::Error grad, function values and gradients alone: 0 +DEAL:3d::Error Lapl, function Laplacians alone: 0 +DEAL:3d:: +DEAL:3d::Testing FE_Q<3>(3) +DEAL:3d::Error val, function values alone: 0 +DEAL:3d::Error grad, function gradients alone: 0 +DEAL:3d::Error val, function values and gradients alone: 0 +DEAL:3d::Error grad, function values and gradients alone: 0 +DEAL:3d::Error Lapl, function Laplacians alone: 0 +DEAL:3d:: +DEAL:3d::Testing FE_Q<3>(4) +DEAL:3d::Error val, function values alone: 0 +DEAL:3d::Error grad, function gradients alone: 0 +DEAL:3d::Error val, function values and gradients alone: 0 +DEAL:3d::Error grad, function values and gradients alone: 0 +DEAL:3d::Error Lapl, function Laplacians alone: 0 +DEAL:3d:: -- 2.39.5