changeset 30781:d041ca628d99 stable

log1p: Correct order of arguments for atan2 (bug #62094). * liboctave/numeric/lo-specfun.cc (log1p): Correct order of arguments for atan2 for calculation of imaginary part of result. * libinterp/corefcn/mappers.cc (Flog1p): Add tests.
author John W. Eaton <jwe@octave.org>
date Mon, 21 Feb 2022 18:50:57 -0500
parents e36380193f2b
children 230c1c39744b a4dbec69e1b5
files libinterp/corefcn/mappers.cc liboctave/numeric/lo-specfun.cc
diffstat 2 files changed, 10 insertions(+), 3 deletions(-) [+]
line wrap: on
line diff
--- a/libinterp/corefcn/mappers.cc	Fri Feb 25 15:37:46 2022 +0100
+++ b/libinterp/corefcn/mappers.cc	Mon Feb 21 18:50:57 2022 -0500
@@ -1763,7 +1763,14 @@
 
 /*
 %!assert (log1p ([0, 2*eps, -2*eps]), [0, 2*eps, -2*eps], 1e-29)
-%!assert (log1p (single ([0, 2*eps, -2*eps])), single ([0, 2*eps, -2*eps]), 1e-29)
+%!assert (log1p (single ([0, 2*eps, -2*eps])), ...
+%!        single ([0, 2*eps, -2*eps]), 1e-29)
+## Compare to result from Wolfram Alpha rounded to 16 significant digits
+%!assert <*62094> (log1p (0.1i), ...
+%!                 0.004975165426584041 + 0.09966865249116203i, eps (0.2))
+%!assert <*62094> (log1p (single (0.1i)), ...
+%!                 single (0.004975165426584041 + 0.09966865249116203i), ...
+%!                 eps (single (0.2)))
 
 %!error log1p ()
 %!error log1p (1, 2)
--- a/liboctave/numeric/lo-specfun.cc	Fri Feb 25 15:37:46 2022 +0100
+++ b/liboctave/numeric/lo-specfun.cc	Mon Feb 21 18:50:57 2022 -0500
@@ -1965,7 +1965,7 @@
         {
           double u = 2*r + r*r + i*i;
           retval = Complex (log1p (u / (1+std::sqrt (u+1))),
-                            atan2 (1 + r, i));
+                            atan2 (i, 1 + r));
         }
       else
         retval = std::log (Complex (1) + x);
@@ -1984,7 +1984,7 @@
         {
           float u = 2*r + r*r + i*i;
           retval = FloatComplex (log1p (u / (1+std::sqrt (u+1))),
-                                 atan2 (1 + r, i));
+                                 atan2 (i, 1 + r));
         }
       else
         retval = std::log (FloatComplex (1) + x);