4
$\begingroup$

Greeting, I'm searching to solve numerically the dispersion relation for Landau Damping, it is: $$ 0=1+\frac{1}{k^2}\int_{-\infty}^{+\infty} \frac{\partial_v f_{0}}{w/k-v}dv$$ Which then can be simplified by using the Principal value.

I know that it can be solved if $f_0$ is a Gaussian and using the plasma dispersion relation $Z(\zeta)$. However i want to do it with an arbitrary function $f_0$.

I know Mathematica can find the roots of equation by using FindRoot[], and i tested it with a dispersion relation with out the integral and it works fine.

So, my problem/question is how can i solve that type of integral?

I tried just to test it with a simple Gaussian Gaussian := Exp[-(((v - V0)/VT)^2)/2]/(Sqrt[2*Pi]*VT)

NIntegrate[D[Gaussian, v]/(w/k - v), {v, -Infinity, Infinity}, PrincipalValue -> True]

and also

NIntegrate[D[Gaussian, v]/(w/k - v), {v, -Infinity, Infinity}]

but the first one says : NIntegrate: Singular points must be specified in the integration range in order to use PrincipalValue.,

and the second one says: NIntegrate: The integrand ... has evaluated to non-numerical values for all sampling points in the region with boundaries {{1.,0.}}.

How can i solve this type of integral and equation? Is there a better way than doing this?

Edit: Second Question: My final objective is to solve the dispersion relation, i found i can use FindRoot[Dispersion relation], but how can i do that? Is it possible numerically?

Thank you.

$\endgroup$
4
$\begingroup$

For NIntegrate all parameters must have numerical values.

Gaussian = Exp[-(((v - V0)/VT)^2)/2]/(Sqrt[2*Pi]*VT);

integrand[v_, V0_, VT_, k_, w_] = D[Gaussian, v]/(w/k - v) //   Simplify

(*   (E^(-((v - V0)^2/(2 VT^2))) k (v - V0))/(Sqrt[2 \[Pi]] VT^3 (k v - w))   *)


nint[V0_, VT_, k_, w_] := 
  NIntegrate[integrand[v, V0, VT, k, w], {v, -Infinity, Infinity}, 
    Method -> "PrincipalValue", Exclusions -> (k v - w) == 0]

nint[1, 2, 3, 4]

(*   0.24312   *)
| improve this answer | |
$\endgroup$
  • $\begingroup$ Oh shit, thanks Second Question: My final objective is to solve the dispersion relation, i found i can use FindRoot[Dispersion relation], but how can i do that?. I will edit the question to add this. Thanks $\endgroup$ – Gundro Mar 30 at 9:57
3
$\begingroup$

Have a look at Normaldistribution and central moment. Than take a look at statistical dispersion. Return to Mathematica documentation: descriptive statistics. They more often towards the reverse situation of given measurement values. This question is targeted towards a fundamental situation for example in plasmas with plane waves. There is not much interest in solving the equation. There is a historical discussion between Landau and Pavlov about the physicality of this dispersion. Pavlov did take a realistics position, Landau an antirealistic one with a preference for pedagogical reason.

There has been a general question group regarding this question. How to deal with complicated Gaussian integrals in Mathematica?.

Form there I took by chance a function for the calculation of this dispersion:

gaussMoment[fPre_, fExp_, vars_] := 
 Module[{coeff, dist, ai, \[Mu], norm}, 
  coeff = CoefficientArrays[fExp, vars, "Symmetric" -> True];
  ai = Inverse[2 coeff[[3]]];
  \[Mu] = -ai.coeff[[2]];
  dist = MultinormalDistribution[\[Mu], -ai];
  norm = 1/PDF[dist, vars] /. Thread[vars -> \[Mu]];
  Simplify[
   norm Exp[1/2 coeff[[2]].\[Mu] + coeff[[1]]] Distribute@
     Expectation[fPre, vars \[Distributed] dist]]]

The moments from there as example for consistency:

RepeatedTiming[gaussMoment[(x^2 + x^4 + x^6), -(x - 1)^2, {x}]]

{0.0028, (223 Sqrt[\[Pi]])/8}

For this question:

RepeatedTiming[gaussMoment[((x - q)^(-1)), -(x - 1)^2, {x}]]

{10.2, Sqrt[\[Pi]]
   Expectation[
   1/(-q + x), {x} \[Distributed] 
    MultinormalDistribution[{1}, {{1/2}}]]}

This is very fast and safe. I think they probably took their formulas from Matlab. In Matlab, more application is in concern for such integrals. As stated earlier in other question this integration can be done nicely numerical for a general parametric solution.

With Your parameters in my suggested integration method:

D[Exp[-(((v - V0)/VT)^2)/2]/(Sqrt[2*Pi]*VT), v]

-((E^(-((v - V0)^2/(2 VT^2))) (v - V0))/(Sqrt[2 \[Pi]] VT^3))

RepeatedTiming[
  gaussMoment[(((v - V0)/(w/k - v))), -(((v - V0)/VT)^2)/
    2, {x}]]/(Sqrt[2 \[Pi]] VT^3)

The result is

{0.00054/VT^3, (E^(
    1/2 (-((v - V0)^2/
        VT^2) + {-((v - V0)^2/(2 VT^2))}[[
         2]].(-Inverse[
            2 {-((v - V0)^2/(2 VT^2))}[[3]]].{-((v - V0)^2/(
              2 VT^2))}[[2]])))
     Expectation[(k (-v + V0))/(
     k v - w), {x} \[Distributed] 
      MultinormalDistribution[-Inverse[
          2 {-((v - V0)^2/(2 VT^2))}[[3]]].{-((v - V0)^2/(2 VT^2))}[[
          2]], -Inverse[2 {-((v - V0)^2/(2 VT^2))}[[3]]]]])/(Sqrt[
    2 \[Pi]] VT^3 PDF[
     MultinormalDistribution[-Inverse[
         2 {-((v - V0)^2/(2 VT^2))}[[3]]].{-((v - V0)^2/(2 VT^2))}[[
         2]], -Inverse[
        2 {-((v - V0)^2/(2 VT^2))}[[3]]]], {-Inverse[
         2 {-((v - V0)^2/(2 VT^2))}[[3]]].{-((v - V0)^2/(2 VT^2))}[[
         2]]}])}

enter image description here

| improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.