分形/复平面迭代/曼德勃罗集/中心

名称 (同义词)

曼德勃罗集的双曲分量的中心
Mu原子的核 ^[1]

定义

中心是一个双曲分量 H 的参数 $c_{0}\in H\,$ (或参数平面的点)，使得对应的周期轨道具有乘子= 0。"^[2]

初始域

"所有 M 的任何周期的中心都包含在 M 的内部，因此圆 C = {z ∈ C: |z − 0.75| = 2} 包围着所有中心。"
"可以利用 M 的实对称性，只使用具有非负虚部的起点；"^[3]

周期

中心的周期 = 临界轨道的周期 = 双曲分量的周期

朱利亚集

它是绘制朱利亚集最简单的情况。

动力学平面由两个超吸引盆组成

外部 (无穷大是一个超吸引不动点)
内部 (z=0 是超吸引周期点之一)

周期 33 中心

参见通用类别

曼德勃罗集组件的数量

方程

$N(p)=2^{p-1}-s$

其中

N(p) 是给定周期 p 的曼德勃罗集的所有组件的数量。
p 是周期
s 是其各个分量在小于周期 p 的周期每个因数的组件数量之和^[4]

所有值都是正整数。它是序列 A000740

Maxima CAS 中的程序

N(p) := block( 
[a:2^(p-1),f:0], 
for f in divisors(p) do 
   if f < p then a : a - N(f), 
return(a)
);

您可以运行它

for p:1 thru 50 step 1 do display(N(p));

N(1)=1
N(2)=1
N(3)=3
N(4)=6
N(5)=15
N(6)=27
N(7)=63
N(8)=120
N(9)=252
N(10)=495
N(11)=1023
N(12)=2010
N(13)=4095
N(14)=8127
N(15)=16365
N(16)=32640
N(17)=65535
N(18)=130788
N(19)=262143
N(20)=523770
N(21)=1048509
N(22)=2096127
N(23)=4194303
N(24)=8386440
N(25)=16777200
N(26)=33550335
N(27)=67108608
N(28)=134209530
N(29)=268435455
N(30)=536854005
N(31)=1073741823
N(32)=2147450880
N(33)=4294966269
N(34)=8589869055
N(35)=17179869105
N(36)=34359605280
N(37)=68719476735
N(38)=137438691327
N(39)=274877902845
N(40)=549755289480
N(41)=1099511627775
N(42)=2199022198821
N(43)=4398046511103
N(44)=8796090925050
N(45)=17592186027780
N(46)=35184367894527
N(47)=70368744177663
N(48)=140737479934080
N(49)=281474976710592
N(50)=562949936643600

或

sum(N(p),p,1,50);

结果

 1125899873221781

还有 c 程序

/*
gcc c.c -lm -Wall
./a.out
Root point between period 1 component and period 987 component  = c = 0.2500101310666710+0.0000000644946597
Internal angle (c) = 1/987
Internal radius (c) = 1.0000000000000000

*/

#include <stdio.h>
#include <math.h>		// M_PI; needs -lm also
#include <complex.h>






// numer of hyberolic components ( and it's centers ) of Mandelbrot set 
int GiveNumberOfCenters(int period){

	//int s = 0;
	int num = 1;
	
	int f; 
	int fMax = period-1; //sqrt(period); // https://stackoverflow.com/questions/11699324/algorithm-to-find-all-the-exact-divisors-of-a-given-integer
	
	
	
	if (period<1 ) {printf("input error: period should be positive integer \n"); return -1;}
	if (period==1) { return 1;}
		
	num = pow(2, period-1);
	
	// minus sum of number of center of all it's divisors (factors)
	for (f=1; f<= fMax; ++f){
	
		if (period % f==0)
			{num = num - GiveNumberOfCenters(f); }
	
	
	
	}
	
	
		
		
	


	return num ;

}


int ListNumberOfCenters(int period){

	int p=1;
	int pMax = period;
	int num=0;
	
	if (period ==1 || period==2) {printf (" for period %d there is only one component\n", period); return 0;}
	
	for (p=1; p<= pMax; ++p){
		num = GiveNumberOfCenters(p);
		printf (" for period %d there are %d components\n", p, num);
		}
		
	return 0;

}


int main (){

	int period = 15;
	
	
	
	 ListNumberOfCenters(period);

	return 0;

}

另请参见

计算主分子中的双曲分量

如何计算中心？

方法

给定周期 $p\,$ $p\,$ 的所有中心
- 求解 $F_{p}(0,c)=0\,$ $F_{p}(0,c)=0\,$ 或 $G_{p}(0,c)=0\,$ $G_{p}(0,c)=0\,$ ^[5]
  - 以图形形式显示结果
  - 使用数值方法^[6]^[7]
    - “快速而肮脏”的算法：检查是否 $abs(z_{n})<eps$ ，则用颜色 n 给 c 点着色。这里 n 是吸引轨道周期，eps 是吸引点周围圆的半径 = 数值计算的精度。
    - 找到多项式的所有根
      - 迭代细化牛顿法^[8]
      - Ehrlich-Aberth 迭代（Mpsolve^[9]）
    - 原子域
    - Myrberg 方法^[10]^[11]
  - 区间算术和分治策略^[12]
- “搜索是使用网格的四叉树细分进行的，这会在根密度更高的区域提高搜索分辨率。在一个搜索区域内，首先以中等分辨率绘制函数幅度的图，然后使用 640 位定点数学细化局部最小值点。由于乘以非常小的数字，根搜索需要目标结果的两倍精度。随着根的添加，四叉树进一步细分，搜索继续。由于对称性，只搜索上半平面。

在更高的周期，四叉树深度被有意地限制，以减少运行时间，因为根的数量变得难以管理。” James Artherton

给定周期 $p\,$ $p\,$ ，给定点 $c\,$ $c\,$ 附近的一个中心 $c_{p}\,$ $c_{p}\,$
- 求解 $F_{p}(0,c_{p})=0\,$ 的牛顿法的一些步骤

定义中心的方程组

双曲分量的中心是一个参数平面的点 $c\,$ ，对于该点，周期性的 z 循环是超吸引的。它给出了一个定义周期为 $p\,$ 的双曲分量的中心的两个方程组

首先定义周期点 $z_{0}\,$ ，
其次使 $z_{0}\,$ 点超吸引。

${\begin{cases}F_{p}(z_{0},c)=z_{0}\\{\frac {df_{c}^{(p)}(z_{0})}{dz}}=0\end{cases}}$

有关详细信息，请参见定义。

“这些多项式具有整数系数。它们可以通过关于度的递推关系获得。令 Pd 为度数为 d 的多项式。我们有” ^[13]

  $P_{0}(z)=1$

  $P_{d+1}(z)=zP_{d}^{2}+1$

求解方程组

当临界点 $z_{cr}=0\,$ 在此循环中时，可以求解第二个方程

$z_{0}=z_{cr}=0\,$

为了解决系统，将 $z_{0}\,$ 代入第一个方程。

定义中心的方程

得到方程：

$F_{p}(0,c)=0\,$

分量的中心是上述方程的根。

因为 $z=0\,$ ，我们可以从这些方程中去掉 z

对于周期 1：z^2+c=z 和 z=0 ，因此 $c=0$
对于周期 2： (z^2+c)^2 +c =z 和 z=0 ，因此 $c^{2}+c=0$
对于周期 3： ((z^2+c)^2 +c)^2 +c = z 和 z=0 ，因此 $(c^{2}+c)^{2}+c=0$

以下是计算上述函数的 Maxima 函数：

P[n]:=if n=1 then c else P[n-1]^2+c;

简化后的定义中心方程

上述方程的根集包含周期为 p 及其因数的中心的集合。这是因为

$F_{p}(0,c)=\prod _{m|p}G_{m}(0,c)\,$

其中

$G\,$ = 格利森周期为 m 的多项式 ^[14] = $F_{p}\,$ 的不可约因式 ^[15] ^[16] ^[1]
使用迭代乘法的大写 Pi 符号
$m|p\,$ 在这里的意思是：对于所有 $p\,$ 的因数（从 1 到 p）。参见因数表。

例如：

$F_{1}=c=G_{1}\,$

$F_{2}=c*(c+1)=G_{1}*G_{2}\,$

$F_{3}=c*(c^{3}+2*c^{2}+c+1)=G_{1}*G_{3}\,$

$F_{4}=c*(c+1)*(c^{6}+3*c^{5}+3*c^{4}+3*c^{3}+2*c^{2}+1)=G_{1}*G_{2}*G_{4}\,$

因此，可以使用以下方法找到不可约多项式

$G_{p}(0,c)={\frac {F_{p}(0,c)}{\prod _{m|p,m<p}G_{m}(0,c)}}\,$

这是用于计算 $G\,$ 的 Maxima 函数

GiveG[p]:=
block(
[f:divisors(p),
t:1], /* t is temporary variable = product of Gn for (divisors of p) other than p */
f:delete(p,f), /* delete p from list of divisors */
if p=1
then return(Fn(p,0,c)),
for i in f do t:t*GiveG[i],
g: Fn(p,0,c)/t,
return(ratsimp(g))
)$

以下是周期为 1 到 10 的 $F\,$ 和 $G\,$ 的阶数表，以及计算这些函数的系数所需的精度（只要使用未展开的形式 $(\ldots (c^{2}+c)^{2}+c\ldots )^{2}+c$ 进行运算，根就可以用更低的精度来计算）。

周期 $p\,$	$\deg(F_{p}(0,c))\,$	$\deg(G_{p}(0,c))\,$	$fpprec\,$
1	1	1	16
2	2	1	16
3	4	3	16
4	8	6	16
5	16	15	16
6	32	27	16
7	64	63	32
8	128	120	64
9	256	252	128
10	512	495	300
11	1024	1023	600

这里

fpprec 是 bigfloat 数值的算术运算中有效的小数位数^[17]

以下是最大系数表。

周期 $p\,$	$log10(gc_{p})\,$	$gc_{p}\,$
1	0	1
2	0	1
3	1	2
4	1	3
5	3	116
6	4	5892
7	11	17999433372
8	21	106250977507865123996
9	44	16793767022807396063419059642469257036652437
10	86	86283684091087378792197424215901018542592662739248420412325877158692469321558575676264
11	179	307954157519416310480198744646044941023074672212201592336666825190665221680585174032224052483643672228833833882969097257874618885560816937172481027939807172551469349507844122611544

精度可以估算为大于最大系数 $log_{2}(gc)\,$ 的二进制形式的大小

 $fpprec>log_{2}(gc)\,$

周期 $p\,$	$\deg(G_{p}(0,c))\,$	$\log _{2}(gc_{p})\,$	$fpprec\,$
1	1	1	16
2	1	1	16
3	3	2	16
4	6	2	16
5	15	7	16
6	27	13	16
7	63	34	32
8	120	67	64
9	252	144	128
10	495	286	300
11	1023	596	600

以下是 Maxima 代码，用于计算 $gc\,$

period_Max:11;
/* ----------------- definitions -------------------------------*/
/* basic function  = monic and centered complex quadratic polynomial 
http://en.wikipedia.org/wiki/Complex_quadratic_polynomial
*/
f(z,c):=z*z+c $
/* iterated function */
fn(n, z, c) :=
 if n=1 	then f(z,c)
 else f(fn(n-1, z, c),c) $
/* roots of Fn are periodic point of  fn function */
Fn(n,z,c):=fn(n, z, c)-z $
/* gives irreducible divisors of polynomial Fn[p,z=0,c] */
GiveG[p]:=
block(
[f:divisors(p),t:1],
g,
f:delete(p,f),
if p=1 
then return(Fn(p,0,c)),
for i in f do t:t*GiveG[i],
g: Fn(p,0,c)/t,  
return(ratsimp(g))
 )$
/* degree of function */
GiveDegree(_function,_var):=hipow(expand(_function),_var);  
log10(x) := log(x)/log(10);
/* ------------------------------*/
file_output_append:true; /* to append new string to file */
grind:true;  
for period:1 thru period_Max step 1 do
(
g[period]:GiveG[period], /* function g */
d[period]:GiveDegree(g[period],c), /* degree of function g */
cf[period]:makelist(coeff(g[period],c,d[period]-i),i,0,d[period]), /* list of coefficients */
cf_max[period]:apply(max, cf[period]), /* max coefficient */
disp(cf_max[period]," ",ceiling(log10(cf_max[period]))),
s:concat("period:",string(period),"  cf_max:",cf_max[period]),
stringout("max_coeff.txt",s)/* save output to file as Maxima expressions */
);

另请参见

math.stackexchange 问题：不可约因子 - 曼德布洛特多项式

寻找中心的图形方法

所有这些方法都显示了周期为 n 及其因子的中心。

颜色与 zn 的大小成正比

颜色与 zn 的大小成正比^[18]
参数平面扫描点 c，使得临界点的轨道消失 ^[19]
YouTube 视频：曼德布洛特振荡 ^[20]

颜色显示 zn 落在哪个象限

这是曼德布洛特集合外部的径向第 n 重分解（与 LSM/M 的第 n 重分解进行比较）

使用 4 种颜色，因为有 4 个象限

re(z_n) > 0 且 im(z_n) > 0（第一象限）
re(z_n) < 0 且 im(z_n) > 0（第二象限）
re(z_n) < 0 且 im(z_n) < 0（第三象限）
re(z_n) > 0 且 im(z_n) < 0（第四象限）。

“... 当四种颜色在某个地方相遇时，你就有了 q_n(c) 的零点，即周期为 n 的因子的中心。此外，浅色或深色表示 c 是否在 M 中。”（Wolf Jung）

以下是用 C 代码计算点 c = cx + cy * i 的 8 位颜色的片段

/* initial value of orbit = critical point Z= 0 */
                       Zx=0.0;
                       Zy=0.0;
                       Zx2=Zx*Zx;
                       Zy2=Zy*Zy;
                       /* the same number of iterations for all points !!! */
                       for (Iteration=0;Iteration<IterationMax; Iteration++)
                       {
                           Zy=2*Zx*Zy + Cy;
                           Zx=Zx2-Zy2 +Cx;
                           Zx2=Zx*Zx;
                           Zy2=Zy*Zy;
                       };
                       /* compute  pixel color (8 bit = 1 byte) */
                       if ((Zx>0) && (Zy>0)) color[0]=0;   /* 1-st quadrant */
                       if ((Zx<0) && (Zy>0)) color[0]=10; /* 2-nd quadrant */
                       if ((Zx<0) && (Zy<0)) color[0]=20; /* 3-rd quadrant */
                       if ((Zx>0) && (Zy<0)) color[0]=30; /* 4-th quadrant */

还可以从 Wolf Jung 的程序 Mandel 源代码的 mndlbrot.cpp 文件中的 mndlbrot 类中的 quadrantqn 函数中找到 cpp 代码 ^[21]

与标准循环的不同之处是

没有溢出测试（例如，对于 abs(zn) 的最大值）。这意味着每个点都有相同的迭代次数。对于较大的 n，可能会发生溢出。此外，不能使用测试 Iteration==IterationMax
迭代限制相对较小（例如，从 IterationMax = 8 开始）

另请参阅 Sage 中的代码 ^[22]

使用高斯波函数法 ^[23]

牛顿法

牛顿法寻找一个根

中心列表

James Artherton 发现的周期为 1-20 的中心（总和 = 1 045 239）。使用自定义软件，采用 320 位定点数学（大约 92 位十进制数）和网格的四叉树细分来提高根密度较高的区域的搜索分辨率。这可以适应天线尖端的非常精细的分辨率。 ^[24]
Robert P Munafo 发现的最大岛屿 ^[25] ^[26]
通过追踪外部射线找到的所有周期不超过 16 的岛屿的数据库（周期、岛屿性、带角度的内部地址、下外部角度分子、分母、上分子、分母、方向、大小、中心实部、虚部） ^[27]
Jay R. Hill 计算的组件的真实中心 ^[28]
中心列表的集合 ^[29]
岛屿
madore.org/~david math/mandpoints.dat

列表

c = 0.339410819995598 -0.050668285162643 i 周期 = 11，扭曲的小曼德布洛特
c = 0.325589509550660 -0.038047880934756 i 周期 = 12，扭曲的小曼德布洛特
c = 0.314559489984000 -0.029273969079000 i 周期 = 13，扭曲的小曼德布洛特
c = 0.272149607027528 +0.005401159465460 i 周期 = 22，扭曲的小曼德布洛特 ^[30]

85 0.355534 -0.337292
134 -1.7492046334590113301 -2.8684660234660531403e-04 
268 -1.7492046334594190961 -2.8684660260955536656e-04 
272 3.060959246472445584e-01 2.374276727158354376e-02 
204 3.06095924633099439e-01 2.3742767284688944e-02 
204 3.06095924629285095e-01 2.37427672645622342e-02 
136 3.06095924643046857e-01 2.374276727237906e-02 
68 3.06095924629536442e-01 2.37427672749394494e-02 
134 -1.7492046334590114 -2.8684660234659111e-04 
267 -1.1822493706369402373694114253571e-01 6.497492977188836930425943612155e-01
268 -1.182402951560276787014475129283e-01 6.4949165134945441813936036487738e-01
18 -0.814158841137593 0.189802029306573 
32 0.2925755 -0.0149977
52 -0.22817920780250860271129306628202459167994 1.11515676722969926888221122588497247465766

程序

MPSolve

mandelbrot-solver

这是一个适用于 Ubuntu 的软件包。一个基于 MPSolve 的曼德布洛特多项式求解器。一个多项式求根器，可以确定具有保证的包含半径的任意精度近似值。它支持利用多项式结构，如稀疏性，并允许多项式隐式定义或以某些非标准基底定义。该二进制文件（作为 MPSolve 软件包中自定义多项式实现的示例提供）利用了此多项式族体的特殊结构来开发一个高效的求解器，该求解器表现出实验性的 O(n^2) 复杂度。

此软件包包含 mpsolve 对曼德布洛特多项式的专门化。

示例

mandelbrot-solver
Usage: mandelbrot-solver n [starting_file]

Parameters: 
 - n is the level of the Mandelbrot polynomial to solve

 - starting_file is an optional file with the approximations that shall be 
                 use as starting points.
a@zelman:~$ mandelbrot-solver 2
 (-1.22561166876654e-01,  7.44861766619744e-01)
 (-1.75487766624669e+00,  0.00000000000000e+00)
 (-1.22561166876654e-01, -7.44861766619744e-01)
a@zelman:~$ mandelbrot-solver 1
 (-1.00000000000000e+00,  0.00000000000000e+00)
a@zelman:~$ mandelbrot-solver 3
 ( 2.82271390766914e-01,  5.30060617578525e-01)
 (-1.56520166833755e-01,  1.03224710892283e+00)
 (-1.94079980652948e+00,  1.26679600613393e-16)
 (-1.00000000000000e+00,  0.00000000000000e+00)
 (-1.31070264133683e+00, -7.22823506473338e-18)
 (-1.56520166833755e-01, -1.03224710892283e+00)
 ( 2.82271390766914e-01, -5.30060617578525e-01)

似乎

 level = period -1

然后结果就是中心。

注释

^ 人们推测这些多项式（用于中心）是不可约的，但这尚未得到证明 - Wolf Jung

参考文献

Eunice Y. S. Chan 的“曼德尔布罗特类多项式求解方法的比较”。2016 年

[1] 核 - 来自 Robert Munafo 的曼德布洛特集合词汇表和百科全书，(c) 1987-2015。

[2] 复杂动力学中的手术，Carsten Lunde Petersen 在线论文

[3] 牛顿法在实践中：有效地找到一百万次多项式的所有根，由 DIERK SCHLEICHER 和 ROBIN STOLL

[4] 维基百科中的除数

[5] FF 主题：曼德布洛特多项式根挑战（阅读 5233 次）

[6] Piers Lawrence 和 Rob Corless：曼德布洛特多项式和矩阵。 1; 048; 575 周期为 21 的根

[7] IAP 海报 2013

[8] Robin Stoll 和 Dierk Schleicher 的“牛顿法在实践中 II：迭代细化牛顿法和寻找某些非常大次数的多项式的所有根的近似最优复杂度”

[9] 维基百科：MPSolve

[10] G Alvarez、M Romera、G Pastor 和 F Montoya 的“通过确定曼德尔布罗特集合的双曲分量中心”。《混沌、孤子与分形》1998 年，第 9 卷第 12 期，第 1997-2005 页

[11] Myrberg，P J，“二次实多项式的迭代 III”，芬兰科学院学报，A. I 编号 348，1964。

[12] 可以计算 100 万个根的 evaldraw 脚本，由 knighty 创建

[13] LEONARDO ROBOL 的“多项式和长期方程的求根算法”

[14] thoverflow 问题：曼德尔布罗特集合的 Gleason 周期 n 多项式的判别式

[15] 维基百科上的不可约因子

[16] V Dolotin 和 A Morozow：关于基本域的形状，或者为什么曼德尔布罗特集合由几乎理想的圆组成？

[17] xima 手册：fpprec

[18] J.D. Long 的“R 中的曼德尔布罗特集合”

[19] Fractal Stream 帮助文件

[20] youtube 视频：曼德尔布罗特振荡

[21] Wolf Jung 的跨平台 cpp 程序 Mandel

[22] Christopher Olah 的“逃逸时间分形的形成”

[23] _3__The_Dark Exploding the Dark Heart of Chaos by Chris King March-Dec 2009

[24] James Artherton 的“分形龙”

[25] Robert P. Munafo 的“最大岛屿”，2010 年 10 月 21 日

[26] rob：largest-islands.txt

[27] Claude 通过跟踪外部射线找到的所有岛屿的数据库，这些岛屿的周期最多为 16

[28] Jay R. Hill 计算的分量的实数中心

[29] set 中心

[30] th.stackexchange 问题：迷你曼德尔布罗特是完全的副本吗？

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]

[9]

[10]

[11]

[12]

[13]

[14]

[15]

[16]

[1]

[17]

[18]

[19]

[20]

[21]

[22]

[23]

[24]

[25]

[26]

[27]

[28]

[29]

[30]