From: Martin Kronbichler Date: Sat, 3 Jun 2017 19:40:11 +0000 (+0200) Subject: Add test where inversion in real_to_unit_cell fails X-Git-Tag: v9.0.0-rc1~1533^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=00f4be9942f8492f0a0663c439911dd100bcd838;p=dealii.git Add test where inversion in real_to_unit_cell fails --- diff --git a/tests/mappings/mapping_real_to_unit_q1_singular.cc b/tests/mappings/mapping_real_to_unit_q1_singular.cc new file mode 100644 index 0000000000..717a3c4d06 --- /dev/null +++ b/tests/mappings/mapping_real_to_unit_q1_singular.cc @@ -0,0 +1,68 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2006 - 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. +// +// --------------------------------------------------------------------- + + + +// check that MappingQ1::transform_real_to_unit_cell can handle the case of a +// singular discriminant. Previously we used to have a division by zero + +#include "../tests.h" + +#include +#include +#include +#include + + +template +void test_real_to_unit_cell() +{ + deallog << "dim=" << dim << std::endl; + + Triangulation triangulation; + GridGenerator::hyper_ball (triangulation); + + Point point; + MappingQGeneric mapping(1); + + point[1] = -1./(1+std::sqrt(2.0))/std::sqrt(2); + + // check on cell 2 + typename Triangulation::cell_iterator cell = triangulation.begin(); + ++cell; + try + { + mapping.transform_real_to_unit_cell(cell, point); + } + catch (typename Mapping::ExcTransformationFailed) + { + deallog << "Transformation for point " << point << " on cell with " + << "center " << cell->center() << " is not invertible" << std::endl; + } + deallog << "OK" << std::endl; +} + + +int +main() +{ + std::ofstream logfile ("output"); + deallog.attach(logfile); + deallog.threshold_double(1.e-10); + + test_real_to_unit_cell<2>(); + + return 0; +} diff --git a/tests/mappings/mapping_real_to_unit_q1_singular.output b/tests/mappings/mapping_real_to_unit_q1_singular.output new file mode 100644 index 0000000000..bdc167b7f5 --- /dev/null +++ b/tests/mappings/mapping_real_to_unit_q1_singular.output @@ -0,0 +1,4 @@ + +DEAL::dim=2 +DEAL::Transformation for point 0.00000 -0.292893 on cell with center -0.500000 -1.38778e-17 is not invertible +DEAL::OK