Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Torbjörn Granlund]

[moved from gmp-devel]


Richard Biener  writes:

  For
>
  #include 
>
  int main()
  {
_Complex double x;
__real x = 3.09126495087690770626068115234375e+8;
__imag x = -4.045689747817175388336181640625e+8;
volatile _Complex double y = ctan (x);
return 0;
  }

If we get a test case somewhat closer to GMP, then it is likely
somebody will look at it.

I expect the maintainers of library called by gcc (mpc?) would be helped
by seeing the actual call to their library.

-- 
Torbjörn
Please encrypt, key id 0xC8601622
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Richard Biener]

On Mon, 19 Nov 2018, Torbjörn Granlund wrote:

> [moved from gmp-devel]
> 
> 
> Richard Biener  writes:
> 
>   For
> >
>   #include 
> >
>   int main()
>   {
> _Complex double x;
> __real x = 3.09126495087690770626068115234375e+8;
> __imag x = -4.045689747817175388336181640625e+8;
> volatile _Complex double y = ctan (x);
> return 0;
>   }
> 
> If we get a test case somewhat closer to GMP, then it is likely
> somebody will look at it.
> 
> I expect the maintainers of library called by gcc (mpc?) would be helped
> by seeing the actual call to their library.

OK, I can reproduce the issue with mpc 1.1.0, mpfr 3.1.4 and gmp 6.1.0
using

#include 
#include 

int main()
{
  mpc_t m;
  mpc_init2 (m, 53);
  mpfr_set_str (mpc_realref (m), "3.09126495087690770626068115234375e+8", 
10, GMP_RNDN);
  mpfr_set_str (mpc_imagref (m), "-4.045689747817175388336181640625e+8", 
10, GMP_RNDN);
  mpfr_clear_flags ();
  mpc_tan (m, m, 0);
  return 0;
}

breaking at tan.c:189 and printing prec and then err I see

Breakpoint 1, mpc_tan (rop=0x7fffdd50, op=0x7fffdd50, rnd=0)
at /space/rguenther/src/svn/trunk2/mpc/src/tan.c:189
189   prec += mpc_ceil_log2 (prec) + err;
$3 = 53
$4 = 7

Breakpoint 1, mpc_tan (rop=0x7fffdd50, op=0x7fffdd50, rnd=0)
at /space/rguenther/src/svn/trunk2/mpc/src/tan.c:189
189   prec += mpc_ceil_log2 (prec) + err;
$5 = 66
$6 = 66

Breakpoint 1, mpc_tan (rop=0x7fffdd50, op=0x7fffdd50, rnd=0)
at /space/rguenther/src/svn/trunk2/mpc/src/tan.c:189
189   prec += mpc_ceil_log2 (prec) + err;
$7 = 139
$8 = 139

...
Breakpoint 1, mpc_tan (rop=0x7fffdd50, op=0x7fffdd50, rnd=0)
at /space/rguenther/src/svn/trunk2/mpc/src/tan.c:189
189   prec += mpc_ceil_log2 (prec) + err;
$29 = 303084
$30 = 303084
...

So somehow err ends up equal to prec and ever increasing?  Which means
we always end up in

  /* OP is no pure real nor pure imaginary, so in theory the real and
 imaginary parts of its tangent cannot be null. However due to
 rounding errors this might happen. Consider for example
 tan(1+14*I) = 1.26e-10 + 1.00*I. For small precision sin(op) and
 cos(op) differ only by a factor I, thus after mpc_div x = I and
 its real part is zero. */
  if (mpfr_zero_p (mpc_realref (x)) || mpfr_zero_p (mpc_imagref (x)))
{
  err = prec; /* double precision */
  continue;
}

As I questioned in the first followup I wonder if GCC can tell mpc
to not exceed 53bits of precision for the result?  Even with denormals
a precision of 303084 bits isn't possible?

Richard.

-- 
Richard Biener 
SUSE LINUX GmbH, GF: Felix Imendoerffer, Jane Smithard, Graham Norton, HRB 
21284 (AG Nuernberg)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Torbjörn Granlund]

Richard Biener  writes:

  > __real x = 3.09126495087690770626068115234375e+8;
  > __imag x = -4.045689747817175388336181640625e+8;
  > volatile _Complex double y = ctan (x);
  > return 0;
  >   }
  > 
  > If we get a test case somewhat closer to GMP, then it is likely
  > somebody will look at it.
  > 
  > I expect the maintainers of library called by gcc (mpc?) would be helped
  > by seeing the actual call to their library.

  OK, I can reproduce the issue with mpc 1.1.0, mpfr 3.1.4 and gmp 6.1.0
  using

  #include 
  #include 

  int main()
  {
mpc_t m;
mpc_init2 (m, 53);
mpfr_set_str (mpc_realref (m), "3.09126495087690770626068115234375e+8", 
  10, GMP_RNDN);
mpfr_set_str (mpc_imagref (m), "-4.045689747817175388336181640625e+8", 
  10, GMP_RNDN);
mpfr_clear_flags ();
mpc_tan (m, m, 0);
return 0;
  }

  As I questioned in the first followup I wonder if GCC can tell mpc
  to not exceed 53bits of precision for the result?  Even with denormals
  a precision of 303084 bits isn't possible?

Use setrlimit(RLIMIT_CPU, ...) and a signal handler?  :-)

I simplified your test case using 1 as the real part and an integer for
the imaginary part.  With the real part set to 1, the computation
finishes in about 10 seconds, with every doubling of it the runtime
almost triples.  (If my observation generalises to an extended real part
range, the compilation for your test case should terminate in about
10*3^log2(4.045689747817175388336181640625e+8/1)/3600/24/365 years.)


-- 
Torbjörn
Please encrypt, key id 0xC8601622
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Richard Biener]

On Mon, 19 Nov 2018, Torbjörn Granlund wrote:

> Richard Biener  writes:
> 
>   > __real x = 3.09126495087690770626068115234375e+8;
>   > __imag x = -4.045689747817175388336181640625e+8;
>   > volatile _Complex double y = ctan (x);
>   > return 0;
>   >   }
>   > 
>   > If we get a test case somewhat closer to GMP, then it is likely
>   > somebody will look at it.
>   > 
>   > I expect the maintainers of library called by gcc (mpc?) would be helped
>   > by seeing the actual call to their library.
> 
>   OK, I can reproduce the issue with mpc 1.1.0, mpfr 3.1.4 and gmp 6.1.0
>   using
> 
>   #include 
>   #include 
> 
>   int main()
>   {
> mpc_t m;
> mpc_init2 (m, 53);
> mpfr_set_str (mpc_realref (m), "3.09126495087690770626068115234375e+8", 
>   10, GMP_RNDN);
> mpfr_set_str (mpc_imagref (m), "-4.045689747817175388336181640625e+8", 
>   10, GMP_RNDN);
> mpfr_clear_flags ();
> mpc_tan (m, m, 0);
> return 0;
>   }
> 
>   As I questioned in the first followup I wonder if GCC can tell mpc
>   to not exceed 53bits of precision for the result?  Even with denormals
>   a precision of 303084 bits isn't possible?
> 
> Use setrlimit(RLIMIT_CPU, ...) and a signal handler?  :-)

:/  But then with GCC we want deterministic behavior -- ISL
added some "computation budget" and error reporting if that was
taken up.  Maybe gmp/mpfr/mpc could consider adding such a thing
(reporting budget overruns via a signal could be acceptable if
the API churn would be too bad for another way).

> I simplified your test case using 1 as the real part and an integer for
> the imaginary part.  With the real part set to 1, the computation
> finishes in about 10 seconds, with every doubling of it the runtime
> almost triples.  (If my observation generalises to an extended real part
> range, the compilation for your test case should terminate in about
> 10*3^log2(4.045689747817175388336181640625e+8/1)/3600/24/365 years.)

Ick.  I wonder if the actual ctan result is infinite or we just
fail to compute it in an efficient manner...

Richard.

-- 
Richard Biener 
SUSE LINUX GmbH, GF: Felix Imendoerffer, Jane Smithard, Graham Norton, HRB 
21284 (AG Nuernberg)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Torbjörn Granlund]

t...@gmplib.org (Torbjörn Granlund) writes:

  I simplified your test case using 1 as the real part and an integer for
  the imaginary part.  With the real part set to 1, the computation
  finishes in about 10 seconds, with every doubling of it the runtime
  almost triples.  (If my observation generalises to an extended real part
  range, the compilation for your test case should terminate in about
  10*3^log2(4.045689747817175388336181640625e+8/1)/3600/24/365 years.)

I meant to write:

I simplified your test case using 1 as the real part and an integer for
the imaginary part.  With the imaginary part set to 1, the computation
finishes in about 10 seconds, with every doubling of it the runtime almost
triples.  (If my observation generalises to an extended imaginary part
range, the compilation for your test case should terminate in about
10*3^log2(4.045689747817175388336181640625e+8/1)/3600/24/365 years.)

-- 
Torbjörn
Please encrypt, key id 0xC8601622
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Torbjörn Granlund]

Richard Biener  writes:

  :/  But then with GCC we want deterministic behavior -- ISL
  added some "computation budget" and error reporting if that was
  taken up.  Maybe gmp/mpfr/mpc could consider adding such a thing
  (reporting budget overruns via a signal could be acceptable if
  the API churn would be too bad for another way).

I'll leave this up to the mpc maintainers.

  > I simplified your test case using 1 as the real part and an integer for
  > the imaginary part.  With the real part set to 1, the computation
  > finishes in about 10 seconds, with every doubling of it the runtime
  > almost triples.  (If my observation generalises to an extended real part
  > range, the compilation for your test case should terminate in about
  > 10*3^log2(4.045689747817175388336181640625e+8/1)/3600/24/365 years.)

  Ick.  I wonder if the actual ctan result is infinite or we just
  fail to compute it in an efficient manner...

The result should be extremely small, so returning Inf would be
quite inaccurate.  :-)

-- 
Torbjörn
Please encrypt, key id 0xC8601622
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Richard Biener]

On Mon, 19 Nov 2018, Torbjörn Granlund wrote:

> Richard Biener  writes:
> 
>   :/  But then with GCC we want deterministic behavior -- ISL
>   added some "computation budget" and error reporting if that was
>   taken up.  Maybe gmp/mpfr/mpc could consider adding such a thing
>   (reporting budget overruns via a signal could be acceptable if
>   the API churn would be too bad for another way).
> 
> I'll leave this up to the mpc maintainers.
> 
>   > I simplified your test case using 1 as the real part and an integer for
>   > the imaginary part.  With the real part set to 1, the computation
>   > finishes in about 10 seconds, with every doubling of it the runtime
>   > almost triples.  (If my observation generalises to an extended real part
>   > range, the compilation for your test case should terminate in about
>   > 10*3^log2(4.045689747817175388336181640625e+8/1)/3600/24/365 years.)
> 
>   Ick.  I wonder if the actual ctan result is infinite or we just
>   fail to compute it in an efficient manner...
> 
> The result should be extremely small, so returning Inf would be
> quite inaccurate.  :-)

I see.  So I guess one issue might be that mpfr numbers have a
precision (mantissa bits?) but no restriction on the exponent?
That is, for ctan eventually zero _is_ a validly rounded result
and the computation should simply terminate...

Richard.

-- 
Richard Biener 
SUSE LINUX GmbH, GF: Felix Imendoerffer, Jane Smithard, Graham Norton, HRB 
21284 (AG Nuernberg)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Vincent Lefevre]

On 2018-11-19 12:49:12 +0100, Richard Biener wrote:
> I see.  So I guess one issue might be that mpfr numbers have a
> precision (mantissa bits?) but no restriction on the exponent?

MPFR has a bounded (user configurable) exponent range, with underflow
and overflow exceptions for out-of-range results, like in IEEE 754.

-- 
Vincent Lefèvre  - Web: 
100% accessible validated (X)HTML - Blog: 
Work: CR INRIA - computer arithmetic / AriC project (LIP, ENS-Lyon)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Vincent Lefevre]

On 2018-11-19 11:14:49 +0100, Richard Biener wrote:
[...]
> Breakpoint 1, mpc_tan (rop=0x7fffdd50, op=0x7fffdd50, rnd=0)
> at /space/rguenther/src/svn/trunk2/mpc/src/tan.c:189
> 189   prec += mpc_ceil_log2 (prec) + err;
> $29 = 303084
> $30 = 303084
> ...
> 
> So somehow err ends up equal to prec and ever increasing?  Which means
> we always end up in
> 
>   /* OP is no pure real nor pure imaginary, so in theory the real and
>  imaginary parts of its tangent cannot be null. However due to
>  rounding errors this might happen. Consider for example
>  tan(1+14*I) = 1.26e-10 + 1.00*I. For small precision sin(op) and
>  cos(op) differ only by a factor I, thus after mpc_div x = I and
>  its real part is zero. */
>   if (mpfr_zero_p (mpc_realref (x)) || mpfr_zero_p (mpc_imagref (x)))
> {
>   err = prec; /* double precision */
>   continue;
> }

A zero value suggests an underflow somewhere, which implies that...

> As I questioned in the first followup I wonder if GCC can tell mpc
> to not exceed 53bits of precision for the result?  Even with denormals
> a precision of 303084 bits isn't possible?

increasing the precision will generally not solve the underflow
issue, which is due to the limited exponent range. The algorithm
should be modified to do some form of scaling (when possible).

MPFR has internally some very partial support for "unbounded floats"
(UBF), i.e. where the exponent is a mpz_t. They are currently used
for mpfr_fma and mpfr_fmma.

-- 
Vincent Lefèvre  - Web: 
100% accessible validated (X)HTML - Blog: 
Work: CR INRIA - computer arithmetic / AriC project (LIP, ENS-Lyon)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Niels Möller]

Vincent Lefevre  writes:

> MPFR has a bounded (user configurable) exponent range, with underflow
> and overflow exceptions for out-of-range results, like in IEEE 754.

Good! I didn't find that in the manual until I knew the phrase to search
for. It's in
https://www.mpfr.org/mpfr-current/mpfr.html#Exception-Related-Functions,
and relevant functions are mpfr_set_emin, mpfr_set_emax and mpfr_subnormalize.

You might want to add an entry for "Exponent range" to the concept index
in the manual.

I guess gcc should set things up similarly to the example on that page,
"This is an example of how to emulate binary double IEEE 754
arithmetic", with the caveat "This emulates a double IEEE 754 arithmetic
with correct rounding in the subnormal range, which may not be the case
for your hardware".

Regards,
/Niels

-- 
Niels Möller. PGP-encrypted email is preferred. Keyid 368C6677.
Internet email is subject to wholesale government surveillance.
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: [Mpc-discuss] Bug in gmp/mpfr/mpc, never ending ctan computation

[Andreas Enge]

Hello all,

On Mon, Nov 19, 2018 at 12:37:16PM +0100, Torbjörn Granlund wrote:
> I'll leave this up to the mpc maintainers.

indeed I suggest to continue the discussion only on the mpc mailing list;
this is clearly not a bug of gmp.

>   > I simplified your test case using 1 as the real part and an integer for
>   > the imaginary part.  With the real part set to 1, the computation
>   > finishes in about 10 seconds, with every doubling of it the runtime
>   > almost triples.  (If my observation generalises to an extended real part
>   > range, the compilation for your test case should terminate in about
>   > 10*3^log2(4.045689747817175388336181640625e+8/1)/3600/24/365 years.)
>   Ick.  I wonder if the actual ctan result is infinite or we just
>   fail to compute it in an efficient manner...
> The result should be extremely small, so returning Inf would be
> quite inaccurate.  :-)

More precisely, the result is essentially epsilon1+(1+epsilon2)*I;
and the mpc function returns round(epsilon1)+1*I with an additional value
indicating the sign of the error made by dropping epsilon2. For this,
an extremely high precision is needed (essentially, as many bits as
the exponent of epsilon2).

Andreas

___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Vincent Lefevre]

On 2018-11-19 14:17:13 +0100, Niels Möller wrote:
> You might want to add an entry for "Exponent range" to the concept index
> in the manual.

The concept index already has "Exponent", and the corresponding
paragraph is:

   An “exponent” is a component of a regular floating-point number.  Its
C data type is ‘mpfr_exp_t’.  Valid exponents are restricted to a subset
of this type, and the exponent range can be changed globally as
described in *note Exception Related Functions::.  Special values do not
have an exponent.

I think that this is quite explicit.

> I guess gcc should set things up similarly to the example on that page,
> "This is an example of how to emulate binary double IEEE 754
> arithmetic", with the caveat "This emulates a double IEEE 754 arithmetic
> with correct rounding in the subnormal range, which may not be the case
> for your hardware".

AFAIK, GCC already changes the exponent range to the binary64 one
(the lack of subnormals is not an issue, because in case of underflow,
GCC ignores the computations done with MPFR).

But for its internal computations, I assume that mpc extends the
exponent range (if it didn't you would have the same kind of problem,
but much more often). Here it seems that the issue is that the maximum
supported exponent range is not large enough for the chosen algorithm
(which may not be the right one for this kind of input).

-- 
Vincent Lefèvre  - Web: 
100% accessible validated (X)HTML - Blog: 
Work: CR INRIA - computer arithmetic / AriC project (LIP, ENS-Lyon)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Richard Biener]

On Mon, 19 Nov 2018, Niels Möller wrote:

> Vincent Lefevre  writes:
> 
> > MPFR has a bounded (user configurable) exponent range, with underflow
> > and overflow exceptions for out-of-range results, like in IEEE 754.
> 
> Good! I didn't find that in the manual until I knew the phrase to search
> for. It's in
> https://www.mpfr.org/mpfr-current/mpfr.html#Exception-Related-Functions,
> and relevant functions are mpfr_set_emin, mpfr_set_emax and mpfr_subnormalize.
> 
> You might want to add an entry for "Exponent range" to the concept index
> in the manual.

Indeed I didn't find these functions either ...

> I guess gcc should set things up similarly to the example on that page,
> "This is an example of how to emulate binary double IEEE 754
> arithmetic", with the caveat "This emulates a double IEEE 754 arithmetic
> with correct rounding in the subnormal range, which may not be the case
> for your hardware".

OK, GCC doesn't use those, for the testcase if I set emax/emin to the
range from 128bit IEEE then the folding finishes quickly.

It looks like emin/emax is global and not per mpfr_t so for GCC we
have to jump through some additional hoops.  I'll see to adjust
GCC to set the limits based on the maximum range ever needed
(largest FP mode supported by the target) and look for fallout.

Richard.

> Regards,
> /Niels
> 
> 

-- 
Richard Biener 
SUSE LINUX GmbH, GF: Felix Imendoerffer, Jane Smithard, Graham Norton, HRB 
21284 (AG Nuernberg)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Richard Biener]

On Mon, 19 Nov 2018, Vincent Lefevre wrote:

> On 2018-11-19 14:17:13 +0100, Niels Möller wrote:
> > You might want to add an entry for "Exponent range" to the concept index
> > in the manual.
> 
> The concept index already has "Exponent", and the corresponding
> paragraph is:
> 
>An “exponent” is a component of a regular floating-point number.  Its
> C data type is ‘mpfr_exp_t’.  Valid exponents are restricted to a subset
> of this type, and the exponent range can be changed globally as
> described in *note Exception Related Functions::.  Special values do not
> have an exponent.
> 
> I think that this is quite explicit.
> 
> > I guess gcc should set things up similarly to the example on that page,
> > "This is an example of how to emulate binary double IEEE 754
> > arithmetic", with the caveat "This emulates a double IEEE 754 arithmetic
> > with correct rounding in the subnormal range, which may not be the case
> > for your hardware".
> 
> AFAIK, GCC already changes the exponent range to the binary64 one
> (the lack of subnormals is not an issue, because in case of underflow,
> GCC ignores the computations done with MPFR).

GCC does ignore results in case of overflow or underflow but it
doesn't currently limit mpfrs exponent range.  Or rather most
languages don't - I do see mpfr_set_e{min,max} calls being
done in the Fortran frontend only, but only temporarily and
seemingly only to call mpfr_subnormalize.

> But for its internal computations, I assume that mpc extends the
> exponent range (if it didn't you would have the same kind of problem,
> but much more often). Here it seems that the issue is that the maximum
> supported exponent range is not large enough for the chosen algorithm
> (which may not be the right one for this kind of input).

It happens that limiting the range makes the CPOW computation fast.

Richard.

-- 
Richard Biener 
SUSE LINUX GmbH, GF: Felix Imendoerffer, Jane Smithard, Graham Norton, HRB 
21284 (AG Nuernberg)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Dennis Clarke]

On 11/19/18 5:40 AM, Torbjörn Granlund wrote:

Richard Biener  writes:

   > __real x = 3.09126495087690770626068115234375e+8;
   > __imag x = -4.045689747817175388336181640625e+8;
   > volatile _Complex double y = ctan (x);
   > return 0;
   >   }
   >
   > If we get a test case somewhat closer to GMP, then it is likely
   > somebody will look at it.
   >
   > I expect the maintainers of library called by gcc (mpc?) would be helped
   > by seeing the actual call to their library.

   OK, I can reproduce the issue with mpc 1.1.0, mpfr 3.1.4 and gmp 6.1.0
   using

   #include 
   #include 

   int main()
   {
 mpc_t m;
 mpc_init2 (m, 53);
 mpfr_set_str (mpc_realref (m), "3.09126495087690770626068115234375e+8",
   10, GMP_RNDN);
 mpfr_set_str (mpc_imagref (m), "-4.045689747817175388336181640625e+8",
   10, GMP_RNDN);
 mpfr_clear_flags ();
 mpc_tan (m, m, 0);
 return 0;
   }


using that test I see 100% cpu and slowly growing RSS with a lot of time 
being spent inside mpn_fft_fft and then mpn_fft_mul_2exp_modF :


Program received signal SIGINT, Interrupt.
0x0008002dcd04 in __gmpn_copyi () from 
/home/dclarke/local/lib/libgmp.so.10

(gdb) where
#0  0x0008002dcd04 in __gmpn_copyi ()
   from /home/dclarke/local/lib/libgmp.so.10
#1  0x00080029a021 in mpn_fft_mul_2exp_modF (r=0x800e47610, 
a=0x80150bcb0,

d=1600, n=40) at mul_fft.c:267
#2  0x00080029a3f0 in mpn_fft_fft (Ap=0x800e8e130, K=32, 
ll=0x800e00098,

omega=160, n=40, inc=16, tp=0x800e47610) at mul_fft.c:391
#3  0x00080029a351 in mpn_fft_fft (Ap=0x800e8dc30, K=64, 
ll=0x800e000a0,

omega=80, n=40, inc=8, tp=0x800e47610) at mul_fft.c:383
#4  0x00080029a3b2 in mpn_fft_fft (Ap=0x800e8dc10, K=128, 
ll=0x800e000a8,

omega=40, n=40, inc=4, tp=0x800e47610) at mul_fft.c:384
#5  0x00080029a351 in mpn_fft_fft (Ap=0x800e8dc10, K=256, 
ll=0x800e000b0,

omega=20, n=40, inc=2, tp=0x800e47610) at mul_fft.c:383
#6  0x00080029a351 in mpn_fft_fft (Ap=0x800e8dc10, K=512, 
ll=0x800e000b8,

omega=10, n=40, inc=1, tp=0x800e47610) at mul_fft.c:383
#7  0x0008002995ff in mpn_mul_fft_internal (op=0x8014d4150, pl=9728, 
k=9,
Ap=0x800e8dc10, Bp=0x800e8b410, A=0x8014fd610, B=0x801060950, 
nprime=40,

l=19, Mp=5, fft_l=0x800e00070, T=0x800e47610, sqr=1) at mul_fft.c:738
#8  0x000800298ef8 in __gmpn_mul_fft (op=0x8014d4150, pl=9728,
n=0x801336980, nl=9475, m=0x801336980, ml=9475, k=9) at mul_fft.c:904
#9  0x0008002cc6d6 in __gmpn_sqrmod_bnm1 (rp=0x8014a8500, rn=19456,
ap=0x801336980, an=9475, tp=0x8014d4150) at sqrmod_bnm1.c:194
#10 0x00080029c568 in __gmpn_nussbaumer_mul (pp=0x8014a8500,
ap=0x801336980, an=9475, bp=0x801336980, bn=9475) at 
nussbaumer_mul.c:61

--Type  for more, q to quit, c to continue without paging--
#11 0x00080029b80f in __gmpn_sqr (p=0x8014a8500, a=0x801336980, n=9475)
at sqr.c:97
#12 0x000800371726 in mpfr_sqrhigh_n (rp=0x801483500, np=0x801324180,
n=18947) at mulders.c:180
#13 0x00080031c0ab in mpfr_mul (a=0x7fffde50, b=0x7fffde50,
c=0x7fffde50, rnd_mode=MPFR_RNDD) at mul.c:968
#14 0x00080036b2c0 in mpfr_sqr (a=0x7fffde50, b=0x7fffde50,
rnd_mode=MPFR_RNDD) at sqr.c:561
#15 0x000800312bd0 in mpfr_exp_3 (y=0x7fffe350, x=0x7fffe380,
rnd_mode=MPFR_RNDD) at exp3.c:232
#16 0x00080031487f in mpfr_exp (y=0x7fffe350, x=0x7fffe380,
rnd_mode=MPFR_RNDD) at exp.c:176
#17 0x00080033a525 in mpfr_sinh_cosh (sh=0x7fffe550,
ch=0x7fffe530, xt=0x7fffe880, rnd_mode=MPFR_RNDN)
at sinh_cosh.c:108
#18 0x0008003bf37b in mpc_sin_cos (rop_sin=0x7fffe7e0,
rop_cos=0x7fffe7a0, op=0x7fffe860, rnd_sin=17, rnd_cos=17)
at sin_cos.c:390
#19 0x0008003c3756 in mpc_tan (rop=0x7fffe860, op=0x7fffe860,
rnd=0) at tan.c:196
#20 0x002013c0 in main () at test_case.c:19
(gdb) quit
A debugging session is active.

Inferior 1 [process 58829] will be killed.

Quit anyway? (y or n) y
titan$

This is gmp 6.1.2 and mpfr-4.0.1-p13 and mpc-1.1.0 on FreeBSD 12.0RC1
 and I also see similar behavior on Solaris SPARC.



I simplified your test case using 1 as the real part and an integer for
the imaginary part.  With the real part set to 1, the computation
finishes in about 10 seconds, with every doubling of it the runtime
almost triples.  (If my observation generalises to an extended real part
range, the compilation for your test case should terminate in about
10*3^log2(4.045689747817175388336181640625e+8/1)/3600/24/365 years.)



Well the test case provided merely packs up and goes away with no
indication that it will return anytime soon.

I am just going to watch this thread closely.

Dennis

___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Vincent Lefevre]

On 2018-11-19 15:26:00 +0100, Richard Biener wrote:
> On Mon, 19 Nov 2018, Vincent Lefevre wrote:
> > AFAIK, GCC already changes the exponent range to the binary64 one
> > (the lack of subnormals is not an issue, because in case of underflow,
> > GCC ignores the computations done with MPFR).
> 
> GCC does ignore results in case of overflow or underflow but it
> doesn't currently limit mpfrs exponent range.  Or rather most
> languages don't - I do see mpfr_set_e{min,max} calls being
> done in the Fortran frontend only, but only temporarily and
> seemingly only to call mpfr_subnormalize.

Thanks for the clarification. FYI, mpfr_subnormalize was added
for gfortran.

> > But for its internal computations, I assume that mpc extends the
> > exponent range (if it didn't you would have the same kind of problem,
> > but much more often). Here it seems that the issue is that the maximum
> > supported exponent range is not large enough for the chosen algorithm
> > (which may not be the right one for this kind of input).
> 
> It happens that limiting the range makes the CPOW computation fast.

Yes, because for some functions, underflow or overflow can be detected
before the actual computation.

-- 
Vincent Lefèvre  - Web: 
100% accessible validated (X)HTML - Blog: 
Work: CR INRIA - computer arithmetic / AriC project (LIP, ENS-Lyon)
___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs

Re: Bug in gmp/mpfr/mpc, never ending ctan computation

[Dennis Clarke]

On 11/19/18 4:59 AM, Torbjörn Granlund wrote:

[moved from gmp-devel]


Richard Biener  writes:

   For



   #include 



   int main()
   {
 _Complex double x;
 __real x = 3.09126495087690770626068115234375e+8;
 __imag x = -4.045689747817175388336181640625e+8;
 volatile _Complex double y = ctan (x);
 return 0;
   }



Unrelated but interesting none the less, even Wolfram Alpha seems to
have troubles :

https://www.wolframalpha.com/input/?i=tan+(+3091264.95087690770626068115234375+-+4045689.747817175388336181640625i+)

That is with the test data divided by ten. Use the actual test data and
the output simply is blank. Interesting.

Dennis


___
gmp-bugs mailing list
gmp-bugs@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-bugs