【入門】対称行列と二次形式(Python)【数値計算】

【入門】対称行列と二次形式(Python)【数値計算】 数値計算
【入門】対称行列と二次形式(Python)【数値計算】

MATLAB、Scilab、Scilab、Julia比較ページはこちら
https://www.simulationroom999.com/blog/comparison-of-matlab-python-scilab/

はじめに

の、
MATLAB,Python,Scilab,Julia比較 第2章 その34【対称行列と二次形式⑥】

を書き直したもの。

正規方程式を導出するための対称行列と二次形式の話。
今回は二次形式の多項式表現と行列表現の計算をPythonで実施する。

ロードマップ【再掲】

ロードマップと今回使用する数式を再掲。

正規方程式に至る道、二次形式、対称行列、二次形式の微分、グラム行列、二乗和誤差、正規方程式

\(
x^TA\circ x^T
\begin{bmatrix}
1\\
\vdots\\
1
\end{bmatrix}
\)

計算する二次形式の数式

\(
3x^2+2y^2+5xy
\)

今回は、Pythonで実施する。

Pythonコード

Pythonコードは以下になる。

import numpy as np
import matplotlib.pyplot as plt

a=3
b=2
c=5
A=np.array([[a,c/2],[c/2,b]])

N=6
ax=np.linspace(0,10,N)
ay=np.linspace(0,10,N)
x,y=np.meshgrid(ax,ay)

fig = plt.figure(figsize=(8, 8))

ax = fig.add_subplot(2, 1, 1, projection='3d')
polyY=a*x**2 +b*y**2 +c*x*y
ax.plot_wireframe(x,y,polyY)

X=np.block([[x.reshape(-1)],[y.reshape(-1)]])
Y=X.T@A*X.T@np.ones((2,1))
matY=Y.reshape(N,N)
ax = fig.add_subplot(2, 1, 2, projection='3d')
ax.plot_wireframe(x, y, matY)

plt.show()

print("polyY=")
print(polyY)
print("matY=")
print(matY)

実行結果

実行結果は以下。

二次形式の多項式表現と行列表現の計算(Python)
polyY=
[[   0.   12.   48.  108.  192.  300.]
 [   8.   40.   96.  176.  280.  408.]
 [  32.   84.  160.  260.  384.  532.]
 [  72.  144.  240.  360.  504.  672.]
 [ 128.  220.  336.  476.  640.  828.]
 [ 200.  312.  448.  608.  792. 1000.]]
matY=
[[   0.   12.   48.  108.  192.  300.]
 [   8.   40.   96.  176.  280.  408.]
 [  32.   84.  160.  260.  384.  532.]
 [  72.  144.  240.  360.  504.  672.]
 [ 128.  220.  336.  476.  640.  828.]
 [ 200.  312.  448.  608.  792. 1000.]]

考察

処理の流れとしてはMATLABと一緒。
二次形式の多項式表現、行列表現で同じ結果になったのも一緒。

気になる点としては、projection=’3d’とplot_wireframe。
projection=’3d’は3Dグラフを使用する場合のオプション。
plot_wireframeはワイヤーフレームをグラフ表示するためのメソッドになる。
plot_wireframeを使用しても。projection=’3d’のオプションが無いと3Dグラフとして表示されない。

まとめ

  • 二次形式の多項式表現と行列表現の計算をPythonで実施。
  • 3Dグラフを表示する際は、projection=’3d’が必要。
  • plot_wireframeでワイヤーフレームでグラフ表示ができる。

MATLAB、Python、Scilab、Julia比較ページはこちら

コメント

タイトルとURLをコピーしました