Package for limits

What is a package that can be used to find the limit of a function? I thought ForwardDiff.jl might have this feature but I couldn’t find it in the documentation.

Here is a “poor man’s” way: stick Inf or -Inf in.

eg,

julia> exp(Inf)
Inf

julia> exp(-Inf)
0.0

EDIT: Sorry, I guess I just assumed you meant limits as the input approached +/- Inf, and that could not be the case. To get more to the point of what you are after, I’m not sure such a thing exists. If the function is continuous, it could be possible.

1 Like

Can you give an example of what you’re looking for?

Something like this? https://github.com/JuliaMath/Richardson.jl

7 Likes

It would be for continuous functions yes, I meant for a tool that substitutes a variable, and rearanges a term such as 1/x or does a limit numerically from each side. Reduce.jl has such a function (but it’s a weird package that won’t load right). I’m thinking a calculus or symbolic algebra package would have these.

Yes, if you want a numerical limit calculation, then Richardson extrapolation would be the first thing that I would try. The Richardson.jl package can handle x → ∞ limits too.

If you want a symbolic limit evaluation, that is a matter for a symbolic-algebra package (e.g. you could try SymPy.jl to call sympy from Julia).

3 Likes

As pointed out, SymPy.jl has the limit function for taking limits. It has a special method for a univariate function. For an example usage, here is an innocuous looking function from a colleague which won’t have a numeric limit:

julia> g(x) = x^2 + 1 + log(abs( 11*x-15 ))/99
g (generic function with 1 method)

julia> limit(g, 15//11)
-∞
1 Like

I’m not sure I know how to read this. What does 15 and 11 mean?

15//11 is a rational type to represent exactly the fraction 15/11, see https://docs.julialang.org/en/v1/manual/complex-and-rational-numbers/#Rational-Numbers

2 Likes

ok thanks.

That would work great, but I do feel silly, I though that subs in SymEngine just changede the value of x, but if I do something like

subs(1x/2x,x=>0)

sure enough, I get 1/2. I’m trying to think of a function with removable discontinuity to test.
I can aslo do one sided limits, using infintesimals

δ=ε=1 * 10^-308 #the smallest number a 64 bit processor can handle
subs(1/x,x=>δ)

and I get a crazy big number, that I can include is infinity.

sub(1/x,x=Inf)

gives me 0.

I was curious, so I gave it a try:

julia> using IntervalArithmetic, Richardson

julia> g(x) = x^2 + 1 + log(abs( 11*x-15 ))/99
g (generic function with 1 method)

julia> extrapolate(g, 15/11)
(1.0273540424380285, 1.4376019930040229e-8)

julia> extrapolate(g, 15//11)
(1.0273540424380285, 1.4376019930040229e-8)

julia> extrapolate(g, big(15)/11)
(-Inf, Inf)

julia> extrapolate(g, interval(15)/11)
([0.846292, ∞], [0, ∞])

julia> extrapolate(g, interval(big(15))/11)
([1.04933, ∞]₂₅₆, [0, ∞]₂₅₆)

julia> g(15/11)
2.5164312852875934

julia> g(15//11)
-Inf

julia> g(big(15)/11)
-Inf

julia> g(interval(15)/11)
[-∞, 2.51644]

Some of them got it, but not all, and the error estimate of Inf from Richardson does not inspire confidence :smile:

3 Likes

Which is from IntervalArithmatic, and which are from Richardson, I like the idea of extrapolate being able to find range, but there is a problem if it makes errors.

extrapolate is from Richardson and interval is from IntervalArithmetic

It is a tricky function as the asymptote is really sharp. But +oo and not -oo?

First, you’re using Richardson.extrapolate incorrectly: The second argument is not the point to extrapolate to, it is the initial step away from the limit point x0 (which defaults to zero).

Second, Richardson extrapolation is based on polynomial extrapolation, so it isn’t going to work if the limit diverges.

Third, this particular extrapolation problem is going to be very difficult to do numerically (as opposed to symbolically) in any case. Because 11/15 is not representable in binary floating point, and the singularity is only logarithmic, in double precision g(11/15) gives a quite small magnitude result of ≈ 1.5573367747276683. (And calling extrapolate(g, 0.1, x0=11/15) “correctly” gives this answer.)

3 Likes

Ah, thanks for the correction!

2 Likes

It is good to have range and domain of many functions. And, it would be good to do limits, I was able to use the subs by an using an infentesimal, which I definite as 10e-302 , or the smallest number a 64 bit machine would do. I’m not sure what is a safe number to put in Julia.

If the limit does not exist, will it give a null or an error?

It will typically halt after a few iterations and give a large error estimate (the second return value):

julia> extrapolate(x -> 1/x, 0.1, x0=0)
(90.0, 80.0)

but log divergences are more difficult to detect with polynomial extrapolation so the error estimates are smaller (but still much larger than the typical error estimates when it succeeds):

julia> extrapolate(log, 0.1, x0=0)
(-744.7275688945767, 1.7346317053750226)

Note that (by design) it computes the one-sided limit (you specify the side using the initial step):

julia> extrapolate(x -> sqrt(x^2 + x^3)/x, 0.1, x0=0)
(1.0000000000000002, 1.3793410857942945e-11)

julia> extrapolate(x -> sqrt(x^2 + x^3)/x, -0.1, x0=0)
(-0.9999999999999998, 1.6191714635738208e-11)
3 Likes