Python SciPy : 連立非線形方程式の求根アルゴリズム

SciPy の連立非線形方程式の求根アルゴリズムについてまとめます。

関数について

SciPy には表 1 の連立非線形方程式の求根アルゴリズムがあります。
表 1 の他に fsolve 等の関数がありますが、全ての計算ルーチンは root の method で選ぶことができるので、実質 root だけ使用します。

Table 1: 連立非線形方程式の求根アルゴリズム
関数 説明 ヤコビ行列
root(method="hybr") 修正 Powell Hybrid 法
root(method="lm") Levenberg-Marquardt 法
root(method="broyden1") good Broyden 法 不要
root(method="broyden2") bad Broyden 法 不要
root(method="krylov") Newton-Krylov 法 不要
root(method="anderson") Anderson 法 不要
root(method="linearmixing") スカラーヤコビ行列近似 不要
root(method="excitingmixing") 調整対角ヤコビ行列近似 不要
root(method="diagbroyden") 対角 Broyden ヤコビ行列近似 不要

hybr は MINPACK の hybrd を使用しています。
収束が高速で、小規模な問題では非常に優秀です。
ヤコビ行列を与えなくても使えますが、ヤコビ行列を与えればより速く収束します。

lm は速さでは hybr に劣りますが、より頑強な方法です。
こちらもヤコビ行列を与えると収束が速くなります。

broyden1 はヤコビ行列の逆行列計算をしない省メモリな方法のため、問題の規模が大きくなって hybr や lm が使えない場合に有用です。

broyden2 は broyden1 よりも計算コストが小さく、大規模な問題でも使える方法です。
また、broyden1 では問題によってゼロ割りが発生するのことがありますが、それを回避できます。
ただし、収束の確実性は broyden1 に劣ります。

krylov は ヤコビ行列の逆行列計算に krylov 部分空間法を用います。
計算コストが小さく、並列計算化しやすいため、大規模問題に広く使われる方法です。
オプションで近似解の決定法を選択することができます。
krylov 部分空間法については この資料 が詳しいと思いました。

anderson はヤコビ行列全体ではなく各反復の方向微分係数を計算することで計算コストを小さくした方法で、大規模問題に適しています。

以下は 単純な反復法で計算負荷が小さく大規模問題に適していますが、適用できる問題が限られます。

linearmixing は一つのスカラーのパラメータで対角行列のヤコビ行列を近似します。

excitingmixing は linearmixing のスカラーパラメータに指定範囲内で調整してヤコビ行列を近似します。

diagbroyden は前回の反復でのヤコビ行列を用いてヤコビ行列を近似します。

最後に、SciPy には既に多くの大規模問題用関数がありますが、ver0.17.0 から DF-SANE (Derivative-Free Spectral Approach for Nonlinear Equations)アルゴリズムが追加されることになっています。

# coding: utf-8
import numpy as np
from scipy import optimize


def fun(x):
    return [x[0]  + 0.5 * (x[0] - x[1])**2 - 1,
	    0.5 * (x[1] + x[0])**3 + x[1]]


def jac(x):
    return np.array([[1 + (x[0] - x[1]),
		      - (x[0] - x[1])],
		     [1.5 * (x[1] + x[0])**2,
		      1 + 1.5 * (x[1] + x[0])**2]])

sol1 = optimize.root(fun, [0, 0], method='broyden1')
sol2 = optimize.root(fun, [0, 0], jac=jac, method='lm')

[sol1.x, sol2.x]
[array([ 0.68832603, -0.10119789]), array([ 0.68832617, -0.10119751])]

コメント

Comments powered by Disqus
書籍更新情報
2016-10-21
Pythonによる科学技術計算 基礎編
PDF版の販売を開始しました。
販売ページはこちら

2016-09-09
Pythonによる科学技術計算 基礎編
1.2版への更新が可能になりました。
サポートページはこちら
電子書籍
Pythonによる科学技術計算 基礎編
Kindle ストア、Leanpubで販売中です
Pythonによる科学技術計算 基礎編
PDF版の販売はこちら
同人誌
技術書典(2016/6/25)
Emacs/org-modeのPDF作成術
電子版をBOOTHで販売中です
Emacs/org-modeのPDF作成術
Share