Pictures of Julia and Mandelbrot Sets/Quaternions

Using More Dimensions
Julia and Mandelbrot sets can, of course, be constructed in spaces of higher dimension than two, but there is not very much gained with this, if we let the picture be a two-dimensional cut coloured in the usual way. We ought to let the cut be three-dimensional and let the set be seen as a body or a surface coloured in the same way as we have coloured a fractal landscape. We will here let the figure be white and let the light be artificial (with parallel rays), so that the picture becomes in grey-tone the nuance growing darker with the distance.

As to the function, we could let it be a mapping from the space into itself, but if it does not satisfy something corresponding to the Cauchy-Riemann differential-equations, the structure will be as in the non-complex case. These equations cannot be generalized to the three-dimensional space, but they can be generalized to the four-dimensional space, because in the four-dimensional space (in contrast to the three-dimensional space) we can introduce operations of calculation that are a natural extension of the operations of the complex numbers.

Properties of Quaternions
The extension of complex numbers to four dimensions was discovered in 1843 by the Irish mathematician W.R. Hamilton (1805-65). Hamilton called his new numbers quaternions. The quaternions are the numbers of the form x + yi + uj + vk, where x, y, u and v are real numbers and where the two new symbols j and k, like i, satisfy $$j^2 = k^2$$ = -1, and where the three symbols are connected by ij = k.

Of these relations it follows that jk = i and ki = j, but also that ji = -k, that is, ji = -ij. This means that the multiplication for the quaternions is not commutative. But apart from this, the quaternions, like the real numbers and the complex number, make up a field: you can operate with them exactly as you operate with real and complex numbers. The skew-field of quaternions is an extension of the field of complex numbers, and the quaternions have the same nice and simple properties as the complex numbers. For instance: for a polynomial p(z) of degree n, the equation p(z) = 0 has exactly n roots (some possibly multiple), and if a quaternion function f(z) is differentiable at every point of an open domain (in the quaternions), then it is differentiable any number of times.

Therefore we can without problems generalize the theory of the Julia and Mandelbrot sets to the quaternions: all we have said (critical points, cycles, potential and distance function, ...) is valid almost without modifications - only the field lines need a comment.

Choices in Rendering Quaternion Fractals
We need a strategy to get from 4 dimensional patterns to 2 dimensions. We assume that our function $$f(z)$$ (z quaternion) is a polynomial with real coefficients, and that it only contains powers of even degree (e.g. $$f(z) = pz^2 + qz^4$$), for then the three-dimensional subspace spanned by 1, i and j is left invariant by the iterations, and therefore we can restrict ourselves to this space. The critical points are complex numbers lying symmetrically around the x-axis, and we construct the Mandelbrot set from a finite critical point and ∞. A Julia set is a fractal surface, and we ought to have two ways of drawing: the filled-in Julia set, which is the complement to the Fatou domain containing ∞, and a Julia set with field lines in the inner Fatou domains (partly because there is nothing to see in an inner Fatou domain, and that the field lines are very decorative, partly because the picture is drawn faster when the domain is filled with field lines, because these lines stop the drawing procedure).

In the space we imagine a plane parallel with the x,y-plane with which we intersect the fractal and look at the part of it lying back (or "below") the plane. The intensity of the light is determined by the distance "down" to the fractal. The intersection plane parallel with the base plane (the x,y-plane) is given by a height h. For the Mandelbrot set the height must be positive, as we are not interested in the part lying under the x,y-plane (because the Mandelbrot set is symmetric around the base plane and abates on the opposite side). The nuance of the colour is calculated by consecutive estimations of the distance to the boundary and, on the basis of these estimations, approximations in steps becoming smaller and smaller. If the boundary is too thin, the sequence of approximation may jump past the fractal (giving a black dot).

Technical Detail of Calculation
We can let the next step down be half of the estimated distance, but we must arrange the program so that we can make the steps smaller if there are faults in the picture. We must choose a little number stepmin, so that if the calculated step u is smaller than this number, we set u = stepmin. For each step, we add the calculated step to a number starting with 0, and when we have reached the boundary (or the interior of the Mandelbrot set or a filled-in Julia set or a field line), we divide the distance by the highest allowable distance, which in case of a Mandelbrot set is the height h of the plane, and in case of a Julia set is h + the maximum distance we can go below the base plane (e.g. 2). The result is a number in the interval [0, 1], and we construct a grey-tone scale indexed by the numbers t in this interval, such that 0 corresponds to white and 1 to black. We can, for instance, let the colour have (equal) RGB-values that are the integral part of the number $$255/(1+turndown+tan(t \pi/2)/intensity)^{exponent}$$, with three parameters to adjust the nuance.

For the Mandelbrot set the iteration stops when the maximum number of iterations is reached or when the iterated point is outside a sphere of a very large radius. In the latter case we calculate the distance to the boundary, and let this number determine the next step down. The same applies to a filled-in Julia set. For an inner Fatou domain, the iteration stops when a point in the sequence is within a given little distance from a (given) point in the cycle. As we do not calculate the real iteration number, we need not calculate the attraction of the cycle, but we must know the order of the cycle and (on account of the field lines) the unit-quaternion corresponding to the angle of rotation in the complex case.

Relation between Quaternion and Complex Fractals
The intersection of the Mandelbrot set with the complex plane (that is, the base plane) is the corresponding complex Mandelbrot set, and this lies symmetrically around the x-axis (because the coefficients in our formula are real). The quaternion-Mandelbrot set (in the 3-dimensional space) is obtained by rotating the complex Mandelbrot set around the x-axis. The boundary of the Mandelbrot set is thus a rotary symmetric fractal surface with the x-axis as generator, and it consists consequently of circles around the x-axis. The Julia set associated to a point in the complex plane (that is, of height 0) also consists of circles around the x-axis. But for a point outside the complex plane the Julia set consists of closed curves or segments of curves lying askew  in relation to the x-axis. In these fractal surfaces the structure is only fractal in one direction: in the direction orthogonal to this it will have the true character of a curve - a phenomenon we have seen in the non-complex fractals.

The Field Lines
The norm of the quaternion z = x + yi + uj + vk is the real number |z| = $$\sqrt{x^2 + y^2 + u^2 + v^2}$$, and z can be written as the product of the norm and a unit-quaternion (quaternion of norm 1), this unit-quaternion is the argument of z.

In our case where the fourth coordinate is 0, the argument is a point on the unit-sphere. And the little circle around one of the points in the cycle, which in the complex case stops the iteration, must now be replaced by a little sphere, and the field lines are constructed on the basis of regularly situated domains on this sphere. As such, we can choose the domains that result from a regularly decomposition of the equator circle and from the projection onto this.

Technical Detail of Calculation
The technique used to draw the field lines is, however, a little complicated. Drawing is done by working our way down in as large steps as possible, but now it must be examined if a field line is hit. As the steps initially are large, we risk jumping past the field line, therefore we must (by a trial iteration) examine if a field line is in our way, and if so, alter the further procedure. If a field line is in our way, we must go forward in smaller steps, and when it is hit, we must go backwards and forwards in smaller and smaller steps, in order to find the precise point of intersection with the field line. The steps we go down are usually (that is, when we disregard field lines) half of the estimated distance, and we may now have to reduce the steps.

If we set $$g(z) = f(z) + c$$, and if the order of the attracting cycle is r and if z* is a point in the cycle, then z* is a fixed point for $$g(g(...g(z)))$$ (the r-fold composition), and near z* this map has (in connection with field lines) character of a rotation with the argument $$\beta$$ of the quaternion $$(d(f(f(...f(z))))/dz)_{z*}$$ = the product of the quaternions $$f'(z_i)$$, for the r points in the cycle.

A given direction $$\theta$$ (= quaternion of norm 1) outgoing from z* determines a field line, and this consists of the points z such that if the argument of $$z_k - z*$$ is $$\psi$$ ($$z_k$$ = the last point in the iteration), then $$\psi \beta^{-k} = \theta$$. Therefore we must know the k'th power of a unit-quaternion, and this can be calculated by this formula:


 * $$\displaystyle (x, y, u, v)^k = (\cos(k \kappa), ty, tu, tv)$$

where $$\kappa$$ is an angle such that $$x = \cos(\kappa)$$ and where t is calculated recursively by the following operation performed k times and starting with t = 0 and i = 0:


 * $$\displaystyle t = \cos(i \kappa) + tx$$ and $$\displaystyle i = i + 1$$.

Let n be the number of field lines and let t be their relative thickness (a number in the interval [0, 1]). For the point z, we have calculated $$\psi \beta^{-k}$$, which is a unit-quaternion and therefore corresponds to a point on the unit sphere. As the domains of this determining the field lines are chosen such that they result from a regular decomposition of the equator circle, the condition depends only on the argument of the projection of this quaternion on the complex plane. If v is this angle (chosen in the interval [0, 2$$\pi$$[) divided by 2$$\pi$$, then z belongs to a field line if |v - i/n| < t/(2n), for one of the integers i = 0, 1, ..., n.

However, as we want to know if a field line is in our way, we do first perform a pre-calculation which calculates the distance to the boundary and determines the field line in question, and this is the one having the number i such that |v - i/n| < 1/(2n). We calculate the number d = 2n|v - i/n| - t. If d is negative, the point z actually belongs to the field line, and the point is coloured white. Otherwise, we regard d as a measure for the distance to the field line, and we multiply the calculated step by d, so that the steps become smaller when we are near a field line.

We must find the precise point in our way down where the field line is hit. Therefore we go backwards and forwards in smaller and smaller steps. We call the actual height h1 and the old height h2. When we for the first time are inside a field line, we set h3 = h1 and h1 = (h1 + h2)/2, and calculate again if we are inside the field line. If no, we set h2 = h1 and h1 = (h1 + h3)/2, if yes we set h3 = h1 and h1 = (h1 + h2)/2. We do this until h2 - h1 < stepmin. The distance to the field line is then h - h1 (h is the height of the plane).