[Python-ideas] Fwd: Trigonometry in degrees

Matt Arcidy marcidy at gmail.com
Tue Jun 12 02:43:33 EDT 2018

```Sorry for top posting, but these aren't really opinions for the
debate, just information.  I haven't seen them mentioned, and none

Number representation::
IEEE-754 doubles cannot represent pi correctly at any bit-depth (i
mean, obviously, but more seriously).
input pi = 3.14159265358979323846264338327
output pi = 3.14159265358979311599796346854

The cutoff is the value everyone uses, which is what is in math.pi
among other places.
3.141592653589793

Any conversion using this value of pi is already wrong, and the input
to the function is already wrong.  The function doesn't matter.  It's
always underestimating if using pi.  I'm clarifying this because it's
separate from any discussion of functions.  More accurate or different
functions cannot use a better input pi.

Calculator: http://www.binaryconvert.com/convert_double.html
input pi source: https://www.piday.org/million/

Libraries:
Boost implements a cos_pi function, it's algorithm will probably be
useful to look at.

glibc implements cos/sin as a look-up table, which is most likely
where any other implementation will end up, as it is a common
endpoint.  I've seen this on different architectures and libraries.
(sincostab.h if you want to google glibc's table).  Maybe there is a
hilarious index lookup off-by-1 there, I didn't look.

Tables are the basis for the calculation in many libraries  in march
architectures, so I wanted to point that out to any function
designers.  A new implementation may come down to calculating right
index locations.
For a degree based implementation, the same algorithm can be used with
a different table input that is calibrated to degrees.
Likewise,  a pi based implementation, the same but for the pi scale
factor input.
Nothing needs to be done designed other than a new lookup table
calibrated for the "base" of a degree, radian, or pi..  it will have
the same output precision issues calculating index values, but at
least the input will be cleaner.

This is not me saying what should be done, just giving information
that may hopefully be useful

The python math library does not implement algorithms, it exposes the
c functions.  You can see C/C++ has this exact issue performing the
calculation there.  As that is "the spec," the spec is defined with
the error.  The math library is technically correct given it's stated
purpose and result.
Technically correct, the best kind of correct.

On Mon, Jun 11, 2018 at 10:53 PM Tim Peters <tim.peters at gmail.com> wrote:
>
> [Steven D'Aprano]
>>
>> ...
>>
>> The initial proposal is fine: a separate set of trig functions that take
>> their arguments in degrees would have no unexpected surprises (only the
>> expected ones). With a decent implementation i.e. not this one:
>>
>>     # don't do this
>>     def sind(angle):
>
>
> But that's good enough for almost all real purposes.  Indeed, except for argument reduction, it's essentially how scipy's sindg and cosdg _are_ implemented:
>
> https://github.com/scipy/scipy/blob/master/scipy/special/cephes/sindg.c
>
>>
>> and equivalent for cos, we ought to be able to get correctly rounded
>> values for nearly all the "interesting" angles of the circle, without
>> those pesky rounding issues caused by π not being exactly representable
>> as a float.
>>
> If people are overly ;-) worried about tiny rounding errors, just compute things with some extra bits of precision to absorb them.  For example, install `mpmath` and use this:
>
>      def sindg(d):
>         import math, mpmath
>         d = math.fmod(d, 360.0)
>         if abs(d) == 180.0:
>             return 0.0
>         with mpmath.extraprec(12):
>
> Then, e.g,
>
> >>> for x in (0, 30, 90, 150, 180, 210, 270, 330, 360):
> ...     print(x, sindg(x))
> 0 0.0
> 30 0.5
> 90 1.0
> 150 0.5
> 180 0.0
> 210 -0.5
> 270 -1.0
> 330 -0.5
> 360 0.0
>
> Notes:
>
> 1. Python's float "%" is unsuitable for argument reduction; e.g.,
>
> >>> -1e-14 % 360.0
> 360.0
>
> `math.fmod` is suitable, because it's exact:
>
> >>> math.fmod(-1e-14, 360.0)
> -1e-14
>
> 2. Using a dozen extra bits of precision make it very likely you'll get the correctly rounded 53-bit result; it will almost certainly (barring bugs in `mpmath`) always be good to less than 1 ULP.
>
> 3. Except for +-180.  No matter how many bits of float precision (including the number of bits used to approximate pi) are used, converting that to radians can never yield the mathematical `pi`; and sin(pi+x) is approximately equal to -x for tiny |x|; e.g., here with a thousand bits:
>
> >>> mpmath.mp.prec = 1000