From e3f3efb3bcfb1c0b89e5ef1ddbf4975f96b0f44a Mon Sep 17 00:00:00 2001 From: David Wells Date: Sat, 28 Jul 2018 13:00:35 -0400 Subject: [PATCH] Use a C++17 function instead of jn. This is not guaranteed to be in the math library, but it is in C++17 (or in boost). --- cmake/checks/check_01_cxx_features.cmake | 17 ++++++ cmake/checks/check_02_system_features.cmake | 18 ------ include/deal.II/base/std_cxx17/cmath.h | 44 ++++++++++++++ source/base/function_lib.cc | 44 +++++--------- tests/base/functions_01.cc | 66 +++++++++++---------- tests/base/functions_01.output | 6 ++ 6 files changed, 116 insertions(+), 79 deletions(-) create mode 100644 include/deal.II/base/std_cxx17/cmath.h diff --git a/cmake/checks/check_01_cxx_features.cmake b/cmake/checks/check_01_cxx_features.cmake index 8ef3144ec6..32d61d3883 100644 --- a/cmake/checks/check_01_cxx_features.cmake +++ b/cmake/checks/check_01_cxx_features.cmake @@ -176,6 +176,22 @@ IF(NOT DEFINED DEAL_II_WITH_CXX17 OR DEAL_II_WITH_CXX17) " DEAL_II_HAVE_CXX17_IF_CONSTEXPR) + # + # Not all compilers with C++17 support include the new special math functions: + # + CHECK_CXX_SOURCE_COMPILES( + " + #include + + int main() + { + std::cyl_bessel_j(1.0, 1.0); + std::cyl_bessel_jf(1.0f, 1.0f); + std::cyl_bessel_jl(1.0, 1.0); + } + " + DEAL_II_HAVE_CXX17_SPECIAL_MATH_FUNCTIONS) + # # Some compilers treat lambdas as constexpr functions when compiling # with C++17 support even if they don't fulfill all the constexpr @@ -211,6 +227,7 @@ IF(NOT DEFINED DEAL_II_WITH_CXX17 OR DEAL_II_WITH_CXX17) IF( DEAL_II_HAVE_CXX17_ATTRIBUTES AND DEAL_II_HAVE_CXX17_IF_CONSTEXPR AND + DEAL_II_HAVE_CXX17_SPECIAL_MATH_FUNCTIONS AND DEAL_II_NON_CONSTEXPR_LAMBDA) SET(DEAL_II_HAVE_CXX17 TRUE) ELSE() diff --git a/cmake/checks/check_02_system_features.cmake b/cmake/checks/check_02_system_features.cmake index 6186752b43..4f0f832786 100644 --- a/cmake/checks/check_02_system_features.cmake +++ b/cmake/checks/check_02_system_features.cmake @@ -18,7 +18,6 @@ # # DEAL_II_HAVE_GETHOSTNAME # DEAL_II_HAVE_GETPID -# DEAL_II_HAVE_JN # DEAL_II_HAVE_SYS_RESOURCE_H # DEAL_II_HAVE_UNISTD_H # DEAL_II_MSVC @@ -40,23 +39,6 @@ CHECK_INCLUDE_FILE_CXX("unistd.h" DEAL_II_HAVE_UNISTD_H) CHECK_CXX_SYMBOL_EXISTS("gethostname" "unistd.h" DEAL_II_HAVE_GETHOSTNAME) CHECK_CXX_SYMBOL_EXISTS("getpid" "unistd.h" DEAL_II_HAVE_GETPID) -# -# Do we have the Bessel function jn? -# -FIND_SYSTEM_LIBRARY(m_LIBRARY NAMES m) -MARK_AS_ADVANCED(m_LIBRARY) - -IF(NOT m_LIBRARY MATCHES "-NOTFOUND") - LIST(APPEND CMAKE_REQUIRED_LIBRARIES ${m_LIBRARY}) - CHECK_CXX_SYMBOL_EXISTS("jn" "math.h" DEAL_II_HAVE_JN) - RESET_CMAKE_REQUIRED() - IF(DEAL_II_HAVE_JN) - LIST(APPEND DEAL_II_LIBRARIES ${m_LIBRARY}) - ENDIF() -ENDIF() - - - ######################################################################## # # # Mac OSX specific setup: # diff --git a/include/deal.II/base/std_cxx17/cmath.h b/include/deal.II/base/std_cxx17/cmath.h new file mode 100644 index 0000000000..1245a0be17 --- /dev/null +++ b/include/deal.II/base/std_cxx17/cmath.h @@ -0,0 +1,44 @@ +// --------------------------------------------------------------------- +// +// 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- +#ifndef dealii_cxx17_cmath_h +#define dealii_cxx17_cmath_h + +#include + +#ifndef DEAL_II_WITH_CXX17 +# include +#endif // DEAL_II_WITH_CXX17 + +#include + +DEAL_II_NAMESPACE_OPEN +namespace std_cxx17 +{ +#ifndef DEAL_II_WITH_CXX17 + double (&cyl_bessel_j)(double, + double) = boost::math::cyl_bessel_j; + float (&cyl_bessel_jf)(float, + float) = boost::math::cyl_bessel_j; + long double (&cyl_bessel_jl)(long double, long double) = + boost::math::cyl_bessel_j; +#else + using std::cyl_bessel_j; + using std::cyl_bessel_jf; + using std::cyl_bessel_jl; +#endif +} // namespace std_cxx17 +DEAL_II_NAMESPACE_CLOSE + +#endif // dealii_cxx17_cmath_h diff --git a/source/base/function_lib.cc b/source/base/function_lib.cc index 02e6d62155..465cde3112 100644 --- a/source/base/function_lib.cc +++ b/source/base/function_lib.cc @@ -17,6 +17,7 @@ #include #include #include +#include #include #include @@ -2283,12 +2284,7 @@ namespace Functions { Assert(dim == 2, ExcNotImplemented()); const double r = p.distance(center); -#ifndef DEAL_II_HAVE_JN - Assert(false, ExcMessage("The Bessel function jn was not found by CMake.")); - return r; -#else - return jn(order, r * wave_number); -#endif + return std_cxx17::cyl_bessel_j(order, r * wave_number); } @@ -2300,17 +2296,11 @@ namespace Functions { Assert(dim == 2, ExcNotImplemented()); AssertDimension(points.size(), values.size()); -#ifndef DEAL_II_HAVE_JN - (void)points; - (void)values; - Assert(false, ExcMessage("The Bessel function jn was not found by CMake.")); -#else for (unsigned int k = 0; k < points.size(); ++k) { const double r = points[k].distance(center); - values[k] = jn(order, r * wave_number); + values[k] = std_cxx17::cyl_bessel_j(order, r * wave_number); } -#endif } @@ -2319,23 +2309,19 @@ namespace Functions Bessel1::gradient(const Point &p, const unsigned int) const { Assert(dim == 2, ExcNotImplemented()); -#ifndef DEAL_II_HAVE_JN - (void)p; - Assert(false, ExcMessage("The Bessel function jn was not found by CMake.")); - return Tensor<1, dim>(); -#else const double r = p.distance(center); const double co = (r == 0.) ? 0. : (p(0) - center(0)) / r; const double si = (r == 0.) ? 0. : (p(1) - center(1)) / r; - const double dJn = (order == 0) ? (-jn(1, r * wave_number)) : - (.5 * (jn(order - 1, wave_number * r) - - jn(order + 1, wave_number * r))); + const double dJn = + (order == 0) ? + (-std_cxx17::cyl_bessel_j(1, r * wave_number)) : + (.5 * (std_cxx17::cyl_bessel_j(order - 1, wave_number * r) - + std_cxx17::cyl_bessel_j(order + 1, wave_number * r))); Tensor<1, dim> result; result[0] = wave_number * co * dJn; result[1] = wave_number * si * dJn; return result; -#endif } @@ -2348,9 +2334,6 @@ namespace Functions { Assert(dim == 2, ExcNotImplemented()); AssertDimension(points.size(), gradients.size()); -#ifndef DEAL_II_HAVE_JN - Assert(false, ExcMessage("The Bessel function jn was not found by CMake.")); -#else for (unsigned int k = 0; k < points.size(); ++k) { const Point &p = points[k]; @@ -2358,15 +2341,15 @@ namespace Functions const double co = (r == 0.) ? 0. : (p(0) - center(0)) / r; const double si = (r == 0.) ? 0. : (p(1) - center(1)) / r; - const double dJn = (order == 0) ? - (-jn(1, r * wave_number)) : - (.5 * (jn(order - 1, wave_number * r) - - jn(order + 1, wave_number * r))); + const double dJn = + (order == 0) ? + (-std_cxx17::cyl_bessel_j(1, r * wave_number)) : + (.5 * (std_cxx17::cyl_bessel_j(order - 1, wave_number * r) - + std_cxx17::cyl_bessel_j(order + 1, wave_number * r))); Tensor<1, dim> &result = gradients[k]; result[0] = wave_number * co * dJn; result[1] = wave_number * si * dJn; } -#endif } @@ -2837,6 +2820,7 @@ namespace Functions template class Monomial<1>; template class Monomial<2>; template class Monomial<3>; + template class Bessel1<1>; template class Bessel1<2>; template class Bessel1<3>; template class InterpolatedTensorProductGridData<1>; diff --git a/tests/base/functions_01.cc b/tests/base/functions_01.cc index a3c269fdff..c0ca0e1dcc 100644 --- a/tests/base/functions_01.cc +++ b/tests/base/functions_01.cc @@ -19,54 +19,58 @@ #include #include #include +#include #include #include #include #include +#include + #include +#include #include #include "../tests.h" #include "functions.h" -#define CHECK(F) \ - { \ - deallog << #F << std::endl; \ - F f; \ - check_function_value_consistency(f, 5); \ - check_function_gradient_consistency(f, 5); \ - check_gradient(f, 5); \ - } - +template +void +check(Args... args) +{ + FunctionType f(args...); + std::string function_name = boost::core::demangle(typeid(f).name()); + // remove the preceding dealii:: to match older output + function_name.erase(0, sizeof("dealii::") - 1 /*remember the NUL*/); + deallog << function_name << std::endl; + check_function_value_consistency(f, 5); + check_function_gradient_consistency(f, 5); + check_gradient(f, 5); +} -#define CHECKN(F, arg) \ - { \ - deallog << #F << '(' << arg << ')' << std::endl; \ - F f(arg); \ - check_function_value_consistency(f, arg + 1); \ - check_function_gradient_consistency(f, arg + 1); \ - } int main() { - std::string logname = "output"; - std::ofstream logfile(logname.c_str()); - deallog.attach(logfile); + initlog(); - CHECK(Functions::SquareFunction<1>); - CHECK(Functions::SquareFunction<2>); - CHECK(Functions::SquareFunction<3>); - CHECK(Functions::CosineFunction<1>); - CHECK(Functions::CosineFunction<2>); - CHECK(Functions::CosineFunction<3>); - CHECK(Functions::CosineGradFunction<1>); - CHECK(Functions::CosineGradFunction<2>); - CHECK(Functions::CosineGradFunction<3>); - CHECK(Functions::ExpFunction<1>); - CHECK(Functions::ExpFunction<2>); - CHECK(Functions::ExpFunction<3>); + check>(); + check>(); + check>(); + check>(); + check>(); + check>(); + check>(); + check>(); + check>(); + check>(); + check>(); + check>(); + // Bessel1 is only implemented in 2D right now + // check>(3, 2.0, Point<1>(1.0)); + check>(3, 2.0, Point<2>(1.0, 2.0)); + // check>(3, 2.0, Point<2>(1.0, 2.0)); + // check(Functions::Bessel1<3>); } diff --git a/tests/base/functions_01.output b/tests/base/functions_01.output index 867ad3cb31..ddca8cb87d 100644 --- a/tests/base/functions_01.output +++ b/tests/base/functions_01.output @@ -71,3 +71,9 @@ DEAL::value list vs vector value list DEAL::gradient vs vector gradient list DEAL::gradient list vs vector gradient list DEAL::gradients vs difference quotients +DEAL::Functions::Bessel1<2> +DEAL::value vs vector value list +DEAL::value list vs vector value list +DEAL::gradient vs vector gradient list +DEAL::gradient list vs vector gradient list +DEAL::gradients vs difference quotients -- 2.39.5