Proof of conjecture concerning numerically satisfactory pairs of solutions to generalized Airy equation?

Question

Question: Given a generalized Airy equation of the form $$ \frac{d^2y}{dx^2}-x^ny=0 $$ I conjecture that the following functions form a pair of linearly independent, numerically satisfactory solutions for $n\in \{0,1,2,3,4,...\}$. How is this proven? $$ A_n(x)=\frac{_0F_1(;\frac{n+1}{n+2};\frac{x^{n+2}}{(n+2)^2})}{(n+2)^{\frac{n+1}{n+2}} \ \ \Gamma(\frac{n+1}{n+2})}-x\frac{_0F_1(;\frac{n+3}{n+2};\frac{x^{n+2}}{(n+2)^2})}{(n+2)^{\frac{n+3}{n+2}-1} \ \ \Gamma(\frac{n+3}{n+2}-1)} \\ $$ $$ B_n(x)=\frac{_0F_1(;\frac{n+1}{n+2};\frac{x^{n+2}}{(n+2)^2})}{(n+2)^{\frac{n+1}{n+2}-\frac{1}{2}} \ \ \Gamma(\frac{n+1}{n+2})}+x\frac{_0F_1(;\frac{n+3}{n+2};\frac{x^{n+2}}{(n+2)^2})}{(n+2)^{\frac{n+3}{n+2}-1-\frac{1}{2}} \ \ \Gamma(\frac{n+3}{n+2}-1)} $$

Context: First, I already know that these functions are solutions because the confluent hypergeometric limit functions within them are simply the series solutions obtained through a straightforward solution attempt with a power series method.

Second, I suspect that they are all linearly independent because these linear combinations appear to have a Wronskian given by $$ W\{A_n(x),B_n(x)\}=\frac{2\sin{(\frac{\pi}{n+2}})}{\pi\sqrt{n+2}} $$ I have verified this Wronskian up to $n=8$ in Mathematica, but I'm not sure how to prove it.

Furthermore, just for clarity, it's worth noting that these solutions reduce for the $n=0$ and $n=1$ cases to $$ A_0(x)=\frac{1}{\sqrt{2\pi}}e^{-x} \\ B_0(x)=\frac{1}{\sqrt{\pi}}e^x \\ A_1(x)=\text{Ai}(x) \\ B_1(x)=\text{Bi}(x) \\ $$ I have tested all of the pairs of solutions possible up to the limits of the numerical precision of Mathematica and found that $A_n(x)$ always converges to $0$ as $x\rightarrow\infty$ while $B_n(x)$ always diverges in the same limit. Likewise, for odd $n$ both solutions appear to oscillate exactly $\frac{\pi}{2}$ out of phase with one another as $x\rightarrow-\infty$. These are (basically) the conditions set forward by J. C. P. Miller in the linked article for numerical satisfaction.

However, I discovered these specific linear combinations through inspection. I have no idea why these constants are the ones that work. If you alter the constants even slightly, the solutions diverge wildly and no longer satisfy any of Miller's criteria.

How could I have deduced the form of these solutions more systematically? Is there a proof to demonstrate that these specific solutions are the preferred ones?

I know that there is a connection to the modified Bessel functions, but I don't really see it. There is clearly a hidden structure here that I'd like to understand and (hopefully) generalize to other, similar ODE families.

Answer 1

Starting from the Bessel ODE https://en.wikipedia.org/wiki/Bessel_function and changing variables $x \rightarrow f(x)$, $y^{'}(x) \rightarrow 1/f^{'}(x) y^{'}(x)$ where $f(x):=B x^n$ we end up with: \begin{eqnarray} y^{''}(x) + \frac{1}{x} y^{'}(x) + \frac{n^2}{x^2} (-\alpha^2 + B^2 x^{2 n}) y(x) = 0 \end{eqnarray} If we now eliminate the coefficient at the first derivative, i.e. write $y(x):=m(x) \cdot v(x)$ where $m(x):= \exp\left( -1/2 \int \frac{dx}{x} \right) = x^{-1/2}$ then we end up with: \begin{eqnarray} v^{''}(x) + \left( \frac{1/4-\alpha^2 n^2}{x^2} + B^2 n^2 x^{-2+2 n} \right) v(x) = 0 \end{eqnarray} The above ODE goes into your ODE on the top of your question if $n \rightarrow (n+2)/2$ and $\alpha=1/(n+2)$ and $B= I 2/(n+2)$. In summary all this means that the ODE: \begin{eqnarray} v^{''}(x) - x^n v(x)=0 \end{eqnarray} is solved by: \begin{eqnarray} v(x) := C_1 \cdot \sqrt{x} J_{\frac{1}{n+2}}(\imath \frac{2}{n+2} x^{n/2+1}) + C_2 \cdot \sqrt{x} Y_{\frac{1}{n+2}}(\imath \frac{2}{n+2} x^{n/2+1}) \end{eqnarray}

In[966]:= n =.; n1 =.; n2 =.; n3 =.; Clear[f]; Clear[y]; Clear[v]; \
Clear[m]; a =.; b =.; c =.; B =.; B1 =.; g =.; d =.; q =.; x =.; k \
=.; mu =.;
f[x_] := B x^n;
subst = {x :> f[x], 
   Derivative[1][y][x] :> 1/f'[x] Derivative[1][y][x], 
   Derivative[2][y][x] :> -f''[x]/(f'[x])^3 Derivative[1][y][x] + 
     1/(f'[x])^2 Derivative[2][y][x]};
myeqn = x^2 y''[x] + x y'[x] + (x^2 - a^2) y[x] /. subst /. 
   y[f[x]] :> y[x] ;


myeqn = Collect[
  myeqn/Coefficient[myeqn, y''[x]], {y[x], y'[x], y''[x]}, Simplify]

m[x_] = Simplify[Exp[-1/2 Integrate[Coefficient[myeqn, y'[x]], x]]]
y[x_] = m[x] v[x];
myeqn = Collect[myeqn, {v[x], v'[x], v''[x]}];
myeqn = Collect[
  myeqn/Coefficient[myeqn, v''[x]], {v[x], v'[x], v''[x], x^_}, 
  FullSimplify]
FullSimplify[(D[#, {x, 2}] - x^n #) & /@ {Sqrt[x] BesselJ[1/(n + 2), 
     I 2/(n + 2) x^((n + 2)/2)], 
   Sqrt[x] BesselY[1/(n + 2), I 2/(n + 2) x^((n + 2)/2)]}]


Out[970]= (n^2 (-a^2 + B^2 x^(2 n)) y[x])/x^2 + 
 Derivative[1][y][x]/x + (y^\[Prime]\[Prime])[x]

Out[971]= 1/Sqrt[x]

Out[974]= ((1/4 - a^2 n^2)/x^2 + B^2 n^2 x^(-2 + 2 n)) v[x] + (
  v^\[Prime]\[Prime])[x]

Out[975]= {0, 0}

Now there are connections between Bessel function and the hypergeometric functions as documented in here https://en.wikipedia.org/wiki/Bessel_function for example. From that you can see that your function are related to the Bessel functions I have given above but are not just equal (i.e. you function are some sort of linear combinations of the Bessel functions I gave above). As a matter of fact the Wronskian of my functions is equal to $(2+n)/\pi$ which is different from your Wronskian.

In[984]:= 
y = {Sqrt[x] BesselJ[1/(n + 2), I 2/(n + 2) x^((n + 2)/2)], 
   Sqrt[x] BesselY[1/(n + 2), I 2/(n + 2) x^((n + 2)/2)]};
Det[{y, D[y, x]}] // FullSimplify

Out[985]= (2 + n)/\[Pi]