Computing the limit of the lognormal density crashes Sage

[7822f248-ba56-45f1-ab3d-4de7482bdf9f]

Define the lognormal density starting from the normal:

var("mu", domain="real")
var("y, m, s, sigma", domain="positive")
dnorm(x, mu, sigma) = e^(-(x-mu)^2/(2*sigma^2))/(sigma*sqrt(2*pi))
dlnorm(y, mu, sigma) = (dnorm(x, mu, sigma).subs(x==log(y))*abs(diff(log(y),y))).simplify()

Let's try to prove that the limit is 0 at 0 and oo:

sage: dlnorm(y, mu, sigma).limit(y=oo)
0

So far so good. But:

dlnorm(y, mu, sigma).limit(y=0, dir="+")
;;;
;;; Detected access to protected memory, also kwown as 'bus or segmentation fault'.
;;; Jumping to the outermost toplevel prompt
;;;

## Numerous repetitions...

Process Sage erreur de segmentation

This seems analogous to but different from #14677...

This limit also seems problematic in other subsystems:

• Sympy returns "Not implemented"
• libgiac returns "Infinity" (wrong)
• Mathematica doesn't return (but (correctly) returns 0 when used directly).
• Used directly, Maxima asks a lot of questions, and fails:

limit(dlnorm(y, mu, sigma), y, 0);
Is sigma^2-mu positive, negative or zero?

p;
Is mu positive, negative or zero?

p;
Is 2*sigma^2-mu positive, negative or zero?

p;
(%o10) ('limit(%e^-((log(y)-mu)^2/(2*sigma^2))/y,y,0))/(sqrt(2)*sqrt(%pi)
                                                               *sigma)

This failure mode is different from the one seen in Sage; we may have a new bug...

Upstream: Not yet reported upstream; Will do shortly.

CC: @slel

Component: symbolics

Keywords: limit, maxima

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

[nbruin]

comment:1

It looks to me like a domain:complex error again. I get:

sage: maxima_calculus.eval("domain:real")
'real'
sage: dlnorm(y,mu,sigma).limit(y=0, dir="+")
0

That's to say, if we set domain to real, maxima in sage seems to do the right thing immediately. So I think it's basically the same problem as #14677.

[slel]

Changed keywords from none to limit, maxima

[slel]

Description changed:

--- 
+++ 
@@ -1,23 +1,23 @@
-dEFINE THE LOGNORMAL DENSITY STARTING FROM THE NORMAL :
+Define the lognormal density starting from the normal:

var("mu", domain="real") -var("y,m,s,sigma", domain="positive") -dnorm(x,mu,sigma)=e^(-(x-mu)^2/(2sigma^2))/(sigmasqrt(2pi)) -dlnorm(u,mu,sigma)=(dnorm(x,mu,sigma).subs(x==log(y))abs(diff(log(y),y))).simplify() +var("y, m, s, sigma", domain="positive") +dnorm(x, mu, sigma) = e^(-(x-mu)^2/(2sigma^2))/(sigmasqrt(2pi)) +dlnorm(y, mu, sigma) = (dnorm(x, mu, sigma).subs(x==log(y))abs(diff(log(y),y))).simplify()

-Let's try to prove that the limits are 0 at 0 and oo :
+Let's try to prove that the limit is 0 at 0 and oo:

-sage: dlnorm(y,mu,sigma).limit(y=oo) +sage: dlnorm(y, mu, sigma).limit(y=oo) 0

-So far so good. But :
+So far so good. But:

-dlnorm(y,mu,sigma).limit(y=0, dir="+") +dlnorm(y, mu, sigma).limit(y=0, dir="+") ;;; ;;; Detected access to protected memory, also kwown as 'bus or segmentation fault'. ;;; Jumping to the outermost toplevel prompt @@ -25,21 +25,20 @@

Numerous repetitions...

- Process Sage erreur de segmentation

-This seems analogous but different of #14677...
+This seems analogous to but different from #14677...

-This limit seems problematic to other subsystems :
+This limit also seems problematic in other subsystems:

 * Sympy returns "Not implemented"
 * libgiac returns "Infinity" (wrong)
 * Mathematica doesn't return (but (correctly) returns 0 when used directly).
-* When used directly, Maxima asks a lot of questions, and fails :
+* Used directly, Maxima asks a lot of questions, and fails:

-limit(dlnorm(y,mu,sigma),y,0); +limit(dlnorm(y, mu, sigma), y, 0); Is sigma^2-mu positive, negative or zero?

p; @@ -53,4 +52,5 @@ *sigma)

-This failure mode is different from the one seen in Sage ; we may have a *new* bug...
+This failure mode is different from the one seen in Sage; we may have a *new* bug...
+