*Article* **Approximate Solution of Nonlinear Time-Fractional Klein-Gordon Equations Using Yang Transform**

**Jinxing Liu 1, Muhammad Nadeem 1,\*, Mustafa Habib <sup>2</sup> and Ali Akgül <sup>3</sup>**


**Abstract:** The algebras of the symmetry operators for the Klein–Gordon equation are important for a charged test particle, moving in an external electromagnetic field in a space time manifold on the isotropic hydrosulphate. In this paper, we develop an analytical and numerical approach for providing the solution to a class of linear and nonlinear fractional Klein–Gordon equations arising in classical relativistic and quantum mechanics. We study the Yang homotopy perturbation transform method (YHPTM),which is associated with the Yang transform (YT) and the homotopy perturbation method (HPM), where the fractional derivative is taken in a Caputo–Fabrizio (CF) sense. This technique provides the solution very accurately and efficiently in the form of a series with easily computable coefficients. The behavior of the approximate series solution for different fractional-order ℘ values has been shown graphically. Our numerical investigations indicate that YHPTM is a simple and powerful mathematical tool to deal with the complexity of such problems.

**Keywords:** fractional Klein–Gordon equation; Yang transform; homotopy perturbation method; series solution

#### **1. Introduction**

Recently, fractional calculus has grown in popularity due to its significant prospective applications in physics and engineering such as biology, mathematics, chemistry, fluid mechanics, physics, and nonlinear optics [1,2]. Fractional partial differential Equations (FPDEs) are a contemporary tool in calculus that can be used to simulate a wide range of classifications in applied sciences and engineering [3–5].

The Klein–Gordon (KG) equation performs a significant role in mathematical physics and many other scientific studies such as quantum field theory, nonlinear optics, and solidstate physic [6–10]. On the other hand, the fractional-order KG equation is derived from the classical KG equation by substituting the time order derivative with the fractional derivative of order ℘. The fractional-order KG equation can be illustrated as below

$$D\_{\theta}^{\phi}\theta(\epsilon,q) - D\_{\epsilon}^{2}\theta(\epsilon,q) + a\_{1}\theta(\epsilon,q) + a\_{2}G(\theta(\epsilon,q)) = f(\epsilon,q),\tag{1}$$

with initial conditions

$$
\theta(\epsilon, 0) = f\_1(\epsilon), \ \theta\_q(\epsilon, 0) = f\_2(\epsilon), \tag{2}
$$

where *D*<sup>℘</sup> *<sup>q</sup>* represents the Caputo fractional time derivative, *a*<sup>1</sup> and *a*<sup>2</sup> are real constants, *f*(, *q*), *f*1() and *f*2() are known as analytical functions, whereas *G*(*ϑ*(, *q*)) is a nonlinear, and *ϑ* is an unknown function of and *q*.

Various authors [11–15] have investigated different analytical and numerical strategies to examine the solution to the KG equation but with some restrictions and lacks. Tamsir and Srivastava [16] used fractional reduced differential transform to obtain the analytical solution of linear and nonlinear KG equation with time-fractional order. Bansu and

**Citation:** Liu, J.; Nadeem, M.; Habib, M.; Akgül, A. Approximate Solution of Nonlinear Time-Fractional Klein-Gordon Equations Using Yang Transform. *Symmetry* **2022**, *14*, 907. https://doi.org/10.3390/sym 14050907

Academic Editor: Eduard Marusic-Paloka

Received: 14 April 2022 Accepted: 26 April 2022 Published: 29 April 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 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/).

Kumar [17] used a radial basis approach, and Kurulay [18] applied the homotopy analysis method to evaluate the numerical solution of the space-time fractional KG equation. Later, Khader and Adel [19] applied a hybridization scheme to achieve the solution of the fractional KG equation. Zhmud and Dimitrov [20] developed the fractional integration method, which is based on extrapolation using a series of integrating and differentiating links with a time constant that changes symmetrically from one step to another. In order to obtain the solution of FPDEs, several valuable strategies have been considered, such as the generalized differential transform method [21], the adomian decomposition method [22], the homotopy analysis method [23], the variational iteration method [24], the homotopy perturbation method [25], the Elzaki transform decomposition method [26], the fractional wavelet method [27,28] and the residual power series method [29,30].

In this paper, we present the Yang homotopy perturbation transform method (YHPTM), which is a composition of YT and HPM. The primary objective of this approach is to investigate the approximate solution of fractional KG equations and minimize the computational work that overcomes nonlinear problems easily. Next, this scheme can promptly deal with the nonlinear KG equation. Finally, this method can reduce the range of the computations and generate an approximate solution with elegantly computed expressions, which is its most impressive advantage. The design of this paper is framed as follows. In Section 2, we start with some primary definitions of Caputo–Fabrizio. In Section 3, we formulate the idea of the Yang homotopy perturbation transform method. In Section 4, we perform this scheme on some illustrative examples to show its capability and efficiency. Concluding remarks are given in Section 5.

#### **2. Preliminaries and Concepts**

**Definition 1.** *The CF derivative is described as [31]*

$${}^{CF}D\_{\eta}^{\wp}\theta(\varepsilon,\eta) = \frac{S(\wp)}{1-\wp} \int\_{0}^{\eta} [Q'(\varrho)K(\varrho,\varrho)]d\varrho, \quad n-1 < \wp \le n \tag{3}$$

*S*(℘) *is the normalization function with S*(0) = *S*(1) = 1*, and then, Equation (3) becomes as*

$${}^{CF}D\_{q}^{\wp}\theta(\varepsilon,\eta) = \frac{\mathcal{S}(\wp)}{1-\wp} \int\_{0}^{q} [\mathcal{Q}(\eta) - \mathcal{Q}(\varrho)] K(\eta,\varrho) d\varrho, \quad n-1 < \wp \le n \tag{4}$$

**Definition 2.** *The fractional CF integral is stated as [32]*

$${}^{CF}I\_{\emptyset}^{\wp}\theta(\varepsilon,q) = \frac{1-\wp}{S(\wp)}\mathcal{Q}(q) + \frac{\wp}{S(\wp)}\int\_{0}^{q} \mathcal{Q}(q)d\varrho, \quad q \ge 0,\ \wp\epsilon(0,1]. \tag{5}$$

**Definition 3.** *For S*(℘) = 1*, the Laplace transform of the CF derivative is [33]*

$$L\left[\prescript{\text{CF}}{}{D\_q^\wp Q}[(q)]\right] = \frac{vL[Q(q) - Q(0)]}{v + \wp(1 - v)}.\tag{6}$$

**Definition 4.** *The* Y*T of Q*(*q*) *is framed as [34]*

$$\mathbb{V}[Q(q)] = \chi(v) = \int\_0^\infty Q(q)e^{-\frac{q}{v}}dq. \ q > 0\tag{7}$$

*Remarks*

The YT of some helpful expressions are as follows:

$$\begin{aligned} \mathbb{Y}[\mathbf{1}] &= v; \\ \mathbb{Y}[q] &= v^2; \\ \mathbb{Y}[q^i] = \Gamma(i+1)v^{i+1}. \end{aligned}$$

**Lemma 1.** *Let the Laplace transform of Q*(*q*) *be F*(*v*)*, and then χ*(*v*) = *F*(1/*v*) *[35].*

**Proof.** From Equation (7), we can obtain the Yang transform by putting *q*/*v* = *ζ* as

$$L[Q(q)] = \int\_0^\infty Q(v\zeta)e^{\zeta}d\zeta,\ \zeta > 0\tag{8}$$

since *L*[*Q*(*q*)] = *F*(*v*), which implies that

$$F(v) = L[Q(q)] = \int\_0^\infty Q(q)e^{-vq}d\eta. \tag{9}$$

Putting *q* = *ζ*/*v* in Equation (9), we obtain

$$F(v) = \frac{1}{v} \int\_0^\infty Q\left(\frac{\zeta}{v}\right) e^{\zeta} d\zeta. \tag{10}$$

Thus, from Equation (8), we obtain:

$$F(\upsilon) = \chi\left(\frac{1}{\upsilon}\right). \tag{11}$$

Furthermore, from Equations (7) and (9), we obtain

$$F\left(\frac{1}{v}\right) = \chi(v). \tag{12}$$

The links between Equations (11) and (12) represent the duality connection among the Laplace and Yang transforms.

**Lemma 2.** *Let Q*(*q*) *be a function, then* Y*T of CF derivatives of Q*(*q*) *is [35]*

$$\mathbb{V}\left[Q(q)\right] = \frac{\mathbb{V}\left[Q(q) - vQ(0)\right]}{v + \wp(v - 1)}.\tag{13}$$

**Proof.** The fractional Laplace transform of CF is defined as in Equation (13)

$$L\left[Q(q)\right] = \frac{L\left[vQ(q) - Q(0)\right]}{v + \wp(1 - v)}.\tag{14}$$

However, we have a correlation among the YT and Laplace properties, namely *χ*(*v*) = *F*(1/*v*), so put 1/*v* for *v* in Equation (14), and we obtain

$$\begin{aligned} \mathbb{V}\left[Q(q)\right] &= \frac{\mathbb{V}\left[\frac{1}{v}Q(q) - Q(0)\right]}{\frac{1}{v} + \wp(1 - \frac{1}{v})},\\ \mathbb{V}\left[Q(q)\right] &= \frac{\mathbb{V}\left[Q(q) - vQ(0)\right]}{1 + \wp(v - 1)}.\end{aligned} \tag{15}$$

Thus, the proof is satisfied.

#### **3. Idea of Yang Homotopy Perturbation Transform Method (**Y**HPTM)**

In this part, we will demonstrate the concept of YHPTM. Let us assume a nonlinear fractional-order PDE, such as

$${}^{CF}D\_q^{\phi} \theta(\varepsilon, q) + R\theta(\varepsilon, q) + N\theta(\varepsilon, q) = \mathfrak{g}(\varepsilon, q), \tag{16}$$

$$
\vartheta(\epsilon, 0) = h(\epsilon),
\tag{17}
$$

where *g*(, *q*) is called the source function. Applying the YT to Equation (16),

$$\frac{1}{\mathcal{U}^{\flat}}\mathbb{Y}\left[\theta(\varepsilon,q)-\upsilon\theta(\varepsilon,0)\right]=-\mathbb{Y}[\mathcal{R}(\theta(\varepsilon,q))+N(\theta(\varepsilon,q))+\mathbb{Y}[\mathcal{g}(\varepsilon,q)]],$$

$$\mathbb{Y}[\theta(\varepsilon,q)]=\upsilon\hbar(\varepsilon)-\upsilon^{\flat}\left[\mathbb{Y}[\mathcal{R}(\theta(\varepsilon,q))+N(\theta(\varepsilon,q))\right]\right]+\mathbb{Y}[\mathcal{g}(\varepsilon,q)].$$

By using inverse YT,

$$\theta(\varepsilon, q) = \theta(\varepsilon, 0) - \mathbb{Y}^{-1} \left[ \upsilon^{\wp} \left[ \mathbb{Y} [\mathbb{R}(\theta(\varepsilon, q)) + \mathbb{N}(\theta(\varepsilon, q))] \right] + \mathbb{Y} [\mathbb{g}(\varepsilon, q)] \right]. \tag{18}$$

However, HPM is stated as

$$\vartheta(\epsilon, q) = \sum\_{i=0}^{\infty} p^i \vartheta\_i(\epsilon, q)\_i \tag{19}$$

and

$$N\theta(\varepsilon,q) = \sum\_{i=0}^{\infty} p^i H\_i \theta(\varepsilon,q). \tag{20}$$

The following strategy can be operated to acquire the He's polynomials,

$$H\_i(\vartheta\_0 + \vartheta\_1 + \cdots + \vartheta\_i) = \frac{1}{n!} \frac{\partial^i}{\partial p^i} \left( N(\sum\_{i=0}^{\infty} p^i \vartheta\_i) \right)\_{p=0} . \quad n = 0, 1, 2, \cdots$$

With the help of Equations (19) and (20), we can obtain Equation (18), such as

$$\sum\_{i=0}^{\infty} p^i \theta\_i(\varepsilon, q) = \vartheta(\varepsilon, 0) - p^\chi \mathbb{I}^{-1} \left[ v^{\psi} \mathbb{I} \left\{ \mathcal{R} \left( \sum\_{i=0}^{\infty} p^i \theta\_i(\varepsilon, q) \right) + \sum\_{i=0}^{\infty} p^i H\_{\mathbb{II}} \theta\_i(\varepsilon, q) \right\} \right].$$

We can obtain the following terms by evaluating the *p* components:

$$\begin{aligned} p^0 &= \theta\_0(\epsilon, q) = \theta(\epsilon, 0), \\ p^1 &= \theta\_1(\epsilon, q) = -\mathbb{Y}^{-1} \left[ v^{\psi} \mathbb{Y} \left\{ R\theta\_0(\epsilon, q) + H\_0(\vartheta) \right\} \right], \\ p^2 &= \theta\_2(\epsilon, q) = -\mathbb{Y}^{-1} \left[ v^{\psi} \mathbb{Y} \left\{ R\theta\_1(\epsilon, q) + H\_1(\vartheta) \right\} \right], \\ p^3 &= \theta\_3(\epsilon, q) = -\mathbb{Y}^{-1} \left[ v^{\psi} \mathbb{Y} \left\{ R\theta\_2(\epsilon, q) + H\_2(\vartheta) \right\} \right], \\ &\vdots \\ p^i &= \theta\_i(\epsilon, q) = -\mathbb{Y}^{-1} \left[ v^{\psi} \mathbb{Y} \left\{ R\theta\_i(\epsilon, q) + H\_i(\vartheta) \right\} \right], \end{aligned} \tag{21}$$

Thus, we can summarize the set of Equations (21) in the series form, such as

$$\begin{aligned} \theta(\varepsilon, q) &= \theta\_0(\varepsilon, q) + \theta\_1(\varepsilon, q) + \theta\_2(\varepsilon, q) + \dotsb \\ \theta(\varepsilon, q) &= \lim\_{N \to \infty} \sum\_{n=0}^{N} \theta\_n(\varepsilon, q) \end{aligned} \tag{22}$$

#### **4. Numerical Applications**

#### *4.1. Example 1*

Consider a linear time-fractional KG problem

$$D\_q^\wp \theta(\epsilon, q) - D\_\epsilon^2 \theta(\epsilon, q) - \theta(\epsilon, q) = 0,\tag{23}$$

with the initial condition

$$
\vartheta(\varepsilon, 0) = 1 + \sin(\varepsilon). \tag{24}
$$

Taking YT of Equation (23), we obtain

$$\mathbb{Y}\left[\frac{\partial^{\wp}\theta}{\partial q^{\wp}}\right] = \mathbb{Y}\left[\frac{\partial^{2}\theta}{\partial \varepsilon^{2}} + \theta\right].$$

Executing the differential property of YT, we obtain

$$\begin{split} \frac{1}{v^{\flat}} \mathbb{Y} \Big[ \theta(\epsilon, q) - v\theta(\epsilon, 0) \Big] &= \mathbb{Y} \Big[ \frac{\partial^2 \theta}{\partial q^2} + \theta \Big], \\\mathbb{Y} \Big[ \theta(\epsilon, q) \Big] &= v\theta(\epsilon, 0) + v^{\flat} \mathbb{Y} \Big[ \frac{\partial^2 \theta}{\partial q^2} + \theta \Big]. \end{split}$$

The inverse YT indicates

$$
\vartheta(\epsilon, q) = \vartheta(\epsilon, 0) + \mathbb{Y}^{-1} \left[ v^{\psi} \left\{ \mathbb{Y} \left( \frac{\partial^2 \vartheta}{\partial q^2} + \vartheta \right) \right\} \right].
$$

Employing HPM such as

$$\theta(\varepsilon, q) = \theta\_0 + p\theta\_1 + p^2\theta\_2 + \cdots$$

$$\sum\_{i=0}^{\infty} p^i \theta\_i(\varepsilon, q) = 1 + \sin(q) + p \left( \mathbb{Y}^{-1} \left[ v^{\psi} \left\{ \mathbb{Y} \left( \sum\_{i=0}^{\infty} p^i \frac{\partial^2 \theta\_i}{\partial q^2} + \sum\_{i=0}^{\infty} p^i \theta\_i \right) \right\} \right] \right)$$

on comparing the identical of *p*, we obtain

$$\begin{split} p^{0} &= \theta\_{0}(\boldsymbol{\epsilon}, q) = \theta(\boldsymbol{\epsilon}, 0), \\ p^{1} &= \theta\_{1}(\boldsymbol{\epsilon}, q) = \mathbb{Y}^{-1} \left[ \boldsymbol{v}^{\psi} \left\{ \mathbb{Y} \left( \frac{\partial^{2} \theta\_{0}}{\partial q^{2}} + \vartheta\_{0} \right) \right\} \right], \\ p^{2} &= \theta\_{2}(\boldsymbol{\epsilon}, q) = \mathbb{Y}^{-1} \left[ \boldsymbol{v}^{\psi} \left\{ \mathbb{Y} \left( \frac{\partial^{2} \theta\_{1}}{\partial q^{2}} + \vartheta\_{1} \right) \right\} \right], \\ p^{3} &= \theta\_{3}(\boldsymbol{\epsilon}, q) = \mathbb{Y}^{-1} \left[ \boldsymbol{v}^{\psi} \left\{ \mathbb{Y} \left( \frac{\partial^{2} \theta\_{2}}{\partial q^{2}} + \vartheta\_{2} \right) \right\} \right], \\ &\vdots \end{split}$$

With help of Equation (24), we gain the iterations successively *ϑi*(), *i* = 1, 2, 3, ··· , as follows:

$$\begin{aligned} \vartheta\_0(\varepsilon, q) &= 1 + \sin(\varepsilon), \\ \vartheta\_1(\varepsilon, q) &= \frac{1}{\Gamma(1 + \wp)} q^{\wp} \end{aligned}$$

$$\begin{aligned} \theta\_2(\varepsilon, q) &= \frac{1}{\Gamma(1 + 2\wp)} q^{2\wp}, \\ \theta\_3(\varepsilon, q) &= \frac{1}{\Gamma(1 + 3\wp)} q^{3\wp}, \\ &\vdots \\ \theta\_i(\varepsilon, q) &= \frac{1}{\Gamma(1 + i\wp)} q^{i\wp}. \end{aligned}$$

Thus, the approximate solution can be obtained by:

$$\begin{split} \theta(\varepsilon, q) &= 1 + \sin(\varepsilon) + \frac{1}{\Gamma(1+\wp)} q^{\wp} + \frac{1}{\Gamma(1+2\wp)} q^{2\wp} + \frac{1}{\Gamma(1+3\wp)} q^{3\wp} + \cdots \\ &= 1 + \sin(\varepsilon) + \sum\_{i=0}^{\infty} p^i \frac{q^{i\wp}}{\Gamma(1+i\wp)} \end{split} \tag{25}$$

which implies the exact solution of Equation (23), In particular, at ℘ = 1, we obtain

$$
\theta(\epsilon, q) = 1 + \sin(\epsilon), \tag{26}
$$

which is in full agreement.

Figure 1a–d indicate the physical behavior of the obtained solution at ∈ [0, 4] and *q* ∈ [0, 0.8]. From these figures, it can be observed that the solution graphs of the problem show the friendly touch with each other. Figure 1a–d demonstrate that the solution achieved by YHPTM approaches the precise solution very rapidly with more iterations. In Figure 2, we have plotted the graph of *ϑ*(, *q*) at different fractional order of ℘ = 0.25, 0.50, 0.75, 1 and ∈ [0, 2*π*] with different values of *q*.

**Figure 1.** The surface solution of *ϑ*(, *q*) with respect to and *q* for distinct values of ℘: (**a**) surface solution of *ϑ*(, *q*) when ℘ = 0.25; (**b**) surface solution of *ϑ*(, *q*) when ℘ = 0.50; (**c**) surface solution of *ϑ*(, *q*) when ℘ = 0.75; (**d**) surface solution of *ϑ*(, *q*) when ℘ = 1.

**Figure 2.** Plot of *ϑ*(, *q*) for different values of ℘.

#### *4.2. Example 2*

Assume a nonlinear time-fractional KG problem

$$D\_q^\wp \theta(\varepsilon, q) - D\_\varepsilon^2 \theta(\varepsilon, q) + \theta^2(\varepsilon, q) = 0,\tag{27}$$

with the initial condition

$$
\vartheta(\varepsilon, 0) = 1 + \sin(\varepsilon). \tag{28}
$$

Taking the Yang transform of Equation (27), we obtain

$$\mathbb{V}\left[\frac{\partial^{\wp}\theta}{\partial q^{\wp}}\right] = \mathbb{V}\left[\frac{\partial^{2}\theta}{\partial \varepsilon^{2}} - \theta^{2}\right].$$

Executing the differential property of YT, we obtain

$$\frac{1}{v^{\wp}}\mathbb{Y}\left[\theta(\varepsilon,q)-v\theta(\varepsilon,0)\right]=\mathbb{Y}\left[\frac{\partial^2\theta}{\partial q^2}-\vartheta^2\right],$$

$$\mathbb{Y}\left[\theta(\varepsilon,q)\right]=v\theta(\varepsilon,0)+v^{\wp}\mathbb{Y}\left[\frac{\partial^2\theta}{\partial q^2}-\vartheta^2\right].$$

The inverse YT indicates

$$
\vartheta(\epsilon, q) = \vartheta(\epsilon, 0) + \mathbb{Y}^{-1} \left[ v^{\phi} \left\{ \mathbb{Y} \left( \frac{\partial^2 \theta}{\partial q^2} - \theta^2 \right) \right\} \right].
$$

Employing HPM such as

$$\sum\_{i=0}^{\infty} p^i \vartheta\_i(\varepsilon, q) = 1 + \sin(q) + p \left( \mathbb{Y}^{-1} \left[ v^{\diamond} \left\{ \mathbb{Y} \left( \sum\_{i=0}^{\infty} p^i \frac{\partial^2 \vartheta\_i}{\partial q^2} - \sum\_{i=0}^{\infty} p^i \vartheta\_i^2 \right) \right\} \right] \right),$$

on comparing the identical of *p*, we obtain

$$\begin{split} p^{0} &= \theta\_{0}(\boldsymbol{\epsilon}, q) = \theta(\boldsymbol{\epsilon}, 0), \\ p^{1} &= \theta\_{1}(\boldsymbol{\epsilon}, q) = \mathbb{Y}^{-1} \left[ \boldsymbol{\upsilon}^{\psi} \left\{ \mathbb{Y} \left( \frac{\partial^{2} \theta\_{0}}{\partial q^{2}} - \theta\_{0}^{2} \right) \right\} \right], \\ p^{2} &= \theta\_{2}(\boldsymbol{\epsilon}, q) = \mathbb{Y}^{-1} \left[ \boldsymbol{\upsilon}^{\psi} \left\{ \mathbb{Y} \left( \frac{\partial^{2} \theta\_{1}}{\partial q^{2}} - 2\theta\_{0} \theta\_{1} \right) \right\} \right], \\ p^{3} &= \theta\_{3}(\boldsymbol{\epsilon}, q) = \mathbb{Y}^{-1} \left[ \boldsymbol{\upsilon}^{\psi} \left\{ \frac{\partial^{2} \theta\_{2}}{\partial q^{2}} - \theta\_{1}^{2} - 2\theta\_{0} \theta\_{2} \right\} \right], \\ &\vdots \end{split}$$

With help of Equation (28), we gain the iterations successively *ϑi*(), *i* = 1, 2, 3, ··· , as follows:

$$\begin{aligned} \theta\_0(\varepsilon, q) &= 1 + \sin(\varepsilon), \\ \theta\_1(\varepsilon, q) &= \frac{-q^{\wp}}{\Gamma(1 + \wp)} \Big( 1 + 3\sin(\varepsilon) + \sin^2(\varepsilon) \Big), \\ \theta\_2(\varepsilon, q) &= \frac{q^{2\wp}}{\Gamma(1 + 2\wp)} \Big( 11\sin(\varepsilon) + 12\sin^2(\varepsilon) + 2\sin^3(\varepsilon) \Big), \\ \theta\_3(\varepsilon, q) &= \frac{q^{3\wp}}{\Gamma(1 + 3\wp)} \Big( 18 - 57\sin(\varepsilon) - 160\sin^2(\varepsilon) - 82\sin^3(\varepsilon) - 10\sin(4\varepsilon) \Big), \\ &\vdots \end{aligned}$$

Thus, the approximate solution can be obtained by:

$$\begin{split} \theta(\varepsilon,q) &= 1 + \sin(\varepsilon) - \frac{q^{\wp}}{\Gamma(1+\wp)} \left( 1 + 3\sin(\varepsilon) + \sin^2(\varepsilon) \right) + \frac{q^{2\wp}}{\Gamma(1+2\wp)} \left( 11\sin(\varepsilon) + 12\sin^2(\varepsilon) + 2\sin^3(\varepsilon) \right) \\ &+ \frac{q^{2\wp}}{\Gamma(1+3\wp)} \left( 18 - 57\sin(\varepsilon) - 160\sin^2(\varepsilon) - 82\sin^3(\varepsilon) - 10\sin(4\varepsilon) \right) + \cdots \end{split}$$

Figure 3a–d indicate the physical behavior of the obtained solution at ∈ [0, 1] and *q* ∈ [0, 1]. From these figures, it can be observed that with the increase in the value of ℘, the approximate solution become close to the exact solution at ℘ = 1. In Figure 4, we have plotted the graph of *ϑ*(, *q*) with different fractional orders of ℘ = 0.25, 0.50, 0.75, 1 at ∈ [0, 2*π*] with different values of *q*. It is obvious that this approximation can only be employed numerically, even though a closed form solution is not accessible. It can be seen that our approximate solution using YHPTM in Table 1 is more significant than that obtained in [36,37].

**Figure 3.** The surface solution of *ϑ*(, *q*) with respect to and *q* for distinct values of ℘: (**a**) surface solution of *ϑ*(, *q*) when ℘ = 0.25 ; (**b**) surface solution of *ϑ*(, *q*) when ℘ = 0.50; (**c**) surface solution of *ϑ*(, *q*) when ℘ = 0.75; (**d**) surface solution of *ϑ*(, *q*) when ℘ = 1.

**Figure 4.** Plot of *ϑ*(, *q*) for different values of ℘

**Table 1.** Comparison between the value *ϑ*(, *q*) for the solution of the KG equation.


#### *4.3. Example 3*

Consider another nonlinear time-fractional KG problem

$$D\_q^\wp \theta(\varepsilon, q) - D\_\epsilon^2 \theta(\varepsilon, q) + \theta(\varepsilon, q) - \theta^\Im(\varepsilon, q) = 0,\tag{29}$$

with the initial condition

$$
\vartheta(\varepsilon, 0) = -\operatorname{sech}(\varepsilon). \tag{30}
$$

According to the idea of YHPTM, we can obtain the following relation

$$\sum\_{i=0}^{\infty} p^i \theta\_i(\epsilon, q) = \text{sech}(\epsilon) + p \left( \mathbb{Y}^{-1} \left[ v^{\psi} \left\{ \mathbb{Y} \left( \sum\_{i=0}^{\infty} p^i \frac{\partial^2 \theta\_i}{\partial q^2} - \sum\_{i=0}^{\infty} p^i \theta\_i + \sum\_{i=0}^{\infty} p^i \theta\_i^2 \right) \right\} \right] \right)$$

when the coefficients of like powers of *p* are compared, we obtain

$$\begin{split} p^{0} &= \theta\_{0}(\boldsymbol{\epsilon}, \boldsymbol{q}) = \theta(\boldsymbol{\epsilon}, \boldsymbol{0}), \\ p^{1} &= \theta\_{1}(\boldsymbol{\epsilon}, \boldsymbol{q}) = \mathbb{Y}^{-1} \left[ \mathbb{V}^{\psi} \left\{ \mathbb{V} \left( \frac{\partial^{2} \theta\_{0}}{\partial q^{2}} - \theta\_{0} + \theta\_{0}^{3} \right) \right\} \right], \\ p^{2} &= \theta\_{2}(\boldsymbol{\epsilon}, \boldsymbol{q}) = \mathbb{Y}^{-1} \left[ \mathbb{V}^{\psi} \left\{ \mathbb{V} \left( \frac{\partial^{2} \theta\_{1}}{\partial q^{2}} - \theta\_{1} + 3\theta\_{0}^{2} \theta\_{1} \right) \right\} \right], \\ p^{3} &= \theta\_{3}(\boldsymbol{\epsilon}, \boldsymbol{q}) = \mathbb{Y}^{-1} \left[ \mathbb{V}^{\psi} \left\{ \mathbb{V} \left( \frac{\partial^{2} \theta\_{2}}{\partial q^{2}} - \theta\_{2} + 3\theta\_{0} \theta\_{1}^{2} + 3\theta\_{0}^{2} \theta\_{2} \right) \right\} \right], \\ &\vdots & \vdots \end{split}$$

with help of Equation (30), we gain the iterations successively *ϑi*(), *i* = 1, 2, 3, ··· , as follows:

$$\begin{aligned} \theta\_0(\varepsilon, q) &= -\operatorname{sech}(\varepsilon), \\ \theta\_1(\varepsilon, q) &= -\frac{q^{\wp}}{\Gamma(1+\wp)} \Big( 2\operatorname{sech}(\varepsilon) - 3\operatorname{sech}^3(\varepsilon) \Big), \\ \theta\_2(\varepsilon, q) &= -\frac{q^{2\wp}}{\Gamma(1+2\wp)} \Big( 3\operatorname{sech}(\varepsilon) - 34\operatorname{sech}^3(\varepsilon) - 18\operatorname{sech}^5(\varepsilon) \Big), \\ \theta\_3(\varepsilon, q) &= -\frac{q^{3\wp}}{\Gamma(1+3\wp)} \Big( 64\operatorname{sech}^3(\varepsilon) - 288\operatorname{sech}^5(\varepsilon) + 240\operatorname{sech}^7(\varepsilon) \Big), \\ &\vdots \end{aligned}$$

Thus, the approximate solution can be obtained by:

$$\begin{split} \theta(\varepsilon,q) &= -\operatorname{sech}(\varepsilon) - \frac{q^{\wp}}{\Gamma(1+\wp)} \Big( 2\operatorname{sech}(\varepsilon) - 3\operatorname{sech}^{\mathfrak{z}}(\varepsilon) \Big) - \frac{q^{2\wp}}{\Gamma(1+2\wp)} \Big( 3\operatorname{sech}(\varepsilon) - 34\operatorname{sech}^{\mathfrak{z}}(\varepsilon) - 18\operatorname{sech}^{\mathfrak{z}}(\varepsilon) \Big) \\ &- \frac{q^{3\wp}}{\Gamma(1+3\wp)} \Big( 64\operatorname{sech}^{\mathfrak{z}}(\varepsilon) - 288\operatorname{sech}^{\mathfrak{z}}(\varepsilon) + 240\operatorname{sech}^{\mathfrak{z}}(\varepsilon) \Big) + \cdots \end{split}$$

Figure 5a–d indicates the physical behavior of the obtained solution at ∈ [−2, 2] and *q* ∈ [0, 0.1]. From these figures, it can be observed that with increase in the value of ℘, the approximate solution graph comes close to the exact exact solution at ℘ = 1. In Figure 6, we have plotted the graph of *ϑ*(, *q*) at different fractional orders of ℘ = 0.25, 0.50, 0.75, 1 and ∈ [0, 2*π*] with different values of *q*. We compared our graphical results obtained by YHPTM, which converges to the exact solution very rapidly with a small amount of *q* compared to [38] at ℘ = 1.

**Figure 5.** The surface solution of *ϑ*(, *q*) with respect to and *q* for distinct values of ℘. (**a**) surface solution of *ϑ*(, *q*) when ℘ = 0.25; (**b**) surface solution of *ϑ*(, *q*) when ℘ = 0.50; (**c**) surface solution of *ϑ*(, *q*) when ℘ = 0.75; (**d**) surface solution of *ϑ*(, *q*) when ℘ = 1.

**Figure 6.** Plot of *ϑ*(, *q*) for different values of ℘.

#### **5. Conclusions**

In this study, YHPTM has been utilized to achieve an approximate solution of nonlinear time-fractional KG equations. We demonstrate some illustrations to show the validity of the method, which reveals that the obtained findings are satisfactory. It is important to note that in order to improve the accuracy, a greater number of iterations with excessive order of *p* are required. The best advantage of YHPTM is that it generates the approximate solution in a quickly convergent power series form. As a result, this strategy can be adopted to elucidate a wide classification of nonlinear challenges in science and engineering with no need for linearization, discretization or perturbation. The proposed strategy has the privilege of being able to tackle linear and nonlinear time-fractional KG problems simultaneously. Mathematica package 11.0.1 has been operated for the graphical analysis as well as for the computations in this paper. We recommend the readers to consider this problem for the Atangana–Baleanu fractional derivative and many others in the place of the Caputo–Fabrizio operator.

**Author Contributions:** Conceptualization, formal analysis and funding acquisition, J.L.; investigation, methodology, software and writing—original draft, M.N.; validation and visualization, M.H.; writing—review and editing and supervision, A.A. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was supported by the Foundation of Yibin University, China (Grant No. 2019QD07).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The authors declare that they have no competing of interest.

#### **Abbreviations**

The following abbreviations are used in this manuscript:


#### **References**


MDPI St. Alban-Anlage 66 4052 Basel Switzerland www.mdpi.com

*Symmetry* Editorial Office E-mail: symmetry@mdpi.com www.mdpi.com/journal/symmetry

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.

Academic Open Access Publishing

mdpi.com ISBN 978-3-0365-9425-5