How does the main cardioid appear in the Mandelbrot set? I also wonder why something "weird" happens at a point with coordinate 2 on the actual coordinate line, I mean why is the point a sort of bounding point? (this question has already been raised, but I have not found a clear answer)
Why the cardioid shows up in the Mandelbrot set
complex-dynamicsfractals
Related Solutions
There are two research papers worth looking at to get a general idea of what's going on:
- The Mandelbrot set is universal (already pointed out by lhf), and
- On the dynamics of polynomial-like mappings.
Both of these papers make it clear that close copies of the Mandelbrot set often appear in the bifurcation locus for a parametrized family of holomorphic functions. While the focus is on rational families of functions, many of the results extend to transcendental families when properly restricted. So, first, let's try to make it clear how it is that a parametrized family arises in your context. You're trying to invert the tangent function. That is, given a complex number $c$, you want to solve $$\tan(z)=c$$ and you want to do so using Newton's method. Put another way, you want the roots of $f(z)=\tan(z)-c$. Thus, the corresponding Newton's method iteration function is $$N_c(z)= z-f(z)/f'(z) = z-\sin (z) \cos (z) + c \cos ^2(z).$$ This is exactly what we mean by a parametrized family; each $N_c$ is a function of $z$ that also depends upon the parameter $c$.
Next, what exactly do we mean by the bifurcation locus $B(F_c)$ of a parametrized family $F_c$? This is defined in the first paper above as the set of all complex parameters $c$ such that the number of attracting cycles of $F_c$ is not locally constant. For the now classic family $F_c(z)=z^2+c$, this results in the boundary of the Mandelbrot set. For points on the boundary of the Mandelbrot set, there are $c$ values nearby but in the exterior such that $F_c$ has only the point at $\infty$ as an attracting cycle. But there are also points nearby in the interior which have a finite attracting cycle. There are two other characterizations of the bifurcation locus given in the paper but I think that this is the simplest.
Now, how do we plot a bifurcation locus for a more general family? A standard technique is to compute the orbit of a critical point under iteration of $F_c$ for each $c$ chosen from a grid of values and then shade the point $c$ according to the type of periodicity found. Often, this involves some escape criteria but, as we move from $z^2+c$ to more general types of functions, a bit of analysis is sometimes necessary.
Unfortunately, your family isn't quite "polynomial-like". I'm not going to worry about the precise definition here as it's easy to see what the problem is by examining some iterates. In fact, to motivate the next step, let's examine the orbit of the point $z_0=0$ under iteration of $N_{\pi}$. As the $c$ and $z_0$ are both real, we can do so with a cobweb plot:
The orbit is clearly unbounded yet a mini-Mandelbrot appears in your picture at $c=\pi$; thus, it seems we should have some sort of attractive behavior. The picture suggests, though, that perhaps we should just mod out by $2\pi$. That is, we define $$NM_c(z) = N_c(z) \: \text{mod} \, 2\pi$$ and we study the dynamics of $NM_c$ on the strip $0\leq \text{Re}(z) <2\pi$. In fact, this works perfectly, because $$N_c(2\pi) = c+2\pi \: \text{ and } \: N_c'(2\pi) = 0 = N_c'(0).$$ That is, $NM_c$ is continuous and differentiable. Technically, we are iterating on a complex manifold and the dynamics of $N_c$ on $\mathbb C$ are conjugate to $NM_c$ on that manifold.
Taking all this into account, I came up with the following image for the bifurcation locus of $NM_c$ on the strip $0\leq \text{Re}(z) <2\pi$:
This image is obtained as follows: We iterate $NM_c$ from the critical point $z_0=0$ until either the imaginary part of the iterate exceeds in 10 in absolute value or we reach 100 iterates. If the iterate did escape, we shade the point blue-ish based on how fast the iterate escaped. If the iterate did not escape, we shade the point yellow if the iterate converged to a root of $\tan(z)-c$ or black otherwise.
The black region is exactly where we see Mandelbrot like sets and we can use it to find fun Julia sets. The red point shown in the figure lies at $c\approx 3.5613938602255226$. The critical orbit of the corresponding function $NM_c$ is super-attracting with period 4 and it Julia set looks like so:
For the promenade around the main cardioid, the point in the Julia set you want is a fixed point $z = f_c(z)$. By tuning (renormalization), the point you want in the "little Julia set" at the center of the Julia set is a periodic point $z = f_c^p(z)$, in this particular case $p = 3$. You can use Newton's method to find the root numerically, the degree is too high for symbolic solution via radicals. You need a good initial guess for Newton's method to converge to the "correct" tuned fixed point (periodic point), as there is another one at the "tip" of the little Julia set (just as there two fixed points for $z = f_c(z)$), and moreover the fixed points $z = f_c(z)$ are also (unwanted) solutions for $z = f_c^p(z)$.
Here is a small C99 program using my mandelbrot-numerics library:
/*
gcc julia-spiral-center.c `PKG_CONFIG_PATH=${HOME}/opt/lib/pkgconfig pkg-config --cflags --libs mandelbrot-numerics` -lm -fopenmp -O3 -std=c99 -Wall -Wextra -pedantic
LD_LIBRARY_PATH=${HOME}/opt/lib/ ./a.out |
ffmpeg -framerate 25 -i - -pix_fmt yuv420p -profile:v high -level:v 4.1 -crf:v 20 -movflags +faststart julia-spiral-center.mp4 -y
*/
#include <complex.h>
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#ifndef M_PI
#define M_PI 3.141592653589793
#endif
// git clone https://code.mathr.co.uk/mandelbrot-numerics.git
// make -C mandelbrot-numerics/c/lib prefix=${HOME}/opt install
#include <mandelbrot-numerics.h>
double cnorm(double _Complex z)
{
double x = creal(z);
double y = cimag(z);
return x * x + y * y;
}
int main(int argc, char **argv)
{
// silence -Wunused-parameter
(void) argc;
(void) argv;
// number of animation frames
const int frames = 25 * 30 * 1;
// parameters for Julia set graphics rendering
const int iterations = 1024;
const double escape_radius_squared = 1 << 24;
const int width = 640;
const int height = 360;
unsigned char *pgm = malloc(width * height);
// how far from the target component to place the moving c
// 1.0 is the boundary of the component
const double radius = 1.25;
// period of the target component
const int period = 3;
// these don't need to be exact, they will be changed by m_d_interior()
double _Complex z = 0;
double _Complex c = -1.754;
// maximum number of steps of Newton's method
const int steps = 64;
// needs a good initial guess to converge to the correct attractor
double _Complex w = -0.1 - 0.01 * I;
for (int frame = 0; frame < frames; ++frame)
{
// path around the exterior of the component
double t = (frame + 0.5) / frames;
double _Complex q = radius * cexp(2 * M_PI * I * t);
// find c,z such that f_c^p(z) = z and d/dz f_c^p(z) = q
m_d_interior(&z, &c, z, c, q, period, steps);
// find w such that f_c^p(w) = w
// use previous frame's w as initial guess for continuity
m_d_attractor(&w, w, c, period, steps);
// scaling for view
double r = cabs(w);
// render Julia set
#pragma omp parallel for
for (int j = 0; j < height; ++j)
{
for (int i = 0; i < width; ++i)
{
// view centered on w
double _Complex z0 =
( ((i + 0.5) / width - 0.5) * 2 * width / height
- I * ((j + 0.5) / height - 0.5) * 2 ) * r + w;
// initialize derivative with pixel spacing
double _Complex dz = 2 * r / height;
// unescaped pixels will be grey
pgm[j * width + i] = 128;
for (int k = 0; k < iterations; ++k)
{
dz = 2 * z0 * dz;
z0 = z0 * z0 + c;
if (cnorm(z0) > escape_radius_squared)
{
// compute exterior distance estimate
double de = 2 * cabs(z0) * log(cabs(z0)) / cabs(dz);
// colour boundary based on distance estimate
pgm[j * width + i] = 255 * tanh(de);
break;
}
}
}
}
// output PGM stream to stdout
fprintf(stdout, "P5\n%d %d\n255\n", width, height);
fwrite(pgm, width * height, 1, stdout);
fflush(stdout);
}
// cleanup
free(pgm);
return 0;
}
The source code for the two functions I used from the library can be browsed at:
They are both applications of Newton's root finding method, in 2 complex variables and 1 complex variables respectively.
You can view the video rendered by this program.
Best Answer
Not responding to the cardioid question. Responding to the $|z| > 2$ implies divergence question.
Let $c = z_0 = x+\mathrm{i} y$ be the starting point for an iteration with the distance between $z_0$ and the origin greater than $2$. That is, $|z_0| = \sqrt{x^2 + y^2} > 2$. It is convenient to also represent this point in polar coordinates, $z_0 = |z_0| \mathrm{e}^{\mathrm{i} \theta_0}$, where $\theta$ is the angle on the plane to the point, measured anticlockwise from the positive real axis.
Applying the iteration for the Mandelbrot set, we get the point $$ z_1 = z_0^2 + c = |z_0|^2 \mathrm{e}^{\mathrm{i} 2\theta_0} + |z_0| \mathrm{e}^{\mathrm{i} \theta_0} $$ This is the sum of two complex numbers. The first has polar angle $2\theta_0$ and the second has polar angle $\theta_0$. Let's be a little lazy and not try to figure out how those two angles are related, but just try to find out the range of things that could happen using the triangle inequality.
If $2\theta_0$ and $\theta_0$ point in opposite directions (are antiparallel), the sum gives the shortest length it can: $$ \text{shortest case}: |z_1| \geq |z_0|^2 - |z_0| = (|z_0| - 1) |z_0| \text{.} $$ Since $|z_0| > 2$, this is the product of a number bigger than $1$ with $|z_0|$, so is a number bigger than $|z_0|$.
If $2\theta_0$ and $\theta_0$ point in the same direction (are parallel), the sum gives the longest length it can: $$ \text{longest case}: |z_1| \leq |z_0|^2 + |z_0| = (|z_0| + 1) |z_0| \text{.} $$ Since $|z_0| > 2$, this is the product of a number bigger than $3$ with $|z_0|$, so is a number bigger than $|z_0|$.
If $2\theta$ and $\theta$ point in less special directions, the length is between the two above, so is still a number bigger than $|z_0|$.
What we have shown is that if your starting point has $|z_0| > 2$, then $|z_1| > |z_0| > 2$. We can continue the analysis and discover $|z_2| > |z_1|$, then $|z_3| > |z_2|$. This tells us we have an increasing sequence of radii. It is not quite enough to show that the sequence of radii diverge. For that we have to show that the amounts of increase don't eventually decrease to zero. (If the radii have to increase by increasing amounts, the radii will ultimately diverge to infinity.) So let's look at that iteration more carefully.
To remind, $|z_{i+1}| > |z_i| > 2$ for $0 \leq i \leq n$. Suppose $|z_n| = 2+\varepsilon$, where $\varepsilon > 0$, so we have a name for how far above $2$ we are at this step. As above, we have $$ (|z_n| - 1)|z_n| < |z_{n+1}| \text{.} $$ We examine \begin{align*} (|z_n| - 1)|z_n| &= (1+\varepsilon)(2+\varepsilon) \\ &= 2 + 3\varepsilon +\varepsilon^2 \text{.} \end{align*} When $\varepsilon$ is small, its square is very small, so when the radius is just above $2$, we triple the distance above $2$ on every iteration. When $\varepsilon$ is large, its square is larger than its triple, so when the radius is far above $2$, we square the distance above $2$ on every iteration. In any case, for a starting point with radius greater than $2$, the sequence of radii of its iterates increases, more and more rapidly, diverging to infinity.