An Efficient Exponential Sum Approximation of Power-Law Kernels for Solving Fractional Differential Equation

In this work, we present a comprehensive framework for approximating the weakly singular power-law kernel $t^{\alpha-1}$ of fractional integral and differential operators, where $\alpha \in (0,1)$ and $t \in [\delta,T]$ with $0<\delta<T<\infty$, using a finite sum of exponentials. This approximation method begins by substituting an exponential function into the Laplace transform of the power function, followed by the application of the trapezoidal rule to approximate the resulting integral. To ensure computational feasibility, the integral limits are truncated, leading to a finite exponential sum representation of the kernel. In contrast to earlier approaches, we pre-specify the admitted computational cost (measured in terms of the number of exponentials) and minimize the approximation error. Furthermore, to reduce the computational cost while maintaining accuracy, we present a two-stage algorithm based on Prony's method that compresses the exponential sum. The compressed kernel is then embedded into the Riemann-Liouville fractional integral and applied to solve fractional differential equations. To this end, we discuss two solution strategies, namely (a) method based on piecewise constant interpolation and (b) a transformation of the original fractional differential equation into a system of first-order ordinary differential equations (ODEs). This reformulation makes the problem solvable by standard ODE solvers with low computational cost while retaining the accuracy benefits of the exponential-sum-approximation. Finally, we apply the proposed strategies to solve some well-known fractional differential equations and demonstrate the advantages, accuracy, and the experimental order of convergence of the methods through numerical results.