分形/彩色Julia集

可视化复函数

HSV投影中的颜色空间

假设我们想要可视化一个复值函数，比如

${\begin{aligned}f:\mathbb {C} &\;\to \mathbb {C} \\z&\;\mapsto w=f(z)\\\end{aligned}}$

为了给$w$着色，我们将其分解为其绝对值 $|w|$ 和其辐角 $\arg(w)$.

然后，我们将颜色

$\mathrm {HSV} \left({\tfrac {1}{2\pi }}\arg(w),\;1-g_{a_{s},b_{s}}(|w|),\;g_{a_{v},b_{v}}(|w|)\right)$

分配给表示$z$ 的点。在这个HSV颜色空间中，所有值都在 0 到 1 之间。HSV 颜色的第一个分量（色调）仅取决于$w$ 的辐角，而第二个和第三个分量（饱和度和明度）仅取决于$w$ 的绝对值。我们在$w$ 上使用变换 $g_{a,b}$ 将其映射到区间 $[0,1]$。关于$g$，参见章节辅助函数.

例子

索引为s（饱和度）的值控制饱和颜色向白色（或灰色）的过渡，即中间值 → 无限大。索引为h（值/色度）的值控制黑色向明亮的过渡，即零 → 非零。参数a控制过渡发生的位置：a只是划分两个区域（暗/亮、饱和/灰等）的圆的半径。参数b控制过渡的锐度：b小 = 平滑，b大 = 锐利。

以下图像均显示范围 [-10,10]×[-10,10]$i$，并使用$a_{s}=5$（彩虹的半径）和$a_{v}=1$（黑色圆盘的半径）。

在 S 和 V 的意义交换的图像中，零用白色打印，而无穷大用黑色打印。

图库

• 一些着色的例子
• 
  $b_{s}=b_{v}=10$
• 
  $b_{s}=b_{v}=1$
• 
  $b_{s}=0.05,b_{v}=1/3$
• 
  $b_{s}=0.05,b_{v}=1/3$
  S 和 V 的角色互换

• 分析 z→z²+c 的临界点
• 
  n=1：起始值（颜色图）
• 
  n=2
• 
  n=3
• 
  n=4
• 
  n=17
• 
  n=18
• 
  n=19
• 
  n=20

可视化 Julia 集

Julia 集 和 Fatou 集 可以很漂亮。它们的计算机图形生成可能很困难。

在本页的剩余部分，我们将研究一种可用于可视化复函数 ƒ 的 Julia 集的方法。它的优点是您不需要知道迭代的吸引子

$z\mapsto f(z)$

生成的图像在 Fatou 部分将是平滑的。

概述：成像方法

请注意，还存在其他方法来对复动力系统进行着色，例如

逆迭代法 (IIM)

计算 ƒ 的原像，即计算反向轨道。由于不动点的稳定性从吸引转变为排斥，反之亦然，人们只需选择一个复数，然后观察它在逆迭代下的变化情况。问题是结果不会均匀分布在$J$中，而且您必须计算 ƒ 的逆。

按吸引速度着色 (CSA)

从用于着色的点阵中取一个值，进行迭代，直到迭代值接近吸引子。根据使点足够接近吸引子所需的迭代次数对该点进行着色。

此方法通常用于可视化多项式的 Julia 集和与牛顿著名方法（用于查找函数的零点）相连的 Julia 集。度数 > 1 的多项式总是具有无穷大作为超吸引不动点。牛顿方法中出现的理性函数总是具有函数的根作为吸引不动点。然而，在这两种情况下，可能存在其他吸引子，而且这些吸引子不一定只包含一个点。

逃逸时间算法 (ETA)

如果 ∞ 是一个吸引子，即过程的一个不动点，那么根据迭代次数（时间）对点进行着色，直到看到点逃逸到 ∞。如果点在最大迭代次数内没有逃逸，那么该点被着色为属于 Julia 集或某个其他吸引子的盆地。此方法适用于多项式。最突出的 Julia 集是z→z²+c 的 Julia 集，其中c 是 Mandelbrot 集的元素或距离 Mandelbrot 集不远。如果您看到这样的 Julia 集的图像，那么很可能使用了 ETA 来获得图像。

柯西收敛算法 (CCA)

在本页的剩余部分，我将介绍一种不同的方法，其基本思想与逃逸时间算法相同。但是，不需要预先知道任何吸引盆地，并且可以将不同的吸引盆地分离并以不同的颜色着色。这种方法使用柯西收敛的概念。这种方法不是观察点的轨道，而是观察两个相邻点z 和 z+ε 的距离在迭代这两个值时的变化情况。如果差值趋于 0，则该点将走向吸引子。如果差值不接近 0，则该点靠近（或部分属于）Julia 集。

度量

令 $\varrho$ 是 紧化 的 复平面 到 黎曼球面 S₂ 上的典范投影。

$\varrho :{\overline {\mathbb {C} }}=\mathbb {C} \cup \{\infty \}\rightarrow S_{2}$

这给了我们一个 度量 d: 作为复平面上两点之间的距离，我们取他们在球面上的距离，即 正形线 的长度。这意味着度量以 π 为界，即使到 ∞（现在是北极）的距离也是有限的。

为了计算两点 z 和 w 之间的距离，我们以这样的方式旋转球体 S₂，使得

1. w 映射到 0
2. z 映射到正实轴

经过这些变换后，距离就可以很容易地计算出来。旋转可以通过一系列 等距 的 莫比乌斯变换 来完成。总之，我们得到

$t_{w}:z\mapsto {\begin{cases}|z|\;,&{\text{ if }}w=0\\{\textstyle {\frac {1}{|z|}}}\;,&{\text{ if }}w=\infty \\{\textstyle {\frac {1}{|w|}}}\;,&{\text{ if }}z=\infty \\&\\{\textstyle {\frac {1}{|w|}}}\cdot \left|{\frac {z{\overline {w}}-|w|^{2}}{z{\overline {w}}+1}}\right|\;,&{\text{ else}}\end{cases}}$

留给读者作为练习。横线表示 复共轭。度量为

$d(z,w)=\,\!2\arctan(t_{w}(z))$

d 引入的一个很好的特性是，以前在无穷远处发散的序列现在收敛到无穷远处，即收敛到 S₂ 的北极。

轨道的稳定性

回顾 Julia 集 的定义，该集合针对 压缩映射 ƒ。该定义暗示了迭代的稳定性的一些事实

$z\mapsto f(z)\mapsto f^{2}(z)\mapsto f^{3}(z)\mapsto \cdots$

其中 ƒⁿ 表示 ƒ 的第 n 次迭代

$f^{\,n}=f\circ f\circ \cdots \circ f$

集合

$\{f^{\,n}(z)\}_{n\in \mathbb {N} _{0}}$

称为z的轨道（在 ƒ 下）。

如果z和w彼此靠近，并且属于Fatou 集 F_ƒ，即 ƒ 的 Julia 集 J_ƒ 的补集，则两个点z和w的轨道在某种意义上表现相似。如果z是J_ƒ的元素，则z和w将表现出相当不同的行为，即使w本身是J_ƒ的元素。

为了获得轨道稳定性的概念，我们设置

$\Sigma _{n}(z)=\Sigma _{n}(z,\varepsilon )=\sum _{k\leqslant n}d\,(f^{k}(z),f^{k}(z+\varepsilon ))$

对于一个小 ε，以及上面的度量d。这意味着我们取两个彼此靠近的点，然后我们总结它们的距离，因为 ƒ 使它们在黎曼球面上跳跃。请注意，对于任何固定的 ε，即使z是 Fatou 点，当n很大时，该和也会发散。

但是，我们可以使用 Σ_n 来衡量一个点离J_ƒ有多近：该和越大，迭代越不稳定，该点离 ƒ 的 Julia 集越近。

为了区分 Julia 集（或接近 Julia 集的点）的点和 Fatou 集中的点，我们需要一个标准。为了获得这个标准，我们计算我们要着色的点的所有 Σ 值，即格点 Γ 中的所有点。计算完这些值后，我们进行一些统计，以获得期望值 E 和标准差 σ 用于 Σ 值集：令

${\begin{aligned}E_{n}&=E_{n}(\Gamma ,\varepsilon )=E\left\{\log(\delta +\Sigma _{n}(z,\varepsilon ))\right\}_{\!z\,\in \,\Gamma }\\\sigma _{n}&={\sqrt {D_{n}}}\\\end{aligned}}$

提醒一下

${\begin{aligned}EX&={\textstyle {\frac {1}{|X|}}}\sum \limits _{x\in X}x\\D\!X&=EX^{2}-(EX)^{2}\\\end{aligned}}$

事实证明，Σ 值广泛分布在多个尺度上。因此，我们不直接使用 Σ。相反，我们使用 Σ 的对数。δ 值只是一个小的常数，以避免对数的输入为零。

着色

现在，我们拥有了为一个点着色所需的一切

1. 选择
   1. 一个点格点 $\Gamma$ 的点 $z$ 进行着色
   2. 要执行的迭代次数 $n$
   3. 一个小 $\varepsilon$
2. 对于所有点，计算 $\Sigma _{n}(z,\varepsilon )$
3. 计算 $E_{n}$ 和 $\sigma _{n}$
4. 计算

$J(z)=\,\ln(\delta +\Sigma _{n}(z))-K$

对于某个常数 $K$。 因为 $K$ 将用于将属于 Julia 集的点 ($J(z)>0$) 与不属于 Julia 集的点 ($J<0$) 分开，$K$ 的合理值大于 $E_{n}$。 尝试设置如 $K=E_{n}+2\sigma _{n}$ 或 $E_{n}+\sigma _{n}$ 等等。 如果 $J(z)>0$，那么我们将 $z$ 着色为属于 Julia 集。 如果 $J(z)<0$，我们可以使用该值来为 Fatou 集着色。 如果我们知道某个吸引子，我们可以检查 ƒ^$n$$(z)$ 是否接近它，并使用该信息。

为了将值映射到饱和度和亮度的有效范围内，我们使用第 辅助函数 h 部分中的函数 $h$。

修改

计算 $E_{n}$ 需要大量时间。 可视化过程需要两次传递

1. 计算来自所有 $\Sigma _{n}(z)$ 的 $E_{n}$
2. 使用 $E_{n}$ 和重新计算的 $\Sigma _{n}(z)$ 对点进行着色

或者

1. 计算并存储所有 $\Sigma _{n}(z)$
2. 计算 $E_{n}$ 并对这些点进行着色

另一种方法看起来像这样

找到最小的 $n$ 使得 $d_{k}(z,z+\varepsilon )$ 在所有 $k\geq n$ 上都低于给定的边界。如果找不到这样的 $n$，则假设 $z$ 是朱利亚集的元素。否则，使用 $n$ 对法图集进行着色。

此算法是逃逸时间算法 (ETA) 的变体。请注意，在 ETA 中，点实际上不会逃逸（至少如果我们在球体上），它只是收敛到 ∞。这种方法类似。但是，我们不需要知道 ƒ 的吸引不动点。

到目前为止，我没有尝试修改后的版本。一个缺点可能是法图集将不再出现平滑的着色。然后我不确定这种修改是否真的有优势，因为迭代必须完成，直到达到给定的最大迭代次数。请注意，即使 $d_{n}$ 在某些 $n$ 上低于边界，距离 $d_{k}$ 可能会再次升高。我无法确定这种影响是否至关重要或可以忽略不计...

图库

• 
  $z\mapsto z^{2}$
• 
  $z\mapsto z^{2}+i$
• 
  $z\mapsto z^{2}-1$
• 
  ${}_{z\mapsto z^{2}-0.742+i/10}$
• 
  $z\mapsto N_{z^{3}-1}$
• 
  $z\mapsto N_{z^{3}-2z+2}$

$N$ 表示牛顿算子

$N:f\mapsto {\mathit {id}}-{\frac {f}{f'}}$

使用临界点吸收 (CPA)

之前的方法可以生成简洁平滑的着色，并且对过程的动力学知识要求最低。然而，它非常耗时。以下方法是多项式逃逸时间算法 (ETA) 的扩展。

令 ƒ 为度数为 d > 1 的复数上的多项式。这样的多项式最多有 d 个临界点：无穷大和 ƒ′ 的最多 d–1 个零点。众所周知，过程 z→ƒ(z) 的每个吸引子至少吸收一个临界点。假设 z_K 是一个临界值，则 ƒⁿ(z_K) 会任意接近 ƒ 的一个吸引周期。

二次多项式 ƒ(z) = z² + c 的过程是最简单的情况：临界值为 0 和  ∞。由于 ∞ 吸收自身，只剩下 0，我们有以下情况。M 表示曼德布罗集合。

$c\notin M$

简单情况。0 被 ∞ 吸收，J(ƒ) 仅由康托尔尘埃组成。可以使用逃逸时间到 ∞ 来对所有点进行着色。

$c\in M\setminus \partial M$

如果 c 是 M 内部的一个元素，则 0 将被一个 (超) 吸引周期吸收。计算

$w=w(c,n)=f^{n}(0)$

对于足够大的 n。由于 w (基本上) 仅依赖于 c，因此可以在开始 J(ƒ) 的可视化之前预先计算它。请注意，w 是可能具有多个元素的周期的元素，即 w 对于该周期来说只是唯一的。

然后，对 C 中的点 z 进行着色与 ETA 一样简单

• 如果 z 被 ∞ 吸收，请使用逃逸时间来对 z 进行着色。
• 如果 ƒⁿ(z) 接近 w，即如果 |ƒⁿ(z) – w| < ε 第一次出现，请使用该 n 来对 z 进行着色。
• 如果 ƒⁿ(z) 既不接近 ∞ 也未接近 w，则将 z 着色为 J(ƒ) 的一部分。只有很少的 z 会属于此类。

辅助函数

辅助函数 t

一些 $t\,$

是一个平滑的、单调的 S形 从 −1 到 1 的转换，它满足 $t'(0)=1$ 和 $t(x)=-t(-x)$。有许多选择。其中一些是

${\begin{aligned}t:\mathbb {R} &\to (-1,1)\\x&\mapsto \tanh \,x\\&\mapsto {\tfrac {2}{\pi }}\arctan \left({\tfrac {\pi }{2}}x\right)\\&\mapsto {\tfrac {2}{\pi }}\operatorname {gd} \left({\tfrac {\pi }{2}}x\right)\\&\mapsto {\frac {x}{\sqrt {1+x^{2}}}}\\&\mapsto {\frac {x}{1+|x|}}\\\end{aligned}}$

gd 表示 古德曼函数。

• 所有函数都在 Desmos 图 中。

辅助函数 h

$h_{a}\,$

辅助函数 $h_{a}$ 与 $t$ 几乎相同，但它映射到 $(0,1)$，我们可以通过参数 $a$ 调整过渡速度。

${\begin{aligned}h_{a}:\mathbb {R} &\to (0,1)\\x\;&\mapsto {\tfrac {1}{2}}+{\tfrac {1}{2}}t(2ax)\end{aligned}}$

其中 $a>0$。那么 $h_{a}$ 关于点 (0, 1/2) 对称，并且

${\begin{alignedat}{2}h_{a}&(0)&&={\tfrac {1}{2}}\\h'_{a}&(0)&&=a\\h_{a}&(-\infty )&&=0^{+}\\h_{a}&(\infty )&&=1^{-}\\\end{alignedat}}$

负值被映射到 0 到 1/2 之间的数值。正值被映射到 1/2 到 1 之间的数值。参数 $a$ 控制过渡的速度。

如果我们想要一个下降函数，我们可以使用对称性

$h_{a}(x)\,=1-h_{a}(-x)=1-h_{-a}(x)$

也就是说，我们取负 $x$ 或 $a$。

辅助函数 g

$g_{a,b}\,$

此函数将正数映射到区间 $[0,1)$。

${\begin{aligned}g_{a,b}:\mathbb {R} ^{\geqslant \,0}&\to (0,1)\\x\;&\mapsto {\tfrac {1}{2}}+{\tfrac {1}{2}}\,t(2b\cdot a\cdot u(x/a))\end{aligned}}$

对于以下定义的某个函数$u$。如果$u$ 选择得当，那么对于$g$ 以下关系成立：

${\begin{array}{lll}g_{a,b}(0)&=&0\\g_{a,b}(\infty )&=&1\\g_{a,b}(a)&=&1/2\\g'_{a,b}(a)&=&b\\g'_{a,b}(x)&>&0{\text{ for }}x>0\end{array}}$

这意味着$g_{a,b}$

• 随着$x$ 从$0$ 到$\infty$ 不断增长，从 0 增长到 1
• 我们可以通过指定参数$a$ 来控制$g$ 穿过 0 到 1 之间中点的位置
• 我们可以通过指定参数$b$ 来控制$g$ 穿过点$(a,1/2)$ 的速度

我们剩下要做的就是确定函数$u$ 在$\mathbb {R} ^{>0}$ 上的形式

$u\,$

$u$ 必须满足以下条件

${\begin{array}{lll}u(0^{+})&=&-\infty \\u(\infty )&=&\infty \\u(1)&=&0\\u'(1)&=&1\end{array}}$

对于 $u$ 我们设置

$u(x)={\tfrac {1}{2}}\left(x-{\tfrac {1}{x}}\right)$

辅助函数 *w*

$w_{a}\,$

此函数将正数映射到区间 $[0,1)$。

${\begin{aligned}w_{a}:\mathbb {R} ^{\geqslant \,0}&\to [0,1)\\x\;&\mapsto t(ax)\end{aligned}}$

通过 $a$ 我们可以控制它在原点的斜率

$w'(0)=\,a$

• 交互式 Desmos 绘制 $w$

圆弧通过 (0,0) 和 (1,1)

一条圆弧通过 (0,0) 和 (1,1)，在 (0,0) 处的斜率为 $a$，在 (1,1) 处的斜率为 $-1/a$

$\operatorname {arc} (x)={\begin{cases}x&;{\text{ if }}a=1\\y_{0}+\operatorname {sign} (\ln a){\sqrt {r^{2}-(x-x_{0})^{2}}}&;{\text{ if }}a\neq 1\\\end{cases}}$

其中

${\begin{aligned}x_{0}&={\frac {a}{a-1}}\\y_{0}&=1-x_{0}\\r&={\sqrt {x_{0}^{2}+y_{0}^{2}}}\\\end{aligned}}$

圆心在直线 $y=1-x$ 上。

• Desmos 绘图

另请参阅

• 对参数平面和曼德勃罗集进行着色
• 通用着色算法： ( 像素的特征 -> 颜色 ): $T:scalarvalue\to colorvalue$