## 3 Discretization in Time on a Non-Uniform Mesh  
The main goal of this work was to construct an unconditionally-stable implicit finite difference method defined on a non-uniform mesh. In general, the existence of a weakly singular kernel \((t-s)^{-\alpha}\), \(0<\alpha<1\), in fractional derivatives makes it more difficult to obtain a higher-order scheme. In particular, when the solutions are not suitably smooth, numerical methods on uniform meshes seem to have a poor convergent rate. For these reasons, numerical schemes on non-uniform meshes have been developed in the last years. In this section, first, we construct the non-uniform mesh and then, in order to approximate the Caputo derivative on the non-uniform grid, we define a suitable discrete fractional derivative formula. We divide the interval, \([0,T]\), into \(N\) subintervals, \([t^{n-1},t^{n}]\), for \(n=1,\cdots,N\) and with \(0=t^{0}<t^{1}\)\(<\cdots t^{N}=T\). We denote the time step \(\Delta t^{n}\) as  
\[\Delta t^{n}=t^{n}-t^{n-1},\qquad 1\leq n\leq N,\]  
and let  
\[\Delta t^{max}=\max_{1\leq n\leq N}\Delta t^{n}\qquad\Delta t^{min}=\min_{1 \leq n\leq N}\Delta t^{n}.\]  
A sequence of mesh is described as quasi-uniform if there a finite constant, \(C\), exists such that  
\[\Delta t^{max}/\Delta t^{min}\leq C.\]  
In this case, it holds that \(\Delta t^{max}\leq\beta TN^{-1}\). When \(C=1\), it holds that \(\Delta t^{n}=TN^{-1}\) for all \(n=1,\cdots,N\), and the mesh is reduced to a uniform mesh. A sequence of meshes is not quasi-uniform if  
\[\Delta t^{max}/\Delta t^{min}\rightarrow+\infty\qquad\text{as}\qquad N \rightarrow+\infty.\]  
The non-uniform mesh and quasi-uniform mesh methods have been used for solving differential equations by several authors. In [39; 40; 41; 42], quasi-uniform meshes were used for solving a class of boundary values problems on infinite domains. In this work, we are interested on the non-uniform mesh [28], which is defined as follows  
\[\Delta t^{n}=(N+1-n)\sigma,\qquad 1\leq n\leq N, \tag{3}\]  
where \(\sigma=\frac{2T}{N(N+1)}\). Note that the time steps, \(\{\Delta t^{n}\}_{n=1}^{N}\), are a monotonically decreasing sequence, with \(\Delta t^{1}=O(N^{-1})\) and \(\Delta t^{N}=O(N^{-2})\). Figure 1 shows a sample of the non-uniform grid (3) obtained for \(N=10\). The fractional derivatives are integrals with weakly singular kernels, and it is well known that the discretization on the uniform mesh may lead to poor accuracy. We decided to use a non-uniform mesh in order to improve the accuracy of the numerical solution. Another interesting choice could be to construct non-uniform meshes that refine the computational domain according to the fractional order of \(\alpha\) or to use adaptive procedures to dynamically choose the size of the time-steps according to the local behavior of the solution. To motivate the choice of the non-uniform mesh (3), we discretized the Caputo derivative (2) of a function \(v(t)\) by means of the following classical approximate formula  
\[\frac{d^{\alpha}}{dt^{\alpha}}v(t^{n}) = \frac{1}{\Gamma(1-\alpha)}\int_{0}^{t^{n}}v^{\prime}(s)(t^{n}-s)^ {-\alpha}ds\] \[= \frac{1}{\Gamma(1-\alpha)}\sum_{k=1}^{n}\frac{v(t^{k})-v(t^{k-1}) }{t^{k}-t^{k-1}}\int_{t^{k-1}}^{t^{k}}(t^{n}-s)^{-\alpha}ds+R^{n},\]  
which is the well-known, so-called L1 formula defined in [2] where \(R^{n}\) is the local truncation error. The L1 formula has been used to solve the fractional differential equations with Caputo derivatives (see [22; 43]). Moreover, using the relationship between the Caputo derivative and the Riemann-Liouville fractional derivative, the L1 formula was also applied to the time fractional diffusion equation with the Riemann-Liouville fractional derivative (see [44; 45]). High-order approximations, such as the compact difference scheme [43; 45; 46; 47] and spectral method [21; 23; 48] were applied to improve the spatial accuracy of fractional diffusion equations.