Excel数式でマンデルブロ集合を描く

:::note warn
以下は訳あって削除済みの過去記事の転載です（published at 2023-05-18 12:56）。
:::

---

Excel数式の`MAKEARRAY`関数を使用すると複数のセルに展開される数式を作成できます。今回はこれでマンデルブロ集合を描画してみます。

https://twitter.com/U_1F992/status/1658843380619571200

# マンデルブロ集合

https://ja.wikipedia.org/wiki/%E3%83%9E%E3%83%B3%E3%83%87%E3%83%AB%E3%83%96%E3%83%AD%E9%9B%86%E5%90%88

こういうやつです。ある複素数 $c$ について

```math
\left\{ \begin{array}{l} z_{0} = 0 \\
z_{n+1} = z_{n} ^ 2 + c \end{array} \right.
```

という数列をしばらく計算しても、 $z_{n}$ の絶対値がめちゃ大きくならない（発散しない） $c$ の集合のことです。$c = x + yi$として平面に表し、めちゃ大きくなっちゃったときは $n$ に応じて色を付ければ、上の画像のようになります。高校数学は万年赤点だったのでこれくらいで許してください。

# `MAKEARRAY`関数

https://support.microsoft.com/ja-jp/office/makearray-関数-b80da5ad-b338-4149-a523-5b221da09097

縦横のサイズと各セルの値を決める関数を渡すと、複数セルに展開される数式です。

```
=MAKEARRAY(2, 2, LAMBDA(r, c,
  r & "," & c
))
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1550728/ca295af8-dbd9-28ef-7fc6-69d43f9fa8fe.png)

今回は、ワークシートを複素平面に見立てます。

# Excel数式でループを表現する

Excel数式にはforやwhileにあたる関数はありません。ループは再帰で表現します。`LET`関数でラムダに名前をつけることができますが、関数定義内ではすでに定義された名前しか使えない、つまり自身の名前を使えないことに注意が必要です。

```
=LET(
cnt, 10,
LOOP, LAMBDA(i,
  IF(i < cnt,
    i & "," & LOOP(i + 1),  # 不可
    ""
  )
),

LOOP(0)
)
```

この制約を乗り越えるには、実行時に引数に自分自身を渡すようにします。

```
=LET(
cnt, 10,
LOOP, LAMBDA(SELF,i,
  IF(i < cnt,
    i & "," & SELF(SELF,i + 1),  # SELFなら使ってもよい
    ""
  )
),

LOOP(LOOP, 0)  # SELFには必ずLOOPを渡す
)
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1550728/e04ddde0-d3c8-7f56-5f59-afe683c43073.png)

これで、Excel数式で再帰を表現できます。

# nの値でぼんやり表現する

あとは実装するだけです。まずは発散しなかった場合は空白、発散した場合はnの値を各セルに書き込むことでぼんやりひょうたん型が見えるようにします。Excel数式には無限大を表す定数がないので、倍精度浮動小数点数型で正確に計算できる最大値の9007199254740991あたりを入れておきます。

```
=LET(
INFINITY, 9007199254740990,
limit, 100,
scale, 10,
size, 41,

offset, (size - 1) / 2 / scale,

COUNTUNTILDIVERGENCE, LAMBDA(SELF, n, c, z,
  IF(AND(IMABS(z) < INFINITY, n <> limit),
    SELF(SELF, n + 1, c, IMSUM(IMPRODUCT(z, z), c)),
    IF(n = limit, "", n)
  )
),

MAKEARRAY(size, size, LAMBDA(row, col, LET(
  rp, -offset + (col - 1) / scale,
  ip, offset - (row - 1) / scale,
  
  COUNTUNTILDIVERGENCE(COUNTUNTILDIVERGENCE, 0, COMPLEX(rp, ip), 0)
)))
)
```

あらかじめ表示位置と倍率は変数にしました。41x41セルを使い、0.1ずつ値を変えることで、各軸-2〜2の範囲を表現できます。

心の清い人にはマンデルブロ集合が浮かび上がって見えます。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1550728/588e29b2-fe5e-9eca-e26d-05935823898f.png)

# nの値に応じて色を付ける

Excelには、セルの値に応じて書式を自動で設定する「条件付き書式」という機能があります。一律で条件を書けるように上の数式を少し変更し、発散しなかった場合は0を書き込み、プリセットにあった赤から白に遷移するグラデーションをかけることにします。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1550728/1d911a9e-c94e-4856-85e2-78c2b56304ed.png)

仕上げに、VBAをちょっとだけ使ってセルを正方形に整形します。セルのサイズはRowHeight, ColumnWidthで設定できますが、単位が異なることに注意してください。RowHeightは素直にポイント単位ですが、ColumnWidthは「標準スタイルの1文字の幅」が単位です。ただしWidthプロパティからはポイント単位の幅が得られるので較正します。

```vb
Cells.RowHeight = 0.75 ' これでだいたい1pxになった
' Dim columnWidth As Double
' columnWidth = 0.75
' Dim diff As Double
' diff = 0.001
' Cells.ColumnWidth = columnWidth
' Do While Cells.RowHeight < Cells(1, 1).Width
'     columnWidth = columnWidth - diff
'     Cells.ColumnWidth = columnWidth
' Loop
' Debug.Print columnWidth ' 0.114999999999999
Cells.ColumnWidth = 0.115
```

801x801、-2〜2の範囲を0.005ずつ区切って計算したものがこちらです。2分ほどで完了したので、Excelいじめとしてはまだ余裕がありそうです。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1550728/fea17fdc-ee2c-49b1-2ef4-5a7ffd2081ec.png)

```
=LET(
INFINITY, 9007199254740990,
limit, 256,
scale, 200,
size, 801,

offset, (size - 1) / 2 / scale,

COUNTUNTILDIVERGENCE, LAMBDA(SELF, n, c, z,
  IF(AND(IMABS(z) < INFINITY, n <> limit),
    SELF(SELF, n + 1, c, IMSUM(IMPRODUCT(z, z), c)),
    IF(n = limit, 0, n)
  )
),

MAKEARRAY(size, size, LAMBDA(row, col, LET(
  rp, -offset + (col - 1) / scale,
  ip, offset - (row - 1) / scale,
  
  COUNTUNTILDIVERGENCE(COUNTUNTILDIVERGENCE, 0, COMPLEX(rp, ip), 0)
)))
)
```