Construction and numerical analysis of a fuzzy non-standard computational method for the solution of an SEIQR model of COVID-19 dynamics

Abbreviations

The following abbreviations are used in this manuscript:

1.   Introduction

This study presents and examines a new nonlinear sequential Hilfer-Hadamard fractional-order integrodifferential equations (HHFIEs) with nonlocal coupled multipoint and Hadamard fractional integrodifferential boundary conditions. The formulation of the problem is as follows:

and it is enhanced by nonlocal coupled multipoint and Hadamard fractional integrodifferential boundary conditions:

Here $\alpha_{1}, \alpha_{2}\in (1, 2]$, $\beta_{1}, \beta_{2}\in [0, 1]$, $\lambda_{1}, \lambda_{2} \in \mathbb{R}_{+}$, $\mathcal{T} > 1$, $\eta_{\mathfrak{j}}, \theta_{\mathfrak{i}}, \lambda_{\mathfrak{k}}, \mathcal{P}_{\mathfrak{u}}, \mathcal{Q}_{\mathfrak{v}}, \mathcal{M}_{\mathfrak{w}} \in \mathbb{R}$, $\xi_{\mathfrak{j}}, \zeta_{\mathfrak{i}}, \mu_{\mathfrak{k}}, \psi_{\mathfrak{u}}, \sigma_{\mathfrak{v}}, \pi_{\mathfrak{w}} \in (1, \mathcal{T})$, $(\mathfrak{j} = 1, 2, ..., m$, $\mathfrak{i} = 1, 2, ..., n$, $\mathfrak{k} = 1, 2, ..., r$, $\mathfrak{u} = 1, 2, ..., a$, $\mathfrak{v} = 1, 2, ..., b$, $\mathfrak{w} = 1, 2, ..., c)$, ${}^{\mathcal{HH}} \mathcal{D}^{\alpha_{i}, \beta_{j}}_{1^{+}}$ denotes the Hilfer-Hadamard fractional derivative (HHFD) operator of order $\alpha_{i}, \beta_{j}; i = 1, 2. j = 1, 2.$ ${}^{\mathcal{H}} \mathcal{I}^{\chi}$ is the Hadamard fractional integral (HFI) of order $\chi \in \{\phi_{\mathfrak{i}}, \delta_{\mathfrak{v}} > 0\} \ \mathfrak{i} = 1, 2, ..., n, \ \mathfrak{v} = 1, 2, ..., b$, ${}^{\mathcal{H}} \mathcal{D}_{1}^{\psi}$ is the Hadamard fractional derivative (HFD) of order $\psi \in \{\omega_{\mathfrak{k}}, \vartheta_{\mathfrak{w}}\} \ \mathfrak{k} = 1, 2, ..., r, \ \mathfrak{w} = 1, 2, ..., c$, $\mathcal{F}, \mathcal{G} : [1, \mathcal{T}] \times \mathbb{R}^{4} \rightarrow \mathbb{R}$ are continuous functions. It's important to highlight and that this study makes a significant contribution to the existing literature by addressing a distinct setup involving sequential HHFIEs along with coupled multipoint and Hadamard fractional integrodifferential boundary conditions. The methodology employed in this study involves using the fixed-point approach to establish both existence and uniqueness results for the problems (1.1) and (1.2). The process involves converting the given problem into an equivalent fixed-point problem, followed by the application of Leray-Schauder alternative and Banach's fixed-point theorem to establish existence and uniqueness results, respectively. Additionally, we investigate the Ulam-Hyers stability of the solution for the proposed system. The findings of this study are innovative and contribute to the existing literature on boundary value problems (BVPs) concerning coupled systems of sequential HHFIEs. In recent decades, fractional calculus has gained considerable attention and become a prominent area of study in mathematical analysis. This growth is largely attributed to the extensive use of fractional calculus techniques in developing innovative mathematical models to represent various phenomena in fields such as economics, mechanics, engineering, science, and others. References ^[1,2,3,4] offer examples and comprehensive discussions on this subject.

In the upcoming section, we will provide an overview of relevant scholarly articles pertaining to the discussed problem. Among various fractional derivatives introduced, the Riemann-Liouville and Caputo fractional derivatives (CFDs) have garnered significant attention due to their practical applications. The Hilfer fractional derivative, introduced by Hilfer in ^[5], incorporates the Riemann-Liouville and CFDs as special cases for certain parameter values. Additional insights into this derivative can be found in ^{[6,7,8,9,10,11,12,13]}. References ^{[14,15,16,17,18]} offer valuable insights into Hilfer-type initial and BVPs. A recent study ^[19] explores the Ulam-Hyers stability and existence of a fully coupled system featuring integro-multistrip-multipoint boundary conditions and nonlinear sequential Hilfer fractional differential equations (HFDEs). Furthermore, ^[20] delves into a hybrid generalized HFDE BVP.

In 1892, Hadamard introduced the HFD, defined by a logarithmic function with an arbitrary exponent in its kernel ^[21]. Subsequent studies, such as those in ^{[22,23,24,25,26]}, have explored variations such as HHFDs and Caputo-Hadamard fractional derivatives (CHFDs). Notably, for specific values of $\beta$-$\beta = 0$ and $\beta = 1$, respectively—HFDs and CHFDs emerge as particular instances of the HHFD.

Stability analysis has been a prominent field of study for fractional differential equations in the last several decades and has drawn a lot of interest from scholars. Numerous stability models, including Lyapunov, exponential, and Mittag-Leffler stability, have been thoroughly examined in the literature. We suggest reviewing publications ^{[27,28,29,30,31]} for historical perspective on Ulam-Hyers stability and current improvements.

The problem of existence and Ulam stability of solutions for the following Hilfer-Hadamard fractional differential equations (HHFDEs) ^[32] is stated as follows:

where $0 < \alpha < 1$, $0 \leq \beta \leq 1$, $\gamma = \alpha + \beta - \alpha\beta$, $\mathcal{T} > 1$, $\varphi \in \mathbb{R}$, and $f : \mathcal{J} \times \mathbb{R} \to \mathbb{R}$ is a given function. ${}^{H}\mathcal{I}^{1-\gamma}_{1}$ denotes the left-sided mixed Hadamard integral of order $1 - \gamma$, and $^{H}\mathcal{D}^{\alpha, \beta}_{1}$ is the HHFD of order $\alpha$ and type $\beta$, introduced by Hilfer. In ^[33], existence results were established for an HHFDE with nonlocal integro-multipoint boundary conditions:

where $\alpha \in (1, 2]$, $\beta\in [0, 1]$, $\theta_{i}, \lambda \in \mathbb{R}$, $\eta, \xi_{i} \in (1, T)$ $(i = 1, 2, ..., m)$, ${^{H}{\mathcal{I}^{\delta}}}$ is the Hadamard fractional integral (HFI) of order $\delta > 0$, and $f : [1, {T}] \times \mathbb{R} \rightarrow \mathbb{R}$ is a continuous function. Problem (1.4) represents a non-coupled system, in contrast to problems (1.1)-(1.2), which is a coupled system. The systems (1.1)-(1.2) presents nonlocal coupled Hadamard fractional integrodifferential and multipoint boundary conditions, whereas the problem (1.4) involves discrete boundary conditions with HFIs. The authors of ^[34] established existence results for nonlocal mixed Hilfer-Hadamard fractional BVPs:

where $\alpha \in (1, 2]$, $\beta\in [0, 1]$, $\eta_{i}, \zeta_{i}, \lambda_{k} \in \mathbb{R}$, $\xi_{i}, \theta_{i}, \mu_{k} \in (1, T)$, $(j = 1, 2, ..., m), (i = 1, 2, ..., n), (k = 1, 2, ..., r)$, ${^{H}{\mathcal{I}^{\phi_{i}}}}$ is the HFI of order $\phi_{i} > 0$, ${_{H}{\mathcal{D}^{\mu_{k}}_{1}}}$ is the HFD of order $\mu_{k} > 0$, and $f : [1, {T}] \times \mathbb{R} \rightarrow \mathbb{R}$ is a continuous function. Equation (1.5) does not represent a coupled system, unlike Eqs (1.1) and (1.2), which does. In the latter, there are nonlocal coupled Hadamard fractional integrodifferential and multipoint boundary conditions. In contrast, Eq (1.5) involves multipoint boundary conditions comprising HFIs and HFDs. Furthermore, in ^[35], investigations were conducted on coupled HHFDEs within generalized Banach spaces. The authors of the aforementioned study ^[36] successfully derived existence results for a coupled system of HHFDEs with nonlocal coupled boundary conditions:

where $\alpha, \gamma\in (1, 2]$, $\beta, \delta\in [0, 1]$, $\mathcal{T} > 1$, ${}^{\mathcal{HH}} \mathcal{D}^{\alpha, \beta}$, ${}^{\mathcal{HH}} \mathcal{D}^{\gamma, \delta}_{1}$ denotes the HHFD operator of order $\alpha, \beta, \gamma, \delta$. ${}^{\mathcal{H}} \mathcal{D}_{1^{+}}^{\chi}$ is the HFD operator of order $\chi \in \{ \varsigma, \vartheta, \varrho_{i}, \eta_{i}, \sigma_{i}, \theta_{i}\}$, $(i = 1, 2, ..., m), (i = 1, 2, ..., n), (i = 1, 2, ..., p), (i = 1, 2, ..., q)$, $f, g : [1, T] \times \mathbb{R} \times \mathbb{R} \rightarrow \mathbb{R}$ are continuous functions. In the boundary conditions, Riemann-Stieltjes integrals with $\mathcal{H}_{i}, \mathcal{K}_{i}, \mathcal{P}_{i}, \mathcal{Q}_{i}$, $(i = 1, 2, ..., m), (i = 1, 2, ..., n), (i = 1, 2, ..., p), (i = 1, 2, ..., q)$, functions of bounded variation. Problem (1.6) involves a coupled system of HHFDEs, while problems (1.1)-(1.2) deals with a coupled system of sequential HHFIEs. In problems (1.1)-(1.2), there is nonlocal coupled multipoint and Hadamard fractional integrodifferential boundary conditions, whereas in problem (1.6), Stieltjes-integral boundary conditions are incorporated, involving HFDs. In problems (1.1)-(1.2), the nonlinearity depends on the unknown function and its fractional integrals at lower orders are included. Conversely, in problem (1.6), the nonlinearity depends on the unknown function, but it does not involve fractional integrals at lower orders. The researchers in ^[37] performed an examination of a coupled system of HHFDEs with nonlocal coupled HFI boundary conditions:

where $\alpha_{1}\in (1, 2]$, $\alpha_{2}\in (2, 3]$, $\beta_{1}, \beta_{2}\in [0, 1]$, $\mathcal{T} > 1$, $\delta_{1}, \delta_{2} > 0$, $\lambda_{1}, \lambda_{2} \in \mathbb{R}$, ${}^{\mathcal{HH}} \mathcal{D}^{\alpha_{i}, \beta_{j}}_{1^{+}}$ denotes the HHFD operator of order $\alpha_{i}, \beta_{j}; i = 1, 2. j = 1, 2.$ ${}^{\mathcal{H}} \mathcal{I}_{1^{+}}^{\chi}$ is the HFI operator of order $\chi \in \{ \delta_{1}, \delta_{2}\}$, and $\varrho_{1}, \varrho_{2} : \mathcal{E} \times \mathbb{R} \times \mathbb{R} \rightarrow \mathbb{R}$ are continuous functions. In Eq (1.7), we encounter a coupled system of HHFDEs, whereas Eqs (1.1) and (1.2) addresses a coupled system of sequential HHFIEs. In the latter, there are nonlocal coupled Hadamard fractional integrodifferential and multipoint boundary conditions, whereas Eq (1.7) involves multipoint and HFI boundary conditions. In Eq (1.7), solutions are obtained for the coupled system of HHFDEs, while in Eqs (1.1) and (1.2), solutions are derived for the coupled system of sequential HHFIEs. In Eqs (1.1) and (1.2), the nonlinearity depends on the unknown function and its fractional integrals at lower orders are included. Conversely, in Eq (1.7), the nonlinearity depends on the unknown function but does not involve fractional integrals at lower orders. Furthermore, in ^[38], a two-point BVP for a system of nonlinear sequential HHFDEs was investigated:

where $\alpha_{1}, \alpha_{2}\in (1, 2]$, $\beta_{1}, \beta_{2}\in [0, 1]$, $\lambda_{1}, \lambda_{2}, \mathcal{A}_{1}, \mathcal{A}_{2} \in \mathbb{R}_{+}$, $f, g : [1, e] \times \mathbb{R} \times \mathbb{R} \rightarrow \mathbb{R}$ are continuous functions. Equation (1.8) features a two-point boundary condition, while problems (1.1)-(1.2) includes multipoint and Hadamard fractional integrodifferential boundary conditions. In Eqs (1.1) and (1.2), the nonlinearity involves the unknown function and its fractional integrals at lower orders. Conversely, in Eq (1.8), the nonlinearity relies on the unknown function but does not incorporate fractional integrals at lower orders.

The document is organized as follows in the following sections: The fundamental ideas of fractional calculus relevant to this research are introduced in Section 2. An auxiliary lemma addressing the linear versions of problems (1.1) and (1.2) is provided in Section 3. The primary findings are presented in Section 4 along with illustrative examples. Finally, Section 5 provides a few recommendations.

2.   Preliminaries

Definition 2.1. The HFI of order $\mathfrak{p} > 0$ for a continuous function $\mathcal{F} : [a, \infty) \rightarrow \mathbb{R}$ is given by

where $\log(\cdot) = \log_{e} (\cdot)$.

Definition 2.2. The HFD of order $\mathfrak{p} > 0$ for a function $\mathcal{F} : [\mathfrak{a}, \infty) \rightarrow \mathbb{R}$ is defined by

where $\delta^{n} = {\varpi}^{n} \frac{d^{n}}{d\omega^{n}}$ and $[\mathfrak{p}]$ denotes the integer part of the real number $\mathfrak{p}$.

Lemma 2.3. If $\mathfrak{p}, \mathfrak{q} > 0$ and $0 < \mathfrak{a} < \mathfrak{b} < \infty$ then

In particular, for $\mathfrak{q} = 1$, we have $\Bigg({}^{\mathcal{H}} \mathcal{D}^{\mathfrak{p}}_{a^{+}} \Bigg) (1) = \frac{1}{\varGamma(1-\mathfrak{p})} \Bigg(\log \frac{\mathfrak{x}}{\mathfrak{a}}\Bigg)^{-\mathfrak{p}} \neq 0, 0 < \mathfrak{p} < 1$.

Definition 2.4. For $n-1 < \mathfrak{p} < n$ and $0 \leq \mathfrak{q}\leq 1$, the HHFD of order $\mathfrak{p}$ and $\mathfrak{q}$ for $\mathcal{F} \in \mathcal{L}'(\mathfrak{a}, \mathfrak{b})$ is defined as

where $^{\mathcal{H}}\mathcal{I}^{(\cdot)}_{\mathfrak{a}^{+}}$ and $^{\mathcal{H}}\mathcal{D}^{(\cdot)}_{\mathfrak{a}^{+}}$ are given and defined by (2.1) and (2.2), respectively.

Theorem 2.5. If $\mathcal{F} \in \mathcal{L}^{1}(\mathfrak{a}, \mathfrak{b}), 0 < \mathfrak{a} < \mathfrak{b} < \infty$, and $\Bigg({}^{\mathcal{H}} \mathcal{I}^{n-\gamma}_{a^{+}} \mathcal{F}\Bigg)({\varpi}) \in \mathcal{AC}^{n}_{\delta}[\mathfrak{a}, \mathfrak{b}]$, then

where $\mathfrak{p} > 0, 0\leq \mathfrak{q} \leq 1$ and $\gamma = \mathfrak{p} +n\mathfrak{q}- \mathfrak{p}\mathfrak{q}, n = [\mathfrak{p}]+1$. Observe that $\varGamma(\gamma-j)$ exists for all $j = 1, 2, \cdots, n-1$ and $\gamma \in [\mathfrak{p}, n]$.

We'll utilize established fixed point theorems in Banach spaces to demonstrate the existence and uniqueness of solutions for Hilfer-Hadamard fractional differential systems.

Lemma 2.6. Let $\mathcal{H}_{1}, \mathcal{H}_{2} \in \mathcal{C}([1, \mathcal{T}], \mathbb{R})$ such that

enhanced by the boundary conditions (1.2) if, and only if,

and

Proof. From the first equation of (2.3), we have

Taking both sides of the HFI of order $\alpha_{1}, \alpha_{2}$ (2.7), we obtain

Equation (2.8) can be written as follows:

and

where $\mathfrak{c}_{0}, \mathfrak{d}_{0}, \mathfrak{c}_{1}$, and $\mathfrak{d}_{1}$ are arbitrary constants. Boundary conditions (1.2) combined with (2.9) and (2.10) produce

from which we have $\mathfrak{c}_{1} = 0$ and $\mathfrak{d}_{1} = 0$. Equations (2.11) and (2.12) can be written as

from which we have

Substitute the values of $\mathfrak{c}_{0}, \mathfrak{c}_{1}, \mathfrak{d}_{0}$, and $\mathfrak{d}_{1}$ in (2.9) and (2.10), and we get solutions (2.4) and (2.5). The converse follows by direct computation. This completes the proof. □

3.   Main results

Let us introduce the Banach space $\mathcal{E} = \mathcal{T}([1, \mathcal{T}], \mathbb{R})$ endowed with the norm defined by $||\mathcal{S}||: = \sup_{{\varpi} \in [1, \mathcal{T}]} |\mathcal{S}({\varpi})|$. Thus, the product space $(\mathcal{E} \times \mathcal{E}, ||\cdot||_{\mathcal{E} \times \mathcal{E}})$ equipped with the norm $||\mathcal{S}, \mathcal{Z}||_{\mathcal{E} \times \mathcal{E}} = ||\mathcal{S}||+||\mathcal{Z}||$ for $(\mathcal{S}\times \mathcal{Z}) \in \mathcal{E} \times \mathcal{E}$ is also a Banach space.

In view of Lemma 2.6, we define as operator $\Upsilon : \mathcal{E} \times \mathcal{E} \rightarrow \mathcal{E} \times \mathcal{E}$ by

where

and

We need the following hypothesis in the sequel:

and

We give now the assumptions we will use in this section.

$(\mathcal{H}_{1})$ Assume that there exist real constants $\mathfrak{M}_{i}, \mathfrak{N}_{i} \geq 0 (i = 1, 2)$ and $\mathfrak{M}_{0} > 0, \mathfrak{N}_{0}$ such that, for all ${\varpi} \in [1, \mathcal{T}], \mathcal{S}_{i} \in \mathbb{R}, i = 1, 2, 3, 4,$

$(\mathcal{H}_{2})$ There exists positive constant $\mathcal{L}, \mathcal{L}_{1}$, such that, for all ${\varpi} \in [1, \mathcal{T}], \mathcal{S}_{i}, \mathcal{Z}_{i} \in \mathbb{R}, i = 1, 2$.

4.   Existence result via Leray-Schauder alternative

The first theorem uses the Leray-Schauder alternative to establish the existence of a result.

Lemma 4.1. ^[1] Let $\theta (\Xi) = \{ \mathcal{S} \in \mathcal{E} : \mathcal{S} = \kappa \Xi (\mathcal{S})$ for some $0 < \kappa < 1\}$, where $\Xi : \mathcal{E} \rightarrow \mathcal{E}$ is a completely continuous operator. Then, either the set $\theta(\Xi)$ is unbounded or there exists at least one fixed for the operator $\Xi$.

Theorem 4.2. Suppose condition $(\mathcal{H}{1})$ is satisfied. Additionally, assume that

Under these conditions, there exists at least one solution to the problems (1.1) and (1.2) on $\mathcal{E}$.

Proof. First, let's establish that the operator $\Upsilon : \mathcal{E} \times \mathcal{E} \rightarrow \mathcal{E} \times \mathcal{E}$, as defined in (3.1), is completely continuous. The continuity of the operator $\Upsilon$ (in terms of $\Upsilon_{1}$ and $\Upsilon_{2}$) is evident from the continuity of $\mathcal{F}$ and $\mathcal{G}$.

Next, we aim to demonstrate that the operator $\Upsilon$ is uniformly bounded. To achieve this, consider a bounded set $\mathcal{B}_{\mathfrak{r}} \subset \mathcal{E} \times \mathcal{E}$. Then, we can find positive constants $\mathcal{N}_{1}$ and $\mathcal{N}_{2}$ satisfying

Consequently, we obtain

This observation, in light of the notation (3.4) and (3.5), yields

Similarly, employing the notation (3.6) and (3.7), we obtain

Then, it follows from (4.4) and (4.5) that

This demonstrates that the operator $\Upsilon$ is uniformly bounded.

To establish the equicontinuity of $\Upsilon$, let ${\varpi}_{1}, {\varpi}_{2} \in \mathcal{E}$ with ${\varpi}_{1} < {\varpi}_{2}$. Then, we find that

independent of $(\mathcal{S}, \mathcal{Z}) \in \mathcal{B}_{\mathfrak{r}}$. Likewise, it can be shown that $|\Upsilon_{1}(\mathcal{S}, \mathcal{Z})({\varpi}_{2})-\Upsilon_{1}(\mathcal{S}, \mathcal{Z})({\varpi}_{1})| \rightarrow 0$ as ${\varpi}_{2}\rightarrow {\varpi}_{1}$ is independent of $(\mathcal{S}, \mathcal{Z}) \in \mathcal{B}_{\mathfrak{r}}$. Consequently, the equicontinuity of $\Upsilon_{1}$ and $\Upsilon_{2}$ implies the equicontinuity of the operator $\Upsilon$. Therefore, by Arzela-Ascoli's theorem, the operator $\Upsilon$ is compact. Finally, we establish the boundedness of the set $\Theta(\Upsilon) = \{\mathcal{S}, \mathcal{Z} \in \mathcal{E} \times \mathcal{E}: \mathcal{S}, \mathcal{Z} = \kappa \Upsilon(\mathcal{S}, \mathcal{Z}); 0 \leq \kappa \leq 1\}$. Let $(\mathcal{S}, \mathcal{Z}) \in \Theta (\Upsilon)$. Then, $(\mathcal{S}, \mathcal{Z}) = \kappa \Upsilon (\mathcal{S}, \mathcal{Z})$, which implies that $\mathcal{S}({\varpi}) = \kappa\Upsilon_{1} (\mathcal{S}, \mathcal{Z})({\varpi}), \mathcal{Z}({\varpi}) = \kappa\Upsilon_{2} (\mathcal{S}, \mathcal{Z})({\varpi})$ for any ${\varpi} \in \mathcal{E}$, and so $|\mathcal{S}({\varpi})| = |\Upsilon_{1} (\mathcal{S}, \mathcal{Z})({\varpi})|, \mathcal{Z}({\varpi}) = |\Upsilon_{2} (\mathcal{S}, \mathcal{Z})({\varpi})| \ \forall \ {\varpi} \in [0, 1]$. Using some inequality proved at the beginning of the proof, we find

which implies that

Similarly, one can find that

From (4.9) and (4.10), we obtained

where $\Psi_{i}$ $i = 0, 1, 2$ are given by (3.8). By (4.1), we deduce that

As a result, $\Theta(\Upsilon)$ is bounded. Consequently, the conclusion of Lemma 2.6 applies, implying that the operator $\Upsilon$ has at least one fixed point, which indeed serves as a solution to the problems (1.1) and (1.2). □

Now, we introduce the constants

where $\Omega_{1}, \Omega_{2}, \bar{\Omega_{1}}, \bar{\Omega_{2}}$ are given by (3.4)–(3.7).

The subsequent result will establish the existence of a unique solution to the problems (1.1) and (1.2) through the application of a fixed point theorem attributed to Banach.

Theorem 4.3. If the assumption ($\mathcal{H}_{2}$) is satisfied and that

where $\Omega_{i}$ and $\bar{\Omega}_{i}$, (i = 1, 2) are given in (3.4)–(3.7), then the problems (1.1) and (1.2) have a unique solution on $\mathcal{E}$.

Proof. By using condition (4.15), we define the positive number

where $\mathscr{G}_{1}, \mathscr{G}_{2}, \mathscr{D}_{1}, \mathscr{D}_{2}$ are given by (4.12). We will prove that $\mathcal{A}(\mathcal{B}_{\mathcal{R}}) \subset \mathcal{B}_{\mathcal{R}}$, where $\mathcal{B}_{\mathcal{R}} = \{(\mathcal{S}, \mathcal{Z}) \in \mathcal{E} \times \mathcal{E} : ||(\mathcal{S}, \mathcal{Z}) \leq {\mathcal{R}}\}$. For $(\mathcal{S}, \mathcal{Z}) \in \mathcal{B}_{\mathcal{R}}$ and ${\varpi} \in [0, 1]$, we obtain

and

Then we deduce that

and

Therefore, by (4.17), (4.18), and the definition of $\mathcal{R}$, we conclude that

which gives us $\Upsilon(\mathcal{B}_{\mathcal{R}})\subset \mathcal{B}_{\mathcal{R}}.$

We will prove next that $\Upsilon$ is a contraction operator. By using ($\mathcal{H}_{2}$), for $(\mathcal{S}_{i}, \mathcal{Z}_{i}) \in \mathcal{B}_{\mathcal{R}}, i = 1, 2,$ and for any ${\varpi} \in [0, 1]$, we find:

Letting $\mathcal{K}_{1} = \sup_{{\varpi} \in [1, \mathcal{T}]} |\mathcal{F}({\varpi}, 0, 0)| < \infty$ and $\mathcal{K}_{2} = \sup_{{\varpi} \in [1, \mathcal{T}]} |\mathcal{G}({\varpi}, 0, 0)| < \infty$, it follows by the assumption ($\mathcal{H}_{1}$) that

and

To begin, we show that $\Upsilon\mathcal{B}_{\rho} \subset \mathcal{B}_{\rho}$, where $\mathcal{B}_{\rho} = \{(\mathcal{S}, \mathcal{Z}) \in \mathcal{E} \times \mathcal{E} : ||(\mathcal{S}, \mathcal{Z}) \leq {\rho}\}$, with

For $(\mathcal{S}, \mathcal{Z}) \in \mathcal{B}_{\rho}$, we have

which yields

Using the notation (3.4)-(3.5), we get

Likewise, we can find that

Then, it follows from (4.21)-(4.22) that

Therefore, $\Upsilon \mathcal{B}_{\rho} \subset \mathcal{B}_{\rho}$ as $(\mathcal{S}, \mathcal{Z}) \in \mathcal{B}_{\rho}$ is an arbitrary element.

In order to verify that the operator $\Upsilon$ is a contraction, let $\mathcal{S}_{i}, \mathcal{Z}_{i} \in \mathcal{B}_{\rho}, i = 1, 2.$ Then, we get

which, by $(\mathcal{H}_{2})$, yields

Similarly, we can observe that

Consequently, it follows from (4.23) and (4.24) that

and by condition (4.15), it follows that $\Upsilon$ is a contraction. Consequently, the operator $\Upsilon$ possesses a unique fixed point as a direct application of the Banach fixed point theorem. Thus, there exists a unique solution for the problems (1.1) and (1.2) on $\mathcal{E}$. □

5.   Hyers-Ulam stability of system

This section is devoted to the investigation of Hyers-Ulam stability for our proposed system. Consider the following inequality:

where $\varepsilon_{1}, \varepsilon_{2}$ are given two positive real numbers.

Definition 5.1. Problem (1.1) is Hyers-Ulam stable if there exist $\Omega_{i} > 0, i = 1, 2, 3, 4$ such that for a given $\varepsilon_{1}, \varepsilon_{2} > 0$ and for each solution $(\mathcal{S}, \mathcal{Z}) \in \mathcal{C}([1, \mathcal{T}], \times \mathbb{R}^{2})$ of inequality (5.1), there exists a solution $(\mathcal{S}^{*}, \mathcal{Z}^{*}) \in \mathcal{C}([1, \mathcal{T}], \times \mathbb{R}^{2})$ of problem (1.1) with

Remark 5.1. $(\mathcal{S}, \mathcal{Z})$ is a solution of inequality (5.1) if there exist functions $\mathcal{Q}_{i} \in \mathcal{C}([1, \mathcal{T}], \mathbb{R}), i = 1, 2$, which depend upon ${\mathcal{S}, \mathcal{Z}}$, respectively, such that

Remark 5.2. If ($\mathcal{S}, \mathcal{Z})$, respectively, is a solution of inequality (5.1), then $(\mathcal{S}, \mathcal{Z})$ is a solution of the following inequality:

As from Remark 5.1, we have

With the help of Definition 5.1 and Remark 5.1, we verified Remark 5.2 in the following lines:

By the same method, we can obtain that

where $\Omega_{1}, \Omega_{2}, \bar{\Omega_{1}}, \bar{\Omega_{2}}$ are given by (3.4)–(3.7). Hence, Remark 5.2 is verified, with the help of (5.7) and (5.8). Thus, the nonlinear sequential coupled system of HHFDEs is Hyers-Ulam stable and, consequently, the system (1.1) is Hyers-Ulam stable.

6.   Examples

Consider the following Hilfer-Hadamard fractional BVP:

with $\alpha_{1} = 5/4$, $\alpha_{2} = 3/2$, $\beta_{1} = 1/2$, $\beta_{2} = 1/4$, $\mathfrak{m} = 2$, $\mathfrak{n} = 2$, $\mathfrak{r} = 2$, $\mathfrak{a} = 2$, $\mathfrak{b} = 2$, $\mathfrak{c} = 2$, $\eta_{1} = 1/5$, $\eta_{2} = 1/10$, $\xi_{1} = 4/3$, $\xi_{2} = 3/2$, $\phi_{1} = 5/3$, $\phi_{2} = 7/3$, $\theta_{1} = 1/2$, $\theta_{2} = 1/2$, $\zeta_{1} = 7/3$, $\zeta_{2} = 5/2$, $\lambda_{1} = 1$, $\lambda_{2} = 1$, $\omega_{1} = 1/4$, $\omega_{2} = 2/3$, $\mu_{1} = 4$, $\mu_{2} = 4/3$, $\mathcal{P}_{1} = 1/18$, $\mathcal{P}_{2} = 1/9$, $\psi_{1} = 4/3$, $\psi_{2} = 5/2$, $\mathcal{Q}_{1} = 1/4$, $\mathcal{Q}_{2} = 1/7$, $\delta_{1} = 1/4$, $\delta_{2} = 3/5$, $\sigma_{1} = 5/3$, $\sigma_{2} = 3/2$, $\mathcal{M}_{1} = 2/3$, $\mathcal{M}_{2} = 2/5$, $\vartheta_{1} = 2/3$, $\vartheta_{2} = 3/5$, $\pi_{1} = 5/2$, $\pi_{2} = 3/2.$ Using the given data, it is found that $\gamma_{1} = 13/8, \gamma_{2} = 13/8$, $\Delta = 0.639100490745, \mathcal{A}_{1} = 0.799441,$$\mathcal{A}_{2} = 0.799441, \mathcal{B}_{1} = 0.1184655, \mathcal{B}_{2} = 0.1251315,$$\Omega_{1} = 1.3283929, \Omega_{2} = 0.718823345, \bar{\Omega_{1}} = 1.028734, \bar{\Omega_{2}} = 0.97432874, \mathcal{T} = 2,$$\mathfrak{p}_{1} = 11/5, \mathfrak{p}_{2} = 25/6, \mathfrak{q}_{1} = 11/5,$$\mathfrak{q}_{2} = 22/7.$

Example 6.1. For illustrating Theorem 4.2, we take

for all ${\varpi} \in [0, 1], \mathcal{S}_{i} \in \mathbb{R}, i = 1, 2, 3, 4.$

We obtained the inequalities

for all ${\varpi} \in [0, 1]$ and $\mathcal{S}_{i} \in \mathbb{R}, i = 1, 2, 3, 4.$ We also have $\mathfrak{M}_{0} = \frac{1}{2}, \mathfrak{M}_{1} = \frac{1}{16}, \mathfrak{M}_{2} = \frac{1}{16}, \mathfrak{M}_{3} = \frac{1}{18}$, $\mathfrak{M}_{4} = \frac{1}{9}, \mathfrak{N}_{0} = \frac{7}{9}, \mathfrak{N}_{1} = \frac{1}{27}$, $\mathfrak{N}_{2} = \frac{4}{27}, \mathfrak{N}_{3} = \frac{4}{5}, \mathfrak{N}_{4} = \frac{4}{5},$

We find here $\Psi_{1} \approx 0.8228588$ and $\Psi_{2} \approx 0.5948088559$. We deduce that the condition $\max\{\Psi_{1}, \Psi_{2}\} = \Psi_{1} < 1$ is satisfied. Then, by Theorem 4.2, we conclude that the problem (6.1) with the nonlinearities (6.2) has at least one solution ${\varpi} \in [0, 1]$.

Example 6.2. For illustrating Theorem 4.3, we take

for all ${\varpi} \in [0, 1], \mathcal{S}_{i} \in \mathbb{R}, i = 1, 2, 3, 4.$

We obtain here the following inequalities

for all ${\varpi} \in [0, 1]$. So, we have $\mathfrak{c}_{0} = 1/4$ and $\mathfrak{d}_{0} = 1/20$. In addition, we find $\rho_{1}\approx 1.4125480, \rho_{2} \approx 1.13889158$, $\mathscr{D}_{1} \approx 0.510067622, \mathscr{D}_{2}\approx 0.43933889$. Then, $\mathscr{D}_{1}+ \mathscr{D}_{2} \approx 0.9410020109 < 1$, that is, the condition (4.15) is satisfied. Therefore, by Theorem 4.3, we conclude that problem (6.1) with the nonlinearities (6.4) has a unique solution ${\varpi} \in [0, 1]$.

7.   Discussion & conclusions

We have presented criteria for the existence, uniqueness, and Ulam-Hyers stability of solutions to a coupled system of nonlinear sequential HHFIEs and nonlocal coupled Hadamard fractional integrodifferential and multipoint boundary conditions. We derive the expected results using a methodology that uses modern analytical tools. It is imperative to emphasize that the results offered in this specific context are novel and contribute to the corpus of existing literature on the topic. Furthermore, our results encompass cases where the system reduces to the boundary conditions of the form:

When $\eta_{\mathfrak{j}} = \mathcal{P}_{\mathfrak{u}} = 0$, then we get

If $\theta_{\mathfrak{i}} = \mathcal{Q}_{\mathfrak{v}} = 0$, we get:

When $\lambda_{\mathfrak{k}} = \mathcal{M}_{\mathfrak{w}} = 0$, the outcome is:

In addition, if $\eta_{\mathfrak{j}} = \mathcal{P}_{\mathfrak{u}} = \lambda_{\mathfrak{k}} = \mathcal{M}_{\mathfrak{w}} = 0$, we obtain:

When $\eta_{\mathfrak{j}} = \mathcal{P}_{\mathfrak{u}} = \theta_{\mathfrak{i}} = \mathcal{Q}_{\mathfrak{v}} = 0$, the boundary condition is:

If $\lambda_{\mathfrak{k}} = \mathcal{M}_{\mathfrak{w}} = \theta_{\mathfrak{i}} = \mathcal{Q}_{\mathfrak{v}} = 0$, we obtain:

These cases represent new findings. Looking ahead, our future plans include extending this work to a coupled system of nonlinear sequential HHFIEs enhanced by the nonlocal coupled mixed integro-differential and discrete type boundary conditions. We also intend to investigate the multivalued analogue of the problem studied in this paper.

Author contributions

Subramanian Muthaiah: Developed the conceptualization and proposed the method, wrote–original draft, reviewed and edited the paper; Manigandan Murugesan: Developed the conceptualization and proposed the method, wrote–original draft, investigated, processed and provided examples; Muath Awadalla: Investigated, processed and provided examples, reviewed and edited the paper; Bundit Unyong: Developed the conceptualization and proposed the method, reviewed and edited the paper; Ria H Egami: Investigated, processed and provided examples, reviewed and edited the paper. All authors have read and agreed to the published version of the manuscript.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

This work was supported by the Deanship of Scientific Research, Vice Presidency for Graduate Studies and Scientific Research, King Faisal University, Saudi Arabia [GrantA113]. This study is supported via funding from Prince Sattam bin Abdulaziz University, project number (PSAU/2024/R/1445). M.Manigandan gratefully acknowledges the Center for Computational Modeling, Chennai Institute of Technology, India, vide funding number CIT/CCM/2023/RP-018.

Conflict of interest

The authors declare no conflicts of interest.