Entropy and fractional calculus are appealing concepts that are increasingly being used to investigate the dynamics of complicated systems. Fractional calculus (FC) has been increasingly used in numerous sectors of research in recent years. Fractional differential equalities (FDEs) efficiently depict the natural evolution associated with viscoelasticity, models of porous electrodes, thermal stresses, electromagnetism, energy transmission in viscous dissipation systems, relaxing oscillations, and thermoelasticity. Most of the mathematical models are obtained through real-world problems which can be modeled via differential equations of the integer order or of the fractional order. The differential equations may arise in diverse areas of technological sciences and biological sciences. In engineering sciences, they may be in the field of fluid dynamics, aerodynamics, the nuclear decay, climate changes, electronic circuits, etc. In biosciences, these models may be of the blood flow, population growth, and decay problems of some kind of species of organisms like bacteria or virus and may be the study of some rate of flow of some gas models or may be the concentration control of some liquids in some other liquids. Similarly, differential equations may also model some problems related to social sciences as in the fields of banking and finance. The customer may be satisfied by preparing a fractional model of interest or premium according to the required efforts of a person. Fractional calculus plays an essential role in these kinds of problems.

Fractional calculation is the generalization of the classical calculus which is an ancient branch of mathematics. The fractional calculus received much more attention of researchers during the last few decades. Fractional calculus has great achievements in the fields of physics, engineering, biology, medicine, hydrology, economics, and finance [1-5]. The models of differential equations may be linear or nonlinear; linear models can be solved easily by different methods and do not require too much difficulty for obtaining the exact solution. But most of the problems of the real world occur non-linearly and cannot be solved easily. They are very hard to solve by simple methods. Most of the non-linear problems do not have exact solutions. Therefore, researchers use different approaches to solve them.

Researchers use numerical methods to solve non-linear problems, but they have discretization issue, are costly, and time-consuming. The famous numerical methods are the following: the collocation method, finite difference technique, finite element procedure, and radial basis function technique [6-8]. Similarly, perturbation methods need small or large parameter assumptions which are very difficult [9, 10]. Non-perturbation methods are the Adomian decomposition methodology (ADM) and differential transformation methodolgy (DTM). These methods work on repetition and that is why these types of problems can be solved with the help of computer software easily. Some well-known iterative methods are the variational iterating methodology (VIM), new iterating methodology, modified variational iteration method (MVIM), etc. [11, 12]. The Zakharov–Kuznetsov (ZK) equation is an enticing modeling formula for studying vortices in geophysics flows. The ZK difficulties appear in many areas of material sciences, implemented arithmetic, and design. It occurs particularly in the realm of physical sciences. The ZK issues govern the behavior of weak non-linear particle acoustic plasma waves, such as cold nanoparticles and hot adiabatic electrons, in the presence of smooth magnetism. The non-direct higher order of the expanded KdV criteria for geometrical removal was used to generate solitary wave configurations. The accurate expository structures of various non-linear advancement equations in numerical materials engineering, namely, space time-fractional Zakharov–Kuznetsov and altered Zakharov–Kuznetsov formulas, were found using a fractional technique. Many approaches, including the new iterating Sumudu understanding of the complex, homotopy perturbation transform method, expanded direct algebra methodology, natural decomposition technique, and q-homotopy analysis transform methodology, have been used to examine it during the last few decades. In this research work, we will find the approximate solution to the fractional order of the Zakharova–Kuznetsova FZK equation (13). The general form of the FZK equation is D t β φ + a ( φ p ) x + b ( φ q ) x x x + c ( φ r ) y y x = 0 ,   0 ≤ t ,   0 < β ≤ 1 , (1) where φ = φ ( x , y , t ) ,   0 < β ≤ 1 signifies an order of the fractional derivative, and a , b and c are optional fixed factors. The integers p , q , r control the behavior of weak non-linear ion acoustical waves in hemoglobin-containing coolant ions and hotness isotherm electrons in the existence of a consistent magneto force. Numerous researchers have tried to solve the FZK equation by using different approaches such as VIM, OHAM, PIA method, and RPSM [14]. We have obtained the solution of the FZK equation by NTIM which is an extension of the natural iterating methodology NIM presented by Gejji and Jafari [15, 16] to obtain the estimated solution of linear and non-linear differential equalities. NTIM was recently applied by Nawaz et al. [17] for solving the fractional order differential equation. In the proposed methodology, NIM is combined with the natural transform for the solution of the FZK equation. We observed that the proposed method was easy to implement and provide an encouraging approximate solution for the linear and non-linear differential equalities of fractional- and integer-order derivatives.

Definition 2.1 [14]: The fractional integral operator I α of order α ≥ 0 in the Riemann–Liouville idea of a function is described as I α f ( ( χ ) ( ) = 1 Γ ( α ) ∫ 0 χ ( χ − s ) α − 1 f ( s ) d s ,   α , χ > 0 , (2) where I 0 f ( χ ) = f ( χ ) and Γ is the well-known function.

Definition 2.2 [14]: Riemann–Liouville fractional derivative can be defined if g ( r ) ∈ C [ a , b ] then I a α = 1 Γ ( α ) ∫ a r g ( χ ) ( r − χ ) 1 − α d χ . (3)

Some properties of the fractional derivative and integral are given as f ∈ C μ , μ ≥ 1 , α , β ≥ 0 and λ > − 1 then • I α I β = I α + β f ( χ ) , • I β I α = I α + β f ( χ ) , • I α χ λ = Γ ( λ + 1 ) Γ ( λ + 1 + α ) χ λ + α .

Definition 2.3: Natural transubstantiate is specified as [18] N + [ ϕ ( t ) ] = R ( s ,   u ) = 1 u ∫ 0 ∞ e − s t u ( ϕ ( t ) ) d t ,   s ,   u   >   0 . (4)

u and s are the transformation variables.

Definition 2.4 [18]: The inverse of the natural transubstantiate R ( s , u ) is defined as N − [ R ( s , u ) ] = ϕ ( t ) = 1 2 π i   ∫ c − i ∞ c + i ∞ e s t u ( R   ( s ,   u ) ) d s . (5)

s = a + b i is the complex plan for executing the integral along s = c , where c ∈ R .

Definition 2.5 [18]: If ϕ ( t ) is a function and ϕ n ( t ) is its nth derivative, then the natural transubstantiate of ϕ ( t ) is N + ( φ n ( t ) ) = R n ( s ,   u ) = s n u n R ( s ,   u ) − ∑ k = 0 n − 1 s n − ( k + 1 ) u n − k   ( φ n (   0   ) ) ,   n ≥ 1 . (6)

Theorem 2.6: If k ( t ) and h ( t ) are defined on a set A and have the natural transforms K ( s , u ) and H ( s , u ) , respectively, then N [ k * h ] = u   K ( s , u ) H ( s , u ) , where [ h * k ] is the convolution of h and k .

Example 3.1

Consider the FZK equation in the following form [13]: D t β φ + ( φ k ) x + 1 8 ( φ k ) x x x + 1 8 ( φ k ) y y x = 0 ,   0 ≤ t ,   0 < β ≤ 1 . (18) Together, the initial condition is φ ( x ,   y ,   0 ) = f ( x ,   y ) . (19)

Eq. (18) is written in the implicit form as D t β φ = − ( φ k ) x − 1 8   ( φ k   ) x x x − 1 8   ( φ k ) y y x . (20)

Using natural conversion to Eq. (20), we get N + [ D β ( φ ) ] = N + [ − ( φ k ) x − 1 8 ( φ k ) x x x − 1 8 ( φ k   ) y y x ] . (21) Utilize the differentiation characteristic of the natural convert as s β u β N + [ φ ] − s β − 1 u β φ ( x , y , 0 ) = N + [ − ( φ k ) x − 1 8   ( φ k ) x x x − 1 8   ( φ k ) y y x ] . (22)

Using the initial condition in Eq. (20) and rearranging, we have N + [ φ ] = f ( x , y ) s + u β s β N + [ − ( φ k ) x − 1 8 ( φ k ) x x x − 1 8 ( φ k ) y y x ] . (23) As φ ( x , y , t ) is the infinite series given as ∑ m = 0 ∞ φ m ( x , y , t ) , (24)

applying natural transform to Eq. (25) and using the idea explained in the method { N + [ φ 0 ( x , y , t ) ] = f ( x , y ) s , N + [ φ 1 ( x , y , t ) ] = u β s β N + [ − ( φ 0 k ) x − 1 8 ( φ 0 k ) x x x − 1 8 ( φ 0 k ) y y x ] , N + [ φ 2 ( x , y , t ) ] = u β s β N + [ − ( ( φ 0 + φ 1 ) k ) x − 1 8 ( ( φ 0 + φ 1 ) k ) x x x − 1 8 ( ( φ 0 + φ 1 ) k ) y y x − ( − ( φ 0 k ) x − 1 8 ( φ 0 k ) x x x − 1 8 ( φ 0 k ) y y x ) ] . ⋮ (25)Now, by taking the inverted natural transmute of Eq. 25, we obtain the solution elements as { φ 0 ( x , y , t ) = N − [ f ( x , y ) s ] , φ 1 ( x , y , t ) = N − [ u β s β N + [ − ( φ 0 k ) x − 1 8 ( φ 0 k ) x x x − 1 8 ( φ 0 k ) y y x ] ] , φ 2 ( x , y , t ) = N − [ u β s β N + [ − ( ( φ 0 + φ 1 ) k ) x − 1 8 ( ( φ 0 + φ 1 ) k ) x x x − 1 8 ( ( φ 0 + φ 1 ) k ) y y x − ( − ( φ 0 k ) x − 1 8 ( φ 0 k ) x x x − 1 8 ( φ 0 k ) y y x ) ] ] . ⋮ (26)The n − t e r m s approximate solution of Eqs. 18 and 19 by NTIM is presented as φ ( x , y ,   t ) = φ 0 + φ 1 + . . . + φ m − 1 . (27)

Example 3.2

Regarding the FZK (2, 2, 2) equality of the structure [13] D t β c φ + ( φ 2 ) x + 1 8 ( φ 2 ) x x x + 1 8 ( φ 2 ) y y x = 0 ,   0 < β ≤ 1 . (28) With preliminary conditions φ ( x , y , 0 ) = 4 3 λ sinh 2 ( x + y ) . (29) Here, λ is an optional fixed value. The accurate solution for β = 1.0 is given by φ ( x , y , t ) = 4 3 λ sinh 2 ( x + y − λ t ) . (30)

Utilizing the procedure of NTIM, we obtain the solution components for Eq. (28) as φ 0 ( x , y , t ) = 4 3 λ sinh 2 ( x + y ) , (31) φ 1 ( x , y , t ) = 8 λ 2 t β ( 4 ⁡ sinh ( 2 ( x + y ) ) − 5 ⁡ sinh ( 4 ( x + y ) ) ) 9 Γ ( β + 1 ) , (32) φ 2 ( x , y , t ) = 64 λ 3 t 2 β ( 3 Γ ( β + 1 ) 2 Γ ( 3 β + 1 ) ( 13 ⁡ cosh ⁡ ( 2 ( x + y ) ) − 70 ⁡ cosh ⁡ ( 4 ( x + y ) ) + 75 ⁡ cosh ⁡ ( 6 ( x + y ) ) ) − 20 λ Γ ( 2 β + 1 ) 2 t β ( 4 ⁡ sinh ⁡ ( 2 ( x + y ) ) + 8 ⁡ sinh ( 4 ( x + y ) ) − 60 ⁡ sinh ( 6 ( x + y ) ) + 85 ⁡ sinh ⁡ ( 8 ( x + y ) ) ) ) 81 Γ ( β + 1 ) 2 Γ ( 2 β + 1 ) Γ ( 3 β + 1 ) . (33)

Adding the elements, the second-order approximated solution can be written as φ ∼ ( x , y , t ) = φ 0 + φ 1 + φ 2 .

Example 3.3

Consider the FZK (3, 3, 3) equation of the structure [13] D t β c φ + ( φ 3 ) x + 2 ( φ 3 ) x x x + 2 ( φ 3 ) y y x = 0 . (34)

Together with initial conditions φ ( x , y , 0 ) = 3 2 λ ⁡ sinh ( x + y 6 ) . (35) Here, λ is an optional fixed amount. The exact solution for b e t a = 1.0 is given by φ ( x , y , t ) = 3 2 λ ⁡ sinh ( 1 6 ( x + y − λ t ) ) . (36)

Using the procedure of the NTIM, the solution elements can be acquired as φ 0 ( x , y , t ) = 3 2 λ ⁡ sinh ( 1 6 ( x + y ) ) , (37) φ 1 ( x , y , t ) = 3 λ 3 t β ⁡ cosh ( x + y 6 ) ( 7 − 9 ⁡ cosh ( x + y 3 ) ) 16 Γ ( β + 1 ) , (38) φ 2 ( x , y , t ) = [ − 3 λ 5 t 2 β 131072 Γ ( β + 1 ) 3 Γ ( 2 β + 1 ) Γ ( 3 β + 1 ) Γ ( 4 β + 1 ) × λ 4 Γ ( 2 β + 1 ) Γ ( 3 β + 1 ) 2 t 2 β ( 9 ( 1700 6 sinh ( 5 ( x + y ) ) + 9135 6 sinh ( 7 ( x + y ) ) − 3412 2 ( sinh ( x + y ) ) − 10935 2 ( sinh ( 3 ( x + y ) ) ) ) + 3550 6 sinh ( x + y ) ) + 16 λ 2 Γ ( β + 1 ) Γ ( 2 β + 1 ) 2 Γ ( 4 β + 1 ) t β ( 9 ( 75 2 cosh ( x + y ) + 1827 6 cosh ( 7 ( x + y ) ) − 1615 6 ( cosh ( 5 ( x + y ) ) ) ) + 1385 6 cosh ( x + y ) ) − 256 Γ ( β + 1 ) 3 Γ ( 3 β + 1 ) Γ ( 4 β + 1 ) ( 70 6 sinh ( x + y ) + 765 6 sinh ( 5 ( x + y ) ) − 621 2 ( sinh ( x + y ) ) ) ] . (39)Adding the components, the second-order series solution can be written as φ ∼ ( x , y , t ) = φ 0 ( x , y , t ) + φ 1 ( x , y , t ) + φ 2 ( x , y , t ) . (40)

In this work, two problems of the FZK equation have been tested by the new developed methodology NTIM. The obtained results are assessed by diverse plots and tabulated data for testing the reliability of the proposed method. Figure 1 shows the 3D surfaces obtained by the NTIM and the accurate result correspondingly for Example 3.2 in the 3D graph by keeping the y parameter constant. By keeping the time parameter constant, the approximate and accurate results are shown in Figure 2 respectively, for problem 3.2. In Figure 3, the absolute error is shown by a 3D plot by the variation of x and y parameters while time is kept constant. A comparison for the variation of the fractional value β is shown by 2D plots for Example 3.2 which shows the consistency of the method by agreeing to the amount of β tactics to the standard amount 1 of a differential equation; the approximated result converges to the accurate solution of the problem. Similarly, Figure 4 shows the approximates and exacted solution by variation of the x and t components, while Figure 5 shows the 3D graphs of the approximated solution and exacted solution by keeping the time parameter constant for problem 3.3. Figure 6 shows the absolute error of the NTIM result and the exacted solution for problem 2. Table 1 shows the comparison of data on the computational amounts of the approximated solution of the PIA and RPS methods for diverse amounts of β , while in Table 2, the absolute errors of our suggested methodology have been matched with the absolute errors of the PIA and RPS methods. Similarly, in Table 3, the fractional value of NTIM has been compared with the third-order RPSM and PIA methods. In Table 4, the absolute errors of NTIM, RPSM, and PIA methods have been compared. The approximate solution in this article is executed up to second order for both problems. The accuracy may be increased by obtaining a higher order of the approximate solution. From the tables and graphs, it is so far clear that the NTIM reveals encouraging approximated results as evaluated with other existing methodologies in the previously published works Table 5.

In the current investigation, the NTIM has been applied successfully to the FZK equations. Two problems have been tested. The proposed results reveal that the method handles the non-linear equations in a good way and provides an efficient approximate solution to non-linear PDEs. The numerical values of approximate and exact solutions through tables show the efficiency and reliability of the method. Also, the graphs verify the efficiency of the proposed method through 3D and 2D plots. The fractional approximation through 2D graphs also shows the consistency of the method by approaching the fractional value β of the equation to the conventional amount 1, so an approximate result converges to the exacted result of the problems. This strategy is also effective when the answer to the integer order model is unknown. As a result, we decided that the current technique is trustworthy and effective in obtaining estimated solutions for various classes of linear and non-linear fractional formulations of ordinary and partial differential equations.