Bessel functions again

[Waldek Hebisch]

I again looked at our Bessel functions. This time I concentrated
on besselI, as it is supposed to be the easiest one. Looking at
old code, besselI(v, z) signals error when either absolute value
of v is bigger than 15, or absolute value of z is bigger than 15.
So, it is quite limited. Within that range testing using v with
nonnegative real part shows a lot of point where relative error
is bigger than 1.0e-11. When real part of v is negative situation
is worse, there are many point with absolute error bigger than
1.0e-7.

I have now a routine that is better, at least for v with positive
real part. Relative error for v and z with absolute values less
than 35 seem to be of order 1.0e-13. In that range I use so
called Miller method. Apparently old code uses that method too
(that is classic method which is recommended in many texts),
but new code has some twists to improve accuracy. For large
arguments I have routine using asymptic expansion, which
in many cases gives similar error.

Now some things that I realized in the process. First, asymptic
expansion was given by Olver in [3] and reformulated by Luke
in [1]. But area of validity of this expansion for besselI is
wrong in both papers. I have modified the formula so that it
is valid in the range we need. Second, Luke [2] claims that
Miller method is convergent and well conditioned numerically.
Convergence claim AFAICS is true, but in numeric context somewhat
misleading, as convergence can be slow. Concerning numeric
conditioning, I had big issue with rounding error: Miller method
using extended precision produced reasonable results, but at
the beginning I was getting significant error when using
double precision arithmetic. It seems that I overcame it by
modifying "normalizing condition": in first logical step Miller
method produces vector of multiples of besselI(v, z), besselI(v + 1, z),
etc up to besselI(v + n, z) for n large enough. Then it uses
known equality for Bessel functions to determine scalar factor.
My twist is to perform the second step as if we wanted to
compute besselI(v + m, z), were m is rougly choosen to make
magintude of v + m and z comparable.

Before my twist, trying to increase n in Miller method led to
decreasing accuracy, and large n is needed for large z. Now
it seems that I can tweak parameters to get reasonable accuracy
from Miller method for large v and z. For large v and z there
is some hard to avoid loss of accuracy. Namely, Bessel
functions behave similarly to exponentials, which means that
for large arguments they tend to expand small errors in
initial arguments or in intermediate calculations. That is
why it is hard to get relative error below 1.0e-13.

All the above was for order v with real part bigger or equal
to 0. When real part of v is negative, then Bessel function
behave differently. Miller method is no longer valid. In
fact, computing besselI for negative v is related to computing
besselK. Namely, there is relation
\[
I_v(z) = \frac{i}{\pi}\left(K_v(-z) - \exp(-\pi vi)K_v(z)\right),
\]
and
\[
I_{-v}(z) = I_v(z) + \frac{2}{\pi}\sin(v\pi)K_v(z).
\]

Current code tries to use the last relation to compute besselK,
but this is problematic for three reasons:
- for z with positive real part $|K_v|$ tends to be much smaller
than $|I_v|$, so we are subtracting two large nubers to obtain
a small one,
- when $v$ is close to an integer $\sin(v\pi)$ is small and
we are dividing be a small number,
- usual methods of computing $I_v$ tend to give large errors when
real part of $v$ is negative (and our curret code has large
errors in such case).

So the next step (harder than besselI case) is to find reasonable
method to compute besselK. The other functions (besselJ, besselY,
Hankel functions) can be reduced to besselK and besselI.

--
Waldek Hebisch

[Waldek Hebisch]

On Tue, Jul 07, 2026 at 06:33:51PM +0200, Waldek Hebisch wrote:
>
> Now some things that I realized in the process. First, asymptic
> expansion was given by Olver in [3] and reformulated by Luke
> in [1]. But area of validity of this expansion for besselI is
> wrong in both papers. I have modified the formula so that it
> is valid in the range we need. Second, Luke [2] claims that
> Miller method is convergent and well conditioned numerically.
> Convergence claim AFAICS is true, but in numeric context somewhat
> misleading, as convergence can be slow. Concerning numeric
> conditioning, I had big issue with rounding error: Miller method
> using extended precision produced reasonable results, but at
> the beginning I was getting significant error when using
> double precision arithmetic.

I forgot to include references:

[1] Y. L. Luke,
Some remarks on uniform asymptotic expansions for Bessel functions,
Comp. & Maths. with Appls. Vol. I, pp. 285--290.
[2] Y. L. Luke,
On generating Bessel functions by use of the backward recurrence formula,
mauscript present online under name AD0743296.pdf
[3] F. W. J. Olver,
The Asymptotic Expansion of Bessel Functions of Large Order,
Philosophical Transactions of the Royal Society of London. Series A,
Mathematical and Physical Sciences, Dec. 28, 1954, Vol. 247,
No. 930 (Dec. 28, 1954), pp. 328-368

--
Waldek Hebisch

[Waldek Hebisch]

I tried now besselJ. It looks worse than I thougt. Coarse
search finds like:

(1) -> besselJ(15.3 + 14.0*%i, 15.3 - 14.0*%i)

(1) 0.0
Type: Complex(DoubleFloat)

Actual value is rather large, that is -12_2930235006.3101744289
-8_7694487856.2102879948*%i, so both absolte error and relative
error are quite large. And this is not isolated thing, notrivial
fraction of test points hits this behaviour.

(2) -> besselJ(1.2999999999999998, 8.3 - 14.0*%i)

(2) 110294.99367399447 - 30403.61231949334 %i
Type: Complex(DoubleFloat)
Here relative error is 1.3787804329893885e-6. And there are several
points like this.

(4) -> besselJ(13.3 + 18.0*%i, 19.3 - 6.0*%i)

(4) - 4.908525004574761e9 - 1.2968865853458834e9 %i
Type: Complex(DoubleFloat)

Here relative error is 0.13231270400584205, so really large.
This one is rather special, seem to be on border of set
exhibiting first behaviour.

For comparison new Miller style routine on the same set of points
has at worst relative error 4.553186381139345e-12.

--
Waldek Hebisch

[Qian Yun]

There are numerous open source numeric scientific programs,
is it a good idea to choose a mature one and model after it?

I mean, there are many many tricks and pitfalls in the numeric
computing world, we don't need to reinvent the wheel here.

- Qian

[Ralf Hemmecke]

On 7/8/26 03:24, Qian Yun wrote:
> There are numerous open source numeric scientific programs,
> is it a good idea to choose a mature one and model after it?

I would even go as far as using such libraries for FriCAS.
Yes, it is new dependency, but why not using open-source effort of
others and continue with other things. We already rely on GMP (if
available), why not doing the same with other libraries, FLINT, for example?

Ralf

[Kurt Pagani]

On 08/07/2026 08:18, 'Ralf Hemmecke' via FriCAS - computer algebra

system wrote:
> On 7/8/26 03:24, Qian Yun wrote:
>> There are numerous open source numeric scientific programs,
>> is it a good idea to choose a mature one and model after it?

Symbolic relations (simpl) of special functions are more important to me
than numeric evaluations, i.e. nice to have but not essential. The
hypergeometricF in EXPR INT, for instance, can be differentiated
correctly (mostly), but has no simplifications like
hypergeometricF([0],[0],z) -> exp(z).

This is is bug by the way.
D(hypergeometricF([0],[0],x),x)

>> Error detected within library code:
catdef: division by zero


D(hypergeometricF([a],[b],z),z) -> ?

>
> I would even go as far as using such libraries for FriCAS.
> Yes, it is new dependency, but why not using open-source effort of
> others and continue with other things. We already rely on GMP (if
> available), why not doing the same with other libraries, FLINT, for
> example?

If at all, then why not use GSL? It's easy to access via CL.
https://www.gnu.org/software/gsl/doc/html/specfunc.html

@Waldek: What do you think about the approach in

Fast and Accurate Bessel Function Computation
John Harrison
Intel Corporation, JF1-13
https://www.cl.cam.ac.uk/~jrh13/papers/bessel.pdf
>
> Ralf
>

[Waldek Hebisch]

On Wed, Jul 08, 2026 at 08:18:40AM +0200, 'Ralf Hemmecke' via FriCAS - computer algebra system wrote:
> On 7/8/26 03:24, Qian Yun wrote:
> > There are numerous open source numeric scientific programs,
> > is it a good idea to choose a mature one and model after it?
>
> I would even go as far as using such libraries for FriCAS.

You guys forgot what I wrote previous time: there are no full
range machine precision open source libraries for Bessel functions.
There are some closed source ones, but it is not clear how they
do this and how good is the result. And asymptotics of Bessel
functions are part of fundamental knowledge needed to write
good library. As I wrote, one piece of this is wrong in
published litarature.

The closest thing to usable library is Mpmath and followers. AFAIK
it uses arbitrary precision arthetic in almost brute force way:
its main method is using power series. In difficult regions
Mpmath uses chain of expansions, which is dramatically better
than using single power series, but nevertheless quite
inefficient in case of Bessel functions.

I consider using arbitrary precision arithmetic to compute double
precision besselK. I already have arbitrary precision power
series, which ATM I use as oracle producing good values and
allowing me to compute error of double precision routines.
Using arbitrary precision arithmetic to compute double
precision result is not entirely satisfactory, but means that
implementation is simpler and easier, so is reasonable. But
for besselI it seems that version computing in double precision
will be OK.

> Yes, it is new dependency, but why not using open-source effort of others
> and continue with other things. We already rely on GMP (if available), why
> not doing the same with other libraries, FLINT, for example?

Already GMP shows problems with libraries. One problem is that
what library authors want to implement/export differs from what
we need. We use 5 functions in GMP, GMP export I think hundreds.
Trouble is, almost all exported GMP functions are of no use to
us. They are "higher level" functions which depend in C memory
management. That make them hard to use for us. And in important
special case of "small bignums" higher level GMP functions are
slow (much slower than low level ones). There are probably about
20 of GMP functions that would be useful to us, but which GMP
does not expose, so w can not use them.

Another problem with libraries is dependency management. GMP is
widely available now, so with GMP this is not a big problem.
But other libraries are much more problematic. And simply number
of needed libraries may be a problem: it is very discouraging
if before building what you want you need to buid first about
50 other things (that count correspond to tex4ht several years
ago).

--
Waldek Hebisch

[Qian Yun]

On 7/8/26 7:44 PM, Waldek Hebisch wrote:
> On Wed, Jul 08, 2026 at 08:18:40AM +0200, 'Ralf Hemmecke' via FriCAS - computer algebra system wrote:
>> On 7/8/26 03:24, Qian Yun wrote:
>>> There are numerous open source numeric scientific programs,
>>> is it a good idea to choose a mature one and model after it?
>>
>> I would even go as far as using such libraries for FriCAS.
>
> You guys forgot what I wrote previous time: there are no full
> range machine precision open source libraries for Bessel functions.
> There are some closed source ones, but it is not clear how they
> do this and how good is the result. And asymptotics of Bessel
> functions are part of fundamental knowledge needed to write
> good library. As I wrote, one piece of this is wrong in
> published litarature.
>

How about https://github.com/flintlib/arb , which is now part of
https://github.com/flintlib/flint/ ?

It is said to support arbitrary-precision with complex arguments.

https://flintlib.org/doc/acb_hypgeom.html#bessel-functions

For numeric computations, if we do not rely on external routines,
then there's too many places that could go wrong -- as I illustrated
in the other thread about complex sqrt.

FriCAS is a symbolic program at heart. And I do believe the open
source world of numeric computation has grown a lot compared
to 10, or 20 years ago.

- Qian

[Waldek Hebisch]

On Wed, Jul 08, 2026 at 01:01:41PM +0200, Kurt Pagani wrote:
>
>
> On 08/07/2026 08:18, 'Ralf Hemmecke' via FriCAS - computer algebra system
> wrote:
> > On 7/8/26 03:24, Qian Yun wrote:
> > > There are numerous open source numeric scientific programs,
> > > is it a good idea to choose a mature one and model after it?
> Symbolic relations (simpl) of special functions are more important to me
> than numeric evaluations, i.e. nice to have but not essential. The
> hypergeometricF in EXPR INT, for instance, can be differentiated correctly
> (mostly), but has no simplifications like hypergeometricF([0],[0],z) ->
> exp(z).

That is somewhat controversial thing. Namely, advanced algorithms
frequently depend on exact form of representation. Unexpected change
to representation can lead to bugs. It is possible that having
by default simplifications like

hypergeometricF([..., a, ...],[..., a, ...], z) ->
hypergeometricF([..., ...],[..., ...], z)

(that is deleting common parameter), would cause no troube.
Similarly

hypergeometricF([], [],z) -> exp(z)

looks relatively safe (and combination of the 2 above would give
what you want). Still, ATM I decided to severly limit default
simplifications of hypergeometricF functions.

There is different thing with transformations done on demand.
We should have various transfomations. Simply there is a lot
of them and there is question how to best organize them.
ATM I am thinking of extending 'normalize' so that it can
detect dependencies between hypergeometricF and meijerG.
And separte "denormalization" routine supposed to detect
cases when hypergeometricF and meijerG can be expressed in
a simpler form.

> This is is bug by the way.
> D(hypergeometricF([0],[0],x),x)
>
> >> Error detected within library code:
> catdef: division by zero

I see.

(24) -> D(hypergeometricF([],[],x), x)

(24) hypergeometricF([],[],x)
Type: Expression(Integer)

works as expected, so probably detecting that 0 cancels is
enough.

> D(hypergeometricF([a],[b],z),z) -> ?

(25) -> D(hypergeometricF([a],[b],z),z)

a hypergeometricF([a + 1],[b + 1],z)
(25) ------------------------------------
b
Type: Expression(Integer)

That is fine. Trouble is that 0 as lower parameter has somewhat
dubious status.

> > I would even go as far as using such libraries for FriCAS.
> > Yes, it is new dependency, but why not using open-source effort of
> > others and continue with other things. We already rely on GMP (if
> > available), why not doing the same with other libraries, FLINT, for
> > example?
>
> If at all, then why not use GSL? It's easy to access via CL.
> https://www.gnu.org/software/gsl/doc/html/specfunc.html
>
> @Waldek: What do you think about the approach in
>
> Fast and Accurate Bessel Function Computation
> John Harrison
> Intel Corporation, JF1-13
> https://www.cl.cam.ac.uk/~jrh13/papers/bessel.pdf

That paper nicely illustrates troubles. First, people do a lot
of work on various special cases, he really considers only real
arguments and oders 0 and 1, that allows using precomputed special
approximations. This does not scale to arbitrary complex
orders and arguments. Second, he admits using Itanium 128 bit
hardware floats. Having extra precision in intermediate
calculations helps a lot. But pupular machines do not have
128 bit hardware floats and software floating point is
_much_ more expensive than hardware one, so once you use
extended precision "fast" part is no longer valid. CMUCL
has thing called double-double, which should be good enough
for our purposes. But IIUC no other Lisp supports it.

--
Waldek Hebisch

[Kurt Pagani]

On 08/07/2026 14:40, Waldek Hebisch wrote:

>> @Waldek: What do you think about the approach in
>>
>> Fast and Accurate Bessel Function Computation
>> John Harrison
>> Intel Corporation, JF1-13
>> https://www.cl.cam.ac.uk/~jrh13/papers/bessel.pdf
>
> That paper nicely illustrates troubles. First, people do a lot
> of work on various special cases, he really considers only real
> arguments and oders 0 and 1, that allows using precomputed special
> approximations. This does not scale to arbitrary complex
> orders and arguments. Second, he admits using Itanium 128 bit
> hardware floats. Having extra precision in intermediate
> calculations helps a lot. But pupular machines do not have
> 128 bit hardware floats and software floating point is
> _much_ more expensive than hardware one, so once you use
> extended precision "fast" part is no longer valid. CMUCL
> has thing called double-double, which should be good enough
> for our purposes. But IIUC no other Lisp supports it.
>

Thanks, I see, actually I also was skeptic.


There is a quite recent paper on the topic

https://arxiv.org/html/2606.14937v1

however, I haven't grasped it yet, though the algorithms seem not to be
very complicated on first sight.

--> 5 Conclusions

This work completes a unified effort aimed at the efficient
multiprecision computation of Bessel functions. In this final part, we
have presented robust, highly accurate, and computationally efficient
algorithms, together with Fortran implementations, for evaluating the
regular Bessel functions Jν​(z) and Yν​(z) for complex arguments and for
both positive and negative real orders.

[Waldek Hebisch]

On Wed, Jul 08, 2026 at 08:23:15PM +0800, Qian Yun wrote:
> On 7/8/26 7:44 PM, Waldek Hebisch wrote:
> > On Wed, Jul 08, 2026 at 08:18:40AM +0200, 'Ralf Hemmecke' via FriCAS - computer algebra system wrote:
> >> On 7/8/26 03:24, Qian Yun wrote:
> >>> There are numerous open source numeric scientific programs,
> >>> is it a good idea to choose a mature one and model after it?
> >>
> >> I would even go as far as using such libraries for FriCAS.
> >
> > You guys forgot what I wrote previous time: there are no full
> > range machine precision open source libraries for Bessel functions.
> > There are some closed source ones, but it is not clear how they
> > do this and how good is the result. And asymptotics of Bessel
> > functions are part of fundamental knowledge needed to write
> > good library. As I wrote, one piece of this is wrong in
> > published litarature.
> >
>
> How about https://github.com/flintlib/arb , which is now part of
> https://github.com/flintlib/flint/ ?
>
> It is said to support arbitrary-precision with complex arguments.
>
> https://flintlib.org/doc/acb_hypgeom.html#bessel-functions

Arb is a follower to Mpmath. As I wrote, using arbitrary precision
to get double precision result is brute force approach and hence
not entirely satisfactory. But for besselK it may be the way
to go.

> For numeric computations, if we do not rely on external routines,
> then there's too many places that could go wrong -- as I illustrated
> in the other thread about complex sqrt.

External routines or not, in numeric computations various things
may go wrong. Still, people frequently want numeric way
because it is faster.

> FriCAS is a symbolic program at heart.

Numerics have some role. There are cases which are easy symbolically,
but hard numerically. And there are cases hard to do in symbolic
way and easy in numeric way. My understending is that Mathematica
in many places uses numeric computations to get final symbolic
result. I do not want to go as far as they, but sometimes numerics
may be the best way to get symbolic result. As a potential
example consider symbolic integration. One of bottlenecks
in Risch-Trager-Bronstein approach is computing order of
torsion divisors. Using known symbolic approaches this has
exponential complexity. It looks like numeric approach could
reduce compexity to polynomial.

> And I do believe the open
> source world of numeric computation has grown a lot compared
> to 10, or 20 years ago.

Sure, but they concentrate on different things.

--
Waldek Hebisch

[Waldek Hebisch]

Thanks for the reference. Abstrat seem to indicate that this is
more or less rework of Amos who also gave algorithms for computing
Bessel function of real orders. Amos considered only positive orders,
but once you have besselK negative orders are an easy extention.

AFAICS Amos algoritm for besselK is justified only for small imaginary
part of order, and arbitrary complex order probably requires a new
algorithm.

--
Waldek Hebisch

[Dima Pasechnik]

On Wed, Jul 8, 2026 at 6:44 AM Waldek Hebisch <de...@fricas.org> wrote:
>
> On Wed, Jul 08, 2026 at 08:18:40AM +0200, 'Ralf Hemmecke' via FriCAS - computer algebra system wrote:
> > On 7/8/26 03:24, Qian Yun wrote:
> > > There are numerous open source numeric scientific programs,
> > > is it a good idea to choose a mature one and model after it?
> >
> > I would even go as far as using such libraries for FriCAS.
>
> You guys forgot what I wrote previous time: there are no full
> range machine precision open source libraries for Bessel functions.
> There are some closed source ones, but it is not clear how they
> do this and how good is the result. And asymptotics of Bessel
> functions are part of fundamental knowledge needed to write
> good library. As I wrote, one piece of this is wrong in
> published litarature.
>
> The closest thing to usable library is Mpmath and followers. AFAIK
> it uses arbitrary precision arthetic in almost brute force way:
> its main method is using power series. In difficult regions
> Mpmath uses chain of expansions, which is dramatically better
> than using single power series, but nevertheless quite
> inefficient in case of Bessel functions.

Mpmath is mostly abandonware (and it's Python, so I don't really see
how FriCAS would be using it in a sane way),
the primary author continues this line of work in Flint/Arb.

>
> I consider using arbitrary precision arithmetic to compute double
> precision besselK. I already have arbitrary precision power
> series, which ATM I use as oracle producing good values and
> allowing me to compute error of double precision routines.
> Using arbitrary precision arithmetic to compute double
> precision result is not entirely satisfactory, but means that
> implementation is simpler and easier, so is reasonable. But
> for besselI it seems that version computing in double precision
> will be OK.
>
> > Yes, it is new dependency, but why not using open-source effort of others
> > and continue with other things. We already rely on GMP (if available), why
> > not doing the same with other libraries, FLINT, for example?
>
> Already GMP shows problems with libraries. One problem is that
> what library authors want to implement/export differs from what
> we need. We use 5 functions in GMP, GMP export I think hundreds.
> Trouble is, almost all exported GMP functions are of no use to
> us. They are "higher level" functions which depend in C memory
> management. That make them hard to use for us.

I don't understand - Lisp systems such as ECL have no problems using
GMP as the main
(and in case of ECL, the only) arbitrary precision arithmetic backend.
Can't you just use them, and not mess around
with low-level C/assembly code?

> And in important
> special case of "small bignums" higher level GMP functions are
> slow (much slower than low level ones). There are probably about
> 20 of GMP functions that would be useful to us, but which GMP
> does not expose, so w can not use them.
>
> Another problem with libraries is dependency management. GMP is
> widely available now, so with GMP this is not a big problem.
> But other libraries are much more problematic.

Flint/Arb is widely available.

> And simply number
> of needed libraries may be a problem: it is very discouraging
> if before building what you want you need to buid first about
> 50 other things (that count correspond to tex4ht several years
> ago).
>
> --
> Waldek Hebisch
>

> --
> You received this message because you are subscribed to the Google Groups "FriCAS - computer algebra system" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to fricas-devel...@googlegroups.com.
> To view this discussion visit https://groups.google.com/d/msgid/fricas-devel/ak44Kec7F2Lhjzqi%40fricas.org.

[Waldek Hebisch]

On Wed, Jul 08, 2026 at 09:47:43AM -0500, Dima Pasechnik wrote:
> On Wed, Jul 8, 2026 at 6:44 AM Waldek Hebisch <de...@fricas.org> wrote:
> >
> > > Yes, it is new dependency, but why not using open-source effort of others
> > > and continue with other things. We already rely on GMP (if available), why
> > > not doing the same with other libraries, FLINT, for example?
> >
> > Already GMP shows problems with libraries. One problem is that
> > what library authors want to implement/export differs from what
> > we need. We use 5 functions in GMP, GMP export I think hundreds.
> > Trouble is, almost all exported GMP functions are of no use to
> > us. They are "higher level" functions which depend in C memory
> > management. That make them hard to use for us.
>
> I don't understand - Lisp systems such as ECL have no problems using
> GMP as the main
> (and in case of ECL, the only) arbitrary precision arithmetic backend.
> Can't you just use them, and not mess around
> with low-level C/assembly code?

In a sense we "just" use GMP. There is low level (partly assembly)
code dealing with interfacing. The point was that we use tiny part
of exported interface. I wrote 5 functions, but actually we use 4:
multiplication, division, gcd and integer square root. C based
Lisp-s can addtionally use addition, subtraction and bit operations
(in our case there would be no/negligible gain so we do not bother).
But, AFAIK GCL for long time it had its own gcd, so in a sense made
less use of GMP than FriCAS. I think that also ECL only uses few
core routines and ignores the rest.

I do not know how much code is needed for "higher level" stuff in
GMP, but dropping it would probably give nontrivial space saving.
OTOH if available we could use Schoenage-Strassen FFT and some
releated routines like coordinatwise operations on vectors with
fixed size components. But AFAIK such routines are present in GMP,
but not part of official interface.

--
Waldek Hebisch

[Waldek Hebisch]

On Wed, Jul 08, 2026 at 01:01:41PM +0200, Kurt Pagani wrote:
>

> The
> hypergeometricF in EXPR INT, for instance, can be differentiated correctly
> (mostly), but has no simplifications like hypergeometricF([0],[0],z) ->
> exp(z).
>
> This is is bug by the way.
> D(hypergeometricF([0],[0],x),x)
>
> >> Error detected within library code:
> catdef: division by zero

Usually zeros (or negative integers) in lower list are considerd
to be invalid. So this is really user error.

We could detect nonpositive integers in lower list when user tries
to create such function and signal error. We could also remove
common element of upper and lower list, but doing this in general
requires some effort and may be undesirable, so it is not now.
Doing this for nopositive integers in the lower list is easier, but
it is not clear to me if this is really useful.

--
Waldek Hebisch

[Kurt Pagani]

On 09/07/2026 00:04, Waldek Hebisch wrote:
> On Wed, Jul 08, 2026 at 01:01:41PM +0200, Kurt Pagani wrote:
>>
>> The
>> hypergeometricF in EXPR INT, for instance, can be differentiated correctly
>> (mostly), but has no simplifications like hypergeometricF([0],[0],z) ->
>> exp(z).
>>
>> This is is bug by the way.
>> D(hypergeometricF([0],[0],x),x)
>>
>> >> Error detected within library code:
>> catdef: division by zero
>
> Usually zeros (or negative integers) in lower list are considerd
> to be invalid. So this is really user error.

You're absolutely right. In pFq we have p,q length of the lists and not
0,0 *in* the lists :( I must have been out of my senses.

>
> We could detect nonpositive integers in lower list when user tries
> to create such function and signal error. We could also remove
> common element of upper and lower list, but doing this in general
> requires some effort and may be undesirable, so it is not now.
> Doing this for nopositive integers in the lower list is easier, but
> it is not clear to me if this is really useful.
>

No, no, let's forget all about it.

[Dima Pasechnik]

On Wed, Jul 8, 2026 at 6:44 AM Waldek Hebisch <de...@fricas.org> wrote:
>

I am not sure what you mean by "C memory management".
Do you mean to say you don't know how to overload GMP's memory management,
as CASes which use GMP normally do?
https://gmplib.org/manual/Custom-Allocation

> That make them hard to use for us. And in important
> special case of "small bignums" higher level GMP functions are
> slow (much slower than low level ones). There are probably about
> 20 of GMP functions that would be useful to us, but which GMP
> does not expose, so w can not use them.
>
> Another problem with libraries is dependency management. GMP is
> widely available now, so with GMP this is not a big problem.
> But other libraries are much more problematic. And simply number
> of needed libraries may be a problem: it is very discouraging
> if before building what you want you need to buid first about
> 50 other things (that count correspond to tex4ht several years
> ago).

mpfr, mpc, flint, gsl are widely available...

[Kurt Pagani]

>
> mpfr, mpc, flint, gsl are widely available...
>

I guess the philosophy of fricas is to stand on its own (two) feet as
far as possible, what's not bad at all. It's at the user's discretion to
interface whatever he/she likes. IMO there's no need for another
composite ;)

When using SBCL I don't see the need of GMP anyway as it features its
own internal, native assembly-optimized implementation of
arbitrary-precision integers (bignums) and rationals. I guess there was
a GSOC 2014 proposal https://www.sbcl.org/gsoc2014/ideas/#sec-1.2 which
eventually produced
ls /usr/lib/sbcl/contrib/sb-gmp.*
/usr/lib/sbcl/contrib/sb-gmp.asd /usr/lib/sbcl/contrib/sb-gmp.fasl
There's /usr/lib/sbcl/contrib/sb-mpfr.asd as well.

I assume SBCL's choice to remain independent from GMP was a deliberate
design choice to reduce C-runtime dependencies.

Your POV seems to be ECL where it's fairly easy to interface C-libs,
however, most users I know prefer SBCL.

Regarding GSL I would make an exception when using GSLL because it is CL
and wasn't touched for many years
(https://gitlab.common-lisp.net/antik/gsll). Unfortunately it lacks just
the functionality of the topic we're talking about here. GSLL can only
compute Bessel functions of real order and argument.

[Qian Yun]

On 7/10/26 7:36 PM, Kurt Pagani wrote:
>
> Regarding GSL I would make an exception when using GSLL because it is CL
> and wasn't touched for many years (https://gitlab.common-lisp.net/antik/
> gsll). Unfortunately it lacks just the functionality of the topic we're
> talking about here. GSLL can only compute Bessel functions of real order
> and argument.
>

Off topic a bit, has the development of GSL been stalled?
Looks like less than a dozen commits per year.

https://cgit.git.savannah.gnu.org/cgit/gsl.git/log/

- Qian

[Dima Pasechnik]

On July 10, 2026 6:36:32 AM CDT, Kurt Pagani <nil...@gmail.com> wrote:
>
>>
>> mpfr, mpc, flint, gsl are widely available...
>>
>
>I guess the philosophy of fricas is to stand on its own (two) feet as far as possible, what's not bad at all.

arguably it's a choice to remain inefficient, and to neglect work of experts in the domain which is far from the core interests of the project. An amateur's choice.

> It's at the user's discretion to interface whatever he/she likes. IMO there's no need for another composite ;)
>
>When using SBCL I don't see the need of GMP anyway as it features its own internal, native assembly-optimized implementation of arbitrary-precision integers (bignums) and rationals. I guess there was a GSOC 2014 proposal https://www.sbcl.org/gsoc2014/ideas/#sec-1.2 which eventually produced
>ls /usr/lib/sbcl/contrib/sb-gmp.*
>/usr/lib/sbcl/contrib/sb-gmp.asd /usr/lib/sbcl/contrib/sb-gmp.fasl
>There's /usr/lib/sbcl/contrib/sb-mpfr.asd as well.
>
>I assume SBCL's choice to remain independent from GMP was a deliberate design choice to reduce C-runtime dependencies.

Perhaps it's more of a license-based decision - unwillingness to use GPL?

Everybody in the Lisp world is certainly entitled to making bad choices, and thus to contribute to demise of Lisp systems as competitive modern environments :-)

[Kurt Pagani]

Regarding GSLL (https://gsll.common-lisp.dev/)

Status

GSLL is largely complete and usable, with functioning interfaces to most
of GSL.


GSL itself:

https://github.com/ampl/GSL

Current State of Development

Latest Release: The current stable version is GSL 2.8, which introduced
up-to-date source refinements and algorithm patches.

Security & Maintenance: The development team actively monitors and
responds to security vulnerabilities. When an integer signedness issue
(CVE-2024-50610) was discovered in late 2024 within the library's solver
mechanics, it was promptly addressed with vendor-provided patches and
mitigations.

Core Philosophy: Because thousands of physics, engineering, and data
science projects depend on GSL as a backend engine, the maintainers
prioritize ensuring that existing APIs do not break. Development focuses
on adding heavily tested, peer-reviewed mathematical algorithms rather
than rewriting the library's core structure.

[Kurt Pagani]

On 10/07/2026 14:37, Dima Pasechnik wrote:
>
>
> On July 10, 2026 6:36:32 AM CDT, Kurt Pagani <nil...@gmail.com> wrote:
>>
>>>
>>> mpfr, mpc, flint, gsl are widely available...
>>>
>>
>> I guess the philosophy of fricas is to stand on its own (two) feet as far as possible, what's not bad at all.
>
> arguably it's a choice to remain inefficient, and to neglect work of experts in the domain which is far from the core interests of the project. An amateur's choice.

I don't think so. You're free to implement your domain knowledge in SPAD
as much as you like, but tell me why one should interface other systems
when I can use them in their native environment? IMO one of the big
advantages of fricas is the ease the algorithms can be reenacted. And do
not forget that fricas has no armada of developers to build the
universal cas ...

Of course, I respect your opinion and know that you are not alone with
it, however, I'm contented with the status quo :)

[Dima Pasechnik]

On Fri, Jul 10, 2026 at 8:27 AM Kurt Pagani <nil...@gmail.com> wrote:
>
>
>
> On 10/07/2026 14:37, Dima Pasechnik wrote:
> >
> >
> > On July 10, 2026 6:36:32 AM CDT, Kurt Pagani <nil...@gmail.com> wrote:
> >>
> >>>
> >>> mpfr, mpc, flint, gsl are widely available...
> >>>
> >>
> >> I guess the philosophy of fricas is to stand on its own (two) feet as far as possible, what's not bad at all.
> >
> > arguably it's a choice to remain inefficient, and to neglect work of experts in the domain which is far from the core interests of the project. An amateur's choice.
>
> I don't think so. You're free to implement your domain knowledge in SPAD
> as much as you like, but tell me why one should interface other systems
> when I can use them in their native environment?

Huh, you never do anything across different maths domains, right?
Otherwise, sure, output intermediate results into text files, edit
them to be able to readable by the other system, read them into the
other system, rinse and repeat.
Very efficient :-)

> IMO one of the big
> advantages of fricas is the ease the algorithms can be reenacted. And do
> not forget that fricas has no armada of developers to build the
> universal cas ...
>
> Of course, I respect your opinion and know that you are not alone with
> it, however, I'm contented with the status quo :)
>
>
>
> >
> >
> >> It's at the user's discretion to interface whatever he/she likes. IMO there's no need for another composite ;)
> >>
> >> When using SBCL I don't see the need of GMP anyway as it features its own internal, native assembly-optimized implementation of arbitrary-precision integers (bignums) and rationals. I guess there was a GSOC 2014 proposal https://www.sbcl.org/gsoc2014/ideas/#sec-1.2 which eventually produced
> >> ls /usr/lib/sbcl/contrib/sb-gmp.*
> >> /usr/lib/sbcl/contrib/sb-gmp.asd /usr/lib/sbcl/contrib/sb-gmp.fasl
> >> There's /usr/lib/sbcl/contrib/sb-mpfr.asd as well.
> >>
> >> I assume SBCL's choice to remain independent from GMP was a deliberate design choice to reduce C-runtime dependencies.
> >
> > Perhaps it's more of a license-based decision - unwillingness to use GPL?
> >
> > Everybody in the Lisp world is certainly entitled to making bad choices, and thus to contribute to demise of Lisp systems as competitive modern environments :-)
> >
> >
> >>
> >> Your POV seems to be ECL where it's fairly easy to interface C-libs, however, most users I know prefer SBCL.
> >>
> >> Regarding GSL I would make an exception when using GSLL because it is CL and wasn't touched for many years (https://gitlab.common-lisp.net/antik/gsll). Unfortunately it lacks just the functionality of the topic we're talking about here. GSLL can only compute Bessel functions of real order and argument.
> >>
> >
>

> --
> You received this message because you are subscribed to the Google Groups "FriCAS - computer algebra system" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to fricas-devel...@googlegroups.com.

> To view this discussion visit https://groups.google.com/d/msgid/fricas-devel/56cda29e-b91b-46ab-8eeb-225078079600%40gmail.com.

[Kurt Pagani]

On 10/07/2026 16:59, Dima Pasechnik wrote:
> On Fri, Jul 10, 2026 at 8:27 AM Kurt Pagani <nil...@gmail.com> wrote:
>>
>>
>>
>> On 10/07/2026 14:37, Dima Pasechnik wrote:
>>>
>>>
>>> On July 10, 2026 6:36:32 AM CDT, Kurt Pagani <nil...@gmail.com> wrote:
>>>>
>>>>>
>>>>> mpfr, mpc, flint, gsl are widely available...
>>>>>
>>>>
>>>> I guess the philosophy of fricas is to stand on its own (two) feet as far as possible, what's not bad at all.
>>>
>>> arguably it's a choice to remain inefficient, and to neglect work of experts in the domain which is far from the core interests of the project. An amateur's choice.
>>
>> I don't think so. You're free to implement your domain knowledge in SPAD
>> as much as you like, but tell me why one should interface other systems
>> when I can use them in their native environment?
>
> Huh, you never do anything across different maths domains, right?
> Otherwise, sure, output intermediate results into text files, edit
> them to be able to readable by the other system, read them into the
> other system, rinse and repeat.

Please, don't understand me wrong, I wrote/use a lot of interfaces to
fricas privately (all on gh, public), but I do not want/urge other
people to depend on / or even use them.

> Very efficient :-)

I'm an amateur and need not be efficient by the way.