Re: [PATCH, Fortran] Extension: COTAN and degree-valued trig intrinsics with -fdec-math

[Fritz Reese <fritzoreese at gmail dot com>]

Tobias,

Many thanks for the comments. I will adjust the patch according to
your advice shortly.

- Fritz

On Mon, Sep 26, 2016, 11:59 Tobias Burnus
<tobias.burnus@physik.fu-berlin.de> wrote:
>
> Fritz Reese wrote:
> > Attached is a patch extending the GNU Fortran front-end to support
> > some additional math intrinsics, enabled with a new compile flag
> > -fdec-math. The flag adds the COTAN intrinsic (cotangent), as well as
> > degree versions of all trigonometric intrinsics (SIND, TAND, ACOSD,
> > etc...).
> >
> > +     /* There is no builtin mpc_cot, so compute x = 1 / tan (x).  */
> > +     val = &result->value.complex;
> > +     mpc_tan (*val, *val, GFC_MPC_RND_MODE);
> > +     mpc_div (*val, one, *val, GFC_MPC_RND_MODE);
>
> The internet remarks: TAN(x) for x -> pi/2 goes to +inf or -inf - while
> COTAN(pi/2) == 0. For values near pi/2, using cos(x)/sin(x) should be
> numerically better than 1/tan(x). [Cf. mpfr_sin_cos  and BUILT_IN_SINCOS.]
>
> I am not sure how big the prevision difference really is. A quick test
> showed results which didn't look to bad:
>
> implicit none
> real(8), volatile  :: x
> x = 3.14159265358d0/2.0d0
> print '(g0)', cotan(x), cos(x)/sin(x), cotan(x)-cos(x)/sin(x)
> end
>
> yields:
>
> 0.48965888601467483E-011
> 0.48965888601467475E-011
> 0.80779356694631609E-027
>
> Usint N[1/Tan[314159265358/10^11/2],200], Mathematica shows
>   4.8966192313216916397514812338... * 10^12.
> I am not quite sure whether I should expect that already the 5th digit
> differs from the results above.
>
>
> > mpfr_set_d (tmp, 180.0l, GFC_RND_MODE);
>
> With some fonts the L after 180.0 looks like a one; I was initially
> puzzled why one uses not 180 but 180.01. Hence, I'd prefer an "L" instead
> of an "l".
>
> However, the second argument of mpfr_set_d is a "double" and 180.0 is
> already a double. Hence, one can even completely get rid of the "l".
>
>
>
> Regarding the decimal trigonometric functions, the internet suggests to
> use
>    sin (fmod ((x), 360) * M_PI / 180)
> instead of
>    sin (x * M_PI / 180)
> to yield better results for angles which are larger than +/-360 degrees.
>
>
> Cheers,
>
> Tobias