跳轉到內容

分形/複平面迭代/拋物線

來自華夏公益教科書

"大多數用於計算 Julia 集的程式在底層動力學為雙曲線時工作良好,但在拋物線情況下會遇到指數級減速。" (Mark Braverman)[1]

換句話說,這意味著使用標準/樸素演算法可能需要幾天時間才能製作出拋物線 Julia 集的優質影像。

這裡有兩個問題

  • 緩慢(= 懶惰)的區域性動力學(在拋物線不動點附近)
  • 某些部分非常薄(使用標準平面掃描很難找到)

動態平面

[編輯 | 編輯原始碼]
在復二次多項式的情況下,離散動態。外部 = 白色,內部 = 灰色,未知 = 紅色;
影像如何隨著不同的最大迭代次數而變化

動態平面 = 復 z 平面

  • Julia 集 是一個公共邊界:
  • Fatou 集
    • Julia 集的外部 = 無窮大的吸引盆地:
    • Julia 集的內部 = 有限的拋物線不動點 p 的吸引盆地:
      • 直接盆地 = 包含拋物線不動點 p 在其邊界上的分量的總和;p 的直接拋物線盆地是拋物線盆地的週期性連通分量的並集。
        • 吸引 Lea-Fatou 花 = n 個吸引花瓣的總和 = 2*n 個吸引萼片的總和
          • 花瓣 = 花的一部分。每個花瓣包含兩個萼片的一部分。
          • 萼片 (設 1 是第一象限中的一條不變曲線,L 1 是 1 ∪ {0} 包圍的區域,稱為萼片。)[2]

另請參閱

  • 填充的Julia集[3]

關鍵詞

[編輯 | 編輯原始碼]
  • 拋物線棋盤或棋盤
  • 拋物線內爆
  • Fatou 座標
  • 夏威夷耳環
  • Gevrey 漸近展開
    • (7.1) 在原點的 Écalle-Voronin 不變數具有 Gevrey- 1/2 漸近展開[5]
  • 胚芽[6][7][8]
    • 函式的胚芽:函式的泰勒展開
  • 重數[9]
  • Julia-Lavaurs 集
  • Leau-Fatou 花定理:[10] 排斥或吸引花。花由花瓣組成
    • Leau-Fatou 花
  • 拋物線線性化定理
  • 喇叭對映
  • Blaschke 積
  • Inou 和 Shishikura 的近拋物線重整化
  • 復多項式向量場[11]
  • 數字
    • “一個正整數 ν,即 **拋物線退化**,具有以下性質:存在 νq 個吸引花瓣和 νq 個排斥花瓣,它們圍繞不動點交替迴圈。”[12]
    • 組合旋轉數
  • 函式 f 在拋物線不動點處的龐加萊線性化器[13]
  • “拋物線鉛筆。這是一個圓族,它們都有一個共同點,因此彼此相切,內部或外部。”[14]

Leau-Fatou 花定理

[編輯 | 編輯原始碼]

Leau-Fatou 花定理指出,如果函式 具有泰勒展開[15]

則複數 在單位圓上 描述單位向量(方向)

  • 如果 ,則 f 在原點處的排斥向量(從原點到 v)
  • 如果 ,則吸引向量(從 v 到原點)

存在 n 個等距的吸引方向,由 n 個等距的排斥方向隔開。整數 n+1 稱為不動點的重數


函式展開:z+z^5

/* Maxima CAS */
display2d:false;
taylor(z+z^5, z,0,5);
                                
z + z^5 

所以

  • a = 1
  • n = 4

計算吸引方向

/* Maxima CAS */
solve (4*v^4 = 1, v);

[v = %i/sqrt(2),v = -1/sqrt(2),v = -%i/sqrt(2),v = 1/sqrt(2)]

函式 m*z+z^2

(%i7) taylor(m*z+z^2, z,0,5);
(%o7) m*z+z^2

所以

  • a = 1
  • n = 1
(%i9) solve (v = 1, v);
(%o9) [v = 1]

Ecalle 柱面

[編輯 | 編輯原始碼]

Ecalle 柱面[16] 或 Ecalle-Voronin 柱面(由 Jean Ecalle[17][18])[19]

“... 在將 z 和 f(z) 等同(如果 z 和 f(z) 都屬於 P)的等價關係下,花瓣 P 的商。這個商流形稱為 Ecalle 柱面,它與無限柱面 C/Z 構形同構。”[20]

攪拌器動力學

[編輯 | 編輯原始碼]


物理模型:使用打蛋器時蛋糕的行為。

數學模型:具有 2 箇中心(二階退化點)的二維向量場[21][22]

該場圍繞中心旋轉,但似乎沒有發散。

也許對拋物線動力學的更好描述是夏威夷耳環

拋物線胚芽

[編輯 | 編輯原始碼]

胚芽:[23][24][25]

  • z+z^2
  • z+z^3
  • z+z^{k+1}
  • z+a_{k+1}z^{k+1}
  • z+a_{k+1}z^{k+1}
  • "胚芽 透過全純共軛到其線性部分 " (Sylvain Bonnot)[26]

向量場的胚芽

喇叭對映

[編輯 | 編輯原始碼]

"喇叭對映 h = Φ ◦ Ψ,其中 Φ 是 Φattr 的簡寫,Ψ 是 Ψrep 的簡寫(擴充套件的 Fatou 座標和引數化)。"[27]

花瓣

  • "花瓣是由 不變的 Jordan 域。" R PEREZ-MARCO[28]
  • 花瓣是一個陷阱,它捕獲任何趨向於拋物線點的軌道
  • 花瓣是花的一部分
萼片和花瓣
內部角為 1/1 的拋物線萼片

定義

所有 花瓣 的總和形成一朵花[30],中心位於拋物線週期點。[31]

花椰菜

[編輯 | 編輯原始碼]

花椰菜或西蘭花[32]

  • 空(其內部為空)對於 Mandelbrot 集之外的 c。Julia 集是完全不連通的。
  • 對於 Mandelbrot 集邊界上的 c=1/4,填充的花椰菜。Julia 集是 Jordan 曲線(準圓)。

請注意

  • 影像大小因 z 平面不同而異。
  • 使用了不同的演算法,因此顏色很難比較。

花椰菜的分岔

[編輯 | 編輯原始碼]

Julia 集如何沿著實軸變化(從 c=0 到 c=1/4,然後更遠)

擾動函式 由複數

當新增一個大於0的ε(沿實軸向正無窮大移動)時,會出現一個拋物型不動點的分叉

  • 吸引不動點 (ε<0)
  • 一個拋物型不動點 (ε = 0)
  • 一個拋物型不動點分裂成兩個共軛排斥不動點 (ε > 0)

“如果我們用ε<0稍微擾動,那麼拋物型不動點將在實軸上分裂成兩個實不動點(一個吸引,一個排斥)。”

參見

  • Wolf Jung 的 Mandel 程式中的演示 2,第 9 頁

拋物線內爆

[edit | edit source]

拋物型內爆

  • 在引數平面中
    • 點 c 從分量的內部,穿過邊界,到曼德爾布羅特集的外部
    • 核 (c=0),沿內部射線 0,拋物點 (c= 0.25),沿外部射線 0
  • 在動力學平面中
    • 連線的朱利亞集(有內部)內爆,並斷開連線(沒有內部)
    • 不動點從內部移動到朱利亞集(拋物型)
    • 一個盆地(內部)消失
沿逃逸路徑 0 的動力學演化(拋物型內爆)
引數 c 位置 朱利亞集 內部 動力學型別 臨界點 不動點
c = 0 中心,內部 連線的 存在 超吸引 被吸引到 alfa 不動點 固定臨界點等於 alfa 不動點,alfa 是超吸引的,beta 是排斥的
0<c<1/4 內部射線 0,內部 連線的 存在 吸引 被吸引到 alfa 不動點 alfa 是吸引的,beta 是排斥的
c = 1/4 尖點,邊界 連線的 存在 拋物型 被吸引到 alfa 不動點 alfa 不動點等於 beta 不動點,兩者都是拋物型的
c>1/4 外部射線 0,外部 斷開連線的 消失 - 排斥到無窮大 兩個有限不動點都是排斥的



YouTube 上的影片[33]

向量場

[edit | edit source]
  • 二維向量場及其

奇點

[edit | edit source]

奇點型別

  • 中心型別:“在這種情況下,可以找到奇點的鄰域,其中所有積分曲線都是封閉的,彼此巢狀,並在其內部包含奇點”[34]
  • 非中心型別:奇點的鄰域由幾個曲線扇區組成:[35]

“曲線扇區定義為由半徑任意小的圓 C 和兩個流線 S 和 S! 圍成的區域,這兩個流線都收斂到奇點。然後考慮透過開扇區 g 的流線,以區分三種可能的曲線扇區型別。”

動力學

[edit | edit source]

動力學

  • 全域性
  • 區域性


區域性動力學

  • 在朱利亞集的外部
  • 在朱利亞集上
  • 靠近拋物型不動點(在朱利亞集內部)


另請參閱

靠近拋物型不動點

[edit | edit source]
靠近拋物型不動點並在朱利亞集內部的軌道

為什麼分析 f^p 而不是 f?

f 靠近拋物型不動點的正向軌道是複合的。它由兩種運動組成

  • 圍繞不動點
  • 向不動點移動/遠離不動點

如何計算拋物型 c 值

[edit | edit source]

拋物型引數的型別

  • 根點
  • 尖點
曼德爾布羅特集(引數平面)的週期 1 分量的拋物點
n 內部角(旋轉數)t = 1/n 根點 c = 拋物型引數 落在根點 c 上的引數射線的兩個外部角(1/(2^n+1); 2/(2^n+1) 不動點 落在不動點
1 1/1 0.25 (0/1 ; 1/1) 0.5 (0/1 = 1/1)
2 1/2 -0.75 (1/3; 2/3) -0.5 (1/3; 2/3)
3 1/3 0.64951905283833*%i-0.125 (1/7; 2/7) 0.43301270189222*%i-0.25 (1/7; 2/7; 3/7)
4 1/4 0.5*%i+0.25 (1/15; 2/15) 0.5*%i (1/15; 2/15; 4/15; 8/15)
5 1/5 0.32858194507446*%i+0.35676274578121 (1/31; 2/31) 0.47552825814758*%i+0.15450849718747 (1/31; 2/31; 4/31; 8/31; 16/31)
6 1/6 0.21650635094611*%i+0.375 (1/63; 2/63) 0.43301270189222*%i+0.25 (1/63; 2/63; 4/63; 8/63; 16/63; 32/63)
7 1/7 0.14718376318856*%i+0.36737513441845 (1/127; 2/127) 0.39091574123401*%i+0.31174490092937 (1/127; 2/127, 4/127; 8/127; 16/127; 32/127, 64/127)
8 1/8 0.10355339059327*%i+0.35355339059327 0.35355339059327*%i+0.35355339059327
9 1/9 0.075191866590218*%i+0.33961017714276 0.32139380484327*%i+0.38302222155949
10 1/10 0.056128497072448*%i+0.32725424859374 0.29389262614624*%i+0.40450849718747

對於內部角 n/p 拋物週期 p 迴圈,它包含一個具有 p 重數的 z 點[36] 和一個乘子 = 1.0。這個點 z 等於不動點

週期 1

[edit | edit source]

可以輕鬆計算邊界點 c

對於給定的內部角(旋轉數)t,使用 Wolf Jung 編寫的這段cpp 程式碼[37]

 t *= (2*PI); // from turns to radians
 cx = 0.5*cos(t) - 0.25*cos(2*t); 
 cy = 0.5*sin(t) - 0.25*sin(2*t);

或這段 Maxima CAS 程式碼

 
/* conformal map from circle to cardioid ( boundary
 of period 1 component of Mandelbrot set */
F(w):=w/2-w*w/4;

/* 
circle D={w:abs(w)=1 } where w=l(t,r) 
t is angle in turns ; 1 turn = 360 degree = 2*Pi radians 
r is a radius 
*/
ToCircle(t,r):=r*%e^(%i*t*2*%pi);

GiveC(angle,radius):=
(
 [w],
 /* point of  unit circle   w:l(internalAngle,internalRadius); */
 w:ToCircle(angle,radius),  /* point of circle */
 float(rectform(F(w)))    /* point on boundary of period 1 component of Mandelbrot set */
)$

compile(all)$

/* ---------- global constants & var ---------------------------*/
Numerator: 1;
DenominatorMax: 10;
InternalRadius: 1;

/* --------- main -------------- */
for Denominator:1 thru DenominatorMax step 1 do
(
 InternalAngle: Numerator/Denominator,
 c: GiveC(InternalAngle,InternalRadius),
 display(Denominator),
 display(c),
  /* compute fixed point */
 alfa:float(rectform((1-sqrt(1-4*c))/2)), /* alfa fixed point */
 display(alfa)
)$

週期 2

[edit | edit source]
// cpp code by W Jung http://www.mndynamics.com
 t *= (2*PI);  // from turns to radians
 cx = 0.25*cos(t) - 1.0;
 cy = 0.25*sin(t);

週期 1-6

[edit | edit source]
/*

batch file for Maxima CAS 
computing bifurcation points for period 1-6

 Formulae for cycles in the Mandelbrot set II
Stephenson, John; Ridgway, Douglas T.
Physica A, Volume 190, Issue 1-2, p. 104-116.
*/

kill(all);
remvalue(all);

start:elapsed_run_time ();

/* ------------ functions ----------------------*/

/* exponential for of complex number with angle in turns */
 /* "exponential form prevents allroots from working", code by Robert P. Munafo */

GivePoint(Radius,t):=rectform(ev(Radius*%e^(%i*t*2*%pi), numer))$ /* gives point of unit circle for angle t in turns */

GiveCirclePoint(t):=rectform(ev(%e^(%i*t*2*%pi), numer))$ /* gives point of unit circle for angle t in turns Radius = 1 */

/* gives a list of iMax points of unit circle */
GiveCirclePoints(iMax):=block(
 [circle_angles,CirclePoints],
 CirclePoints:[],
 circle_angles:makelist(i/iMax,i,0,iMax),
 for t in circle_angles do CirclePoints:cons(GivePoint(1,t),CirclePoints),
 return(CirclePoints) /* multipliers */
)$

/* http://commons.wikimedia.org/wiki/File:Mandelbrot_set_Components.jpg 
Boundary equation  b_n(c,P)=0 
    defines relations between hyperbolic components and unit circle for given period n ,
    allows computation of exact coordinates of hyperbolic componenets.

b_n(w,c), is boundary polynomial (implicit function of 2 variables).

*/

GiveBoundaryEq(P,n):=
block(
 if n=1 then return(c + P^2 - P),
 if n=2 then return(- c + P - 1),
 if n=3 then return(c^3 + 2*c^2 - (P-1)*c + (P-1)^2),
 if n=4 then return(c^6 + 3*c^5 + (P+3)* c^4 + (P+3)* c^3 - (P+2)*(P-1)*c^2 - (P-1)^3),
 if n=5 then return(c^15 + 8*c^14 + 28*c^13 + (P + 60)*c^12 + (7*P + 94)*c^11 + 
  (3*P^2 + 20*P + 116)*c^10 + (11*P^2 + 33*P + 114)*c^9 + (6*P^2 + 40*P + 94)*c^8 + 
  (2*P^3 - 20*P^2 + 37*P + 69)*c^7 + (3*P - 11)*(3*P^2 - 3*P - 4)*c^6 + (P - 1)*(3*P^3 + 20*P^2 - 33*P - 26)*c^5 +
  (3*P^2 + 27*P + 14)*(P - 1)^2*c^4 - (6*P + 5)*(P - 1)^3*c^3 + (P + 2)*(P - 1)^4*c^2 - c*(P - 1)^5  + (P - 1)^6),
if n=6 then return (c^27+
13*c^26+
78*c^25+
(293 - P)*c^24+
(792 - 10*P)*c^23+
(1672 - 41*P)*c^22+
(2892 - 84*P - 4*P^2)*c^21+
(4219 - 60*P - 30*P^2)*c^20+
(5313 + 155*P - 80*P^2)*c^19+
(5892 + 642*P - 57*P^2 + 4*P^3)*c^18+
(5843 + 1347*P + 195*P^2 + 22*P^3)*c^17+
(5258 + 2036*P + 734*P^2 + 22*P^3)*c^16+
(4346 + 2455*P + 1441*P^2 - 112*P^3 + 6*P^4)*c^15 + 
(3310 + 2522*P + 1941*P^2 - 441*P^3 + 20*P^4)*c^14 + 
(2331 + 2272*P + 1881*P^2 - 853*P^3 - 15*P^4)*c^13 + 
(1525 + 1842*P + 1344*P^2 - 1157*P^3 - 124*P^4 - 6*P^5)*c^12 + 
(927 + 1385*P + 570*P^2 - 1143*P^3 - 189*P^4 - 14*P^5)*c^11 + 
(536 + 923*P - 126*P^2 - 774*P^3 - 186*P^4 + 11*P^5)*c^10 + 
(298 + 834*P + 367*P^2 + 45*P^3 - 4*P^4 + 4*P^5)*(1-P)*c^9 + 
(155 + 445*P - 148*P^2 - 109*P^3 + 103*P^4 + 2*P^5)*(1-P)*c^8 + 
2*(38 + 142*P - 37*P^2 - 62*P^3 + 17*P^4)*(1-P)^2*c^7 + 
(35 + 166*P + 18*P^2 - 75*P^3 - 4*P^4)*((1-P)^3)*c^6 + 
(17 + 94*P + 62*P^2 + 2*P^3)*((1-P)^4)*c^5 + 
(7 + 34*P + 8*P^2)*((1-P)^5)*c^4 + 
(3 + 10*P + P^2)*((1-P)^6)*c^3 + 
(1 + P)*((1-P)^7)*c^2 +
-c*((1-P)^8) + (1-P)^9)
)$

/* gives a list of points c on boundaries on all components for give period */
GiveBoundaryPoints(period,Circle_Points):=block(
 [Boundary,P,eq,roots],
  Boundary:[],
 for m in Circle_Points do (/* map from reference plane to parameter plane */
  P:m/2^period,
  eq:GiveBoundaryEq(P,period), /* Boundary equation  b_n(c,P)=0  */
  roots:bfallroots(%i*eq),
  roots:map(rhs,roots),
  for root in roots do Boundary:cons(root,Boundary)),
  return(Boundary)
)$

/* divide llist of roots to 3 sublists = 3  components */
/* ---- boundaries of period 3 components 
period:3$
Boundary3Left:[]$
Boundary3Up:[]$
Boundary3Down:[]$

Radius:1;

 for m in CirclePoints do (
  P:m/2^period,
  eq:GiveBoundaryEq(P,period),
  roots:bfallroots(%i*eq),
  roots:map(rhs,roots),
  for root in roots do 
     (
       if realpart(root)<-1  then Boundary3Left:cons(root,Boundary3Left),
       if (realpart(root)>-1 and imagpart(root)>0.5) 
            then Boundary3Up:cons(root,Boundary3Up),
       if (realpart(root)>-1 and imagpart(root)<0.5) 
            then Boundary3Down:cons(root,Boundary3Down)
               
     )

)$
--------- */

/* gives a list of parabolic points for given: period and internal angle */
GiveParabolicPoints(period,t):=block
(
 [m,ParabolicPoints,P,eq,roots],
 m: GiveCirclePoint(t), /* root of unit circle, Radius=1, angle t=0 */
 ParabolicPoints:[],
 /* map from reference plane to parameter plane */
 P:m/2^period,
 eq:GiveBoundaryEq(P,period), /* Boundary equation  b_n(c,P)=0  */
 roots:bfallroots(%i*eq),
 roots:map(rhs,roots),
 for root in roots do ParabolicPoints:cons(float(root),ParabolicPoints),
 return(ParabolicPoints)

)$

compile(all)$

/* ------------- constant values ----------------------*/

fpprec:16;

/* ------------unit circle on a w-plane -----------------------------------------*/
a:GiveParabolicPoints(6,1/3);
a$

/*
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>



/* 
   c functions using complex type numbers
   computes c from  component  of Mandelbrot set */
complex double Give_c( int Period,  int n, int d , double InternalRadius )
{
  complex double c; 
  complex double w;  // point of reference plane  where image of the component is a unit disk 
 // alfa = ax +ay*i = (1-sqrt(d))/2 ; // result 
  double t; // InternalAngleInTurns
  
  t  = (double) n/d; 
  t = t * M_PI * 2.0; // from turns to radians
  
  w = InternalRadius*cexp(I*t); // map to the unit disk 
  
  switch ( Period ) // of component 
    {
    case 1: // main cardioid = only one period 1 component
      c = w/2 - w*w/4; // https://wikibook.tw/wiki/Fractals/Iterations_in_the_complex_plane/Mandelbrot_set/boundary#Solving_system_of_equation_for_period_1
      break;
    case 2: // only one period 2 component 
      c = (w-4)/4 ; // https://wikibook.tw/wiki/Fractals/Iterations_in_the_complex_plane/Mandelbrot_set/boundary#Solving_system_of_equation_for_period_2
      break;
      // period > 2 
    default: 
      printf("higher periods : to do, use newton method \n");
      printf("for each q = Period of the Child component  there are 2^(q-1) roots \n");
      c = 10000.0; // bad value 
       
      break; }
      
      
  return c;
}



void PrintAndDescribe_c( int period,  int n, int d , double InternalRadius ){


	complex double c = Give_c(period, n, d, InternalRadius);
	
	printf("Root point between period %d component and period %d component  = c = %.16f%+.16f*I\t",period, d, creal(c), cimag(c));
	printf("Internal angle (c) = %d/%d\n",n, d);
	//printf("Internal radius (c) = %.16f\n",InternalRadius);



}

/*
https://stackoverflow.com/questions/19738919/gcd-function-for-c
The GCD function uses Euclid's Algorithm. 
It computes A mod B, then swaps A and B with an XOR swap.
*/

int gcd(int a, int b)
{
    int temp;
    while (b != 0)
    {
        temp = a % b;

        a = b;
        b = temp;
    }
    return a;
}





int main (){


	int period = 1;
	double InternalRadius = 1.0;
	// internal angle in turns as a ratio = p/q
	int n = 1;
	int d = 987;
	
	
	
	
	
	// n/d = local angle in turns
 	
   		for (n = 1; n < d; ++n ){
     			if (gcd(n,d)==1 )// irreducible fraction
       				{ PrintAndDescribe_c(period, n,d,InternalRadius); }
       				}
       		
	
	



	return 0;

}

如何繪製拋物線朱利亞集

[編輯 | 編輯原始碼]

填充朱利亞集內部的所有點都趨向於一個週期軌道(或不動點)。此點位於朱利亞集中,並且是弱吸引的。[38] 您可以分析拋物線不動點附近的行為。可以使用 臨界軌道 來完成此操作。

這裡有兩種情況:簡單和困難。

如果拋物線不動點附近的朱利亞集類似於 n 臂星形(未扭曲),則您可以簡單地檢查 zn 相對於不動點的幅角。例如,請參見 z+z^5。這是一個簡單的例子。

在困難的情況下,朱利亞集圍繞不動點扭曲。

根點的朱利亞集在拓撲上與子週期中心的朱利亞集相同,但

  • 中心(核)朱利亞集非常容易繪製(超吸引盆 = 非常快的動態,因為臨界點也是週期點)
  • 而根朱利亞集(拋物線)很難繪製(拋物線盆和懶惰動態)

示例

  • t = 1/2
    • 根點的朱利亞集 = 胖巴塞利卡朱利亞集:c = -3/4 = - 0.75
    • 週期 2 中心的朱利亞集 = (細長)巴塞利卡朱利亞集:c = -1
  • t = 1/3
    • 根點的朱利亞集 = 胖 杜阿迪兔子:c = -0.125000000000000 +0.649519052838329i
    • 週期 3 中心的朱利亞集 = (細長)杜阿迪兔子 朱利亞集:c = -0.122561166876654 +0.744861766619744i 週期 = 3

從外部估計

[編輯 | 編輯原始碼]

逃逸時間

[編輯 | 編輯原始碼]

描述

長迭代或大步法

[編輯 | 編輯原始碼]

描述

步驟

  • 選擇函式 ,其在原點(z=0)處有一個固定拋物線點
  • 選擇內部角度
  • 計算 ,它是函式 的較高迭代次數的近似值,z 接近於零,使用以零為中心的冪級數(泰勒級數 = 麥克勞林級數)
  • 實驗性地找出使用冪級數的項數以及使用特定 的環帶


每個 將被用於一個環帶

其中 K 為固定值

具有在原點處具有乘子 的無差異不動點的復二次多項式的 Lambda 形式[41]

其中

  • 定點的乘數
  • 內角是一個有理數和真分數


選擇

  • 因此
  • 冪級數的30項
  • 環k=4,5,...,10的近似函式,對於z的較大值(在環外),使用預設函式f^n
  • 函式的delta等於10^-5


(* code by Professor: Mark McClure from https://marksmath.org/classes/Spring2019ComplexDynamics/ *)
n = 7;
f[z_] = Exp[2 Pi*I/n] z + z^2;

Remove[F];
F[0][z_]    = N[Normal[Series[f[z], {z, 0, 30}]]];
Do[F[0][z_] = Chop[N[Normal[Series[F[0][f[z]], {z, 0, 30}]]], 10^-5],  {n - 1}];
Do[F[k][z_] = Chop[N[Normal[Series[F[k - 1][F[k - 1][z]], {z, 0, 30}]]], 10^-5], {k, 1, 10}]

(* define and compile function FF *)
FF = With[{
    n = n,
    f4 = Function[z, Evaluate[F[4][z]]],
    f5 = Function[z, Evaluate[F[5][z]]],
    f6 = Function[z, Evaluate[F[6][z]]],
    f7 = Function[z, Evaluate[F[7][z]]],
    f8 = Function[z, Evaluate[F[8][z]]],
    f9 = Function[z, Evaluate[F[9][z]]],
    f10 = Function[z, Evaluate[F[10][z]]]
    },
   Compile[{{z, _Complex}},
    Which[
     Abs[z] > 1/2^3, 
     Nest[Function[zz, N[Exp[2 Pi*I/n]] zz + zz^2], z, n],
     Abs[z] <= 1/2^9, f10[z],
     Abs[z] <= 1/2^8, f9[z],
     Abs[z] <= 1/2^7, f8[z],
     Abs[z] <= 1/2^6, f7[z],
     Abs[z] <= 1/2^5, f6[z],
     Abs[z] <= 1/2^4, f5[z],
     Abs[z] <= 1/2^3, f4[z],
     True, 0]]];
     
(* iterate 1000 times and then see what happens *)
iterate = With[{FF = FF, n = n},
   Compile[{{z0, _Complex}},
    Module[{z, i},
     z = z0;
     i = 0;
     While[1/2^9 < Abs[z] <= 2 && i++ < 1000 n,
      z = FF[z]];
     z],
    RuntimeOptions -> "Speed", CompilationTarget -> "C"]];

(* now compute some iteration data *)
data = Monitor[
    Table[iterate[x + I*y], {y, Im[center] + 1.2, Im[center], -0.0025},
     {x, Re[center] - 1.2, Re[center] + 1.2, 0.0025}], 
    y]; 


(* use some symmetry to cut computation time in half *)
center = First[Select[z /. NSolve[f[z] == 0, z], Im[#] < 0 &]]/2 (* center = -0.311745 - 0.390916*I *)
data = Join[data, Reverse[Rest[Reverse /@ data]]];


(* plot it  *)
kernel = {
   {1, 1, 1},
   {1, -8, 1},
   {1, 1, 1}
   };
   
 (*   
classifyArg = Compile[
   {{z, _Complex}, {z0, _Complex}, {v, _Complex}, {n, _Integer}},
   Module[{check, check2},
    check = n (Arg[(z0 - z)/v] + Pi)/(2 Pi);
    check2 = Ceiling[check];
    If[check == check2, 0, check2]]];   
   
   
   
classified = Map[classify, data, {2}];

convolvedData = ListConvolve[kernel, classified];

ArrayPlot[Sign[Abs[convolvedData]]]


Wolfram 語言的數學函式( 

  • Series[f,{x,x0,n}] 生成關於點的f的冪級數展開,階數為,其中n是一個明確的整數
  • N[expr] 給出expr的數值
  • Chop[expr,delta] 將絕對值小於delta的數字替換為0
  • Normal[expr] 將expr從各種特殊形式轉換為正常表示式。
  • Do[expr,n] 評估expr n次
  • Do[expr,{i,imin,imax}] 評估expr,從i=imin開始
  • With[{x=x0,y=y0,…},expr] 指定expr中所有x、y等的出現應該被替換為x0、y0等。
  • Compile[{{x1,t1},…},expr] 假設xi的型別與ti匹配。
  • Which[test1,value1,test2,value2,…] 依次評估每個testi,返回第一個產生True的valuei的值
  • Nest[f,expr,n] 給出對expr應用n次f的表示式。

動態射線

[edit | edit source]
內角為1/15的拋物線朱利亞集 - 使用外部射線作為朱利亞集在alfa不動點附近的近似值

可以使用落在拋物線不動點的週期性動態射線來尋找外部的狹窄部分。

讓我們檢查外部半徑=4的週期性射線上的點需要多少次反向迭代才能到達距離拋物線不動點0.003的位置

週期 反向迭代 時間
1 340 0m0.021s
2 55 573 0m5.517s
3 8 084 815 13m13.800s
4 1 059 839 105 1724m28.990s

可以使用外部射線的點z的幅角及其到alfa不動點的距離(參見影像中的程式碼)。它適用於週期數不超過15(也許更多……)。

從內部估計

[edit | edit source]

朱利亞集是填充的朱利亞集Kc的邊界。

  • 找到Kc內部的點
  • 使用邊緣檢測找到Kc內部的邊界

如果內部的元件彼此非常接近,則使用以下方法找到元件:[42]

color = LastIteration % period

對於父元件和子元件之間的拋物線元件:[43]

periodOfChild = denominator*periodOfParent  
color = iLastIteration % periodOfChild

其中分母是曼德勃羅集的父元件的內角的分母。

角度

[edit | edit source]

"如果zn的迭代趨於一個固定的拋物點,那麼初始種子z0根據zn−z0的幅角進行分類,分類由花定理提供"(Mark McClure[44])

吸引時間

[edit | edit source]
各種型別的動力學

填充的朱利亞集的內部由元件組成。所有元件都是前週期的,其中一些是週期的(吸引力的直接盆地)。

換句話說

  • 一次迭代將z移動到另一個元件(以及整個元件移動到另一個元件)
  • 所有元件的點都具有相同的吸引時間(到達吸引子周圍的目標集所需的迭代次數)

可以使用它來對元件進行著色。因為在拋物線情況下吸引子是弱的(弱吸引),所以某些點需要很多次迭代才能到達它。

 // i = number of iteration
 // iPeriodChild = period of child component of Mandelbrot set ( parabolic c value is a root point between parant and child component
 /* distance from z to Alpha  */
 Zxt=Zx-dAlfaX;
 Zyt=Zy-dAlfaY;
 d2=Zxt*Zxt +Zyt*Zyt;
 // interior: check if fall into internal target set (circle around alfa fixed point)  
 if (d2<dMaxDistance2Alfa2) return  iColorsOfInterior[i % iPeriodChild];

這裡有一些示例值

 iWidth  = 1001 // width of image in pixels
 PixelWidth  = 0.003996  
 AR  = 0.003996 // Radius around attractor
 denominator  = 1 ; Cx  = 0.250000000000000; Cy  = 0.000000000000000 ax  = 0.500000000000000; ay  = 0.000000000000000   
 denominator  = 2 ; Cx  = -0.750000000000000; Cy  = 0.000000000000000 ax  = -0.500000000000000; ay  = 0.000000000000000   
 denominator  = 3 ; Cx  = -0.125000000000000; Cy  = 0.649519052838329 ax  = -0.250000000000000; ay  = 0.433012701892219  
 denominator  = 4 ; Cx  = 0.250000000000000; Cy  = 0.500000000000000 ax  = 0.000000000000000; ay  = 0.500000000000000   
 denominator  = 5 ; Cx  = 0.356762745781211; Cy  = 0.328581945074458 ax  = 0.154508497187474; ay  = 0.475528258147577   
 denominator  = 6 ; Cx  = 0.375000000000000; Cy  = 0.216506350946110 ax  = 0.250000000000000; ay  = 0.433012701892219    
 
 denominator  = 1 ;   i =               243.000000 
 denominator  = 2 ;   i =            31 171.000000 
 denominator  = 3 ;   i =         3 400 099.000000 
 denominator  = 4 ;   i =       333 293 206.000000 
 denominator  = 5 ;   i =    29 519 565 177.000000 
 denominator  = 6 ;   i = 2 384 557 783 634.000000 

其中

C = Cx + Cy*i 
a = ax + ay*i // fixed point alpha
i // number of iterations after which critical point z=0.0 reaches disc around fixed point alpha with radius AR
denominator of internal angle (in turns)
internal angle = 1/denominator

請注意,吸引時間i與分母成正比。

各種分母的吸引時間

現在您看到了弱吸引意味著什麼。

你可以

  • 使用蠻力方法(吸引半徑<畫素大小;迭代最大值足夠大,讓內部的所有點都到達目標集;長時間或快速計算機
  • 找到更好的方法(:-)),如果您覺得時間太長

內部距離估計

[edit | edit source]

陷阱

[edit | edit source]
陷阱 = 目標集

從內部和外部估計

[編輯 | 編輯原始碼]

Julia 集是填充的 Julia 集和無窮大吸引域的共同邊界。

  • 找到 Kc 內部/分量的點
  • 找到逃逸點
  • 使用 Sobel 濾波器找到邊界點

它適用於分母不超過 4。

排斥點的逆迭代

[編輯 | 編輯原始碼]

逆迭代 阿爾法不動點。它只對切割點(外部射線落入的地方)有效。其他點仍然沒有被擊中。

外部示例


有關其他多項式對映,請參閱此處

另請參閱

[編輯 | 編輯原始碼]

參考文獻

[編輯 | 編輯原始碼]
  1. Mark Braverman:關於拋物 Julia 集的有效計算
  2. 關於拋物線的動態穩定擾動的註記 作者:河原智己
  3. 維基百科中的填充的 Julia 集
  4. Barile, Margherita. "夏威夷耳環。" 來自 MathWorld - Wolfram 網頁資源,由 Eric W. Weisstein 建立。https://mathworld.tw/HawaiianEarring.html
  5. Augustin Fruchard, Reinhard Sch¨afke. 複合漸近展開和差分方程。應用資訊學與數學非洲研究評論,INRIA,2015,20,pp.63-93。<hal-01320625>
  6. 維基百科:胚芽 (數學)
  7. 微分同胚的不動點、向量場的奇點及其軌道的 epsilon 鄰域 作者:Maja Resman
  8. 展開拋物不動點的解析微分同胚胚芽族模空間 作者:Colin Christopher,Christiane Rousseau
  9. 維基百科:重數 (數學)
  10. 表面同胚的動力學 Leau-Fatou 花定理和穩定流形定理的拓撲版本 作者:Le Roux,F
  11. C 中復多項式向量場的動力學 作者:Kealey Dias
  12. 退化拋物二次有理對映的極限 作者:XAVIER BUFF、JEAN ECALLE 和 ADAM EPSTEIN
  13. 更高維的龐加萊線性化 作者:Alastair Fletcher
  14. 圓的鉛筆 作者:James King
  15. math.stackexchange 問題:拋物臨界軌道的形狀是什麼
  16. 全純不變數理論。國家論文,奧賽,1974 年 3 月
  17. 法語維基百科中的讓·埃卡勒
  18. 讓·埃卡勒主頁
  19. 盧卡斯·蓋耶 - 透過一致化實現的標準形式 (2016 年 10 月 28 日)。這是使用一致化定理證明吸引、排斥和拋物不動點區域性標準形式的草案,在課堂上分發。最終它將被納入講義。
  20. 對映 作者:Luna Lomonaco
  21. 通用拋物微分同胚展開的解析分類模 作者:P. Mardesic、R. Roussarie¤ 和 C. Rousseau
  22. mathoverflow 問題:函式方程 ffxxfx2
  23. 維基百科中的胚芽
  24. 通用拋物微分同胚展開的解析分類模 作者:P. Mardesic、R. Roussarie 和 C. Rousseau
  25. 展開拋物不動點的解析微分同胚胚芽族模空間 作者:Colin Christopher,Christiane Rousseau
  26. Mathoverflow:在固定點附近函式的無窮小分類,直到共軛
  27. 單奇異全純對映的近拋物重整化 作者:Arnaud Cheritat
  28. Ricardo Pérez-Marco. "不動點和圓周對映." Acta Math. 179 (2) 243 - 294, 1997. https://doi.org/10.1007/BF02392745
  29. 關於拋物線的動態穩定擾動的註記 由 川平知己
  30. 拋物不動點:Leau-Fatou 花 由 Davide Legacci 2021 年 3 月 18 日
  31. 維基百科:玫瑰 (拓撲)
  32. 花椰菜 在 MuEncy 由 Robert Munafo
  33. 圓圈內爆成火焰 - 影片由 sinflrobot
  34. 二維向量場的拓撲簡化方法 由 Xavier Tricoche Gerik Scheuermann Hans Hagen
  35. 數學百科:常微分方程理論中的扇形
  36. 維基百科:數學中的重數
  37. Mandel:用於實數和複數動力學的軟體 由 Wolf Jung
  38. 不動點處的區域性動力學 由 Evgeny Demidov
  39. 拋物型朱利亞集是多項式時間可計算的 由 Mark Braverman
  40. Mark McClure 春季課程 2019:複數動力學,參見程式碼“拋物點附近的迭代”(Mathematica 筆記本)
  41. Michael Yampolsky, Saeed Zakeri:配對西格爾二次多項式。
  42. 不動點和週期軌道 由 Evgeny Demidov
  43. 繪製拋物型朱利亞集的 C 程式原始碼
  44. Stack Exchange 問題:拋物型臨界軌道的形狀是什麼?
  45. PlanetMath:聖馬可分形
  46. 維基百科:杜阿迪兔子
  47. PlanetMath:聖馬可分形
  48. 圖片:非標準拋物線 由 Cheritat
  49. Maxima CAS 中的拋物型情況下的朱利亞集
華夏公益教科書