Anomalous diffusion of particles has been described by the time-fractional reaction–diffusion equation. A hybrid formulation of numerical technique is proposed to solve the time-fractional-order reaction–diffusion (FRD) equation numerically. The technique comprises the semi-discretization of the time variable using an L1 finite-difference scheme and space discretization using the quintic Hermite spline collocation method. The hybrid technique reduces the problem to an iterative scheme of an algebraic system of equations. The stability analysis of the proposed numerical scheme and the optimal error bounds for the approximate solution are also studied. A comparative study of the obtained results and an error analysis of approximation show the efficiency, accuracy, and effectiveness of the technique.

The generalized Burgers–Huxley (gBH) equation arises as a canonical model in nerve-pulse propagation (generalizing the Hodgkin–Huxley/FitzHugh–Nagumo excitable-media framework), in population dynamics with Allee-threshold reaction kinetics, and in nonlinear wave propagation in dispersive media; accurate benchmark solutions are essential for quantitative predictions in these domains. We couple the fourth-order exponential time differencing scheme ETDRK4 with a Chebyshev collocation spatial discretization and a linear boundary-lifting procedure to solve the gBH equation on a bounded interval with non-homogeneous Dirichlet data. On the canonical Ismail–Raslan–Rabboh travelling-wave benchmark the scheme attains $L^{\infty }$ errors at the level of floating-point round-off (∼10⁻¹⁹ absolute, ∼10⁻¹⁵ relative) with as few as $N=2$ collocation points and a single time step of size $\Delta t=1.0$—that is, three total nodes and one ETDRK4 advance. In strongly nonlinear regimes ($\gamma =0.1, 0.3, 0.5, 0.9$) the scheme exhibits approximately $\mathcal{O}(\Delta t^{2.45})$ temporal convergence across all four parameter values, consistent with the classical Hochbruck–Ostermann order reduction for exponential integrators on parabolic PDEs with non-homogeneous Dirichlet data. Used as a high-accuracy probe, the scheme provides a diagnostic of independent interest: the wave-speed formula of Wang, Zhu and Lu, still appearing as the exact-solution benchmark in numerical studies as recently as 2020, does not satisfy the partial differential equation. The corrected formula stated by Deng and verified symbolically by Appadu and Tijani is the unique value that makes the travelling-wave ansatz a genuine solution. We derive the residual associated with Wang’s formula in closed form, $R=\frac{\gamma A_1}{2}(A_2-A_2^{\mathrm{W}})(1-v^2)$, and show both analytically and numerically that reported errors for schemes benchmarked against Wang’s formula coincide with the analytical wave-profile gap $\frac{\gamma A_1}{2}|A_2-A_2^{\mathrm{W}}|$ rather than with true scheme accuracy. At the Ismail benchmark this gap equals $3.748\times {10}^{-7}$, which matches the N- and $\Delta t$-independent plateau observed when the scheme is measured against Wang’s profile. In the nerve-pulse and excitable-media interpretation, the two formulas correspond to action-potential propagation speeds of opposite sign at the Ismail benchmark, underscoring that the correction is not a mere algebraic curiosity but changes the qualitative physical prediction of the model.

The generalized Burgers–Huxley (gBH) equation arises as a canonical model in nerve-pulse propagation (generalizing the Hodgkin–Huxley/FitzHugh–Nagumo excitable-media framework), in population dynamics with Allee-threshold reaction kinetics, and in nonlinear wave propagation in dispersive media; accurate benchmark solutions are essential for quantitative predictions in these domains. We couple the fourth-order exponential time differencing scheme ETDRK4 with a Chebyshev collocation spatial discretization and a linear boundary-lifting procedure to solve the gBH equation on a bounded interval with non-homogeneous Dirichlet data. On the canonical Ismail–Raslan–Rabboh travelling-wave benchmark the scheme attains \( L^\infty \) errors at the level of floating-point round-off (\( \sim 10^{-19} \) absolute, \( \sim 10^{-15} \) relative) with as few as N=2 collocation points and a single time step of size Δt = 1.0—that is, three total nodes and one ETDRK4 advance. In strongly nonlinear regimes (γ = 0.1, 0.3, 0.5, 0.9) the scheme exhibits approximately \( \mathcal{O}(\Delta t^{2.45}) \) temporal convergence across all four parameter values, consistent with the classical Hochbruck–Ostermann order reduction for exponential integrators on parabolic PDEs with non-homogeneous Dirichlet data. Used as a high-accuracy probe, the scheme provides a diagnostic of independent interest: the wave-speed formula of Wang, Zhu and Lu (1990), still appearing as the exact-solution benchmark in numerical studies as recently as 2020, does not satisfy the partial differential equation. The corrected formula stated by Deng (2008) and verified symbolically by Appadu and Tijani (2022) is the unique value that makes the travelling-wave ansatz a genuine solution. We derive the residual associated with Wang's formula in closed form, \( R=\tfrac{\gamma A_1}{2}(A_2-A_2^{\mathrm{W}})(1-v^2) \), and show both analytically and numerically that reported errors for schemes benchmarked against Wang's formula coincide with the analytical wave-profile gap \( \tfrac{\gamma A_1}{2}|A_2-A_2^{\mathrm{W}}| \) rather than with true scheme accuracy. At the Ismail benchmark this gap equals 3.748 × 10^-7, which matches the N- and Δt-independent plateau observed when the scheme is measured against Wang's profile. In the nerve-pulse and excitable-media interpretation, the two formulas correspond to action-potential propagation speeds of opposite sign at the Ismail benchmark, underscoring that the correction is not a mere algebraic curiosity but changes the qualitative physical prediction of the model.