Improve support for Jacobi elliptic functions

[eviatarbach]

Currently, all evaluation (numeric and symbolic) of the Jacobi elliptic functions is done through Maxima. No derivatives or arbitrary-precision numeric evaluation are defined. Worse still, the Maxima numerical evaluation is wrong for some of the inverse Jacobi functions; see https://sourceforge.net/p/maxima/bugs/2615/.

CC: @burcin @kcrisman @benjaminfjones @fredrik-johansson

Component: symbolics

Author: Eviatar Bach

Branch/Commit: 6df0d14

Reviewer: Ralf Stephan, Travis Scrimshaw

Issue created by migration from https://trac.sagemath.org/ticket/14996

[eviatarbach]

comment:1

Attachment: trac14996.patch.gz

The patch rewrites the Jacobi functions as standard symbolic functions. I wrote a new numerical routine for the inverse Jacobi functions, since it didn't seem to be available elsewhere. I've also added the Jacobi amplitude function.

[eviatarbach]

comment:4

I hope you don't mind me CC'ing you.

In this patch I add the functions inverse_jacobi_f and jacobi_amplitude_f which use mpmath. Do you think these could be added to a future version? I know I'm not smart about increasing the context precision for successive computations but it seems to work.

[fredrik-johansson]

comment:5

Sure, they could be added to a future version. I never even tried to implement them because it seemed like too much work to figure out all the case distinctions and writing test code for branch cuts. You seem to restrict the code to real arguments, which makes things easier and probably makes a good enough start.

[eviatarbach]

comment:6

Ah, good. Should I submit a pull request on GitHub or something?

The reason I used real numbers for the inverse Jacobi is because it was a pain to get continuity otherwise. For the amplitude function, it's defined in the DLMF for real numbers (http://dlmf.nist.gov/22.16).

[fredrik-johansson]

comment:7

Sure, you're very welcome to supply a patch if you have time.

Better put this code in Sage now though, and then perhaps replace it with calls to mpmath in the future (can't guarantee when the next release of mpmath will happen, regrettably).

[eviatarbach]

comment:8

Thanks, I will do so when I have time.

New patch fixes two bugs.

[eviatarbach]

Attachment: trac14996_2.patch.gz

[eviatarbach]

comment:9

New patch rebases to Sage 5.11, gets coverage to 100%, adds tons of new tests, and makes all Python ints in simplifications and derivatives Sage integers.

Please don't be worried about the length of the patch, much of it derives from identities and tests. Since there are 12 Jacobi elliptic functions and all their inverses, this adds up to a lot.

Patchbot apply trac14996_3.patch

[sagetrac-mjd]

comment:10

Attachment: trac14996_3.patch.gz

I did some tests, and the doctest for plot/line.py failed (it tries to plot a Jacobi elliptic function at line 413 and this goes awry in trying to convert a complex number which is very close into a real number to a real number).

A small issue is that perhaps jacobi should complain if one passes to it as the first argument a string which is not of the form 'pq' with p, q one of s, c, n, d. Moreover, if the argument is 'pp' then it could perhaps use the convention used in the documentation and return 1?

I'm new to sage, apologies for the inanity of my comments.

[eviatarbach]

comment:11

No, they're not inane at all! Thanks for the comments, they bring up important issues to fix.

You're right about the doctest and the need to raise an error; the Jacobi class raises an error but I forgot to add it for the jacobi function. Fix coming up shortly.

As for the pp=1 convention, it could be done but I don't see how it would be useful. In all sources I've seen they're called the "12 Jacobi functions". Also, it might give users the impression that Sage is actually calculating the ratios of the functions and simplifying. To remain consistent with the other functions we would also have to add the aliases jacobi_nn, jacobi_ss, etc., which amounts to just adding four constants that equal 1 to the global namespace.