mmcdara

4029 Reputation

17 Badges

6 years, 186 days

MaplePrimes Activity


These are questions asked by mmcdara

@ASHAN's question was about a linestyle which wasn't displayed as expected.

Maybe a handling error on my part?

To ASHAN:

Does anyone have any idea why Maple can obtain a closed form of 

Int(f, t=1..3);

but doesn't for 

Int(g, t=1..3);

?

The Int forms are quite close and I don't understand what makes Maple's task that difficult in the second case (the issue seems to come from the conditions in the piecewise function).
Can we force Maple to perform the second integration ?

TIA

restart:

interface(Version)

`Standard Worksheet Interface, Maple 2015.2, Mac OS X, December 21 2015 Build ID 1097895`

(1)

f := piecewise(z/t < 0, 0, z/t < 1, 630*z^8/t^8-2520*z^7/t^7+3780*z^6/t^6-2520*z^5/t^5+630*z^4/t^4, 0)/t:

Int(f, t=1..3);

value(%);  # returns a closed form of the integral

g := piecewise(z/t < 1, 0, z/t < 2, 10080-60480*z/t+156240*z^2/t^2-226800*z^3/t^3+202230*z^4/t^4-113400*z^5/t^5+39060*z^6/t^6-7560*z^7/t^7+630*z^8/t^8, 0)/t:

Int(g, t=1..3);

value(%): # unable to return a closed form of the integral

Download MyIntegral.mw

I'm working on finding the analytic expression of the PDF of a sum of abstract Uniform Random Variables URV).
Here "abstract" means that the supports are not numeric but litteral.

Maple is capable to find such a PDF for numeric supports but unable to determine the PDF of U1+U2 where 

U1 := RandomVariable(Uniform(a1, b1)):
U2 := RandomVariable(Uniform(a2, b2)):

A way to deal with abstract URV is to complute explicitely the convolution product of th PDFs.
Ir seems that this fails (MAPLE 2015.2) for these PDF are piecewise functions.
A workaround is to convert them first into Heaviside(s).

One done the explicit expression of the convolution product can be obtained for a sum of 2 abstract URVs, but not for a sum of a larger number of abstract URVs.

Thus the second workaround which consists in using direct and inverse Fourier transform.

The question is :
obviously, the PDF f(t ; a1...aN, b1...bN) of  U1 + ... + UN is a continuous function of t: why does discont(f(t ; ...), t) returns the non empty set of the values where the Heaviside functions are undefined ?

Here is a very simple result

u1 := RandomVariable(Uniform(-1, 1)): 
u2 := RandomVariable(Uniform(-1, 1)):
p := PDF(u1+u2, t):
discont(p, t);

print("-------------------------------------");

f1 := convert(PDF(u1, t), Heaviside):
f2 := convert(PDF(u2, t), Heaviside):
g1 := fourier(f1, t, xi):
g2 := fourier(f2, t, xi):
g  := g1*g2:
f  := invfourier(g, xi, t):
discont(f, t);

                               {}
            "-------------------------------------"
                           {-2, 0, 2}

The problem is (IMO) that discont(f, t) should return { }, but that some function (does it exists) should say that d is undefined at points t=-2, t=0, t=2.
The output of discont(f, t) doesn't seem consistent with the definition of the continuity

limit(f, t=-2, left);
limit(f, t=-2, right);
eval(f, t=-2)
                               0
                               0
                           undefined

restart:

with(inttrans);
with(Statistics):

[addtable, fourier, fouriercos, fouriersin, hankel, hilbert, invfourier, invhilbert, invlaplace, invmellin, laplace, mellin, savetable]

(1)

N := 3:
for n from 1 to N do
  U||n := RandomVariable(Uniform(a__||n, b__||n)):
end do;

_R

 

_R0

 

_R1

(2)

# Maple fails to compute the PDF of a sum of abstract (meaning with symbolic support) uniform RVS
# PDF(U1+U2, t)

# For N <=3 the computation of the convolution product is possible.
#
# For N > 4 (Maple 2015) it seems this is no longer the case. The trick used here is based on
# the fact that Fourier(Conv(f, g)) : Fourier(f)*Fourier(g)
# Thus PDF(U1+U2) = conv(PDF(U1), PDF(U2)) = invFourier(Fourier(PDF(U1)).Fourier(PDF(U1)))


for n from 1 to N do
  f||n := convert(PDF(U||n, t), Heaviside):
end do:

for n from 1 to N do
  g||n := fourier(f||n, t, xi):
end do:
g := mul(g||n, n=1..N):

hyp := seq(b__||n > a__||n, n=1..N);
f   := invfourier(g, xi, t) assuming hyp;

a__1 < b__1, a__2 < b__2, a__3 < b__3

 

((1/2)*(t-a__1-a__2-a__3)^2*Heaviside(-t+a__1+a__2+a__3)-(1/2)*(t-a__1-a__2-b__3)^2*Heaviside(-t+a__1+a__2+b__3)-(1/2)*(t-a__1-a__3-b__2)^2*Heaviside(-t+a__1+a__3+b__2)+(1/2)*(t-a__1-b__2-b__3)^2*Heaviside(-t+a__1+b__2+b__3)-(1/2)*(t-a__2-a__3-b__1)^2*Heaviside(-t+a__2+a__3+b__1)+(1/2)*(t-a__2-b__1-b__3)^2*Heaviside(-t+a__2+b__1+b__3)+(1/2)*(t-b__1-a__3-b__2)^2*Heaviside(-t+b__1+a__3+b__2)-(1/2)*(t-b__1-b__2-b__3)^2*Heaviside(-t+b__1+b__2+b__3))/((-b__1+a__1)*(-b__2+a__2)*(-b__3+a__3))

(3)

# f is obviously a continuous function of t, but I get this strange result

discont(f, t);

{a__1+a__2+a__3, a__1+a__2+b__3, a__1+a__3+b__2, a__1+b__2+b__3, a__2+a__3+b__1, a__2+b__1+b__3, b__1+a__3+b__2, b__1+b__2+b__3}

(4)

# note that this "error" also appears if the a__n's and b__n's are numeric

r := rand(0. .. 1.):
P := convert(indets(g, name) minus{xi}, list):
E := NULL:
for n from 1 to N do
  a := r():
  b := a + r():
  E := E, a__||n = a, b__||n = b
end do:
E := [E];

G := eval(g, E):
f := invfourier(G, xi, t);
discont(f, t);

[a__1 = .3055679837, b__1 = .7906643786, a__2 = .8311025583, b__2 = .9857257095, a__3 = .4223879539, b__3 = .7034757826]

 

0.1333206533e-8*(-0.1000000000e11*t+0.1136670542e11)*Heaviside(-t+1.136670542)+0.1333206533e-8*(0.1000000000e11*t-0.1291293693e11)*Heaviside(-t+1.291293693)+0.1333206533e-8*(0.1000000000e11*t-0.1621766937e11)*Heaviside(-t+1.621766937)+0.1333206533e-8*(-0.1000000000e11*t+0.1776390088e11)*Heaviside(-t+1.776390088)

 

{1.136670542, 1.291293693, 1.621766937, 1.776390088}

(5)

u1 := RandomVariable(Uniform(-1, 1)):
u2 := RandomVariable(Uniform(-1, 1)):
p := PDF(u1+u2, t):
discont(p, t);

print("-------------------------------------");

f1 := convert(PDF(u1, t), Heaviside):
f2 := convert(PDF(u2, t), Heaviside):
g1 := fourier(f1, t, xi):
g2 := fourier(f2, t, xi):
g  := g1*g2:
f  := invfourier(g, xi, t):
discont(f, t);
 

{}

 

"-------------------------------------"

 

{-2, 0, 2}

(6)

limit(f, t=-2, left);
limit(f, t=-2, right);
eval(f, t=-2);

print("-------------------------------------");

limit(f, t=0, left);
limit(f, t=0, right);
eval(f, t=0)

0

 

0

 

undefined

 

"-------------------------------------"

 

1/2

 

1/2

 

undefined

(7)

 

Download Sum_of_Uniform_RVs.mw

I use to write modules with option package and store them in a .mla files.
Let's assume M is the name of the module and P1, P2, ..., PN are the procedures M contains.
Each time one of these procedures hass to call another one, I write 

Pm := proc(...)
 ...
 a := M:-Pn(...)
 ...
end proc:

instead of simply

Pm := proc(...)
 ...
 a := Pn(...)
 ...
end proc:

It seemed to me that the first method sometimes saved me a lot of trouble

Is it a common/necessary/safer/useless practice to write a:=M:Pn(...) instead of a:=Pn(...) ?

Consider this piece of program

restart:
interface(version)
Standard Worksheet Interface, Maple 2015.2, Mac OS X, December 21 2015 Build ID 1097895
with(Units):
a := 3*Unit('m');
                       3 Units:-Unit('m')

I would like to define a new quantity b which is the dimensionless variant of a.
I browsed the Units package to look for a function that would "remove" the dimension of a dimensional quantity like a above (that is to get '3' alone).

As I couldn't find such a function I use to use this workaround.

b := remove(has, a, Unit);
                               3

Is this a robust strategy?
What would you propose to "remove" the dimension of a?

Thanks in advance

1 2 3 4 5 6 7 Last Page 1 of 34