Advanced Search

mmcdara

7022 Reputation

18 Badges

8 years, 230 days
Contact mmcdara

MaplePrimes Activity


These are questions asked by mmcdara

Why am I not able to get the desired int...

mmcdara 7022 Maple 2015
integral integration int
January 02 2024
0 2

There are things that seem simple but rapidly turn into a nightmare.

Here is an example: what I want is to the expression given at equation (4) in the attached file.

Using Int gives a wrong result.
Using int gives a right one but not of the desired form (some double integrals are nested while others are not).

I've been stuck on this problem for hours, can you please help me to fix it?

TIA

> 

restart
> 

use Statistics in
  # For more generality defina an abstract probability distribution.
  AbstractDistribution := proc(N)
    Distribution(
      PDF = (x -> varphi(seq(x[n], n=1..N)))
    )
  end proc:

  # Define two random variables pf AbstractDistribution type.
  X__1 := RandomVariable(AbstractDistribution(2)):
  X__2 := RandomVariable(AbstractDistribution(2)):

end use;

proc (N) Statistics:-Distribution(Statistics:-PDF = (proc (x) options operator, arrow; varphi(seq(x[n], n = 1 .. N)) end proc)) end proc

  

_R

  

_R0

 (1)
> 

F := (U1, U2) -> U1/(U1+U2);
T := mtaylor(F(X__1, X__2), [X__1=1, X__2=1], 2):

proc (U1, U2) options operator, arrow; U1/(U1+U2) end proc

 (2)


Error: x[2] is droped out of the double integral in the rightmost term

> 

use IntegrationTools in

J := eval([op(expand(T))], [seq(X__||i=x[i], i=1..2)]);
L := add(
       map(
         j ->  
         if j::numeric then
           j
         else
           (Expand@CollapseNested)(
             Int(
               j * Statistics:-PDF(X__1, x)
               , seq(x[i]=-infinity..+infinity, i=1..2)
             )
           )
         end if
         , J
       )  
     ):
ET := %
end use;

[1/2, (1/4)*x[1], -(1/4)*x[2]]

  

1/2+(1/4)*(Int(x[1]*varphi(x[1], x[2]), [x[1] = -infinity .. infinity, x[2] = -infinity .. infinity]))-(1/4)*x[2]*(Int(varphi(x[1], x[2]), [x[1] = -infinity .. infinity, x[2] = -infinity .. infinity]))

  

1/2+(1/4)*(Int(x[1]*varphi(x[1], x[2]), [x[1] = -infinity .. infinity, x[2] = -infinity .. infinity]))-(1/4)*x[2]*(Int(varphi(x[1], x[2]), [x[1] = -infinity .. infinity, x[2] = -infinity .. infinity]))

 (3)


I want this

> 

'ET' = 1/2
       +
       (1/4)*(Int(Int(x[1]*varphi(x[1], x[2]), x[1] = -infinity .. infinity), x[2] = -infinity .. infinity))
       -(1/4)*(Int(Int(x[2]*varphi(x[1], x[2]), x[1] = -infinity .. infinity), x[2] = -infinity .. infinity))

ET = 1/2+(1/4)*(Int(Int(x[1]*varphi(x[1], x[2]), x[1] = -infinity .. infinity), x[2] = -infinity .. infinity))-(1/4)*(Int(Int(x[2]*varphi(x[1], x[2]), x[1] = -infinity .. infinity), x[2] = -infinity .. infinity))

 (4)


With int instead of Int one integral is double the other is double-nested

> 

L := add(
       map(
         j ->  
         if j::numeric then
           j
         else
             int(
               j * Statistics:-PDF(X__1, x)
               , seq(x[i]=-infinity..+infinity, i=1..2)
             )
         end if
         , J
       )  
     ):
ET := %

1/2+int(int((1/4)*x[1]*varphi(x[1], x[2]), x[1] = -infinity .. infinity), x[2] = -infinity .. infinity)+int(-(1/4)*x[2]*(int(varphi(x[1], x[2]), x[1] = -infinity .. infinity)), x[2] = -infinity .. infinity)

 (5)


As the expression of ET is now correct, I tried to use IntegrationTools to get the
form I want (equation (4)).

But as soon as I replace int by Int x[2] is again droped out.

So it's not even worth thinking about using CollapseNested!

 

> 

use IntegrationTools in
  eval(ET, int=Int);  
end use;

1/2+Int(Int((1/4)*x[1]*varphi(x[1], x[2]), x[1] = -infinity .. infinity), x[2] = -infinity .. infinity)+Int(-(1/4)*x[2]*(Int(varphi(x[1], x[2]), x[1] = -infinity .. infinity)), x[2] = -infinity .. infinity)

 (6)
> 

 

Download Int_int.mw

Should this be considered a simplify bug...

mmcdara 7022 Maple 2015
December 25 2023
1 6

A case where simplify(...) and simplify~(...) both return the wrong result.
Should we consider this a simplify bug?

> 

restart:


A simple case

> 

J := Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity)
     *
     Int(r[2]^2*varphi[2](r[2]), r[2] = -infinity .. infinity)

(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[2]^2*varphi[2](r[2]), r[2] = -infinity .. infinity))

 (1)
> 

# OK

op(1, J) = simplify(op(1, J))

Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity) = Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity)

 (2)
> 

# OK

op(2, J) = simplify(op(2, J))

Int(r[2]^2*varphi[2](r[2]), r[2] = -infinity .. infinity) = Int(r[2]^2*varphi[2](r[2]), r[2] = -infinity .. infinity)

 (3)
> 

# But...
#
# Not OK

simplify(J)

(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[1]^2*varphi[2](r[1]), r[1] = -infinity .. infinity))

 (4)
> 

# Not OK

simplify~(J)

(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[1]^2*varphi[2](r[1]), r[1] = -infinity .. infinity))

 (5)
> 

# OK

map(simplify, J)

(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[2]^2*varphi[2](r[2]), r[2] = -infinity .. infinity))

 (6)


A slightly more complex case

> 

J := (Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[2]^2*varphi[2](r[2]), r[2] = -infinity .. infinity))-(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[2]*varphi[2](r[2]), r[2] = -infinity .. infinity))^2;

(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[2]^2*varphi[2](r[2]), r[2] = -infinity .. infinity))-(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[2]*varphi[2](r[2]), r[2] = -infinity .. infinity))^2

 (7)
> 

is(J=simplify(J))

false

 (8)
> 

is(J=simplify~(J))

false

 (9)
> 

is(J=map(simplify, J));
map(simplify, J);

false

  

(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[1]^2*varphi[2](r[1]), r[1] = -infinity .. infinity))-(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[1]*varphi[2](r[1]), r[1] = -infinity .. infinity))^2

 (10)
> 

add(map(u -> map(simplify, u), [op(J)]));

is(J=%);

(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[2]^2*varphi[2](r[2]), r[2] = -infinity .. infinity))-(Int(r[1]^2*varphi[1](r[1]), r[1] = -infinity .. infinity))*(Int(r[2]*varphi[2](r[2]), r[2] = -infinity .. infinity))^2

  

true

 (11)
> 

 

Download Simplify_is_wrong.mw

How to define an ordering on partial de...

mmcdara 7022 Maple
ordering sort mtaylor
December 21 2023
1 1

Notional example:
I use mtaylor compute the Taylor expansion of a function f (U) of several variables U1, .., UN.
In the resul the terms are ordered this way:

  • the leftmost term is f (P) where P denotes the point where the expansion is done
  • followed by a succession of terms :
    • firstly ranked according to the total order of derivation of  f.
    • and among terms of same derivation order, ranked in some kind of lexicographic order

For instance

> 

Vars    := [seq(U[i], i=1..2)]:
AtPoint := [seq(P[i], i=1..2)]:
mt      := mtaylor(f(Vars[]), Vars =~ AtPoint, 3)
 

f(P[1], P[2])+(D[1](f))(P[1], P[2])*(U[1]-P[1])+(D[2](f))(P[1], P[2])*(U[2]-P[2])+(1/2)*(D[1, 1](f))(P[1], P[2])*(U[1]-P[1])^2+(D[1, 2](f))(P[1], P[2])*(U[1]-P[1])*(U[2]-P[2])+(1/2)*(D[2, 2](f))(P[1], P[2])*(U[2]-P[2])^2

 (1)
> 

map(t -> op([0, 0], select(has, [op(t)], D)[]), [op(mt)][2..-1])

[D[1], D[2], D[1, 1], D[1, 2], D[2, 2]]

 (2)
 

 

Download mtaylor.mw

How could I define another ordering of the terms in this mtaylor expansion?
For instance 1 being identified to some letter and 2 to another one such that 1 <  2 a lexicographic order would be 

D[1], D[1, 1], D[1, 2], D[2], D[2, 2]

Thanks in advance

Is it possible to go from eq (2) to eq (...

mmcdara 7022 Maple
sum product
December 04 2023
1 1

Is it possible to transform relation (2) into relation (7) without using the hand-made relation (3) and the  Sum -> Int -> Sum trick?

> 

restart
> 

f := Product(x[i]^a*(1-x[i])^b, i)

Product(x[i]^a*(1-x[i])^b, i)

 (1)
> 

Lf := ln(f);

ln(Product(x[i]^a*(1-x[i])^b, i))

 (2)
> 

Sum(ln(x[i]^a*(1-x[i])^b), i)

Sum(ln(x[i]^a*(1-x[i])^b), i)

 (3)
> 

expand(%) assuming x[i] > 0, x[i] < 1, a > 0, b > 0

Sum(a*ln(x[i])+b*ln(1-x[i]), i)

 (4)
> 

eval(%, Sum=Int)

Int(a*ln(x[i])+b*ln(1-x[i]), i)

 (5)
> 

IntegrationTools:-Expand(%);

a*(Int(ln(x[i]), i))+b*(Int(ln(1-x[i]), i))

 (6)
> 

Lf = eval(%, Int=Sum)

ln(Product(x[i]^a*(1-x[i])^b, i)) = a*(Sum(ln(x[i]), i))+b*(Sum(ln(1-x[i]), i))

 (7)
> 

 

Download From2to7.mw

TIA

Can Maple dsolve formally this ode?...

mmcdara 7022 Maple
integration ode dsolve
November 10 2023
0 0

I recently answered a question concerning the Lane-Emden equation (see here LaneEmden) where the main topic was about finding its numerical solution.

The generic form of the Lane-Emden equation with parameter n is

LaneEmden := n -> (Diff(xi^2*(Diff(theta(xi), xi)), xi)) = -theta(xi)^n * xi^2

      d   /  2 / d            \\             n   2
n -> ---- |xi  |---- theta(xi)|| = -theta(xi)  xi 
      dxi \    \ dxi          //


I have just realized that I missed a "small" point in the original question: the OP ( @shashi598 ) wrote
"[...] Maple never comes out of evaluating [the] analytical solution when n=5 [...] ".
The important point here is that this solution (at least for some initial conditions) is known and simple (in the sense it doen't involve any special function).

So I tried for a few hours to verify this claim, and ended wondering myself if it might not be right?

Could you please tell me (I guess @shashi598 would be interested too in your return) if the differential equation LaneEmden(5) can be solved formally?
TIA.

Emden_equation.mw


EDITED:
After a little research it seems that very specigic method are used to build the analytic solution of the LaneEmden(n) (n not equal to 0, 1 and 5): serie expansions, homotopy, Adomian decomposition for instance.
I wasn't capable to find how the solution for LaneEmden(5) have been got for the first time (iseems to be atthe end of the 19th century).

5 6 7 8 9 10 11 Last Page 7 of 45