SpecialFunctions.jl Successor

I would personally really like efforts to be directed towards SIMD implementations of special functions.
Optimal SIMD and serial definitions are likely to be different, as the former ought to be branch-free.


A SimdSpecialFunctions.jl would indeed be quite nice. But in many ways, I see it as a higher priority that we at least have a package where serial special functions can live given that SpecialFunctions.jl appears to not be accepting any PRs anymore. There’s a bunch of work just sitting in the SpecialFunctions.jl PR queue just waiting to bit-rot.

Unfortunately, serial special functions aren’t a very sexy topic, they’re more of a utility that users expect to just exist.


Okay. Despite the name, SLEEFPirates is open for PRs that include implementations other than SLEEF ports, so perhaps a new library for SIMD implementations isn’t as needed.

EDIT: I’d also like to point out that @musm deserves all the credit for the the real work behind SLEEFPirates.

@viralbshah just showed me Libm.jl.


Note that Numerical Recipes has a very restrictive software license which means that no code from the book can be used. And basically you shouldn’t even look at the code. (Yes, that is ridiculous. But the code is pretty unreadable anyway.)

+1 for implementing special functions in Julia, but it’s (much) harder than it appears.


This seems to be a good, recent reference:

The Mathematical Function Computation Handbook, by Nelson Beebe


But the actual code for the library he describes is LGPL.

Implementing various special functions under an MIT license is a perfect idea for a GSOC (JSOC) project: it is well-defined, has immediate results, and valuable to the community.

I would recommend just starting from the literature, instead of looking at code with a license that is not compatible with a permissive one like MIT. Most libraries cite the relevant papers anyway.

Some knowledge of numerical analysis is helpful, but testing can be done by comparing with BigFloat calculations.


Love the idea of SIMD optimizations for special functions, but we can’t omit generic fallbacks. A big advantage of implementing these in Julia is automatic differentiation, and SIMD optimizations are a disadvantage for such purposes. Now, a SIMD-optimized dual numbers implementation…


SpecialFunctions is pretty much already covered by AD rules, so I don’t think it’ll do much here.


Good point. It’s worth keeping mind, though, for the day when someone implements pFq(a1,...,ap; b1,...,bp; x) (generalized hypergeometric function).


I looked at Gil, Segura, and Temme (2007) and it looks like another useful reference.

Some general thoughts about a native Julia package:

  1. I don’t think that fragmentation of “special” functions into smaller packages, depending on some natural constraints and details of implementation & testing, is necessarily a bad idea. Specifically, orthogonal families using recurrence relations and serial implementations of special functions like Bessel, Airy, and the error function just use different techniques and approaches.

    I would focus on (mostly) univariate (but possibly parametric) functions that use a clever approximation over \mathbb{R} to achieve a small relative error. In other words, the yucky formulas involving branches and tables of constants that come from continued fractions, Padé approximation, etc.

  2. Julia’s number types are generic, so it’s not clear what this package should target. I would suggest implementing methods for Float64, with the understanding that less and more precise types just fall back to this, with the implied extra cost or lack of possible precision. Interested contributors can extend these.

  3. The package should establish coding standard from the very beginning, including being explicit about errors we target (do we want custom variants, eg less precise functions?), testing against tables, BigFloat comparisons, random testing, etc; and importantly citation of sources and explicitly discouraging even looking at non-MIT-license compatible code. The primary sources should be papers.

  4. Ideally the package would work the following way: whenever someone needs a special function that is not available, they would just locate the relevant paper, code it, and make a PR, which is reviewed and merged in a reasonable timeframe. This way it would just channel contributions that may happen anyway. The benefit for contributors is that they get a code review (presumably they are implementing these functions because they are needed for their own work).

  5. Implementing functions within GSOC and similar should be considered — these are well-defined tasks and may serve as a great introduction to Julia but may require expertise that few students have.

  6. There should be 2–3 core contributors who have some numerical expertise and are willing to review and merge PRs.


One of my tasks is to study the fully generalized hypergeometric functions using geometry. I don’t claim to know every single thing necessary for this, but I am studying the general theory behind it and also pursuing this problem using some tricks. Fully generalized hypergeometric function is one of my long term goals.


I’m not sure why a package can’t be both a wrapper for openspecfun and also have native julia implementations. There are already a bunch of native implementations in there, and it seems simpler/less confusing to keep things together.


As far as my research on this goes, I will be creating my own repositories and not do GSOC or contribute to existing repositories. I work from first principles and do a completely original implementation from scratch, and will not want the hassle of dealing with a committee (as an unpaid volunteer, I want 100% full control of my repo). Also, my default license will probably be AGPL.

1 Like

Mmmm, if I may be so bold, let me say one thing… If some of the people around here were willing to act as a committee in any endeavour of mine and offer their expertise to me for free I would be on cloud nine, to be quite honest :smiley: . I used to think a bit like you in this respect, but then I discovered open source. I’m not interested in being in full control anymore, just want to share.

But back on topic, I agree a SpecialFunctions replacement is a key piece in the ecosystem, and deserves a good deal of attention to make it top notch and remain maintained. I’ve used it more that once in the past and found it very good. I didn’t even realize it was not pure Julia.


I’m not opposed to receiving guidance and free help from the community, I just dont like being required to answer to a committee when I’m a volunteer, as an obligation. Service obligations like that require payment.

I don’t think anyone expects that from you so no need to worry :). And if they do, what leverage do they have heh.


Apparently the wrapper for openspecfun (SpecialFunctions.jl) just wants to remain that.

Whether to separate wrapper and native libraries is a matter of style/preference. Personally I find it cleaner for many reasons, one of which is versioning: for a wrapper package you mostly follow changes in the library that is being wrapped, and aim for feature completeness, while a native library is more open.

Given how lightweight Julia packages are, I usually keep things modular unless there is a reason not to.

1 Like

Starting from first principles requires starting from scratch.

SpecialFunctions already is more than just a wrapper. There are a couple of native implementations in there.

openspecfun itself doesn’t really change anymore, but I suppose it would be fine to make an OpenSpecFun.jl package and have SpecialFunctions depend on that. To me it doesn’t seem worth the effort though.

My point was mostly about the segregation of the functions. It’s nice having a single package where people can go look for these things.