Two New Iterated Maps for Numerical Nth Root Evaluation

In this paper we propose two original iterated maps to numerically approximate the nth root of a real number. Comparisons between the new maps and the famous Newton-Raphson method are carried out, including fixed point determination, stability analysis and measure of the mean convergence time, which is confirmed by our analytical convergence time model. Stability of solutions is confirmed by measuring the Lyapunov exponent over the parameter space of each map. A generalization of the second map is proposed, giving rise to a family of new maps to address the same problem. This work is developed within the language of discrete dynamical systems.


Introduction
Recent applications of iterated maps in numerical analysis have been found in literature, using and extending the techniques of dynamical systems to the study of numerical algorithms and number theory [1]- [3].Application in technology and hardware devices are also frequent nowadays [4]- [7].
We propose and study in this work two new methods for numerical root approximations, both of which based on iterated maps.In the following sections we present a detailed study of each map, their fixed points and stability, the occurrence of bifurcations and chaotic behavior.
Some common tools of nonlinear dynamics [8] [9] are used to the study of the orbits, i.e., the numerical time series obtained for each map are investigated.The numerical approximations , 0,1, 2, i x i =  to solve the gen- eral equation n x k = , where n ∈  and k ∈ + are obtained, but the validity of the results can be extended to n ∈ + .
In Section 2 we present a new map proposed by one of us (C.C. Dias), named as First Dias Map (FDM), showing the existence of a fixed point for roots in the range [ ] 1, 4 n ∈ , and that this fixed point corresponds exactly to n k .
In Section 3, we generalize the FDM by adding a new parameter for studying the stability of its fixed point by defining a new class of maps called Weighted Average Map (WAM).For this class of maps, we investigate the dependence of the fixed point corresponding to the nth root of k over the parameter space ( ) , n p .Finally, in Section 4, we measure and compare the Mean Convergence Time (MCT) for all the studied maps, for 2 n = and 3 n = , varying k over an uniform grid of initial conditions and computing the average amount of iterations to converge to the root n k within the standard numerical double precision.An analytical model is proposed and used to confirm the numerical results of MCT with the analytical convergence time (ACT) for the WAM.

The First Dias Map (FDM)
The map which we will study now was created by Charles C. Dias to extract real roots of numbers numbers, by solving the equation n x k = .The proposed map is one-dimensional and is defined as Comparing this with the Newton-Raphson Method (NRM) equation and Babylonian Method (BABM) [10] noticed that this statement is a mixture of both, and the FDM is an arithmetic average between the linear and nonlinear terms in Equation (1.1), and can be used to approximate the nth root of 0 k > for ( ) , as we shall see.Outside this interval of n , this map presents chaotic dynamics through after entering a bifurcation cascade, whose roots having no longer relationship to the nth root of k .
The base function that appears in the iterated map defined by Equation (1.1) can be derived dividing n and adding x at both sides, and dividing it by 2, we recover functional form of Equation (1.1).

Geometrical Construction
To construct geometrically the FDM time series, the first step is to find the auxiliary equations of the lines i A B (see Figure 1), writing their slopes i i i m y x = ∆ ∆ .From this figure, ( ) + and i i x x ∆ =, and the slopes ( ) , 0,1, 2,3, Knowing that their linear coefficients are all k − , then all the auxiliary lines pass through the point ( ) we obtain the working lines ( ) − and the auxiliary points i s , the intersections points of the work- ing lines with the X -axis, are 1 , 0,1, 2,3, , and taking the arithmetic mean between the auxiliary points i s and the i x points we recover the original FDM equation (Equation (1.1)).
Figure 2(a) shows the cobweb for the FDM time series for 2 n = and 2 k = , and Table 1 (top) shows time series used in this figure.In this example, the convergence to the root is achieved after only five steps, considering the standard double precision, and is exactly the same time series of NRM for these parameters.shows a numerical development of the FDM series for the parameters 3 n = and 2 k = , based on the time se- ries shown in Table 1 (bottom), where the convergence to the root occurs after 27 steps.Some intermediary time steps are omitted in this table.

Fixed Point and Stability Analysis
Solving ( ) x * * = we find the FDM fixed point to be . Applying the stability criterion [11], to the map function ( ) ( ) and solving the last equation we have the range of parameter where the fixed point of the map is stable.

Numerical Results
The FDM time series have different dynamics depending on the parameters ( ) , k n , presenting a fixed point, periodicity or chaos, as occurs to the logistic map [12].
To measure the rate of divergent orbits, i.e., the sensitive dependence on initial conditions, we can use is cha- . From Figure 3(a) we see that, at 4 n = , FDM enters a bifurcation cascade, therefore its fixed point is no longer stable.In the white to gray regions the exponent is negative indicating that for this region of parameters the FDM not is chaotic, and at the black stripes the Lyapunov exponent  3(b).The values of n studied are uniformly distributed in a grid of 600 points in the interval [ ] 1,10 .Also plotted over the bifurcation diagram is the exact root n k , the fixed point of the map, plotted in black.
We also study numerically the FDM return diagrams for different values of the parameter n , for 2 k = and

The Weighted Average Map (WAM)
Instead of adding x , if a more general term px is added in Equation (1.2), we get a new map (WAM) that depends on parameters n , k and p .The new parameter p is a positive real number and corresponds to the weight of the linear term of the map.This term is directly linked to the parameter n , since for each value of n there is a minimum value of p for the fixed point to be stable, as we shall see.
Adding px to both sides of Equation (1.2) we gain and after collecting x and dividing by ( ) and solving its fixed point equation we obtain the expected value n x k * = .

Fixed Point and Stability Analysis
Applying the stability criterion [11], i.e., ( ) ) and solving this inequality we obtain 2 1 p n > − to guarantee fixed point stability, and Figure 4 shows the line corresponding to this condition, below which the fixed point is unstable.As n is increased the value of p should also be increased to avoid the unstable region, where the time series do not converges to the fixed point.

WAM Subclasses and Hierarchy
A special subclass of WAM is FDM, when 1 p = , so that the fixed point on the map according to Figure 4, is stable in the range [1,4] for n , thus in accordance with the stability analysis the fixed point n k loses stability    4 n = .Other very important subclass is NRM, for ( ) − , so that the weight p is chosen within the stable region.
For NRM, the derivative of the mapping function at the fixed point is null, satisfying the stability criterion and resulting in the most efficient rate of convergence of the time series near the fixed point.When the starting point 0 x is chosen far away from the fixed point, we observe numerically that the initial rate of convergence is greater for FDM, in general.Finally, there is a third subclass of WAM, the Babylonian square root method (BABM), for 1 p = and 2 n = , the oldest and perhaps the most efficient known method to solve the equation 2 x k = .At this parameters, WAM reduces itself to BABM, therefore the same occurs with FDM and NRM.The hierarchy of WAM subclasses are shown in Figure 5.
From the definition of the Lyapunov characteristic exponent for a unidimensional map we conclude that the derivative of the mapping function at the fixed point ( ) ′ defines the rate of convergence of its time series, discarded the transient.For WAM, it is easy to show that ( ) ( ) ( )

Mean Convergence Time (MCT)
This section reports the numerical results for the mean convergence time (MCT) for NRM, FDM and WAM, based on the average number of iterates to converge within different precisions ε , from single ( ) , varying the initial condition 0 x on a second uniform grid with 3 10 points, whose limits are given by a maximum relative difference of 25% around the exact value of the root of k at each point.Using this schema, the MCT is computed for cubic roots ( ) 3 n = , and for WAM we set 3 p = .Figures 6(a)-(c) show the numerical results.In Figure 6(a) we see that the NRM MCT is close to 4, which means that after 4 iterations, on average, there has been convergence to the root.From this figure, we conclude that FDM is around 10 times slower than NRM, and WAM is around has twice the speed of FDM.In this test, the most efficient is NRM, with the lowest MCT.
Both NRM and FDM belong to the same WAM family, as discussed in Section 3, and the stability of the fixed point n k of WAM depends on the parameters n and p .Changing the parameter p of WAM we get a new map subclass, for example, for 2 n = have the FDM.From these fact, we tried to detect numerically the optimal value of the parameter p , to minimize the ACT over the whole WAM family.For this we used a FORTRAN program to measure extensively the WAM MCT varying parameters ( ) , n p on a uniform grid of 500 500 × points, for a radicand k .The result for 2 k = is shown in Figure 4(a).The region in gray corresponds to unstable fixed point n k map WAM, as found in Section 1.3.The other colors seen in the graph are the regions of stability of the fixed point.For best visualization the MCT scale of this figure is truncated at a maximum value of 20, and higher values as inked light gray.
We can see in Figure that, as we approach the line that corresponds to the weighting term, NRM the minimum MCT over this line and therefore the most efficient of all studied maps is the NRM.
Summarizing the key information about NRM, FDM and WAM, with the numerical results for the MCT within double precision for these maps, for the parameters 3 n = and  2.

Analytical Convergence Time Model
The Lyapunov characteristic exponent for a unidimensional map ( ) usually defined by ( ) ( ) ( ) since the derivative of the mapping function at the fixed point ( ) f x * ′ defines the rate of convergence of its time series, after discarded the transient.
For WAM, it is easy to show that ( ) ( ) + , so that this derivative assumes the value 1 2 n − when 1 p = (FDM) and is zero only if 1 p n = − (NRM or BABM, if 2 n = ).In Figure 3(c) we observe that ( ) 0 f x * ′ = , for 2 n = FDM reduces to NRM, and in Figure 3(d), ( ) Using the original Lyapunov's idea, the characteristic exponent λ measures the average rate of convergence between two solutions separated by an initial distance 0 ε , that is the case of time series dominated by a fixed point.For this orbits, the distance after i iterates is 0 e i i λ ε ε = so that, if we assume that one orbit is initialized at 0 x x * = and other at 0 1 x x * = + , i.e., the initial distance is unitary between orbits, we can use last equation to measure the error found in the second orbit, the root to be approximated.Within the standard double precision, the maximum error is of order 10 D − , where D is the number of decimal significants, typically 15 16 D = ∼ places.
Applying the natural logarithm to both sides of the above equation we have ( ) for the number of iterations needed to reduces the error in the second orbit to ε .In the same manner we defined MCT, we define now the analytical convergence time (ACT), estimated by valid for any fixed point of a unidimensional map, where the approximated λ * was used.
Applying this model to our more general map (WAM), we have

Conclusions
In the study of iterated maps to extract the real root of real numbers we have applied some common tools from nonlinear dynamics that allowed us to predict the fixed point of the studied maps associated with the nth root of k , and their stabilities could be analyzed in details.
We conclude, through the geometric argument used to recover the original analytical form of FDM, that both NRM and FDM can be reduced to averages between terms, one linear and other nonlinear.From this observation, we generalize the original FDM idea to a new family of maps on which we add a new parameter p , whose value defines the stability of the map, the WAM.We show that FDM and NRM belong to this family of maps, being FDM recover when 1 p = and NRM when 1 p n = − .The mean convergence time (MCT) numerical results indicate that NRM is the most efficient subclass of the more general weighted average map (WAM) proposed in this work, as pointed out in Figure 4(a), over the line 1 p n = − .The analytical model for ACT is in complete agreement with the numerical results for WAM, the most general class of map studied.The model presented in Equation 1.7 is general, and can be adapted to any unidimensional map to study its fixed point attractor.
The main results of this work are obtained for x ∈  , but their generalization is straightforward over the complex set  .

Figure 4 .
Figure 4.For 2 k = , (a) the numerical results of WAM MCT (n, p) and (b) the analytical model for WAM ACT (n, p). at

Figure 3
FDM reduces to NRM, and in Figure3(d), we observe that

10 −
. For this, we varied k on a uniform grid with 10 3 points in the interval ]
shown in Table

1 D
n and p .To double precision this approximated model function is plotted in Figure4(b), that is remarkably very close to the numerical version plotted in Figure4(a).Both figures uses the same color palette and truncated maximum, for better comparisons.

Table 1 .
FDM signing the period bifurcations.The yellow to red regions indicate the a positive Lyapunov exponents, the signature of chaos.The FDM bifurcation diagram, discarded a transient of 10 3 iterations and plotted the next 500 values of x , is depicted in Figure

Table 2 .
MCT numerical results for NRM, FDM and WAM maps.