ダウンロード

ジュリア集合(Julia set)

ジュリア集合とは
集合の計算手順
終了条件
色の設定
エクセルでの限界
シート構成
サンプル画像2010
VBAコード
サンプルファイルのダウンロード
リンクリスト
VBAコードを 全て展開 折り畳む




ジュリア集合(Julia set)とは


計算する漸化式はマンデルブロ集合と同じですが、定数部を固定し、
座標の初期値を変化させます。

ジュリア集合
ジュリア集合


0.原図
1.左のピンク枠内を拡大




集合の計算手順

漸化式

漸化式

で定義される複素数列において C を固定し、Zn の初期値を変化させ n → ∞ の極限で無限大に発散しないという条件を満たす複素数 Z 全体が作る集合、これがジュリア集合となります。

図中の黒い部分は収束域、 周囲の色は無限大に発散する速さを表します。


ここでは、ジュリア集合の zn を点 (xn, yn) に、c を点 (a, b) にそれぞれ置き代えて、

漸化式

として計算します。

   Zn = Xn + i Yn とすると
   Zn+1 = Zn2 + C
       = Xn2 - Yn2 + i 2XnYn + a + ib

   Zn+1 = Xn+1 + iYn+1 

   なので、実数部、虚数部は、それぞれ
   Xn+1 = Xn2 - Yn2 + a
   Yn+1 = 2XnYn + b


 下図は、a=-0.3 b=-0.63 とした場合の集合

座標


マンデルブロ集合とジュリア集合の比較

マンデルブロ集合 ジュリア集合
漸化式 Z=Z2+C Z=Z2+C
Cの値 複素平面上を変化させる 固定(外部から設定)
Zの初期値Z0 0 + 0i に固定 複素平面上を変化させる
集合の特徴 フラクタル性、一個の塊 フラクタル性、一個の塊の場合とバラバラの場合がある


終了条件


下記の場合に計算を打ち切る

 ・計算回数が所定の回数に達する
 ・|Z| > 2 になった場合(発散が証明済み)
 ・Zn+1 = Zn となる場合(収束)




発散回数を求める「CalDivCount」関数は、マンデルブロで用いたものと同 じ。

計算座標はDouble型(倍精度浮動小数点型)としているので有効数字は15桁となり、1兆倍以上の拡大可能。

使用したVBAコード

Function CalDivCout(x0 As Double, y0 As Double, a As Double, b As Double, _
                    Optional MaxCnt As Long = 100) As Long
'
' 発散するまでの反復回数を計算
'
    Dim Cnt As Long
    Dim x As Double, y As Double, xn As Double, yn As Double
    Dim S As Double, S0 As Double
        
    xn = x0: yn = y0
    Cnt = 0
    Do
    
        x = xn * xn - yn * yn + a
        y = 2 * xn * yn + b
       
        S = x * x + y * y

        
        xn = x: yn = y
        Cnt = Cnt + 1
        
        If S = S0 Then Cnt = MaxCnt
        
        S0 = S
    
    Loop While ((S < 4#) And (Cnt < MaxCnt))

    CalDivCout = Cnt

End Function



VBAを使って色で表すと、不思議で複雑な美しい画像が現れて来ます。

色の設定

・白黒の場合

 計算回数に応じて、255〜1までの輝度を割り当てる。
 発散する場合は、0(黒)

・カラーの場合

 収束色(黒)と区別できるように、RGB(100,100,100)を基礎色とし、
 それより大きい値〜白までの値を計算回数に応じて割り当てる。
 割り当てる値は、最大計算回数に対する比とする。

 具体的には、

  色 = ColorBase + ( &HFFFFFF-ColorBase ) * v / MaxLimit
  としている。 発散する場合は、0(黒)

  計算回数:v
  基礎色:ColorBase=RGB(100,100,100)
  最大計算回数:MaxLimit

 
            '-------------------------------------
            ' 着色処理
            '-------------------------------------
            Color = 0   '黒
            If v < max_cnt Then '反復計算回数最大値より小さい場合
            
                If blColor Then '着色要?
                    Dim ColorBase As Long
                    ColorBase = RGB(100, 100, 100)
                    Color = ColorBase + (&HFFFFFF - ColorBase) * v / MAX_CNT_LIMIT   'カラー
                Else
                    v = 255 * (1 - v / max_cnt)             '白黒
                    Color = RGB(v, v, v)
                End If
            
            End If

    

エクセルでの限界


解像度

エクセルワークシートを使用しているため、エクセル2003以前では解像度の限界があります。
列数は255なので、ワークシート列を最大限近くまで使用しています。
ここでは、列数と行数を切りよく4:3の比で、240列x180行とします。
(列幅:0.08、行高:0.75)

列幅を大きくするとギザギザが目立ってしまいますので、列幅を最小幅に縮小しています。

かなり描画エリアが小さくなりますが、鑑賞に堪えられる画質にはなっています。

エクセル2010で列幅を拡張したものを同梱しています。
エクセル2007以降では列数は6万5千以上なので、理論的には十分大きな描画エリアをとれますが、
計算時間の壁がありますので、ここでは控えめに、720列x580行としています。

拡大倍率

倍精度(Double)型( -1.79769313486231570E+308 〜 -4.94065645841246544E-324)
を使用しているため、1兆倍くらいまでは拡大可能。
(限界に近づくとモザイク状になり、最後は単色表示)



シート構成



フォーム構成  A:定数値−入力可
   a - 実数部 b - 虚数部

@:初期状態(左画面)に戻す

A選択域拡大
  拡大したい部分を選択すると有効になる

B選択域細密
  解像度が 「x 1」以外の場合、
  選択域のみを「x 1」に細密表示

  画像表示に時間がかかるエクセル2010ブックの場合、
  解像度を「x 1」以外で粗く表示しておき、
  必要な部分だけを細密化すると効率的に画像探索ができる。


C再計算
  同じ描画範囲で再計算する。
  粗い画像で確認後、解像度を上げて再描画したい場合、
  または、計算回数変更後の再描画時など。

D戻る
  ひとつ前の画像に戻る
  拡大画像をやり直したい場合など。

ESTOP
  計算停止

F解像度
  x 1:最大解像度(描画時間大)
  x 7:最小解像度(描画時間小)

G最大計算回数
  拡大を繰り返していくと、発散までの計算回数が多い領域が出てきます。
  最大計算回数を超えると収束と見なして黒く塗り潰すため、
  輪郭がぼやけてくることがあります。
  不自然な黒の領域が現れたら、最大計算回数を上げて下さい。

HColor
  オンの場合、カラー表示。オフは白黒画像。

Iセル座標、前図からの倍率、原図からの倍率

J計算座標



サンプル画像


それぞれ、 a + b i

0.32 + 0.043 i
サンプル画像

0.27334 + 0.00742 i
サンプル画像

-0.15652 + 1.03225 i
サンプル画像

-0.11 + 0.67 i
サンプル画像

-0.12 + 0.74 i
サンプル画像

-0.39054 -0.58679 i
サンプル画像

0.11031 - 0.67037 i
サンプル画像

0.11031 - 0.67037 i
サンプル画像

-0.74543 + 0.11301 i
サンプル画像

VBAコード

変数定義部 VBAコード

Option Explicit

'-------------------------------------
' ジュリア集合
'
' 複素数の数列 Zn+1=Zn^2+C
'
'  が発散する回数を複素平面での色として求める
'
'  Z = X+iY  C=a+ib (i 虚数単位)として
'
'  実数 Xn+1 = Xn^2 - Yn^2 + a
'  虚数 Yn+1 = 2*Xn*Yn + b
'
'  の数列を計算する
'
'--------------------------------------

Dim Rng As Range        'セル範囲オブジェクト
Dim Selected As Range   '選択されている範囲

Dim 中止 As Boolean     '中止フラグ

'---------------------
' 初期値
'---------------------

'計算範囲

Dim C0x As Double ' < Julia >
Dim C0y As Double ' < Julia >

Const INIT_AMAX = 2#     '実数軸最大
Const INIT_AMIN = -2#   '  最小
Const INIT_BMAX = 1.5   '虚数軸最大
Const INIT_BMIN = -1.5   '  最小

'セル範囲
Dim CMax As Integer        '列最大
Dim CMin As Integer        ' 最小
Dim RMax As Integer        '行最大
Dim RMin As Integer        ' 最小
Dim RC_ratio As Double     '比

Dim xmax As Double, xmin As Double
Dim ymax As Double, ymin As Double

Dim cs As Integer, ce As Integer
Dim rs As Integer, re As Integer


Dim max_cnt As Long         '指定反復計算最大回数
Const MAX_CNT_LIMIT = 1000  '同上限

'------------------------------------
' 履歴
'------------------------------------
Private Type TYPE_HISTORY
    xmin As Double
    xmax As Double
    ymin As Double
    ymax As Double
    max_cnt As Long
    rs As Integer
    re As Integer
    cs As Integer
    ce As Integer
End Type

Dim History(1) As TYPE_HISTORY  '履歴


イベント処理 VBAコード

'------------------------------------
' イベント処理
'------------------------------------

Private Sub cmdStop_Click()     '中止
    中止 = True
End Sub


Private Sub cmdInit_Click() ' 初期状態に戻す
    Call Init                       '初期化
    Rng.Interior.Pattern = xlNone   'クリア
    Set Selected = Nothing              '未選択状態
    Range("XY").Value = GetXY       '計算範囲
    Range("MAG").Value = GetMag     '拡大率表示
    Call Go                         '開始
End Sub


Private Sub cmdReCal_Click()  '再計算
    Call Go
End Sub


Private Sub cmdBack_Click()
'
' 戻る
'
    Set Rng = Range("DATA")                 '描画域
    With History(1)
        xmin = .xmin: xmax = .xmax: ymin = .ymin: ymax = .ymax
        max_cnt = .max_cnt
        rs = .rs: re = .re: cs = .cs: ce = .ce
    End With
    cmbCalCount.Value = max_cnt
    cmdBack.Enabled = False

    
    Call Go '作成開始
    
    Range("DATA").Range("A1").Select 'セル選択

End Sub


Private Sub cmdMagnify_Click()
'
' 選択範囲を拡大
'
    Dim xs As Double, xe As Double
    Dim ys As Double, ye As Double
    Dim Lc As Double, Lr As Double
    Dim Lx As Double, Ly As Double

    Set Rng = Range("DATA")                 '描画域
    Rng.Interior.Pattern = xlNone   'クリア
    
    '-------------------------
    ' 選択範囲チェック
    '-------------------------
    Dim blErr As Boolean
    blErr = False
    If Selected Is Nothing Then '未選択
        blErr = True
    ElseIf Selected.Rows.Count + Selected.Columns.Count < 4 Then    '範囲不足
        blErr = True
    End If
    If blErr Then
        MsgBox "拡大範囲を選択して下さい!"
        Exit Sub
    End If
    '-------------------------
   
    If xmax = 0 And xmin = 0 Then Call Init    '範囲選択−未

    Lc = CMax - CMin + 1
    Lr = RMax - RMin + 1
    Lx = xmax - xmin
    Ly = ymax - ymin
    
    cs = Selected.Column - CMin             '選択セル 列開始
    ce = cs + Selected.Columns.Count - 1    '         列終わり
    
    rs = Selected.Row - RMin                '選択セル 行開始
    re = rs + (ce - cs) * RC_ratio          '     列幅にあわせる
    
    '------------------------------------
    ' 履歴
    '------------------------------------
    With History(1)
        .xmin = xmin: .xmax = xmax: .ymin = ymin: .ymax = ymax
        .max_cnt = max_cnt
        .rs = rs: .re = re: .cs = cs: .ce = ce
    End With
    '------------------------------------

    xmax = xmin + Lx * (ce - 1) / (Lc - 1)
    xmin = xmin + Lx * (cs - 1) / (Lc - 1)
    
    ymin = ymin + Ly * (rs - 1) / (Lr - 1)
    ymax = ymin + (xmax - xmin) * (INIT_BMAX - INIT_BMIN) / (INIT_AMAX - INIT_AMIN)
    

    cmdBack.Enabled = False  '「戻る」ボタン無効

    Call Go '作成開始
    
    cmdBack.Enabled = True  '「戻る」ボタン有効
    
End Sub


Private Sub cmdDetail_Click()
'
' 選択域を細密表示
'
    中止 = False

    With Selected
        Go_Sub .Row, .Row + .Rows.Count - 1, .Column, .Column + .Columns.Count - 1, 1 '解像度に応じて
        .Select
    End With
    
End Sub


Private Sub Worksheet_Activate()
    Call Init   '初期値
End Sub


Private Sub Worksheet_SelectionChange(ByVal Target As Range)    '選択範囲が変更された
    
    Dim Data As Range

    Set Data = Range("DATA")
    Set Selected = Intersect(Target, Data)

    Dim Rng As Range
    Set Rng = Intersect(Target, Data)
    Application.StatusBar = GetMag(Rng) & "    " & _
                            GetSel(Rng) & "     " & _
                            GetXY(Rng)   '計算範囲、拡大率表示

End Sub




内部サブルーチン/関数 VBAコード

'////////////////////////////////////////////
'-----------------------------------------
' サブルーチン/関数
'-----------------------------------------

Private Sub Init()          '初期値設定
    Set Rng = Range("DATA")                 '描画域
    RMin = 1: RMax = Rng.Rows.Count         'セル行−開始/終了
    CMin = 1: CMax = Rng.Columns.Count      '  列−開始/終了
    RC_ratio = (RMax - RMin) / (CMax - CMin) '行列数比
    xmax = INIT_AMAX: xmin = INIT_AMIN      '描画範囲 x軸 開始/終了
    ymax = INIT_BMAX: ymin = INIT_BMIN      '     y軸
End Sub

Private Sub Go(Optional Stp As Integer = 1)

    Set Rng = Range("DATA")                 '描画域
    Range("SELCEL").Value = GetSel(Selected)  '選択範囲
    Range("A2").Select      '描画エリアを未選択状態にする
    
    中止 = False

    Go_Sub RMin, RMax, CMin, CMax, 2 * lstResolution.ListIndex + 1 '解像度に応じて

End Sub

Private Sub Go_Sub(rrmin As Integer, rrmax As Integer, _
                   ccmin As Integer, ccmax As Integer, _
                   Optional Stp As Integer = 1)
'
' 作成開始
'
    Dim v As Double
    
    Dim r As Integer, c As Integer  '行、列
    Dim rrs As Integer, rre As Integer
    Dim ccs As Integer, cce As Integer
    
    Dim x As Double, y As Double

    Dim blColor As Boolean  '着色−描画時に有効


    Dim Color As Long

    
    cmdStop.Enabled = True  '中止ボタンを有効

    If xmax = 0 And xmin = 0 Then Call Init     '計算範囲未設定→初期化
    
    '---------------------------
    ' 可変値
    '---------------------------
    
    C0x = Val(Range("C0x")) ' < Julia >
    C0y = Val(Range("C0y")) ' < Julia >
    
    blColor = CBool(chkColor.Value)     '着色
    max_cnt = Val(cmbCalCount.Value)    '反復計算回数

    Range("MAG0").Value = GetMag0       '拡大率(前回比)
    Range("MAG").Value = GetMag         '拡大率

    Range("XY").Value = GetXY                   '計算範囲
    
    MyWait 100, 中止
    
    '--------------------------------------
    ' 解像度に応じて計算中心位置を決定
    '--------------------------------------
    Dim Span As Integer
    Span = (Stp + 1) / 2 - 1
    rrs = rrmin + Span: rre = rrmax - Span
    ccs = ccmin + Span: cce = ccmax - Span
    
    For r = rrs To rre - Span Step Stp

        For c = ccs To cce - Span Step Stp
        
            x = xmin + (xmax - xmin) * (c - CMin + 1) / (CMax - CMin + 1)
            y = ymin + (ymax - ymin) * (r - RMin + 1) / (RMax - RMin + 1)
            v = CalDivCout(x, y, C0x, C0y, max_cnt)     ' < Julia >

            '-------------------------------------
            ' 着色処理
            '-------------------------------------
            Color = 0   '黒
            If v < max_cnt Then '反復計算回数最大値より小さい場合
            
                If blColor Then '着色要?
                    Dim ColorBase As Long
                    ColorBase = RGB(100, 100, 100)
                    Color = ColorBase + (&HFFFFFF - ColorBase) * v / MAX_CNT_LIMIT   'カラー
                Else
                    v = 255 * (1 - v / max_cnt)             '白黒
                    Color = RGB(v, v, v)
                End If
            
            End If
            
            Rng.Range(Cells(r - Span, c - Span), Cells(r + Span, c + Span)).Interior.Color = Color
            '--------------------------------------

        Next c

        If 中止 Then Exit For
        Application.StatusBar = (r & "/" & RMax)    'ステータスバーに進行状況表示
        
        DoEvents

    Next r
    
    cmdStop.Enabled = False '中止ボタンを無効

End Sub


Private Function GetSel(Optional Target As Range) As String
'
' 選択セル
'
    Dim r0 As Integer, c0 As Integer, r As Integer, c As Integer
    Dim Data As Range
    
    Set Data = Range("DATA")
    
    With Target
        
        If Target Is Nothing Then
            
            cmdMagnify.Enabled = False  '選択域拡大 有効
            cmdDetail.Enabled = False   '   細密
            
            GetSel = RMin & "," & CMin & " - " & RMax & "," & CMax
        Else
            
            cmdMagnify.Enabled = True
            cmdDetail.Enabled = True
            
            r0 = Data.Row - 1: c0 = Data.Column - 1
            r = .Row - r0: c = .Column - c0
            GetSel = r & "," & c & " - " & _
                               (r + .Rows.Count - 1) & "," & (c + .Columns.Count - 1)   '選択セル
        End If
        
    End With
    
End Function


Private Function GetMag0(Optional Target As Range) As String
'
' 拡大率(前回比)
'
    Dim v As Double
    v = 0
    With History(1)
        If (.xmax <> .xmin) And (xmax <> xmin) Then
            v = (.xmax - .xmin) / (xmax - xmin)                '拡大率
        End If
    End With
    GetMag0 = ""
    If v >= 1 Then GetMag0 = Format(v, "x#,###")
    
End Function

Private Function GetMag(Optional Target As Range) As String
'
' 拡大率
'
    Dim XY As Range
    Dim x1 As Double, y1 As Double, x2 As Double, y2 As Double
    
    Set XY = Range("XY")
    
    If Target Is Nothing Then
        x1 = xmin: x2 = xmax: y1 = ymin: y2 = ymax
    Else
        With Target
            Call CellToXY(.Row, .Column, x1, y1)
            Call CellToXY(.Row + .Rows.Count - 1, .Column + .Columns.Count - 1, x2, y2)
        End With
    End If

    If x1 <> x2 Then GetMag = _
            Format((INIT_AMAX - INIT_AMIN) / (x2 - x1), "x#,###")    '拡大率
End Function


Private Function GetXY(Optional Target As Range) As String
'
' 計算範囲
'
    Dim x1 As Double, y1 As Double, x2 As Double, y2 As Double
    
    If Target Is Nothing Then
        x1 = xmin: x2 = xmax: y1 = ymin: y2 = ymax
    Else
        With Target
            Call CellToXY(.Row, .Column, x1, y1)
            Call CellToXY(.Row + .Rows.Count - 1, .Column + .Columns.Count - 1, x2, y2)
        End With
    End If
    
    GetXY = x1 & " , " & y1 & " : " & x2 & " , " & y2            '計算範囲

End Function


Private Function CellToXY(r As Integer, c As Integer, ByRef x As Double, ByRef y As Double)
'
' セル座標から計算座標を求める
'
'  出力:x,y
        If CMax = 0 Then Exit Function
        x = xmin + (xmax - xmin) * (c - CMin) / (CMax - CMin)
        y = ymin + (ymax - ymin) * (r - RMin) / (RMax - RMin)
End Function





×
PageTop