HP 17B Solver - another ARCTAN(Y/X) approximation

[Albert Chan]

(05-16-2021 05:04 PM)Albert Chan Wrote:  We could do 1 more, with θ = pi/4:

atan(x) = pi/4 + atan((x-1)/(x+1))

However, this is more messy with unknown sign of x (we may need θ = -pi/4 instead)

Updated my original version (1st post) with θ = ±pi/4 reduction.
This may not be optimized (more like a patch to original).

(tan(pi/8) = √(2)-1) < |z| ≤ (tan(pi/4) = 1):

atan(z)
= s * atan(|z|)                               // where s = sign(z)
= s * (pi/4 + atan((|z|-1)/(|z|+1)))
= s*pi/4 + atan((z-s)/(z*s+1))

Code:

function my_atan2(y, x)
    local z, offset, p = y/x, 0, (y >= 0 and pi or -pi)
    if abs(z)>1 then z, offset = -x/y, p/2 end
    if x < 0 then offset = p - offset end -- Q2, Q3
    if abs(z) > sqrt(2)-1 then
        local s = (z >=0 and 1 or -1)
        z, offset = (z-s)/(z*s+1), offset + s*pi/4
    end
    local s = z*z
    z = z * (((15159/35*s + 4213)*s + 9867)*s + 6435) /
      ((((35*s + 1260)*s + 6930)*s + 12012)*s + 6435)
    return z + offset
end

Both this and previous post version had about same max rel. error ≈ 1e-11
This version reduced atan(z), with |z| ≤ tan(pi/8) = √(2)-1 = 0.41421356 ...

(05-16-2021 05:04 PM)Albert Chan Wrote:  atan(1) = 2*atan(29/70) + atan(-1/8119)

29/70 = 0.41428571 ...

We can use another Pade approximation, to force the split to below tan(pi/8).
(it does not reduce max. rel. error by much, thus not used in code)
Let s = x*x

atan(x) = 2*atan(x*(s+4)*(6*s+8) / horner([1,24,80,64],s))
            + atan(x*s^6 / horner([13,364,2912,9984,16640,13312,4096],s))

atan(1) = 2*atan(70/169) + atan(1/47321)

70/169 = 0.41420118 ...

Trivia: both 29/70 and 70/169 are convergents of √(2)-1