Ascertaining when a basin is Wada: the merging method

Trying to imagine three regions separated by a unique boundary seems a difficult task. However, this is exactly what happens in many dynamical systems showing Wada basins. Here, we present a new perspective on the Wada property: A Wada boundary is the only one that remains unaltered under the action of merging the basins. This observation allows to develop a new method to test the Wada property, which is much faster than the previous ones. Furthermore, another major advantage of the merging method is that a detailed knowledge of the dynamical system is not required.

Trying to imagine three regions separated by a unique boundary seems a difficult task. However, this is exactly what happens in many dynamical systems showing Wada basins. Here, we present a new perspective on the Wada property: A Wada boundary is the only one that remains unaltered under the action of merging the basins. This observation allows to develop a new method to test the Wada property, which is much faster than the previous ones. Furthermore, another major advantage of the merging method is that a detailed knowledge of the dynamical system is not required.
Wada basins are one of those unexpected encounters that often happen in science. The story begins when a topologist named Takeo Wada tried to answer the following question: Can three or more open regions be separated by a single boundary? Our daily experience makes us think that this is impossible. It suffices to look at a common political map to realize that the boundaries separate two different regions, except perhaps some isolated points that separate three or more regions (think about the Four Corners in the USA, for example). However, Takeo Wada devised an iterative process to make this counter-intuitive situation possible, as reported by his student Kunizo Yoneyama 1,2 . The Wada lakes were conceived in a topological context as a way to separate three connected regions in a plane by means of a continuous boundary 1 . From a topological point of view, Wada boundaries have intriguing properties. For example, the Polish topologist Kazimierz Kuratowski showed that in the plane, continuous Wada boundaries must be indecomposable continua 3 (though the situation in higher dimensions is quite different).
This discovery remained as a mathematical curiosity until James Yorke and his collaborators found that the basins of attraction of some dynamical systems presented the Wada property 4,5 . From the dynamical point of view, the most interesting feature of Wada basins is the fact that an arbitrarily small perturbation of a system with initial conditions lying in a Wada boundary can drive it to any of the possible attractors, which implies a special kind of unpredictability 6 . Therefore, in this context, Wada boundaries are usually referred to as those that separate three or more basins at a time, but the basins need not to be connected. Since the earliest references to the Wada property in dynamical systems, many authors claim that the boundaries have the Wada property for disconnected basins [7][8][9][10][11][12] . In this work, we adopt this latter definition: Wada boundaries are those that separate three or more basins, no matter whether the basins are connected or not. Therefore, using this definition, we believe that the methodology and the results presented in this work are valid for any number of dimensions.
Despite our primary intuition, Wada basins are a common feature appearing in many dynamical systems. Since its first report, Wada basins have been found in open Hamiltonian systems 10 , ecological models 11 , delayed differential equations 12 , hydrodynamical systems 13 , and many engineering problems [14][15][16] . This is possible because Wada boundaries are related to iterative processes and fractal structures, which are a common feature in the basins of nonlinear dynamical systems 17 .
So far, two methods have been proposed to determine when the basins of a system possess the Wada property. The first one was developed by Nusse and Yorke 18,19 , and involved the computation of the unstable manifold of a saddle point of the basin boundary. This method requires a detailed knowledge of the system and the computation of unstable trajectories, which can be cumbersome in many cases. Indeed, many papers 10,13,14,20 are devoted just to determine whether the Nusse-Yorke condition is fulfilled in a particular dynamical system and for a certain set of parameters. Years after the original works by Yorke and collaborators, Daza et al. 21 developed a grid method based on the successive refinement of the grid in order to determine whether all the boundary points were Wada points (points that separate three or more basins at a time). This latter method can be automated and used for every dynamical system. As a drawback, it requires the precise computation of new trajectories at very high SCIeNTIFIC REPORTS | (2018) 8:9954 | DOI:10.1038/s41598-018-28119-0 resolutions. Although it supposes a qualitative and quantitative improvement with respect to the Nusse-Yorke method, the grid method needs several hours or even days of parallel computation to check the Wada property in a given dynamical system. In this paper, we present a new method to determine when a basin is Wada, which is founded on the observation that: A Wada boundary is the only one that remains unaltered under the action of merging the basins. This new method, that we call the merging method, can test the Wada property in a few seconds, and furthermore it does not require a detailed knowledge of the system. Thus, the merging method supposes a new quantitative and qualitative leap with respect to the previous available methods to test the Wada property.
The paper is organized as follows. In Sec. 2, we explain how Wada boundaries can be defined as the only ones that remain unaltered after the action of merging the basins. In Sec. 3.1, based on the previous definition, a new method to test the Wada property is presented. Section 3.2 is devoted to the detailed analysis of the merging method using different model examples. Finally, we discuss the results and present the main conclusions of the paper in Sec. 4.

Merging Basins
The set of all initial conditions leading to a particular attractor is called the basin of attraction of a dynamical system. We will focus on a very special set of initial conditions called the boundary. A point p is in the boundary If the point satisfies the previous condition for all the basins B i with N a ≥ 3 basins of attraction, we call it a Wada point. If all the boundary points are Wada points, then the basin of attraction has the Wada property, and we call it a Wada basin.
However, we can formulate the Wada property in slightly different terms. Assume we have N a ≥ 3 basins of attraction. Now, we want to determine the boundary ∂B i of each basin B i , but instead of using its complement B i  to determine which points belong to the boundary, we will say that a point p is in the boundary if it is arbitrarily close to B i and also arbitrarily close to at least one of the other basins Then, we determine each boundary ∂B i as the boundary between a basin B i and all the other merged basins ∪ ≠ B j i j , so that we end with as many different boundaries as different possible attractors, i = 1, …, N a . Thus, all the boundaries created following this previous procedure and the boundary of the original basins corresponding to the N a attractors are exactly the same ∂B i = ∂B j for ∀i ≠ j, i = 1,…, N a , if and only if the system is Wada.
The two previous definitions of Wada basins are completely equivalent. However, the second definition emphasizes the striking idea that Wada basins can be merged and the boundary will still remain unaltered. To be more precise, given N a ≥ 3 Wada basins, it is possible to merge up to N a − 1 without any change in the boundary (if we merge the N a basins then there would be only one basin and the boundary would be lost). This notable effect is better appreciated when Wada boundaries are compared to non-Wada boundaries. The time-2π (Poincaré) map of the forced damped pendulum defined by + .
sin 166 cos possesses three attractors, and consequently its phase space  x x ( , ) contains three basins. This is a paradigmatic system showing Wada basins 4 . In the top-left panel of Fig. 1(a), we display the original three-colored Wada basins of the forced damped pendulum described in 4 . The other three plots show the result of merging the basins according to the color code sketched in Fig. 1(c) (yellow = red + green, magenta = blue + red, cyan = blue + green). It is important to notice that each color represents a different basin, being impossible to establish a one-to-one correspondence between basins of different colors. However, even though the basins are different, the boundaries are the same for the four panels of Fig. 1(a).
If we look at the colored disks in Fig. 1(b), we can see how the action of merging affects usual (non-Wada) basins. Here we can clearly notice that the boundaries change under the action of merging the basins. In fact, the center of the disk is the only point that is in the boundary of the four panels, so that it is a Wada point.
Despite the abundant research devoted to Wada basins, the effect that the Wada boundaries remain unaltered after the action of merging the basins has been unnoticed. In the next sections, we will use it to develop a new way of testing Wada basins in dynamical systems.

Merging Basins to Test the Wada Property
Description of the merging method. The property that we have just described, that is, that Wada basins can be merged without any change in their boundary, can be used to build a new method to test the Wada property. From a purely mathematical point of view, given the basins of a system, it suffices to check that the boundary remains unaltered under the merging of the basins. However, it would require an arbitrarily high resolution of the basins to guarantee that the boundaries of the merged basins are exactly the same.
Usually, the basins are computed by means of a regular grid of finite size. In this approach, every pixel of the grid has a linear size ε and contains only one corresponding initial condition, in such a way that the fate of this initial condition determines the color of the pixel. Therefore, the computation of the boundaries is limited by the size of the pixel ε. In Fig. 2(a,b), we can see that the computed boundaries of the merged basins, which we call slim boundaries, are not exactly the same, even though they are Wada boundaries. It can be observed at naked eye that although the boundaries seem similar, they are not strictly identical. It is noticeable that the boundary depicted in Fig. 2(a) is thicker than the boundary depicted in Fig. 2(b).
To overcome this issue, we can try to fatten the boundaries for their subsequent comparison (see the fattened boundary of Fig. 2(c)). In the fattening procedure, we replace each pixel belonging to the boundary by a fat pixel defined by the fattening parameter r. This fattening parameter is the radius of the fat pixel according to the Chebyshev metric or maximum distance metric. This metric preserves the square shape of the pixels and it is defined in the plane as r = max(|x 2 − x 1 |, |y 2 − y 1 |), where (x, y) are the usual Cartesian coordinates. Sometimes this metric is also called the chessboard distance, since it represents the number of moves that a king would have to make to go from one position to another (see Fig. 2(d)). The way the fattening is made can be changed according to different metrics, this is not crucial for the method. In the next section, we will analyze how the method that we are describing depends on the fattening parameter r. Now, let us move forward to the last part of the procedure.
For the moment, we have the original boundaries of the merged basins, the slim boundaries ∂B i , and their fattened versions, the fat boundaries ∂B i . The final step of the procedure is to compare all the slim boundaries with all the fat boundaries. If all the slim boundaries fit in all the fat boundaries ∂B i ⊂ ∂B j ∀i,j = 1, …, N a ; then we will say that the basin is Wada. Otherwise, we will say that the system is not Wada, and the method will determine which points are Wada and which ones are not. This last step verifies if each pixel of the slim boundaries ∂B i lies in the set ∂B j . To connect with our previous definition of a basin with the Wada property, the algorithm checks if the points p i in the boundaries B i are within a ball b(p j , r) of radius r (r is the fattening parameter) around the points p j of the boundary B j .
In the case of partially Wada basins 16 , where Wada and non-Wada boundaries coexist, we can characterize them by the Wada parameter W N a defined in the grid method of Daza et al. 21 . This parameter W N a provides the ratio of Wada points to boundary points (Wada and non-Wada), in such a way that = W 1 N a means that the system has the full Wada property, whereas indicates only partially Wada basins. In the merging method, given a basin, we can compute the pixels lying in the boundary of that basin n b , and we can also register the number of boundary points which are not Wada n NW . Then, the Wada parameter for a fixed resolution can be calculated simply as Again, for a better understanding of the comparison between slim and fat boundaries, it is convenient to observe an example of non-Wada basins, such as the disks of Fig. 1(b). We would have to fatten the boundaries by a very large amount (comparable to the size of the disks) in order to make the fat boundaries able to contain the slim ones. We can conclude, as mentioned before, that the only Wada point is the center of the disk.
The whole procedure described before can be fully automated and the only input needed is a finite resolution basin. For basins with a resolution of 1000 × 1000 and three different colors (attractors), the merging method takes around two seconds to determine whether a basin is Wada running in a personal computer. This contrasts with previous methods to test Wada basins. The grid method 21 needs to compute new trajectories at finer resolutions, which can take several hours or even days of parallel computation in a cluster with one hundred cores. The Nusse-Yorke method 18,19 requires detailed knowledge of the dynamics of the system and, in general, it cannot be automated. In fact, many works 10,13,14,20 are exclusively devoted to the application of the Nusse-Yorke method to one particular system and one specific set of parameters due to the difficulty of the task. In comparison, the merging method is incredibly fast and easy, since it only requires a finite resolution basin to be applied. Furthermore, since the merging method does not need any further computation of new points, we do not even need to know the underlying dynamical system nor its parameters. Next we summarize the steps that this merging algorithm takes, which can also be visualized in the flowchart of Fig. 3. (a) At first, we have a picture of the basins at a given resolution ε. As we will discuss later, the higher the resolution the more reliable the determination of the Wada property will be.   Analysis of the merging method. The whole method described above relies on a single parameter: the fattening parameter r. This parameter determines the confidence that we have in the result of the algorithm, since we will be able to say that the basin is Wada up to the resolution defined by the fat pixels that we use. Therefore, it is natural to analyze the behavior of the method for different values of r in dynamical systems with different features. This is exactly our aim in this section. In order to examine the behavior of the procedure with respect to the fattening parameter r, we can apply it to different Wada boundaries. We can characterize fractal boundaries by their fractal dimension and by the number of basins that they separate at the same time. Here, we examine two dynamical systems with Wada boundaries where we can easily vary these two quantities. Namely, the two paradigmatic dynamical systems under study are the Hénon-Heiles Hamiltonian 10 and the Newton method to find complex roots 8  . For values of the energy above the critical one, the escape basins of this open Hamiltonian system show the Wada property 10 . The fractal dimension of the boundaries diminishes as the energy increases, but the Wada property is preserved, as reported in 22,23 . We have used the merging method described in the previous section with different values of the energy E and of fattening parameter r. The results are plotted in Fig. 4(a). We have plotted only three different values of the energy for clarity, but we can observe that the algorithm correctly determines that the basins are Wada for r ≥ 4 at every tested value of the energy. It can also be noticed that there is no relation between the number of non-Wada points detected by the algorithm and the fractal dimension of the boundaries. Furthermore, we have carried out similar computations for increasing resolutions with analogous results. Thus, from these numerical experiments, we conclude that there is no relation between the value of the fattening parameter r needed to correctly predict the Wada property and the fractal dimension of the boundaries.
The second system where we have tested the merging method described before is the Newton method to find complex roots. This method can be described by the discrete complex variable map z n+1 = z n − (z R − 1)/(R⋅z R−1 ), where the parameter R determines the number of roots and therefore the number of attractors N a = R. It has been reported that the basins produced by this complex variable map show the disconnected Wada property (the basins are disconnected and also Wada) no matter the number of attractors determined by R 8,24,25 . Thus, we ran the merging algorithm for an increasing number of roots R, and consequently of attractors N a . As shown in Fig. 4(b), the merging method correctly classifies the basins as Wada for all r > 4 even for a large number of attractors N a , which seldom appears in typical dynamical systems. Moreover, we have found no trivial relation between the number of attractors N a and the percentage of non-Wada points. Again, we have performed the computations at different resolutions (up to 5000 × 5000) with consistent results. Hence, we conclude that the value of the fattening parameter r needed for a correct classification of the basins does not directly depend on the fractal dimension of the boundaries nor on the number of attractors. This also proves that the method works correctly for disconnected Wada basins.
Finally, we would like to mention another adjustment that could be added to the merging algorithm in case of need. Just as described before, the merging algorithm is an all or nothing test. If there is a single pixel of a slim boundary that does not fit into a fat boundary, then the basin is labeled as non-Wada. However, it is clear that this can be too restrictive in some cases. For instance, if the basin is obtained by experimental procedures, it is very likely to have some wrong pixels. In these cases, we could complement the merging algorithm with the measure of the fractal dimension of the non-Wada boundary, using a box-counting algorithm on the resulting image of the non-Wada points, for example. If the fractal dimension of these non-Wada boundary points is close to zero, then we can admit that the basins have the Wada property, despite the misbehaved pixels. In any case, the merging method is able to determine whether a basin is Wada or not up to a given resolution using minimum requirements.

Conclusions
In the study of the asymptotic behavior of nonlinear dynamical systems, Wada basins appear frequently. Initial conditions lying in the boundary of Wada basins can suffer an arbitrarily small perturbation leading the trajectory to any of the possible attractors of the system. This supposes a special kind of unpredictability different not only from basins with smooth boundaries, but also from other fractal basins 6,17 .
In this paper, we have seen how the action of merging the basins reveals a new aspect of Wada basins. Actually, Wada boundaries are those that remain unaltered under the action of merging the basins. This perspective provides a new way to test Wada basins, that is faster than previous methods by orders of magnitude, and also much easier to use. Given a basin with three attractors with 1000 × 1000 initial conditions, it takes around two seconds to test the Wada property in a personal computer. Furthermore, no knowledge of the underlying dynamical system is needed. This means that this method is especially suitable for cases in which the exact equations and parameters governing the phenomenon are unknown. Besides, this method can be easily automated. Previous methods 19,21 required a detailed knowledge of the dynamical system and important computational efforts. The only possible black spot of the merging method is that it tests the Wada property up to a given resolution (this is also true for the Wada test proposed in 21 ). Nevertheless, the merging method is the best option to check the Wada property with minimum requirements. This is why we believe that the merging method will become a fundamental tool in the study of the Wada property with applications to many scientific and engineering contexts.