#### Abstract

The appearance of nonlinear equations in science, engineering, economics, and medicine cannot be denied. Solving such equations requires numerical methods having higher-order convergence with cost-effectiveness, for the equations do not have exact solutions. In the pursuit of efficient numerical methods, an attempt is made to devise a modified strategy for approximating the solution of nonlinear models in either scalar or vector versions. Two numerical methods of second-and sixth-order convergence are carefully merged to obtain a hybrid multi-step numerical method with twelfth-order convergence while using seven function evaluations per iteration, resulting in the efficiency index of about 1.4262. The convergence is also ascertained theoretically, and the asymptotic error constant is computed. Furthermore, various numerical methods of varying orders are used to compare the numerical results obtained with the proposed hybrid method in approximate solution, number of iterations, absolute error, absolute functional value, and the machine time measured in seconds. The obtained results outperformed the chosen methods when applied models arising from physical and natural fields were solved. Finally, to observe the convergence graphically, some complex polynomials are plotted as polynomiographs, wherein the rapid convergence of the proposed method is guaranteed.

#### 1. Introduction

Computing the approximate zeros of the nonlinear scalar and vector functions is one of the most important and interesting research areas in the modern age. There are many applications of the root-finding methods in different disciplines of science as well as in arts and economics. With the help of several mathematical techniques, a variety of complex problems in different applied sciences can be modulated in the form of nonlinear equations and then can be readily solved via different root-finding techniques. A root-finding method in mathematics and computer technology is a method for finding zeros, commonly known as ”roots” of continuous functions. A zero of a real-valued or a complex-valued function , is a value such that . Mostly, the root-finding techniques give approximations to zeros, expressed either as floating-point integers or as tiny isolating intervals, or discs for real or complex roots, because the zeros of a function cannot be calculated accurately with available analytical techniques.

Most numerical approaches for root-finding rely on the iteration process, which generates a series of discrete points that should converge towards the root as a limit. These iteration schemes start with one or more estimations of the root as initial inputs, and every iteration of the process generates a more accurate estimation of the exact root [1–4].

Since iterations must be ended at some point, these approaches yield estimation to the root rather than a precise solution. Many approaches compute successive approximations by considering an auxiliary function on the values that came before them. Therefore, the limit is a fixed point of the function , which is selected to have the solution of the original equation as fixed points and to converge to these fixed points quickly. In pursuing accurate and efficient roof-finding techniques, several techniques have been proposed in the past and current literature from the Newton method through the techniques proposed in the ongoing year.

Perhaps, the most commonly used root-finding method is the Newton Rahpson method [5] with quadratic convergence. Its computational step is shown below that uses two function evaluations: 1 for the function itself and 1 for the first-order derivative :where .

The researchers in [6] devised a modified version of an existing algorithm aiming at the removal of first-order derivatives. They came up with a two-step method with fourth-order convergence abbreviated by . One of the advantages of the algorithm was the use of only four function evaluations per iteration, as depicted in the following computational scheme:where

The researchers in [7] devised a three-step method having eighth-order convergence as denoted by . Despite being eighth-order convergent, one of the advantages of the algorithm was the use of only four function evaluations per iteration, as depicted in the following computational scheme:

A ninth-order convergence for a method as denoted by was proved in [8]. The algorithm merely needs 5 function evaluations per iteration, as depicted in the following computational scheme:

Another higher-order three-step iterative method as denoted by is in [9]. The method is shown to be three-step that requires only 6 function evaluations per iteration, as depicted in the following computational scheme:for

We have also chosen an iterative process with eleventh-order convergence as denoted by in [10]. The process consists of four steps and requires 7 function evaluations per iteration, as depicted in the following computational scheme:

Finally, a three-step iterative method with twelfth-order convergence as denoted by is in [11, 12]. The method consists of 5 function evaluations per iteration, as depicted in the following computational scheme:

#### 2. Formulation of the Proposed Method

It has been observed in the current literature that new modified root-finding techniques are being proposed because of increasing the efficiency of the existing ones. In search of such algorithms, some researchers have merged two existing iterative methods of convergence order and to obtain a hybrid method with convergence order . In this respect, the convergence order is improved. Nonetheless, the computational aspect was ignored, resulting in an increased number of function evaluations in most newly modified hybrid approaches. For example, authors in [13] proposed an iterative third-order method with five function evaluations required per iteration, including another algorithm presented in the same research work with a fourth-order three-step method that requires eight function evaluations. Likewise, authors in [14,15] have employed an excessive number of first-order derivatives, leading to high computational effort and machine time. The primary concern of the present work is to propose a hybrid method with possible higher-order convergence with the minimum number of function evaluations so that the computational cost in terms of arithmetic operations and CPU time be reduced. The proposed method comes from Newton’s method and a three-step method in [16,17], leading to produce the proposed approach with convergence order while using just seven function evaluations per iteration. The four-step proposed method results in the following computational steps, whose flowchart is depicted in Figure 1:

for .

The methods as mentioned earlier, including the one proposed as the four-step hybrid method in the present article, are compared in Figure 2 with each other based on efficiency index , order of convergence , and the number of function evaluations used by each method per iteration.

#### 3. Convergence Analysis of the Proposed Method

This section has been divided into two sections wherein the convergence of the proposed four-step hybrid method in both scalar and vector form has been discussed in detail. It is noted that the twelfth-order convergence is theoretically verified in each case with the aid of Taylor’s expansion.

##### 3.1. Convergence Theory with Scalar Form

In this subsection, we theoretically prove the local order of convergence for the proposed method given in (9).

Theorem 1. *Suppose that be the required simple root for a sufficiently differentiable function within an open real interval . Then, the proposed four-step numerical method (9) possesses twelfth-order convergence with the error equation given by:where and *

*Proof. *Suppose be a simple root of , where be the th approximation for the root by the proposed method (9), and be the error term in variable at the th iteration step. Employing the single real variable Taylor’s series in [9] for around , we obtainSimilarly, using the Taylor’s series for around , we obtainMultiplying (11) and (12), we obtainSubstituting (13) in the first step of (9), we obtainwhere . Using the Taylor’s series for around , we obtainSimilarly, using the Taylor’s series for around , we obtainMultiplying (15) and (16), we obtainSubstituting (17) in the second step of (9), we obtainwhere . Using the Taylor’s series for around , we obtainSimilarly, using the Taylor’s series for around , we obtainMultiplying (19) and (20), we obtainSubstituting (21) in the first step of (9), we obtainwhere . Using the Taylor’s series for around , we obtainSubstituting (19), (20), and (23) in the fourth step of (9), one obtainsFinally, using , and (14) and (18), (22) for the above equation, we obtainHence, the twelfth-order convergence of the proposed method P12 given by (9) for the nonlinear functions in single variable is proved.

##### 3.2. Convergence Theory with Vector Form

.

*Proof. *Suppose be a simple root of , where be the th approximation for the root by the proposed method (9), and be the error term in variable at the th iteration step. Employing the multi variable Taylor’s series given in the theorem [9] for around , we obtainAgain, using the Taylor’s expansion for the inverted Jacobian matrix around , we obtainMultiplying (27) and (28), we obtainSubstituting (29) in the first step of (9), we obtainwhere . Using the Taylor’s series for around , we obtainAgain, using the Taylor’s expansion for the inverted Jacobian matrix around , we obtainMultiplying (31) and (32), we obtainSubstituting (33) in the second step of (9), we obtainwhere . Using the Taylor’s expansion for around , we obtainAgain, using the Taylor’s expansion for the inverted Jacobian matrix around , we obtainMultiplying (35) and (36), we obtainSubstituting (37) in the first step of (9), we obtainwhere . Using the Taylor’s series for around , we obtainSubstituting (35), (36), and (39) in the fourth step of (9), one obtainsFinally, using , and (30) and (34), (38) for the above equation, we obtainHence, the twelfth-order convergence of the proposed multi-step (four-step) hybrid method P12 mentioned by (9) for the nonlinear functions in multi-variable case is proved. □

#### 4. Polynomiography

Polynomiography is a process that integrates mathematics and art to create a new type of visual art. The produced graphics result from algorithmic visualization of iterative approaches for solving a polynomial equation. This term was first introduced by Dr. Bahman Kalantari at the start of the 21st century [18]. Dr. Bahman Kalantari’s study on polynomial root-finding, which is an old and traditional discipline that continues to find new implications with each generation of mathematicians and scientists, inspired the concepts of Polynomiography. Dr. Kalantari invented the term ”polynomiography,” which is a mixture of the word ”polynomial” with the suffix ”graphy.” A ”polynomiograph” is a separately produced picture resulting from Polynomiography. It is defined as ”An iterative procedure for producing two-dimensional colored pictures (polynomiographs) that represent polynomials.”

In recent years, researchers worked in the field of Polynomiography along with its implementations in other fields. In [19], the authors introduced a new mathematical art with the help of Newton–Ellipsoid method. Gdawiec et al. in [20], presented the visual analysis of Newton’s method with fractional-order derivatives. The authors employed the techniques of coloring by roots and coloring via iterations to study the convergence and dynamical aspects of the processes visualized by polynomiographs.

Naseem et al. in [21] presented some new graphical art with the help of newly suggested ninth-order iteration schemes. Scot et al. in [22], presented the basin of attraction for several methods and examined its dependence on their convergence orders. In [23], the authors introduced a new family of eighth-order methods and then drew their basins of attractions by assuming different polynomials. Finally, in [24], the authors generated some new fractal patterns by combing two root-finding methods. The obtained fractal patterns were diverse and had many applications in the textile and ceramic industries.

We use a rectangle along with the dimension , accuracy and the max. no. of iterations to create the polynomiographs over the complex plane C through the computer software by taking multiple complex polynomials. The color black is allocated to the spots where the method failed to converge. The partitioning of determines the pixel density of the created visual representations; for example, if we partition the rectangle into a grid of , the plotting polynomiographs will then have better resolution.

For drawing graphical objects in the complex plane, we use the four complex polynomials listed below:

For coloring the iterations, we employ the colormap given in Figure 3.

*Problem 1. *Polynomiographs for the Polynomial Through Various Methods

In this example, we investigate and compare the dynamical results obtained through different iteration schemes with our presented method by considering the cubic polynomial which possesses three distinct simple zeros:1, , . We executed all the methods to achieve the simple zeros of the considered polynomials and the results *can be visualized in* Figure 4.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(g)**

**(h)**

*Problem 2. *Polynomiographs for the Polynomial Through Various Methods

In this example, we investigate and compare the dynamical results obtained through different iteration schemes with our presented method by considering the 3rd-degree polynomial which possesses three distinct simple zero*s:*We executed all methods to achieve the simple zeros of the considered polynomials and the results can be seen in Figure 5.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(g)**

**(h)**

*Problem 3. *Polynomiographs for the Polynomial Through Various Methods

In this example, we investigate and compare the dynamical results obtained through different methods with our presented method by considering the 4th-degree polynomial which possesses the following simple zeros:1, , , and . We executed all the methods to achieve the simple zeros of the considered polynomials and the results are given in Figure 6.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(g)**

**(h)**

*Problem 4. *Polynomiographs for the Polynomial Through Various Methods

In this example, we investigate and compare the dynamical results obtained through different iteration methods with our presented method by considering 4th-degree polynomial which possesses the simple zeros:1, , 2, and . We executed all methods to achieve the simple zeros of the considered polynomials and the results can be visualized in Figure 7.

In all given examples, a detailed graphical analysis of the designed algorithm has been provided via polynomigraphs. For plotting polynomiographs on the complex plane, we take two cubic-degree polynomials namely: , and and two quatric-degree polynomials represented by , and respectively. The plotted graphs tell us about the convergence speed and the iterations performed by the method for drawing these objects. The second characteristic is the dynamics of the iteration scheme. In each polynomiographs, the individual root has been denoted by the blue colored dot. The black colored zones denote the divergence area or deficiency of the method through which the polynomiographs has been plotted. The darker or brighter zones in the provided polynomiographs showing less iterations performed to approximate the solution. One can easily observe the superiority of the proposed method over the others by examining the more darker or brighter zones of the polynomiographs drawn by the suggested method.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(g)**

**(h)**

#### 5. Numerical Simulations: Real-world Scenarios

This section of the paper discusses the real-life applications by applying the newly proposed hybrid method. We also present a numerical comparison with other existing most frequently used methods: N2, N4, W8, N9, P10, N11, and N12, whose computational steps are shown in the introductory section above. In each applied model, we set the tolerance to be as the stopping criterion of the iterative process of every method under consideration: . Two additional methods with fifth- and sixth-order of convergence are also included for the simulations of a six-dimensional model chosen from the field of neurophysiology based on the reason that some of the methods under discussion in the above introductory section did not prove to be valuable candidates when it comes to the simulations of the nonlinear models presented in the system or vector form.

*Problem 5. *The Plank’s radiation law in physics explains the spectral density of radiation emitted by a black body in thermal equilibrium at the temperature and the condition that there must not be a flow of energy between the body and its surroundings. In other words, the law is introduced to determine the amount of energy density in a black body based on isothermal properties. Moreover, it is sometimes used to estimate the maximal radiation’s wavelength. As described in [25], the maximal wavelength of the radiation may be written in the form of the following nonlinear equation in scalar version:where stands for the maximal wavelength. The exact solution of the above equation is as follows: 0.0.

The Plank’s radiation model described in (44) is simulated with several numerical algorithms while assuming two different initial guesses. The maximum number of iterations in each case is set to be . It can be observed in Table 1 that the accuracy is maximum for P12 at the cost of a slightly higher amount of CPU time, regardless of the initial estimate’s location. The eleventh-order method N11 did not converge with the initial estimate taken to be while it converges second to P12 when the initial guess is taken. With this second initial guess, the CPU time consumption also slightly increases for N4 and N11.

*Problem 6. *Fraction Conversion of Nitrogen-Hydrogen to Ammonia [26]. This nonlinear scalar problem depicts the fractional conversion of nitrogen-hydrogen to ammonia and has appeared in several research works conducted in the past and recent literature. In this experiment, we set the pressure value to be 250 atmospheric pressure while the temperature is set to 500 degrees Celsius. In terms of a nonlinear function, the model mentioned above has the following polynomial form:It has been identified in the recent work [27] that one of its positive real roots lies in the open unit interval which is estimated to be 0.2777595428.

For this nonlinear model, the numerical simulations are shown in Table 2 while the number of iterations for each method under consideration is set to be 7. Two different initial guesses, that is, (near to the exact solution) and (away from the exact solution), are chosen. It is seen that the fourth-order method N4 converges towards some other solution for the first initial guess while the method abbreviated as N11 failed after three iterations while the same is the case for the W8 method, but the method W8 produced the correct approximate solution till four iterations and failed after that. The Newton method N2 has the comparatively most significant absolute error at the fourth iteration compared to other methods. Nonetheless, the proposed hybrid method, in addition to a few more methods, always converged towards the required solution. More so, the proposed method has achieved the minor error tolerance with a reasonable amount of time. Hence, it can be concluded that the initial location of the estimate does not matter much when it comes to the proposed hybrid method devised in this article.

*Problem 7. *Application of mechanical engineering [28]: There are several fields wherein the use of thermodynamics is extensively required. The particular areas include mechanical, civil, mechatronics, electronic, chemical engineering, and many others. The fourth-order polynomial is used to show a relation between the zero-pressure specific heat of dry air (kJ/kgK) to temperature :where is used.

As described in [28], the above nonlinear model given in terms of fourth-order polynomial has two real distinct roots given as: and . It can be observed in Table 3 that each method converges for the initial guesses chosen to determine the approximate solution of the above model. The eleventh-order method N11 performed better from an accuracy viewpoint, followed by the proposed hybrid method. Looking at the CPU values, it is clear that the methods N11 and P12 take the same amount of time, thereby being equally time-efficient for this particular fourth-order polynomial. Moreover, the absolute functional values are identical for both methods, including others such as N9 and P10.

*Problem 8. *Neurophysiology Application: As a final experiment, we consider a six-dimensional nonlinear system first proposed in [29] and later was used by several researchers for the simulation purpose of their newly developed algorithms. See, for example, [30], and some cited references therein. The nonlinear model consists of the following six equations:The constants in the above model can be randomly chosen. In our experiment, we co*nsidered**.*

The simulations for the above neurophysiology application system are shown in Table 4 wherein the two columns represent the absolute error at the last iteration and the CPU time consumption. Two more methods, including fifth-order Halley’s (HM5) in [6] and a sixth-order Hameer-Mujtaba method (HM6) in [17], is used for the simulations of the above system. It is evident from Table 4 that the accuracy is much higher for the proposed approach in comparison to other competitive methods, while the CPU time is also reasonable.

#### 6. Concluding Remarks and Future Directions

A new four-step nonlinear method for solving type models is introduced in this research work with twelfth-order convergence, and seven function evaluations are required per iteration. The theoretical order of convergence for the proposed hybrid method is proved under both scalar and vector cases, along with an asymptotic error constant. Comparison with various existing numerical methods discloses the better performance of the proposed approach when the absolute errors, the absolute functional value computed at the last iteration, and the time of machine in seconds are taken into consideration. The proposed method brings out the slightest absolute error regardless of initial conditions chosen for the simulations of the nonlinear models that belong to real-world scenarios from science and engineering. The rapid convergence of the proposed hybrid method is confirmed with the aid of polynomiography when the method is applied to some complex-valued polynomials.

Declarations

#### Data Availability

All the data required for this paper is included within the paper.

#### Ethical Approval

We with this affirm that the contents of this article are original. Furthermore, it has been neither published elsewhere fully or partially in any language nor submitted for publication (fully or partially) elsewhere simultaneously. It contains no matter that is scandalous, obscene, fraud, plagiarism, libelous, or otherwise contrary to law.

#### Conflicts of Interest

The authors do not have any conflicts of interest.

#### Authors’ Contributions

All authors contributed equally in this paper.