Abstract

Translate

Two nonlinear Duffing equations are numerically treated in this article. The nonlinear fractional-order Duffing equations and the second-order nonlinear Duffing equations are handled. Based on the collocation technique, we provide two numerical algorithms. To achieve this goal, a new family of basis functions is built by combining the sets of Fibonacci and Lucas polynomials. Several new formulae for these polynomials are developed. The operational matrices of integer and fractional derivatives of these polynomials, as well as some new theoretical results of these polynomials, are presented and used in conjunction with the collocation method to convert nonlinear Duffing equations into algebraic systems of equations by forcing the equation to hold at certain collocation points. To numerically handle the resultant nonlinear systems, one can use symbolic algebra solvers or Newton’s approach. Some particular inequalities are proved to investigate the convergence analysis. Some numerical examples show that our suggested strategy is effective and accurate. The numerical results demonstrate that the suggested collocation approach yields accurate solutions by utilizing Fibonacci–Lucas polynomials as basis functions.

Full text

Turn on search term navigation

Translate

1. Introduction

For a wide variety of applications, special functions are essential. One can consult [1,2,3] for specific findings and uses of special functions. The Fibonacci and Lucas polynomials are important special functions [4]. A number of authors were keen on presenting and exploring various extensions and variations of these polynomials.

Research on generalized Lucas polynomials and their connections to Fibonacci and Lucas polynomials was covered in [5], whereas [6] focused on Gauss–Fibonacci and Gauss–Lucas polynomials and their uses. Additional contributions on these polynomials and their uses may be found in [7,8,9,10].

Finding different formulae for special functions is something that a number of mathematicians are interested in. The development of numerical approaches to solve different types of differential equations (DEs) can greatly benefit from these equations. Many approaches for solving different DEs using spectral methods require expressing the derivatives of special polynomials as combinations of their original ones. This approach has been used in previous works (see [11,12]). Furthermore, for any orthogonal or nonorthogonal polynomials, it is an important goal to construct the operational matrices for the derivatives of these polynomials. These operational matrices aid in transforming the DEs using an appropriate spectral approach into an algebraic system of equations that can be solved using appropriate linear algebra procedures; see [13,14,15,16,17].

Nonlinear DEs are indispensable in many scientific fields, including physics, biology, chemistry, economics, and engineering. These equations describe systems that are too complex to be described by linear equations. The study of nonlinear differential equations (DEs)—which are utilized to depict phenomena such as chaos, turbulence, and nonlinear waves—provides valuable insights that find widespread application in fields ranging from climate modeling and fluid dynamics to telecommunications (e.g., [18,19]).

Fractional differential equations (FDEs) are fundamental in a number of areas of the applied sciences. Standard DEs are unable to catch some events which they explain. This is due to their exceptional ability to mimic genetic and memory functions. For instance, as mentioned in [20], they mimic a variety of physiological and biological processes, such as neuronal activity and tumor formation. Since these equations cannot be solved analytically, numerical analysis typically comes into play while solving them. For instance, as demonstrated in [21], a collocation approach might be useful when dealing with several equations. Some of these approaches are referenced in [22,23,24].

Among the important nonlinear DEs are the different Duffing equations introduced by the engineer Georg Duffing in 1918. They arise in physics and engineering to model a variety of physical phenomena. Examples of these uses include the description of a system’s chaotic behavior [25] and the control of a chaotic system’s mobility around less complicated attractors by the injection of modest dampening signals [26]. Many authors were interested in handling the different Duffing equations. For example, the authors of [27] found analytical solutions for some Duffing equations. The author in [28] applied the Pell–Lucas approach to treating the Duffing equation. A machine learning was used in [29]. A certain Runge–Kutta approach was utilized in [30]. A general solution of the Duffing equation of third-order nonlinearity was proposed in [31]. In [32], the authors used an adapted block hybrid method to handle Duffing equations. To access many contributions about diverse forms of Duffing equations, one may refer to [33,34,35,36,37,38].

Spectral methods have become widely recognized as a significant class of numerical techniques for addressing various problems in different disciplines. These methods have many advantages comparable to various numerical methods (see [39]). The approximate solutions obtained by these methods are highly accurate and provide exponential convergence rates. In addition, unlike the finite element or finite difference method, they provide global solutions, not local ones. These methods are adaptable in treating different types of differentiable equations. We can choose the suitable method that we can use according to the type of differential equation and the type of underlying conditions. There are three main spectral methods. Every method has its advantages and uses. The Galerkin and Petrov–Galerkin methods can be applied successfully for linear problems and some specific nonlinear problems; see, for example, [40,41,42]. The Tau method has a wider range of applications than the Galerkin method due to its ability to handle more complex boundary conditions; for example, see [43,44,45]. The collocation method is advantageous, since it can handle any type of differential equation regulated with any type of underlying conditions, for example, [46,47,48,49,50].

In this paper, we are devoted to proposing two numerical algorithms to solve the second-order and fractional-order Duffing equations. The spectral collocation algorithm is applied for such a purpose. Moreover, a new set of basis functions that generalize the Fibonacci and Lucas basis are introduced and employed. To our knowledge, the basis functions introduced for our algorithm’s derivation are new. In what follows, we summarize the main points, including the novelty of our contribution in this paper:

Introducing a new type of generalized Fibonacci and Lucas polynomials.

Establishing some theoretical results concerning these polynomials that will be the backbone of our numerical results.

Designing a numerical algorithm for treating the nonlinear second-order Duffing equation.

Designing a numerical algorithm for treating the nonlinear fractional Duffing equation.

Discussing the error analysis of the proposed method.

Testing our algorithms numerically by presenting some numerical examples with some comparisons.

The advantages of our proposed technique can be summarized as follows:

By choosing combined Fibonacci-Lucas polynomials as basis functions, a few retained modes produce highly accurate approximations.

The approach requires fewer computations to achieve the desired precision.

We can obtain several approximate solutions based on the presence of two free parameters, a and b.

Our technique can treat both linear and non-linear equations.

Here is the outline of the paper: An overview of Fibonacci and Lucas polynomials is given in Section 2. In addition, a combined Fibonacci–Lucas class of polynomials is presented in this section. The nonlinear second-order Duffing problem is solved using a matrix collocation method in Section 3. In Section 4, the fractional-order Duffing equation is addressed using a collocation method. Section 5 discusses the expansion’s convergence and truncation error bound. Section 6 provides a few comparisons and examples. Some last thoughts are presented in Section 7.

2. Introducing a Unified Sequence of Fibonacci and Lucas Polynomials

This section discusses some essential features of Fibonacci and Lucas polynomials and introduces a unified sequence of Fibonacci and Lucas sequences.

2.1. Fibonacci and Lucas Polynomial Sequences

The Fibonacci and Lucas polynomial sequences can be generated using, respectively, the following two recursive formulas:

(1) $F_{i} (t) - t F_{i - 1} (t) - F_{i - 2} (t) = 0, F_{0} (t) = 0, F_{1} (t) = 1,$

(2) $L_{i} (t) - t L_{i - 1} (t) - L_{i - 2} (t) = 0, L_{0} (t) = 2, L_{1} (t) = t .$

The generating function for

F_{i} (t)

is given by

$G (t, z) = \frac{z}{1 - z^{2} - z t} = \sum_{i = 0}^{\infty} F_{i} (t) z^{i},$

while the generating function for

L_{i} (t)

is given by

$\bar{G} (t, z) = \frac{1 + z^{2}}{1 - z^{2} - z t} = \sum_{i = 0}^{\infty} L_{i} (t) z^{i} .$

For an approach for generating a function, one can refer to [51].

From (1) and (2), it is evident that both Fibonacci and Lucas polynomials satisfy the same recursive formula but with different initials; thus, it is clear that the following recurrence relation

(3) $ϕ_{i} (t) - t ϕ_{i - 1} (t) - ϕ_{i - 2} (t) = 0, ϕ_{0} (t) = a, ϕ_{1} (t) = b t,$

generalizes the two sequences in (1) and (2), and we will denote

ϕ_{i} (t) = F L_{i}^{a, b} (t)

, that is,

(4) $F L_{i}^{a, b} (t) - t F L_{i - 1}^{a, b} (t) - F L_{i - 2}^{a, b} (t) = 0, F L_{0}^{a, b} (t) = a, F L_{1}^{a, b} (t) = b t .$

It is also clear that

$F_{i + 1} (t) = F L_{i}^{1, 1} (t), L_{i} (t) = F L_{i}^{2, 1} (t) .$

We also have the following expressions

(5) $\begin{matrix} F_{m + 1} (t) = & \sum_{s = 0}^{⌊\frac{m}{2}⌋} (\binom{m - s}{s}) t^{m - 2 s}, \end{matrix}$

(6) $\begin{matrix} L_{m} (t) = & m \sum_{s = 0}^{⌊\frac{m}{2}⌋} \frac{(\binom{m - s}{s})}{m - s} t^{m - 2 s}, \end{matrix}$

and their inverse expressions

(7) $\begin{matrix} t^{m} = & \sum_{s = 0}^{⌊\frac{m}{2}⌋} \frac{{(- 1)}^{s} (m - 2 s + 1) {(m - s + 2)}_{s - 1}}{s!} F_{m - 2 s + 1} (t), \end{matrix}$

(8) $\begin{matrix} t^{m} = & \sum_{s = 0}^{⌊\frac{m}{2}⌋} \frac{c_{m - 2 s} {(- 1)}^{s} {(m - s + 1)}_{i}}{s!} L_{m - 2 s} (t), \end{matrix}$

where

$c_{m} = \{\begin{matrix} \frac{1}{2}, & m = 0, \\ 1, & m \geq 1 . \end{matrix}$

Remark 1.

The key idea to develop the formulas concerned with the polynomials $F L_{i}^{a, b} (t)$ that satisfy (3) is the following theorem, in which we will show that $F L_{i}^{a, b} (t)$ may be expressed as a combination of two Fibonacci polynomials.

Theorem 1.

Consider any non-negative integer j. The polynomials $F L_{j}^{a, b} (t)$ can be represented as

(9) $F L_{j}^{a, b} (t) = b F_{j + 1} (t) + (a - b) F_{j - 1} (t) .$

Proof.

Consider the following polynomial:

(10) $ξ_{j} (t) = b F_{j + 1} (t) + (a - b) F_{j - 1} (t) .$

It is clear that

ξ_{0} (t) = F L_{0}^{a, b} (t) = a, ξ_{1} (t) = F L_{1}^{a, b} (t) = b t

, and hence, it is sufficient to prove that

ξ_{j} (t)

satisfies the same recurrence relation of

F L_{j}^{a, b} (t)

, for

j \geq 2

, that is, we are going to prove that

$ξ_{j + 2} (t) - t ξ_{j + 1} (t) - ξ_{j} (t) = 0 .$

Using the recursive formula of the Fibonacci polynomials (1) in the form

$t F_{j} (t) = F_{j + 1} (t) - F_{j - 1} (t),$

along with the definition in (10), it can be shown that

$ξ_{j + 2} (t) - t ξ_{j + 1} (t) - ξ_{j} (t) = 0 .$

This proves the theorem. □

The inverse connection formula of (9) is also interesting. The following theorem exhibits this result.

Theorem 2.

The Fibonacci polynomials are linked by Fibonacci–Lucas polynomials by the following two formulas:

(11) $\begin{matrix} F_{2 i + 1} (t) = & \frac{{(b - a)}^{i}}{b^{i}} + \sum_{r = 0}^{i - 1} \frac{{(b - a)}^{r}}{b^{r + 1}} F L_{2 i - 2 r}^{a, b} (t), \end{matrix}$

(12) $\begin{matrix} F_{2 i + 2} (t) = & \sum_{r = 0}^{i} \frac{{(b - a)}^{r}}{b^{r + 1}} F L_{2 i - 2 r + 1}^{a, b} (t) . \end{matrix}$

Proof.

The proofs of (11) and (12) are similar. We prove now (12). We will show that the right-hand side of (12) equals its left-hand side. Now, let

$η_{i} (x) = \sum_{r = 0}^{i} \frac{{(b - a)}^{r}}{b^{r + 1}} F L_{2 i - 2 r + 1}^{a, b} (t),$

and we will prove that

η_{i} (t) = F_{2 i + 2} (t)

Now, making use of Formula (9), we can write

$η_{j} (t) = \sum_{r = 0}^{i} {(b - a)}^{r} b^{- 1 - r} (b F_{2 i - 2 r + 2} (t) + (a - b) F_{2 i - 2 r} (t)),$

and therefore, we obtain

$η_{i} (t) = \sum_{r = - 1}^{i - 1} {(b - a)}^{r + 1} b^{- r - 1} F_{2 i - 2 r} (x) + \sum_{r = 0}^{i} {(- 1)}^{r} {(a - b)}^{r + 1} b^{- r - 1} F_{2 i - 2 r} (t),$

which is equivalent to

(13) $η_{i} (t) = F_{2 + 2 i} (x) + \sum_{r = 0}^{i - 1} {(b - a)}^{r + 1} b^{- r - 1} F_{2 i - 2 r} (t) + \sum_{r = 0}^{i} {(- 1)}^{r} {(a - b)}^{r + 1} b^{- r - 1} F_{2 i - 2 r} (t) .$

Noting that

F_{0} (t) = 0

, it is easy to show that the sum of the last two sums in (13) is zero, and accordingly,

$η_{i} (t) = F_{2 i + 2} (t) .$

This proves (12). □

Remark 2.

The two connection Formulas (11) and (12) can be merged to give

(14) $F_{i + 1} (t) = \sum_{r = 0}^{⌊\frac{i}{2}⌋} ϵ_{i - 2 r} \frac{{(b - a)}^{r}}{b^{r + 1}} F L_{i - 2 r}^{a, b} (t),$

where

$ϵ_{k} = \{\begin{matrix} \frac{b}{a}, & i f k = 0, \\ 1, & o t h e r w i s e . \end{matrix}$

Theorem 3.

Consider a positive integer k. The power for representation of $F L_{k}^{a, b} (t)$ is

(15) $\begin{matrix} F L_{k}^{a, b} (t) = & \sum_{m = 0}^{⌊\frac{k}{2}⌋} \frac{{(k - 2 m + 1)}_{m - 1} ((k - 2 m) b + m a)}{m!} t^{k - 2 m} . \end{matrix}$

Proof.

In virtue of the combination (9), together with (5), we obtain the following expression:

$F L_{k}^{a, b} (t) = b \sum_{m = 0}^{⌊\frac{k}{2}⌋} (\binom{k - m}{m}) t^{k - 2 m} + (a - b) \sum_{m = 0}^{⌊\frac{k}{2}⌋ - 1} (\binom{k - m - 2}{m}) t^{k - 2 m - 2},$

which can be written as

$F L_{k}^{a, b} (t) = b \sum_{m = 0}^{⌊\frac{k}{2}⌋} (\binom{k - m}{m}) t^{k - 2 m} + (a - b) \sum_{m = 1}^{⌊\frac{k}{2}⌋} (\binom{k - m - 1}{m - 1}) t^{k - 2 m} .$

The last formula is equivalent to

$F L_{k}^{a, b} (t) = \sum_{m = 0}^{⌊\frac{k}{2}⌋} (b (\binom{k - m}{m}) + (a - b) (\binom{k - m - 1}{m - 1})) t^{k - 2 m} .$

Noting the identity

$b (\binom{k - m}{m}) + (a - b) (\binom{k - m - 1}{m - 1}) = \frac{{(k - 2 m + 1)}_{m - 1} ((k - 2 m) b + m a)}{m!},$

then the expression in (15) can be obtained. □

Remark 3.

It is useful to express Formula (15) in the following alternative form:

$F L_{k}^{a, b} (t) = \sum_{m = 0}^{k} \frac{(\frac{a (k - m)}{2} + b m) δ_{k + m} {(m + 1)}_{\frac{k - m}{2} - 1}}{(\frac{k - m}{2})!} t^{m},$

where

(16) $δ_{r} = \{\begin{matrix} 1, & i f r even, \\ 0, & o t h e r w i s e . \end{matrix}$

Lemma 1.

The inversion formula of $F L_{i}^{a, b} (t)$ is

(17) $t^{m} = \sum_{r = 0}^{⌊\frac{m}{2}⌋} ϵ_{m - 2 r} S_{r, m} F L_{m - 2 r}^{a, b} (t),$

where

$S_{r, m} = \sum_{j = 0}^{r} \frac{{(- 1)}^{r} a^{r - j} (m - r - j) {(m - r + 1)}_{j - 1}}{j! b^{r + 1 - j}} .$

Proof.

Starting from the inversion formula of Fibonacci polynomials, one may write

(18) $t^{m} = \sum_{r = 0}^{⌊\frac{m}{2}⌋} \frac{{(- 1)}^{r} (m - 2 r + 1) {(m - r + 2)}_{r - 1}}{r!} F_{m - 2 r + 1} (t) .$

inserting the connection formula (14) into (18) yields

$t^{m} = \sum_{r = 0}^{⌊\frac{m}{2}⌋} \frac{{(- 1)}^{r} (m - 2 r + 1) {(m - r + 2)}_{r - 1}}{r!} \sum_{k = 0}^{⌊\frac{1}{2} (m - 2 r)⌋} ϵ_{m - 2 k - 2 r} {(b - a)}^{k} b^{- 1 - k} F L_{m - 2 r - 2 k}^{a, b} (t) .$

After some algebraic computations, the last formula leads to Formula (17). □

2.2. Derivatives and Operational Matrices of the Fibonacci–Lucas Polynomials

In this part, we will develop the high-order derivatives of the Fibonacci–Lucas polynomials and, after that, establish their operational matrices of integer derivatives, which will be pivotal in designing our numerical algorithm.

Theorem 4.

Consider two positive integers q and j with $j \geq q$ . The qth derivative of $F L_{i}^{a, b} (t)$ takes the form

(19) $\frac{d^{q} F L_{j}^{a, b} (t)}{d t^{q}} = \sum_{p = 0}^{j - q} χ_{j, p}^{q} F L_{p}^{a, b} (t),$

where

$χ_{j, p}^{q} = δ_{j + p - q} \sum_{r = 0}^{\frac{1}{2} (j - p - q)} \frac{(j - r - 1)! (a r + b (j - 2 r)) G_{\frac{j - p - q}{2} - r, j - q - 2 r}}{r! (j - q - 2 r)!},$

and

$G_{m, r} = \{\begin{matrix} \frac{1}{a}, & i f m = 0, \\ ϵ_{r - 2 m} \sum_{j = 0}^{m} \frac{{(- 1)}^{m + 1} a^{m - j} b^{j - m - 1} (j + m - r) {(- m + r + 1)}_{j - 1}}{j!}, & i f m > 0, \end{matrix}$

and $δ_{k}$ is defined as in (16).

Proof.

The analytic formula in (15) yields the following formula:

$\frac{d^{q} F L_{j}^{a, b} (t)}{d t^{q}} = \sum_{r = 0}^{⌊\frac{j - q}{2}⌋} \frac{(b (j - 2 r) + a r) (j - r - 1)!}{r! (j - q - 2 r)!} t^{j - 2 r - q} .$

By virtue of (17), the last formula turns into

$\begin{matrix} \frac{d^{q} F L_{j}^{a, b} (t)}{d t^{q}} = \sum_{r = 0}^{⌊\frac{j - q}{2}⌋} \frac{(b (j - 2 r) + a r) (j - r - 1)!}{r! (j - q - 2 r)!} \sum_{s = 0}^{⌊\frac{j - q}{2}⌋ - r} G_{s, j - 2 r - q} F L_{j - 2 r - q - 2 s}^{a, b} (t), \end{matrix}$

where

$G_{m, r} = \{\begin{matrix} \frac{1}{a}, & if m = 0, \\ ϵ_{r - 2 m} \sum_{j = 0}^{m} \frac{{(- 1)}^{m + 1} a^{m - j} b^{j - m - 1} (j + m - r) {(- m + r + 1)}_{j - 1}}{j!}, & if m > 0, \end{matrix}$

and hence, some computations lead to

(20) $\frac{d^{q} F L_{j}^{a, b} (t)}{d t^{q}} = \sum_{p = 0}^{⌊\frac{j - q}{2}⌋} (\sum_{r = 0}^{p} \frac{((b (j - 2 r) + a r) (j - r - 1)!) G_{p - r, j - 2 r - q}}{r! (j - q - 2 r)!}) F L_{j - q - 2 p}^{a, b} (t) .$

It is easy to see that (20) is equivalent to (19). □

Now, if we consider the vector $U (t)$ defined as

(21) $U (t) = {[F L_{0}^{a, b} (t), F L_{1}^{a, b} (t), \dots, F L_{N}^{a, b} (t)]}^{T},$

then based on Formula (19), we can write the following general derivative expression:

(22) $\frac{d^{q} U (t)}{d t^{q}} = M_{q} U (t),$

where

M_{q} = (M_{i, p, q})

is the general operational matrix of derivatives of order

{(N + 1)}^{2}

whose elements can be written in the following form:

$M_{i, p, q} = \{\begin{matrix} χ_{i, p}^{q}, & if i > p + q - 1, \\ 0, & o t h e r w i s e . \end{matrix}$

For our subsequent purposes, it is necessary to compute the two operational matrices of derivatives for the two cases corresponding to

q = 1

and

q = 2

. The following corollary presents these results.

Corollary 1.

For $q = 1, 2$ , Formula (22) gives, respectively, the following two derivative expressions:

(23) $\begin{matrix} \frac{d U (t)}{d t} = H U (t), \end{matrix}$

(24) $\begin{matrix} \frac{d^{2} U (t)}{d t^{2}} = F U (t), \end{matrix}$

where $H = (H_{i, p})$ and $F = (F_{i, p})$ are operational matrices of derivatives of order ${(N + 1)}^{2}$ whose elements can be written in the following form:

$\begin{matrix} H_{i, p} = \{\begin{matrix} χ_{i, p}^{1}, & i f i > p, \\ 0, & o t h e r w i s e, \end{matrix} \\ F_{i, p} = \{\begin{matrix} χ_{i, p}^{2}, & i f i > p + 1, \\ 0, & o t h e r w i s e . \end{matrix} \end{matrix}$

As an example, for $N = 5$ , $H$ and $F$ take the following forms:

$H = (\begin{matrix} 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{b}{a} & 0 & 0 & 0 & 0 & 0 \\ 0 & 2 & 0 & 0 & 0 & 0 \\ \frac{b}{a} - 2 & 0 & 3 & 0 & 0 & 0 \\ 0 & - \frac{2 a}{b} & 0 & 4 & 0 & 0 \\ \frac{2 a}{b} - 2 + \frac{b}{a} & 0 & - \frac{2 a}{b} - 1 & 0 & 5 & 0 \end{matrix}),$

$F = (\begin{matrix} 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 \\ \frac{2 b}{a} & 0 & 0 & 0 & 0 & 0 \\ 0 & 6 & 0 & 0 & 0 & 0 \\ \frac{4 b}{a} - 10 & 0 & 12 & 0 & 0 & 0 \\ 0 & - \frac{2 (7 a + b)}{b} & 0 & 20 & 0 & 0 \end{matrix}) .$

Remark 4.

The operational matrices of integer derivatives in (22) will play an essential role in deriving the proposed algorithm.

Remark 5.

After establishing the fundamental background for the combined Fibonacci–Lucas polynomials, they may be utilized in solving other types of differential equations, both linear and nonlinear, using the matrix approach.

Remark 6.

Although the Fibonacci and Lucas polynomials were utilized in several publications to act as basis functions in spectral methods, it is worth mentioning that our Fibonacci–Lucas polynomial basis has the advantage of merging both the Fibonacci and Lucas polynomial bases to obtain several approximate solutions.

3. A Matrix Collocation Approach for the Nonlinear Second-Order Duffing Equation

In this section, we consider the following nonlinear second-order Duffing equation (NSDE) [28,52]:

(25) $v^{″} (t) + f v^{'} (t) + h v (t) + c v^{3} (t) + d v^{5} (t) + e v^{7} (t) = g (t), t \in [0, T],$

which is subject to the conditions

(26) $v (0) = g_{1}, v^{'} (0) = g_{2} .$

The main idea to solve (25) and (26) is to employ the operational matrices of the derivatives of

F L_{i}^{a, b} (t)

, together with applying the collocation method.

Now, let us define the following space function:

$Δ_{N} = span {F L_{i}^{a, b} (t) : i = 0, 1, \dots, N} .$

Then, any function

v_{N} (t) \in Δ_{N}

may be approximated as

(27) $v_{N} (t) \approx \sum_{i = 0}^{N} θ_{i} F L_{i}^{a, b} (t) = θ U (t),$

where

U (t)

is defined in (21), and

$θ = [θ_{0}, θ_{1}, \dots, θ_{N}] .$

By virtue of (23), (24), and (27), the residual

R (t)

of Equation (25) is given by

(28) $\begin{matrix} R (t) & = v_{N}^{″} (t) + f v_{N}^{'} (t) + h v_{N} (t) + c {v_{N}}^{3} (t) + d {v_{N}}^{5} (t) + e {v_{N}}^{7} (t) - g (t) \\ = θ F U (t) + f θ H U (t) + h θ U (t) + c {(θ U (t))}^{3} + d {(θ U (t))}^{5} + e {(θ U (t))}^{7} - g (t) . \end{matrix}$

Now, by using the collocation method, we can obtain the following

(N + 1)

equation algebraic system in the unknown expansion coefficients of

θ_{i}

(29) $\begin{matrix} R (t_{i}) = 0, i = 1, 2, \dots, N - 1, \\ θ U (0) = g_{1}, θ \frac{d U (0)}{d t} = g_{2}, \end{matrix}$

where

{t_{i} = \frac{i}{N + 1} : i = 1, 2, \dots, N - 1}

are some collocation points. Therefore, a system in (29) can be solved to obtain

θ_{i}

with the aid of the well-known Newton’s iterative method.

4. A Matrix Collocation Approach for the Nonlinear Fractional-Order Duffing Equation

This section is confined to presenting a numerical algorithm for treating the nonlinear fractional Duffing equation. First, some fundamental properties regarding the fractional calculus are mentioned below.

Definition 1

([53,54,55,56]). The Gerasimov–Caputo fractional derivative of order μ is defined as

$D^{μ} h (z) = \frac{1}{Γ (m - μ)} \int_{0}^{z} {(z - t)}^{m - μ - 1} h^{(m)} (t) d t, μ > 0, z > 0,$

where $m - 1 < μ \leq m, m \in N$ .

The operator $D^{μ}$ satisfies the following properties for all $m - 1 < μ \leq m, m \in N$ :

(30) $\begin{matrix} D^{μ} b = & 0, b is a constant, \end{matrix}$

(31) $\begin{matrix} D^{μ} z^{m} = & \{\begin{matrix} 0, & if m \in N_{0} a n d m < ⌈ μ ⌉, \\ \frac{Γ (m + 1)}{Γ (m - μ + 1)} z^{m - μ}, & if m \in N_{0} a n d m \geq ⌈ μ ⌉, \end{matrix} \end{matrix}$

where

N = {1, 2, 3, \dots},

N_{0} = {0} \cup N

, and the notation

⌈ μ ⌉

denotes the ceiling function.

Now, consider the following nonlinear fractional-order Duffing equation (NFDE) [28,52]:

(32) $D^{α} v (t) + f D^{β} v (t) + h v (t) + c v^{3} (t) + d v^{5} (t) + e v^{7} (t) = g (t), t \in [0, T],$

which is directed to the constraints

(33) $v (0) = g_{1}, v^{'} (0) = g_{2},$

where

α \in (1, 2), β \in (0, 1) .

4.1. The Operational Matrix of Fractional Derivatives for $F L_{i}^{a, b} (t)$

Theorem 5.

Let $U (t)$ be the vector defined in (21). The following formula holds for all $μ > 0$ and $t \in [0, T]$ .

(34) $D^{μ} U (t) = t^{- μ} λ^{(μ)} U (t),$

where

(35) $λ^{(μ)} = {[\begin{matrix} 0 & 0 & 0 & \dots & 0 \\ ⋮ & ⋮ & ⋮ & \dots & ⋮ \\ γ_{⌈ μ ⌉, 0}^{μ} & γ_{⌈ μ ⌉, ⌈ μ ⌉}^{μ} & 0 & 0 \\ ⋮ & ⋮ & ⋮ & ⋮ \\ γ_{i, 0}^{μ} & \dots & γ_{i, 0}^{μ} & \dots & 0 \\ ⋮ & ⋮ & ⋮ & ⋮ \\ γ_{N, 0}^{μ} & γ_{N, 1}^{μ} & γ_{N, 2}^{μ} & \dots & γ_{N, N}^{μ} \end{matrix}]}_{(N + 1) \times (N + 1)} .$

Proof.

The proof of this theorem can be divided into two cases corresponding to the value of k and $⌈ μ ⌉$ :

Case 1: If $⌈ μ ⌉ \leq k \leq N + 1$ .

Applying $D^{μ}$ to the series of of $F L_{k}^{a, b} (t)$ in (3) yields

(36) $D^{μ} F L_{k}^{a, b} (t) = \sum_{m = ⌈ μ ⌉}^{k} \frac{(\frac{a (k - m)}{2} + b m) δ_{k + m} {(m + 1)}_{\frac{k - m}{2} - 1} m!}{(\frac{k - m}{2})! Γ (m - μ + 1)} t^{m - μ} .$

As a result of (17), the last relation can be written as

(37) $D^{μ} F L_{k}^{a, b} (t) = t^{- μ} \sum_{m = ⌈ μ ⌉}^{k} \sum_{r = 0}^{⌊\frac{m}{2}⌋} \frac{(\frac{a (k - m)}{2} + b m) δ_{k + m} {(m + 1)}_{\frac{k - m}{2} - 1} m! ϵ_{m - 2 r} S_{r, m}}{(\frac{k - m}{2})! Γ (m - μ + 1)} F L_{m - 2 r}^{a, b} (t),$

which can be written again in the form

(38) $D^{μ} F L_{k}^{a, b} (t) = t^{- μ} \sum_{n = 0}^{k} \sum_{m = ⌈ μ ⌉}^{k} \frac{ϵ_{n} δ_{k + m} δ_{n + m} Γ (\frac{k + m}{2}) S_{\frac{m - n}{2}, m} (a (k - m) + 2 b m)}{2 a (\frac{k - m}{2})! Γ (m - μ + 1)} F L_{n}^{a, b} (t),$

which can be rewritten as

(39) $D^{μ} F L_{k}^{a, b} (t) = t^{- μ} \sum_{n = 0}^{k} γ_{k, n}^{μ} F L_{n}^{a, b} (t),$

where

$γ_{k, n}^{μ} = \sum_{m = ⌈ μ ⌉}^{k} \frac{ϵ_{n} δ_{k + m} δ_{n + m} Γ (\frac{k + m}{2}) S_{\frac{m - n}{2}, m} (a (k - m) + 2 b m)}{2 a (\frac{k - m}{2})! Γ (m - μ + 1)} .$

Case 2: If $0 \leq k \leq ⌈ μ ⌉ - 1$ .

Here, we find that

(40) $D^{μ} F L_{k}^{a, b} (t) = 0 .$

Finally, Cases 1 and 2 may be joined in matrix form as

$D^{μ} U (t) = t^{- μ} λ^{(μ)} U (t),$

where the elements of the matrix

λ^{(μ)}

are given in the following form:

$λ_{k, n}^{(μ)} = \{\begin{matrix} γ_{k, n}^{μ}, & i f k \geq ⌈ μ ⌉, k \geq n, \\ 0, & otherwise . \end{matrix}$

This ends the proof of this theorem. □

4.2. Collocation Algorithm for the NFDE

Using similar procedures as in the preceding section, we can obtain the following residual based on Theorem 5:

(41) $\begin{matrix} \bar{R} (t) & = D^{α} v_{N} (t) + f D^{β} v_{N} (t) + h v_{N} (t) + c {v_{N}}^{3} (t) + d {v_{N}}^{5} (t) + e {v_{N}}^{7} (t) - g (t) \\ = t^{- α} θ λ^{(α)} U (t) + f t^{- β} θ λ^{(β)} U (t) + h θ U (t) + c {(θ U (t))}^{3} + d {(θ U (t))}^{5} \\ + e {(θ U (t))}^{7} - g (t) . \end{matrix}$

As a result, we may obtain the following $(N + 1)$ system of equations by using the collocation method:

$\begin{matrix} \bar{R} (t_{i}) = 0, i = 1, 2, \dots, N - 1, \\ θ U (0) = g_{1}, θ \frac{d U (0)}{d t} = g_{2}, \end{matrix}$

where

{t_{i} = \frac{i}{N + 1} : i = 1, 2, \dots, N - 1} .

The above system can be solved to obtain $θ_{i}$ with the aid of the well-known Newton’s iterative method.

5. Error Bound

In this section, we aim to demonstrate that when N approaches infinity, $R (t)$ and $\bar{R} (t)$ converge to zero. For the unknown function $v (t)$ , we derive various error bounds and derivatives of this function.

Theorem 6.

Assume that $\frac{d^{i} v (t)}{d t^{i}} \in C ([0, T])$ , $i = 0, 1, 2, \dots, N$ , let $v_{N} (t)$ be the proposed approximate solution belonging to $Δ_{N}$ , and define

$ℓ_{N} = sup_{t \in [0, T]} |\frac{d^{N + 1} v (t)}{d t^{N + 1}}| .$

Consequently, this estimate holds:

${∥ v (t) - v_{N} (t) ∥}_{2} \leq \frac{ℓ_{N} T^{N + \frac{3}{2}}}{\sqrt{2 N + 3} (N + 1)!} .$

Proof.

Consider the following Taylor expansion of $v (t)$ about the point $t = 0$ :

(42) $χ_{N} (t) = \sum_{i = 0}^{N} {(\frac{d^{i} v (t)}{d t^{i}})}_{t = 0} \frac{t^{i}}{i!},$

$v (t) - χ_{N} (t) = \frac{t^{N + 1}}{(N + 1)!} {(\frac{d^{N + 1} v (t)}{d t^{N + 1}})}_{t = c}, c \in [0, T] .$

Because

v_{N} (t)

is the best approximation solution of

v (t)

, we have, using the idea of best approximation,

(43) $\begin{matrix} {∥ v (t) - v_{N} (t) ∥}_{2}^{2} & \leq {∥ v (t) - χ_{N} (t) ∥}_{2}^{2} \\ \leq \int_{0}^{T} \frac{{ℓ_{N}}^{2} t^{2 (N + 1)}}{{((N + 1)!)}^{2}} d t \\ = \frac{{ℓ_{N}}^{2} T^{2 N + 3}}{(2 N + 3) {((N + 1)!)}^{2}}, \end{matrix}$

and therefore, we have

${∥ v (t) - v_{N} (t) ∥}_{2} \leq \frac{ℓ_{N} T^{N + \frac{3}{2}}}{\sqrt{2 N + 3} (N + 1)!} .$

□

Theorem 7.

Suppose that $v (t)$ and $v_{N} (t)$ meet the assumption of Theorem 6 and

$σ_{N, m} = sup_{t \in [0, T]} |{(\frac{d^{N + 1} v (t)}{d t^{N + 1}})}^{m}|, m = 2, 3, 4, \dots .$

The following estimation holds:

${∥v^{m} (t) - v_{N}^{m} (t)∥}_{2} \leq \frac{σ_{N, m} T^{m (N + 1) + \frac{1}{2}}}{\sqrt{(2 m (N + 1) + 1)} {((N + 1)!)}^{m}} .$

Proof.

At the point $t = 0$ , we can use the Taylor expansion of $v (t)$ in (42) to write

(44) $v^{m} (t) - χ_{N}^{m} (t) = \frac{t^{m (N + 1)}}{{((N + 1)!)}^{m}} {(\frac{d^{N + 1} v (t)}{d t^{N + 1}})}_{t = c}^{m}, c \in [0, T] .$

Imitating similar steps as in Theorem 6 in accordance with the best approximation’s concept, one has

$\begin{matrix} {∥ v^{m} (t) - v_{N}^{m} (t) ∥}_{2}^{2} & \leq {∥ v^{m} (t) - χ_{N}^{m} (t) ∥}_{2}^{2} \\ \leq \int_{0}^{T} \frac{{σ_{N, m}}^{2} t^{2 m (N + 1)}}{{((N + 1)!)}^{2 m}} d t \\ = \frac{{σ_{N, m}}^{2} T^{2 m (N + 1) + 1}}{(2 m (N + 1) + 1) {((N + 1)!)}^{2 m}}, \end{matrix}$

and this leads to

${∥v^{m} (t) - v_{N}^{m} (t)∥}_{2} \leq \frac{σ_{N, m} T^{m (N + 1) + \frac{1}{2}}}{\sqrt{(2 m (N + 1) + 1)} {((N + 1)!)}^{m}} .$

□

Theorem 8.

Suppose that $v (t)$ , $v_{N} (t)$ , and $\frac{d^{i} v (t)}{d t^{i}}$ meet the assumption of Theorem 6 and

$τ_{N, m} = sup_{t \in [0, T]} |\frac{d^{N - m + 1} v (t)}{d t^{N - m + 1}}|, m \in N .$

Then, the following estimation holds:

${∥\frac{d^{m}}{d t^{m}} [v (t) - v_{N} (t)]∥}_{2} \leq \frac{τ_{N, m} T^{N - m + \frac{3}{2}}}{\sqrt{2 (N - m) + 3} (N - m + 1)!} .$

Proof.

Assume that $\frac{d^{m} χ_{N} (t)}{d t^{m}}$ is the Taylor expansion of $\frac{d^{m} v (t)}{d t^{m}}$ about the point $t = 0$ ; then, the residual between $\frac{d^{m} v (t)}{d t^{m}}$ and $\frac{d^{m} χ_{N} (t)}{d t^{m}}$ can be written as

$\frac{d^{m}}{d t^{m}} [v (t) - χ_{N} (t)] = \frac{t^{N - m + 1}}{(N - m + 1)!} {(\frac{d^{N - m + 1} v (t)}{d t^{N - m + 1}})}_{t = n}, n \in [0, T] .$

Since

\frac{d^{m} v_{N} (t)}{d t^{m}}

is the best approximate solution of

\frac{d^{m} v (t)}{d t^{m}}

, then according to the definition of the best approximation, we obtain

${∥\frac{d^{m}}{d t^{m}} [v (t) - v_{N} (t)]∥}_{2} \leq {∥\frac{d^{m}}{d t^{m}} [v (t) - χ_{N} (t)]∥}_{2} .$

We obtain the desired result by performing steps as in Theorem 6. □

Theorem 9.

Assume that the Gerasimov–Caputo operator $D^{μ} v (t) \in C ([0, T])$ and the conditions of Theorem 6 hold. Then,

${∥D^{μ} [v (t) - v_{N} (t)]∥}_{2} \leq \frac{ℓ_{N} T^{N - μ + \frac{3}{2}}}{\sqrt{2 (N - μ) + 3} Γ (N - μ + 2)} .$

Proof.

Using Equation (43) along with the application of the operator $D^{μ}$ , one obtains

$\begin{matrix} {∥ D^{μ} [v (t) - v_{N} (t)] ∥}_{2}^{2} & \leq \int_{0}^{T} \frac{{ℓ_{N}}^{2} t^{2 (N - μ + 1)}}{{(Γ (N - μ + 2))}^{2}} d t \\ = \frac{{ℓ_{N}}^{2} T^{2 (N - μ) + 3}}{(2 (N - μ) + 3) {(Γ (N - μ + 2))}^{2}} . \end{matrix}$

Therefore, we obtain

${∥D^{μ} [v (t) - v_{N} (t)]∥}_{2} \leq \frac{ℓ_{N} T^{N - μ + \frac{3}{2}}}{\sqrt{2 (N - μ) + 3} Γ (N - μ + 2)} .$

□

Theorem 10.

Let $R (t)$ be the residual of Equation (25) given by (28). Then, ${∥ R (t) ∥}_{2}$ will be sufficiently small for the sufficiently large values of N.

Proof.

$R (t)$ of Equation (28) can be written as

(45) $\begin{matrix} R (t) & = v_{N}^{″} (t) + f v_{N}^{'} (t) + h v_{N} (t) + c {v_{N}}^{3} (t) + d {v_{N}}^{5} (t) + e {v_{N}}^{7} (t) - g (t) \\ = \frac{d^{2}}{d t^{2}} [v_{N} (t) - v (t)] + f \frac{d}{d t} [v_{N} (t) - v (t)] + h [v_{N} (t) - v (t)] + c [{v_{N}}^{3} (t) - v^{3} (t)] \\ + d [{v_{N}}^{5} (t) - v^{5} (t)] + e [{v_{N}}^{7} (t) - v^{7} (t)] . \end{matrix}$

Taking

{∥ . ∥}_{2}

and using Theorems 6–8, we obtain

(46) $\begin{matrix} {∥ R (t) ∥}_{2} & \leq \frac{τ_{N, 2} T^{N - \frac{1}{2}}}{\sqrt{2 N - 1} (N - 1)!} + f \frac{τ_{N, 1} T^{N + \frac{1}{2}}}{\sqrt{2 N + 1} (N)!} + h \frac{ℓ_{N} T^{N + \frac{3}{2}}}{\sqrt{2 N + 3} (N + 1)!} \\ + c \frac{σ_{N, 3} T^{3 N + \frac{7}{2}}}{\sqrt{(6 N + 7)} {((N + 1)!)}^{3}} + d \frac{σ_{N, 5} T^{5 N + \frac{11}{2}}}{\sqrt{(10 N + 11)} {((N + 1)!)}^{5}} \\ + e \frac{σ_{N, 7} T^{7 N + \frac{15}{2}}}{\sqrt{(14 N + 15)} {((N + 1)!)}^{7}} . \end{matrix}$

At last, it is clear from Equation (46) that

{∥ R (t) ∥}_{2}

will be small enough for suitably large values of N. This concludes the proof of the theorem. □

Theorem 11.

Let $\bar{R} (t)$ be the residual of Equation (32) given by (41). Then, ${∥ \bar{R} (t) ∥}_{2}$ will be sufficiently small for the sufficiently large values of N.

Proof.

$\bar{R} (t)$ of Equation (41) can be written as

(47) $\begin{matrix} \bar{R} (t) & = D^{α} v_{N} (t) + f D^{β} v_{N} (t) + h v_{N} (t) + c {v_{N}}^{3} (t) + d {v_{N}}^{5} (t) + e {v_{N}}^{7} (t) - g (t) \\ = D^{α} [v_{N} (t) - v (t)] + f D^{β} [v_{N} (t) - v (t)] + h [v_{N} (t) - v (t)] + c [{v_{N}}^{3} (t) - v^{3} (t)] \\ + d [{v_{N}}^{5} (t) - v^{5} (t)] + e [{v_{N}}^{7} (t) - v^{7} (t)] . \end{matrix}$

Taking

{∥ . ∥}_{2}

and using Theorems 6–9, we obtain

(48) $\begin{matrix} {∥ \bar{R} (t) ∥}_{2} & \leq \frac{ℓ_{N} T^{N - α + \frac{3}{2}}}{\sqrt{2 (N - α) + 3} Γ (N - α + 2)} + f \frac{ℓ_{N} T^{N - β + \frac{3}{2}}}{\sqrt{2 (N - β) + 3} Γ (N - β + 2)} \\ + h \frac{ℓ_{N} T^{N + \frac{3}{2}}}{\sqrt{2 N + 3} (N + 1)!} + c \frac{σ_{N, 3} T^{3 N + \frac{7}{2}}}{\sqrt{(6 N + 7)} {((N + 1)!)}^{3}} + d \frac{σ_{N, 5} T^{5 N + \frac{11}{2}}}{\sqrt{(10 N + 11)} {((N + 1)!)}^{5}} \\ + e \frac{σ_{N, 7} T^{7 N + \frac{15}{2}}}{\sqrt{(14 N + 15)} {((N + 1)!)}^{7}} . \end{matrix}$

Finally, it is clear from Equation (48) that

{∥ \bar{R} (t) ∥}_{2}

will be small enough for sufficiently high N values. So, this theorem’s proof is complete. □

6. Illustrative Examples

Evaluation of our suggested collocation methods is the focus of this section. We solve a few test problems and present a few comparisons to ensure that our suggested methods are applicable and accurate.

The absolute errors (AEs) in the given tables are

$A E = | v (t) - v_{N} (t) | .$

Example 1

([28,52]). Consider the following NFDE:

$D^{α} v (t) + f D^{β} v (t) + h v (t) + c v^{3} (t) + d v^{5} (t) + e v^{7} (t) = g (t), t \in [0, T],$

which is subject to the conditions

$v (0) = \frac{1}{2}, v^{'} (0) = \frac{- 1}{2},$

and $g (t)$ is selected in a way that makes the exact solution become $v (t) = \frac{e^{- t}}{2} .$

Equation (1) is solved using our proposed algorithm for $T = 1$ and $T = 12$ :

Case 1: For $T = 1$ and $f = 4, h = - 2, c = 3, d = 0, e = 0$ , Table 1 presents the AEs at different values of $a, b$ at $N = 13$ when $α = 2, β = 1$ . Furthermore, Figure 1 shows the AEs at $a = 2, b = 3$ at different values of N when $α = 2, β = 1$ . Figure 2 shows that the approximate solutions have smaller variations for values of α and β near the values $α = 2$ and $β = 1$ when $a = 2, b = 3 .$ Table 2 presents the absolute errors (AEs) at different values of $a, b$ at $N = 12$ when $α = 1.9$ and $β = 0.9$ .

Case 2: For $T = 12$ , $α = 2, β = 1$ and $f = 2, h = 1, c = 8, d = 0, e = 0$ , Table 3 presents a comparison between our method at $N = 14$ and method in [52]. Figure 3 shows the AE (left) and exact, approximate solutions (right) of the example at $N = 17$ and $a = 2, b = 3,$ which demonstrates that the results of our method are extremely close to the exact solution.

Example 2

([52]). Consider the following NFDE:

$D^{α} v (t) + f D^{β} v (t) + h v (t) + c v^{3} (t) + d v^{5} (t) + e v^{7} (t) = g (t), t \in [0, T],$

which is subject to the conditions

$v (0) = 1, v^{'} (0) = 0,$

and $g (t)$ is selected in a way that makes the exact solution become $v (t) = cos (t)$ at $f = h = c = 1$ , where $d = 0$ .

Equation (2) is solved using our algorithm for $T = 1$ and $T = 5$ :

Case 1: For $T = 1$ , Table 4 and Table 5 present the AEs at different values of $a, b$ at, respectively, $N = 2$ and $N = 3$ when $α = 1.9, β = 0.9$ . Table 6 presents the AEs at different values of $a, b$ at $N = 12$ when $α = 2, β = 1$ . Figure 4 shows that the approximate solutions have smaller variations for values of α and β near the values $α = 2$ and $β = 1$ when $a = 1, b = 3 .$ Table 7 presents the absolute errors (AEs) at different values of $a, b$ at $N = 13$ when $α = 1.8$ and $β = 0.8$ .

Case 2: For $T = 5$ and $α = 2, β = 1$ , Figure 5 shows the AEs at different values of N when $a = b = 1$ . Also, Figure 6 illustrates the AEs at different values of N when $a = b = 2$ . These figures show the accuracy of our method.

Remark 7.

The results of Table 4 demonstrate that the small values of N cause clear variation in the error values for different choices of the parameters a and b. For instance, for $N = 2$ , the error changes from $2.8 \times 10^{- 2}$ to $1.05 \times 10^{- 4}$ at $N = 2$ due to the change of the involved parameters. Especially at larger time steps, the error differences brought on by variations in a and b become slightly less obvious.

Example 3.

Consider the following NFDE:

$D^{α} v (t) + f D^{β} v (t) + h v (t) + c v^{3} (t) + d v^{5} (t) + e v^{7} (t) = 0, t \in [0, T],$

which is governed by

$χ (0) = 0, χ^{'} (0) = 1 .$

Due to the nonavailability of the exact solution, let us define the following absolute residual error norm at $T = 1$ :

$R E = max_{t \in (0, 1)} |D^{α} v_{N} (t) + f D^{β} v_{N} (t) + h v_{N} (t) + c {v_{N}}^{3} (t) + d {v_{N}}^{5} (t) + e {v_{N}}^{7} (t)|,$

and apply our method at $f = h = c = d = e = T = 1$ when $α = 2, β = 1$ .

Figure 7 illustrates the RE (left) and approximate solution (right) at $a = b = 1$ and $N = 16$ . Also, Figure 8 illustrates the RE at different values of N when $a = 1$ and $b = 2 .$

Remark 8.

The numerical results obtained in this section show that we have received several highly accurate approximate solutions using the combined Fibonacci–Lucas polynomials. This gives us an advantage in introducing these generalized polynomials.

Remark 9.

We comment that the approximations resulting from utilizing other generalized polynomials, such as ultraspherical and Jacobi polynomials, do not change significantly due to the change of their parameters, especially for large values of the retained modes; see, for example, [57].

Remark 10.

The combined Fibonacci–Lucas polynomials provide excellent approximations, since an order error $10^{- 16}$ is sometimes reached for certain choices of $a, b$ , and N.

Remark 11.

Combining Fibonacci–Lucas polynomials leads to little improvement in numerical results, since the change of the two parameters a and b leads to small changes in the resulting errors.

7. Concluding Remarks

This paper established a generalized sequence of polynomials, namely, unified Fibonacci–Lucas polynomials. The well-known polynomial sequences of Fibonacci and Lucas are particular types of these polynomials. These polynomials have two parameters, yielding various solutions for every choice of them. Some theoretical results concerned with these polynomials were the keys to implementing our numerical algorithms for solving the second-order and the fractional-order Duffing nonlinear DEs via the celebrated collocation method. The operational matrices of derivatives of the Fibonacci–Lucas polynomials that are derived using the derivative formula of these polynomials were employed to design the proposed numerical algorithm. We comment here that For every choice of the two parameters a and b, a numerical solution was obtained. The numerical results show that the change in the absolute errors caused by variations in the two parameters of the combined Fibonacci–Lucas polynomials is minimal when choosing large values of the retained modes; however, these variations become larger for small values of the retained modes. We aim to investigate the impact of these parameters when solving other types of differential equations. We believe that this is the first time these polynomials have been employed in applications. Future research directions may involve employing these polynomials to solve other types of DEs.

Author Contributions

Conceptualization, W.M.A.-E. and A.G.A.; Methodology, W.M.A.-E., O.M.A. and A.G.A.; Software, W.M.A.-E. and A.G.A.; Validation, W.M.A.-E., O.M.A., A.K.A. and A.G.A.; Formal analysis, W.M.A.-E. and A.G.A.; Investigation, W.M.A.-E., O.M.A., A.K.A. and A.G.A.; Writing—original draft, W.M.A.-E. and A.G.A.; Writing—review & editing, W.M.A.-E. and A.G.A.; Supervision, W.M.A.-E.; Funding acquisition, A.K.A. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data are contained within the article.

Acknowledgments

The authors extend their appreciation to Umm Al-Qura University, Saudi Arabia, for funding this research work through grant number 25UQU4331287GSSR03.

Conflicts of Interest

The authors declare no conflicts of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1 The AEs of Example 1 at $a = 2, b = 3$ .

Figure 2 Different solutions of Example 1 at $N = 12$ and different values of $α, β$ .

Figure 3 The AEs (left) and exact, approximate solutions (right) of Example 1 at $N = 17$ and $a = 2, b = 3 .$

Figure 4 Different solutions of Example 2 at $N = 13$ and different values of $α, β$ .

Figure 5 The AEs of Example 2 at $a = b = 1$ .

Figure 6 The AEs of Example 2 at $a = b = 2$ .

Figure 7 The RE (left) and approximate solution (right) of Example 3 at $a = b = 1$ and $N = 16$ .

Figure 8 The RE of Example 3 at $a = 1, b = 2$ .

Table 1

AEs of Example 1 at $N = 13$ .

t	$a = 1, b = 1$	$a = 1, b = 3$	$a = 2, b = 1$	$a = 2, b = 2$
0.1	1.77636 × $10^{- 15}$	1.9984 × $10^{- 15}$	3.88578 × $10^{- 16}$	1.16573 × $10^{- 15}$
0.2	3.44169 × $10^{- 15}$	3.66374 × $10^{- 15}$	6.66134 × $10^{- 16}$	1.9984 × $10^{- 15}$
0.3	4.4964 × $10^{- 15}$	4.82947 × $10^{- 15}$	6.66134 × $10^{- 16}$	2.66454 × $10^{- 15}$
0.4	5.21805 × $10^{- 15}$	5.71765 × $10^{- 15}$	7.77156 × $10^{- 16}$	3.10862 × $10^{- 15}$
0.5	5.77316 × $10^{- 15}$	6.21725 × $10^{- 15}$	8.88178 × $10^{- 16}$	3.38618 × $10^{- 15}$
0.6	6.10623 × $10^{- 15}$	6.66134 × $10^{- 15}$	1.05471 × $10^{- 16}$	3.71925 × $10^{- 15}$
0.7	6.46705 × $10^{- 15}$	6.93889 × $10^{- 15}$	1.08247 × $10^{- 15}$	3.85803 × $10^{- 15}$
0.8	6.74461 × $10^{- 15}$	7.35523 × $10^{- 15}$	1.19349 × $10^{- 15}$	4.05231 × $10^{- 15}$
0.9	7.13318 × $10^{- 15}$	7.60503 × $10^{- 15}$	1.02696 × $10^{- 15}$	4.16334 × $10^{- 15}$

Table 2

AEs of Example 1 at $N = 12$ .

t	$a = 1, b = 3$	$a = 2, b = 4$	$a = 3, b = 5$
0.1	1.11022 × $10^{- 16}$	3.33067 × $10^{- 16}$	1.11022 × $10^{- 16}$
0.2	2.22045 × $10^{- 16}$	4.44089 × $10^{- 16}$	2.77556 × $10^{- 16}$
0.3	2.22045 × $10^{- 16}$	7.21645 × $10^{- 16}$	3.33067 × $10^{- 16}$
0.4	2.22045 × $10^{- 16}$	8.32667 × $10^{- 16}$	3.33067 × $10^{- 16}$
0.5	1.66533 × $10^{- 16}$	8.32667 × $10^{- 16}$	3.88578 × $10^{- 16}$
0.6	3.33067 × $10^{- 16}$	8.32667 × $10^{- 16}$	3.88578 × $10^{- 16}$
0.7	2.49805 × $10^{- 16}$	8.88178 × $10^{- 16}$	4.16334 × $10^{- 16}$
0.8	3.60822 × $10^{- 16}$	1.08247 × $10^{- 15}$	4.99649 × $10^{- 16}$
0.9	4.44089 × $10^{- 16}$	1.05471 × $10^{- 15}$	6.38378 × $10^{- 16}$

Table 3

Comparison of AEs of Example 1 at $a = 1, b = 1$ .

t	Method in [52]	Our Method at $N = 14$
0	2.9762 × $10^{- 3}$	0
1.008	1.0924 × $10^{- 3}$	9.86874 × $10^{- 8}$
2.016	3.9923 × $10^{- 4}$	5.27874 × $10^{- 8}$
3.012	1.4726 × $10^{- 4}$	2.30097 × $10^{- 8}$
4.008	5.4206 × $10^{- 5}$	9.63606 × $10^{- 9}$
5.004	1.9913 × $10^{- 5}$	3.97005 × $10^{- 9}$
6	7.3004 × $10^{- 6}$	1.61747 × $10^{- 9}$
7.008	2.6387 × $10^{- 6}$	6.45425 × $10^{- 10}$
8.004	9.6327 × $10^{- 7}$	2.57282 × $10^{- 10}$
9	3.5087 × $10^{- 7}$	1.08199 × $10^{- 10}$
10.008	1.2596 × $10^{- 7}$	1.90535 × $10^{- 10}$
11.004	4.5672 × $10^{- 8}$	2.68969 × $10^{- 9}$
11.988	1.6717 × $10^{- 8}$	3.06842 × $10^{- 6}$

Table 4

AEs of Example 2 at $N = 2$ .

t	$a = 0.5, b = 4$	$a = - 2, b = 8$	$a = 1, b = 3$
0.1	1.05497 × $10^{- 4}$	3.47492 × $10^{- 4}$	6.66711 × $10^{- 4}$
0.2	3.72070 × $10^{- 4}$	1.34005 × $10^{- 3}$	2.61693 × $10^{- 3}$
0.3	6.50468 × $10^{- 4}$	2.82842 × $10^{- 3}$	5.70139 × $10^{- 3}$
0.4	6.93596 × $10^{- 4}$	4.56552 × $10^{- 3}$	9.67302 × $10^{- 3}$
0.5	1.58985 × $10^{- 4}$	6.20886 × $10^{- 3}$	1.41893 × $10^{- 2}$
0.6	1.38779 × $10^{- 3}$	7.32403 × $10^{- 3}$	1.88159 × $10^{- 2}$
0.7	4.46876 × $10^{- 3}$	7.38901 × $10^{- 3}$	2.30307 × $10^{- 2}$
0.8	9.68835 × $10^{- 3}$	5.79933 × $10^{- 3}$	2.62293 × $10^{- 2}$
0.9	1.77274 × $10^{- 2}$	7.87423 × $10^{- 3}$	2.77310 × $10^{- 2}$

Table 5

AEs of Example 2 at $N = 3$ .

t	$a = 0.5, b = 4$	$a = - 2, b = 8$	$a = 1, b = 3$
0.1	1.08735 × $10^{- 4}$	3.14522 × $10^{- 4}$	1.70807 × $10^{- 4}$
0.2	2.99588 × $10^{- 4}$	1.01358 × $10^{- 3}$	5.05742 × $10^{- 4}$
0.3	4.43908 × $10^{- 4}$	1.80481 × $10^{- 3}$	8.12949 × $10^{- 4}$
0.4	5.10886 × $10^{- 4}$	2.49366 × $10^{- 3}$	9.98418 × $10^{- 4}$
0.5	5.65089 × $10^{- 4}$	2.98097 × $10^{- 3}$	1.06351 × $10^{- 3}$
0.6	7.63035 × $10^{- 4}$	3.25955 × $10^{- 3}$	1.10154 × $10^{- 3}$
0.7	1.34886 × $10^{- 3}$	3.40979 × $10^{- 3}$	1.29345 × $10^{- 3}$
0.8	2.64908 × $10^{- 3}$	3.59449 × $10^{- 3}$	1.90254 × $10^{- 3}$
0.9	5.06659 × $10^{- 3}$	4.05282 × $10^{- 3}$	3.26851 × $10^{- 3}$

Table 6

AEs of Example 2 at $N = 12$ .

t	$a = 1, b = 1$	$a = 1, b = 3$	$a = 2, b = 1$	$a = 2, b = 2$
0.1	1.9984 × $10^{- 15}$	6.66134 × $10^{- 16}$	8.88178 × $10^{- 16}$	1.55431 × $10^{- 15}$
0.2	3.88578 × $10^{- 15}$	8.88178 × $10^{- 16}$	9.99201 × $10^{- 16}$	2.66454 × $10^{- 15}$
0.3	5.32907 × $10^{- 15}$	1.33227 × $10^{- 15}$	2.22045 × $10^{- 15}$	3.55271 × $10^{- 15}$
0.4	6.99441 × $10^{- 15}$	1.9984 × $10^{- 15}$	3.21965 × $10^{- 15}$	5.10703 × $10^{- 15}$
0.5	7.99361 × $10^{- 15}$	2.22045 × $10^{- 15}$	2.88658 × $10^{- 15}$	5.66214 × $10^{- 15}$
0.6	8.65974 × $10^{- 15}$	2.66454 × $10^{- 15}$	3.44169 × $10^{- 15}$	5.88418 × $10^{- 15}$
0.7	8.21565 × $10^{- 15}$	2.33147 × $10^{- 15}$	3.55271 × $10^{- 15}$	5.66214 × $10^{- 15}$
0.8	8.32667 × $10^{- 15}$	2.33147 × $10^{- 15}$	2.55351 × $10^{- 15}$	5.88418 × $10^{- 15}$
0.9	8.65974 × $10^{- 15}$	2.44249 × $10^{- 15}$	3.77476 × $10^{- 15}$	6.10623 × $10^{- 15}$

Table 7

AEs of Example 2 at $N = 13$ .

t	$a = 2, b = 3$	$a = 3, b = 3$	$a = 3, b = 4$
0.1	1.11022 × $10^{- 16}$	1.11022 × $10^{- 16}$	1.11022 × $10^{- 16}$
0.2	0	0	0
0.3	1.11022 × $10^{- 16}$	1.11022 × $10^{- 16}$	1.11022 × $10^{- 16}$
0.4	4.44089 × $10^{- 16}$	4.44089 × $10^{- 16}$	4.44089 × $10^{- 16}$
0.5	3.33067 × $10^{- 16}$	3.33067 × $10^{- 16}$	3.33067 × $10^{- 16}$
0.6	2.22045 × $10^{- 16}$	2.22045 × $10^{- 16}$	2.22045 × $10^{- 16}$
0.7	1.11022 × $10^{- 16}$	1.11022 × $10^{- 16}$	1.11022 × $10^{- 16}$
0.8	0	0	0
0.9	2.22045 × $10^{- 16}$	2.22045 × $10^{- 15}$	2.22045 × $10^{- 16}$

References

1. Boyd, J.P. Chebyshev and Fourier Spectral Methods; Courier Corporation: North Chelmsford, MA, USA, 2001.

2. Hesthaven, J.; Gottlieb, S.; Gottlieb, D. Spectral Methods for Time-Dependent Problems; Cambridge University Press: Cambridge, UK, 2007; Volume 21.

3. Trefethen, L.N. Spectral Methods in MATLAB; SIAM: Philadelphia, PA, USA, 2000; Volume 10.

4. Koshy, T. Fibonacci and Lucas Numbers with Applications; John Wiley & Sons: Hoboken, NJ, USA, 2019.

5. Ozkan, E.; Altun, I. Generalized Lucas polynomials and relationships between the Fibonacci polynomials and Lucas polynomials. Comm. Algebra; 2019; 47, pp. 4020-4030. [DOI: https://dx.doi.org/10.1080/00927872.2019.1576186]

6. Özkan, E.; Taştan, M. On Gauss Fibonacci polynomials, on Gauss Lucas polynomials and their applications. Comm. Algebra; 2020; 48, pp. 952-960. [DOI: https://dx.doi.org/10.1080/00927872.2019.1670193]

7. Du, T.; Wu, Z. Some identities involving the bi-periodic Fibonacci and Lucas polynomials. AIMS Math.; 2023; 8, pp. 5838-5846. [DOI: https://dx.doi.org/10.3934/math.2023294]

8. Haq, S.; Ali, I. Approximate solution of two-dimensional Sobolev equation using a mixed Lucas and Fibonacci polynomials. Eng. Comput.; 2022; 38, pp. 2059-2068. [DOI: https://dx.doi.org/10.1007/s00366-021-01327-5]

9. Mohamed, A.S. Fibonacci collocation pseudo-spectral method of variable-order space-fractional diffusion equations with error analysis. AIMS Math.; 2022; 7, pp. 14323-14337. [DOI: https://dx.doi.org/10.3934/math.2022789]

10. Manohara, G.; Kumbinarasaiah, S. An Innovative Fibonacci Wavelet Collocation Method for the Numerical Approximation of Emden-Fowler Equations. Appl. Numer. Math.; 2024; 201, pp. 347-369.

11. Abd-Elhameed, W.M.; Al-Harbi, A.K.; Alqubori, O.M.; Alharbi, M.H.; Atta, A.G. Collocation method for the time-fractional generalized Kawahara equation using a certain Lucas polynomial sequence. Axioms; 2025; 14, 114. [DOI: https://dx.doi.org/10.3390/axioms14020114]

12. Abd-Elhameed, W.M. Novel expressions for the derivatives of sixth-kind Chebyshev polynomials: Spectral solution of the non-linear one-dimensional Burgers’ equation. Fractal Fract.; 2021; 5, 74. [DOI: https://dx.doi.org/10.3390/fractalfract5020053]

13. Youssri, Y.H.; Abd-Elhameed, W.M.; Abdelhakem, M. A robust spectral treatment of a class of initial value problems using modified Chebyshev polynomials. Math. Methods Appl. Sci.; 2021; 44, pp. 9224-9236. [DOI: https://dx.doi.org/10.1002/mma.7347]

14. Ahmed, H.M. A new first finite class of classical orthogonal polynomials operational matrices: An application for solving fractional differential equations. Contemp. Math.; 2023; 4, pp. 974-994. [DOI: https://dx.doi.org/10.37256/cm.4420232716]

15. Tohidi, E.; Bhrawy, A.H.; Erfani, K. A collocation method based on Bernoulli operational matrix for numerical solution of generalized pantograph equation. Appl. Math. Model.; 2013; 37, pp. 4283-4294. [DOI: https://dx.doi.org/10.1016/j.apm.2012.09.032]

16. Abd-Elhameed, W.M.; Alsuyuti, M.M. New spectral algorithm for fractional delay pantograph equation using certain orthogonal generalized Chebyshev polynomials. Commun. Nonlinear Sci. Numer. Simul.; 2025; 141, 108479. [DOI: https://dx.doi.org/10.1016/j.cnsns.2024.108479]

17. Ahmed, H.M. New generalized Jacobi Galerkin operational matrices of derivatives: An algorithm for solving multi-term variable-order time-fractional diffusion-wave equations. Fractal Fract.; 2024; 8, 68. [DOI: https://dx.doi.org/10.3390/fractalfract8010068]

18. Enns, R.M. Nonlinear Phenomena in Physics and Biology; Springer Science & Business Media: New York, NY, USA, 2012; Volume 75.

19. Hilborn, R.C. Chaos and Nonlinear Dynamics: An Introduction for Scientists and Engineers; Oxford University Press: Oxford, UK, 2000.

20. Magin, R. Fractional calculus in bioengineering, part 1. Crit. Rev. Biomed. Eng.; 2004; 32, pp. 1-104. [DOI: https://dx.doi.org/10.1615/CritRevBiomedEng.v32.10]

21. Abd-Elhameed, W.M.; Ahmed, H.M. Spectral solutions for the time-fractional heat differential equation through a novel unified sequence of Chebyshev polynomials. AIMS Math.; 2024; 9, pp. 2137-2166. [DOI: https://dx.doi.org/10.3934/math.2024107]

22. Wang, F.; Hou, E.; Salama, S.A.; Khater, M.M.A. Numerical investigation of the nonlinear fractional Ostrovsky equation. Fractals; 2022; 30, 2240142. [DOI: https://dx.doi.org/10.1142/S0218348X22401429]

23. Amin, A.Z.; Abdelkawy, M.A.; Solouma, E.; Al-Dayel, I. A spectral collocation method for solving the non-linear distributed-order fractional Bagley–Torvik differential equation. Fractal Fract.; 2023; 7, 780. [DOI: https://dx.doi.org/10.3390/fractalfract7110780]

24. Heydari, M.H.; Razzaghi, M.; Baleanu, D. A numerical method based on the piecewise Jacobi functions for distributed-order fractional Schrödinger equation. Commun. Nonlinear Sci. Numer. Simul.; 2023; 116, 106873. [DOI: https://dx.doi.org/10.1016/j.cnsns.2022.106873]

25. Kovacic, I.; Brennan, M.J. The Duffing Equation: Nonlinear Oscillators and Their Behaviour; John Wiley & Sons: Hoboken, NJ, USA, 2011.

26. Alvarez-Ramirez, J.; Espinosa-Paredes, G.; Puebla, H. Chaos control using small-amplitude damping signals. Phys. Lett. A; 2003; 316, pp. 196-205. [DOI: https://dx.doi.org/10.1016/S0375-9601(03)01147-2]

27. Aghdam, M.M.; Fallah, A. Analytical Solutions for Generalized Duffing Equation. Nonlinear Approaches in Engineering Applications: Advanced Analysis of Vehicle Related Technologies; Springer: Cham, Switzerland, 2016; pp. 263-278.

28. El-Sayed, A.A. Pell-Lucas polynomials for numerical treatment of the nonlinear fractional-order Duffing equation. Demonstr. Math.; 2023; 56, 20220220. [DOI: https://dx.doi.org/10.1515/dema-2022-0220]

29. Wang, Y.R.; Chen, G.W. Predicting multiple numerical solutions to the Duffing equation using machine learning. Appl. Sci.; 2023; 13, 10359. [DOI: https://dx.doi.org/10.3390/app131810359]

30. Kamiński, M.; Corigliano, A. Numerical solution of the Duffing equation with random coefficients. Meccanica; 2015; 50, pp. 1841-1853. [DOI: https://dx.doi.org/10.1007/s11012-015-0133-0]

31. Elías-Zúñiga, A. A general solution of the Duffing equation. Nonlinear Dynam.; 2006; 45, pp. 227-235. [DOI: https://dx.doi.org/10.1007/s11071-006-1858-z]

32. Abdulganiy, R.I.; Wen, S.; Feng, Y.; Zhang, W.; Tang, N. Adapted block hybrid method for the numerical solution of Duffing equations and related problems. AIMS Math.; 2021; 6, pp. 14013-14034. [DOI: https://dx.doi.org/10.3934/math.2021810]

33. Geng, F. Numerical solutions of Duffing equations involving both integral and non-integral forcing terms. Comput. Math. Appl.; 2011; 61, pp. 1935-1938. [DOI: https://dx.doi.org/10.1016/j.camwa.2010.07.053]

34. Tabatabaei, K.; Gunerhan, E. Numerical solution of Duffing equation by the differential transform method. Appl. Math. Inf. Sci. Lett.; 2014; 2, pp. 1-6.

35. Salas, A.H.; Castillo, J.E. Exact solutions to cubic Duffing equation for a nonlinear electrical circuit. Visión Electrónica; 2014; 8, pp. 46-53.

36. Singh, H.; Srivastava, H.M. Numerical investigation of the fractional-order Liénard and Duffing equations arising in oscillating circuit theory. Front. Phys.; 2020; 8, 120. [DOI: https://dx.doi.org/10.3389/fphy.2020.00120]

37. Kim, V.A.; Parovik, R.I.; Rakhmonov, Z.R. Implicit finite-difference scheme for a Duffing oscillator with a derivative of variable fractional order of the Riemann-Liouville Type. Mathematics; 2023; 11, 558. [DOI: https://dx.doi.org/10.3390/math11030558]

38. Kim, V.A.; Parovik, R.I. Some aspects of the numerical analysis of a fractional duffing oscillator with a fractional variable order derivative of the Riemann-Liouville type. AIP Conf. Proc.; 2022; 2467, 060014.

39. Canuto, C.; Hussaini, M.Y.; Quarteroni, A.; Zang, T.A. Spectral Methods in Fluid Dynamics; Springer: Berlin/Heidelberg, Germany, 1988.

40. Alsuyuti, M.M.; Doha, E.H.; Ezz-Eldien, S.S. Galerkin operational approach for multi-dimensions fractional differential equations. Commun. Nonlinear Sci. Numer. Simul.; 2022; 114, 106608. [DOI: https://dx.doi.org/10.1016/j.cnsns.2022.106608]

41. Atta, A.G.; Abd-Elhameed, W.M.; Moatimid, G.M.; Youssri, Y.H. Shifted fifth-kind Chebyshev Galerkin treatment for linear hyperbolic first-order partial differential equations. Appl. Numer. Math.; 2021; 167, pp. 237-256. [DOI: https://dx.doi.org/10.1016/j.apnum.2021.05.010]

42. Hafez, R.M.; Youssri, Y.H. Fully Jacobi–Galerkin algorithm for two-dimensional time-dependent PDEs arising in physics. Int. J. Mod. Phys. C; 2024; 35, 2450034. [DOI: https://dx.doi.org/10.1142/S0129183124500347]

43. Atta, A.G.; Abd-Elhameed, W.M.; Moatimid, G.M.; Youssri, Y.H. Modal shifted fifth-kind Chebyshev tau integral approach for solving heat conduction equation. Fractal Fract.; 2022; 6, 619. [DOI: https://dx.doi.org/10.3390/fractalfract6110619]

44. El-Sayed, A.A.; Boulaaras, S.; Sweilam, N.H. Numerical solution of the fractional-order logistic equation via the first-kind Dickson polynomials and spectral tau method. Math. Meth. Appl. Sci.; 2023; 46, pp. 8004-8017. [DOI: https://dx.doi.org/10.1002/mma.7345]

45. Tu, H.; Wang, Y.; Yang, C.; Liu, W.; Wang, X. A Chebyshev–Tau spectral method for coupled modes of underwater sound propagation in range-dependent ocean environments. Phys. Fluids; 2023; 35, 037113. [DOI: https://dx.doi.org/10.1063/5.0138012]

46. Mostafa, D.; Zaky, M.A.; Hafez, R.M.; Hendy, A.S.; Abdelkawy, M.A.; Aldraiweesh, A.A. Tanh Jacobi spectral collocation method for the numerical simulation of nonlinear Schrödinger equations on unbounded domain. Math. Meth. Appl. Sci.; 2023; 46, pp. 656-674. [DOI: https://dx.doi.org/10.1002/mma.8538]

47. Weera, W.; Kumar, R.S.V.; Sowmya, G.; Khan, U.; Prasannakumara, B.C.; Mahmoud, E.E.; Yahia, I.S. Convective-radiative thermal investigation of a porous dovetail fin using spectral collocation method. Ain Shams Eng. J.; 2023; 14, 101811. [DOI: https://dx.doi.org/10.1016/j.asej.2022.101811]

48. Atta, A.G. Two spectral Gegenbauer methods for solving linear and nonlinear time fractional Cable problems. Int. J. Mod. Phys. C; 2023; 35, 2450070. [DOI: https://dx.doi.org/10.1142/S0129183124500700]

49. Abd-Elhameed, W.M.; Alqubori, O.M.; Atta, A.G. A collocation procedure for treating the time-fractional FitzHugh–Nagumo differential equation using shifted Lucas polynomials. Mathematics; 2024; 12, 3672. [DOI: https://dx.doi.org/10.3390/math12233672]

50. Abd-Elhameed, W.M.; Alqubori, O.M.; Atta, A.G. A collocation procedure for the numerical treatment of FitzHugh–Nagumo equation using a kind of Chebyshev polynomials. AIMS Math.; 2025; 10, pp. 1201-1223. [DOI: https://dx.doi.org/10.3934/math.2025057]

51. Lyapin, A.P.; Akhtamova, S.S. Recurrence relations for the sections of the generating series of the solution to the multidimensional difference equation. Vestnik Udmurtskogo Universiteta Matematika Mekhanika Komp’yuternye Nauki; 2021; 31, pp. 414-423. [DOI: https://dx.doi.org/10.35634/vm210305]

52. Pirmohabbati, P.; Sheikhani, A.H.R.; Najafi, H.S.; Ziabari, A.A. Numerical solution of full fractional Duffing equations with Cubic-Quintic-Heptic nonlinearities. AIMS Math.; 2020; 5, pp. 1621-1641. [DOI: https://dx.doi.org/10.3934/math.2020110]

53. Novozhenova, O.G. Life and science of Alexey Gerasimov, one of the pioneers of fractional calculus in Soviet Union. Fract. Calc. Appl. Anal.; 2017; 20, pp. 790-809. [DOI: https://dx.doi.org/10.1515/fca-2017-0040]

54. Caputo, M.; Fabrizio, M. On the notion of fractional derivative and applications to the hysteresis phenomena. Meccanica; 2017; 52, pp. 3043-3052. [DOI: https://dx.doi.org/10.1007/s11012-017-0652-y]

55. Caputo, M. Linear models of dissipation whose Q is almost frequency independent—II. Geophys. J. Int.; 1967; 13, pp. 529-539. [DOI: https://dx.doi.org/10.1111/j.1365-246X.1967.tb02303.x]

56. Gerasimov, A.N. Generalization of linear deformation laws and their application to internal friction problems. Appl. Math. Mech.; 1948; 12, pp. 529-539.

57. Hafez, R.M.; Zaky, M.A.; Abdelkawy, M.A. Jacobi spectral Galerkin method for distributed-order fractional Rayleigh–Stokes problem for a generalized second grade fluid. Front. Phys.; 2020; 7, 240. [DOI: https://dx.doi.org/10.3389/fphy.2019.00240]

© 2025 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.

Numerical Solutions for Nonlinear Ordinary and Fractional Duffing Equations Using Combined Fibonacci–Lucas Polynomials

Content area

Abstract

Full text

1. Introduction

2. Introducing a Unified Sequence of Fibonacci and Lucas Polynomials

2.1. Fibonacci and Lucas Polynomial Sequences

2.2. Derivatives and Operational Matrices of the Fibonacci–Lucas Polynomials

3. A Matrix Collocation Approach for the Nonlinear Second-Order Duffing Equation

4. A Matrix Collocation Approach for the Nonlinear Fractional-Order Duffing Equation

4.1. The Operational Matrix of Fractional Derivatives for FLia,b(t)

4.2. Collocation Algorithm for the NFDE

5. Error Bound

6. Illustrative Examples

7. Concluding Remarks

4.1. The Operational Matrix of Fractional Derivatives for $F L_{i}^{a, b} (t)$