An incorrect graph (computation ?) I do not undertstand

[Emmanuel Charpentier]

I have this (very naïvely defined) Python function :

def dunif(x,a,b):
    if(x<a):return(0)
    elif(x>b):return(0)
    else:return(1/(b-a))

I can plot that easily

plot(lambda x:dunif(x,-1,1), (x,-10,10), figsize=3)

But I can't seem to plot *correctly* its numerical integral :

plot(lambda x:numerical_integral(lambda t:dunif(t,-1,1),-oo,x)[0], (x,-10,10), figsize=3)

The problem is not a marginal display : I can zoom on (one of) the problematic zone(s) :

lot(lambda x:numerical_integral(lambda t:dunif(t,-1,1),-oo,x)[0], (x,4,6), figsize=3)

It is not a problem of "too much lambda" either : I have the same problem after defining :

def punif(x,a,b):return(numerical_integral(lambda t:dunif(t,a,b),-oo,x)[0])

plot(lambda x:punif(x,-1,1), (x,4,6), figsize=3)

Same difference... I do not see an "obvious" mistake in what I coded ; I'm aware that it's extremely ineffcient, that type conversions will occur constantly and return a variety of types, etc... But this was only a prototype for something more complicated, and I do not understand the absurdity that sage spews on my screen in this absurdly simple case ("Hey, sage : spewing absurdities is my job" :).

Note : what I pasted here comes from a 6.7beta0 ipython notebook. But I reproduced it with 6.7beta0 notebook, 6.7beta0 command line, 6.7beta0 via emacs, and with 6.6 (final) on notebook, ipython-notebook and command line.

Can some kind soul point me in the direction of my error(s), possibly with an ndication of its (their) magnituds(s) ?

Sincerely,

--
Emmanuel Charpentier

[Emmanuel Charpentier]

Never mind the plots (you can paste plots in the Web editor, but they do not record...).

Instead, see this minimal example on Sagemath cloud.

Any idea ?

--
Emmanuel Charpentier

[William Stein]

On Thu, Apr 16, 2015 at 3:27 AM, Emmanuel Charpentier
<emanuel.c...@gmail.com> wrote:
> Never mind the plots (you can paste plots in the Web editor, but they do not
> record...).
>

> Instead, see this minimal example on Sagemath cloud:

Hi -- content is private by default on SageMathCloud (unlike github).
I've taken the liberty to make just that one worksheet public

https://cloud.sagemath.com/projects/ec50c702-c618-4167-be42-ec4381d49144/files/Minimal%20numerical%20integration%20problem%20demonstration.sagews

by clicking "i", then "Share publicly".

>
> Any idea ?
>
> --
> Emmanuel Charpentier
>

> --
> You received this message because you are subscribed to the Google Groups
> "sage-support" group.
> To unsubscribe from this group and stop receiving emails from it, send an
> email to sage-support...@googlegroups.com.
> To post to this group, send email to sage-s...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sage-support.
> For more options, visit https://groups.google.com/d/optout.



--
William (http://wstein.org)

[Emmanuel Charpentier]

Thank you very much, William. The interface of Sagemath cloud is a bit ... cloudy... for me.

I'm pursuing my problem. The numerical_integral docstring shows that this function has arguments looking suspiciously like Maxima's quad_qag*, interface to quadpack. And indeed, Maxima has the same problem (documentation on the above worksheet to come).

Now sympy (part of sage) has a nice "quad" function for (low-dimension) quadrature. If I could figure a way to use it from sage...

Are you aware of other numerical integration/quadrature function in Sagemath's components (allowing for infinite bounds) that could be used for this purpose ?

To be continued,

--
Emmanuel Charpentier

[William Stein]

On Thu, Apr 16, 2015 at 5:45 AM, Emmanuel Charpentier
<emanuel.c...@gmail.com> wrote:
> Thank you very much, William. The interface of Sagemath cloud is a bit ...
> cloudy... for me.
>
> I'm pursuing my problem. The numerical_integral docstring shows that this
> function has arguments looking suspiciously like Maxima's quad_qag*,
> interface to quadpack. And indeed, Maxima has the same problem
> (documentation on the above worksheet to come).
>
> Now sympy (part of sage) has a nice "quad" function for (low-dimension)
> quadrature. If I could figure a way to use it from sage...

You should be able to just use sympy from Sage, following any of their
examples/documentation:

http://docs.sympy.org/dev/modules/mpmath/calculus/integration.html

??

[Emmanuel Charpentier]

Indeed, using quad (hich is from mpbath, by the way (by bad !)) is easy. But it doesn't solve my problem : it worsens it. I must catch special cases of conversions from rationals to reals. See the updated worksheet (still public).

--
Emmanuel Charpentier
Contemplating reverting to dental surgery (simpler...)

[Emmanuel Charpentier]

FWIW, I found variations of the same problematic behaviour, so far in sahe(-maxima), mpmath and R. See this worksheet for examples...

--
Emmanuel Charpentier
Where did I put those syndesmotomes ?