Because infinite series are expensive for numerical simulations, I am trying to approximate the following time and radius dependent function:
$$ f(t,r)={\sum_{n=1}^\infty}{{J_0(r\cdot a_n)}\over {{(a_n)^3e^{(a_n)^2t}J_1(a_n)}}}\tag1 $$
where:
$ t$ is the dimensionless time
$ r$ is the dimensionless radius
$ J_0$ and $J_1$ are Bessel functions of the first kind of order zero and one.
$ a_n$ is $n$-th zero of Bessel function of the first kind of order zero.
The boundaries are $t\in[0,t_{\text{finite}}] $ and $r\in[0,1] $.
I noticed a linear characteristic by plotting the function $a_n$. It appears that a good analytic expression of $a_n$ is $$ a_n = \pi n-{\pi\over4}\tag2 $$
By using the above expression, I plotted equation $(1)$ for $n_{\max} =500$. A $2$D and $3$D plot is shown below.
Please tell me with what numerically cheap function can I approximate $f(t,r)$. The summation is probably too hard to solve considering it has two Bessel functions with different arguments.
Thank you for your time.
Quite funny, I was running into exactly the same problem the other day.
I am assuming that the source of this sum is from Batchelor's solution to the time-varying flow of laminar pressure driven flow in a pipe. In dimensionless form, it can be expressed as
$$u(t,r)=2\left(1-r^2\right)-\sum_{n=1}^\infty \frac{16}{{s_n}^3}\frac{1}{J_1(s_n)}J_0\left(s_nr\right)\exp\left(-{s_n}^2~t\right)$$ The symbols used have the meanings
$u=\frac{8\nu}{ga^2}u_{\text{real}}$ is the dimensionless velocity where $\nu$ is the kinematic viscosity, $g$ is the kinematic pressure gradient, and $a$ is the radius of the pipe
$J_\alpha$ is a Bessel function of the first kind
$s_n$ is the $n$th positive real root of $J_0$
$t=\frac{\nu t_{\text{real}}}{a^2}\in[0,\infty)$ is the dimensionless time
$r=\frac{r_{\text{real}}}{a}\in[0,1]$ is the dimensionless radius
In fact asymptotic expressions for $s_n$ are known. In particular it is known that $$s_n=-\frac{\pi}{4}+n\pi+\frac{1}{8n\pi-\pi/4}+\mathrm O(n^{-3})\\ \text{as}~n\to\infty$$
We can combine this with known asymptotic formulas for Bessel functions (given $x\in \Bbb R$):
$$J_1(x)=\sqrt{\frac{2}{\pi x}}\cos\left(x-\frac{3\pi}{4}\right)+\mathrm O\left(|x|^{-3/2}\right) \\ \text{as}~|x|\to\infty$$
Which means,
$$J_1(s_n)\asymp \sqrt{\frac{2}{\pi}}~\frac{4\pi(32n-1)}{16+\pi^2-36\pi^2 n+128\pi^2n^2}\cos\left(\frac{4+\pi^2-33\pi^2 n+32\pi^2n^2}{\pi(32n-1)}\right) \\ \text{as}~n\to\infty$$
A crude, but nonetheless very good approximation of this is
$$J_1(s_n)\approx 2^{1/2}\pi^{-3/2}\frac{1}{9/32-n}\cos(\pi n)\sim -2^{1/2}\pi^{-3/2}\frac{1}{n}$$
Also, from our asymptotic expansion for $s_n$, we can deduce $$\frac{1}{{s_n}^3}\approx\frac{1}{\pi^3}\frac{1}{\big(n-1/4\big)^3}\sim\frac{1}{\pi^3n^3}$$
Putting these together, this means that $$\boxed{\left|\frac{1}{{s_n}^3}\frac{1}{J_1(s_n)}\right|\sim2^{-1/2}\pi^{-3/2}n^{-2}} $$ (check for yourself)
Finally, we approximate the exponential bit. By the same reasoning as before, we have $${s_n}^2\approx \pi^2 \big(n-1/4\big)^2$$ So $$\exp(-{s_n}^2 t)\sim \frac{1}{e^{\pi^2n^2 t}}$$
So in all we have the very crude but nonetheless enlightening approximation
$$\boxed{\left|\frac{16}{{s_n}^3}\frac{1}{J_1(s_n)}\exp(-{s_n}^2 t)\right|\sim 16\cdot 2^{-1/2}\pi^{-3/2}\frac{1}{n^2e^{\pi^2n^2 t}}}$$
Note that this expression shrinks w.r.t $n$ very fast, and does so even faster as $t$ increases.
Keeping in mind that the maximum value that $u$ takes in the channel is $2$, once the terms in the sum become smaller than, say $\frac{2}{1000}$, they are probably ok to neglect. As far as I can tell from my numerical investigation, for all $t>0.01$, this happens in the worst case after $6$ terms. So, all together, my suggestion to you is to use the expression
$$\boxed{\small{u(t,r)\approx 2(1-r^2)-\sum_{n=1}^{\left\lceil\frac{6}{t+1}\right\rceil}\frac{16}{\sqrt{2\pi^3}}\frac{\frac{9}{32}-n}{\big(n-1/4\big)^3}~\sec(\pi n)~J_0\left(\left(\frac{-\pi}{4}+n\pi\right)r\right)\exp\left(-\pi^2\big(n-1/4\big)^2 t\right)}}$$
And hopefully you should have pretty good results :)
EDIT:
As an example, the exact expression for $u(t,r)$ gives the velocity after half a second at the center of the pipe as $u(0.5,0)=1.82032$ whereas the approximate expression gives $u(0.5,0)\approx 1.78435$. Not bad. The Mathematica code is shown below: