3
$\begingroup$

Consider the Sturm-Liouville problem

$$y''(x) + \lambda x^2 y=0, \ y(0)=0,\ y(1)=0$$

The analytical solution is given by

$$\lambda_n=4\alpha_n^2, \ y_n(x)=\sqrt{x}J_\frac14(\alpha_nx^2)$$

where $\alpha_n$is the $n$-th zero of the Bessel function of the first kind of order 1/4. Thus, $\lambda_1\approx30.93$ and $y_1(x) \approx \sqrt{x} \,J_\frac14(2.78x^2)$

{vals, funs} = 
  NDEigensystem[
    {-(y''[x]/x^2), 
     DirichletCondition[y[x] == 0, x == 0], 
     DirichletCondition[y[x] == 0, x == 1]}, 
    y[x], {x, 0, 1}, 1]

yields vals = 16.1035 which is clearly incorrect.

$\endgroup$

2 Answers 2

5
$\begingroup$

That is an interesting case I have not seen before. You can work around this by using a time dependent formulation like so:

{vals, funs} = 
 NDEigensystem[{x^2*D[y[t, x], {t, 1}] + Laplacian[y[t, x], {x}] == 0,
    DirichletCondition[y[t, x] == 0, x == 0 || x == 1]}, {y[t, x]}, 
  t, {x, 0, 1}, 1];

vals
{30.933443250231416`}

I suspect that the original formulation has numerical issues close to x=0. Might be worth putting in the possible issues section, what do you think?

$\endgroup$
2
$\begingroup$

This is not an answer, just long comment. This verifies that eigenvalues from NDEigensystem are not correct.

ClearAll[y,x,lam];
sol=y[x]/.First@DSolve[y''[x]+lam*x^2*y[x]==0,y[x],x];
eq1=(sol/.x->0)==0;
eq2=(sol/.x->1)==0;
c=C[1]/.First@Solve[eq1,C[1]]; (*Solve for C[1]*)
eq2=eq2/.C[1]->c (*plug in second equation*)

Mathematica graphics

eq2=eq2/.C[2]->1 (*divide by C[2]*)

Mathematica graphics

Plot[eq2,{lam,0,200}]

Mathematica graphics

There are two zeros, one near 30 and one near 140

FindRoot[eq2,{lam,30}]//Chop

Mathematica graphics

4*(N@BesselJZero[1/4,1])^2

Mathematica graphics

FindRoot[eq2,{lam,140}]//Chop

Mathematica graphics

4*(N@BesselJZero[1/4,2])^2

Mathematica graphics

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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