thomasokken / free42

Free42 : An HP-42S Calculator Simulator
https://thomasokken.com/free42/
GNU General Public License v2.0
280 stars 54 forks source link

Kahan complex asinh #69

Closed achan001 closed 5 months ago

achan001 commented 5 months ago

Free42 already have Kahan's accurate acos/acosh, atan/atanh This issue to to complete Kahan's complex arc trig/hyp list.

asin(z) + acos(z) = pi/2

static Complex Casin(Complex z)
{
  Complex a = csqrt(1+z), b = csqrt(1-z);
  double x2 = atan(Real(z) / (Real(a)*Real(b) - Imag(a)*Imag(b)));
  double y2 = asinh(Real(a)*Imag(b) - Imag(a)*Real(b));
  return CMPLX(copysign(x2, Real(z)), copysign(y2, Imag(z)));
}

asinh(z) = swap(asin(swap(z)))

static int Lasinh(lua_State *L) {
  return pushcomplex(L, cswap(Casin(cswap(Z(1)))));
}

https://www.hpmuseum.org/forum/thread-131-post-151459.html#pid151459 https://opensource.apple.com/source/Libm/Libm-47.1/complex.c.auto.html

achan001 commented 5 months ago

Reuse issue #62, #63 test case for asin

p2> from mpmath import *
p2> mp.pretty = 1
p2> mp.dps = 34
p2> x = '3.141592653589793238462643383279503e-10'
p2> z = mpc(x,x)

p2> t = z + z**3/6   # asin(z) taylor series
p2> re(t)
0.0000000003141592653589793238359289127678504
p2> im(t)
0.0000000003141592653589793238565997638880502

p2> a, b = sqrt(1+z), sqrt(1-z)
p2> x2 = atan(re(z) / (re(a)*re(b) - im(a)*im(b)))
p2> y2 = asinh(re(a)*im(b) - im(a)*re(b))
p2> sign(re(z)) * abs(x2)
0.0000000003141592653589793238359289127678504
p2> sign(im(z)) * abs(y2)
0.0000000003141592653589793238565997638880502

Kahan's asin(z) matched taylor series.

Free42-Decimal asin(z) result (bad digits in bold):

3.141592653589793238359289127678505e-10 + 3.141592653589793238565994941221248e-10 i

thomasokken commented 5 months ago

https://github.com/thomasokken/free42/commit/5aefe1ce67ef629e04e622043f496cc0eeb341a3

BTW, Albert, do you like me to create test builds for you (and if so, for which platform), or do you build from source yourself?

achan001 commented 5 months ago

Wow, this is quick!

Test build for windows please.

thomasokken commented 5 months ago

https://thomasokken.com/free42/download/test/

I prefer to respond to your contributions more quickly, but the previous two got held up for months while I was working off a backlog of features and UI improvements that I had originally planned for Plus42 1.0. My backlog is finally clear now, so I expect to be able to be more responsive from here on out.

achan001 commented 5 months ago

Free42-Decimal tested very well. Again, bad digits bolded.

z = (PI, PI) / 1e10

asin(z) 3.141592653589793238359289127678505e-10 + 3.141592653589793238565997638880502e-10i

sin(asin(z)) --> z 3.141592653589793238462643383279505e-10 + 3.141592653589793238462643383279503e-10i