From: Giovanni Alzetta Date: Mon, 23 Oct 2017 13:58:31 +0000 (+0000) Subject: Added test for Functions::FEFieldFunction::compute_point_locations X-Git-Tag: v9.0.0-rc1~892^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b3cc35702f147e8746ba338c51bfdd16c166b12c;p=dealii.git Added test for Functions::FEFieldFunction::compute_point_locations --- diff --git a/doc/news/changes/minor/20171023GiovanniAlzetta b/doc/news/changes/minor/20171023GiovanniAlzetta new file mode 100644 index 0000000000..fdbf497def --- /dev/null +++ b/doc/news/changes/minor/20171023GiovanniAlzetta @@ -0,0 +1,3 @@ +Added: Test for Functions::FEFieldFunction::compute_point_locations +
+(Giovanni Alzetta, 2017/10/23) diff --git a/tests/fe/fe_compute_point_locations_01.cc b/tests/fe/fe_compute_point_locations_01.cc new file mode 100644 index 0000000000..4de0d5dfa0 --- /dev/null +++ b/tests/fe/fe_compute_point_locations_01.cc @@ -0,0 +1,113 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 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 Functions::FEFieldFunction::compute_point_locations + +#include "../tests.h" +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include +#include +#include + +#include + +using namespace dealii; + +template +void test_compute_pt_loc(unsigned int n_points) +{ + deallog << "Testing for dim = " << dim << std::endl; + deallog << "Testing on: " << n_points << " points." << std::endl; + + // Creating a grid in the square [0,1]x[0,1] + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(std::max(6-dim,2)); + + //Creating the finite elements needed: + FE_Q fe(1); + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + + // Creating the random points + std::vector> points; + + for (size_t i=0; i p; + for (unsigned int d=0; d::active_cell_iterator> cells; + std::vector > > qpoints; + std::vector > maps; + + // Creating a dummy vector/fe_field_function in order to use FEFieldFunction + Vector dummy; + Functions::FEFieldFunction fe_function(dof_handler, dummy, StaticMappingQ1::mapping); + size_t n_cells = fe_function.compute_point_locations(points, cells, qpoints, maps); + + deallog << "Points found in " << n_cells << " cells" << std::endl; + + // testing if the transformation is correct: + // For each cell check if the quadrature points in the i-th FE + // are the same as maps[i] + for (unsigned int i=0; i fev(fe, quad, update_quadrature_points); + fev.reinit(cell); + const auto &real_quad = fev.get_quadrature_points(); + + for (unsigned int q=0; q 1e-10) + deallog << "Error on cell : " << cell + << " at local point " << i + << ", corresponding to real point " + << points[local_map[q]] + << ", that got transformed to " + << real_quad[q] + << " instead." << std::endl; + } + } + deallog << "Test finished" << std::endl; +} + +int main() +{ + initlog(); + + deallog << "Deal.II compute_pt_loc:" << std::endl; + test_compute_pt_loc<2>(100); + test_compute_pt_loc<3>(200); +} diff --git a/tests/fe/fe_compute_point_locations_01.output b/tests/fe/fe_compute_point_locations_01.output new file mode 100644 index 0000000000..0a6e0e1554 --- /dev/null +++ b/tests/fe/fe_compute_point_locations_01.output @@ -0,0 +1,10 @@ + +DEAL::Deal.II compute_pt_loc: +DEAL::Testing for dim = 2 +DEAL::Testing on: 100 points. +DEAL::Points found in 81 cells +DEAL::Test finished +DEAL::Testing for dim = 3 +DEAL::Testing on: 200 points. +DEAL::Points found in 170 cells +DEAL::Test finished