# Homework Help: Solving for complex phase speed eigenvalue

1. Aug 9, 2011

### primalfido

Hey there, first timer poster here

1. The problem statement, all variables and given/known data
I'm working on a barotropic linear instability analysis and I've been having trouble getting an expression for the complex phase speed eigenvalue $C = C_r + i*C_i$ for the purpose of plotting a dispersion diagram ($C_i$ vs $k$ or $C_r$ vs $k$)

2. Relevant equations

The Normal mode solution:

$\psi(x,y,t) = \Psi(y)e^{ik(x - Ct)}$

where k is the zonal wavenumber

And the eigenvalue problem I'm trying to solve for C

$\frac{\partial^{2} \Psi}{\partial y^{2}} \Psi + (\frac{\beta - \frac{\partial^{2} U}{\partial y^2}}{U-C} -k^{2} )\Psi = 0$

where U is a velocity profile thats given in gaussian form \e^{-y^2} and \beta is a given constant.

3. The attempt at a solution

The only things undetermined are C,k and y. The purpose is to find C for a given k and y. However I'm not too sure on the procedure in solving an eigenvalue expression like this so I tried to follow the Chebyshev method in Mathematica to solve for the eigenvalue C and I've ended up with a large output at the
Eigenvalues[result] part of the code.
that is too large to post full so here is a small part of it showing values and the general form:
$Root[......] \ , \ Root[-600. C^5 E^((3 y^2)/2000000000000) k^12 + 150000. C^4 E^((3 y^2)/2500000000000) k^12 + ......] \ ,Root[....]$

But I'm not sure whether I'm on the right track and if I am, how to go about plotting C vs k.

Thanks

Sorry if I'm not sounding clear

Last edited: Aug 9, 2011
2. Aug 17, 2011

### primalfido

Just thought I'd post the Mathematica code I'm working with:

Code (Text):

Uo = 100
b = 6*10^(-11)
R = 60268000
r = R*Sin[76]
U = Uo*Exp[-b*(y^2)/(2*Uo)]
Ld = N[Sqrt[2*Uo/b]]
\[CapitalOmega] = N[(2*\[Pi]/(39 + 60*39 + 10*60*60))]
\[Beta] = 2*\[CapitalOmega]*Cos[(76) Degree]/R
Clear[y]

Uyy2 =  (b^{2} y^{2}/Uo)*Exp[-b*(y^2)/(2*Uo)] - b*Exp[-b*(y^2)/(2*Uo)]

DE = \[Phi]''[
y] + (- k^2 +  (\[Beta] - D[U, {y, 2}] )/(U - g)) \[Phi][y]

t1 = DE /. {\[Phi][y] -> a[i]*ChebyshevT[i, y],
D[\[Phi][y], {y, a1_}] -> D[a[i]*ChebyshevT[i, y], {y, a1}]}

t2 = ExpandAll[t1]

n = 12

t3 = Sum[t2, {i, 0, n}]
t4 = Simplify[t3]
t5 = (t4 *(1/Sqrt[1 - y^2])*ChebyshevT[i, y] 2)/(\[Pi] c[i])
t6 = Table[t5, {i, 0, n, 2}];
t7 = (t6 /. c[0] -> 2) /. Table[c[i] -> 1, {i, 2, n}];
t8 = Collect[Expand[Numerator[t7]], y];
t9 = Transpose[Table[Coefficient[t8, y, i], {i, 0, 2 n, 2}]]
t10 = Table[y^i/(\[Pi] Sqrt[1 - y^2]), {i, 0, 2 n, 2}]
Timing[t11 = \!$$\*SubsuperscriptBox[\(\[Integral]$$, $$-1$$, $$1$$]$$t10 \ \[DifferentialD]y$$\)]
t12 = Simplify[Expand[t9.t11]]
t13 = Transpose[Table[Coefficient[t12, a[i]], {i, 0, n, 2}]]
b1 = Table[1, {i, 0, n, 2}]
t14 = Append[Table[t13[[i]], {i, 1, n/2}], b1]
s1 = Solve[t14[[n/2 + 1]].Table[a[i], {i, 0, n, 2}] == 0, a[n]]
t15 = Table[t14[[i]].Table[a[i], {i, 0, n, 2}], {i, 1, n/2}]
t16 = Expand[t15 /. s1]
t17 = Transpose[
Table[Table[Coefficient[t16[[1, i]], a[j]], {i, 1, n/2}], {j, 0,
n - 2, 2}]]
(*How would I go about geting the complex phase speed g from this?*)
(*Like this? Am I approaching it wrong*)
vals = Eigenvalues[-t17 /. g -> 0]

I'm just not very sure on how to go about solving this complex eigenvalue problem numerically.