Expanding definite integrals in a power series about a parameter
- To: mathgroup at smc.vnet.net
- Subject: [mg35444] Expanding definite integrals in a power series about a parameter
- From: "Carl K. Woll" <carlw at u.washington.edu>
- Date: Fri, 12 Jul 2002 04:29:17 -0400 (EDT)
- Sender: owner-wri-mathgroup at wolfram.com
Hi newsgroup,
I have developed a package called SeriesSolve that can do two neat
things (I think). There are a few bugs in the package, so that it is
not quite ready for MathSource. The package works sufficiently well
for me as it is, so I am checking with you guys to see if the
functionality provided by the package would interest you.
First, given a polynomial with coefficients which are themselve power
series in another variable, the function SeriesSolve will find the
roots of the polynomial in a series in that variable. As a very simple
example, suppose you have the polynomial (a+a^2)+x+a x^2+(1+a^2)x^3 and
are interested in the roots of this polynomial. Using SeriesSolve, this
is simple (you will need to view all of the following formulas in a
fixed font to understand what's going on):
rt=SeriesSolve[(a+a^2)+x+a x^2+(1+a^2)x^3,x,{a,0,5}]//Chop
2 4 5 6
{-1. a - 1. a + 1. a + 3. a + O[a] ,
2 3
1. I + (0.5 - 0.5 I) a + (0.5 + 0.5 I) a -
4 5 6
(0.5 - 1.25 I) a - 2. a + O[a] ,
2 3
-1. I + (0.5 + 0.5 I) a + (0.5 - 0.5 I) a -
4 5 6
(0.5 + 1.25 I) a - 2. a + O[a] }
Let's check that the putative solutions are indeed correct
In[4]:=
(1+a^2)Times@@(x-rt)//ExpandAll//Chop
Collect[%,x]
Out[4]=
3 2 3 2 6
(1. x + x ) + (1. + 1. x ) a + (1. + x ) a + O[a]
Out[5]=
2 2 2 3
1. a + 1. a + 1. x + 1. a x + (1 + a ) x
Admittedly, it is not difficult to produce the above output for low order
power series in a. However, for high order power series, for power series
with multiple roots (at lowest order) or for power series
where the coefficient of the highest order term in the polynomial vanishes
at leading order in the power series, it does become difficult to produce
the roots of the polynomials. My package can handle all of these
situations. For a slightly more difficult example, consider
In[11]:=
rt=SeriesSolve[(1+a)x^3+a x^2+a^3 x+a^4,x,{a,0,5}]//Chop
Out[11]=
3/2 3 7/2 4
{1. I a + 0.5 a - 0.5 I a - 0.5 a -
9/2 5 11/2 6
0.125 I a - 1. a + 2.125 I a + O[a] ,
3/2 3 7/2 4
-1. I a + 0.5 a + 0.5 I a - 0.5 a +
9/2 5 11/2 6
0.125 I a - 1. a - 2.125 I a + O[a] ,
2 3 4 5 6
-1. a + 1. a - 2. a + 2. a + 1. a + O[a] }
and checking the solution:
In[12]:=
(1+a)Times@@(x-rt)//ExpandAll//Chop
Out[12]=
3 2 3 3 4 6
x + (1. x + x ) a + 1. x a + 1. a + O[a]
I hope that gives you a sufficient idea about how SeriesSolve works.
In particular, when the coefficient of the highest order term vanishes
at leading order in the expansion, it becomes quite difficult to find
the series form of the roots.
You might wonder why I developed this package. Consider evaluating an
integral of the form
Integrate[ q^m / poly^n, {q,0,Infinity}]
where poly stands for some polynomial in q with coefficients which are
themselves polynomials in another variable, say a. Then, using the
calculus of residues it is possible to evaluate these integrals, and in
particular one can find the series expansion in a of the integral. This
is exactly what the function SeriesRationalIntegrate does. For example,
if the polynomial is a+a q+q^2+q^3, then let
In[13]:=
int=SeriesRationalIntegrate[a+a q+q^2+q^3,q,{a,0,5}];
where I put a semicolon at the end as the output of SeriesRationalIntegrate
is rather unappetizing. Now, int[n,m] will return the series expansion of
the integral
Integrate[q^m/(a+a q+q^2+q^3)^n,{q,0,Infinity}]
Note that when a=0, the above integral is divergent. Hence, you can't simply
expand
the integrand in a series in a and then do term by term integration. At any
rate, we
find for example,
In[14]:=
int[1,1/2]
Out[14]=
2.22144 1/4 3/4
------- - 3.14159 + 2.22144 a - 2.22144 a +
1/4
a
5/4 7/4 2
3.14159 a - 2.22144 a + 2.22144 a - 3.14159 a +
9/4 11/4 3
2.22144 a - 2.22144 a + 3.14159 a -
13/4 15/4 4
2.22144 a + 2.22144 a - 3.14159 a +
17/4 19/4 5 21/4
2.22144 a - 2.22144 a + 3.14159 a + O[a]
Deriving this series expansion is not at all trivial.
Of course, there may not be much interest in doing a series on a root
object (SeriesSolve) or doing a series on a particular kind of definite
integral (SeriesRationalIntegrate). If there is, I will polish up the
package and post it sometime.
Carl Woll
Physics Dept
U of Washington