Fix #5663: Inaccurate Jacobi sn, cn, and dn evaluation
I did apply the new condition replacing the old |m| < eps. Without that, jacobi_sn(1+100*%i, 1e-5) is way off, as mentioned below. The new elliptic-sn-descending should use keyword args. I also dislike returning different number of values based on some parameter. It might be better to have elliptic-cn-descending and elliptic-dn-descending. Common code could be factored out. I haven't looked at that yet. (Or maybe made into a macro.) All the comments I had previously made are gone, so it's harder...
Refine thresholds for new tests
Fix #5663: Fix inaccurate jacobi_sn and friends
Assuming I implemented the suggested changes correctly, jacobi_sn(20.0,2.5) is still a bit inaccurate: 0.2282498060629115. To get the correct answer you show, I need to transform m to 1/m. Then I get the right answer. In addition, these tests fail now. The accuracy has been reduced. ***************** rtest_elliptic.mac: Problem 236 (line 1083) ***************** Input: closeto(jacobi_am(5.0b-1, 1.5b0) - 4.7071978970469916313225786086966b-1, 7.7038b-34) Result: 3.0814879110195773648895647081359b-33...
Assuming I implemented the suggested changes correctly, jacobi_sn(20.0,2.5) is still a bit inaccurate: 0.2282498060629115. To get the correct answer you show, I need to transform m to 1/m. Then I get the right answer.
That looks good, but I haven't tested it yet. Still has a problem when m is near 1, but I guess there's not much we can do about that.
For the example jacobi_sn(20,2.5), it looks like we should apply the transformation http://functions.wolfram.com/09.36.16.0024.01 first to get the modulus into the range [-1,1]: jacobi_sn(sqrt(m)*z, 1/m) = sqrt(m)*jacobi_sn(z,m) Or 1 jacobi_sn(sqrt(m) z, ─) m (%o38) ─────────────────────── = jacobi_sn(z, m) sqrt(m) Then subst([m=2.5,z=20],%o38); (%o42) 0.2282498061811312 = 0.2282499166210549 The LHS if very close to the expected 0.2282498061811287.