Working with integration of rational functions I came to the

conclusion

that when RootSum is generated in the course of the computation

strange behavior of Mathematica (at least for me) takes place.

In a recent thread I show one such example.

As another example consider the integral of (x^3 + 6*x + 7)/(x^6 +

4*x^4 + 3*x^3 + 6)

in [0,Infinity).

The indefinite integral obtained by Mathematica is

F=Integrate[(x^3 + 6*x + 7)/(x^6 + 4*x^4 + 3*x^3 + 6), x]

RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[x - #1] + 6*Log[x -

#1]*#1 + Log[x - #1]*#1^3)/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

which is of course correct.

Together[D[F, x]]

(7 + 6*x + x^3)/(6 + 3*x^3 + 4*x^4 + x^6)

Application of the Newton-Leibniz formula to get from the

andiderivative the requested

definite integral fails.

Limit[RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[x - #1] + 6*Log[x

- #1]*#1 + Log[x - #1]*#1^3)/

(9*#1^2 + 16*#1^3 + 6*#1^5) & ], x -> Infinity] -

Limit[RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[x - #1] +

6*Log[x - #1]*#1 + Log[x - #1]*#1^3)/

(9*#1^2 + 16*#1^3 + 6*#1^5) & ], x -> 0]

Limit[RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[x - #1] + 6*Log[x

- #1]*#1 + Log[x - #1]*#1^3)/

(9*#1^2 + 16*#1^3 + 6*#1^5) & ], x -> Infinity] - RootSum[6 +

3*#1^3 + 4*#1^4 + #1^6 & ,

(7*Log[-#1] + 6*Log[-#1]*#1 + Log[-#1]*#1^3)/(9*#1^2 + 16*#1^3 +

6*#1^5) & ]

My question comes now...

Since for integrals like this Mathematica uses (as it appears in a

recent forum) the Newton-Leibniz formula

how it can evaluate the definite integral?

Integrate[(x^3 + 6*x + 7)/(x^6 + 4*x^4 + 3*x^3 + 6), {x, 0, Infinity}]

-((Log[-Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 1]]*(7 + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 1] +

Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 1]^3))/(Root[6 + 3*#1^3 +

4*#1^4 + #1^6 & , 1]^2*

(9 + 16*Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 1] + 6*Root[6 +

3*#1^3 + 4*#1^4 + #1^6 & , 1]^3))) -

(Log[-Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 2]]*(7 + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 2] +

Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 2]^3))/(Root[6 + 3*#1^3 +

4*#1^4 + #1^6 & , 2]^2*

(9 + 16*Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 2] + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 2]^3)) -

(Log[-Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 3]]*(7 + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 3] +

Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 3]^3))/(Root[6 + 3*#1^3 +

4*#1^4 + #1^6 & , 3]^2*

(9 + 16*Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 3] + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 3]^3)) -

(Log[-Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 4]]*(7 + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 4] +

Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 4]^3))/(Root[6 + 3*#1^3 +

4*#1^4 + #1^6 & , 4]^2*

(9 + 16*Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 4] + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 4]^3)) -

(Log[-Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 5]]*(7 + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 5] +

Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 5]^3))/(Root[6 + 3*#1^3 +

4*#1^4 + #1^6 & , 5]^2*

(9 + 16*Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 5] + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 5]^3)) -

(Log[-Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 6]]*(7 + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 6] +

Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 6]^3))/(Root[6 + 3*#1^3 +

4*#1^4 + #1^6 & , 6]^2*

(9 + 16*Root[6 + 3*#1^3 + 4*#1^4 + #1^6 & , 6] + 6*Root[6 + 3*#1^3

+ 4*#1^4 + #1^6 & , 6]^3))

What am I missing here?

Dimitris

ello,

I think, here is nothing very mysterious...

The antiderivative is

F = Integrate[(x^3 + 6*x + 7)/(x^6 + 4*x^4 + 3*x^3 + 6), x]

RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[x - #1] + 6*Log[x - #1]*#1 + Log[x - #1]*#1^3)/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

The limit for x->Infinity is 0 as is easily obtained

RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[x] + 6*Log[x]*#1 + Log[x]*#1^3)/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

0

Therefore, the integral F is given by

-RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[-#1] + 6*Log[-#1]*#1 + Log[-#1]*#1^3)/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

according to Newton-Leibnitz.

That this is the correct result is easily checked numerically.

Directly calculating the definite integral yields

F = Integrate[(x^3 + 6*x + 7)/(x^6 + 4*x^4 + 3*x^3 + 6), {x, 0, Infinity}]

-RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (Log[-#1]*#1)/(9 + 16*#1 + 6*#1^3) & ] -

6*RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , Log[-#1]/(9*#1 + 16*#1^2 + 6*#1^4) & ] -

7*RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , Log[-#1]/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

which is equivalent to the pedestrian Newton-Leibnitz result given above.

The only surprising thing may be the fact that the integrator in this example gets the limit x-> Infinity right, while in some other

examples it does not do that. And also your way to calculate this limit using the build in command Limit[] does not seem to work

properly... But this in not a problem of RootSum, but of Limit, I think...

Regards Michael

"dimitris" < XXXX@XXXXX.COM > schrieb im Newsbeitrag news:evfkhu$73g$ XXXX@XXXXX.COM ...

I think, here is nothing very mysterious...

The antiderivative is

F = Integrate[(x^3 + 6*x + 7)/(x^6 + 4*x^4 + 3*x^3 + 6), x]

RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[x - #1] + 6*Log[x - #1]*#1 + Log[x - #1]*#1^3)/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

The limit for x->Infinity is 0 as is easily obtained

RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[x] + 6*Log[x]*#1 + Log[x]*#1^3)/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

0

Therefore, the integral F is given by

-RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (7*Log[-#1] + 6*Log[-#1]*#1 + Log[-#1]*#1^3)/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

according to Newton-Leibnitz.

That this is the correct result is easily checked numerically.

Directly calculating the definite integral yields

F = Integrate[(x^3 + 6*x + 7)/(x^6 + 4*x^4 + 3*x^3 + 6), {x, 0, Infinity}]

-RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , (Log[-#1]*#1)/(9 + 16*#1 + 6*#1^3) & ] -

6*RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , Log[-#1]/(9*#1 + 16*#1^2 + 6*#1^4) & ] -

7*RootSum[6 + 3*#1^3 + 4*#1^4 + #1^6 & , Log[-#1]/(9*#1^2 + 16*#1^3 + 6*#1^5) & ]

which is equivalent to the pedestrian Newton-Leibnitz result given above.

The only surprising thing may be the fact that the integrator in this example gets the limit x-> Infinity right, while in some other

examples it does not do that. And also your way to calculate this limit using the build in command Limit[] does not seem to work

properly... But this in not a problem of RootSum, but of Limit, I think...

Regards Michael

"dimitris" < XXXX@XXXXX.COM > schrieb im Newsbeitrag news:evfkhu$73g$ XXXX@XXXXX.COM ...

lthough I understand the presence of RootSum in cases of rational

functions P(x)/Q(x) where Q(x) is a polynomial of degree 5 or higher

(recall the Abel's impisibility theorem) I believe that for degree

less than

5 (in particular for 3 and 4) Integrate could not generate RootSum.

In[307];

Quit[]

In[1]:=

$Version

Out[1]=

"5.2 for Microsoft Windows (June 20, 2005)"

Consider the function

In[2]:=

g[x_] := (2*x)/((x + 1)*(x^3 + 3*x^2 + 2*x + 1))

The following analysis show that the function is integrable

in the range [0,Infinity).

In[3]:=

Plot[g[x], {x, 0, 10}]

Out[3]=

Graphics[]

In[4]:=

(g[x] + O[x, #1]^2 & ) /@ {0, Infinity}

Out[4]=

{SeriesData[x, 0, {2}, 1, 2, 1], SeriesData[x, Infinity, {}, 2, 2, 1]}

In[5]:=

Reduce[Denominator[g[x]] == 0 && x > 0, x, Reals]

Out[5]=

False

Mathematica tries really hard but it returns the definite integral

back.

In[6]:=

Timing[Integrate[g[x], {x, 0, Infinity}]]

Out[6]=

{142.078*Second, Integrate[(2*x)/((1 + x)*(1 + 2*x + 3*x^2 + x^3)),

{x, 0, Infinity}]}

Let's obtain the (Mathematica's) antiderivative of g(x)

In[6]:=

G[x_] = Integrate[g[x], x]

Out[6]=

2*(-Log[1 + x] + RootSum[1 + 2*#1 + 3*#1^2 + #1^3 & , (Log[x - #1] +

2*Log[x - #1]*#1 + Log[x - #1]*#1^2)/

(2 + 6*#1 + 3*#1^2) & ])

Note at this point how Integrate works. First Apart is called in order

to write the integrand

as follows

In[8]:=

Apart[g[x]]

Out[8]=

-(2/(1 + x)) + (2*(1 + 2*x + x^2))/(1 + 2*x + 3*x^2 + x^3)

and then Integrate...integrates each term of the sum.

In[9]:=

(Integrate[#1, x] & ) /@ %

Out[9]=

-2*Log[1 + x] + 2*RootSum[1 - #1 + #1^3 & , (Log[1 + x - #1]*#1^2)/(-1

+ 3*#1^2) & ]

Note also that version 5.2 is more careful and doesn't return Infinity

for the definite

integral as version 4.0 does (I think because it considers the

integrals of the parts in the last output, which

both are divergent at Infinity).

So, the generation of the RootSum object in the antiderivative is due

to the failure if the Apart function to do the partial fraction

decomposition of (2*(1 + 2*x + x^2))/(1 + 2*x + 3*x^2 + x^3).

Note that I lack serious knowledge on symbolics algebra and so, I may

miss something fundamental but (although a few months ago I have a

query about this behavior of Apart) I still don't understand why there

are cases that Apart returns its argument back.

Here are two examples

In[12]:=

Apart /@ {x^2/(x^2 + 3*x + 2), x/(x^3 + 2*x^2 + x + 1)}

Out[12]=

{1 + 1/(1 + x) - 4/(2 + x), x/(1 + x + 2*x^2 + x^3)}

It seems that Apart does the PFD provided it can evaluate the roots of

the denominator of the

rational function. In the second case it can't evaluate and that's why

it fails to do the PFD.

Could/Should Apart have an option to make it do something like?

In[13]:=

Simplify[x /. Solve[1 + x + 2*x^2 + x^3 == 0, x]]

Simplify /@ Factor[1 + x + 2*x^2 + x^3, Extension -> %]

Simplify /@ Apart[x/%]

(*output is ommited*)

I think if Apart could do this, then the generation of RootSum objects

could be avoided

(especially for definite integrals where one or both limits of

integration is infinity,

and undoubtfully Integrate faces problems. Of course may be what I ask

is too difficult

and in view of the output of In[13] it seems quite complicated, and in

the very next version

RootSum at infinity does not face problems; If I knew that I would

strop right now the message! If the other CAS I use can already do it

pr

functions P(x)/Q(x) where Q(x) is a polynomial of degree 5 or higher

(recall the Abel's impisibility theorem) I believe that for degree

less than

5 (in particular for 3 and 4) Integrate could not generate RootSum.

In[307];

Quit[]

In[1]:=

$Version

Out[1]=

"5.2 for Microsoft Windows (June 20, 2005)"

Consider the function

In[2]:=

g[x_] := (2*x)/((x + 1)*(x^3 + 3*x^2 + 2*x + 1))

The following analysis show that the function is integrable

in the range [0,Infinity).

In[3]:=

Plot[g[x], {x, 0, 10}]

Out[3]=

Graphics[]

In[4]:=

(g[x] + O[x, #1]^2 & ) /@ {0, Infinity}

Out[4]=

{SeriesData[x, 0, {2}, 1, 2, 1], SeriesData[x, Infinity, {}, 2, 2, 1]}

In[5]:=

Reduce[Denominator[g[x]] == 0 && x > 0, x, Reals]

Out[5]=

False

Mathematica tries really hard but it returns the definite integral

back.

In[6]:=

Timing[Integrate[g[x], {x, 0, Infinity}]]

Out[6]=

{142.078*Second, Integrate[(2*x)/((1 + x)*(1 + 2*x + 3*x^2 + x^3)),

{x, 0, Infinity}]}

Let's obtain the (Mathematica's) antiderivative of g(x)

In[6]:=

G[x_] = Integrate[g[x], x]

Out[6]=

2*(-Log[1 + x] + RootSum[1 + 2*#1 + 3*#1^2 + #1^3 & , (Log[x - #1] +

2*Log[x - #1]*#1 + Log[x - #1]*#1^2)/

(2 + 6*#1 + 3*#1^2) & ])

Note at this point how Integrate works. First Apart is called in order

to write the integrand

as follows

In[8]:=

Apart[g[x]]

Out[8]=

-(2/(1 + x)) + (2*(1 + 2*x + x^2))/(1 + 2*x + 3*x^2 + x^3)

and then Integrate...integrates each term of the sum.

In[9]:=

(Integrate[#1, x] & ) /@ %

Out[9]=

-2*Log[1 + x] + 2*RootSum[1 - #1 + #1^3 & , (Log[1 + x - #1]*#1^2)/(-1

+ 3*#1^2) & ]

Note also that version 5.2 is more careful and doesn't return Infinity

for the definite

integral as version 4.0 does (I think because it considers the

integrals of the parts in the last output, which

both are divergent at Infinity).

So, the generation of the RootSum object in the antiderivative is due

to the failure if the Apart function to do the partial fraction

decomposition of (2*(1 + 2*x + x^2))/(1 + 2*x + 3*x^2 + x^3).

Note that I lack serious knowledge on symbolics algebra and so, I may

miss something fundamental but (although a few months ago I have a

query about this behavior of Apart) I still don't understand why there

are cases that Apart returns its argument back.

Here are two examples

In[12]:=

Apart /@ {x^2/(x^2 + 3*x + 2), x/(x^3 + 2*x^2 + x + 1)}

Out[12]=

{1 + 1/(1 + x) - 4/(2 + x), x/(1 + x + 2*x^2 + x^3)}

It seems that Apart does the PFD provided it can evaluate the roots of

the denominator of the

rational function. In the second case it can't evaluate and that's why

it fails to do the PFD.

Could/Should Apart have an option to make it do something like?

In[13]:=

Simplify[x /. Solve[1 + x + 2*x^2 + x^3 == 0, x]]

Simplify /@ Factor[1 + x + 2*x^2 + x^3, Extension -> %]

Simplify /@ Apart[x/%]

(*output is ommited*)

I think if Apart could do this, then the generation of RootSum objects

could be avoided

(especially for definite integrals where one or both limits of

integration is infinity,

and undoubtfully Integrate faces problems. Of course may be what I ask

is too difficult

and in view of the output of In[13] it seems quite complicated, and in

the very next version

RootSum at infinity does not face problems; If I knew that I would

strop right now the message! If the other CAS I use can already do it

pr

3 post • Page:**1** of **1**