Dynamical Analysis of Rubella Disease Model in the Context of Fractional Piecewise Derivative: Simulations with Real Statistical Data

Abstract

Rubella is a viral disease that can lead to severe health complications, especially in pregnant women and their unborn babies. Understanding the dynamics of the Rubella disease model is crucial for developing effective strategies to control its spread. This paper introduces a major innovation by employing a novel piecewise approach that incorporates two different kernels. This innovative approach significantly enhances the accuracy of modeling Rubella disease dynamics. In the first interval, the Caputo operator is employed to address initial conditions, while the Atangana–Baleanu derivative is utilized in the second interval to account for anomalous diffusion processes. A thorough theoretical analysis of the piecewise derivative for the problem is provided, discussing mathematical properties, stability, and convergence. To solve the proposed problem effectively, the piecewise numerical Newton polynomial technique is employed and the numerical scheme for both kernels is established. Through extensive numerical simulations with various fractional orders, the paper demonstrates the approach’s effectiveness and flexibility in modeling the spread of the Rubella virus. Furthermore, to validate the findings, the simulated results are compared with real data obtained from Rubella outbreaks in Uganda and Tanzania, confirming the practical relevance and accuracy of this innovative model.

1. Introduction

A wide range of infectious diseases exist in human society, some of which have little effect on our health, while others are fatal. Rubella was originally recorded in the mid-eighteenth century. Friedrich Hoffmann originally reported rubella in 1740, and it was later verified by de-Bergen (1752) and Orlow (1758) [1]. Rubella, sometimes called German measles, is a contagious and potentially life-threatening disease that mostly affects children and young people. The disease is characterized by the appearance of a rash on the face, which spreads to the chest and limbs and usually disappears after a period of 3 days. This rash is the most common sign of rubella infection, according to [2]. People’s lymph nodes and skin are particularly sensitive to this viral infection. Congenital rubella syndrome (CRS) and fetal death and are the most common complications in pregnant women. Rubella is a viral infection, and vaccination is only option to protect oneself (MMR vaccine). Because the illness is so mild, half of individuals are unaware they are afflicted, according to [3]. Rubella is a common virus in various countries of the globe, with about 100,000 infections of the CRS documented each year [4]. It spreads through a variety of routes, which includes direct contact to infected people and airborne droplets of infected persons during sneezing, coughing, or speaking. It is possible that an individual who comes into contract with it will be unaware of it for about a week or two. In contrast, the outbreak only lasts 3 to 5 days.

Mathematical models for infections can be useful in understanding how these diseases behave. Using methods from differential calculus, numerous models have been developed to better explain infectious illnesses such as in [5], where Jajarmi et al. investigated a new fractional mathematical model involving a non-singular derivative operator to discuss the clinical implications of diabetes and tuberculosis coexistence. A space–time spectral order sinc-collocation technique is used to investigate the heat model in viscoelasticity [6]. Furthermore, the predictor–corrector compact difference approach for nonliear fractional differential equations is presented in [7]. The integer order nonlinear HIV/AIDS infection model is extended to the non- integer nonlinear model [8]. Further, a general formulation for the SIRV epidemiological model is presented as a system of fractional order derivatives with respect to time to characterize some infectious diseases alongside the proportion of $u_1$ and $u_2$, which describe of vaccination and treatment, respectively [9]. For the deterministic modeling of epidemics, which consists of first-order ordinary differential equations with an autonomous nature, classical (integer-order) derivatives have been utilized [10]. An example of such a model is the rubella epidemic presented below.

$$\begin{array}{ccc}\hfill \frac{d\mathbb{S}}{d\mathsf{t}}& =& \mathcal{A}-[\nu +Q+\omega ]\mathbb{S}\left(\mathsf{t}\right),\hfill \\ \hfill \frac{d\mathbb{E}}{d\mathsf{t}}& =& \nu \mathbb{S}\left(\mathsf{t}\right)-(\xi +\omega )\mathbb{E}\left(\mathsf{t}\right),\hfill \\ \hfill \frac{d\mathbb{I}}{d\mathsf{t}}& =& \xi \mathbb{E}\left(\mathsf{t}\right)-(\epsilon +\omega )\mathbb{I}\left(\mathsf{t}\right),\hfill \\ \hfill \frac{d\mathbb{R}}{d\mathsf{t}}& =& \epsilon \mathbb{I}\left(\mathsf{t}\right)-\omega \mathbb{R}\left(\mathsf{t}\right),\hfill \\ \hfill \frac{d\mathbb{V}}{d\mathsf{t}}& =& \mathcal{D}\mathbb{V}\left(\mathsf{t}\right)-\omega \mathbb{V}\left(\mathsf{t}\right).\hfill \end{array}$$

(1)

With initial conditions

$$\begin{array}{c}\hfill \mathbb{S}\left(0\right)={\mathbb{S}}_0>0,\mathbb{E}\left(0\right)={\mathbb{E}}_0\ge 0,\mathbb{I}\left(0\right)={\mathbb{I}}_0\ge 0,\mathbb{R}\left(0\right)={\mathbb{R}}_0\ge 0,\mathbb{V}\left(0\right)={\mathbb{V}}_0\ge 0.\end{array}$$

In the above model, $\mathbb{S}\left(\mathsf{t}\right)$ stands for susceptible, $\mathbb{E}\left(\mathsf{t}\right)$ stands for latent, $\mathbb{I}\left(\mathsf{t}\right)$ shows infected, $\mathbb{R}\left(\mathsf{t}\right)$ shows recovered, and $\mathbb{V}\left(\mathsf{t}\right)$ represents vaccinated individuals. $\mathcal{A}$ shows those individuals that have been immunized with vaccination, and $\nu$ stands for force of infection at the age $\eta$ at $\mathsf{t}$. Further, $\xi$ and $\epsilon$ stand for latent rate and infection rate, respectively.

Classical models, on the other hand, rely on local differential and integral operators that lack the ability to preserve the memory of the epidemic under investigation. As a result, the memory characteristics of the underlying epidemic are ignored in classical calculus. Nonlocal operators must be utilized to preserve memory effects in the deterministic epidemic model due to their superiority over classical operators, as demonstrated in a number of recent research investigations. For example, in [11], Atangana and Qureshi have modeled attractors of chaotics systems. The applications of fractional calculus to epidemiology are discussed in [12]. The reaction disffusion model of dengue disease is studied in [13]. Further, the dengue model with the novel fractional operator is studied in [14].

Fractional calculus, notably introduced by Riemann–Liouville in 1832, deals with generalizing the concepts of integer-order derivatives and integrals to non-integer orders. In 1967, Caputo provided a practical definition of fractional operators, which has found applications in various fields, describing global phenomena in a wide range of disciplines. The Caputo indicator, also known as the Caputo derivative, is introduced in research contexts to capture complex dynamics, enhance modeling accuracy, address real-world scenarios, and overcome limitations associated with traditional derivatives. This is especially important when dealing with phenomena that exhibit non-standard, fractional-order behavior, making the Caputo derivative a valuable tool for more accurately representing and understanding such complex systems. In 2015, the Caputo–Fabrizio derivative operator was introduced as an update to the traditional Caputo operator. This new fractional derivative eliminated singular kernels, leading to more stable and well-behaved calculations. Its benefits include improved accuracy and the ability to analyze complex systems and phenomena more effectively [15]. In 2016, Atangana and Baleanu introduced the $\mathrm{ABC}$ operator [16], providing a valuable tool for solving global problems that require modern calculus [17]. The Atangana–Baleanu derivative, a key component of the ABC operator, plays a crucial role in capturing memory effects in various systems by incorporating the power-law kernel. Memory effects refer to the phenomenon where a system’s current state depends not only on its immediate past state but also on its past states over an extended period. Such memory effects can arise in a wide range of physical, biological, and engineering systems. The Atangana–Baleanu derivative enables the modeling of systems with long-term memory and anomalous diffusion, providing a more accurate description of complex behaviors that cannot be adequately captured by traditional integer-order calculus. Several researchers have used different types of operators and applied them to numerous types of mathematical models such as the finiancial bubble mathematical model [18], the dicrete time SIR model [19], COVID-19 pandemic [20], energy prices [21], system risks and spillover networks [22], and the prey–predator model [23]. Additionally, Qurashi [24] used the fractal-fractional operator to design a mathematical model for the rubella epidemic while taking care of dimensional consistency among the model equations.

Atangana and Araz [25] recently constructed a novel sorts of operators known as piecewise integrals and derivatives. These operators were formulated to address the limitations of exponential or Mittag–Leffler kernels in specifying the timing of crossover. Their unique piecewise derivative approach was presented in [25] as a solution to these challenges. Subsequently, researchers explored the crossover behaviors via these operators in a new way, for instance, Huanglomngbing disease [26], COVID-19 transmission [27], rompurs spreading models [28], the tumor immune interaction model [29], and the Rota virus disease [30]. Inspired by the advantages of these operators, we will apply the piecewise Caputo and ABC piecewise operators to investigate the model (1) in the following manner:

$$\begin{array}{c}\hfill \begin{array}{cc}& {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{S}\left(\mathsf{t}\right)\right)=\mathcal{A}-[\nu +Q+\omega ]\mathbb{S}\left(\mathsf{t}\right),\hfill \\ & {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{E}\left(\mathsf{t}\right)\right)=\nu \mathbb{S}\left(\mathsf{t}\right)-(\xi +\omega )\mathbb{E}\left(\mathsf{t}\right),\hfill \\ & {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{I}\left(\mathsf{t}\right)\right)=\xi \mathbb{E}\left(\mathsf{t}\right)-(\epsilon +\omega )\mathbb{I}\left(\mathsf{t}\right),\hfill \\ & {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{R}\left(\mathsf{t}\right)\right)=\epsilon \mathbb{I}\left(\mathsf{t}\right)-\omega \mathbb{R}\left(\mathsf{t}\right),\hfill \\ & {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{V}\left(\mathsf{t}\right)\right)=\mathcal{D}\mathbb{V}\left(\mathsf{t}\right)-\omega \mathbb{V}\left(\mathsf{t}\right),\hfill \end{array}\end{array}$$

(2)

where, $\eta \in (0,1],\mathsf{t}\in [0,\mathsf{T}]$. A more detailed form of the model (2) is given as:

$$\begin{array}{c}\hfill {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{S}\left(\mathsf{t}\right)\right)=\left\{\begin{array}{l}{}_0^C{D}_{\mathsf{t}}^{\eta }\left(\mathbb{S}\left(\mathsf{t}\right)\right){=}^C{\mathbf{Q}}_1(\mathbb{S},\mathsf{t}),0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ {}_0^{ABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{S}\left(\mathsf{t}\right)\right){=}^{ABC}{\mathbf{Q}}_1(\mathbb{S},\mathsf{t}),{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.,\\ \hfill {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{E}\left(\mathsf{t}\right)\right)=\left\{\begin{array}{c}{}_0^C{D}_{\mathsf{t}}^{\eta }\left(\mathbb{E}\left(\mathsf{t}\right)\right){=}^C{\mathbf{Q}}_2(\mathbb{E},\mathsf{t}),0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ {}_0^{ABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{E}\left(\mathsf{t}\right)\right){=}^{ABC}{\mathbf{Q}}_2(\mathbb{E},\mathsf{t}),{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.,\\ \hfill {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{I}\left(\mathsf{t}\right)\right)=\left\{\begin{array}{c}{}_0^C{D}_{\mathsf{t}}^{\eta }\left(\mathbb{I}\left(\mathsf{t}\right)\right){=}^C{\mathbf{Q}}_3(\mathbb{I},\mathsf{t}),0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ {}_0^{ABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{I}\left(\mathsf{t}\right)\right){=}^{ABC}{\mathbf{Q}}_3(\mathbb{I},\mathsf{t}),{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.,\\ \hfill {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{R}\left(\mathsf{t}\right)\right)=\left\{\begin{array}{c}{}_0^C{D}_{\mathsf{t}}^{\eta }\left(\mathbb{R}\left(\mathsf{t}\right)\right){=}^C{\mathbf{Q}}_4(\mathbb{R},\mathsf{t}),0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ {}_0^{ABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{R}\left(\mathsf{t}\right)\right){=}^{ABC}{\mathbf{Q}}_4(\mathbb{R},\mathsf{t}),{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.,\\ \hfill {}_0^{CABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{V}\left(\mathsf{t}\right)\right)=\left\{\begin{array}{c}{}_0^C{D}_{\mathsf{t}}^{\eta }\left(\mathbb{V}\left(\mathsf{t}\right)\right){=}^C{\mathbf{Q}}_5(\mathbb{V},\mathsf{t}),0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ {}_0^{ABC}{D}_{\mathsf{t}}^{\eta }\left(\mathbb{V}\left(\mathsf{t}\right)\right){=}^{ABC}{\mathbf{Q}}_5(\mathbb{V},\mathsf{t}),{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.,\end{array}$$

(3)

where, ${}_0^C{D}_{\mathsf{t}}^{\eta }$ is the Caputo and ${}_0^{ABC}{D}_{\mathsf{t}}^{\eta }$ is the $ABC$ derivative operator.

This article applies the piecewise operator with two different kernels: the $\mathrm{ABC}$ and Caputo fractional operators. The Caputo operator is used for the first interval, while the $\mathrm{ABC}$ derivative is applied for the second interval. The existence and uniqueness of the solution are examined within the framework of the piecewise derivative. The dual-kernel structure of our model requires careful consideration of different intervals, making the Newton polynomial technique [31] a suitable choice due to its ability to approximate functions over separate intervals while maintaining smooth transitions. Therefore, the approximate solution is evaluated using the technique of Newton polynomials, establishing the solutions for both the $\mathrm{ABC}$ and Caputo fractional operators. The numerical simulation section compares the real data from Uganda and Tanzania.

The paper’s structure unfolds as follows: Section 2 is dedicated to an in-depth examination of preliminaries, establishing the foundational concepts. Transitioning to Section 3, we delve into the exploration of existence and uniqueness results, shedding light on the core theoretical aspects. In Section 4, our attention shifts towards the development and discussion of the numerical scheme, a critical component of our methodology. Finally, Section 5 serves as a concluding reflection, where we provide a concise summary encapsulating the achievements, insights, and key findings uncovered in this study.

2. Preliminaries

Some fundamental definitions are provided in this section.

Definition 1 

([16]). Let $\mathcal{M}\left(\mathsf{t}\right)\in {\mathcal{H}}^1(0,\mathsf{t})$, with fractional order $0<\eta <1$; then,

$$\begin{array}{c}\hfill {}^{ABC}{{}_{0}D}_{\mathsf{t}}^{\eta }\left(\mathcal{M}\left(\mathsf{t}\right)\right)=\frac{ABC\left(\eta \right)}{1-\eta }{\int }_0^{\mathsf{t}}\frac{d}{d\varsigma }\mathcal{M}\varsigma E_{\eta }\left[\frac{-\eta }{1-\eta }{\left(\mathsf{t}-\varsigma \right)}^{\eta }\right]d\varsigma ,\end{array}$$

(4)

is called $ABC$ derivative of a function $\mathcal{M}\left(\mathsf{t}\right)$ and the generalization function defined as $ABC\left(\eta \right)=1-\eta +\frac{\eta }{\mathsf{\Gamma }\left(\eta \right)}$. For the fractional $ABC$ derivative, we have $ABC\left(0\right)=ABC\left(1\right)=1$.

While the integration can be written as:

$$\begin{array}{c}\hfill {}^{ABC}{{}_{0}I}_{\mathsf{t}}^{\eta }\mathcal{M}\left(\mathsf{t}\right)=\frac{1-\eta }{ABC\left(\eta \right)}\mathcal{M}\left(\varsigma \right)\left(\mathsf{t}\right)+\frac{\eta }{\mathsf{\Gamma }\left(\eta \right)ABC\left(\eta \right)}{\int }_0^{\mathsf{t}}\mathcal{M}{(\mathsf{t}-\varsigma )}^{\eta -1}d\varsigma .\end{array}$$

(5)

Definition 2 

([16]). For $\mathcal{M}\left(\mathsf{t}\right)\in C[0,T]$, the fractional operator in Caputo sense is

$$\begin{array}{c}\hfill {}_0^C{D}_{\mathsf{t}}^{\eta }\mathcal{M}\left(\mathsf{t}\right)=\frac{1}{\mathsf{\Gamma }(1-\eta )}{\int }_0^{\mathsf{t}}\frac{d}{d\varsigma }M\varsigma {(\mathsf{t}-\varsigma )}^{n-\eta -1}d\varsigma .\end{array}$$

Definition 3 

([25]). The piecewise differential operator in the sense of Caputo and ABC for a piecewise differentiable $\mathcal{M}\left(\mathsf{t}\right)$ is:

$$\begin{array}{c}\hfill {}_0^{PCABC}{D}_{\mathsf{t}}^{\eta }\mathcal{M}\left(\mathsf{t}\right)=\left\{\begin{array}{c}{}_0^C{D}_{\mathsf{t}}^{\eta }\mathcal{M}\left(\mathsf{t}\right),0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ {}_0^{ABC}{D}_{\mathsf{t}}^{\eta }\mathcal{M}\left(\mathsf{t}\right){\mathsf{t}}_1<\mathsf{t}\le \mathsf{T}.\hfill \end{array}\right.\end{array}$$

Definition 4 

([25]). The corresponding piecewise integral for $\mathcal{M}\left(\mathsf{t}\right)$ is

$$\begin{array}{c}\hfill {}_0^{PCABC}{I}_{\mathsf{t}}\mathcal{M}\left(\mathsf{t}\right)=\left\{\begin{array}{c}\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}\mathcal{M}\left(\varsigma \right){(\mathsf{t}-\varsigma )}^{\eta -1}d,0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \frac{1-\eta }{ABC\left(\eta \right)}\mathcal{M}\left(\mathsf{t}\right)+\frac{\eta }{ABC\eta \mathsf{\Gamma }\eta }{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}\mathcal{M}\left(\varsigma \right){(\mathsf{t}-\varsigma )}^{\eta -1}d{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.\end{array}$$

here, ${}_0^{PCABC}{I}_{\mathsf{t}}^{\eta }$ represents the piecewise integral operator.

3. Existence and Uniqueness Results

The existence and uniqueness result of the suggested model (2) in the piecewise notion are found in this part. In order to do this, the system (3) can be written as follows.

$$\begin{array}{c}\hfill {}_0^{PCABC}{D}_{\mathsf{t}}^{\eta }\mathcal{L}\left(\mathsf{t}\right)=\mathbf{Q}(\mathsf{t},\mathcal{L}\left(\mathsf{t}\right)),0<\eta \le 1,\end{array}$$

(6)

is

$$\begin{array}{c}\hfill \mathcal{L}\left(\mathsf{t}\right)=\left\{\begin{array}{c}{\mathcal{L}}_0+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{\mathsf{t}}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right)){(\mathsf{t}-\varsigma )}^{\eta -1}d\varsigma ,0<\mathsf{t}\le {\mathsf{t}}_1\hfill \\ \mathcal{L}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{ABC\left(\eta \right)}\mathbf{Q}(\mathsf{t},\mathcal{L}\left(\mathsf{t}\right))+\frac{\eta }{ABC\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))d\left(\varsigma \right),{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.\end{array}$$

(7)

where

$$\mathcal{L}\left(\mathsf{t}\right)=\left\{\begin{array}{c}\mathbb{S}\left(\mathsf{t}\right)\hfill \\ \mathbb{E}\left(\mathsf{t}\right)\hfill \\ \mathbb{I}\left(\mathsf{t}\right)\hfill \\ \mathbb{R}\left(\mathsf{t}\right)\hfill \\ \mathbb{V}\left(\mathsf{t}\right)\hfill \end{array}\right.,{\mathcal{L}}_0=\left\{\begin{array}{c}{\mathbb{S}}_0\hfill \\ {\mathbb{E}}_0\hfill \\ {\mathbb{I}}_0\hfill \\ {\mathbb{R}}_0\hfill \\ {\mathbb{V}}_0\hfill \end{array}\right.,{\mathcal{L}}_{\left({\mathsf{t}}_1\right)}=\left\{\begin{array}{c}{\mathbb{S}}_{\left({\mathsf{t}}_1\right)}\hfill \\ {\mathbb{E}}_{\left({\mathsf{t}}_1\right)}\hfill \\ {\mathbb{I}}_{\left({\mathsf{t}}_1\right)}\hfill \\ {\mathbb{R}}_{\left({\mathsf{t}}_1\right)}\hfill \\ {\mathbb{V}}_{\left({\mathsf{t}}_1\right)}\hfill \end{array}\right.,\mathbf{Q}(\mathsf{t},\mathcal{L}\left(\mathsf{t}\right))=\left\{\begin{array}{c}\hfill {\mathbf{Q}}_1=\left\{\begin{array}{c}{}^C{\mathbf{Q}}_1(\mathbb{S},\mathsf{t})\hfill \\ {}^{ABC}{\mathbf{Q}}_1(\mathbb{S},\mathsf{t}),\hfill \end{array}\right.\\ \hfill {\mathbf{Q}}_2=\left\{\begin{array}{c}{}^C{\mathbf{Q}}_2(\mathbb{E},\mathsf{t})\hfill \\ {}^{ABC}{\mathbf{Q}}_2(\mathbb{E},\mathsf{t}),\hfill \end{array}\right.\\ \hfill {\mathbf{Q}}_3=\left\{\begin{array}{c}{}^C{\mathbf{Q}}_3(\mathbb{I},\mathsf{t})\hfill \\ {}^{ABC}{\mathbf{Q}}_3(\mathbb{I},\mathsf{t}),\hfill \end{array}\right.\\ \hfill {\mathbf{Q}}_4=\left\{\begin{array}{c}{}^C{\mathbf{Q}}_4(\mathbb{R},\mathsf{t})\hfill \\ {}^{ABC}{\mathbf{Q}}_4(\mathbb{R},\mathsf{t}),\hfill \end{array}\right.\\ \hfill {\mathbf{Q}}_5=\left\{\begin{array}{c}{}^C{\mathbf{Q}}_5(\mathbb{V},\mathsf{t})\hfill \\ {}^{ABC}{\mathbf{Q}}_5(\mathbb{V},\mathsf{t}).\hfill \end{array}\right.\end{array}\right.$$

(8)

Taking $\infty >{\mathsf{t}}_2\ge \mathsf{t}>{\mathsf{t}}_1>0$ and the Banach space $\mathsf{\Lambda }=C[0,\mathsf{T}]$ with

$$\parallel \mathcal{L}\parallel =\underset{\mathsf{t}\in [0,\mathsf{T}]}{max}\left|\mathcal{L}\left(\mathsf{t}\right)\right|.$$

Further,

(C1) 

∃$L_{\mathcal{L}}>0$ such that ∀$\mathcal{L},\overline{\mathcal{L}}\in \mathsf{\Lambda }$, one may have

$$|\mathbf{Q}(\mathsf{t},\mathcal{L})-\mathbf{Q}(\mathsf{t},\overline{\mathcal{L}})|\le L_{\mathcal{L}}|\mathcal{L}-\overline{\mathcal{L}}|.$$

(C2) 

∃$C_{\mathbf{Q}}>0$ & $M_{\mathbf{Q}}>0,$ one have

$$\left|\mathbf{Q}(\mathsf{t},\mathcal{L}\left(\mathsf{t}\right))\right|\le C_{\mathbf{Q}}\left|\mathcal{L}\right|+M_{\mathbf{Q}}.$$

Theorem 1. 

If $\mathbf{Q}$ is piecewise continuous on $0<\mathsf{t}\le {\mathsf{t}}_1$ and ${\mathsf{t}}_1<\mathsf{t}\le {\mathsf{t}}_2$ on $[0,\mathsf{T}]$, also satisfying $\left(C2\right)$, then (6) has at least one solution.

Proof. 

Consider $\mathsf{\Theta }\subseteq \mathsf{\Lambda }$, which is given by

$$\mathsf{\Theta }=\{\mathcal{L}\in \mathsf{\Lambda }:\parallel \mathcal{L}\parallel \le R_{1,2},R>0\},$$

Define $\mathbb{H}:\mathsf{\Theta }\to \mathsf{\Theta }$ by

$$\begin{array}{c}\hfill \mathbb{H}\left(\mathcal{L}\right)=\left\{\begin{array}{c}{\mathcal{L}}_0+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right)){(t-\varsigma )}^{\eta -1}d\varsigma ,0<\mathsf{t}\le {\mathsf{t}}_1\hfill \\ \mathcal{L}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{ABC\left(\eta \right)}\mathbf{Q}(\mathsf{t},\mathcal{L}\left(\mathsf{t}\right))+\frac{\eta }{ABC\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^t{(\mathsf{t}-\varsigma )}^{\eta -1}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))d\left(\varsigma \right),{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.\end{array}$$

(9)

for $\mathcal{L}\in \mathsf{\Theta }$, one can obtain

$$\begin{array}{ccc}\hfill \left|\mathbb{H}\right(\mathcal{L}\left)\right(\mathsf{t}\left)\right|& \le & \left\{\begin{array}{c}|{\mathcal{L}}_0|+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1}\left|\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))\right|d\varsigma ,\hfill \\ |{\mathcal{L}}_{(}{\mathsf{t}}_1\left)\right|+\frac{1-\eta }{ABC\left(\eta \right)}\left|\mathbf{Q}(\mathsf{t},\mathcal{L}\left(\mathsf{t}\right))\right|+\frac{\eta }{ABC\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}\left|\mathbf{Q}\left(\varsigma \mathcal{L}\left(\varsigma \right)\right)\right|d\left(\varsigma \right),\hfill \end{array}\right.\hfill \\ & \le & \left\{\begin{array}{c}|{\mathcal{L}}_0|+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1}[C_{\mathbf{Q}}\left|\mathcal{L}\right|+M_{\mathbf{Q}}]d\varsigma ,\hfill \\ |{\mathcal{L}}_{(}{\mathsf{t}}_1\left)\right|+\frac{1-\eta }{ABC\left(\eta \right)}[C_{\mathbf{Q}}\left|\mathcal{L}\right|+M_{\mathbf{Q}}]+\frac{\eta }{ABC\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}[C_{\mathbf{Q}}\left|\mathcal{L}\right|+M_{\mathbf{Q}}]d\left(\varsigma \right),\hfill \end{array}\right.\hfill \\ & \le & \left\{\begin{array}{c}|{\mathcal{L}}_0|+\frac{{\mathbf{T}}^{\eta }}{\mathsf{\Gamma }(\eta +1)}[C_H\left|\mathbb{U}\right|+M_{\mathbf{Q}}]=R_1,0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ |{\mathcal{L}}_{(}{\mathsf{t}}_1\left)\right|+\frac{1-\eta }{ABC\left(\eta \right)}[C_{\mathbf{Q}}\left|\mathcal{L}\right|+M_{\mathbf{Q}}]+\frac{\eta {(T-\mathbf{T})}^{\eta }}{ABC\left(\eta \right)\mathsf{\Gamma }\eta +1}[C_{\mathbf{Q}}\left|\mathcal{L}\right|+M_{\mathbf{Q}}]d\left(\varsigma \right)=R_2,{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T},\hfill \end{array}\right.\hfill \\ & \le & \left\{\begin{array}{c}{R}_1,0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ R_2,{\mathsf{t}}_1<\mathsf{t}\le T.\hfill \end{array}\right.\hfill \end{array}$$

Since $\mathcal{L}\in \mathsf{\Theta }$. Hence, $\mathbb{H}\left(\mathbf{\Theta }\right)\subset \mathbf{\Theta }$. So, $\mathbb{H}$ is closed and complete. Next, taking ${\mathsf{t}}_i<{\mathsf{t}}_j\in [0,{\mathsf{t}}_1]$ gives

$$\begin{array}{ccc}\hfill |\mathbb{H}\left(\mathcal{L}\right)\left({\mathsf{t}}_j\right)-\mathbb{H}\left(\mathcal{L}\right)\left({\mathsf{t}}_i\right)|& =& |\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_j}{({\mathsf{t}}_j-\varsigma )}^{\eta -1}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))d\varsigma \hfill \\ & -& \frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_i}{({\mathsf{t}}_i-\varsigma )}^{\eta -1}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))d\varsigma |\hfill \\ & \le & \frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_i}[{({\mathsf{t}}_i-\varsigma )}^{\eta -1}-{({\mathsf{t}}_j-\varsigma )}^{\eta -1}]\left|\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))\right|d\varsigma \hfill \\ & +& \frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_i}^{{\mathsf{t}}_j}{({\mathsf{t}}_j-\varsigma )}^{\eta -1}\left|\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))\right|d\varsigma \hfill \\ & \le & \frac{1}{\mathsf{\Gamma }\left(\eta \right)}[{\int }_0^{{\mathsf{t}}_i}[{({\mathsf{t}}_i-\varsigma )}^{\eta -1}-{({\mathsf{t}}_j-\varsigma )}^{\eta -1}]d\varsigma \hfill \\ & +& {\int }_{{\mathsf{t}}_i}^{{\mathsf{t}}_j}{({\mathsf{t}}_j-\varsigma )}^{\eta -1}d\varsigma \left]\right(C_H\left|\mathcal{L}\right|+M_{\mathbf{Q}})\hfill \\ & \le & \frac{(C_{\mathbf{Q}}\mathcal{L}+M_{\mathbf{Q}})}{\mathsf{\Gamma }(\eta +1)}[{\mathsf{t}}_j^{\eta }-{\mathsf{t}}_i^{\eta }+2{({\mathsf{t}}_j-{\mathsf{t}}_i)}^{\eta }].\hfill \end{array}$$

(10)

This implies

$$|\mathbb{H}\left(\mathcal{L}\right)\left({\mathsf{t}}_j\right)-\mathbb{H}\left(\mathcal{L}\right)\left({\mathsf{t}}_i\right)|\to 0,\mathrm{as}{\mathsf{t}}_i\to {\mathsf{t}}_j.$$

Hence, $\mathbb{H}$ is equi continuous in $[0,{\mathsf{t}}_1]$. Now, for ${\mathsf{t}}_i,{\mathsf{t}}_j\in [{\mathsf{t}}_1,\mathsf{T}]$, one can have

$$\begin{array}{ccc}\hfill |\mathbb{H}\left(\mathcal{L}\right)\left({\mathsf{t}}_j\right)-\mathbb{H}\left(\mathcal{L}\right)\left({\mathsf{t}}_i\right)|& =& |\frac{1-\eta }{ABC\left(\eta \right)}\mathbf{Q}(\mathsf{t},\mathcal{L}\left(\mathsf{t}\right))+\frac{\eta }{\mathsf{\Gamma }\left(\eta \right)ABC\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_j}{({\mathsf{t}}_j-\varsigma )}^{\eta -1}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))d\varsigma ,\hfill \\ & -& \frac{1-\eta }{ABC\left(\eta \right)}\mathbf{Q}(\mathsf{t},\mathcal{L}\left(\mathsf{t}\right))+\frac{\left(\eta \right)}{ABC\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_i}{({\mathsf{t}}_i-\varsigma )}^{\eta -1}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))d\varsigma |\hfill \\ & \le & \frac{\eta }{ABC\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_i}[{({\mathsf{t}}_i-\varsigma )}^{\eta -1}-{({\mathsf{t}}_j-\varsigma )}^{\eta -1}]\left|\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))\right|d\varsigma \hfill \\ & +& \frac{\eta }{ABC\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_i}^{{\mathsf{t}}_j}{({\mathsf{t}}_j-\varsigma )}^{\eta -1}\left|\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))\right|d\varsigma \hfill \\ & \le & \frac{\eta }{ABC\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}[{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_i}[{({\mathsf{t}}_i-\varsigma )}^{\eta -1}-{({\mathsf{t}}_j-\varsigma )}^{\eta -1}]d\varsigma \hfill \\ & +& {\int }_{{\mathsf{t}}_i}^{{\mathsf{t}}_j}{({\mathsf{t}}_j-\varsigma )}^{\eta -1}d\varsigma \left]\right(C_{\mathbf{Q}}\left|\mathcal{L}\right|+M_{\mathbf{Q}})\hfill \\ & \le & \frac{\eta (C_{\mathbf{Q}}\mathcal{L}+M_{\mathbf{Q}})}{ABC\left(\eta \right)\mathsf{\Gamma }(\eta +1)}[{\mathsf{t}}_j^{\eta }-{\mathsf{t}}_i^{\eta }+2{({\mathsf{t}}_j-{\mathsf{t}}_i)}^{\eta }].\hfill \end{array}$$

(11)

As a result, the map $\mathbb{H}$ is equi-continuous. The operator $\mathbb{H}$ satisfies the properties of being completely continuous, uniformly continuous, and bounded by the Arzela’-Ascoli result. Additionally, the piecewise problem (3) has at least one solution within each sub-interval, by the Schauder result. □

Theorem 2. 

With assumption $\left(C1\right)$, and $\frac{{\mathbf{T}}^{\eta }}{\mathsf{\Gamma }(\eta +1)}\le 1$, the considered model has a unique solution.

Proof. 

Since $\mathbb{H}:\mathbf{\Theta }\to \mathbf{\Theta }$ is continuous, taking $\mathcal{L}$ and $\overline{\mathcal{L}}\in \mathsf{\Theta }$, for $\mathsf{t}\in [0,{\mathsf{t}}_1]$ one has

$$\begin{array}{ccc}\hfill \parallel \mathbb{H}\left(\mathcal{L}\right)-\mathbb{H}\left(\overline{\mathcal{L}}\right)\parallel & =& \underset{\mathsf{t}\in [0,{\mathsf{t}}_1]}{max}|\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}\mathbf{Q}(\varsigma ,\mathcal{L}\left(\varsigma \right))d\varsigma -\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}\mathbf{Q}(\varsigma ,\overline{\mathcal{W}}\left(\varsigma \right))d\varsigma |\hfill \\ & \le & \frac{{\mathbf{T}}^{\eta }}{\mathsf{\Gamma }(\eta +1)}{L}_{\mathbf{Q}}\parallel \mathcal{L}-\overline{\mathcal{L}}\parallel .\hfill \end{array}$$

(12)

From (12), it follows that

$$\begin{array}{c}\hfill \parallel \mathbb{H}\left(\mathcal{L}\right)-\mathbb{H}\left(\overline{\mathcal{L}}\right)\parallel \le \frac{{\mathbf{T}}^{\eta }}{\mathsf{\Gamma }(\eta +1)}{L}_{\mathbf{Q}}\parallel \mathcal{L}-\overline{\mathcal{L}}\parallel .\end{array}$$

(13)

Therefore, from the Banach theorem the solution of the suggested problem is unique in $\mathsf{t}\in [0,{\mathsf{t}}_1]$. Now, for $\mathsf{t}\in [{\mathsf{t}}_1,\mathsf{T}]$ with the $ABC$ operator, one can obtain

$$\begin{array}{c}\hfill \parallel \mathbb{H}\left(\mathcal{L}\right)-\mathbb{H}\left(\overline{\mathcal{L}}\right)\parallel \le \frac{1-\eta }{ABC\left(\eta \right)}{L}_{\mathbf{Q}}\parallel \mathcal{L}-\overline{\mathcal{L}}\parallel +\frac{\eta \left({\mathbf{T}-T}^{\eta }\right)}{ABC\left(\eta \right)\mathsf{\Gamma }(\eta +1)}{L}_{\mathbf{Q}}\parallel \mathcal{L}-\overline{\mathcal{L}}\parallel ,\end{array}$$

(14)

or

$$\begin{array}{c}\hfill \parallel \mathbb{H}\left(\mathcal{L}\right)-\mathbb{H}\left(\overline{\mathcal{L}}\right)\parallel \le L_{\mathbf{Q}}\left[\frac{1-\eta }{ABC\left(\eta \right)}+\frac{\eta {(T-\mathbf{T})}^{\eta }}{ABC\left(\eta \right)\mathsf{\Gamma }(\eta +1)}\right]\parallel \mathcal{L}-\overline{\mathcal{L}}\parallel .\end{array}$$

(15)

Hence, from the Banach theorem solution the considered problem is unique. So, from Equations (13) and (15) the model has a unique solution on considered sub-intervals. □

4. Numerical Scheme

Here, the numerical algorithm for the Rubella model (2) in the sense of the Caputo and ABC piecewise derivatives is presented. The Rubella disease model (2) in the context of the fractional derivative is studied in [32]; further in [33], this model is analyzed by the application of the Adam Bashforth scheme in the fractal fractional scenario. Here, by applying the piece-wise integral, one can obtain

$$\begin{array}{ccc}\hfill \mathbb{S}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{S}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_1(\mathsf{t},\mathbb{S})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{S}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_1(\mathsf{t},\mathbb{S})d\varsigma +\frac{\eta }{AB\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_1(\mathsf{t},\mathbb{S})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T,\hfill \end{array}\right.\hfill \\ \hfill \mathbb{E}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{E}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_2(\mathsf{t},\mathbb{E})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{E}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_2(\mathsf{t},\mathbb{E})d\varsigma +\frac{\eta }{AB\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_2(\mathsf{t},\mathbb{E})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T,\hfill \end{array}\right.\hfill \\ \hfill \mathbb{I}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{I}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_3(\mathsf{t},\mathbb{I})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{I}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_3(\mathsf{t},\mathbb{I})d\varsigma +\frac{\eta }{AB\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_3(\mathsf{t},\mathbb{I})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T,\hfill \end{array}\right.\hfill \\ \hfill \mathbb{R}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{R}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_4(\mathsf{t},\mathbb{R})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{R}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_4(\mathsf{t},\mathbb{R})d\varsigma +\frac{\eta }{AB\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_4(\mathsf{t},\mathbb{R})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T,\hfill \end{array}\right.\hfill \\ \hfill \mathbb{V}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{V}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_5(\mathsf{t},\mathbb{V})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{V}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_5(\mathsf{t},\mathbb{V})d\varsigma +\frac{\eta }{AB\left(\eta \right)\mathsf{\Gamma }\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{\mathsf{t}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_5(\mathsf{t},\mathbb{V})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T.\hfill \end{array}\right.\hfill \\ \hfill \end{array}$$

Considering $\mathsf{t}={\mathsf{t}}_{n+1}$, this gives:

$$\begin{array}{ccc}\hfill \mathbb{S}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{S}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_1(\mathsf{t},\mathbb{S})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{S}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_1(\mathsf{t},\mathbb{S})d\varsigma +\frac{\eta }{\mathsf{\Gamma }\left(\eta \right)AB\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_{n+1}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_1(\mathsf{t},\mathbb{S})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T,\hfill \end{array}\right.\hfill \\ \hfill \mathbb{E}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{E}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_2(\mathsf{t},\mathbb{E})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{E}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_2(\mathsf{t},\mathbb{E})d\varsigma +\frac{\eta }{\mathsf{\Gamma }\left(\eta \right)AB\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_{n+1}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_2(\mathsf{t},\mathbb{E})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T,\hfill \end{array}\right.\hfill \\ \hfill \mathbb{I}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{I}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_3(\mathsf{t},\mathbb{I})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{I}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_3(\mathsf{t},\mathbb{I})d\varsigma +\frac{\eta }{\mathsf{\Gamma }\left(\eta \right)AB\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_{n+1}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_3(\mathsf{t},\mathbb{I})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T,\hfill \end{array}\right.\hfill \\ \hfill \mathbb{R}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{R}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_4(\mathsf{t},\mathbb{R})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{R}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_4(\mathsf{t},\mathbb{R})d\varsigma +\frac{\eta }{\mathsf{\Gamma }\left(\eta \right)AB\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_{n+1}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_4(\mathsf{t},\mathbb{R})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le T,\hfill \end{array}\right.\hfill \\ \hfill \mathbb{V}\left(\mathsf{t}\right)& =& \left\{\begin{array}{c}\mathbb{V}\left(0\right)+\frac{1}{\mathsf{\Gamma }\left(\eta \right)}{\int }_0^{{\mathsf{t}}_1}{(\mathsf{t}-\varsigma )}^{\eta -1c}{\mathbf{Q}}_5(\mathsf{t},\mathbb{V})d\varsigma 0<\mathsf{t}\le {\mathsf{t}}_1,\hfill \\ \mathbb{V}\left({\mathsf{t}}_1\right)+\frac{1-\eta }{AB\left(\eta \right)}{\mathbf{Q}}_5(\mathsf{t},\mathbb{V})d\varsigma +\frac{\eta }{\mathsf{\Gamma }\left(\eta \right)AB\left(\eta \right)}{\int }_{{\mathsf{t}}_1}^{{\mathsf{t}}_{n+1}}{(\mathsf{t}-\varsigma )}^{\eta -1}{\mathbf{Q}}_5(\mathsf{t},\mathbb{V})d\varsigma ,{\mathsf{t}}_1<\mathsf{t}\le \mathsf{T}.\hfill \end{array}\right.\hfill \\ \hfill \end{array}$$

After applying the Newton interpolation formula, one can obtain

$$\begin{array}{c}\hfill \mathbb{S}\left({\mathsf{t}}_{n+1}\right))=\left\{\begin{array}{c}{\mathbb{S}}_0+\left\{\begin{array}{c}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =2}^i\left[{}^C{H}_1({\mathbb{S}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =2}^i{[}^C{H}_1({\mathbb{S}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {-}^C{H}_1({\mathbb{S}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{2\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =2}^i{[}^C{H}_1({\mathbb{S}}^{\kappa },{\mathsf{t}}_{\kappa })-2^C{H}_1({\mathbb{S}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^C{H}_1({\mathbb{S}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\aleph \hfill \\ \end{array}\right\}\hfill \\ \mathbb{S}\left({\mathsf{t}}_1\right)+\left\{\begin{array}{c}\frac{1-\eta }{ABC\left(\eta \right)}{}^{ABC}{\mathbf{Q}}_1({\mathbb{S}}^n,{\mathsf{t}}_n)\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =i+3}^n\left[{}^{ABC}{\mathbf{Q}}_1({\mathbb{S}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_1({\mathbb{S}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_1({\mathbb{S}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_1({\mathbb{S}}^{\kappa },{\mathsf{t}}_{\kappa })\hfill \\ -2^{ABC}{\mathbf{Q}}_1({\mathbb{S}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_1({\mathbb{S}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mathsf{\Delta }.\hfill \end{array}\right\},\hfill \end{array}\right.\end{array}$$

(16)

$$\begin{array}{c}\hfill \mathbb{E}\left({\mathsf{t}}_{n+1}\right))=\left\{\begin{array}{c}{\mathbb{E}}_0+\left\{\begin{array}{c}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =2}^i\left[{}^C{H}_2({\mathbb{E}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =2}^i{[}^C{H}_2({\mathbb{E}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {-}^C{H}_2({\mathbb{E}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{2\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =2}^i{[}^C{H}_2({\mathbb{E}}^{\kappa },{\mathsf{t}}_{\kappa })-2^C{H}_2({\mathbb{E}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^C{H}_2({\mathbb{E}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\aleph \hfill \\ \end{array}\right\}\hfill \\ \mathbb{E}\left({\mathsf{t}}_1\right)+\left\{\begin{array}{c}\frac{1-\eta }{ABC\left(\eta \right)}{}^{ABC}{\mathbf{Q}}_2({\mathbb{E}}^n,{\mathsf{t}}_n)\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =i+3}^n\left[{}^{ABC}{\mathbf{Q}}_2({\mathbb{E}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_2({\mathbb{E}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_2({\mathbb{E}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_2({\mathbb{E}}^{\kappa },{\mathsf{t}}_{\kappa })\hfill \\ -2^{ABC}{\mathbf{Q}}_2({\mathbb{E}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_2({\mathbb{E}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mathsf{\Delta }.\hfill \end{array}\right\},\hfill \end{array}\right.\end{array}$$

(17)

$$\begin{array}{c}\hfill \mathbb{I}\left({\mathsf{t}}_{n+1}\right))=\left\{\begin{array}{c}{\mathbb{I}}_0+\left\{\begin{array}{c}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =2}^i\left[{}^C{H}_2({\mathbb{I}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =2}^i{[}^C{H}_2({\mathbb{I}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {-}^C{H}_2({\mathbb{I}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{2\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =2}^i{[}^C{H}_2({\mathbb{I}}^{\kappa },{\mathsf{t}}_{\kappa })-2^C{H}_2({\mathbb{I}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^C{H}_2({\mathbb{I}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\aleph \hfill \\ \end{array}\right\}\hfill \\ \mathbb{I}\left({\mathsf{t}}_1\right)+\left\{\begin{array}{c}\frac{1-\eta }{ABC\left(\eta \right)}{}^{ABC}{\mathbf{Q}}_3({\mathbb{I}}^n,{\mathsf{t}}_n)\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =i+3}^n\left[{}^{ABC}{\mathbf{Q}}_3({\mathbb{I}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_3({\mathbb{I}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_3({\mathbb{I}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_3({\mathbb{I}}^{\kappa },{\mathsf{t}}_{\kappa })\hfill \\ -2^{ABC}{\mathbf{Q}}_3({\mathbb{I}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_3({\mathbb{I}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mathsf{\Delta }\hfill \end{array}\right\},\hfill \end{array}\right.\end{array}$$

(18)

$$\begin{array}{c}\hfill \mathbb{R}\left({\mathsf{t}}_{n+1}\right))=\left\{\begin{array}{c}{\mathbb{R}}_0+\left\{\begin{array}{c}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =2}^i\left[{}^C{H}_2({\mathbb{R}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =2}^i{[}^C{H}_2({\mathbb{R}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {-}^C{H}_2({\mathbb{R}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{2\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =2}^i{[}^C{H}_2({\mathbb{I}}^{\kappa },{\mathsf{t}}_{\kappa })-2^C{H}_2({\mathbb{R}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^C{H}_2({\mathbb{R}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\aleph \hfill \\ \end{array}\right\}\hfill \\ \mathbb{R}\left({\mathsf{t}}_1\right)+\left\{\begin{array}{c}\frac{1-\eta }{ABC\left(\eta \right)}{}^{ABC}{\mathbf{Q}}_4({\mathbb{I}}^n,{\mathsf{t}}_n)\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =i+3}^n\left[{}^{ABC}{\mathbf{Q}}_4({\mathbb{R}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_4({\mathbb{R}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_4({\mathbb{R}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_4({\mathbb{R}}^{\kappa },{\mathsf{t}}_{\kappa })\hfill \\ -2^{ABC}{\mathbf{Q}}_4({\mathbb{R}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_4({\mathbb{R}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mathsf{\Delta }\hfill \end{array}\right\},\hfill \end{array}\right.\end{array}$$

(19)

$$\begin{array}{c}\hfill \mathbb{V}\left({\mathsf{t}}_{n+1}\right))=\left\{\begin{array}{c}{\mathbb{V}}_0+\left\{\begin{array}{c}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =2}^i\left[{}^C{H}_2({\mathbb{V}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =2}^i{[}^C{H}_2({\mathbb{V}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {-}^C{H}_2({\mathbb{V}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{2\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =2}^i{[}^C{H}_2({\mathbb{I}}^{\kappa },{\mathsf{t}}_{\kappa })-2^C{H}_2({\mathbb{V}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^C{H}_2({\mathbb{V}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\aleph \hfill \\ \end{array}\right\}\hfill \\ \mathbb{V}\left({\mathsf{t}}_1\right)+\left\{\begin{array}{c}\frac{1-\eta }{ABC\left(\eta \right)}{}^{ABC}{\mathbf{Q}}_5({\mathbb{I}}^n,{\mathsf{t}}_n)\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +1)}\sum _{\kappa =i+3}^n\left[{}^{ABC}{\mathbf{Q}}_5({\mathbb{V}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})\right]\mathsf{\Pi }\hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{{\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +2)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_5({\mathbb{V}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_5({\mathbb{V}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mho \hfill \\ +\frac{\eta }{ABC\left(\eta \right)}\frac{\eta {\left(\mathrm{h}\right)}^{\eta -1}}{\mathsf{\Gamma }(\eta +3)}\sum _{\kappa =i+3}^n{[}^{ABC}{\mathbf{Q}}_5({\mathbb{V}}^{\kappa },{\mathsf{t}}_{\kappa })\hfill \\ -2^{ABC}{\mathbf{Q}}_5({\mathbb{V}}^{\kappa -1},{\mathsf{t}}_{\kappa -1})\hfill \\ {+}^{ABC}{\mathbf{Q}}_5({\mathbb{V}}^{\kappa -2},{\mathsf{t}}_{\kappa -2})]\mathsf{\Delta }\hfill \end{array}\right\}.\hfill \end{array}\right.\end{array}$$

(20)

Here,

$$\begin{array}{ccc}\hfill \mathsf{\Pi }& =& {(1+n+\kappa )}^{\eta }\left(2{(n-\kappa )}^2+(3\eta +10)(n-\kappa )+2{\eta }^2+9\eta +12\right)\hfill \\ & & -(n-\kappa )\left(2{(n-\kappa )}^2+(5\eta +10)(n-\kappa )+6{\eta }^2+18\eta +12\right),\hfill \end{array}$$

$$\begin{array}{ccc}\hfill \mho & =& {(1+\kappa +n)}^{\eta }\left(3+2\eta +n-\kappa \right)-(-\kappa +n)\left(n-\kappa +3\eta +3\right),\hfill \end{array}$$

$$\begin{array}{c}\hfill \aleph ={(1+\kappa +n)}^{\eta }-{(-\kappa +n)}^{\eta },\end{array}$$

and

$$\begin{array}{c}\hfill \begin{array}{c}{}^C{\mathbf{Q}}_1(\mathsf{t},\mathbb{S})={}^{ABC}{\mathbf{Q}}_1(\mathsf{t},\mathbb{S})=\mathcal{A}-[\nu (\xi ,\mathsf{t})+Q+\omega ]\mathbb{S}\left(\mathsf{t}\right),\hfill \\ {}^C{\mathbf{Q}}_2(\mathsf{t},\mathbb{E})={}^{ABC}{\mathbf{Q}}_2(\mathsf{t},\mathbb{E})=\nu \mathbb{S}\left(\mathsf{t}\right)-(\xi +\omega )\mathbb{E}\left(\mathsf{t}\right),\hfill \\ {}^C{\mathbf{Q}}_3(\mathsf{t},\mathbb{I})={}^{ABC}{\mathbf{Q}}_3(\mathsf{t},\mathbb{I})=\xi \mathbb{E}\left(\mathsf{t}\right)-(\epsilon +\omega )\mathbb{I}\left(\mathsf{t}\right)\hfill \\ {}^C{\mathbf{Q}}_4(\mathsf{t},\mathbb{R})={}^{ABC}{\mathbf{Q}}_4(\mathsf{t},\mathbb{R})=\epsilon \mathbb{I}\left(\mathsf{t}\right)-\omega \left(xi\right)\mathbb{R}\left(\mathsf{t}\right),\hfill \\ {}^C{\mathbf{Q}}_5(\mathsf{t},\mathbb{V})={}^{ABC}{\mathbf{Q}}_5(\mathsf{t},\mathbb{V})=\mathfrak{D}\mathbb{V}\left(\mathsf{t}\right)-\omega \mathbb{V}\left(\mathsf{t}\right).\hfill \end{array}\end{array}$$

5. Graphical Analysis

In this section, we present a graphical exploration of the system (2) using MATLAB-2020 for simulation. We conducted simulations using two distinct sets of parameters: set 1, with parameter values as detailed in

Table 1. Parameters and their values of model (2).

Parameter	Description	Value
$\mathcal{A}$	Transmission rate	0.8
$\nu$	Rate of infection	0.4
Q	Immunized with vaccines	0.3
$\omega$	Natural death-rate	0.1
$\xi$	Rate of recovery	0.4
$\epsilon$	Rate of exposure to disease	0.4
$\mathfrak{D}$	Rate of the second dose of vaccine	0.1

, and set 2, characterized by $[A=0.2,$ $\nu =0.5,Q=0.3,\omega =0.4,\xi =0.2,\epsilon =0.4,\mathfrak{D}=0.1]$. The initial values were set to $[\mathbb{S}0,\mathbb{E}0,\mathbb{I}0,\mathbb{R}0,{\mathbb{V}}_0]=[200,0,0,0,0]$.

For simulations with the parameter values from set 1, we showcase the results in Figure 1a–e. In contrast, simulations for set 2 parameters are presented in Figure 2a–e. To facilitate simulations, we divided the time interval into two sub-intervals: $(0,{\mathsf{t}}_1]=(0,5]$ and $({\mathsf{t}}_1,\mathsf{T}]=(5,20]$. We applied the ABC operator in the second half of the interval after initially utilizing the Caputo derivative in the first half. During simulations, we employed a step size of $dt=0.01$. Consequently, the curves in the second half of the interval demonstrate the behavior of the suggested problem with varying orders in the ABC sense, while the first half of the interval depicts the dynamics of the proposed model (2) under the Caputo operator.

In Figure 1a–e, the values of $\eta$ are represented using different colors and values: (blue, $0.99)$, (green, $0.96)$, (magenta, $0.95)$, (black, $0.94)$, and (red, $0.93)$. The plots were generated using MATLAB-2020.

Figure 1a, depicts the dynamics of susceptible people. Figure 1b,c show the evolution of latent and infected populations, respectively. Similarly, Figure 1d,e are projected to see the graphs of recovered and vaccinated individuals. In Figure 1a, one may notice a continuous decrease in the susceptible population $\mathbb{S}$ over time. Notably, a significant rapid decline occurs as we progress towards the second interval where the ABC operator is applied. Similarly, as depicted in Figure 1b, the latent population $\mathbb{E}$ displays an initial increase followed by a gradual decline. This behavior indicates a crossover at $\mathsf{t}={\mathsf{t}}_1$ and stabilizes after 8 days at a fractional order of 0.96. On the top of that, in Figure 1c,d one sees the epidemic’s infected individuals $\mathbb{I}$ and recovered individuals $\mathbb{R}$ reach their highest values around $\mathsf{t}=4$. After this peak, both populations decrease quickly and stabilize in the second sub-interval. Finally, it is observed that Figure 1e shows that at the start more individuals are vaccinated and that after 7 to 10 days there is no disease so that the number of vaccinations also becomes stable.

Next, in Figure 2a the dynamics in susceptible population are presented. Similarly, Figure 2b,c demonstrate the behavior of latent and infected populations, respectively. Further, Figure 2d,e are presented to see the dynamics in the recovered and vaccinated population of the considered epidemic model (2). Here, one can observe that at the start the susceptible population is high, which decreases with time and becomes stable after $\mathsf{t}=18$. Further, in Figure 2b,c the latent and infected populous shows that in these two classes the number of individuals gradually increases, which shows a rapid decrease as the system advances from the first to the second sub-interval where the ABC operator is applied and these two classes becomes stable, which shows that the disease dies out. Moreover, Figure 2d,e depict that the recovered population and those who are vaccinated show decreased and rapid stability at a lower fractional order after $\mathsf{t}=15$. This shows that the recovery rate increases and the rate of infection decreases as the rate of vaccinations increases.

To ensure the effectiveness of our proposed method, it is crucial to compare the theoretical results with real data in biological models. Thus, we present the comparisons between our simulated results and actual data in Figure 3a–c [34,35].

In Figure 3a, we demonstrate the comparison between real data from Uganda and our simulated data using the parameter values of $\mathcal{A}=1.16$, $\nu =0.28$, $\mathcal{Q}=0.5$, $\omega =0.395$, $\xi =0.2$, $ϵ=0.4$, and $\mathfrak{D}=0.1$. Moving on to Figure 3b, we estimate the parameter values to be $\mathcal{A}=1.16$, $\nu =0.28$, $\mathcal{Q}=0.5$, $\omega =0.4$, $\xi =0.2$, $ϵ=0.4$, and $\mathfrak{D}=0.1$, and we compare the simulated results with data from Tanzania. In this case, we consider fractional orders of $(\mathrm{blue},0.985)$, $(\mathrm{cyan},0.980)$, $(\mathrm{black},0.975)$, $(\mathrm{green},0.970)$, and $(\mathrm{magenta},0.960)$, highlighting the best-fitted dynamics.

Similarly, in Figure 3c we compare our simulated data with real data from Ethiopia using estimated parameter values of $\mathcal{A}=1.16$, $\nu =0.28$, $\mathcal{Q}=0.5$, $\omega =0.384$, $\xi =0.26$, $ϵ=0.4$, and $\mathfrak{D}=0.1$.

Through these comparisons, we observe that our proposed technique is highly effective for studying biological models as it accurately captures the best-fitted dynamics when compared to real data. This showcases the proficiency of our approach in understanding and analyzing these models. Overall, the medical importance of comparing simulated data with real data lies in its potential to improve disease modeling, treatment optimization, drug development, and public health planning. These applications have far-reaching implications, ultimately leading to better healthcare outcomes and advancements in medical research and practice.

The comparisons between the simulated results and real data have significant implications for disease prevention, vaccination campaigns, and policy-making decisions. By accurately capturing the dynamics of disease models and validating them against real data, we gain valuable insights that can contribute to various aspects of public health. These comparisons enable the assessment of intervention strategies, such as vaccination campaigns, by analyzing the alignment between simulated results and real data. Policymakers can make informed decisions regarding the implementation and optimization of vaccination programs based on this information. Additionally, the simulations provide insights into disease spread and transmission dynamics, allowing for the identification of critical periods or regions of high transmission. This knowledge enables the prioritization of targeted interventions, resource allocation, and the early detection of outbreaks, which are crucial for effective public health efforts.

Furthermore, the comparisons between simulated and real data contribute to the optimization of treatment strategies. By understanding the dynamics of disease models and their alignment with real-world outcomes, healthcare professionals can refine and tailor treatment protocols. This optimization of treatment strategies leads to improved patient outcomes and reduced disease burden. Moreover, the insights gained from the simulations can inform policy-making decisions in public health planning and response. Policymakers can assess the effectiveness of various control measures, such as quarantine measures or social distancing policies, based on the accurate representation of disease dynamics. This information guides the development of evidence-based policies aimed at mitigating the impact of diseases on populations, ultimately improving public health outcomes.

6. Conclusions

In this study, the behavior of the Rubella disease model is analyzed using a piecewise operator approach. Specifically, this manuscript investigated the Caputo and ABC operators in relation to the model. The study examines the existence and uniqueness of a solution with a piecewise derivative for this particular illness model. The piecewise Newton polynomial approach is utilized to obtain a rough estimate of the suggested problem’s solution. Additionally, simulations are conducted for the Rubella model considering different fractional orders.

Throughout the simulations, a decrease in the susceptible individuals ($\mathbb{S}$) over time has been notices, with a notable rapid decrease as the model progresses into the second interval where the ABC operator is used. It is also noticed that the dynamics of the recovered individuals and vaccinated individuals show a decrease and rapid stabilization at lower fractional orders. Furthermore, by comparing the last three figures, we conclude that the fractional piecewise operator accurately captures the dynamics of the considered disease. This is evident from the simulation results of the infected class, which closely match the real data on Rubella infections in Uganda and Tanzania, considering the fractional order dynamics.

Possible future work includes conducting a sensitivity analysis to assess the robustness of the Rubella disease model, validating the model using additional real-world data, exploring different intervention strategies, comparing the proposed piecewise operator approach with alternative models, and developing parameter estimation methods to improve model accuracy and data fit. These avenues of research would contribute to a deeper understanding of Rubella dynamics and aid in the development of effective control strategies for infectious diseases. Nowadays, the neural network [36,37,38,39,40,41,42] and integer and fractional order delay differential [43,44,45,46] are being used for investigation of different mathematical models in applied sciences. One can use these concepts to study the considered model in the future.