Mandelbrot set as a particular Julia set of Fractional Order, echipotential lines and external rays of Mandelbrot and Julia sets of Fractional Order

In this paper it is shown analytically and computationally that the Mandelbrot set of integer order are particular cases of Julia sets of Caputo s like fractional order. Also the differences between the fractional-order Mandelbrot set and Julia sets from their integer-order variants are revealed. Equipotential lines and external rays of Mandelbrot set and Julia sets of fractional order are determined.


Introduction
The fractional-order (FO) Mandelbrot and Julia sets in the sense of q-th Caputolike discrete fractional differences, for q ∈ (0, 1), generated by the quadratic complex (Mandelbrot) map f c (z) = z 2 + c, with z and c complex and starting from the initial value z 0 = 0 (the critical point), are introduced in [2] (see also [37,4,5,6]).The algorithms for generating FO sets, base on the known Mandelbrot set and Julia sets of Integer Order (IO) which, after they were discovered, still represent a huge source of inspiration for computer graphics programmers but also for mathematicians.The first who drawn the Mandelbrot set of IO are Robert W. Brooks and Peter Matelski in 1978 [8] before the American-French-Polish mathematician Benoit B. Mandelbrot made it famous and given him the importance and the place in the chaos theory [9].Moreover, the Mandelbrot set, an invariant universal set, is considered to play a similar role as π and e have in mathematics.History of the origin of the Mandelbrot can be found at [27].Before the born of the Mandelbrot set, the study of the dynamics of complex maps was initiated by P. Fatou and G. Julia in the early of twentieth century [1,3].For fractal structures see e.g.[11,12,13,14], while for details and background on Mandelbrot set and Julia sets see [10,15,16,11,17].Mandelbrot set can be considered as a book with an infinity of pages, each page being a Julia set.
While to generate the Mandelbrot set, c is considered variable within a lattice in the parametric plane C, the Julia sets are obtained with fixed c, the origin of iterations of f c (z) being variable in the considered lattice.
The infinite beauty of these fractal sets, generated by the quadratic map, does not represent the subject of this paper, interested readers being directed to, e.g., [11,13,14].
It was found that the fractional calculus represents more accurate the natural behavior in areas bioengineering, recurrent neural network, of engineering, image encryption, electronics, viscoelasticity, robotics, control theory and so on (see e.g.[33,34,35,36,38,39,40]).First definitions of a fractional difference operator are proposed in 1974 [18].Aspects related to Caputo fractional sums and differences can be found in [19,20], while Initial Value Problems (IVPs) in fractional differences are studied in [30].The stability of fractional differences is snalyzed in [21,22], and weakly fractional difference equations, symmetry-breaking of fractional maps can be found in [23,24,25].For nonexistence of periodic solutions see [26].
In this paper the fractional Mandelbrot and Julia sets are introduced and some important properties are analytically and computationally presented.The discrete fractional calculus is used as a natural extension of difference calculus and Mandelbrot's idea of creating fractals from iterative of complex mappings to study iterations of complex fractional difference equations is extended.
Notations utilized in this paper: • IOM: Mandelbrot of IO; • IOK c : Filled Julia set of IO; • IOJ c : Julia set of IO; • F OM: Mandelbrot set of FO; • F OK c : Filled Julia set of FO;

Mandelbrot set and Julia sets of FO
Next a briefly recall of some elementary notions about Mandelbrot set and Julia sets of IO required by the FO counterparts is presented (see [11,13] or [7] for more information).
The iteration of f c with z(0) = 0 generates the sequence which will be used to generate Mandelbrot sets, while for z(0) ̸ = 0, the sequence of iterates becomes used to generate Julia sets.The Mandelbrot set of IO, IOM, is the set of complex values c for which the absolute value of z n remains bounded and does not tend to infinite, for all n ≥ 1, |z(n)| < M , usually r taken as r = M [13], but could be taken even of order of thousands.
To define, for a fixed c, the Julia sets of IO, IOJ c , let consider the attraction basin of ∞, A c (∞), the set of points z 0 which tend to ∞ through the iteration (1) The boundary of A c (∞), which depends on c, represents the Julia set of IO, IOJ c [13] Another notion related to Julia sets of IO considered in this paper, is the filled Julia set of IO, IOK c , which for a fixed c, is the set of all points z 0 ∈ C for which the orbit (2) remains bounded The IOJ c set is contained in the IOK c set and is the boundary of the IOK c set [13] ∂IOK c = IOJ c = ∂A c (∞).
In this paper for computer graphic reasons, without loss of generality, one consider the analysis and computations of the IOK c sets.
To draw the IOM and IOK c the so called escape-time algorithm (direct algorithm) is used [7].For Mandelbrot set, z 0 = 0, and c is varied within a finite complex parametric domain (lattice), while for Julia sets, z 0 is varied within a finite lattice for fixed c.If after a finite number of iterations N of f c the modulus |z n |, n = 1, 2, ..., N , remains bounded (|z N | ≤ 2), then c in the case of IOM set, or z 0 in the case of IOK c sets, belongs to the IOM or IOK c , respectively.Otherwise, c or z 0 does not belongs to IOM or IOK c , respectively.
In the graphical representations of this paper, the complex plane of Mandelbrot sets is the plane (c x , c y ) (the parameter space), where c x and c y are the coordinates of c, while the Julia sets are drawn in the complex plane (x, y) of the initial point z(0) of coordinates x and y.The sets IOM and IOK c are usually plotted as color, most often black, while outside points of these sets can be plotted colored using smooth coloring schemes or black-white [11] i.e.IOM set is the set of all parameters c for which IOK c is a connected set; P3.For c chosen at the boundary of the Mandelbrot set the related IOK c set is a "dendrite", while for c situated outside of IOM the corresponding IOK c is a Cantor sets, or "dust"-like, composed of infinitely many disjoint points [12] (see Fig. 1 (c) and (d) where the IOK c sets corresponding to points A and B are presented).

Mandelbrot and Julia maps of FO
Let the time scale N a = {a, a + 1, a + 2, ...}.The q-th Caputo-like discrete fractional difference of a function u : N a → R, for q > 0 and q ̸ ∈ N, is defined as is defined [31] as where t ∈ N a+n−q and n = [q] + 1.
∆ n is the n-th order forward difference operator, while ∆ −q a represents the fractional sum of order q of u, namely, The falling factorial t (q) , is defined as follows: .
Note that the fractional operator ∆ −q a maps functions on N a to functions on N a+q (for time scales see, e.g., [32]).
For q ∈ (0, 1), when ∆u(s) = u(s + 1) − u(s), n = 1, and starting point a = 0, case considered in this paper, q-th Caputo's difference, ∆ q , becomes Then, the real FO autonomous Initial Value Problem (IVP) in the sense of Caputo ∆ q u(t) = f (u(t + q − 1)), t ∈ N 1−q , u(0) = u 0 , with f a continuous real valued map and q ∈ (0, 1), has the numerical solution with the commonly used form in numerical applications Consider next the complex variant of the IVP of FO with q ∈ (0, 1), z = x + ıy ∈ C, scaled c within a parametric complex domain, and z 0 ∈ C.Then, the numerical integral (4) becomes [7] z which represents the mathematical description of the complex logistic (Mandelbrot) map of FO used to generate the F OM set (with z(0) = 0) or F OK c sets (with z(0) variable).
Remark 1.To note that whole for the IOM (IOK c ) set one iterates f c , usually few dozens iterations, N , to generate the F OM (F OK c ) set, every iteration of (5) necessary to obtain the set, the expression of z(n) requires the calculation of the sum in the right hand side of (5) for each n = 1, 2, ..., N .Only small values of N , of order of few tens (e.g.N = 20 as for the IOM set) does not provide good accuracy in the calculation of z(n).Therefore, a compromise between a relatively higher N (time consuming), e.g.N = 70 − 100, and computing time is desirable.
To obtain the F OM set, one iterates (5), for c scanning a complex domain (usually a rectangular lattice [7]).The set of points c for which the sequence of modules |z(n)| remains bounded after a finite number of iterations N , forms the F OM set.
To obtain the F OK c sets one fix c (see P2 and P3) and one iterates (5) with z(0) variable within a complex domain.Like for the F OM set, if after N iterations |z(n)| remains bounded, z(0) belongs to the underlying F OK c set.In Fig. 2 (b) is drawn the F OM for q = 0.5, while in Figs. 2 (c),(d) and (e) several F OK c sets are drawn for q = 0.5 and c = −0.0781+ 0.6694i, c = −1.0516+ 0.0913i, c = −0.5, respectively.Details on FO algorithms can be found in [7], while a matlab code for F OM can be found at [28].

Properties of the F OM set
Several properties of the F OM map are analyzed in [7], especially for real c.In this paper some properties are analytically proved and verified with aided of the scientific computation.
Contrary to expectations, the IOM set is not a particular case of the F OM for q = 1, but only for q ↓ 0 as proved bellow Proposition 2. The IOM set is the F OM set for q ↓ 0.
Proof.Consider the limit of z(n) for q ↓ 0, with z(0) = 0 in (5) rearranged as follows and for a finite times of iterations, n−1 i=1 and, therefore, lim or, in simplified form i.e. the IOM map.
Because at q = 0, Γ(q) has a simple pole singularity, in simulations q cannot be set 0 and, therefore, it is considered numerically, e.g., q = 10 −m with m some positive integer.In Fig. 2 (a) the F OM set is presented for m = 6.
An animation showing the metamorphosis of the F OM sets for q varying from q = 1 to q = 0 is presented as a supplementary video.
For c ̸ = 0, Proposition 2 no longer takes place for the F OK c map, because in this case one obtains: lim which by iteration, generates the sequence (6) different from the sequence generating the IOK c sets (see (3)) For example, compare the F OK c set and the IOK c obtained for c = −0.5 in Figs 2 (e) and (f), respectively.However, for c = 0, probably the most important property is the following Proposition 3. The IOM set is the IOK c set for q ↓ 0 and c = 0 Proof.Consider, as required by the F OK c set, z 0 variable within a complex lattice, and q ↓ 0. For c = 0, the sequence (6) generating the F OK 0 , becomes To note that, while the F OK 0 set with q ↓ 0 for c = 0 is identic with the IOM set (Fig. 3 (a)), the IOK 0 set for c = 0 is, as known [11], a filled disc (Fig. 3 (b)) While generally it is considered that a FO continuous or discrete system for q = 1 should identify with his own IO variant, the next result shows another surprising property of F OM map (see also.[25,42,43,44] for differences between real FO systems and their IO counterparts) Proposition 4. For q ↑ 1 the F OM set differs from the IOM set.
Proof.Consider the limit of z(n) for q ↑ 1, with z(0) = 0 in (5) rearranged as follows or, in simplified form i.e., for q ↑ 1, the F OM map generates a set different from the IOM set.
The difference revealed by Proposition 4 can be viewed in Figs. 4, where the IOM set and F OM set, for q = 1, are presented.

Basic notions on echipotential lines and external rays for IOM and IOK c sets
The external arguments theory of the IOM set has been developed in [45,46] and popularized in [13] and makes use of an analogy to electrodynamics.It was shown that the exterior of the IOM set can be viewed as an electrostatic field.Consider, as described in [52,45,46], a capacitor made of a hollow metallic cylinder with a great diameter inside of which an axis of aluminum shaped in such a way that its cross-section is the Mandelbrot set.The ensembles of cylinder and axial bar are supposed to be infinitely long.In other words, one has an aluminum bar with cross section the Mandelbrot set, situated in the middle of a large hollow metallic cylinder.
If the interior bar is set at potential 0 and the exterior cylinder at a high potential, between the two metallic pieces appears an electric field which creates surfaces in the surrounding space.If one considers an orthogonal section through this ensemble of metallic corps and echipotential surfaces, one obtains the IOM set, surrounded by echipotential curves (lemniscates), sections through the equipotential surfaces with constant potential.It has been proved that the equipotential lines are also lines of equal escape time in the time escape algorithm to generate the IOM or IOK c sets [52].
A particle starting from the frontier of the IOM set will reach the great circle surrounding the IOM set by following the external rays, perpendicular curve on the equipotential lines, being gradient lines of potential.
The equipotential curves are given by while the external rays are defined as follows where, for the Mandelbrot map, ) 2 + c,...,n = 1, 2, ... Using the Green's function, the Douady-Hubbard potential U (c) of a point c situated between the cylinder and the outside of the Mandelbrot set can be written [48,49] (see (7)) U (c) is zero at points c belonging to the boundary of the IOM set, while for large c U (c) is approximated by log(c).
An equipotential line defined by a constant c ∈ [2, ∞), is a closed curve surrounding the IOM set and is defined as the set of points c with property: {c |U (c) = c}, i.e. the set The external argument θ(c) of an external ray that pass through a point c, with a large |c|, and which determines the point where it reaches the great circle, is the argument of the function Φ(c) given as follows [50] The external argument of Φ(c) is (see (8)) and, therefore, the external ray for a fixed angle θ ∈ [0, 2π) is the locus of points c in the complex plane that have all a same external argument with property: {c |θ(c) = θ}, i.e. the set By arg(c) ∈ [0, 2π) is denoted the principal value of the argument of a complex number.To every point on the frontier of the IOM set there could exist several external rays.

Approximations of equipotential lines and external rays
To generate computationally equipotential curves and external rays using the above relations is a quite hard task.However, there exist several other simpler methods which can be applied both to the IOM set and IOK c sets. 1.In [13] the potential U (c) is approximated by the value Ū (c) defined as following: if |z n | > M , while f c is iterated, where M is e.g.10000 [13], the potential can be approximated by Ū (c) ≈ log |z n |/2 n , otherwise, i.e. |z n | remains smaller than M , the potential is set 0. If, for a considered constant c, Ū (c) is close or equal to c, then c, or z 0 in the case of Julia sets, belongs to the equipotential line and the point c, or z 0 , respectively, is plotted.An even more simple method is the Level Set Method (LSM [13]) which, to a point c, or z 0 , respectively, within a complex lattice, attributes a color (e.g.black), depending to the number of iterations of f c for which |z n | remains bounded.Therefore, for each n ∈ [1, N ] one obtains a level set which, approximatively is identic to an equipotential line.2. On the other side, the external rays (9) can be approximated by the Binary Decomposition Method (BDM [13]) with respect to fixed angle α 0 ∈ [0, 2π).Thus to a point c, or z 0 , respectively, within a complex lattice one attributes a color (e.g.black or white) if the argument of z n , arg(z n ), belongs or not to the intervals 2 n α 0 ≤ arg(z n ) ≤ 2 n α 0 + π( mod 2π) Drawing precisely of an external ray, e.g.inside a detail of IOM or F OM set is strongly restricted by the number of bits of the floating-point arithmetic used by the computer program [48].Therefore, the external rays cannot be drawn in certain details with computer programs (see e.g.[51]) using the double 64 bits format.Also, deep studies of external arguments require external arguments measured not in radians, but as fractions of complete turns.Using this unit, most of the notable points of the Mandelbrot set boundary have rational external arguments [48,49].

Echipotential lines and external rays of F OM and F OK c sets
The external rays for rational angles land to the frontier of the IOM set and connected IOK c sets and, as verified numerically in this paper, this property holds for FO sets too.
Because writing analogue expressions of equipotential lines (7) and external rays (8) for the F OM or F OK c maps, called hereafter Fractional Equipotential Lines (FEL)s and Fractional External Rays (FER)s, are a difficult task, a possible overcome is to adapt the approximations given in Subsection 5.1.1where z n is used instead f n c .In Fig. 5 (a) are presented the FELs of the IOM set and F OM set for q = 10 −6 (see Property 2), and in Figs 5 (b) and (c) are presented the FELs for F OM sets for q = 1 and q = 0.5 respectively.
As known, the following property of equipotential hods Proposition 5. Equipotential lines can not intersect each other.
While for the IOM set and the F OM set, for q = 1, this property is obviously verified computationally (Fig. 5 (a)), for the F OM set with q < 1 this property seems to be no more verified and the FELs intersect (see Figs 5 (b) and (c)).
Similarly, the FELs of F OK c sets seem to intersect (see Figs 5 (d),(e),(f) where the F OK c sets are determined for q = 0.5 and correspond to c, chosen at the points denoted A, B and C in the F OM set in Fig. 5 (c)) To draw FERs one divided the complex plane C into sectors where we set the same color if θ is within some certain interval [47].Another simple way to identify the external rays is to use the ordinary escape-iterations algorithm with large escape radius and plot those points c for which ℑ(z n ) > 0. In Figs 6 (a), are presented as overplot the FERs for the IOM set and F OM set for q = 10 −6 , in Fig 6 (b) and (c) the FERs for F OM sets with q = 1 and q = 0.5 respectively and Figs 6 (d), (e) and f are drawn the FERs for F OK c sets corresponding to the points A, B and C taken in the F OM set in Fig. 5 (c).
Conflict of Interest: The authors declare that they have no conflict of interest    Figure 5: a), b), c) FELs of three F OM sets for q = 10 −6 , q = 1 and q = 0.5, respectively; d), e), f) FELs of F OK c for q = 0.5, corresponding to points A, B, C from the F OM, respectively.
, [41] (see Fig. 1 (b) for the IOM and Fig. 1 (c) and (d) for IOK c sets generated for c considered as points A and B in the IOM set).Some of the most important properties of the IOM set and IOK c sets analyzed in this paper are : P1.B. B. Mandelbrot found empirically that some isolated islands out of the body of the IOM set were actually connected to the mainland by very thin filaments [29] (Fig. 1 (b).Now, it is conjectured that the Mandelbrot set is locally connected [15, 27]; P2.IOK c sets are connected if the underlying c belong to the interior of the IOM, If one denote c = z 0 , one obtains the sequence defining the IOM set c, c2 + c, (c 2 + c) 2 + c, ((c 2 + c) 2 + c) 2 + c, ... i.e. the IOM obtained with the map f c(z) = z 2 + c.

Figure 1 :
Figure 1: a) IOM set; b) Zoomed detail of zone C of the IOM set; c) IOK c set corresponding to c = 0.3596 + 0.6294i (point A exterior to the IOM set); d) IOK c set corresponding to c = 0.3902 + 0.3591i (point B from the interior of the IOM set).

Figure 2 :
Figure 2: a) F OM set for q = 10 −6 ; b) F OM set for q = 0.5; c) F OK c set for q = 0.5 and c = −0.0781+ 0.6694i (see point A in Fig. 2 (b)); d) F OK c set for q = 0.5 corresponding to c = −1.0516+ 0.0913i (see point B in Fig. 2 (b)); e) F OK c set for q = 0.5 corresponding to c = −0.5 (see point C in Fig. 2 (b)); f) IOK c set for c = −0.5.

Figure 3 :
Figure 3: a) The F OK c set for q ↓ 0 and c = 0 is identic to the IOM set; b) The IOK c set for c = 0.

Figure 4 :
Figure 4: Differences between the IOM set (a) and the F OM set for q = 1 (b); c) Overplot of the IOM and F OM set for q = 1 (image from animated video).