In this page you can see and download some Mathematica notebooks that show how to solve FDEs by means of the methods discussed in:
S. B. Yuste and L. Acedo, On an explicit finite difference method for fractional diffusion equations, arXiv:cs/0311011v1 [cs.NA] (2003); S. B. Yuste and L. Acedo, An explicit finite difference method and a new von Neumann-type stability analysis for fractional diffusion equations, SIAM Journal of Numerical Analysis, 52, 1862-74 (2005).
S. B. Yuste, Weighted average finite difference methods for fractional diffusion equations, Proceedings of FDA'04, pages 335-340 (2004) (see also arXiv:cs/0408053v1 [cs.NA] (2004) ) and S. B. Yuste, Weighted average finite difference methods for fractional diffusion equations, Journal of Computational Physics 216 (2006) 264-274. [DOI: 10.1016/j.jcp.2005.12.006]. One gets a fractional version of the Crank-Nicolson method if one chooses the weight factor λ=1/2.
Here you can find:
A Mathematica notebook that uses this method: FractionalWeightedAveragedMethod.nb