integrate

47 views
Skip to first unread message

Kjetil brinchmann Halvorsen

unread,
Jun 23, 2012, 8:47:57 PM6/23/12
to sy...@googlegroups.com
Here is one other integral sympy cannot do, but which shouls not be to
difficult, it should be possible to do it symbolically,
the last integral in the following output:


In [1]: t1,t2,t3 = symbols('t1 t2 t3', real=True)

In [2]: beta=symbols('beta', real=True)

In [3]: integrate( (t1*t2*t3)**beta * (t1-t2)*(t1-t3)*(t2-t3),
(t1,0,1),(t2,0,1),(t3,0,1))
Out[3]: 0

In [4]: integrate( (t1*t2*t3)**beta * abs((t1-t2)*(t1-t3)*(t2-t3)),
(t1,0,1),(t2,0,1),(t3,0,1))
Out[4]:
1 1 1
⌠ ⌠ ⌠
⎮ ⎮ ⎮ β
⎮ ⎮ ⎮ (t₁⋅t₂⋅t₃) ⋅│(t₁ - t₂)⋅(t₁ - t₃)⋅(t₂ - t₃)│ d(t₁) d(t₂) d(t₃)
⌡ ⌡ ⌡
0 0 0


Kjetil

--
"If you want a picture of the future - imagine a boot stamping on the
human face - forever."

George Orwell (1984)

Aaron Meurer

unread,
Jun 24, 2012, 1:10:05 AM6/24/12
to sy...@googlegroups.com
Yes, unfortunately SymPy currently does not do so well with integrals
of absolute values.

Aaron Meurer

On Jun 23, 2012, at 6:47 PM, Kjetil brinchmann Halvorsen
> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To post to this group, send email to sy...@googlegroups.com.
> To unsubscribe from this group, send email to sympy+un...@googlegroups.com.
> For more options, visit this group at http://groups.google.com/group/sympy?hl=en.
>

Vinzent Steinberg

unread,
Jun 25, 2012, 3:11:08 PM6/25/12
to sy...@googlegroups.com
Am Sonntag, 24. Juni 2012 07:10:05 UTC+2 schrieb Aaron Meurer:
Yes, unfortunately SymPy currently does not do so well with integrals
of absolute values.

Actually you can work around this limitation by using piecewise functions:

In [12]: myabs = lambda x: Piecewise((x, x>=0), (-x, x<0))

In [13]: integrate( (t1*t2*t3)**beta * myabs((t1-t2)*(t1-t3)*(t2-t3)), 
(t1,0,1),(t2,0,1),(t3,0,1)).doit()
Out[13]: 0

Maybe we should implement something like abs(x).rewrite(Piecewise).

Vinzent

Vinzent Steinberg

unread,
Jun 25, 2012, 3:12:24 PM6/25/12
to sy...@googlegroups.com
Am Montag, 25. Juni 2012 21:11:08 UTC+2 schrieb Vinzent Steinberg:
Am Sonntag, 24. Juni 2012 07:10:05 UTC+2 schrieb Aaron Meurer:
Yes, unfortunately SymPy currently does not do so well with integrals
of absolute values.

Actually you can work around this limitation by using piecewise functions:

In [12]: myabs = lambda x: Piecewise((x, x>=0), (-x, x<0))

In [13]: integrate( (t1*t2*t3)**beta * myabs((t1-t2)*(t1-t3)*(t2-t3)), 
(t1,0,1),(t2,0,1),(t3,0,1)).doit()
Out[13]: 0

(The .doit() is not necessary, sorry.)

Vinzent 

Kjetil brinchmann Halvorsen

unread,
Jun 25, 2012, 3:14:38 PM6/25/12
to sy...@googlegroups.com
see inline!

On Mon, Jun 25, 2012 at 3:11 PM, Vinzent Steinberg
<vinzent....@googlemail.com> wrote:
> Am Sonntag, 24. Juni 2012 07:10:05 UTC+2 schrieb Aaron Meurer:
>>
>> Yes, unfortunately SymPy currently does not do so well with integrals
>> of absolute values.
>
>
> Actually you can work around this limitation by using piecewise functions:
>
> In [12]: myabs = lambda x: Piecewise((x, x>=0), (-x, x<0))
>
> In [13]: integrate( (t1*t2*t3)**beta * myabs((t1-t2)*(t1-t3)*(t2-t3)),
> (t1,0,1),(t2,0,1),(t3,0,1)).doit()
> Out[13]: 0

But that answer 0 is wrong!

Kjetil


>
> Maybe we should implement something like abs(x).rewrite(Piecewise).
>
> Vinzent
>
> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To view this discussion on the web visit
> https://groups.google.com/d/msg/sympy/-/FZ1u8fmryOMJ.
>
> To post to this group, send email to sy...@googlegroups.com.
> To unsubscribe from this group, send email to
> sympy+un...@googlegroups.com.
> For more options, visit this group at
> http://groups.google.com/group/sympy?hl=en.



Vinzent Steinberg

unread,
Jun 29, 2012, 5:04:00 AM6/29/12
to sy...@googlegroups.com
Am Montag, 25. Juni 2012 21:14:38 UTC+2 schrieb kjetil1001:
> Actually you can work around this limitation by using piecewise functions:
>
> In [12]: myabs = lambda x: Piecewise((x, x>=0), (-x, x<0))
>
> In [13]: integrate( (t1*t2*t3)**beta * myabs((t1-t2)*(t1-t3)*(t2-t3)),
> (t1,0,1),(t2,0,1),(t3,0,1)).doit()
> Out[13]: 0

But that answer 0 is wrong!

You are right, it should not be 0, this is a bug. Do you know the correct value?

Vinzent

Kjetil brinchmann Halvorsen

unread,
Jun 29, 2012, 10:11:01 AM6/29/12
to sy...@googlegroups.com
see inline.
Yes. This is a special case of the so-called Selberg integral (named
for Atle Selberg) , see
http://en.wikipedia.org/wiki/Selberg_integral

The special case above has the value: (G is the gamma function)

S_3(beta+1,1,gamma=1/2) = \prod_{j=1}^3 \frac{ G(beta+1+(j-1)gamma)
G(1+(j-1)gamma) G(1+j gamma) }

{ G(beta+2+(3+j-2)gamma) G(1+gamma) }


Kjetil



> Vinzent
>
> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To view this discussion on the web visit
> https://groups.google.com/d/msg/sympy/-/V4hsllBr76sJ.

Vinzent Steinberg

unread,
Jul 1, 2012, 8:28:04 AM7/1/12
to sy...@googlegroups.com
Am Freitag, 29. Juni 2012 16:11:01 UTC+2 schrieb kjetil1001:
> You are right, it should not be 0, this is a bug. Do you know the correct
> value?
>

Yes. This is a special case of the so-called Selberg integral (named
for Atle Selberg) ,  see
http://en.wikipedia.org/wiki/Selberg_integral

The special case above has the value: (G is the gamma function)

S_3(beta+1,1,gamma=1/2) = \prod_{j=1}^3  \frac{ G(beta+1+(j-1)gamma)
G(1+(j-1)gamma) G(1+j gamma) }

{ G(beta+2+(3+j-2)gamma) G(1+gamma) }

Thanks, I made this issue 3317 [1].

Vinzent


Reply all
Reply to author
Forward
0 new messages