Plotting $\omega(k_x, k_z)$ in $(\omega, k)$ plane with assumption of $k^2 = k_x^2 + k_z^2$

I'm currently studying atmospheric-gravity waves and their dispersion relations.

For analysis I need to plot the dispersion function:

$$\omega^4 - \omega^2\cdot C_s^2(k_x^2 + k_z^2 + \dfrac{1}{4H^2}) + (\gamma - 1)g^2\cdot k_x^2 = 0$$ where, $C_s^2$, $H$, $\gamma$ and $g$ are some constant values.

I used Solve to get $\omega$ as a function of $k_x^2$ and $k_z^2$ - $\omega = \omega(k_x^2, k_z^2)$.

Code for solving:

   Solve[\[Omega]^4 - \[Omega]^2*(Subscript[k, x]^2 +
Subscript[k, z]^2 + 1/(4*H^2)) + (\[Gamma] - 1)*g^2*Subscript[k, x]^2 == 0, \[Omega]]


What I want:

Plot[f, {k, 0, 20}]


where f is any solution of Solve. It plot nothing, because Mathematica didn't know about $k^2 = k_x^2 + k_z^2$ relation.

-
Could you please share your Mathematica code and specify what are the problems you get ? –  b.gatessucks Oct 29 '12 at 9:51
And how should I share code, when copying from Mathematica I get plenty of some strange additions like \FractioBox and \SuperscriptBox? –  m0nhawk Oct 29 '12 at 10:04
If you start with the equation I'll try to edit it so you can see how it's done. –  b.gatessucks Oct 29 '12 at 10:06
so you want say contours of iso kx in the (omega,k) plane which are solution to F(omega,k,kx)=0 above? –  chris Oct 29 '12 at 10:16
I've updated the question, if it isn't clear enough now - I will update further. –  m0nhawk Oct 29 '12 at 10:20

I'm pretty unsure about this, because I don't know much about manifolds and other difficult words. So I hope I learn something here too.

I use your equation with all parameters set to simple values as Chris did

eqn = w^4 - w^2 ((kx^2 + kz^2) + 1/4) + 2 kx^2 == 0;


Your second equation k^2==kx^2+kz^2 defines a tube with different radii k. We could for instance look at both equations and choose k==3

contPlot = ContourPlot3D[{w^4 - w^2 ((kx^2 + kz^2) + 1/4) + 2 kx^2 == 0,
kx^2 + kz^2 == 9}, {kx, -5, 5}, {w, -5, 5}, {kz, -5, 5},
ContourStyle -> {Directive[Opacity[0.3], Red], Automatic}]


We are interested in the intersection of both contours. What we can do here, is to parametrize k^2==kx^2+kz^2 differently. A tube along w with radius k can expressed as the parametric equation

$$f(w,\phi)=\{k \cos(\phi), w,k \sin(\phi)\}$$

therefore we define a transformation rule

rule = Thread[{kx, w, kz} :> {k*Cos[phi], w, k*Sin[phi]}]

(* {kx :> k Cos[phi], w :> w, kz :> k Sin[phi]} *)


We can now apply this transformation and solve your initial equation for w. With this we get solutions for w which only depend on phi and k. With those solutions, we have an explicit parametrization of the curve in 3d

sol = Solve[eqn /. rule, w];
paramPlot =
ParametricPlot3D[{kx, w, kz} /. rule /. sol /. k :> 3, {phi, 0,
2 Pi}, PlotStyle -> Red] /. Line[pts_] :> Tube[pts, 0.1]


And we can of course combine them to see whether it fits with our imagination

Show[{contPlot, paramPlot}]


What you want to have now is a plot where on the first axis is k and on the second w

ParametricPlot[Evaluate[{k, w} /. sol], {k, 0, 9}, {phi, -Pi, Pi}]


-
Its only a curve because I fixed arbitrarily one degree of freedom (I think anyway...) –  chris Oct 29 '12 at 15:45
and stop playing the trick of throwing lots of nice looking 3D curves: I invented it ;o) (just kidding nice plots and good answer!) –  chris Oct 29 '12 at 15:47
@chris Thanks. I see now that you set one dimension. I was just looking at your images and started trying it by myself. –  halirutan Oct 29 '12 at 15:53

eqn = ω^4 - ω^2 ((kx^2 + kz^2) + 1/4) + 2 kx^2


Let's now apply the extra condition on the modulus of $k$

eqn= eqn /. kz -> Sqrt[k^2 - kx^2];


then we can find zeros of the dispersion relation at fixed kx

 Table[ContourPlot[
eqn == 0 /. kx -> i/4 // Release, {ω, 0, 2}, {k, 0, 2},
ContourStyle -> ColorData[10][i],
RegionFunction -> Function[{k, ω}, k > i/4],
FrameLabel -> {k, ω}], {i, 1, 4}] // Show


which produces

where I have added a cut to impose k>kx as suggested by Rahul Narain.

Alternatively, following more closely his strategy (i.e. proceed at fixed angle between 'kx' and 'ky')

 eqn = ω^4 - ω^2 ((kx^2 + kz^2) + 1/4) + 2 kx^2 /.
kx -> k Cos[θ] /. kz -> k Sin[θ]
Table[ContourPlot[eqn == 0 /.θ-> Pi i/8 // Release,
{ω, 0, 2}, {k, 0,2}, ContourStyle -> ColorData[10][i],
FrameLabel -> {k, ω}], {i, 1, 8}] // Show


-

Using the same choice as @chris (and assuming you only want $\omega>0$, you can get the solution (in this case 2 solutions actually) as :

eqn[cs2_, oofh2_, gg_, kx_, kz_] = \[Omega]^4 - \[Omega]^2 cs2 ((kx^2 + kz^2) + oofh2) + gg kx^2 ;

sol[cs2_?NumericQ, oofh2_?NumericQ, gg_?NumericQ, kx_?NumericQ, kz_?NumericQ] :=
Solve[{eqn[cs2, oofh2, gg, kx, kz] == 0, \[Omega] >= 0}, {\[Omega]}, Reals][[All, 1, 2]]


Then you can simply plot it as a function of the parameters :

ParametricPlot[Thread[{Sqrt[kx^2 + kz^2], sol[1., 1/4., 2, kx, kz]}], {kx, -2, 2}, {kz, -2, 2},
AxesLabel -> {" k ", " \[Omega] "}]


However you get possibly better information with a more canonical plot :

Plot3D[sol[1., 1/4., 2, kx, kz], {kx, 0, 5}, {kz, 0, 5}, AxesLabel -> {" kx ", " kz ", " \[Omega] "}]


-