# Optimal Homotopy Asymptotic Method-Least Square for Solving Nonlinear Fractional-Order Gradient-Based Dynamic System from an Optimization Problem

Authors: Oluwaseun Olumide Okundalaye; Wan Ainun Mior Othman; Nallasamy Kumaresan
Journal: Advances in Mathematical Physics
 (2020)

---

## Abstract

In this paper, we consider an approximate analytical method of optimal homotopy asymptotic method-least square (OHAM-LS) to obtain a solution of nonlinear fractional-order gradient-based dynamic system (FOGBDS) generated from nonlinear programming (NLP) optimization problems. The problem is formulated in a class of nonlinear fractional differential equations, (FDEs) and the solutions of the equations, modelled with a conformable fractional derivative (CFD) of the steepest descent approach, are considered to find the minimizing point of the problem. The formulation extends the integer solution of optimization problems to an arbitrary-order solution. We exhibit that OHAM-LS enables us to determine the convergence domain of the series solution obtained by initiating convergence-control parameterCj′s. Three illustrative examples were included to show the effectiveness and importance of the proposed techniques.


---

## Body

## 1. Introduction

Consider a nonlinear programming-constrained optimization problems (NLPCOPs) of the form(1)minx∈Rnfxsubjecttogkx≤0andhkx=0∀k∈I=1,2..m,where f:Rn⟶R,hk:Rn⟶R, and gk:Rn⟶R, k, are C2 functions. Let X0=x∈Rn∣hk=0,gk≤0,i∈I be the feasible set of Equation (1), and we assume that X0 is not empty. The general idea of obtaining an approximate analytical solution to Equation (1) is to transform to an unconstrained nonlinear programming problem by any suitable technique such as augmented Lagrange method, barrier method, and penalty method [1, 2]; it can then be solved by any unconstrained optimization numerical method like the steepest descent method, conjugate gradient method, and Newton method. In optimization, the penalty method is the most efficient method to transform a constrained optimization problem into an unconstrained optimization problem [3–5]. An efficient penalty function for equality and inequality problem Equation (1) is given below
(2)Ppenaltyhkx=μ1σ∑i=1phkxσ,(3)Ppenaltygkx=μ∑i=1pmax0,gkxσ,where σ=2. It can be seen that under some conditions, the solutions to Equation (1) are solutions of the unconstrained below [6],
(4)minFx,μ=fx+μ1σ∑iphkxσ+∑i=1pmax0,gkxσ,subjecttox∈Rn,where μ>0 is an auxiliary penalty variable. The corollary connecting the minimizer of the constraint problem in Equation (1) and unconstrained problem in Equation (4) is seen in [7]. The gradient descent method as a standard optimization algorithm has been widely applied in many engineering applications, such as optimization machine learning and image [8–10]. Through diverse research and studies, it is established that the gradient method is one of the most reliable and efficient ways to find the optimal solution of optimization problems [11]. Nowadays, one of the critical points of the gradient method is how to improve the performance further. As an important area of mathematics, fractional calculus is believed to be an excellent tool to enhance the old gradient descent method, mainly because of its special long memory characteristics and nonlocality [12–14]. In the past decade, several methods have been considered to solve unconstrained nonlinear optimization in the form of ordinary differential equation (ODE) dynamic system of which the gradient-based method is one of the approaches. The technique transforms the nonlinear optimization problem to an ODE dynamic system with some optimality conditions, to obtain optimal solutions to the optimization problem. The gradient-based method was first proposed by [15], was developed by [16, 17], and was later extended to solve differential nonlinear programming problems [18]. However, the studies of nonlinear fractional-order gradient-based dynamic systems are still in the infant stage and are considered further in this paper.Arbitrary-order ODEs, which are the generalizations of integer-order ODEs, are mostly used to model problems in applied sciences. Several numerical methods had been used to solve linear and nonlinear problems of FDEs, such as the Adomian decomposition method (ADM) [19], variational iteration method (VIM) [20], homotopy perturbation method for solving fractional Zakharov-Kuznetsov equation [21], a numerical method for FDEs [22], and multivariate padé approximation (MPA) [23]. The usefulness of an arbitrary-order started receiving tremendous attention of researchers in the field of applied science and engineering in the last two decades where some authors in the area of optimization focused on developing approximate analytical methods for different types of nonlinear constrained optimization problems in the form of IVPs of nonlinear FDE systems including multistage ADM for NLP [24], a fractional dynamics trajectory approach [25], the convergence of HAM and application [26], fractional steepest descent approach [27], studied optimal solution of fractional gradient [28], gradient descent direction with Caputo derivative sense for BP neural networks [29], fractional-order gradient methods [30], and conformable fractional gradient-based system [31]. In 2008, Marinca and Herisanu [32] introduced a numerical method called OHAM to solve a nonlinear problem, later extended by Azimi et al. [33] for strong nonlinear differential equations (NLDEs). This powerful tool called OHAM has not been applied in the area of FOGBDS, which motivates this work.So, in this paper, we showed that the steady-state solutionsxt of the proposed system can be approximated analytically to the expected exact optimal solution x∗ of the nonlinear programming constrained optimization problem by OHAM-LS as t⟶∞. The significant contribution is summarized as follows:
(1)
The reason why OHAM-LS is preferable to be the method used [25, 31] to solve FOGBDS(2)
The reason why some existing approximate analytical method cannot guarantee the convergence of the series solution is discussed(3)
From the previous approximation analytical method of solving FOGBDS, accurate optimal values control-convergence parameter had been a little bit difficult to achieve which is easily address with least square optimization techniques(4)
OHAM-LS with guaranteed convergence ability is proposed with conformable fractional derivative sense to solve FOGBDS. The fastest convergence ability of the proposed compared with fourth-order Runge-Kutta is also shownWe arrange the paper as follows: a brief introduction to the fractional calculus and OHAM-LS derivation is given in Section2. Section 3 is devoted to problem formulation of OHAM-LS with FOGBDS and the key contributions. In Section 4, we solved some NLP constrained optimization problems to show the effectiveness of the proposed method. The results obtained from OHAM-LS are plotted in several figures with numerical method comparisons to confirm the validity and ability of the method to solve the problem. In the last section are the conclusions.


## 2. Preliminaries

### 2.1. Fractional Calculus

The most common arbitrary-order in literature is the Riemann-Liouville’s and the Caputo fractional derivative. The arbitrary-order definitions are generally used for mathematical modelling within many areas, especially when the classical-order derivative operator fails or additional memory effect is required. However, the limitation of these two definitions is that they do not provide some of the features that the classical derivative provides, such as chain rule, quotient rule, product rule, and derivative of constant. Recently, Khalil et al. [34] have characterized a new fractional derivative operator, which is an extension of the usual conformable fractional derivative, to overcome these deficiencies. Besides these advantages, the conformable fractional derivative does not show the memory effect, which is inherent for the other classical fractional derivatives.Definition 1.
Letf:0,∞⟶R be a given function. The αth order CFD of f given by
(5)Tαfx=limε→0fx+εx1−α−fxε,
∀x>0 and α∈0,1
This new definition preserves many properties of the classical derivatives refer to [34, 35]. Some features that we will adopt are as follows:Theorem 2.
. Let0<α≤1 and f,g be α-differentiable at a point x>0; if f is a differentiable function, then dαf/dxα=x1−αdf/dx.Definition 3.
Iaαfx=Ia1xα−1f=∫axft/t1−αdt, where the integral is the regular Riemann improper integral, and α∈0,1.Theorem 4.
Letf be any continuous function in the domain of Iα, then TαIaαfx=fx∀x≥a.Theorem 5.
Letf:a,b⟶R be differentiable and 0<α≤1. Then, for all x>a, we have IaαTαfx=fx−fa.


### 2.2. The Elementary Concepts of OHAM-LS

We start from the fundamental principle of OHAM as described in [36–38]. Consider the IVPs
(6)Lizit+Nizit+git=0t∈φi=1,2,⋯m,with initial conditions
(7)zib=ai,where Li is a linear operator, Ni is a nonlinear operator, t is an independent variable, zit is an unknown function, φ is the problem domain, and git is a known function. According to OHAM, one can construct an homotopy map Hiϕit,p:φ×0,1⟶φ which satisfies
(8)1−pLiϕit,p+git=HipLiϕit,p+Niϕit,p+git,where p∈0,1 is an embedding parameter, Hip is a nonzero auxiliary function for p≠0, H0=0, and ϕit,p is an unknown function. Obviously, when p=0 and p=1, it holds that ϕit,0=zi,0t and ϕit,1=zit, respectively. Thus, as p varies from 0 to 1, the solution ϕit,p approaches from zi,0t to zit where zi,0t is the initial guess that satisfies the linear operator which is obtained from Equation (8) for p=0 as
(9)Lizi,0t+git=0.zi,0b=0.Hip is chosen in the form
(10)Hip=pC1+p2C2+p3C3+⋯j=1,2,⋯n,where Cj would be determined in the last part of this work. We consider Equation (8) in the form
(11)ϕit,p,Cj=zi,0t+∑k≥1zi,kt,Cjpkj=1,2,⋯n.Now substituting Equation (11) in Equation (8) and equating the coefficient of like power of p, we obtain the governing equation of zi,0t in a linear form, given in Equation (9). The first- and second-order problems are given by
(12)Lizi,1t+git=C1N0zi,0t,zi,1b=0,(13)Lizi,2t−Lizi,1t=C2Ni,0zi,0t+C1Lizi,1t+Ni,1zi,1t,zi,2b=0,and the general governing equations for zi,kt are given by
(14)Lizi,kt−Lizi,k−1t=CkNi,0zi,0t+∑m=1k−1Cj,mLizi,k−mt+Ni,k−mzi,k−1t,(15)zi,kb=0,k=2,3,⋯,where Ni,mz0t,zi,1t⋯,zi,mt is the coefficient of pm, obtained by expanding Niϕit,p,Cj in series with respect to the embedding parameter p(16)Niϕit,p,Cj=Ni,0zi,0t+∑m=1∞Ni,mzi,mtpm,where ϕit,p,Cj is obtained from Equation (11). It should noted that zi,k for k≥0 is governed by the linear Equations (9), (12), and (14) with linear initial conditions that come from the original problem, which can be easily solved.It has been shown that the convergence of the series Equation (16) depends upon the Cj. If it is convergent at p=1, we have
(17)zit,Cj=zi,0t+∑k≥1zi,kt,Cj,The result of themth-order approximation is given as
(18)z~it,Cj=zi,0t+∑k=1mzi,kt,Cj,Substituting Equation (18) in Equation (6), we get the following expression for the residual
(19)Rit,Cj=Liz~it,Cj+Niz~it,Cj+git,i=1,2,⋯m,IfRit,Cj=0, then z~it,Cj is the exact solution. Usually, such a case does not arise for nonlinear problems. Several methods [39, 40] can be used to find the optimal values of convergence-control parameters Cj′s like the method of the least square method, collocation method, Ritz method, and Galerkin’s method. By applying the least square method, we have minimized the functional
(20)JkC1,C2,C3,⋯Cm=∫abRk2t,C1,C2,C3,⋯Cmdt,where the value a and b depends on the given problem.
(21)ψk=∂JkCk∂Ck=0,k=1,2,⋯m.With these knownCj′s, the approximate solution (of mth-order) is well determined.The correctness of the method by(1)
Error NormL2.
(22)L2=Zexact−ZN≈b−aN∑i=0Nziexact−zNi2(2)
Error NormL∞.
(23)L∞=Zexact−ZN∞≈maxiZiexact−ZNiThe OHAM-LS is based on hybridization of OHAM with the least square method of optimization technique. The OHAM enable us to determine the convergence domain of the series solution, and the least square method allows us to obtain the optimal values of theCks.Remark 6.
OHAM-LS is preferable because VIM, HPM, and HAM are just a case as proved by [41–44].Remark 7.
The existing approximate analytical for FOGBDS cannot guarantee convergence mainly because they possess no criteria for the establishment for convergence of the series solution Equations (20) and (21).


## 2.1. Fractional Calculus

The most common arbitrary-order in literature is the Riemann-Liouville’s and the Caputo fractional derivative. The arbitrary-order definitions are generally used for mathematical modelling within many areas, especially when the classical-order derivative operator fails or additional memory effect is required. However, the limitation of these two definitions is that they do not provide some of the features that the classical derivative provides, such as chain rule, quotient rule, product rule, and derivative of constant. Recently, Khalil et al. [34] have characterized a new fractional derivative operator, which is an extension of the usual conformable fractional derivative, to overcome these deficiencies. Besides these advantages, the conformable fractional derivative does not show the memory effect, which is inherent for the other classical fractional derivatives.Definition 1.
Letf:0,∞⟶R be a given function. The αth order CFD of f given by
(5)Tαfx=limε→0fx+εx1−α−fxε,
∀x>0 and α∈0,1
This new definition preserves many properties of the classical derivatives refer to [34, 35]. Some features that we will adopt are as follows:Theorem 2.
. Let0<α≤1 and f,g be α-differentiable at a point x>0; if f is a differentiable function, then dαf/dxα=x1−αdf/dx.Definition 3.
Iaαfx=Ia1xα−1f=∫axft/t1−αdt, where the integral is the regular Riemann improper integral, and α∈0,1.Theorem 4.
Letf be any continuous function in the domain of Iα, then TαIaαfx=fx∀x≥a.Theorem 5.
Letf:a,b⟶R be differentiable and 0<α≤1. Then, for all x>a, we have IaαTαfx=fx−fa.


## 2.2. The Elementary Concepts of OHAM-LS

We start from the fundamental principle of OHAM as described in [36–38]. Consider the IVPs
(6)Lizit+Nizit+git=0t∈φi=1,2,⋯m,with initial conditions
(7)zib=ai,where Li is a linear operator, Ni is a nonlinear operator, t is an independent variable, zit is an unknown function, φ is the problem domain, and git is a known function. According to OHAM, one can construct an homotopy map Hiϕit,p:φ×0,1⟶φ which satisfies
(8)1−pLiϕit,p+git=HipLiϕit,p+Niϕit,p+git,where p∈0,1 is an embedding parameter, Hip is a nonzero auxiliary function for p≠0, H0=0, and ϕit,p is an unknown function. Obviously, when p=0 and p=1, it holds that ϕit,0=zi,0t and ϕit,1=zit, respectively. Thus, as p varies from 0 to 1, the solution ϕit,p approaches from zi,0t to zit where zi,0t is the initial guess that satisfies the linear operator which is obtained from Equation (8) for p=0 as
(9)Lizi,0t+git=0.zi,0b=0.Hip is chosen in the form
(10)Hip=pC1+p2C2+p3C3+⋯j=1,2,⋯n,where Cj would be determined in the last part of this work. We consider Equation (8) in the form
(11)ϕit,p,Cj=zi,0t+∑k≥1zi,kt,Cjpkj=1,2,⋯n.Now substituting Equation (11) in Equation (8) and equating the coefficient of like power of p, we obtain the governing equation of zi,0t in a linear form, given in Equation (9). The first- and second-order problems are given by
(12)Lizi,1t+git=C1N0zi,0t,zi,1b=0,(13)Lizi,2t−Lizi,1t=C2Ni,0zi,0t+C1Lizi,1t+Ni,1zi,1t,zi,2b=0,and the general governing equations for zi,kt are given by
(14)Lizi,kt−Lizi,k−1t=CkNi,0zi,0t+∑m=1k−1Cj,mLizi,k−mt+Ni,k−mzi,k−1t,(15)zi,kb=0,k=2,3,⋯,where Ni,mz0t,zi,1t⋯,zi,mt is the coefficient of pm, obtained by expanding Niϕit,p,Cj in series with respect to the embedding parameter p(16)Niϕit,p,Cj=Ni,0zi,0t+∑m=1∞Ni,mzi,mtpm,where ϕit,p,Cj is obtained from Equation (11). It should noted that zi,k for k≥0 is governed by the linear Equations (9), (12), and (14) with linear initial conditions that come from the original problem, which can be easily solved.It has been shown that the convergence of the series Equation (16) depends upon the Cj. If it is convergent at p=1, we have
(17)zit,Cj=zi,0t+∑k≥1zi,kt,Cj,The result of themth-order approximation is given as
(18)z~it,Cj=zi,0t+∑k=1mzi,kt,Cj,Substituting Equation (18) in Equation (6), we get the following expression for the residual
(19)Rit,Cj=Liz~it,Cj+Niz~it,Cj+git,i=1,2,⋯m,IfRit,Cj=0, then z~it,Cj is the exact solution. Usually, such a case does not arise for nonlinear problems. Several methods [39, 40] can be used to find the optimal values of convergence-control parameters Cj′s like the method of the least square method, collocation method, Ritz method, and Galerkin’s method. By applying the least square method, we have minimized the functional
(20)JkC1,C2,C3,⋯Cm=∫abRk2t,C1,C2,C3,⋯Cmdt,where the value a and b depends on the given problem.
(21)ψk=∂JkCk∂Ck=0,k=1,2,⋯m.With these knownCj′s, the approximate solution (of mth-order) is well determined.The correctness of the method by(1)
Error NormL2.
(22)L2=Zexact−ZN≈b−aN∑i=0Nziexact−zNi2(2)
Error NormL∞.
(23)L∞=Zexact−ZN∞≈maxiZiexact−ZNiThe OHAM-LS is based on hybridization of OHAM with the least square method of optimization technique. The OHAM enable us to determine the convergence domain of the series solution, and the least square method allows us to obtain the optimal values of theCks.Remark 6.
OHAM-LS is preferable because VIM, HPM, and HAM are just a case as proved by [41–44].Remark 7.
The existing approximate analytical for FOGBDS cannot guarantee convergence mainly because they possess no criteria for the establishment for convergence of the series solution Equations (20) and (21).


## 3. Construction of OHAM-LS with FOGBDS Generated by NLPCOPs

We begin by considering a NLP constrained in the form(24)minx∈Rnfxsubjecttogkx≤0andhkx=0∀k∈I=1,2..m,where f:Rn⟶R is the objective function, hkx:Rn⟶R are equality constraint functions, gkx:Rn⟶R are inequality constraint functions, and C2 are continuous differentiable functions. One of the main ideas of solving unconstrained NLP is by searching for the next point by choosing proper search direction dk and the stepsize αk as in the Newton direction [45], trust-region algorithm for unconstrained optimization [46]; the descent method [47], conjugate gradient method [48], three-term conjugate gradient method [49], and subspace method for nonlinear optimization [50]; the hybrid method for convex NLP [51]; CCM for optimization problem and application [52]; and descent direction stochastic approximation for optimization problem [53]. But there are studies for other approaches. In this paper, we obtain the minimizing point of the problem by solving a certain initial-value system of FDEs. This kind of FOGBDS was first proposed by Evirgen and Özdemir [24].Using the penalty function Equation (2) and (3) for Equation (24) with ρ=2, the conformable FOGBDS model can be constructed as
(25)Tαxt=−∇xFx,μ,subject to the initial conditions
(26)xk0=xk0,k=1⋯.,m.where ∇xFx,μ is the gradient vector of Equation (25) with respect to xk∈Rn and Tα is the CFD of 0<α≤1.Note that a pointxe is called an equilibrium point of Equation (25) if it satisfies the RHS of Equation (25). We reformulate fractional dynamic system Equation (25) as
(27)Tαxkt=gkt,μ,x1,x2⋯.xn,k=1,2⋯,m.We used OHAM-LS to obtain the solution of system Equation (27) by constructing the following homotopy
(28)Tαxkt=pgkt,μ,x1,x2⋯.xn,where k=1,2⋯,n and p∈0,1. If p=0, Equation (28) becomes
(29)Tαxkt=0,and when p=1, the homotopy Equation (28) becomes
(30)Tαxkt=gkt,μ,x1,x2⋯.xn,k=1,2⋯.m,t∈0,1,0<α≤1,subject to the initial conditions,
(31)xkb=ak,k=1,2⋯,m.The correction functional for the system of conformable fractional nonlinear differential equation Equation (30), according to OHAM-LS, can be constructed as
(32)1−pTαφkt,p=HkpTαφkt,q+Nφkt,q+gkt,μ,φ1t,q,φ2t,q..⋯φnt,q,Thus asp varies from 0 to 1, the solution φkt,p approaches from xk,0t to xkt where xk,0t is the initial guess that satisfies the linear operator which is obtained from Equation (32) for p=0 as
(33)Tαxk,0t=0,xk,0b=0.Hkp is chosen in the form
(34)Hkp=pC1+p2C2+p3C3⋯,where Cj can be determined later. We get an approximate solution by expanding φkt,p,Cj in Taylor’s series with respect to p; we have
(35)φkt,p,Cj=xi,0t+∑i≥1xk,it,Cjpi,j=1,2,⋯,n.Now using Equation (35) in Equation (32) and equating the coefficient of like power of p, we obtain the governing equation of xi,0t in a linear form, given in Equation (33). The 1st- and 2nd-order problems are given by
(36)Tαxk,1t+gkt=C1N0xk,0t,xk,1b=0,Tαxk,2t−Tαxk,1t=C2Nk,0xk,0t+C1Tαxk,1t+Nk,1xk,1t,xk,2b=0and the general governing equations for xk,it are given by
(37)Tαxk,it−Tαxk,i−1t=CiNk,0xk,0t+∑m=1i−1Cj,mTαxk,i−mt+Nk,i−mxk,i−1t,xk,ib=0,i=2,3,⋯m,where Nk,mx0t,xk,1t⋯,xk,mt is the coefficient of pm, obtained by expandingNkφkt,p,Cj in series with respect to p.
(38)Nkφkt,p,Cj=Nk,0xk,0t+∑m≥1Nk,mx0,x1,⋯xmpm.It has been shown that the convergence of the series Equation (38) depends upon the Cj. If it is convergent at p=1, one has
(39)xkt,Cj=xk,0t+∑i≥1xk,it,Cj.The solution of Equation (30) is determined approximately in the form,
(40)x~kt,Cj=xk,0t+∑i=1mxk,it,Cj,j=1,2,⋯.,n.Substituting Equation (40) in Equation (30), we get the following expression for the residual error
(41)Rkt,Cj=Tαx~kt,Cj+Nx~kt,Cj+gkx~kt,Cj.IfRkt,Cj=0, then x~kt,Cj is the exact solution. Usually, such a case does not arise for nonlinear problems. Using the least square method as below minimizes the functional
(42)JkC1,C2,C3,⋯Cm=∫abRk2t,C1,C2,C3,⋯Cmdt,where the value of a and b depends on the given problem.
(43)ψk=∂JkCk∂Ck=0,k=1,2,⋯m.With these knownCk, the analytical approximate solution (of mth-order) is well determined.The steps for optimal homotopy asymptotic method-least square (OHAM-LS) are as follows:Step 1.
We transform the nonlinear constrained optimization problem to the unconstrained optimization problem by a penalty method.Step 2.
We find the gradient of the unconstrained optimization problem, with given initial conditions.Step 3.
We choose the linear and nonlinear operators for OHAM-LS.Step 4.
We construct homotopy for the conformable fractional nonlinear differential equation which includes embedding parameter, auxiliary function, and the unknown function.Step 5.
We substitute the series solution results into the governing equation and equate to zero for an exact solution. Usually, such case a does not arise in nonlinear problems.Step 6.
We find the optimal values forCj′s by using the optimization method called least square method, for good analytical approximate solution.


### 3.1. Convergence Analysis of OHAM-LS with FOGBDS

Theorem 8.
As long as the seriesx~kt,Cj=xk,0t+∑i=1mxk,it,Cj,j=1,2,⋯.,n converges where x~kt,Cj is governed by Equation (40) under the definitions Equations (37) and (38), it must be the solution of Equations (25) and (26).Proof.
If we assume∑m=1∞x~k,mt,Cj, k=1,2⋯.n, converges to x~kt,Cj, then
(44)limm→∞x~k,mt,Cj=0∀k=1,2⋯n.
From Equation (37), we can write
(45)∑i=1∞CiNk,0xk,0t+∑m=1i−1Cj,mTαxk,i−mt+Nk,i−mxk,i−1t=∑i=1∞Tαxk,it−Tαxk,i−1t=limn→∞∑i=1nTαxk,it−Tαxk,i−1t=Tαx11t+Tαx22t−Tαx21t+..+Tαxnnt−Tαxnn−1t=Tαlimn→∞∑m=1nxnnt=Tαlimn→∞xnnt=0,∑i=1∞CiNk,0xk,0t+∑m=1i−1Cj,mTαxk,i−mt+Nk,i−mxk,i−1t=0.
So, using above gives(46)=∑m=1∞Tαxkm−1+Nxkm−1+gkt,μ,x1m−1,x2m−1⋯..xnm−1,(47)=∑m=1∞Tαxkm−1+∑m=1∞Nxkm−1+∑m=1∞gkt,μ,xkm−1,(48)=Tαx~kt,Cj+Nx~kt,Cj+gkx~kt,Cj.
From Equation (48), we have
(49)Tαx~kt,Cj+Nx~kt,Cj+gkx~kt,Cj=0∀k=1,2⋯m,


## 3.1. Convergence Analysis of OHAM-LS with FOGBDS

Theorem 8.
As long as the seriesx~kt,Cj=xk,0t+∑i=1mxk,it,Cj,j=1,2,⋯.,n converges where x~kt,Cj is governed by Equation (40) under the definitions Equations (37) and (38), it must be the solution of Equations (25) and (26).Proof.
If we assume∑m=1∞x~k,mt,Cj, k=1,2⋯.n, converges to x~kt,Cj, then
(44)limm→∞x~k,mt,Cj=0∀k=1,2⋯n.
From Equation (37), we can write
(45)∑i=1∞CiNk,0xk,0t+∑m=1i−1Cj,mTαxk,i−mt+Nk,i−mxk,i−1t=∑i=1∞Tαxk,it−Tαxk,i−1t=limn→∞∑i=1nTαxk,it−Tαxk,i−1t=Tαx11t+Tαx22t−Tαx21t+..+Tαxnnt−Tαxnn−1t=Tαlimn→∞∑m=1nxnnt=Tαlimn→∞xnnt=0,∑i=1∞CiNk,0xk,0t+∑m=1i−1Cj,mTαxk,i−mt+Nk,i−mxk,i−1t=0.
So, using above gives(46)=∑m=1∞Tαxkm−1+Nxkm−1+gkt,μ,x1m−1,x2m−1⋯..xnm−1,(47)=∑m=1∞Tαxkm−1+∑m=1∞Nxkm−1+∑m=1∞gkt,μ,xkm−1,(48)=Tαx~kt,Cj+Nx~kt,Cj+gkx~kt,Cj.
From Equation (48), we have
(49)Tαx~kt,Cj+Nx~kt,Cj+gkx~kt,Cj=0∀k=1,2⋯m,


## 4. Numerical Examples and Results

In this section, three examples are presented to illustrate the efficiency of the new method for solving NLPCOPs. The calculations are performed using maple software 2018, HP ENVY laptop 13 corei7 8th Gen 16GB.Example 1.
Consider the NLPCOP test problem from Schittkowski [54] (No. 216).
(50)Minimizefx=100x12−x22+x1−12,subjecttohx=x1x1−4−2x2+12=0,whose exact solution is not known, but expected optimal solution is x1∗=1.9993,x2∗=3.9998. First, we transform the constraint problem to an unconstrained problem by quadratic penalty function for σ=2; then, we have
(51)fx,μ=100x12−x22+x1−12+12μx1x1−4−2x2+122,where μ∈R+, and so that the nonlinear FOGBDS can be given as
(52)Tαx1t=−400x12−x2x1−2x1−1−μ2x1−4x12−4x1−2x2+12,Tαx2t=200x12−x2+2μx12−4x1−2x2+12,x10=0,x20=0,where 0<α≤1. By using OHAM-LS with auxiliary penalty variable μ=200, the terms of the OHAM-LS solutions for fractional order are acquired by using the concept of homotopy. According to Equation (6)), we choose the linear and nonlinear operators in the following forms:
(53)L1φ1t,p=Tαφ1t,p,L2φ2t,p=Tαφ2t,p,N1φ1t,p=Tαφ1t,p+400φ1t,p2−φ2t,pφ1t,p+2φ1t,p−1+2002φ1t,p−4φ1t,p2−4φ1t,p−2φ2t,p+12,N2φ2t,p=Tαφ2t,p−200φ1t,p2−φ2t,p−400φ1t,p2−4φ1t,p−2φ2t,p+12.We can construct the following homotopy(54)1−pTαφ1t,p=HpTαφ1t,p+400φ1t,p2−φ2t,pφ1t,p+2φ1t,p−1+2002φ1t,p−4φ1t,p2−4φ1t,p−2φ2t,p+12,(55)1−pTαφ2t,p=HpTαφ2t,p−200φ1t,p2−φ2t,p−400φ1t,p2−4φ1t,p−2φ2t,p+12,where
(56)φ1t,p=x1,0t+∑j≤1x1,jtpj,(57)φ2t,p=x2,0t+∑j≤1x2,jtpj,(58)Hkp=pC1+p2C2+p3+C3+⋯,k=1,2⋯m.Substituting Equations (56)-(58) into Equations (54) and (55) and equating the coefficient of the same powers of p result to the following set of linear FDEs.
(59)p0:Tαx1,0t=0,(60)p0:Tαx2,0t=0,(61)p1:Tαx1,1t=2000x1,03C1+Tαx1,0C1−3600x2,0x1,0C1−9600x1,02C1−Tαx1,0+6400x2,0C1+32002x1,0C1−38402C1=0,(62)p1:Tαx2,1t=Tαx2,0C1−1800x1,02C1−Tαx2,0+6400x1,0C1+3400x2,0C1−19200C1=0,(63)p2:Tαx1,2t=2000x1,03C2+6000x1,02x1,1C1+Tαx1,0C2+Tαx1,1C1−3600x2,0x1,0C2−3600x2,0x1,1C1−9600x1,02C2−19200x1,0x1,1C1−3600x1,0x2,1C1−Tαx1,1+6400x2,0C2+32002x1,0C2+32002x1,1C1+6400x2,1C1−38402C2=0,(64)p2:Tαx2,2t=Tαx2,0C2+Tαx2,1C1−1800x1,02C2−3600x1,0x1,1C1−Tαx2,1+6400x1,0C2+6400x1,1C1+3400x2,0C2+3400x2,1C1−19200C2=0.Applying the operatorIα to both sides of Equations (59)-(64) with initial conditions given in Equation (5.6), we obtain
(65)x1,0t=0,(66)x2,0t=0,(67)x1,1t,C1=384020t1/10C1,(68)x2,1t,C1=192000t1/10C1,(69)x1,2t,C1,C2=−6.759104020×1010t1/5C12−384020t1/10+384020t1/10C1+384020C2t1/10,(70)x2,2t,C1,C2=−1.555264000×1010t1/5C12−192000t1/10C12+192000t1/10C1+192000C2t1/10.Adding up the solution components Equations (65)-(70), the 2nd-order approximate solution obtained by OHAM-LS at α=0.9, for p=1, are
(71)x12t,C1,C2=768040C1−384020C12+384020C2t1/10−6.759104020×1010t1/5C12,(72)x22t,C1,C2=384000C1−192000C12+192000C2t1/10−1.555264000×1010t1/5C12.For the calculations ofC1 and C2 in x12t and x22t given in Equations (71) and (72), we apply the procedure mentioned in Equations (19)-(21); we obtain, for x12t,
(73)c1=1.800506863×10−6,c2=6.594892833×10−6,and for x22t,
(74)c1=0.111906918×10−4,c2=0.2190543167×10−4.Substituting these optimal values into Equations (71) and (72) becomes
(75)x12t=4.196444315t1/10−1.084631569t1/5,x22t=7.546421106t1/10−0.7996784175t1/5.Table1 shows the Ck at different values of α for Example 1. Table 2 shows the comparisons and the absolute error between OHAM-LS and RK4 at different values of α=1. Figure 1 shows the analytical approximate solutions obtained by OHAM-LS for α=1,0.9,0.8, and 0.7 with RK4 at α=1.Example 2.
Consider the NLPCOPs test problem from Schittkowski [54] [No 320].
(76)Minimizefx=x1−202+x2+202,subjecttohx=x12100+x224−1=0.Table 1
Control-convergence parametersCk at different values of α.
Variablex1tx1tx2tx2tαC1C2C1C211.912514527×10–67.294797236×10–60.113698734×10–40.229361274×10–40.91.800506863×10–66.594892833×10–60.111906918×10–40.2190543167×10–40.81.714313871×10–65.524430129×10–60.10992623×10–40.201017632×10–40.71.593611093×10–65.294592861×10–60.107191284×10–40.197911283×10–4Table 2
Comparisons and absolute error between OHAM-LS and RK4,α=1.
tkOHAM-LSx1tOHAM-LSx2tRK4x1tRK4x2tErrorx1tErrorx2t0.00000.0000000.0000000.0000000.0000000.0000000.0000000.00051.9702113.8717211.9708993.8718870.0006880.0001660.00101.9771343.9078231.9782743.9079930.0011360.000170.00131.9811023.9224031.9813843.9225540.0002820.0001510.00151.9822113.9305341.9831323.9305780.0009214.4E-050.00201.9832143.9356531.9842523.9356540.0010381E-06Figure 1
(a) Different values ofα (OHAM-LS; α=1, dot; α=0.9, dash; α=0.8, dash dot; and α=0.7, long dash) and Rk4 (α=1, solid) at x1. (b) Different values of α (OHAM-LS; α=1, dot; α=0.9, dash; α=0.8, dash dot; and β=0.7, long dash) and RK4 (α=1, solid) at x2.
(a)(b)This is a practical problem, and the exact solution is not known, but the expected optimal solution isx1∗=9.395, x2∗=−0.6846. First, the quadratic penalty function is used to get the unconstrained optimization problem as follows:
(77)Fx,μ=x1−202+x2+202+12μx12100+x224−12,where μ∈R+ and so that the nonlinear FOGBDS be given as
(78)Tαx1t=2x1−40+μ15000x13+1200x1x22−150x1,Tαx2t=2x2+40+μ1200x2x12+18x23−12x2,0<α≤1,x10=0,x20=0.By using OHAM-LS withμ=106, the terms of the OHAM-LS solutions for fractional order are acquired by using the concept of homotopy. According to Equation (6), we choose the linear and nonlinear operators in the following forms:
(79)L1φ1t,p=Tαφ1t,p,(80)L2φ2t,p=Tαφ2t,p,(81)N1φ1t,p=Tαφ1t,p−2φ1t,p+40−10615000φ1t,p3+1200φ1t,pφ2t,p2−150φ1t,p,(82)N2φ2t,p=Tαφ2t,p−2φ2t,p−40−1061200φ2t,p×φ1t,p2−18φ2t,p3+12φ2t,p,(83)1−pTαφ1t,p=HpTαφ1t,p−2φ1t,p+40−10615000φ1t,p3+1200φ1t,p×φ2t,p2−150φ1t,p,(84)1−pTαφ2t,p=HpTαφ2t,p−2φ2t,p−40−1061200φ2t,p×φ1t,p2−18φ2t,p3+12φ2t,p,where
(85)φ1t,p=x1,0t+∑j≤1x1,jtpj,(86)φ2t,p=x2,0t+∑j≤1x2,jtpj,(87)Hkp=pC1+p2C2+p3+C3+⋯,k=1,2.Substituting Equations (85)-(87) into Equations (83) and (84) and equating the coefficient of the same powers of p yields the following set of linear FDEs:
(88)p0:Tαx1,0t=0,(89)p0:Tαx2,0t=0,(90)p1:Tαx1,1t=−200x1,03C1−5000x1,0x2,02C1+Tαx1,0C1−Tαx1,0+19998x1,0C1+40C1,(91)p1:Tαx2,1t=−125000x2,03C1−5000x2,0x1,02C1−Tαx2,0C1−2x2,0C1−Tβx2,0+500000x1,0C1−40C1=0,(92)p2:Tαx1,2t=−200x1,03C2−6000x1,02x1,1C1−5000x1,0x2,02C2−10000x1,0x2,0x2,1C1−5000x2,02x1,1C1+Tαx1,0C2+Tαx1,1C1−Tαx1,1+19998x1,0C2+19998x1,1C1+40C2=0,(93)p2:Tαx2,2t=−125000x2,03C2−375000x2,02x2,1C1−10000x2,0x1,1x1,0C1−5000x2,0x1,02C2−5000x2,1x1,02C1+Tαx2,1C1+Tαx2,0C2−2x2,0C2−2x2,1C1+500000x1,1C1−Tαx2,1+500000x1,0C2−40C2=0.Applying the operatorIα on both sides of Equations (88)-(93) with initial conditions given in Equation (4.38), we obtain
(94)x1,0t=0,(95)x2,0t=0,(96)x1,1t,C1=−200t1/5C1,(97)x2,1t,C1=200t1/5C1,(98)x1,2t,C1,C2=9.999000×106t2/5C12+200t1/5C12−200t1/5C1−200C2t1/5,(99)x2,2t,C1,C2=2.50001000×108t2/5C12−200t1/5C12+200t1/5C1+200C2t1/10.Adding up the solution components Equations (94)-(99), the 2nd-order approximate solution obtained by OHAM-LS at α=0.9, for p=1, is
(100)x12t=t1/5−400C1−9.999000×106t1/5C12+200C12−200C2,(101)x22t=t1/5400C1−2.50001000×108t1/5C12−200C12+200C2.For the calculations ofC1 and C2 in x12t and x22t given in Equations (100) and (101), we apply the procedure mentioned in Equations (19)-(21), we obtain for x12t,
(102)c1=−0.1208162856×10−3,c2=−0.02826592550.And forx22t,
(103)c1=−1.343994006×10−5,c2=−2.536649902×10−3.Substituting these optimal values into Equations (100) and (101), we have
(104)x12t=5.701514534+0.1459511521t1/5t1/5,x22t=−0.5167059926+0.04515817783t1/5t1/5.Table3 shows the Ck at different values of α for example 2. Table 4 show the comparisons and the absolute error between OHAM-LS and RK4 at different values of α=1. Figure 2 show the analytical approximate solutions obtained by OHAM-LS for α=1,0.9,0.8, and 0.7 with RK4 at α=1.Example 3.
Consider the NLPCOP test problem from Schittkowski [54] (No. 300).
(105)Minimizefx=x12+x22+2x32+x42−5x1−5x2−21x3+7x4,subjectto8−x12−x22−x32−x42−x1+x2−x3+x4≤0,10−x12−2x22−x32+x1+x4≤0,5−2x12−x22−x32−2x1+x2+x4≤0.Table 3
Control-convergence parametersCk at different values of α.
Variablex1tx1tx2tx2tαC1C2C1C21−0.1198434251×10–3-0.02645325610−1.3256727843×10–5−2.527402984×10–30.9−0.1208162856×10–3-0.02826592550−1.343994006×10–5−2.536649902×10–30.8−0.1487623674×10–3-0.02983123651−1.3619012564×10–5−2.550122356×10–30.7−0.1598723560×10–3-0.03154109428−1.3801234527×10–5−2.573641295×10–3Table 4
Comparison and absolute error between OHAM-LS and RK4,α=1.
tkOHAM-LSx1tOHAM-LSx2tRK4x1tRK4x2tErrorx1tErrorx2t0.00.0000000.0000000.0000000.0000000.0000000.0000001.05.847225-0.4735475.847465-0.4775470.00024-0.0042.06.741648-0.5422756.741904-0.5425530.000256-0.0002784.08.455454-0.5841108.457551-0.5841330.002097-2.3E-056.08.731031-0.6152058.732010-0.6153050.000979-1E-047.08.977101-0.6403248.977187-0.6404668.6E-05-0.0001428.09.090736-0.6612019.090838-0.6616730.000102-0.00047210.09.303112-0.6800329.303232-0.6800660.00012-3.4E-05Figure 2
(a) Different values ofα (OHAM-LS; α=1 dot, α=0.9 dash, α=0.8 dash dot, and α=0.7 long dash) and Rk4 (α=1, solid) at x1. (b) Different values of α (OHAM-LS; α=1 dot, α=0.9 dash, α=0.8 dash dot, and α=0.7 long dash) and RK4 (α=1, solid) at x2.
(a)(b)This is a practical problem, and the exact solution is not known, but the expected optimal solution isx1∗=0, x2∗=1, x3∗=2, and x4∗=−1. From the above procedure, the second-order approximate solution obtained by OHAM-LS at α=0.9, for p=1, is
(106)x12t=16100C1−8050C12+8050C2t1/10+1.54569500×108t1/5C12,(107)x22t=−25900C1+12950C12−12950C2t1/10−3.75020500×108t1/5C12,(108)x32t=16420C1−8210C12+8210C2t1/10+1.62525800×108t1/5C12,(109)x42t=−46140C1+23070C12−23070C2t1/10−1.24619300×108t1/5C12.For the calculations ofC1 and C2 in x12t, x22t, x32t, and x42t given in Equations (4.79)-(4.82), we apply the procedure mentioned in Equations (19)-.(21); we obtain for x12t,
(110)c1=0,c2=0,and for x22t,
(111)c1=−4.494712729×10−13,c2=−0.1618198317×10−4,and for x32t,
(112)c1=−1.096696787×10−14,c2=0.599293243×10−3,and for x42t,
(113)c1=−5.935109529×10−12,c2=0.8773054262×10−4.Substituting these optimal values into Equations (106)-(109), we have
(114)x~1t=0,x~2t=0.095566833×t1/10−7.576330095×10−17t1/5,x~3t=3.119322656×t1/10+1.954769053×t1/5,x~4t=−2.023943344×t1/10−4.389780283×10−15t1/5.Tables5 and 6 show the Ck at different values of α for Example 3. Tables 7 and 8 show the comparisons and the absolute error between OHAM-LS and RK4 at α=1. Also, Figure 3 shows the comparisons of OHAM-LS at α=1,0.9,0.8, and 0.7 with RK4 at α=1, which verifies the performance of the present method as an excellent tool for NLPCOPs. For α=1, it can be seen that the approximate analytical solution agrees with the ideal solution. Thus, as α approaches 1, the classical solution for the system is recovered.Table 5
Control-convergence parametersCk at different values of α.
Variablex2tx2tx3tαC1C2C11−4.470470112×10–13−0.1542801253×10–3−1.0784243190×10–40.9−4.494712729×10–13−0.1618198317×10–3−1.096696787×10–40.8−4.5167327196×10–13−0.1832920121×10–3−1.1079094521×10–40.7−4.5371220162×10–13−0.2087212810×10–3−1.1261409123×10–4Table 6
Control-convergence parametersCk at different values of α.
Variablex3tx4tx4tαC2C1C21−0.5711237191×10–3−5.921274832×10–120.8973526178×10–40.90.599293243×10–3−5.935109529×10–120.8773054262×10–40.8−0.6190253261×10–3−5.953158452×10–120.8696526178×10–40.7−0.62131534211×10–3−5.975232801×10–120.8572034710×10–4Table 7
Comparisons and absolute error between OHAM-LS and RK4,α=1.
tkOHAM-LSx2tOHAM-LSx3ttOHAM-LSx4tRK4x2t0.0000.0000000.0000000.0000000.0000000.0010.8341651.835129-0.8055090.8342600.0020.8941011.886099-0.8592110.8941370.0030.9311121.928753-0.9637650.9311370.0040.9582111.965657-0.9786250.9583130.0050.9796231.998334-0.9918990.979937Table 8
Comparison and absolute error between OHAM-LS and RK4,α=1.
tkRK4x3tRK4x4tErrorx2tErrorx3tErrorx4t0.0000.0000000.0000000.0000000.0000000.0000000.0011.835278-0.8057469.5E-050.000149-0.0002370.0021.886124-0.8593103.6E-052.5E-05-9.9E-050.0031.928838-0.9638582.5E-058.5E-05-9.3E-050.0041.965795-0.9788310.0013720.000138-0.0002060.0051.998449-0.9919900.0003140.000115-9.1E-05Figure 3
(a) Different values ofα (OHAM-LS; α=1 dot, α=0.9 dash, α=0.8 dash dot, and α=0.7 long dash) and Rk4 (α=1, solid) at x2. (b) Different values of α (OHAM-LS; α=1 dot, α=0.9 dash, α=0.8 dash dot, and α=0.7 long dash), and RK4 (α=1, solid) at x3. (c) Different values of α (OHAM-LS; α=1 dot, α=0.9 dash, α=0.8 dash dot, and α=0.7 long dash), and RK4 (α=1, solid) at x4.
(a)(b)(c)


## 5. Conclusions

In this paper, we implemented OHAM-LS for solving nonlinear FOGBDS from the optimization problem. The fractional derivative is considered in a new conformable fractional derivative sense. The optimization minimization approach of the least square method helps to obtain optimal values of theCsj for accurate approximate analytical solutions. The comparisons between the fourth-order Runge-Kutta (α=1) and OHAM-LS show that our present method performs rapid convergence to the expected optimal solutions of the optimization problem. The results obtained are in close agreement with the exact solution, and those from the RK4 and OHAM-LS are reliable, dependable, and efficient for finding an approximate analytical solution for nonlinear FOGBDS optimization problem.



---

*Source: 8049397-2020-07-26.xml*