Solutions of fractional order differential equations modeling temperature distribution in convective straight fins design

*Correspondence: msulaiman@awkum.edu.pk; poom.kum@kmutt.ac.th 1Department of Mathematics, Abdul Wali Khan University, Mardan 23200, KP, Pakistan 2KMUTT Fixed Point Research Laboratory, KMUTT-Fixed Point Theory and Applications Research Group, SCL 802 Fixed Point Laboratory, Department of Mathematics, Faculty of Science, King Mongkut’s University of Technology Thonburi (KMUTT), 126 Pracha Uthit Rd., Bang Mod, Thung Khru, Bangkok 10140, Thailand Full list of author information is available at the end of the article Abstract


Introduction
A vast number of problems which model physical phenomena, for example heat transfer, involve the nonlinear function [1]. In mechanical engineering, heat transfer is a very common science because in various objects it can be required. The problems of improvement of heat transfer are solved on an extended surface which is known as fin. The heat transfer mechanism of fin is to conduct heat through its thermal conduction from the source of heat to the fin's surface, and then heat is dissipated into the air through the effects of ther-fins in [5]. The heat transfer dynamics of fins in space radiator and one-dimensional radiation fins are further investigated in [6][7][8]. The system of heat-rejecting comprising parallel tubes connected by web plates has been studied by Bartas and Sellers [9]. In [10], Hug and Aziz used a perturbation-based technique to find the closed form solutions for straight convective fin with thermal conductivity which is dependent on temperature. To find analytical solution for dimensionless temperature and investigate the efficiency of the fin having thermal conductivity dependent on temperature, Arslanturk [11,12] used Adomian decomposition technique.
In recent years, studying the heat transfer on extended surfaces has become very pivotal with growing importance of performance of heat transfer fins having lower volumes, weights, initial and operating costs of the systems [11]. In addition to the developments in standard methods of numerical computation, a new methodology, which is known as homotopy perturbation Sumudu transformation technique, was proposed to analyze less or strongly nonlinear systems. HPSTM has recently been used for solution of nonlinear fractional equation of gas dynamics and some other physical phenomena [13]. Sumudu transform method, homotopy perturbation method (HPM), and He's polynomials are combined to design the HPSTM technique [14]. Fractional calculus [15,16] is a branch of applied mathematics which deals with arbitrary order differentiation and integration. It has found many applications in different areas of science and engineering over the last three decades [17][18][19][20]. The HPSTM is used for the solution of energy balance equation of fractional order [21]. More works on analytical and numerical techniques for the solving integer and fractional differential equations are available in .
In this work, we have designed a hybrid of ANNs and biogeography-based heterogeneous cuckoo search algorithm (BHCS) for the solution of integer and fractional order energy balance equation in order to analyze the temperature distribution in convective straight fins. We have named our algorithm the ANN-BHCS algorithm. We have analyzed seven cases of the problems and the results are compared with other techniques such as HPM, VIM, and HPSTM [6,21,52]. The results show that the ANN-BHCS algorithm is better than other techniques in terms of obtaining the solutions with high accuracy. We have further tested the efficiency of the ANN-BHCS algorithm by utilizing the performance metrics MAD, TIC, and ENSE.  Consider a straight convective fin having thermal conductivity dependent on temperature, and it has an arbitrary constant area of cross section A c , length b, and perimeter P, and the heat transfer coefficient is denoted by h. The fin is associated with the temperature T at the base surface, and extends into temperature T a of the fluid. The geometry of the straight fin is given in Fig. 2. The one-dimensional energy balance equation is given by [52,53] The thermal conductivity of fin's material is considered as a linear function according to Eq. (2), where k is the parameter defining the variation of the thermal conductivity and k a is the thermal conductivity at the ambient fluid temperature of the fin. Introducing the dimensionless parameters: Now Eq. (1) reduces to the following equation: with the following boundary conditions: The computational domain 0 ≤ x ≤ b is transformed to 0 ≤ ξ ≤ 1 by introducing the dimensionless parameters given in Eq. (3).
To understand the anomalous behavior of this system, we fractionalize the energy balance Eq. (4) into fractional order (v > 0) as follows in order to find fin temperature in straight fins: with the following boundary conditions as in Eq. (5).

Basic definitions
This section consists of some definitions and important relations from fractional calculus that have been used in the construction of ANN for FDEs. Fractional derivatives and integrals have been expressed in different ways in literature, i.e., Riemann-Liouville, Caputo, Erdélyi-Kober, Hadamard, Grünwald-Letnikov, and Riesz type etc. In standard fractional calculus, equivalence of these definitions for some functions has been given [15,54,55]. All of these definitions have their own importance and advantages in different kinds of problems in mathematics. Definitions of Riemann-Liouville and Caputo fractional derivatives are given below.
Definition 1 (The fractional order Riemann-Liouville integral and derivative) The integral of fractional order v > 0 can be written as [16] Here, I v shows the fractional integral of order v. The fractional derivative of order v > 0 is normally given as Here, D v represents the fractional derivative of order v and n is an integer.
Definition 2 (Caputo fractional derivatives) There are some limitations of the definition of fractional derivatives given by Riemann-Liouville, when it is used for modeling of some real world phenomena related to differential equations of fractional order. Therefore, a modified definition for fractional differential operator D v is introduced by Caputo [16,56]: where I v is given in Eq. (3). Caputo integral operator is given by The ordinary derivative followed by a fractional integral gives the Caputo fractional derivative, while the calculation in reverse order gives the Riemann-Liouville derivative. Using Caputo fractional derivative, we can use the traditional homogeneous and nonhomogeneous initial/boundary conditions occurring in general applications. However, for homogeneous initial conditions, Riemann-Liouville and Caputo formulations coincide [16,57].
Definition 3 (Mittag-Leffler function (MLF)) The Mittag-Leffler function (MLF) is one of the most important functions having widespread applications in fractional calculus. It plays an important role in the solution of differential equations of integer and fractional orders because of its exponential nature.
The classical MLF has the definition as given below [58]: It becomes the exponential function when α = 1. The MLF function with two parameters α and β is as follows: For β = 1, it becomes a standard MLF function.

ANN modeling
This section presents the mathematical modeling of artificial neural networks (ANN) for differential equations of fractional order. Neural networks modeling has already been implemented for solving integer order differential equations. Now, we model the ANN to find the solution for fractional order differential equations. The exponential function is used as an activation function for ANN. It has the capability to approximate the functions and its fractional derivative is also calculated with terms represented by classical MLF. The fractional derivative of an exponential function can be written as Approximate solution for the problem considered in this research and its v order derivative is given bŷ Equations (15) and (16) are used to approximate a solution of the fractional order differential equation given in Eq. (6). The neural networks architecture for fractional differential equations is given in Fig. 3. The objective function for the problem considered in this paper is given by where E 1 and E 2 are given by Here, E 1 is related to the differential equation and E 2 is related to the initial and boundary conditions. We try to find the weights α i , ω i , and β i in Eq. (15) such that E 1 and E 2 approach zero, then E will also approach zero. Hence the approximate solutionθ (ξ ) will approach the exact solution θ (ξ ).

Cuckoo search
Inspired by the cuckoo bird's breeding behavior, a metaheuristic algorithm was developed which is called the cuckoo search algorithm [59]. The female bird lays eggs in other host birds' nests and they unintentionally raise her brood. When the host bird finds the egg of the cuckoo bird in her nest, it either throws it out of the nest or starts making her own brood elsewhere [60].
In the cuckoo search algorithm, the solution is represented by the egg of the host bird and the new candidate solution is represented by cuckoo's egg. There are three rules that are described for cuckoo search and those are [61]: (1) the cuckoo lays a single egg at a time and puts it in the host's nest; (2) the nests that have a high quality egg, i.e., a better solution will go to the next generation; and (3) there is a fixed number of host nests, and the host bird can find an alien egg with certain probability. Assuming . , x iD ) as the position for the ith egg (solution) then updated solution x i new is generated by Levy flights as given below: where the product ⊕ is entry-wise multiplication; the Levy flight exponent is denoted by β; the step size for a cuckoo is determined by a positive parameter α; the best solution within the current population is denoted by x g ; u and v are random numbers: where is used for gamma function, and β controls the value of σ u . There is a discovery operator in CS which is used to replace the discovered nests with a probability (pa). The equation that is used to update the solution is given as follows: where x ij new is the jth element of the ith solution x i new ; x r1,j and x r2,j are the jth elements of the two solutions x r1 and x r2 , where r1 and r2 are two different integers in interval [1,NP], where NP represents size of population, pa represents the discovery probability, P and rand are some random numbers that belong to the interval [0, 1].

Biogeography-based optimization
Biogeography-based optimization (BBO) is an evolutionary algorithm which is inspired by different characteristics of species living in the islands [62]. In BBO, each habitat is considered as a candidate solution having some habitat's suitability index (HSI), which is employed for measurement of the quality of a habitat. A habitat (solution) is represented by some suitability index variables (SIV). Two types of operators, i.e., migration and mutation, are used in BBO that are employed for the evolution of the population. In migration process, the solutions with high HSI share their characteristics with the solutions having low HSI and the solutions with low HSI accept new characteristics from the solutions with high HSI. In BBO, population is randomly initialized with NP habitats (solutons). Each generation sorts the population from the best to the worst and immigration and emigration rates λ and μ respectively are assigned to each habitat: where immigration (I) and emigration (E) rates are such that I = E = 1; S i represents the number of species of the habitats and S i = NPi. Accordingly, for the best solution the S i value is NP -1, and for the second best solution the S i value is NP -2, and for the worst solution the S i value is 0. The migration mixes the features within the population that modifies the solutions. After migration, to modify the solutions, BBO also uses the mutation operator.

Heterogeneous cuckoo search algorithm based on BBO
CS and BBO are hybridized because CS uses the Levy flights to modify the solutions as it is good at exploration, and BBO modifies the solutions using the migration operator as it is good at exploitation. Combining the exploration and exploitation, a hybrid metaheuristic algorithm is developed which is known as BBO-based heterogeneous cuckoo search (BHCS) algorithm. The proposed BHCS algorithm has two main stages that are the heterogeneous cuckoo search and the discovery based on biogeography. The details of these two stages are explained in the next section.

Heterogeneous cuckoo search strategy
At first stage, the BHCS algorithm uses the Levy flights and quantum mechanism based heterogeneous cuckoo search. This strategy is inspired by quantum mechanism and was first presented in [60,63]. The rules to update the solutions by heterogeneous cuckoo search are given as follows [60,63]: where L = δ ln(1/η), ε = δ exp(η), x g is used for the best solution at the current iteration; x i represents the mean of all solutions; sr and η are random numbers in the interval [0, 1]. Equation (25) shows that heterogeneous cuckoo search employs three equations to update the solutions with the same probabilities. The first equation is related to Levy flights in original cuckoo search and the second and third equations to update the solutions are based on quantum mechanism. Updating the solutions using heterogeneous rules diversifies the search and follows the direction towards the real global region.

Biogeography-based discovery operator
At the second stage, new solutions are generated using a discovery operator. When the host bird finds an alien egg with probability pa, it abandons the old nest and starts making a new nest based on the migration operator.
Initially, solutions are listed from the best to worst, and an immigration rate μ is assigned to each solution: where E = 1 represents the maximum emigration rate; S i = NPi represents the number of species in solutions. In biogeography-based discovery operator, those solutions having best fitness value share their characteristics with other solutions, which helps to enhance the exploitation.

Overall BHCS algorithm
The BHCS algorithm uses a cascading structure for the implementation of its two steps. The cooperation between heterogeneous search strategy and biogeography-based discovery operator can efficiently balance the exploitation and exploration.

Performance metrics
We have performed 100 simulations on all four problems to establish the stability, adaptability, and certainty of the BHCS algorithm. For this purpose, we have determined the mean absolute deviation (MAD) in solutions, root-mean-square error (RMSE), error in Nash-Sutcliffe efficiency (ENSE), Theil's inequality coefficient (TIC), and Nash-Sutcliffe efficiency (NSE). The analytical definition of these indexes are provided in Eqs. (27)-(30),

Results and discussion
In this paper, the biogeography-based heterogeneous cuckoo search (BHCS) algorithm is used to train the ANN model for the solution of fractional differential equations. We have considered the fractional form of the energy balance equation (6) in order to find the temperature in straight fins. We have considered the cases with integer and fractional orders. The problems are solved using the ANN modeling given in Eqs. (15) and (16). In our ANN network, we have taken 10 neurons with 30 unknown weights α, ω, and β. The training is performed over the interval [0, 1] with a step size of 0.1. The domain has 11 grid points.
Numerical solutions for case 2 with different values of thermal conductivity β are given in Table 7. Solutions are also plotted in Fig. 9(a). In this case, the thermo-geometric fin parameter is taken as ψ = 0.5 and the thermal conductivity β is varied from -0.5 to 0.5 with a step size of 0.2. From Table 7 and Fig. 9(a), we see that the dimensionless temperature θ increases with increase in the value of β. The accuracy of the ANN-BHCS algorithm can be seen from convergence of the fitness values and histograms of fitness values in Fig. 9.
Numerical solutions of case 4 for different values of β are given in Table 9 and Fig. 13(a). From Fig. 13(a), it is clear that the dimensionless temperature θ increases as the value of thermal conductivity β goes from -0.5 to 0.5. Convergence of the fitness values for case 4 is given in Fig. 13(b). Histograms of the fitness values for different values of β are given in Figs. 13(c)-13(h). The figures show that more than 90% of the fitness values are very close to zero, which shows the efficiency of the ANN-BHCS algorithm.
Numerical solutions of case 5 for different values of β are given in Table 10 and Fig. 15(a). From Fig. 15(a), it is clear that the dimensionless temperature θ increases as the value of thermal conductivity β goes from -0.5 to 0.5. Convergence of the fitness values for case 5 is given in Fig. 15
Numerical solutions of case 6 for different values of β are given in Table 11 and Fig. 17(a). From Fig. 17(a), it is clear that the dimensionless temperature θ increases as the value of thermal conductivity β goes from -0.5 to 0.5. Convergence of the fitness values for case 6 is given in Fig. 17
Numerical solutions for case 7 are given in Table 12. Solutions are also plotted in Fig. 19(a) for different values of β. The solution figures show that dimensionless temperature θ increases with increasing the values of thermal conductivity from -0.5 to 0.5. Convergence of the fitness values for case 7 is given in Fig. 19

Sensitivity analysis of parameters
In this section, we analyze the sensitivity of different parameters which are number of neurons, population size, and discovery probability pa. The ODE in case 1 with ψ = 0.2 is solved for sensitivity analysis of parameters. Solutions and absolute errors for different values of parameters are given in Tables 13, 14

Conclusion
In this paper, we have used ANN-based biogeography-based heterogeneous cuckoo search algorithm (ANN-BHCS) to analyze the problem of temperature distribution for convective straight fins. The ANN-BHCS algorithm is an efficient technique for the solution of fractional differential equations. We have considered seven cases of the problem, in which three cases are of integer order and four cases are of fractional order. In the first case, we have solved integer order energy balance equation for different values of thermogeometric fin parameter ψ. The series solutions for case 1 are given in Eqs. (34), (35) and (36). The comparison of exact and numerical solutions obtained by HPM, VIM, HPSTM, and ANN-BHCS for case 1 is given in Tables 1-6. Solutions obtained for case 1 are also plotted in Fig. 5(a). The results show that the ANN-BHCS algorithm gives better solutions        The results show that the ANN-BHCS algorithm can efficiently solve the integer and fractional order differential equations.