Evolutionary algorithm based heuristic scheme for nonlinear heat transfer equations

In this paper, a hybrid heuristic scheme based on two different basis functions i.e. Log Sigmoid and Bernstein Polynomial with unknown parameters is used for solving the nonlinear heat transfer equations efficiently. The proposed technique transforms the given nonlinear ordinary differential equation into an equivalent global error minimization problem. Trial solution for the given nonlinear differential equation is formulated using a fitness function with unknown parameters. The proposed hybrid scheme of Genetic Algorithm (GA) with Interior Point Algorithm (IPA) is opted to solve the minimization problem and to achieve the optimal values of unknown parameters. The effectiveness of the proposed scheme is validated by solving nonlinear heat transfer equations. The results obtained by the proposed scheme are compared and found in sharp agreement with both the exact solution and solution obtained by Haar Wavelet-Quasilinearization technique which witnesses the effectiveness and viability of the suggested scheme. Moreover, the statistical analysis is also conducted for investigating the stability and reliability of the presented scheme.


Introduction
Most of the real world problems especially heat transfer equations which emerge in many scientific and engineering fields possess nonlinear behavior, are modeled by nonlinear ordinary differential equations and hence attained an ever increasing attention of scientists and engineers towards the solution of these nonlinear problems. Consequently several analytical and numerical techniques have been suggested by researchers for solving nonlinear heat transfer equations as quoted under references [1][2][3][4][5][6][7][8][9][10].
The heat transfer equations have major practical importance in cooling of electronic components, cooling of heated stirred vessels and heated parts of the space vehicles etc [4].
This work is mainly concerned to present a heuristic scheme which is stochastic in nature for the study of two nonlinear heat transfer equations. First equation refers to a well known boundary value temperature distribution equation in a lumped system of combined convection-radiation that can be described by the following mathematical model: With the boundary conditions: x 0 ð0Þ ¼ a and xð1Þ ¼ b While second equation is regarded as an initial value cooling equation of a lumped system by combined convection and radiation whose governing equation is given by: With the initial condition: where a, b, c and δ are real valued constants. Researchers put forward an extensive collection of techniques for obtaining the numerical solution of such nonlinear heat transfer equations. To highlight a few, Ganji et al. [1] tailored nonlinear heat transfer equations by using homotopy perturbation method (HPM) and variational iteration method (VIM). A. Atangana in [2] adopted a new iterative analytical technique along with the stability, convergence and the uniqueness analysis of adopted technique for dealing with nonlinear fractional partial differential equations arising in biological population dynamics system. Yang Juan-Cheng et al. in [3] employed Direct Numerical Simulation (DNS) to look at the mechanism of heat transfer enhancement (HTE) with the findings that DNS procedure is trustworthy. Abbasbandy [4] encountered nonlinear heat transfer equations by using homotopy analysis method (HAM). Kolade M. Owolabi et al. in [5] estimated the fractional Schrodinger equation with the Riesz fractional derivative by making use of the improved exponential time differencing Runge-Kutta scheme with fourth-order accuracy and realized a different distribution of the complex wave functions both for the focusing and defocusing cases. R. A. Khan in [6] adopted generalized approximation method (GAM) and found that GAM perform quite effectively irrespective of dependency on small parameter as in case of HPM and PM. Dehghan et al. [7] employed semi-analytical methods i.e. homotopy perturbation method (HPM) and finite difference method (FDM) which enhances heat transfer and validates semianalytical methods for study of heat transfer rates while Saeed et al. [8] obtained the solutions of nonlinear heat transfer equations by using Haar Wavelet-Quasilinearization Technique.
Because of limitations in accuracy and efficiency of these classical analytical and numerical methods, the evolutionary computing techniques (ECT) have demonstrated their importance, investigated thoroughly and applied successfully to solve nonlinear problems in engineering and applied science domain [11][12]. For instance, Arqub et al. [13,14] considered continuous GA based technique for dealing with boundary value problem, Troesch's and Bratu's Problems. Malik et al. [15,16] taken into account an evolutionary computing scheme of hybrid genetic algorithm (HGA) for solving biochemical reaction and Singular Boundary Value Problems Arising in Physiology. Kadri et al. [17] solved nonlinear heat conduction problems by using genetic algorithm (GA). Khan et al. [18,19] solved fractional order system of Bagley-Torvik equation and differential equations of first orders by using HGA and swarm intelligence approaches with the findings that the proposed evolutionary computing techniques are reliable and effective. Arqub et al. in [20] investigated the efficiency, accuracy and convergence analysis of the continuous genetic algorithm by solving a class of nonlinear systems of secondorder boundary value problems.
In this paper, nonlinear heat convection-radiation equations have been treated numerically by using a heuristic scheme which is stochastic in nature and hybridized with Log Sigmoid and Bernstein Polynomial basis functions with unknown coefficients. The equation representing nonlinear heat transfer phenomena is converted into an error minimization problem. The hybrid approach of genetic algorithm is exploited for solving the error minimization problem and to obtain the unknown coefficients that further gives the numerical solution of the problem under consideration.

Methodology overview
In this section, the research methodology for the hybridization of Evolutionary Algorithm (EA) with two different basis functions i.e. Log Sigmoid and Bernstein polynomials for treating nonlinear heat transfer equations is presented in detail as under.

Methodology overview for EA hybridization with Log Sigmoid Basis Function
For the hybridization of EA with Log Sigmoid Basis Function, we assume that the approximate numerical solution x(t) of Eq (1) and its first and second derivatives that is x 0 (t) and x@(t), is linear Combination of some log sigmoid basis functions that can be represented by the following equations:- x 0 ðtÞ ¼ Where ϕ(t) is log sigmoid function defined by (α i , ω i , β i ) are real valued unknown parameters and k is the number of basis functions. The given nonlinear ODE of Eq (1) is converted into an error minimization problem that provides the required unknown parameters (α i , ω i , β i ) as follows: Where N is the total number of steps taken in the solution range [0, 1]. ε 1 is the mean of the sum of square error of the given system (1), ε 2 is the mean of the sum of square error due to boundary conditions of Eq (1) and j is the number of generations executed. The minimization of Eq (9) is carried out by using Genetic Algorithm. The optimal values of unknown parameters (α i , ω i , β i ) resulting in minimum ε j are used in Eq (3) that gives the approximate numerical solution of Eq (1). The same methodology applies for obtaining the approximate numerical solution of Eq (2) using Log Sigmoid Basis Function.

Methodology overview for EA hybridization with Bernstein Polynomial Basis Function
For the hybridization of EA with Bernstein Polynomial Basis Function, again we assume that the approximate numerical solution x(t) of Eq (1) and its first and second derivatives that is x 0 (t) and x@(t), is linear combination of Bernstein Polynomial (B-Polynomials) basis functions of degree k and can be represented by the following equations:- Where (α 0 , α 1 ,. . ...α k ) are unknown parameters and k is the degree of B-polynomials. The given nonlinear ODE of Eq (1) is then converted into an error minimization problem that provides the required unknown parameters (α 0 , α 1 ,. . ...α k ) by using the same Eqs (7-9). The minimization of fitness function ε j as defined by Eq (9) for B-Polynomials basis functions is carried out by using Genetic Algorithm. The optimal values of unknown parameters (α 0 , α 1 ,. . ...α k ) resulting in minimum ε j achieved are then used in Eq (10) that gives the approximate numerical solution of Eq (1). The same methodology applies for obtaining the approximate numerical solution of Eq (2) using Bernstein Polynomial Basis Function.

Heuristic optimization technique
Genetic algorithm is a widely known global heuristic search optimization technique in evolutionary computing field. Genetic algorithm is inspired by the Darwinian principle of evolution and natural selection in which stronger individuals are likely to be the winners in a competing environment.Genetic algorithm operates on a population of individuals referred to as chromosomes.Each chromosome represents a possible solution to the problem and assigned a real numbered fitness value, which is a measure of the excellence of the solution to the particular problem. Genetic algorithm starts with a randomly generated population of chromosomes, carries out a process of fitness based selection and recombination to produce the next generation. During recombination, two or more selecting parent chromosomes recombine to produce child chromosomes. This process is repeated iteratively and recursively by means of genetic operators; selection, crossover and mutation with a hope that the average fitness of the chromosomes increase until some defined stopping criteria is fulfilled. In this way, Genetic algorithm searches for the best optimum possible solution to the given problem [11][12].
Interior-Point Algorithm (IPA) is a local search optimizer which is used extensively in variety of optimization problems. IPA solves Karush-Kuhn-Tucker (KKT) equations in the system by applying either newton step or conjugate gradient (CG) step iteratively to optimize problem's defined merit function [11,16].
Hybrid approach of evolutionary algorithms is efficient in the sense that they are accurate and fast convergent. In hybrid approach, the optimum chromosomes found by global optimizer (GA) are given to the local optimizer (IPA) as starting point which improves the solution provided by global optimizer (GA) [11,16].
The steps for the hybridized scheme of GA with IPA used in this work are outlined below in Table 1 [16].

Implementation of proposed scheme
In this section, the proposed scheme is implemented on nonlinear heat transfer equations for the validity of the suggested approach by using both Log Sigmoid and Bernstein Polynomial basis functions. MATLAB tool is used for systems simulations.

Implementation of proposed scheme using Log Sigmoid basis function
Problem 01. Let us consider temperature distribution equation in combined convection-radiation lumped system: Let the system has volume υ, surface area C, density ρ, specific heat η, η c is specific heat at temperature T c , initial temperature T 0 , temperature of the convection environment T c and heat transfer coefficient ℏ. The temperature distribution equation of the system is given by nonlinear boundary value problem of (1) as [8]: The boundary conditions are: x 0 ð0Þ ¼ 0 and xð1Þ ¼ 1 The approximate solution of Eq (13) is obtained in the range [0, 0.8] with increment of 0.2 and basis function k = 10. The results obtained are restricted within the bounds [-25, +25]. Two cases for Eq (13) have been treated by considering δ = 0.6 and δ = 2.0.
The fitness functions formulated for Eq (13) for both the cases i.e. for δ = 0.6 and δ = 2.0 are as follows: Where x(t), x 0 (t), x@(t) are same as defined by Eqs (3-5) respectively. Table 1. Steps for hybridization of GA with IPA.

Step 1 (Population Initialization)
Random population of N individuals or chromosomes having M genes per chromosome is generated in a bounded limit.

Step 2 (Evaluation & Ranking)
Evaluate each individual using problem specific fitness function and rank them proportionate to their fitness value.

Step 3 (Stoppage Criteria)
The algorithm keeps executing until and unless some user defined stoppage criteria is met. If the stopping criterion is satisfied then go to step 6, else repeat steps 2 to 5.

Step4 (Selection & Reproduction)
Based on fitness value, the chromosomes from current population are chosen as parents for new generation. These parents then produce further offspring as a result of crossover operation which became parents for next generation.

Step 5 (Mutation)
Mutation operation is an optional operator and execute only if no improvement in the fitness value of the generation is seen. It randomly changes the offspring resulted from crossover to find a good solution.
Step 6 (Local Search Improvement) The optimum chromosomes found by GA fed to local optimizer IPA as starting point which tends to optimize the results further. https://doi.org/10.1371/journal.pone.0191103.t001 The fitness functions in Eqs (14 and 15) are functions of unknown parameters (α i , ω i , β i ) and these unknown parameters are achieved by employing GA, IPA and hybrid scheme of GA-IPA for minimizing Eqs (14 and 15) and accordingly the approximate solution of Eq (13) is obtained for both the values of δ. The parameter settings used for the execution of GA and IPA corresponding to minimum fitness ε j1 , ε j2 are given in Table 2. The fitness functions have been minimized by fixing the size of chromosome i.e. total number of unknown parameters equal to 30. The values of unknown parameters obtained by GA, IPA and GA-IPA are provided in Tables 3 and 4 for case 1 and case 2 respectively. The approximate solutions achieved by GA, IPA and GA-IPA for both the cases are provided in Tables 5 and 6. The comparisons of absolute errors necessary to demonstrate the applicability of the suggested approach are provided in Tables 7 and 8.
The comparisons of absolute errors given in Tables 7 and 8 for both the cases i.e. δ = 0.6 and δ = 2.0 reveals that the proposed scheme dealt the given nonlinear problem with quite excellent accuracy than GA and HPM while the results are quite compromising with the exact results and results obtained by Haar Wavelet method which can be made more competitive by increasing the number of basis functions.
Problem 02. Let us consider case of cooling of a lumped system by combined convection and radiation: Let the system has volume υ, surface area C, density ρ, specific heat η, η c is specific heat at temperature T c , initial Temperature T 0 , temperature of the convection environment T c , emissivity μ, heat transfer coefficient ℏ and sink temperature is represented by T α as the system loses heat through radiation. The mathematical model of the system is given by the  Table 3. Unknown parameters achieved by GA, IPA and GA-IPA for δ = 0.6. nonlinear initial value problem of Eq (2) as [8]: The approximate solution of Eq (16)  The fitness function is formulated as given by: GA and IPA are implemented using the settings as given in Table 9. The values of unknown parameters corresponding to minimum ε j for different values of δ are omitted in this case. These unknown parameters are necessary to known for obtaining the approximate solutionxðtÞ of Eq (17). The results achieved by GA, IPA and GA-IPA are shown in Table 10. The comparison of absolute errors provided by our method is made with available methods and is given in Table 11. The comparison of absolute errors given in Table 11 above shows that the absolute errors generated by proposed method are quite smaller than errors generated by available numerical methods i.e. VIM and HPM while the generated errors are in quite good comparison with the errors generated by Haar Wavelet method [8]. However, by increasing the number of basis

Implementation of proposed scheme using Bernstein Polynomial basis function
Problem 01. The same system of nonlinear boundary value problem along with boundary value conditions is considered as for Log Sigmoid basis function which is: The boundary conditions are: x 0 ð0Þ ¼ 0 and xð1Þ ¼ 1 All the settings and parameters remains the same as used for Log Sigmoid basis function for tuning the algorithms. B-Polynomial degree of 8 i.e. k = 8 is opted to use. Two cases for Eq (18) have been treated by considering δ = 0.6 and δ = 2.0. The fitness functions formulated for Eq (18) for both the cases i.e. for δ = 0.6 and δ = 2.0 are as follows: ðx 00 ðt i Þ À 0:6x 4 ðt i Þ 2 Þ þ 1 2 ððx 0 ð0ÞÞ 2 þ ðxð1Þ À 1Þ 2 Þ ð19Þ Where x(t), x 0 (t), x@(t) are same as defined by Eqs (10-12) respectively. The fitness functions in Eqs (19 and 20) are functions of unknown parameters (α 0 , α 1 ,. . ... α k ) that are necessary to be known for obtaining the optimum numerical solution to given nonlinear ODE. The values of unknown parameters obtained by GA, IPA and GA-IPA are  Tables 12 and 13 for case 1 and case 2 respectively. The approximate solutions achieved by GA, IPA and GA-IPA for both the cases are provided in Tables 14 and 15. The comparisons of absolute errors necessary to demonstrate the applicability of the approach are provided in Tables 16 and 17. From the comparison of results shown in Table 17, it is clear that the proposed method has provided optimum results as compared to numerical methods GA, HPM and Haar Wavelet method. Problem 02. Again the same system, considered for Log Sigmoid basis function, which is represented by the nonlinear initial value problem of Eq (2) is taken into account as: The initial condition is For tuning the algorithms using Bernstein Polynomial basis function, the settings and parameters used are same as for Log Sigmoid basis function shown under section 4.1. The fitness function formulated for Eq (21) is given by: Where x(t) and x 0 (t) are same as defined by Eqs (10 and 11) respectively. The values of unknown parameters corresponding to minimum ε j for different values of δ are omitted in this case. The results provided by GA, IPA and GA-IPA are shown in Table 18. The comparison of absolute errors generated by proposed method is made with available numerical methods and is provided in Table 19.
It is evident from the comparison of absolute errors given above in Table 19 that the absolute errors generated by proposed method are quite smaller than errors generated by available numerical methods i.e. VIM and HPM which testifies the viability and affirms the accuracy of proposed approach while the generated absolute errors are in good comparison with the errors generated by Haar Wavelet method [8] which can be further improved by increasing the degree k of basis function.

Statistical analysis for proposed approach
In this section, statistical analysis has been performed for nonlinear heat transfer problems considered in this paper both with log sigmoid and B-Polynomial approaches to demonstrate the stability and authenticity of the proposed methodology. For this sake, 10 independent runs of the proposed methodology are performed without any deviation in the parametric settings of the problem as shown in Tables 2 and 9 Table 20.
From Table 20, It can be seen that the mean of the average absolute error lies within the range 10 −05 to 10 −03 while the STD of the absolute error lies within the range 10 −06 to 10 −03 for Log Sigmoid based approach. Similarly, the mean of the average absolute error lies within the range 10 −06 to 10 −05 while the STD of the absolute error lies within the range 10 −12 to 10 −04 for B-Polynomial based approach. The very close ranges of mean and STD determines the close

Conclusion
On the basis of simulations and results provided by the proposed approach, it is observed that the proposed hybrid heuristic approach of GA-IPA out performs the other existing solutions for the nonlinear heat transfer equations. Further, the proposed scheme shows the supremacy over numerical methods i.e. Generalized Approximation, Homotopy Perturbation Method and Variational Iteration Method. In some cases the results are even more satisfactory than results provided by the numerical method Haar Wavelet Quasilinearization Technique. In this connection, it is established that the proposed approach is a good and trusty alternative approach for researchers for solving nonlinear problems in engineering and applied sciences domain.