ENH: special: add logiv and logive

[sethtroisi]

Fixes #12607

I'm not sure if this is the right structure, the right place, or the right name.

1. logiv, log_iv, lniv, ivln
   context: scipy.special contains ndtr(x) and log_ndtr(x) | gamma(z), gammaln(x), and loggamma(z) | beta(a, b), betaln(a, b)

2. I'm not sure what precision goal I should aim for;

   1. percent error seem to be very small when v is large but can be as high as 1e-9 when v is just greater than the threshold.

3. I'm not sure if I need to support complex numbers for these functions.

4. I don't yet support array/arraylike input.

   1. Is there a guide or example PR someone can point me at?

[sethtroisi]

@WarrenWeckesser I'd like to hear your opinion if these are worth including.

[ev-br]

A couple of drive-by comments, not a definitive review:

• it'd be nice to deduplicate with _hyp0f1, if possible

• the test_mpmath machinery is very useful for smoking out corner cases. Would be nice to include this

• best ensure that what you add here is available from the cython_special interface.

• these functions should be ufuncs, so cython kernels need to be hooked up into the ufunc-generation machinery

[ev-br]

Provably want to add a See Also section

[ev-br]

Maybe move the dispatch to Cython, so that cython_special interface is consistent.

[sethtroisi]

A couple of drive-by comments, not a definitive review:

@ev-br These are very useful.

• it'd be nice to deduplicate with _hyp0f1, if possible

I keep this in mind but won't look at it till we're closer to the end.

• the test_mpmath machinery is very useful for smoking out corner cases. Would be nice to include this

I wasn't aware of this test; I added a test_logiv method which confirms accurate to 3e-14

• best ensure that what you add here is available from the cython_special interface.
• these functions should be ufuncs, so cython kernels need to be hooked up into the ufunc-generation machinery

Can you point me at some links for doing this? I'm trying to follow 4e86352 but I'm still missing something.

[ev-br]

I keep this in mind but won't look at it till we're closer to the end.

Makes sense!

Re ufuncs: you'll need to add an entry to functions.json, the docstrings go to add_newdocs.py.

I'm not sure about the moving parts now: It used to be necessary to run generate_ufuncs.py manually, not sure if it's still needed (if it is, please make a separate GEN commit).

The most authoritative review should be from Josh or Pauli.

[sethtroisi]

Re ufuncs: you'll need to add an entry to functions.json, the docstrings go to add_newdocs.py.

I'm not sure about the moving parts now: It used to be necessary to run generate_ufuncs.py manually, not sure if it's still needed (if it is, please make a separate GEN commit).

I don't believe this is needed anymore; _generate_pyx.py seems to be run during setup.py

I defined my function in _logiv.pxd as cdef double log_iv_asym(double v, double z) nogil: and added that definition (dd->d) to functions.json. I can import my function as from ._ufuncs import _log_iv_asym but I'm not able to pass an ndarray as parameter. I suspect and I'm missing a ufuncs wrapper function somewhere (_ufuncs.pyi maybe?)

[ev-br]

Seems to work?

In [8]: from scipy.special._ufuncs import _log_iv_asym                          

In [9]: _log_iv_asym([1, 2, 3], [2, 3, 5])                                      
Out[9]: array([0.46121925, 0.80862499, 2.33515425])

However, I'm not sure what is the end goal here: is it really useful to add the asymptotic expansion ufunc to the API (albeit a private part of it). I suspect you'll change it later so that from scipy.special import log_ive works without a python wrapper.

[dschmitz89]

CC @fancidev : these functions would come in very handy for statistical computations, e. g. the Von Mises-Fisher distribution, in case you are interested ;)

[fancidev]

@dschmitz89 The log BesselI function would be very useful. I also wanted it (as well as log BesselK function) when I worked on a compound Poisson process. The ive and kve functions were helpful but they still overflow / underflow for large $x$ if I remember correctly.

However, I guess it’s quite challenging to write a solid logiv. I’ll probably do some research after I finish some of the other easier ones, especially the univariate “all roots” finder which I’m quite interested in.