Let $P = a_n(x) D_x^n + a_{n-1}(x) D_x^{n-1} + \ldots + a_0(x)$ be a linear ordinary differential operator with polynomial (or real analytic) coefficients $a_j(x)$. Suppose that $a_n(x)$ doesn't vanish on an interval $(a,b)$ and that $u$ is a weak solution of $P u = 0$ on $(a,b)$. It can be concluded from the elliptic regularity theorem that $u$ is in fact real analytic and a classical solution to the differential equation.

In other words, the operator $P$ is analytically hypoelliptic.

However, appealing to the elliptic regularity theorem seems a bit much for the above case of an ordinary differential operator. Does anyone know a (canonical) reference for this (probably much more classical) case? I browsed through several books on differential equations but the closest I could find was in G. Folland's Fourier analysis and its applications where he mentions this fact (with just $C^\infty$ smoothness) on the top of page 344 without reference or proof.

Thank you!