From: Matthias Maier Date: Fri, 5 Apr 2024 03:19:25 +0000 (-0500) Subject: Test fe/fe_values_view_{07|10|13|16}: make Assert more robust X-Git-Tag: v9.6.0-rc1~406^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F16849%2Fhead;p=dealii.git Test fe/fe_values_view_{07|10|13|16}: make Assert more robust The Assert checks that we computed with a relative error of less than 1e-12. This fails however, if the vector in question is 0. Thus, also check for an absolute error. --- diff --git a/tests/fe/fe_values_view_07.cc b/tests/fe/fe_values_view_07.cc index 044bf2257c..217d646dc8 100644 --- a/tests/fe/fe_values_view_07.cc +++ b/tests/fe/fe_values_view_07.cc @@ -72,9 +72,10 @@ test(const Triangulation &tr, const FiniteElement &fe) for (const auto q : fe_values.quadrature_point_indices()) { deallog << scalar_values[q] << std::endl; - Assert((scalar_values[q] - vector_values[q][c]).norm() <= - 1e-12 * scalar_values[q].norm(), - ExcInternalError()); + const auto norm = (scalar_values[q] - vector_values[q][c]).norm(); + const auto tolerance = + std::max(1.e-11, 1.e-12 * scalar_values[q].norm()); + Assert(norm <= tolerance, ExcInternalError()); } } } diff --git a/tests/fe/fe_values_view_10.cc b/tests/fe/fe_values_view_10.cc index 68e1df2e96..d0592ae015 100644 --- a/tests/fe/fe_values_view_10.cc +++ b/tests/fe/fe_values_view_10.cc @@ -78,9 +78,10 @@ test(const Triangulation &tr, const FiniteElement &fe) deallog << scalar_values[q][d][e] << (d < dim - 1 || e < dim - 1 ? " " : ""); deallog << std::endl; - Assert((scalar_values[q] - vector_values[q][c]).norm() <= - 1e-12 * scalar_values[q].norm(), - ExcInternalError()); + const auto norm = (scalar_values[q] - vector_values[q][c]).norm(); + const auto tolerance = + std::max(1.e-11, 1.e-12 * scalar_values[q].norm()); + Assert(norm <= tolerance, ExcInternalError()); } } } diff --git a/tests/fe/fe_values_view_13.cc b/tests/fe/fe_values_view_13.cc index 26ee8fbbe4..50d9bafc07 100644 --- a/tests/fe/fe_values_view_13.cc +++ b/tests/fe/fe_values_view_13.cc @@ -78,9 +78,11 @@ test(const Triangulation &tr, const FiniteElement &fe) for (unsigned int d = 0; d < dim; ++d) { deallog << selected_vector_values[q][d] << std::endl; - Assert((selected_vector_values[q][d] - vector_values[q][c + d]) - .norm() <= 1e-12 * selected_vector_values[q][d].norm(), - ExcInternalError()); + const auto norm = + (selected_vector_values[q][d] - vector_values[q][c + d]).norm(); + const auto tolerance = + std::max(1.e-11, 1.e-12 * selected_vector_values[q][d].norm()); + Assert(norm <= tolerance, ExcInternalError()); } } } diff --git a/tests/fe/fe_values_view_16.cc b/tests/fe/fe_values_view_16.cc index 43ed89355c..fe281b40e3 100644 --- a/tests/fe/fe_values_view_16.cc +++ b/tests/fe/fe_values_view_16.cc @@ -83,9 +83,11 @@ test(const Triangulation &tr, const FiniteElement &fe) for (unsigned int f = 0; f < dim; ++f) deallog << selected_vector_values[q][d][e][f] << ' '; deallog << std::endl; - Assert((selected_vector_values[q][d] - vector_values[q][c + d]) - .norm() <= 1e-12 * selected_vector_values[q][d].norm(), - ExcInternalError()); + const auto norm = + (selected_vector_values[q][d] - vector_values[q][c + d]).norm(); + const auto tolerance = + std::max(1.e-11, 1.e-12 * selected_vector_values[q][d].norm()); + Assert(norm <= tolerance, ExcInternalError()); } } }