Fractional Modified Euler Method (FMEM)

• FMEM is a family of discretization schemes for fractional differential and stochastic equations that incorporate correction terms to address memory effects and nonlocal behavior.
• It extends traditional Euler methods by achieving high-order convergence rates and robust error analysis in both deterministic and stochastic fractional models.
• Fast implementations, such as the SOE-based approach, significantly reduce computational cost while maintaining the theoretical accuracy of the method.

The Fractional Modified Euler Method (FMEM) refers to a comprehensive family of strong-order time-discretization schemes for differential and stochastic differential equations with fractional or nonlocal time structure. FMEM algorithms are distinguished from standard Euler and Euler–Maruyama (EM) methods by the presence of correction terms that account for the memory effects and non-Markovian increments characteristic of fractional operators, especially those driven by fractional Brownian motion (fBm) or Riemann–Liouville type derivatives. FMEMs achieve optimal or near-optimal convergence orders in both deterministic and stochastic fractional models, provide mathematically natural extensions to Stratonovich calculus for $H \to 1/2$, and enable efficient high-order numerical solutions for multi-term and time-variable fractional evolution equations (Zhang et al., 2022, Hu et al., 2014, Liu et al., 2017, Ngondiep, 2022, Hu et al., 2013).

1. Fundamental Structure of Fractional Modified Euler Methods

FMEMs are defined for models in which either the underlying dynamics incorporates fractional time derivatives (e.g., Riemann–Liouville, Caputo), or the driving noise is an fBm with Hurst parameter $H \in (0,1)$. In both classes, a prototypical FMEM one-step update takes the form

$$\begin{aligned} Y_{n+1} &= y_0 - \sum_{i=1}^m \frac{h^{1-\alpha_i}}{\Gamma(1-\alpha_i)} \sum_{j=0}^n (n+1-j)^{-\alpha_i} Y_j \ &\quad + h \sum_{j=0}^n f(t_j,Y_j) + \sum_{j=0}^n g(t_j,Y_j) \Delta W_j, \end{aligned}$$

for multi-term Riemann–Liouville SFDEs with $0 < \alpha_1 < ... < \alpha_m < 1$ (Zhang et al., 2022). For SDEs driven by fBm ($H > 1/3$), the FMEM corrects the classical Euler increment by a deterministic or random term dependent on the local regularity: $$y_{k+1} = y_k + b(y_k)\,h + V(y_k)\, (B_{t_{k+1}} - B_{t_k}) + \tfrac12\,(\partial V V)(y_k) \,h^{2H}$$ (Liu et al., 2017, Hu et al., 2014). In weakly singular (e.g., Caputo, variable-order) deterministic problems, multi-step and predictor–corrector FMEMs achieve high stability and accuracy by tailored convolution quadrature and spatial discretizations (Ngondiep, 2022).

2. Convergence Properties and Error Analysis

The strong convergence rate of FMEMs critically depends on the memory exponent or the regularity of the stochastic driver. For stochastic systems with Riemann–Liouville fractional derivatives,

$$\mathbb{E}|Y_n - y(t_n)|^2 = O\bigl(h^{\min\{2(1-\alpha_m),1\}}\bigr), \qquad \text{order } \min\{1-\alpha_m, 1/2\},$$

where $\alpha_m$ is the largest fractional order in the multi-term equation (Zhang et al., 2022). For SDEs driven by fBm, the strong order for $H \in (1/2,1)$ is $2H-1/2$, improving upon the naive Euler's $H$-order and converging to $1/2$ as $H \to 1/2$ (Hu et al., 2014, Hu et al., 2013). In the rough case ($H \in (1/3, 1/2)$), the order reduces to $n^{1/2-2H}$ (Liu et al., 2017). This improvement derives from the inclusion of correction terms compensating for local increments' nontrivial second-order moments.

For the two-step fourth-order FMEM/Crank–Nicolson scheme in deterministic advection–dispersion models with Caputo derivatives, the method achieves temporal order $O(k)$ and spatial order $O(h^4)$, with unconditional $L^\infty(0,T; L^2)$ stability for $0 < \bar\beta < 2/3$ (Ngondiep, 2022). The fourth-order spatial accuracy arises via compact finite-difference operators consistent with the theoretical smoothness requirements.

3. Algorithmic Formulations and Fast Implementation

Typical FMEMs require the evaluation and storage of convolution sums or iterated increments with nonlocal weights. For multi-term SFDEs, the direct implementation incurs $O(N^2)$ cost, motivating a fast FMEM based on the sum-of-exponentials (SOE) approximation of power-law kernels: $t^{-\alpha} \approx \sum_{q=1}^{N_\exp} \omega_q e^{-s_q t}, \quad | \cdot | \le \varepsilon,$ with $$N_{\exp}^{(\alpha)} = O((\ln \varepsilon^{-1})(\ln\ln\varepsilon^{-1} + \ln(T/\delta)))$$. This allows each mode to be updated recursively in $O(1)$ time, reducing overall cost to $O(N \log \varepsilon^{-1})$ without loss of theoretical accuracy (Zhang et al., 2022). For fBm-driven schemes, efficient simulation strategies for correlated increments are required (e.g., circulant embedding, Cholesky).

The predictor–corrector FMEM for time-fractional PDEs employs two interlinked sub-steps advancing from $t_n$ to $t_{n+1}$, using banded pentadiagonal matrices and convolution weights for the Caputo operator (Ngondiep, 2022). The careful composition of explicit and Crank–Nicolson updates, together with high-order finite-differences, is key to its unconditional stability and efficiency.

4. Theoretical Techniques and Analytical Tools

FMEM analysis relies on fractional calculus (Weyl–Malliavin integration by parts, bounds for Riemann–Liouville/Caputo integrals), discrete Grönwall-type inequalities adapted to fractional integrals, and Malliavin calculus for quantifying strong and weak errors in SDEs. Error decompositions typically separate the stochastic drift, memory/convolution history, and high-order Taylor remainders. The dominant error term in the stochastic setting is often a weighted Lévy area, whose scaling in $L^p$ determines the rate.

The FMEM is shown to be closer to Stratonovich discretizations for $H \to 1/2$, recapturing the classical EM method for Brownian-driven SDEs. Central limit theorems (CLT) for the error process in rough SDEs ($H \leq 3/4$) reveal convergence to solutions of linear SDEs driven by matrix-valued Brownian motion, with a transition to non-Gaussian limits (generalized Rosenblatt processes) for $H > 3/4$ (Hu et al., 2014).

In deterministic fractional PDEs, energy methods, summation-by-parts for high-order finite-difference operators, and discrete convolution inequalities yield unconditional stability and optimal convergence bounds.

5. Numerical Performance and Practical Considerations

Empirical results demonstrate that FMEMs achieve their predicted convergence rates. For multi-term stochastic fractional systems, the observed root mean square error is $O(h^{\min\{1-\alpha_m,0.5\}})$; fast FMEM implementations using SOE exhibit more than an order of magnitude CPU-time savings over direct $O(N^2)$ methods, particularly as $N$ increases (e.g., $n=1024$: 500s for direct, 50s for fast FMEM) (Zhang et al., 2022).

For time-variable advection–dispersion problems, the two-step FMEM/Crank–Nicolson method is unconditionally stable, achieving $O(k + h^4)$ error, and outperforms standard explicit-Euler and Crank–Nicolson schemes in both accuracy and computational efficiency (Ngondiep, 2022).

Implementation of FMEM for fBm-driven systems requires the reproducible generation of fBm increments, deterministic correction terms based on $H$, and high-order derivatives of drift/diffusion coefficients. For high-dimensional problems, efficient storage and update of correction terms and convolution sums is critical.

6. Connections, Extensions, and Scope

FMEM forms the foundation for higher-order fractional schemes (e.g., fractional Milstein), robust simulation of option pricing in rough volatility models, stochastic control with fractional noise, parameter estimation, and numerical solutions of subdiffusive or anomalous transport phenomena.

The connection between FMEM and rough paths theory is particularly pronounced for SDEs with rough drivers ($H \leq 1/2$), where lifting the process and its error into geometric rough paths spaces enables sharp error bounds and limit theorems via the sewing lemma and weighted Young integrals (Liu et al., 2017).

A plausible implication is that FMEMs, by virtue of their natural extension of classic Euler/EM methods to the fractional and rough settings, serve as a reference discretization for future advances in high-order and adaptively stabilized schemes for fractional models across stochastic and deterministic domains.

7. Summary Table: FMEM Core Properties Across Model Classes

The following table summarizes key properties of the FMEM variants documented in the cited works.

Application Domain	Model Class	FMEM Strong Order	Notable Features	Source
Stochastic SFDEs	Multi-term RL, $0<\alpha_m<1$	$\min\{1-\alpha_m, 0.5\}$	SOE-based acceleration ($O(N)$ ops)	(Zhang et al., 2022)
SDEs with fBm, $H>1/2$	Additive/multiplicative fBm	$2H-1/2$	Correction for optimal strong order	(Hu et al., 2014, Hu et al., 2013)
SDEs with fBm, $1/3	“Rough” regime	$n^{1/2-2H}$	Error CLT, rough paths structure	(Liu et al., 2017)
Fractional PDEs	Caputo/variable-order	$O(k+h^4)$	Fourth-order compact differences	(Ngondiep, 2022)

Each FMEM variant is tailored to the structure of its respective equation class, employing bespoke correction or convolution terms, and is supported by rigorous error analysis, often using advanced tools from fractional and Malliavin calculus.

---

References: (Zhang et al., 2022): "A modified EM method and its fast implementation for multi-term Riemann-Liouville stochastic fractional differential equations" (Hu et al., 2014): "Rate of convergence and asymptotic error distribution of Euler approximation schemes for fractional diffusions" (Liu et al., 2017): "First-order Euler scheme for SDEs driven by fractional Brownian motions: the rough case" (Ngondiep, 2022): "Unconditional Stability Of A Two-Step Fourth-Order Modified Explicit Euler/Crank-Nicolson Approach For Solving Time-Variable Fractional Mobile-Immobile Advection-Dispersion Equation" (Hu et al., 2013): "Modified Euler approximation scheme for stochastic differential equations driven by fractional Brownian motions"