Solving Fuzzy Volterra Fractional Integro-Differential Equations with a New Least Squares Technique

This part provides an overview of the essential definitions and notations of fractional calculus that are necessary to understand the content in the next sections. Definitions 2.1 and 2.2, extracted from references [14, 22], provide precise explanations of the fractional derivative and integral. Additionally, the property that the Caputo derivative of a constant is equal to zero makes them commonly used, especially in the representation of physical and engineering problems. In addition, the Caputo fractional derivative enables the use of beginning or boundary conditions in the problem formulation. Also, this section will provide precise definitions for important terminology in fuzzy set theory, including fuzzy functions and fuzzy integers. The range of numbers is from references [23, 24].

Definition 2.1 [14, 24] (Riemann-Liouville): The R-L fractional integral of order $\alpha>0$ is define as:

$\mathrm{I}_{\mathrm{t}}^\alpha \mathfrak{f}(\mathrm{t})=\frac{1}{\Gamma(\alpha)} \int_0^{\mathrm{t}}(t-\tau)^{\alpha-1} \mathfrak{f}(t) \mathrm{d} \tau, t>0, \quad \alpha \in \mathbb{R}^{+}$

where, $\Gamma$ represents the Gamma function.

Definition 2.2 [11, 24] (Caputo fractional order derivative): The Caputo fractional order derivative of order $\alpha>0$ is define as follows:

${ }^{\mathrm{c}} \mathfrak{D}^\alpha \mathfrak{f}(\mathrm{t})=\left\{\begin{array}{c}\frac{1}{\Gamma(\mathrm{~m}-\alpha)} \int_0^{\mathrm{t}} \frac{\mathrm{f}^{(\mathrm{m})}(\tau)}{(t-\tau)^{\alpha+1-\mathrm{m}}} \mathrm{d} \tau, \mathrm{m}-1<\alpha<\mathrm{m} \\ \frac{\mathrm{d}^{\mathrm{m}}}{\mathrm{dt}^{\mathrm{m}}} \mathfrak{f}(t), \alpha=\mathrm{m}\end{array}\right.$

For $\alpha>0$, the Caputo fractional derivative has the following properties:

(1) ${ }^{\mathrm{c}} \mathfrak{D}^\alpha\left(\mathrm{I}_{\mathrm{t}}^\alpha \mathfrak{f}(t)\right)=\mathrm{f}(t)$

(2) $\mathrm{I}_{\mathrm{t}}^\alpha\left({ }^{\mathrm{c}} \mathfrak{D}^\alpha \mathfrak{f}(t)\right)=\mathrm{f}(\mathrm{t})-\sum_{\mathrm{k}=0}^{\mathrm{n}-1} \mathfrak{f}^{(\mathrm{k})}\left(0^{+}\right) \frac{\mathrm{x}^{\mathrm{k}}}{\mathrm{k}!}$

(3) ${ }^{\mathrm{c}} \mathfrak{D}^\alpha(\mathrm{C})=0, \mathrm{C} \in \mathbb{R}$

(4) $I_{\mathrm{t}}^\alpha t^{\mathfrak{B}}=\frac{\Gamma(\mathfrak{B}+1)}{\Gamma(\alpha+\mathfrak{B}+1)} t^{+\alpha}$

(5) ${ }^{\mathrm{c}} \mathfrak{D}^\alpha\left(\mathrm{t}^{\mathrm{k}}\right)=\left\{\begin{array}{c}\frac{\Gamma(\mathrm{k}+1)}{\Gamma(\mathrm{k}-\alpha+1)} \mathrm{t}^{\mathrm{k}-\alpha}, \mathrm{k} \in\{1,2,3, \ldots\}, \mathrm{k} \geq[\alpha] \\ 0, \mathrm{k} \in\{1,2,3, \ldots\}, \mathrm{k}<[\alpha]\end{array}\right.$

where, $[\alpha]$ is the floor function of $\alpha$.

Definition 2.3 [25] Let X be any of elements, a fuzzy set $\widetilde{\mathfrak{B}}$ is characterized by a membership function $\mu_{\mathfrak{B}}(x): \mathrm{X} \rightarrow[0,1]$, and may be written as the set of points:

$\widetilde{\mathfrak{B}}=\left\{\left(x, \mu_{\mathfrak{B}}(\mathrm{x})\right) x \in \mathrm{X}, 0 \leq \mu_{\mathfrak{B}}(x) \leq 1\right\}$

Definition 2.4 [26]. The crisp set of components that make up the set $\widetilde{\mathfrak{B}}$ at least to the degree $\alpha$ is called the weak $\alpha$-level set (or weak $\alpha$-cut), and is defined by $\mathfrak{B}_b=\{x \in \mathrm{X}: \left.\mu_{\mathfrak{B}}(\mathrm{x}) \geq \alpha\right\}$.

While the strong $\alpha$-level set (or strong $\alpha$-cut) is defined by: $\mathfrak{B}_\alpha^{\prime}=\left\{x \in \mathrm{X}: \mu_{\mathfrak{B}}(x)>\alpha\right\}$.

Definition 2.5 [25, 27] (Triangular Fuzzy Number)

A fuzzy number with the following function of membership:

$\mu(x ; c, b, a)=\left\{\begin{array}{rl}0, & x<c \\ \frac{x-c}{b-c}, & c \leq x \leq b \\ \frac{a-x}{a-b}, & b \leq x \leq a \\ 0, & x>a\end{array}\right.$,

which is indicated by the fuzzy number of triangular kinds, where $a, b$ and $c \in \mathbb{R}$ with $a>b>c$.

In Figure 1, we can observe that the graph illustrates the symmetry of the membership function $\mu(x)$ takes the triangular shape with the base over the interval $[c, b]$ and vertex at $x=b$.

ping_mu_jie_tu_2025-08-20_113234.png

Figure 1. Triangular fuzzy number

The related $\alpha$-cut is $[\tilde{\mu}]_\alpha=[c+\alpha(b-c), a-\alpha(a- b)], \alpha \in[0,1]$.

Definition 2.6 [28, 29]: Given two fuzzy numbers (two fuzzy sets) $\widetilde{a}, \tilde{b} \in \widetilde{U}$, then the Hausdorff distance between $\widetilde{a}$ and $\tilde{b}$ symbolized by $D_H\left([\tilde{a}, \tilde{b}]_\alpha\right)$ is such that $D([\tilde{a}, \tilde{b}])= \sup \left\{D_H\left([\tilde{a}, \tilde{b}]_\alpha\right) \mid \alpha \in[0,1]\right\}$, where ( $\widetilde{U}, D$ ) indicates the Cauchy space (complete metric space).

Suppose that $\widetilde{U}$ is the set of all semi-continuous upper normal convex fuzzy numbers with bounded $\alpha$-level sets. The $\alpha$-cuts of fuzzy numbers are always bounded and closed such that the intervals are $[\tilde{\mu}(x)]_\alpha=[\underline{\mu}(x), \bar{\mu}(x)]_\alpha, x \in \mathbb{R}, \forall \alpha \in [0,1]$.

Let $[\tilde{b}(\alpha)]=[\underline{b}(\alpha), \bar{b}(\alpha)],[\tilde{a}(\alpha)]=[\underline{a}(\alpha), \bar{a}(\alpha)]$ represent two fuzzy numbers [Definition 2.5], for $\gamma \geq 0$. We can define the binary operations on fuzzy numbers $\tilde{b}$ and $\tilde{a}$ of addition and scalar multiplication described as follows:

1. $(\underline{b+a})(\alpha)=(\underline{b}(\alpha)+\underline{a}(\alpha))$

2. $(\overline{b+a})(\alpha)=(\bar{b}(\alpha)+\bar{a}(\alpha))$

3. $(\underline{\gamma b})(\alpha)=\gamma \cdot \bar{b}(\alpha)$

Definition 2.7 [29]: Let $\mathrm{f} \in \mathrm{C}(\mathrm{J} ; \mathrm{E}) \backslash \mathrm{L}^1(\mathrm{~J} ; \mathrm{E})$ exhibit characteristics of a fuzzy set-valued function. Function $f$ is considered Caputo's fuzzy differentiable at point $x$ when

${ }^c \mathfrak{D}^\alpha \mathfrak{f}(t)=\frac{1}{\Gamma(1-\alpha)} \int_0^t \frac{f^{(m)}(\tau)}{(t-\tau)^\alpha} d \tau, 0<\alpha \leq m$

Definition 2.8 [30]: Let $\tilde{f}(x)$ be a fuzzy-valued approximate solution to a fuzzy differential equation, and let $\widetilde{R}\left(\widetilde{a}_i, x, \alpha\right)$ expressed using $\alpha$-cuts, is defined as:

$[\tilde{R}(\alpha)]=[[\underline{b}(t, \alpha)]-[f(\alpha)]]_{\boldsymbol{\alpha}} \mid=[\underline{R}(t, \alpha), \bar{R}(t, \alpha)]$

In this section, we will present a new form of standard LST in order to solve FV-FIDEs, consider the following F-FIDE of order $\alpha \in[0,1]$ for $w \in[0,1]$:

$\begin{gathered}\mathrm{c} \widetilde{\mathfrak{D}}^\alpha \tilde{y}(x ; w)=\tilde{\mathfrak{f}}(x ; w) +\int_0^t k(x, t) F(\tilde{y}(t ; w)) d t\end{gathered}$                      (1)

$0 \leq x, t \leq 1,0<\alpha \leq 1$

Under fuzzy initial conditions:

$\tilde{y}(0 ; w)=\tilde{y}_0(w)$                       (2)

where, Eq. (1) can be restated with a broader generalization by applying the concepts of defuzzification, where $\tilde{\mathfrak{f}}(x ; w)=[\overline{\mathfrak{f}}, \underline{\mathfrak{f}}]^\omega$ is a fuzzy function hence the approximate solution will be a fuzzy function which can be written as $\tilde{y}(t ; w)= [\bar{y}, \underline{y}]^\omega$, where, ${ }^{\mathrm{c}} \widetilde{\mathfrak{D}}^\alpha$ is the fuzzy Caputo fractional derivative of order $\alpha$ is a fuzzy function of the crisp variable $x$, and $\tilde{y}(x ; w)$ is a fuzzy function. Thus, Eq. (1) can be written in terms level sets $w \in[0,1]$ as:

$\begin{gathered}{ }^c \underline{\mathfrak{D}}^\alpha \underline{y}(x ; w)=\underline{\mathfrak{f}}(x ; w)+\int_0^t k(x, t) \underline{F}(\underline{y}(t ; w)) d t, \\ 0 \leq x, t \leq 1,0<\alpha \leq 1\end{gathered}$                    (3)

$\begin{aligned} c \overline{\mathfrak{D}}^\alpha \bar{y}(x ; w) & =\overline{\mathfrak{f}}(v ; w)+\int_0^t k(x, t) \bar{F}(\bar{y}(t ; w)) d t,\\ & 0 \leq x, t \leq 1,0<\alpha \leq 1 \\ \bar{y}(0 ; w) & =\bar{y}_0(w) \text { and } \underline{y}(0 ; w)=\underline{y}_0(w)\end{aligned}$

From Section 3 application, the approximate solution of solve $\tilde{y}(x ; w)$ in terms of $\alpha$ and $\omega$ in the following LST polynomial approximation solution of LST such that

$\tilde{y}(x ; w)=\sum_{i=0}^n\left[\tilde{a}_i\right]^\omega \tilde{P}_{l_{i+1}}(t ; w)$

then substitute the values of $\left[\tilde{a}_i\right]^\omega=\left[\bar{a}_i, \underline{a}_i\right]^\omega$ into Eq. (3) yields:

$\begin{aligned} & \sum_{i=0}^n \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w)=I\left(\underline{f}(x ; w)+\int_0^t k(x, t) \underline{F}\left(\sum_{i=0}^n \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t\right) \\ & \sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)=I\left(\bar{f}(x ; w)+\int_0^t k(x, t) \bar{F}\left(\sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t\right)\end{aligned}$                     (4)

Next, Eq. (4) can be expanded into the following forms

$\begin{aligned} & \sum_{i=0}^n \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w)-I\left(\underline{f}(x ; w)+\int_0^t k(x, t) \underline{F}\left(\sum_{i=0}^n \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t\right)=0 \\ & \sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)-I\left(\bar{f}(x ; w)+\int_0^t k(x, t) \bar{F}\left(\sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t\right)=0\end{aligned}$                     (5)

Then, by minimizing Eq. (5) such that

$\begin{gathered}\underline{E}\left(\left[\underline{a}_i\right]^w\right)=\binom{\sum_{i=0}^n \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w)-I(\underline{\mathfrak{f}}(x ; w)+}{\int_0^t k(x, t) \underline{F}\left(\sum_{i=0}^n \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t}^2 \\ \bar{E}\left(\left[\bar{a}_i\right]^w\right)=\binom{\sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)-I(\overline{\mathfrak{f}}(x ; w)+}{\int_0^t k(x, t) \bar{F}\left(\sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t}^2\end{gathered}$                    (6)

where, $\tilde{E}\left(\left[\tilde{a}_i\right]^\omega\right)$ is minimizing requires finding the values of $\left[\tilde{a}_i\right]^w$, to determining the best approximate solution of Eq. (1). To get the minimum possible value of $\tilde{E}\left(\left[\tilde{a}_i\right]^w\right)$, then the following partial derivatives with respect to ${ }^c \widetilde{\mathrm{D}}^\alpha \tilde{y}(x ; w)$ as follows:

$\begin{aligned} & \frac{\partial \underline{E}\left(\left[\underline{a}_i\right]^w\right)}{\partial\left[\underline{a}_i\right]^w}=0, i=0,1, \ldots, n \\ & \frac{\partial \bar{E}\left(\left[\bar{a}_i\right]^w\right)}{\partial\left[\bar{a}_i\right]^w}=0, i=0,1, \ldots, n\end{aligned}$                       (7)

Such that:

$\begin{gathered}\frac{\partial \underline{E}\left(\left[\underline{a}_i\right]^w\right)}{\partial\left[\underline{a}_i\right]^w}=2\left(\sum_{i=0}^n \underline{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)-I\left(\overline{\mathfrak{f}}(x ; w)+\int_0^t k(x, t) \bar{F}\left(\sum_{i=0}^n \underline{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t\right)\right) \\ \frac{\partial \bar{E}\left(\left[\bar{a}_i\right]^w\right)}{\partial\left[\bar{a}_i\right]^w}=2\left(\sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)-I\left(\overline{\mathfrak{f}}(x ; w)+\int_0^t k(x, t) \bar{F}\left(\sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t\right)\right)\end{gathered}$                 (8)

$\frac{\partial \underline{E}\left(\left[\underline{a}_i\right]^w\right)}{\partial\left[\underline{a}_i\right]^w}$ and $\frac{\partial \bar{E}\left(\left[\bar{a}_i\right]^w\right)}{\partial\left[\bar{a}_i\right]^w}$ for all i, causing $\tilde{E}\left(\left[\tilde{a}_i\right]^w\right)$ to be reduced with each i. By rearranging Eq. (8), we can derive the following normal equations, to find coefficients $\left[\tilde{a}_i\right]^w$ such that

$\begin{gathered}\frac{\partial E\left(\left[\underline{a}_i\right]^w\right)}{\partial\left[\underline{a}_i\right]^w}=\sum_{i=0}^n \underline{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)-I\left(\overline{\mathfrak{f}}(x ; w)+\int_0^t k(x, t) \bar{F}\left(\sum_{i=0}^n \underline{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t\right) \\ \frac{\partial \bar{E}\left(\left[\bar{a}_i\right]^w\right)}{\partial\left[\bar{a}_i\right]^w}=\sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)-I\left(\overline{\mathfrak{f}}(x ; w)+\int_0^t k(x, t) \bar{F}\left(\sum_{i=0}^n \bar{a}_i \overline{\mathcal{P}}_{l_{i+1}}(t ; w)\right) d t\right)\end{gathered}$                       (9)

Thus, the derivative of the sum represented by Eq. (9) with respect to each of the coefficients $\left[\tilde{a}_i\right]^w$. The F-FLST polynomial approximation of Eq. (9) is defined as $\widetilde{\mathrm{X}}\left[\tilde{a}_i\right]^w= f(x ; w)$.

We will use the methodology from the previous section to test the effectiveness of the FV-FIDEs. Each example showcases the efficacy of the proposed method by contrasting the estimated results with the precise solutions.

Example 6.1 [35]: Given the linear V-FFIDEs

$\begin{gathered}\mathfrak{D}^\alpha \tilde{y}(x ; w)=[w-1,1-w] +\int_0^x \tilde{y}(t ; w) d t, 0<w \leq 1,0 \leq x, t \leq 1\end{gathered}$                      (13)

with fuzzy initial condition $\tilde{y}(0, w)=0$.

${ }^c \mathfrak{D}^\alpha$ represents the $\alpha$-th Caputo fractional derivative of y(x) and $0<\alpha \leq 1,[\omega-1,1-\omega]$ are fuzzy triangular numbers, $y\left(x_0\right)$ represents the fuzzy initial condition, and $y$ ($x$) represent a crisp function of non-fuzzy independent variable $x$. The exact solution of Eq. (13) when $\alpha$ = 1:

$[y(x)]^w=[w-1,1-\omega] \sinh (x)$                         (14)

By utilizing the FVLST algorithm as presented previously, the $k^{t h}$ order deformation equation of Eq. (11) and operating $I^\alpha$ on the both sides of Eq. (13) and plugging property (4) of defined (2).

$\begin{gathered}\underline{y}(x ; w)=I^\alpha\left((w-1)+\int_0^x \underline{y}(t ; w) d t\right), \quad 0<\alpha \leq 1,0 \leq x, t \leq 1 \\ \bar{y}(x ; w)=I^\alpha\left((1-w)+\int_0^x \bar{y}(t ; w) d t\right), \quad 0<\alpha \leq 1,0 \leq x, t \leq 1\end{gathered}$                      (15)

By applying the Shifted Legendre Polynomials from section (3), yields.

$\left\{\begin{array}{l}\sum_{i=0}^2 \underline{a_i} \underline{\mathcal{P}}_{l_{i+1}}(x ; w)-I^\alpha\left((w-1)+\int_0^x \sum_{i=0}^2 \underline{a_i} \underline{\mathcal{P}}_{l_{i+1}}(t ; w) d t\right)=0 \\ \sum_{i=0}^2 \overline{b_i} \overline{\mathcal{P}}_{l_{i+1}}(x ; w)-I^\alpha\left((1-w)+\int_0^x \sum_{i=0}^2 \overline{b_i} \overline{\mathcal{P}}_{l_{i+1}}(t ; w) d t\right)=0\end{array}\right.$                         (16)

Thus, we minimized Eq. (16) as

$\begin{aligned} & E 1\left(\underline{a}_0, \underline{a}_1, \underline{a}_2\right)^w=\sum_{i=0}^2 \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(x ; w)-I^\alpha\left((w-1)+\int_0^x \sum_{i=0}^2 \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w) d t\right) \\ & E 2\left(\overline{b_0}, \overline{b_1}, \overline{b_2}\right)^w=\sum_{i=0}^2 \overline{b_i} \overline{\mathcal{P}}_{l_{i+1}}(x ; w)-I^\alpha\left((1-w)+\int_0^x \sum_{i=0}^2 \overline{b_i} \overline{\mathcal{P}}_{l_{i+1}}(t ; w) d t\right)\end{aligned}$                     (17)

Therefore, minimizing E required for finding the values of $a_i$, as well as to determining the best estimate for the solution of the fuzzy fractional integro differential Eq. (17). To get the minimum possible value of E, then the following partial derivatives with respect to ${ }^c D^\alpha y(x)$ and $0<\alpha \leq 1$ for $i$ =0,1,2.

$\left\{\begin{array}{l}\left(\frac{\partial E}{\partial \underline{a_i}}\right)^w=0 \\ \left(\frac{\partial E}{\partial \overline{b_i}}\right)^w=0\end{array}\right.$                  (18)

$\begin{aligned} & E 1\left(\underline{a}_0, \underline{a}_1, \underline{a}_2\right)^w=\binom{\sum_{i=0}^2 a_i \underline{\mathcal{P}}_{l_{i+1}}(x ; w)}{-I^\alpha\left((w-1)+\int_0^x \sum_{i=0}^2 \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w) d t\right)}^2 \\ & E 2\left(\overline{b_0}, \overline{b_1}, \overline{b_2}\right)^w=\binom{\sum_{i=0}^2 \overline{b_i} \overline{\mathcal{P}}_{l_{i+1}}(x ; w)}{-I^\alpha\left((1-w)+\int_0^x \sum_{i=0}^2 \overline{b_i} \overline{\mathcal{P}}_{l_{i+1}}(t ; w) d t\right)}^2\end{aligned}$                       (19)

Applying Eq. (18) in Eq. (19) we get:

$\begin{aligned} & \left(\frac{\partial E 1}{\partial\left(\underline{a}_0, \underline{a}_1, \underline{a}_2\right)}\right)^w=2\left[\sum_{i=0}^2 \underline{a_i} \underline{\mathcal{P}}_{l_{i+1}}(x ; w)-I^\alpha\left((w-1)+\int_0^x \sum_{i=0}^2 \underline{a}_i \underline{\mathcal{P}}_{l_{i+1}}(t ; w) d t\right)\right] \\ & \left(\frac{\partial E 2}{\partial\left(\overline{b_0}, \overline{b_1}, \overline{b_2}\right)}\right)^w=2\left[\sum_{i=0}^2 \overline{b_i} \overline{\mathcal{P}}_{l_{i+1}}(x ; w)-I^\alpha\left((1-w)+\int_0^x \sum_{i=0}^2 \overline{b_i} \overline{\mathcal{P}}_{l_{i+1}}(t ; w) d t\right)\right]\end{aligned}$                      (20)

$\frac{\partial E 1}{\partial a_i}, \frac{\partial E 2}{\partial b_i}$ for all i=0, 1, 2, causing E to be reduced. with each i=0, 1, 2. By rearranging Eq. (20), one can derive the following normal equations, to find fuzzy coefficients $\tilde{a}_0, \tilde{a}_1, \tilde{a}_2$ and $\tilde{b}_0, \tilde{b}_1, \tilde{b}_2$ using Eq. (11) to verify the accuracy of our proposed method as shown in the tables below.

Table 1 displays the fuzzy approximation solutions derived from the Laplace–Sumudu Transform (LST) method for fractional orders $\alpha$ = 0.5, 0.7, and 0.9, as applied to Eq. (13).

Table 1. Numerical results for Example 1, with different values of $\alpha$ at $w=0.5,0.75$

The preliminary estimate employed for this equation is $\tilde{y}_0(x)=[\underline{y}(x), \bar{y}(x)]$, as elaborated in Section 5. The values of $\alpha$ are derived from Eq. (15), and by adjusting $I^\alpha$ for $\alpha \in[0,1]$, the convergence domain of the solution is investigated. To assess the convergence behavior, we quantitatively examined the absolute and relative errors between the approximate and exact fuzzy solutions over the domain. For example, for $x$ = 1, the average absolute error across $\alpha$-levels was determined to be less than $10^{-3}$, signifying exceptional precision. The relative error consistently stayed below 1% for all evaluated values of $\alpha$, with the lowest error recorded at $\alpha$ = 0.7, indicating optimal convergence at this parameter. Figure 2 depicts the lower and upper bound solutions $\bar{y}(x)$ and $\underline{y}(x)$ for the eighth-order fuzzy Volterra linear fractional integro-differential equation associated with Eq. (15), displayed at $w$ = 1. The proximity and consistency of these bounds across various $\alpha$-levels validate the method's robustness and define the effective convergence zone of the LST methodology. The results indicate that the LST method yields dependable and precise approximations for FFIDEs when suitable fractional orders are used. The precision is enhanced with smooth kernel functions and uniform $\alpha$-level decomposition.

ping_mu_jie_tu_2025-08-20_151443.png

Figure 2. Exact and approximate fuzzy solutions of V-FIDEs of Eq. (13) when $w$ = 1

It also shows the valid region of the fuzzy convergence parameters $\tilde{y}_0(x)=[\underline{y}(x), \bar{y}(x)]$. Hence, the order of the fuzzy fractional LST solution can be seen clearly, as can the region of valid convergence control parameter ($\alpha$) values and the line segment at $\alpha=1$ for $\alpha[0,1]$.

Furthermore, using a distinct value for $\omega$ and various $\alpha$ would thus result in a triangle shape for the graph. When the value of $w=0,0.3,0.6,0.8$, $\alpha=0.75$ we would have created a single shape in Figure 3.

ping_mu_jie_tu_2025-08-20_151502.png

Figure 3. Fuzzy approximate solution where $w=0,0.3,0.6,0.8$

Eq. (5) delineates a triangular fuzzy number, which underpins the fuzzy beginning conditions employed in this work. The fuzzy approximate solutions are shown graphically in Figure 3, with special attention to the triangle form that is depicted there. This form directly embodies the theoretical characteristics of the fuzzy solution framework, particularly the application of triangular fuzzy numbers and $\alpha$-cuts. The triangle arrangement verifies that the fuzziness disseminates in accordance with the established features of fuzzy sets. A comprehensive examination provides insights into the exact and approximate solutions for various values of www and $\alpha$, emphasizing their respective similarities and distinctions. The graphical findings validate that the suggested method produces fuzzy approximate solutions that nearly correspond to the exact solution. The strong similarity supports the accuracy of the triangular fuzzy representation and shows that the chosen method works well within the fuzzy theoretical framework.

Example (6.2): [36] Given the non-linear V-FIDEs

$\begin{gathered}\mathfrak{D}^\alpha \tilde{y}(x)=\frac{25 x^{6 / 5}}{3 \Gamma(1 / 5)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma(1 / 5)} +\int_0^t(t-\mathrm{s})[\tilde{y}(t)]^2 d \mathrm{~s}, 0<\alpha \leq 1, \mathrm{x} \in[0,1]\end{gathered}$                (21)

From the Definition 2.6 of the fuzzy number, we can defuzzify the initial condition let $[\tilde{0}]_\alpha$ and by Definition 2.5 triangular fuzzy number such that $\forall \alpha \in[0,1]$, we have non-linear V-FFIDEs:

$\begin{gathered}\mathfrak{D}^\alpha \tilde{y}(x ; \omega)=\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)} -[w-1,1-w]+\int_0^t(t-\mathrm{s})[\tilde{y}(t ; \omega)]^2 d \mathrm{~s} ， 0<\alpha \leq 1, \mathrm{x} \in[0,1]\end{gathered}$                    (22)

with fuzzy initial condition: $(\tilde{y}(0))^w=[w-1,1-w]$.

The exact solution of Eq. (21) when $\alpha$ = 1:

$[y(x)]^\omega=[w-1,1-w]\left(x^2-x\right)$

From Section 5, the $k^{\text {th }}$ order deformation equation of Eq. (11) and operating $\mathrm{I}^\alpha$ on the both sides of Eq. (22) and plugging property Eq. (4) of defined (2).

$\begin{gathered}\overline{\mathrm{y}}(x ; w)= \mathrm{I}^\alpha\binom{\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}-}{(w-1)+\int_0^t(t-\mathrm{s})[\bar{y}(t ; \omega)]^2 d \mathrm{~s}} \\ \underline{y}(x ; w)= \mathrm{I}^\alpha\binom{\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}-}{(1-w)+\int_0^t(t-\mathrm{s})[\underline{y}(t ; \omega)]^2 d \mathrm{~s}}\end{gathered}$                             (23)

By applied Shifted Legendre Polynomials in Section 3 yields.

Hence, in order to minimize $E 1, E 2$, it is necessary to determine the values of $a_i, b_i$ that will yield the optimal estimate for the solution of the fuzzy fractional integro differential Eq. (25). In order to obtain the least value of $E 1, E 2$, we need to calculate the partial derivatives of ${ }^{\mathrm{c}} \mathfrak{D}^\alpha y(x)$ with respect to $\alpha$, where $\alpha$ is a value between 0 and 1.

$\begin{aligned} & \sum_{\mathrm{i}=0}^2 \mathrm{a}_{\mathrm{i}} \overline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{x} ; w)-\mathrm{I}^\alpha\left(\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}\right. \left.+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}-(w-1)+\int_0^t(t-\mathrm{s})[\overline{\mathrm{y}}(t ; \omega)]^2 d \mathrm{~s}\right)=0 \\ & \sum_{\mathrm{i}=0}^2 b_{\mathrm{i}} \mathcal{P}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{x} ; w)-\mathrm{I}^\alpha\left(\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}\right. \left.+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}-(1-w)+\int_0^t(t-\mathrm{s})[\underline{y}(t ; \omega)]^2 d \mathrm{~s}\right)=0\end{aligned}$                          (24)

Thus, we minimized Eq. (24) as

$\begin{gathered}E 1\left(\underline{a}_0, \underline{a}_1, \underline{a}_2\right)^\omega=\sum_{\mathrm{i}=0}^2 \underline{a}_i \underline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{t} ; w) -I^\alpha\binom{\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}-(w-1)}{+\int_0^1(t-\mathrm{s})\left(\sum_{\mathrm{i}=0}^2 \mathrm{a}_{\mathrm{i}} \overline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(t ; w)\right)^2 \mathrm{~d} t} \\ E 2\left(\overline{\mathrm{~b}_0}, \overline{\mathrm{~b}_1}, \overline{\mathrm{~b}_2}\right)^\omega=\sum_{\mathrm{i}=0}^2 \overline{\mathrm{~b}_{\mathrm{i}}} \overline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{x} ; w) -I^\alpha\binom{\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}-(1-w)}{+\int_0^1(t-\mathrm{s})\left(\sum_{\mathrm{i}=0}^2 b_{\mathrm{i}} \underline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(t ; w)\right)^2 \mathrm{~d} t}\end{gathered}$                         (25)

$\left\{\begin{array}{l}\frac{\partial E}{\partial a_i}=0, i=0,1,2 \\ \frac{\partial E}{\partial b_i}=0, i=0,1,2\end{array}\right.$                    (26)

$\begin{gathered}\frac{\partial E}{\partial {a_i}}=\sum_{\mathrm{i}=0}^2 \underline{a}_i \underline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{t} ; w) \mathrm{I}^\alpha -\binom{\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}}{+\int_0^1(t-\mathrm{s})\left(\sum_{\mathrm{i}=0}^2 \underline{a}_i \underline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{t} ; w)\right)^2 \mathrm{~d} t} \\ \frac{\partial E}{\partial \overline{\mathrm{~b}}_{\mathrm{i}}}==\sum_{\mathrm{i}=0}^2 \overline{\mathrm{~b}}_{\mathrm{i}} \overline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{x} ; w)-\mathrm{I}^\alpha\binom{\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}}{+\int_0^1(t-\mathrm{s})\left(\sum_{\mathrm{i}=0}^2 \overline{\mathrm{~b}}_{\mathrm{i}} \overline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{x} ; w)\right)^2 \mathrm{~d} t}\end{gathered}$                                (27)

Applying Eq. (26) in Eq. (27) we get $\frac{\partial \mathrm{E} 1}{\partial \mathrm{a}_{\mathrm{i}}}$, $\frac{\partial \mathrm{E} 2}{\partial \mathrm{~b}_{\mathrm{i}}}$ for all $i=0,1,2$ causing E to be reduced.

By manipulating Eq. (28), we may deduce the subsequent normal equations, which allow us to determine the coefficients $\tilde{a}_0, \tilde{a}_1, \tilde{a}_2$ and $\tilde{b}_0, \tilde{b}_1, \tilde{b}_2$. To assess the precision of our proposed approach, we will employ Eq. (22) and present the results in the figures and tables provided below. Therefore, we can identify the convergence region for $I^\alpha$ by plotting the curves of the eighth order fuzzy Volterra linear fractional integro-differential equation LST lower bound solution $\underline{y}(x)$ and upper bound solution $\bar{y}(x)$ for Eq. (22) at $w=1$ in Figure 4.

$\begin{gathered}\left(\frac{\partial \mathrm{E} 1}{\partial\left(\mathrm{a}_0, \mathrm{a}_1, \mathrm{a}_2\right)}\right)^\omega= 2\left[\sum_{\mathrm{i}=0}^2 \underline{a_i} \underline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{t} ; w)-\mathrm{I}^\alpha\left(\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}+\right.\right. \\ \left.\left.\int_0^1(t-\mathrm{s})\left(\sum_{\mathrm{i}=0}^2 \mathrm{a}_{\mathrm{i}} \overline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(t ; w)\right)^2 \mathrm{~d} t\right)\right] \\ \left(\frac{\partial E 2}{\partial\left(\mathrm{~b}_0, \mathrm{~b}_1, \mathrm{a}_2\right)}\right)^\omega= 2\left[\sum_{\mathrm{i}=0}^2 \overline{\mathrm{~b}}_{\mathrm{i}} \overline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(\mathrm{x} ; w)-\mathrm{I}^\alpha\left(\frac{25 x^{\frac{6}{5}}}{3 \Gamma\left(\frac{1}{5}\right)}+\frac{x^6}{30}-\frac{x^5}{10}+\frac{x^4}{12}+\frac{5 \sqrt[5]{x}}{3 \Gamma\left(\frac{1}{5}\right)}+\right.\right. \\ \left.\left.\int_0^1(t-\mathrm{s})\left(\sum_{\mathrm{i}=0}^2 b_{\mathrm{i}} \underline{\mathcal{P}}_{\mathrm{l}_{\mathrm{i}+1}}(t ; w)\right)^2 \mathrm{~d} t\right)\right]\end{gathered}$                          (28)

As stated in Section 5, the chosen starting approximation for Eq. (22) is $\tilde{y}_0(x)=[\underline{y}(x), \bar{y}(x)]$. The values of ($\alpha$) can be determined from Eq. (24) in Section 5 by testing all the values of $I^\alpha$ at $\alpha$ = 0.5 for $\alpha \in[0,1]$. This allows us to identify the convergence region for $I^{0.5}$. Figure 4 shows the curves of the lower bound solution $\underline{y}(x)$ and the upper bound solution $\bar{y}(x)$ for Eq. (23) of the fuzzy Volterra linear fractional integro-differential equation LST with eight roots.

ping_mu_jie_tu_2025-08-20_151513.png

Figure 4. Exact and first order solution of V-FFIDEs of Eq. (22) when $w$ = 1 and $\mathrm{t} \in[0,1]$

Figure 4 illustrates the graphical representation of the fuzzy solution bounds for Eq. (23), which is a fuzzy linear fractional Volterra integro-differential equation solved using the (LST). The two curves represent the lower bound solution $\bar{y}(x)$ the upper bound solution $\underline{y}(x)$. These bounds encapsulate the fuzzy solution interval and are obtained under the condition of using eight distinct roots within the LST-based approach.

Table 2 shows the fuzzy approximate solutions by using LST for fractional orders of $\alpha=0.5,0.7,0.9$ for Eq. (22).

Table 2. Numerical results for Example 2, with different values of $\alpha$ at $w$

According to Section 5, the selected initial guess of Eq. (23) is $\tilde{y}_0(x)=[\underline{y}(x), \bar{y}(x)]$. From Section 5, the values of ($\alpha$) can be obtained from Eq. (24) and by testing all the values of $\mathrm{I}^\alpha$ for $\alpha \in[0,1]$.

Furthermore, if we vary the value of w, the resulting drawing will have a triangular shape, as depicted in Figure 5.

ping_mu_jie_tu_2025-08-20_151523.png

Figure 5. Lower and upper approximate solution were $w$=0, 0.3, 0.6, 0.8

Figure 5 illustrates the fuzzy approximate solutions represented as Triangular Fuzzy Numbers, as defined in Definition 5. The results are shown for different values of the parameters $w$ and $\alpha$. The comparison clearly indicates that the proposed method provides approximations that closely match the exact solution. Moreover, the similarity between the approximate solutions for different parameter values confirms the stability and accuracy of the suggested approach.