Series Sol'n to a DE with trig coefficient?

[outhsakotad]

Homework Statement

I'm trying to solve this DE: 4xy''+2y'+(cosx)y=0 using a series solution.

Homework Equations

The Attempt at a Solution

See attached PDF. I tried writing the cosine term as an infinite series, but that gives me a messy double summation. And in the end, I get an index that is not an integer... Any hints on how to deal with the cosine term?

 

[jackmell]

My computer is not letting me have access to that document but I got a suggestion for you: muscle through it. But start at the end and work backwards. Since it's singular at x=0, then the solution would be valid for x>0 so say work towards solving:

4xy''+2y'+\cos(x)y=0,\quad y(1)=0,\quad y'(1)=1

and since the roots of the indical equation are 0 and 1/2 and the difference non-integer, the solution is:

y(x)=C_1 \sum_{n=0}^{\infty} a_n x^n+C_2 \sum_{n=0}^{\infty} b_n x^{n+1/2}

and we obtain that solution by substituting each power series into the DE, and letting the first coefficient equal to 1 and compute all the rest recursively even though there's a double sum in there and you just have to learn how to work with a "Cauchy product of double sums".

I'll do c=0 for you, you do the messy one with all the fractional powers of x:

\sum_{n=0}^{\infty} a_n(4n^2-2n)x^{n-1}+\sum_{n=0}^{\infty} \sum_{k=0}^{n} \frac{(-1)^k}{(2n)!} a_{n-k} x^{n+k}=0

Ok, that's it. You're really not expected to reduce that further manually are you? That's just not reasonable to do in my opinion in the 21's century. It's just so simple with a CAS like Mathematica to get all the terms recursively: set a_0=1, extract all the coefficients of powers of x, solve for a_n:

nmax = 50; 
Remove[a, b]

myalist = Sum[Subscript[a, n]*(4*n^2 - 2*n)*x^(n - 1), {n, 0, nmax}] + Sum[((-1)^k/(2*k)!)*Subscript[a, n - k]*x^(n + k), 
     {n, 0, nmax}, {k, 0, n}]; 

Subscript[a, 0] = 1; 

myatable = Table[Subscript[a, n + 1] = First[Subscript[a, n + 1] /. N[Solve[Coefficient[myalist, x, n] == 0, 
         Subscript[a, n + 1]]]], {n, 0, nmax - 1}]; 

myatable = Prepend[myatable, 1]; 

f1[x_] := c1*Sum[Subscript[a, n]*x^n, {n, 0, nmax - 1}]

Once you do that for c=1/2 and find:

f2(x)=C_2 \sum_{n=0}^{\infty}b_n x^{n+1/2}

Then solve for c1 and c2 using the initial conditions:

they[x_] := f1[x] + f2[x]
thec = NSolve[{they[1] == 0, they'[1] == 1}, {c1, c2}] // First
pa = Plot[they[x] /. thec, {x, 1, 2}]

 

[outhsakotad]

Thanks... I'll look at it. According to my prof, we are supposed to get the recursion relation by hand. I'm just confused by the fractional indices I'm getting... :/

 

[jackmell]

Ok. I'm sorry. I shouldn't have told you to cheat. You need the practice learning how to do it manually. I mean really, I found a_1, a_2, a_3 by hand initially but got dizzy tryin' to find a_4. Tell you what, how about I just use Mathematica to find the recursive expressions for the first six (this is in Mathematica Latex so I don't have to type it in):
<br /> \left\{\left\{a_1\to -\frac{a_0}{2}\right\},\left\{a_2\to -\frac{a_1}{12}\right\},\left\{a_3\to \frac{1}{60} \left(a_0-2 a_2\right)\right\},\left\{a_4\to \frac{1}{112} \left(a_1-2 a_3\right)\right\},\left\{a_5\to \frac{-a_0+12 a_2-24 a_4}{2160}\right\},\left\{a_6\to \frac{-a_1+12 a_3-24 a_5}{3168}\right\}\right\}

See what I mean. I don't think that's happening for me.

Also, this is what I get for the fractional part.\text{myexp2}=\left(\sum _{n=0}^{\text{nmax}} b_n(4(n+c)(n+c-1)+2(n+c))x^{n+c-1}+\sum _{n=0}^{\text{nmax}} \sum _{k=0}^n \frac{(-1)^kx^{2k}}{(2k)!} b_{n-k}x^{n+c-k}\right)\text{/.}c\to 1/2;

and here's what the first 4 terms look like:

\sqrt{x} b_0-\frac{1}{2} x^{5/2} b_0+\frac{1}{24} x^{9/2} b_0-\frac{1}{720} x^{13/2} b_0+6 \sqrt{x} b_1+x^{3/2} b_1-\frac{1}{2} x^{7/2} b_1+\frac{1}{24} x^{11/2} b_1+20 x^{3/2} b_2+x^{5/2} b_2-\frac{1}{2} x^{9/2} b_2+42 x^{5/2} b_3+x^{7/2} b_3

. . . worst

Also keep in mind I compared this analytic solution with the numerical solution of the IVP above and got excellent agreement in the interval [1,5] using the first 50 terms and so that gives me high confidence the series expressions above are correct.

 

[outhsakotad]

Okay, it seems I'm having a conceptual problem here. When I try to find the recursion relation, I usually try to bring all the terms to have the same power of x, but I'm not sure how to deal with the "n+k" power that comes from the Cauchy Product...

 

[jackmell]

That's ok for regular easy problems. When you have a double sum, or something else complicated, need to just learn to use what you have. So don't try to bring everything to the same power. I don't think you can do that with this one but maybe I'm wrong. Also, I'm approaching it the way I was taught in Rainville and Bedient. Perhaps there is another way to approach it that would make it easier to determine the recursion relations.

 

[outhsakotad]

I honestly don't know how to do it another way... Could you perhaps give me a general idea of the method you are using to get the coefficients?

 

[outhsakotad]

Ahhh... I think I'm starting to see it.

 

[jackmell]

I honestly don't know how to do it another way... Could you perhaps give me a general idea of the method you are using to get the coefficients?

I used Mathematica:

\text{myatable}=\text{Table}\left[a_{n+1}=a_{n+1}\text{/.}\text{Solve}\left[\text{Coefficient}[\text{myaexp},x,n]==0,a_{n+1}\right]\text{//}\text{First},\{n,0,\text{nmax}-1\}\right]

That is, find the coefficient of x^n, set the resulting expression equal to zero, solve for a_(n+1).

Really, I believe in doing things by hand to understand them. That's important. But I've already done a bunch of these that way and understand the concept well enough to do that. Also I'm practical: if I can't figure it out, understand or no understantd I'm still going to use Mathematica to help me find the answer.

 

[outhsakotad]

I'm still having a lot of trouble understanding how to do it by hand. On page 1202 near the top of this document, there appears to be a general formula they derive for the coefficients. But I'm don't understand how to apply this formula. http://www.cacr.caltech.edu/~sean/applied_math.pdf And I don't see where the a_n+1 term comes from... From the Cauchy product, we have a_n's and a_n-k's.

 

[jackmell]

Ok, I looked at it. Suppose I write the DE then as:

y&#039;&#039;+\frac{1/2}{x}y&#039;+\frac{1}{x^2}\sum_{n=0}^{\infty} \frac{(-1)^n x^{2n+1}}{4(2n)!}y=0

or now:

y&#039;&#039;+\frac{P(x)}{x} y&#039;+\frac{Q(x)}{x^2}y=0

can we then now use that formula for z_n to compute the coefficients? Keep in mind that:

Q(x)=\sum_{n=0}^{\infty} \frac{(-1)^n x^{2n+1}}{4(2n)!}=\sum_{n=0}^{\infty}a_n x^n with a_n=0 for even n though.

 

[outhsakotad]

Ah, I get it! And I'm getting the same coefficients as you for the 0 indicial root. Thanks so much!

 

[outhsakotad]

Okay... so our s=0 coefficients match, but I'm getting y2=b0*x^(1/2)*{1-(1/6)x+(1/120)*x^2+(59/5040)*x^3+...) for the s=1/2 ones using the formula... that doesn't seem to be what you got...