diff mbox series

[libstdc++;,pr66689;,pr68397] Backport special function fixes to 7.

Message ID 0b3064f7-99c3-41e0-26c6-1c4fb3391d34@verizon.net
State New
Headers show
Series [libstdc++;,pr66689;,pr68397] Backport special function fixes to 7. | expand

Commit Message

Ed Smith-Rowland April 27, 2018, 3:17 a.m. UTC
I'm thinking of going back on my choice not to fix special function bugs 
on 7.

These build and test clean on gcc-7.

OK?
2018-04-30  Edward Smith-Rowland  <3dw4rd@verizon.net>

	PR libstdc++/pr66689 - comp_ellint_3 and ellint_3 return garbage values
	* include/tr1/ell_integral.tcc: Correct the nu sign convention
	in ellint_3 and comp_ellint_3.
	* testsuite/tr1/5_numerical_facilities/special_functions/
	06_comp_ellint_3/check_value.cc: Regen with correct values.
	* testsuite/tr1/5_numerical_facilities/special_functions/
	14_ellint_3/check_value.cc: Ditto.
	* testsuite/special_functions/06_comp_ellint_3/check_value.cc: Ditto.
	* testsuite/special_functions/13_ellint_3/check_value.cc: Ditto.
	* testsuite/special_functions/06_comp_ellint_3/pr66689.cc: New.
	* testsuite/special_functions/13_ellint_3/pr66689.cc: New.
	* testsuite/tr1/5_numerical_facilities/special_functions/
	06_comp_ellint_3/pr66689.cc: New.
	* testsuite/tr1/5_numerical_facilities/special_functions/
	14_ellint_3/pr66689.cc: New.

2018-04-30  Edward Smith-Rowland  <3dw4rd@verizon.net>

	PR libstdc++/68397 std::tr1::expint fails ... long double arguments.
	* include/tr1/exp_integral.tcc: Increase iteration limits.
	* testsuite/tr1/5_numerical_facilities/special_functions/15_expint/
	pr68397.cc: New test.
	* testsuite/special_functions/14_expint/pr68397.cc: New test.
Index: include/tr1/exp_integral.tcc
===================================================================
--- include/tr1/exp_integral.tcc	(revision 246797)
+++ include/tr1/exp_integral.tcc	(working copy)
@@ -86,7 +86,7 @@
       _Tp __term = _Tp(1);
       _Tp __esum = _Tp(0);
       _Tp __osum = _Tp(0);
-      const unsigned int __max_iter = 100;
+      const unsigned int __max_iter = 1000;
       for (unsigned int __i = 1; __i < __max_iter; ++__i)
         {
           __term *= - __x / __i;
@@ -156,7 +156,7 @@
     _Tp
     __expint_En_series(unsigned int __n, _Tp __x)
     {
-      const unsigned int __max_iter = 100;
+      const unsigned int __max_iter = 1000;
       const _Tp __eps = std::numeric_limits<_Tp>::epsilon();
       const int __nm1 = __n - 1;
       _Tp __ans = (__nm1 != 0
@@ -202,7 +202,7 @@
     _Tp
     __expint_En_cont_frac(unsigned int __n, _Tp __x)
     {
-      const unsigned int __max_iter = 100;
+      const unsigned int __max_iter = 1000;
       const _Tp __eps = std::numeric_limits<_Tp>::epsilon();
       const _Tp __fp_min = std::numeric_limits<_Tp>::min();
       const int __nm1 = __n - 1;
Index: testsuite/tr1/5_numerical_facilities/special_functions/15_expint/pr68397.cc
===================================================================
--- testsuite/tr1/5_numerical_facilities/special_functions/15_expint/pr68397.cc	(nonexistent)
+++ testsuite/tr1/5_numerical_facilities/special_functions/15_expint/pr68397.cc	(working copy)
@@ -0,0 +1,46 @@
+// Copyright (C) 2017 Free Software Foundation, Inc.
+//
+// This file is part of the GNU ISO C++ Library.  This library is free
+// software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the
+// Free Software Foundation; either version 3, or (at your option)
+// any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License along
+// with this library; see the file COPYING3.  If not see
+// <http://www.gnu.org/licenses/>.
+
+// PR libstdc++/68397 -  std::tr1::expint fails in __expint_En_cont_frac
+// for some long double arguments due to low __max_iter value
+
+#include <tr1/cmath>
+#include <testsuite_hooks.h>
+
+void
+test01()
+{
+  // Answers from Wolfram Alpha.
+  long double ans_ok = -0.10001943365331651406888645149537315243646135979573L;
+  long double ans_bomb = -0.10777727809650077516264612749163100483995270163783L;
+
+  long double Ei_ok = std::tr1::expint(-1.500001L);
+  long double diff_ok = std::abs(Ei_ok - ans_ok);
+  VERIFY(diff_ok < 1.0e-15L);
+
+  long double Ei_bomb = std::tr1::expint(-1.450001L);
+  long double diff_bomb = std::abs(Ei_bomb - ans_bomb);
+  VERIFY(diff_bomb < 1.0e-15L);
+}
+
+int
+main()
+{
+  test01();
+  return 0;
+}
+
Index: testsuite/special_functions/14_expint/pr68397.cc
===================================================================
--- testsuite/special_functions/14_expint/pr68397.cc	(nonexistent)
+++ testsuite/special_functions/14_expint/pr68397.cc	(working copy)
@@ -0,0 +1,47 @@
+// { dg-do run { target c++11 } }
+// { dg-options "-D__STDCPP_WANT_MATH_SPEC_FUNCS__" }
+// Copyright (C) 2017 Free Software Foundation, Inc.
+//
+// This file is part of the GNU ISO C++ Library.  This library is free
+// software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the
+// Free Software Foundation; either version 3, or (at your option)
+// any later version.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License along
+// with this library; see the file COPYING3.  If not see
+// <http://www.gnu.org/licenses/>.
+
+// PR libstdc++/68397 -  std::tr1::expint fails in __expint_En_cont_frac
+// for some long double arguments due to low __max_iter value
+
+#include <cmath>
+#include <testsuite_hooks.h>
+
+int
+test01()
+{
+  // Answers from Wolfram Alpha.
+  long double ans_ok = -0.10001943365331651406888645149537315243646135979573L;
+  long double ans_bomb = -0.10777727809650077516264612749163100483995270163783L;
+
+  auto Ei_ok = std::expint(-1.500001L);
+  auto diff_ok = Ei_ok - ans_ok;
+  VERIFY(std::abs(diff_ok) < 1.0e-15);
+
+  auto Ei_bomb = std::expint(-1.450001L);
+  auto diff_bomb = Ei_bomb - ans_bomb;
+  VERIFY(std::abs(diff_bomb) < 1.0e-15);
+}
+
+int
+main()
+{
+  test01();
+  return 0;
+}

Comments

Jonathan Wakely April 27, 2018, 12:51 p.m. UTC | #1
On 26/04/18 23:17 -0400, Ed Smith-Rowland wrote:
>I'm thinking of going back on my choice not to fix special function 
>bugs on 7.
>
>These build and test clean on gcc-7.
>
>OK?

OK for gcc-7-branch, thanks.
diff mbox series

Patch

Index: include/tr1/ell_integral.tcc
===================================================================
--- include/tr1/ell_integral.tcc	(revision 259666)
+++ include/tr1/ell_integral.tcc	(working copy)
@@ -685,8 +685,8 @@ 
           const _Tp __kk = __k * __k;
 
           return __ellint_rf(_Tp(0), _Tp(1) - __kk, _Tp(1))
-               - __nu
-               * __ellint_rj(_Tp(0), _Tp(1) - __kk, _Tp(1), _Tp(1) + __nu)
+               + __nu
+               * __ellint_rj(_Tp(0), _Tp(1) - __kk, _Tp(1), _Tp(1) - __nu)
                / _Tp(3);
         }
     }
@@ -735,9 +735,9 @@ 
 
           const _Tp __Pi = __s
                          * __ellint_rf(__cc, _Tp(1) - __kk * __ss, _Tp(1))
-                         - __nu * __sss
+                         + __nu * __sss
                          * __ellint_rj(__cc, _Tp(1) - __kk * __ss, _Tp(1),
-                                       _Tp(1) + __nu * __ss) / _Tp(3);
+                                       _Tp(1) - __nu * __ss) / _Tp(3);
 
           if (__n == 0)
             return __Pi;
Index: testsuite/special_functions/06_comp_ellint_3/pr66689.cc
===================================================================
--- testsuite/special_functions/06_comp_ellint_3/pr66689.cc	(nonexistent)
+++ testsuite/special_functions/06_comp_ellint_3/pr66689.cc	(working copy)
@@ -0,0 +1,24 @@ 
+// { dg-do run { target c++11 } }
+// { dg-require-c-std "" }
+// { dg-options "-D__STDCPP_WANT_MATH_SPEC_FUNCS__" }
+// { dg-add-options ieee }
+
+#include <cmath>
+#include <testsuite_hooks.h>
+
+void
+test01()
+{
+  double Pi1 = std::comp_ellint_3(0.75, 0.0);
+  VERIFY(std::abs(Pi1 - 1.91099) < 0.00001);
+
+  double Pi2 = std::comp_ellint_3(0.75, 0.5);
+  VERIFY(std::abs(Pi2 - 2.80011) < 0.00001);
+}
+
+int
+main()
+{
+  test01();
+  return 0;
+}
Index: testsuite/special_functions/13_ellint_3/pr66689.cc
===================================================================
--- testsuite/special_functions/13_ellint_3/pr66689.cc	(nonexistent)
+++ testsuite/special_functions/13_ellint_3/pr66689.cc	(working copy)
@@ -0,0 +1,26 @@ 
+// { dg-do run { target c++11 } }
+// { dg-require-c-std "" }
+// { dg-options "-D__STDCPP_WANT_MATH_SPEC_FUNCS__" }
+// { dg-add-options ieee }
+
+#include <cmath>
+#include <testsuite_hooks.h>
+
+void
+test01()
+{
+  const double pi = 3.141592654;
+
+  double Pi1 = std::ellint_3(0.75, 0.0, pi / 2.0);
+  VERIFY(std::abs(Pi1 - 1.91099) < 0.00001);
+
+  double Pi2 = std::ellint_3(0.75, 0.5, pi / 2.0);
+  VERIFY(std::abs(Pi2 - 2.80011) < 0.00001);
+}
+
+int
+main()
+{
+  test01();
+  return 0;
+}
Index: testsuite/tr1/5_numerical_facilities/special_functions/06_comp_ellint_3/pr66689.cc
===================================================================
--- testsuite/tr1/5_numerical_facilities/special_functions/06_comp_ellint_3/pr66689.cc	(nonexistent)
+++ testsuite/tr1/5_numerical_facilities/special_functions/06_comp_ellint_3/pr66689.cc	(working copy)
@@ -0,0 +1,20 @@ 
+
+#include <tr1/cmath>
+#include <testsuite_hooks.h>
+
+void
+test01()
+{
+  double Pi1 = std::tr1::comp_ellint_3(0.75, 0.0);
+  VERIFY(std::abs(Pi1 - 1.91099) < 0.00001);
+
+  double Pi2 = std::tr1::comp_ellint_3(0.75, 0.5);
+  VERIFY(std::abs(Pi2 - 2.80011) < 0.00001);
+}
+
+int
+main()
+{
+  test01();
+  return 0;
+}
Index: testsuite/tr1/5_numerical_facilities/special_functions/14_ellint_3/pr66689.cc
===================================================================
--- testsuite/tr1/5_numerical_facilities/special_functions/14_ellint_3/pr66689.cc	(nonexistent)
+++ testsuite/tr1/5_numerical_facilities/special_functions/14_ellint_3/pr66689.cc	(working copy)
@@ -0,0 +1,22 @@ 
+
+#include <tr1/cmath>
+#include <testsuite_hooks.h>
+
+void
+test01()
+{
+  const double pi = 3.141592654;
+
+  double Pi1 = std::tr1::ellint_3(0.75, 0.0, pi / 2.0);
+  VERIFY(std::abs(Pi1 - 1.91099) < 0.00001);
+
+  double Pi2 = std::tr1::ellint_3(0.75, 0.5, pi / 2.0);
+  VERIFY(std::abs(Pi2 - 2.80011) < 0.00001);
+}
+
+int
+main()
+{
+  test01();
+  return 0;
+}