their forecasting capabilities [3; 6; 7]. Nevertheless, the singular nature of the RLFI's kernel, \((t-\tau)^{\alpha-1}\), presents substantial computational hurdles, as conventional numerical methods frequently encounter limitations in accuracy or incur high computational costs. Consequently, the development of efficient and high-precision approximation techniques for the RLFI is of paramount importance for advancing computational modeling across physics, engineering, and biology, where fractional calculus is increasingly employed to tackle real-world problems involving non-local behavior or fractal structures. Existing approaches to RLFI approximation include wavelet-based methods [7; 8; 9; 10; 11; 12], polynomial and orthogonal function techniques [13; 14; 15; 16; 17], finite difference and quadrature schemes [18; 19; 20; 21], operational matrix methods [22; 23; 24; 25], local meshless techniques [26], alternating direction implicit schemes [27], and radial basis functions [28]. While these methods have shown promise in specific contexts, they often struggle with trade-offs between accuracy, computational cost, and flexibility, particularly when adapting to diverse problem characteristics or varying fractional orders. This study introduces the GBFA method, which overcomes these challenges through three key innovations: (i) _Parameter adaptability_, utilizing the tunable SG parameters \(\lambda\) (for interpolation) and \(\lambda_{q}\) (for quadrature approximation) to optimize performance across a wide range of problems; (ii) _Super-exponential convergence_, achieving rapid error decay for smooth functions, often reaching near machine-precision with modest node counts; and (iii) _Computational efficiency_, enabled by precomputable FSGIMs that minimize runtime costs. The proposed GBFA method addresses these challenges by using the orthogonality and flexibility of SG polynomials to achieve super-exponential convergence. This approach offers near machine-precision accuracy with minimal computational effort, particularly for systems necessitating repeated fractional integrations or displaying symmetric patterns. Numerical experiments demonstrate that the GBFA method significantly outperforms established tools like MATLAB's integral function and MATHEMATICA's NIntegrate, achieving up to two orders of magnitude higher accuracy in certain cases. It also surpasses prior methods, such as the trapezoidal approach of [PERSON] [29], the spline-based techniques of [PERSON] and [PERSON] [30], and the neural network method [31], in both precision and efficiency. Sensitivity analysis reveals that setting \(\lambda_{q}<\lambda\) often accelerates convergence, while rigorous error bounds confirm super-exponential decay under optimal parameter choices. The FSGIMs's invariance for fixed points and parameters enables precomputation, making the GBFA method ideal for problems requiring repeated fractional integrations with varying \(\alpha\). While the method itself exploits the mathematical properties of orthogonal polynomials, its application can be crucial in analyzing systems where symmetry or repeating patterns are fundamental, as fractional calculus is used to model phenomena with memory where past states influence future behavior, and identifying symmetries in such systems can simplify analysis and prediction.