7
$\begingroup$

I am trying to solve the eigensystem of a 1st-order linear ODE system in the region $(-\infty,\infty)$ and with Dirichlet boundary condition at the infinities $$ -i\partial_xu(x)+f^*(x)v(x)=\lambda u(x)\\ f(x)u(x)+i\partial_xv(x)=\lambda v(x) $$ where $f(x)$ is a complex-valued function mainly varying around $x=0$ and $f^*(x)$ is its complex conjugate. I define $f(x)$ by its norm and phase in the code, where I give two forms of $f(x)$ (using Tanh or Lorentzian) of the similar profile.
I tried reducing it to a 2nd-order ODE. But unfortunately, it messes up the eigenstructure, i.e., $\lambda$ appears in several terms. I have no idea if it can be solved analytically.

It is known that the system will have a few (at least one) discrete real eigenvalues in $(-1,1)$ and the eigenfunction is more or less localized around $x=0$. Outside $(-1,1)$, there will be a continuous spectrum with probably not quite confined eigenfunctions. I am interested in the eigenvalue(s) in $(-1,1)$ only.

The code below gives 4 solutions but some eigenfunctions are very saw-tooth like. Therefore, I want to make sure if these solutions are true or spurious and if any other way to solve it reliably.

{0.0856771, -0.169612, -0.881437, 0.94686}

a = 1; b = 0.8 a; c = 1; Nless = 8; cutoff = 20;
xR = cutoff; xL = -cutoff;
f[x_, pm_] = (a + b (Tanh[x/c - 1] - Tanh[x/c + 1])/(2 Tanh[1])) Exp[
    I pm \[Pi] (Tanh[x/c] + 1)/2];
(*f[x_, pm_] = (a + b (-c^2)/(x^2 + c^2)) Exp[
    I pm \[Pi] (Tanh[x/c] + 1)/2];*)
Fop1[F_] := D[F, x];
variables = {u, v};
lhs = {-I Fop1[u[x]] + f[x, 1] v[x], I Fop1[v[x]] + f[x, -1] u[x]};
bc = DirichletCondition[
   Table[component@x == 0, {component, variables}], True];
{vals, funcs} = 
  NDEigensystem[{lhs, bc}, {u[x], v[x]}, {x, xL, xR}, Nless, 
   Method -> {"PDEDiscretization" -> {"FiniteElement", {"MeshOptions" \
-> {"MaxCellMeasure" -> 0.001, "MeshOrder" -> 2}}}}];
vals
Table[Plot[
  If[iband <= Nless, Norm, _]@funcs[[iband, ab]], {x, xL, xR}, 
  PlotRange -> All, ImageSize -> Small], {ab, {1, 2}}, {iband, 1, 
  6}(*,{part,{Norm(*,Re,Im*)}}*)]
$\endgroup$
  • $\begingroup$ The short answer is yes, i think the 0.0856771 value is spurious and the other three are fine. I'll hopefully post an answer using my package later. $\endgroup$ – KraZug Jul 12 at 16:46
  • $\begingroup$ I don't really know why NDEigensystem gives these spurious roots. The Evans function has a double root at zero here, in a similar way to your previous question, that could be part of the reason. $\endgroup$ – KraZug Jul 12 at 19:25
  • $\begingroup$ @KraZug I notice a weird behaviour. Your package doesn’t give any solution if I change $\pi$ even by a tiny bit in the phase factor of $f$ while NDEigen still produces some solutions. Any remedy or reason? Thanks. $\endgroup$ – xiaohuamao Jul 12 at 23:07
  • $\begingroup$ The Evans function becomes complex for real values of $\lambda$, so Plot won't show anything unless you use ReIm or Abs. Or use ContourPlot or just FindRoot. $\endgroup$ – KraZug Jul 13 at 7:23
  • $\begingroup$ @KraZug I just use FindRoot and it gives nothing when I set $\pi$ to $1.1\pi$ in the phase factor. Does your answer give any solution for this case? $\endgroup$ – xiaohuamao Jul 13 at 7:33
7
$\begingroup$

I have a package for solving eigenvalue boundary value problems using the Compound Matrix Method with the Evans function, which I'll use here. The package is available on my GitHub (which has a notebook with examples), more details are in my other answers to questions here, and a good introduction to the method is in this pdf.

Install and load the package:

PacletInstall["CompoundMatrixMethod", 
     "Site" -> "http://raw.githubusercontent.com/paclets/Repository/master"];
Needs["CompoundMatrixMethod`"]

First write the system in the form $\mathbf{y}' = \mathbf{A} \cdot \mathbf{y}$ using the function ToMatrixSystem:

sys1 = ToMatrixSystem[Thread[lhs == λ {u[x], v[x]}], {}, {u, v}, {x, -L, L}, λ]

Note by not giving boundary conditions the requirement is that $\mathbf{A}$ tends to the decaying behaviour of $\mathbf{A}_{\pm \infty} =\lim_{x\to\pm\infty} \mathbf{A}(x)$. This will only work when there is a negative (resp. positive) eigenvalue of $\mathbf{A}_{\infty} (\mathbf{A}_{- \infty})$. Here the eigenvalues become imaginary for real (not complex) values outside (-1,1):

Plot[Evaluate[Eigenvalues[Limit[sys1[[1]], x -> ∞]]], {\[FormalLambda], -1, 1}]

enter image description here

and Evans will fail without explicit boundary conditions at those points. This method may be wrong to apply here, I think I need to consider further whether that is only applicable for higher order differential equations (e.g. a second order ODE in one variable).

The Evans function is an analytic function whose zeroes correspond to eigenvalues, so we can plot that and look for roots:

Plot[Evans[λ, sys1 /. L -> 10], {λ, -0.99, 0.99}, 
 Epilog -> Point[{#, 0} & /@ vals]]

enter image description here

It looks initially like all 4 of the NDEigensystem values are roots, but if we zoom in we can see that there is a double root at zero, so the value at $0.0856771$ is spurious:

enter image description here

The winding number of the circle with radius 0.99 is 5, showing there are no further eigenvalues with $|\lambda|<=0.99$. The double root at zero is contributing two to the count here.

PlotEvansCircle[sys1 /. L -> 10, ContourRadius -> 0.99, nPoints -> 1000, Joined -> True]

enter image description here

Note that providing explicit boundary conditions in gives very similar results, except for the eigenvalue at zero.

$\endgroup$
  • $\begingroup$ Sorry, but it seems that your links point to wrong content. Could you please check them again? $\endgroup$ – Pinti Jul 15 at 6:59
  • $\begingroup$ @Pinti Ah sorry, I copied those from a previous post but forgot to include the links. Duh. Edited them. $\endgroup$ – KraZug Jul 15 at 8:59
  • $\begingroup$ @Pinti, if you have any questions regarding my package I'm more than happy to discuss them. $\endgroup$ – KraZug Jul 15 at 9:56

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.