# How to solve Sturm-Liouville problems in SageMath?

I've conducted DuckDuckGo searches for "Sturm-Liouville sagemath" and received no usable results, and as a result I am here to ask, how can one solve Sturm-Liouville problems with SageMath? I know how to do so numerically with Julia, GNU Octave and Scilab, using matrices, but is there any way to solve at least some of them analytically in SageMath (e.g. $\frac{d^2 y}{dx^2}-xy = \lambda y$ has solutions $y_n=C\cdot \textit{Ai}(x+\lambda_n)$, where $\lambda_n $ are the zeros of the Airy $\textit{Ai}(x)$ function, if $y(0)=y(\infty)=0$)? Failing this, is there a way to numerically solve them in SageMath (without resorting to its interfaces with the languages I just mentioned)?

If it helps, here's the code I use in GNU Octave to solve the aforementioned SLE:

```
# Standard values are for transformation x=(y-L)/(y+L), y in 0, inf
# *lin are for x=2/L*y-1 (linear transformation)
clear all
format long e;
N=5000;
perc=0.266
Nfrag=round(perc*N);
n=0:N;
# 75 at N=1k; 162 for N=3k; 225 for N=5k;
L=226;
t=pi+n'*pi/N;
x=cos(t);
y=L*(cot(t/2)).^2;
ysub=y(2:N);
ylintrans = L/2*(x+1);
# Chebyshev polys of the first kind
T = cos(t*n);
# Now for arrays that do not include the endpoints
xsub = x(2:N);
Tsub = T(2:N,:);
Usub = diag(1./sqrt(1-xsub.^2))*sin(acos(xsub)*n);
dTsub = Usub*diag(n);
dT = [-(n.^2).*(-1).^(n); dTsub ; (n).^2];
dTL = diag(2*L./((y+L).^2))*dT;
dTlin = 2/L*dT;
D1 = dTL/T;
D1lin = dTlin/T;
D2 = D1^2;
D2lin = D1lin^2;
H = D2-diag(y);
Hlin = D2lin-diag(ylintrans);
H = H(2:N,2:N);
Hlin = Hlin(2:N,2:N);
[Y, Lam] = eig(H);
[Ylin, Lamlin] = eig(Hlin);
Lam = diag(Lam);
Lamlin = diag(Lamlin);
[Lam, IX] = sort(Lam, 'ascend');
[Lamlin, IXlin] = sort(Lamlin, 'ascend');
Y = Y(:,IX);
Ylin = Ylin(:,IXlin);
aizerochk = abs(airy(0,Lam));
aizerochklin = abs(airy(0,Lamlin));
L
errrms = sqrt(aizerochk(1:Nfrag)'*aizerochk(1:Nfrag)/(Nfrag+1))
errrmslin = sqrt(aizerochklin(1:Nfrag)'*aizerochklin(1:Nfrag)/(Nfrag+1))
ain=quad("airynorm",0,100);
figure(1)
plot(ysub(1:3*N/5), Y(1:3*N/5,1));
figure(2)
plot(ysub(1:3*N/5), Y(1:3*N/5,2));
figure(3)
plot(ysub(1:3*N/5), Y(1:3*N/5,Nfrag));
```

. It uses a Chebyshev spectral method on the extrema grid, and then analyses the results.

I think that anything which can be computed in Julia, Octave or Scilab can also be computed in SageMath. So I would suggest to post the code you use to solve the problem in one of these languages, e.g. in Octave. In this way, you could get help to translate it to SageMath.

Thanks, shall do.

From

`airy_ai?`

:But Sage doesn't use it when solving the ODE...