gijyutsu-keisan.com

2.二変数のニュートン法

2.1.二変数のニュートン法手順

二つの関数が与えられているとします。
\[ \begin{eqnarray} & f(x,y) & \ =0 \\ & g(x,y) & \ =0 \end{eqnarray} \tag{2.1-1} \]
今、\( f(x,y)=0, g(x,y)=0 \)を同時に満足する点\( x_t,y_t \)を求めます。 これは\( f(x,y)=0 \)と\( g(x,y)=0 \)の交点を求めることに他なりません。

二変数のニュートン法は次のStepを踏みます。

Step1:

真の解\( ( x_t, y_t) \)の近傍点\( ( x^{(0)} \, y^{(0)} \) \)を初期解として、適当に決めます。

Step2:

\( f(x,y)=0, g(x,y)=0 \)をそれぞれ点\( ( x^{(0)} \ , y^{(0)} \ ) \)でテイラー展開して、その一次近似をとります。
\[ \begin{eqnarray} f^{(0)} + \frac{\partial f^{(0)}}{\partial x} \ ( x - x^{(0)} ) + \frac{\partial f^{(0)}}{\partial y} \ ( y - y^{(0)} ) = 0 \\ g^{(0)} + \frac{\partial g^{(0)}}{\partial x} \ ( x - x^{(0)} ) + \frac{\partial g^{(0)}}{\partial y} \ ( y - y^{(0)} ) = 0 \end{eqnarray} \tag{2.1-2} \] \[ ( f^{(0)} = f(x^{(0)} \ , y^{(0)} \ ), g^{(0)} = g(x^{(0)} \ , y^{(0)} \ ) ) \]
この式は行列を使って表現できます。
\[ \boldsymbol{f^{(0)}} = \begin{pmatrix} f^{(0)} \\ g^{(0)} \end{pmatrix} , \ \boldsymbol{x^{(0)}} = \begin{pmatrix} x^{(0)} \\ y^{(0)} \end{pmatrix} , \ \boldsymbol{x} = \begin{pmatrix} x \\ y \end{pmatrix} \] \[ J^{(0)} = \begin{pmatrix} \partial f^{(0)} / \partial x & \ \partial f^{(0)} / \partial y & \\ \partial g^{(0)} / \partial x & \ \partial g^{(0)} / \partial y & \end{pmatrix} \]
特に、偏微分係数で構成される行列\(J\)をヤコビ行列と呼びます。
\[ J(\boldsymbol{x}) = \begin{pmatrix} \partial f / \partial x & \ \partial f / \partial y & \\ \partial g / \partial x & \ \partial g / \partial y & \end{pmatrix} \tag{2.1-3} \]
さて、元に戻って(2.1-2)式行列表現すると、
\[ \boldsymbol{f^{(0)}} + J^{(0)} ( \boldsymbol{x} - \boldsymbol{x^{(0)}} ) = \boldsymbol{o} \\ \\ \leftrightarrow \boldsymbol{x} = \boldsymbol{x^{(0)}} - \left(J^{(0)} \right)^{-1} \boldsymbol{f^{(0)}} \tag{2.1-4} \]
となって、\(J^{(0)}\)の逆行列が求まれば\( x \)が求まります。 逆行列は存在するものとして、求まった解を\( x^{(1)} \)とします。

なお、通常数値計算では\( H^{-1} \)は求めず、 ガウスの消去法等を用いて\( \boldsymbol{x^{(1)}} \)を求めます。

Step3:

\( \boldsymbol{x^{(0)}} \ \)と\( \boldsymbol{x^{(1)}} \ \)の大きさ (ノルム:\( \| \boldsymbol{x^{(0)}} \|, \| \boldsymbol{x^{(1)}} \| \))を求め、 その差分、つまり前Stepとの解の変化が微小値\( \epsilon \)(収束判定値)より小さければ、\( x^{(1)} \)を\( ( x_t, y_t) \)の近似解とします。
\[ | \| \boldsymbol{x^{(1)}} \| - \| \boldsymbol{x^{(0)}} \| | < \epsilon \tag{2.1-5} \] \[ \| \boldsymbol{x} \| = \sqrt{ x^2 + y^2 } \]
上式を満たさない場合は、解の変化が\( \epsilon \)より小さくなるまで、\( \boldsymbol{x^{(2)}},\boldsymbol{x^{(3)}},\cdots \)を求めて解の差分を評価する操作、Step2、Step3、・・・を繰り返します 。
以上の手順から、二変数のニュートン法は次の計算を繰り返すだけです。 ただし、一変数とは違い、行列計算が必要になるのが厄介です。
\[ \begin{eqnarray} & & J(\boldsymbol{x}) = \begin{pmatrix} \partial f / \partial x & \ \partial f / \partial y & \\ \partial g / \partial x & \ \partial g / \partial y & \end{pmatrix} \\ \\ & & J^{(i)} = J(\boldsymbol{x^{(i)}}) , \ \boldsymbol{f^{(i)}} = \begin{pmatrix} f^{(i)} \\ g^{(i)} \end{pmatrix} , \ \boldsymbol{x^{(i)}} = \begin{pmatrix} x^{(i)} \\ y^{(i)} \end{pmatrix} \\ \\ & & J ( \boldsymbol{x^{(i+1)}} - \boldsymbol{x^{(i)}} ) = -\boldsymbol{f^{(i)}} \\ \\ & & \| \boldsymbol{x} \| = \sqrt{ x^2 + y^2 } \\ \\ & & | \| \boldsymbol{x^{(i+1)}} \ \| - \| \boldsymbol{x^{(i)}} \ \| | < \epsilon \end{eqnarray} \tag{2.1-6} \]

2.2.二変数のニュートン法の例

例として、\( x^2 + y^2 =1 \)と\( y = x^3 \)の交点を求めてみます。
まずは関数を次のように書き換えます。
\[ \begin{eqnarray} & f(x,y) & \ = x^2+y^2-1=0 \\ & g(x,y) & \ = x^3-y=0 \end{eqnarray} \tag{2.2-1} \]
次に、(2.1-6)式の\( J(\boldsymbol{x}) \)を計算するために偏微分計算を行います。
\[ \begin{eqnarray} & \frac{\partial f}{\partial x} & \ = 2x , \quad & \frac{\partial f}{\partial y} & \ = 2y \\ & \frac{\partial g}{\partial x} & \ = 3x , \quad & \frac{\partial g}{\partial y} & \ = -1 \end{eqnarray} \tag{2.2-2} \]
ここで、初期解を\( (x^{(0)} \ , y^{(0)} \ ) = ( 1, 1) \)として(2.2-1)、(2.2-2)式に代入し繰り返し計算すると、次の結果が得られます。
\[ ( x_t, y_t ) \simeq ( 0.826031357794845 , 0.563624161955114 ) \]
このときの繰り返し計算数は9回です。

上記例の計算は、本サイトで公開しているExcelツールNewtonMethod.xlsmで計算しています。
その他、Matlabサイトの関数等を代入して計算してみてください。

3.多変数のニュートン法

ニュートン法は多変数に対しても拡張できます。 計算方法は、二変数の場合と同じです。
n個のn変数関数について、
\[ \begin{eqnarray} f_1(x_1,x_2, \cdots, x_n)=0 \\ f_2(x_1,x_2, \cdots, x_n)=0 \\ \cdots \\ f_n(x_1,x_2, \cdots, x_n)=0 \end{eqnarray} \quad または \quad \boldsymbol{f}(\boldsymbol{x})=\boldsymbol{o} \tag{3-1} \]
をすべて満たす解\( (x_{t1},x_{t2},\cdots,x_{tn}) \)を求めます。

Step1:

真の解解\( (x_{t1},x_{t2},\cdots,x_{tn}) \)の近傍点\( (x_1^{(0)} \ ,x_2^{(0)} \ ,\cdots \ , x_n^{(0)} \ ) \)を初期解として適当に決めます。

Step2:

二変数関数のときと同様に、以下のヤコビ行列を求めます。
\[ J( \boldsymbol{x} ) = \begin{pmatrix} \partial f_1 / \partial x_1 & \ \partial f_1 / \partial x_2 & \cdots & \partial f_1 / \partial x_n & \\ \partial f_2 / \partial x_2 & \ \partial f_1 / \partial x_2 & \cdots & \partial f_2 / \partial x_n & \\ \vdots & \vdots & \ddots & \vdots & \\ \partial f_n / \partial x_1 & \ \partial f_n / \partial x_2 & \cdots & \partial f_n / \partial x_n & \end{pmatrix} \tag{3-2} \]
\( \boldsymbol{x} \)に初期解\( \boldsymbol{x^{(0)}} \ \)を代入して\( \boldsymbol{x^{(1)}} \ \)を計算します。
\[ \boldsymbol{f^{(0)}} = J^{(0)}(\boldsymbol{x^{(1)}} - \boldsymbol{x^{(0)}} ) \\ \\ \leftrightarrow \boldsymbol{x} = \boldsymbol{x^{(0)}} - \left(J^{(0)} \right)^{-1} \boldsymbol{f^{(0)}} \tag{3-3} \]

>\( ( J^{(0)}=J(\boldsymbol{x^{(0)}}) ) \)

なお、(3-3)式は二変数のときの式と全く同じで、変数の次元が違うだけです。

Step3:

\( \boldsymbol{x^{(0)}} \ \)と\( \boldsymbol{x^{(1)}} \ \)の大きさ (ノルム:\( \| \boldsymbol{x^{(0)}} \|, \| \boldsymbol{x^{(1)}} \| \))を求め、 その差分、つまり前Stepとの解の変化が微小値\( \epsilon \)(収束判定値)より小さければ、\( x^{(1)} \)を\( ( x_t, y_t) \)の近似解とします。
\[ | \| \boldsymbol{x^{(1)}} \| - \| \boldsymbol{x^{(0)}} \| | < \epsilon \tag{3-4} \] \[ \| \boldsymbol{x} \| = \sqrt{ x_1^2 + x_2^2 + \cdots + x_n^2 } \]
上式を満たさない場合は、解の変化が\( \epsilon \)より小さくなるまで、\( \boldsymbol{x^{(2)}},\boldsymbol{x^{(3)}},\cdots \)を求めて解の差分を評価する操作、Step2、Step3、・・・を繰り返します 。
このように行列表現をとることで、多変数のニュートン法は二変数のニュートン法と全く同じ形で表せます。 多変数で厄介なのは、偏微分係数の数が多く、それらを求める(設定する)のが大変なことです。

参考文献