HP-35’s x^y Why?

[Albert Chan]

(11-03-2021 03:21 PM)Albert Chan Wrote:  Does HP71B has XROOT (or equivalent ?)
I am curious of how to efficiently implement this ...

First try ...

y^(1/n) = y^float(1/n) * y^ε       , where ε = 1/n - float(1/n) = (1 - n*float(1/n)) / n

Correction y^ε = exp(ε*ln(y)) ≈ 1 + ε*ln(y)

Code:

function surd(y,n)      -- nthroot of y
    local z, sgn = 1/n, 1
    if y<0 and n%2==1 then y=-y; sgn=-1 end
    local x = y^z * sgn
    z = fma(n,-z,1) * x
    return z==0 and x or x + z*log(y)/n
end

lua> surd(-8, 3)
-2
lua> surd(-8, 3.1)
-NaN

lua> surd(12345,3)
23.111618749807267
lua> 12345 ^ (1/3)
23.111618749807263