多変数関数のニュートンフラクタル

2024年2月18日 19:00

どうも、108Hassiumです。

以前の記事で、「ニュートンフラクタル」というフラクタル図形を紹介しました。

ニュートンフラクタルは、ニュートン法というアルゴリズムにおいてどの初期値がどの値に収束するかで複素平面を彩色した際に現れるフラクタル図形です。

そしてニュートン法は $f(z)=0$ という形の方程式の解( $f(z)$ の根)を見つけるアルゴリズムで、その内容は $z_{n+1}=z_n-\frac{f(z_n)}{f'(z_n)}$ という数列を計算するというものである、というのがニュートンフラクタルの概要です。

ニュートン法のような「 $f(z)$ の根を求めるアルゴリズム」を「求根アルゴリズム」と呼び、以前の記事では他の求根アルゴリズムに対応するニュートンフラクタルの類似物を紹介したりもしました。

ところで、wikipediaのニュートン法の記事には他の求根アルゴリズムだけではなく、「多変数関数に対するニュートン法」についても書かれています。

というわけで、この記事では多変数関数のニュートン法とそれが生成するフラクタル図形を紹介したいと思います。

※多変数関数と言いつつ、今回紹介するものは2変数のものに限定します。

定義

まず、多変数関数のニュートン法で扱う関数は「 $m$ 次元のベクトルを変数として $m$ 次元のベクトルを返す関数」です。

今までの記事では多変数関数は

$f(x,y)=(x^3-3xy^2-x^2-y^2-1,3x^2y-y^3)$

・・・という書き方をしていましたが、この記事では以下のように書くことにします。

$f(\mathbf{x})=\begin{pmatrix}x^3-3xy^2-x^2-y^2-1\\3x^2y-y^3\end{pmatrix}$

ただし、 $\mathbf x=\begin{pmatrix}x\\y\end{pmatrix}$ です。

要するに関数 $f$ は $\mathbf x=\begin{pmatrix}x\\y\end{pmatrix}$ という2次元のベクトルを受け取って $\begin{pmatrix}x^3-3xy^2-x^2-y^2-1\\3x^2y-y^3\end{pmatrix}$ という2次元ベクトルを返す関数である、ということです。

さて、多変数関数のニュートン法の漸化式は以下の通りです。

$\mathbf{x_{n+1}}=\mathbf{x_n}-\partial f(\mathbf{x_n})^{-1}f(\mathbf{x_n})$

ここで、

$\partial f(\mathbf{x_n})$ は $f(\mathbf{x_n})$ のヤコビ行列
$\partial f(\mathbf{x_n})^{-1}$ は $\partial f(\mathbf{x_n})$ の逆行列
$\partial f(\mathbf{x_n})^{-1}f(\mathbf{x_n})$ は $\partial f(\mathbf{x_n})^{-1}$ と $f(\mathbf{x_n})$ の積

・・・です。

私の記事は基本的に高校レベルまでの知識があれば読めるように書いているのですが、今回は予備知識が多すぎて解説しきれないので、分からない部分は各自で調べるか理解を諦めて画像だけ流し見していってください。

さて、先程の式のままでは扱いにくいので少し式変形をします。

まず、 $f(\mathbf x)=\begin{pmatrix}X\\Y\end{pmatrix}$ 、 $\partial f(\mathbf{x_n})=\begin{pmatrix}X_x,X_y\\Y_x,Y_y\end{pmatrix}$ と略記します。

これを使って $\mathbf{x}-\partial f(\mathbf{x})^{-1}f(\mathbf{x})$ を計算すると、以下のようになります。

$\mathbf{x}-\partial f(\mathbf{x})^{-1}f(\mathbf{x})\\=\mathbf{x}-\begin{pmatrix}X_x,X_y\\Y_x,Y_y\end{pmatrix}^{-1}\begin{pmatrix}X\\Y\end{pmatrix}\\=\mathbf{x}-\frac{1}{X_xY_y-X_yY_x}\begin{pmatrix}Y_y,-X_y\\-Y_x,X_x\end{pmatrix}\begin{pmatrix}X\\Y\end{pmatrix}\\=\mathbf{x}-\frac{1}{X_xY_y-X_yY_x}\begin{pmatrix}XY_y-YX_y\\YX_x-XY_x\end{pmatrix}\\=\begin{pmatrix}x-\frac{XY_y-YX_y}{X_xY_y-X_yY_x}\\y-\frac{YX_x-XY_x}{X_xY_y-X_yY_x}\end{pmatrix}$

この式を利用することで、多変数関数のニュートン法における漸化式を作ることができるようになります。

例えば $f(\mathbf{x})=\begin{pmatrix}x^3-3xy^2-x^2-y^2-1\\3x^2y-y^3\end{pmatrix}$ の場合、 $\begin{pmatrix}X_x,X_y\\Y_x,Y_y\end{pmatrix}=\begin{pmatrix}3x^2-3y^2-2x,-6xy-2y\\6xy,3x^2-3y^2\end{pmatrix}$ になるので漸化式は以下のようになります。

$\begin{cases}x_{n+1}=x_n-\frac{(x_n^3-3x_ny_n^2-x_n^2-y_n^2-1)(3x_n^2-3y_n^2)-(3x_n^2y_n-y_n^3)(-6x_ny_n-2y_n)}{(3x_n^2-3y_n^2-2x_n)(3x_n^2-3y_n^2)-(-6x_ny_n-2y_n)(6x_ny_n)}\\y_{n+1}=y_n-\frac{(3x_n^2y_n-y_n^3)(3x_n^2-3y_n^2-2x_n)-(x_n^3-3x_ny_n^2-x_n^2-y_n^2-1)(6x_ny_n)}{(3x_n^2-3y_n^2-2x_n)(3x_n^2-3y_n^2)-(-6x_ny_n-2y_n)(6x_ny_n)}\end{cases}$

実例

☝(x^3-3xy^2-x^2-y^2-1,3x^2y-y^3)のニュートンフラクタル

☝(x^3-3xy^2-x^2-y^2-0.1,3x^2y-y^3+0.1)のニュートンフラクタル

$\begin{pmatrix}x^3-3xy^2-x^2-y^2+a\\3x^2y-y^3+b\end{pmatrix}$ という関数は、実は以前紹介した $z^3-x^2-y^2+c$ を2変数実関数として表したものです。

$z^3-x^2-y^2+c$ はジュリア集合が3回回転対称になるのが特徴的でしたが、ニュートンフラクタルにも同じ性質が受け継がれるようです。

☝(x^3-3xy^2-x^2-y^2+0.1,3x^2y-y^3+0.03)のニュートンフラクタル

☝(x^3-3xy^2-x^2-y^2+0.15,3x^2y-y^3+0.05)のニュートンフラクタル

$\begin{pmatrix}x^3-3xy^2-x^2-y^2+a\\3x^2y-y^3+b\end{pmatrix}=\begin{pmatrix}0\\0\end{pmatrix}$ の解は $a$ と $b$ の組み合わせによっては9個になることもあるのですが、その場合もちゃんと9種類の解が求まっているようです。