Infertility affects millions of couples worldwide, and in vitro fertilization (IVF) is the foremost assisted reproductive technology. However, IVF success rates vary considerably due to the memory-dependent nature of hormonal regulation, embryo development, and repeated treatment cycles. Clinical data are typically collected at discrete monthly intervals, yet few mathematical models simultaneously capture both the discrete time structure and the memory effects inherent in IVF. This work develops and rigorously analyses a discrete fractional-order (FO) IVF model using the dual Caputo nabla fractional difference (CNFD) operator, which naturally incorporates long-range memory while aligning with the cycle-based format of medical records. The model stratifies the IVF process into six compartments: infertile couples, patients under treatment, high-, medium-, and low-quality embryos, and positive pregnancy outcomes. Existence, uniqueness, and non-negativity of solutions are proved. Equilibrium analysis yields a unique positive endemic equilibrium, and sufficient conditions for global asymptotic stability (GAS) and Mittag-Leffler stability (MLS) are established via a novel Volterra-type discrete Lyapunov function (LF). Numerical simulations, performed with representative parameter values, confirm the theoretical results and show that lower FOs introduce stronger memory effects and slower convergence, thereby reproducing realistic, protracted IVF dynamics. The theoretical findings are not tied to a specific experimental dataset and thus provide a general framework that could assist clinicians in optimizing treatment decisions and personalizing patient management.

This work introduces a novel spectral collocation scheme based on ultraspherical (Gegenbauer) polynomials evaluated at Chebyshev–Gauss–Lobatto nodes for the numerical treatment of the nonlinear inhomogeneous time-fractional FitzHugh–Nagumo differential problem. The proposed methodology exploits the flexibility of the ultraspherical parameter $\lambda >-\frac{1}{2}$ to achieve enhanced accuracy and stability. We derived new operational matrices for both integer-order and Caputo fractional derivatives of the shifted ultraspherical basis, accompanied by rigorous proofs. A comprehensive convergence analysis in the $L^2$ norm was established, demonstrating spectral accuracy. Extensive numerical experiments confirmed that the proposed method outperforms the classical Legendre-based approach for optimal choices of $\lambda$, with errors reduced by several orders of magnitude. The superiority of optimized $\lambda$ over the Legendre case ( $\lambda =\frac{1}{2}$) was demonstrated through both numerical benchmarks and theoretical error bounds that explicitly depend on $\lambda$, showing that the Legendre choice is not universally optimal for problems with boundary layers or specific regularity properties. An efficient algorithmic implementation was provided, and comparative tables illustrate the superiority of the ultraspherical framework across various fractional orders and parameter settings.

Motivated by the challenges faced by standard methods in solving nonlinear fractional and q-fractional models with strong memory effects, this study develops numerical approaches capable of handling these complex behaviors more effectively. The proposed techniques are tested on four representative fractional and q-fractional initial value problems for several values of the order $\alpha \in$ $[0.5,1]$ and $q\in (0,1)$. In particular, the major aim of this work is to propose two generalizations of the higher-order Taylor method: The first one is the fractional Taylor method, and the second one is the q-fractional Taylor method. These methods will then be used to find approximate solutions for several fractional and q-fractional initial value problems. Numerous numerical comparisons will be performed to verify the effectiveness of our proposed generalizations.

Several researchers have dealt with the one-dimensional fractional heat conduction equation in the last decades, but as far as we know, no one has investigated such a problem from the perspective of developing suitable fractional-order methods. This has actually motivated us to address this problem by the way of establishing a proper fractional approach that involves employing a combination of a novel fractional difference formula to approximate the Caputo differentiator of order $\alpha$ coupled with the modified three-point fractional formula to approximate the Caputo differentiator of order $2\alpha$, where $0<\alpha \le 1$. As a result, the fractional heat conduction equation is then reexpressed numerically using the aforementioned formulas, and by dividing the considered mesh into multiple nodes, a system is generated and algebraically solved with the aid of MATLAB. This would allow us to obtain the desired approximate solution for the problem at hand.