On a novel iterative method to compute polynomial approximations to Bessel functions of the first kind and its connection to the solution of fractional diffusion/diffusion-wave problems

We present an iterative method to obtain approximations to Bessel functions of the first kind $J_p(x)$ ($p>-1$) via the repeated application of an integral operator to an initial seed function $f_0(x)$. The class of seed functions $f_0(x)$ leading to sets of increasingly accurate approximations $f_n(x)$ is considerably large and includes any polynomial. When the operator is applied once to a polynomial of degree $s$, it yields a polynomial of degree $s+2$, and so the iteration of this operator generates sets of increasingly better polynomial approximations of increasing degree. We focus on the set of polynomial approximations generated from the seed function $f_0(x)=1$. This set of polynomials is not only useful for the computation of $J_p(x)$, but also from a physical point of view, as it describes the long-time decay modes of certain fractional diffusion and diffusion-wave problems.