はじめに
みなさんこんにちは．豊田工業大学Kaggleサークルアドベントカレンダー20日目担当の@Kofukuです．
現在は学部4年生でPythonは研究で行った数値計算のデータ処理に使っています．
就職をしない学部四年生は修士課程に進学します．修士学生として必要な能力は様々ですが，修士学生になるには大学院入試に合格するだけでいいのです．
ただ，大学院入試合格に向けて困難なことの一つとしては解答がないことです．答えがなければ自分の計算ミスに気づくことも難しくなります．
今回はPythonを使い線形代数と微分積分の問題を解き，解答がない問題に取り組む計算力がない学生でも困りにくいようにします．
実行環境
-
OS
- ProductName: macOS
- ProductVersion:11.6
- BuildVersion: 20G165
-
Python
- anacond3-2020.7
- sympy 1.6.1
- numpy 1.21.4
線形代数
大学院入試で頻出の線形代数の分野として行列式や逆行列やKernelを求める問題があげられます．
詳しい行列式と逆行列の求め方は,豊田工業大学Kaggleサークルアドベントカレンダー1日目の記事に書いてあります．
Numpy入門
名古屋大学機械航空系2018基礎の2-2も上記のアドベントカレンダーの内容を元にすれば解けるはずです．
名古屋大学機械航空系2018基礎問題
それでは2-2をPythonで解いていきます．
import numpy as np
from scipy import linalg
# 行列Cを定義する
C = np.array([[1, 0, -2], [ 0, 1, 1],[2, 1, -2 ]])
print("C\n",C)
# C^2を求める
square_C = np.dot(C,C)
print("Ans1 C^2\n", square_C)
#2)の固有方程式を解く代わりに，固有値と固有ベクトルを求める．
eig_value, eig_vector = linalg.eig(C)
print("Cの固有値", eig_value)
print("Cの固有ベクトル\n", eig_vector)
#3)逆行列を求める.
inverse_C = linalg.inv(C)
print("Cの逆行列\n",inverse_C)
# おまけCのカーネルを求める．
kernel_C = linalg.null_space(C)
print("Cのカーネル", kernel_C)
これを実行すると
C
[[ 1 0 -2]
[ 0 1 1]
[ 2 1 -2]]
Ans1 C^2
[[-3 -2 2]
[ 2 2 -1]
[-2 -1 1]]
Ans2 Cの固有値 [-0.5+0.8660254j -0.5-0.8660254j 1. +0.j ]
Ans2 Cの固有ベクトル
[[-7.07106781e-01+0.00000000e+00j -7.07106781e-01-0.00000000e+00j
4.47213595e-01+0.00000000e+00j]
[ 3.53553391e-01-8.32667268e-17j 3.53553391e-01+8.32667268e-17j
-8.94427191e-01+0.00000000e+00j]
[-5.30330086e-01+3.06186218e-01j -5.30330086e-01-3.06186218e-01j
-3.53471931e-16+0.00000000e+00j]]
Ans3 Cの逆行列
[[-3. -2. 2.]
[ 2. 2. -1.]
[-2. -1. 1.]]
Cのカーネル []
出力はこのようになります．
Cのカーネルとしては何も出力がされていないので，Cのカーネルは0ベクトルであることがわかります．
##微分積分
Sympyによる積分方法は豊田工業大学アドベントカレンダー9日目のScipyによるシンボリック変数での微分積分に詳しく書いてあるので参考にしてください．
そのほかにも積分法についてまとめた記事は豊田工業大学アドベントカレンダー19日目のPythonを用いた数値積分の方法に詳しく書いてあるので参考にしてください．
カージオイドの面積積分と線積分の値を求めよ(名古屋大学機械航空系2021基礎)
カージオイドは次の式で表される.
r=a(1+\cos \theta)
解析的な計算方法は高校数学の美しい物語にあるのでそちらを参照してください．
面積分の値は$\frac{3}{2}\pi a^2$,線積分の値は$8a$です．
媒介変数表示を使って積分値を求める手法
媒介変数表示面積
import sympy as sym
import matplotlib.pyplot as plt
import numpy as np
t = sym.Symbol("t")
t = sym.Symbol("t")
a = sym.Symbol("a")
S = 2 * abs(sym.integrate(y * sym.diff(x, t), (t, 0, np.pi)))
print("面積",S)
面積 4.71238898038469*Abs(a**2)
媒介変数表示線積分
import sympy as sym
import numpy as np
t = sym.Symbol("t")
x = a * (1 + sym.cos(t)) * sym.cos(t)
y = a * (1 + sym.cos(t)) * sym.sin(t)
l = 2 * sym.integrate((((sym.diff(x,t) ** 2) + (sym.diff(y,t) **2))**0.5), (t, 0, np.pi))
print("線積分",l)
線積分 2*(a**2)**0.5*Integral((sin(t)**2 + cos(t)**2)**0.5*(sin(t)**2 + cos(t)**2 + 2*cos(t) + 1)**0.5, (t, 0, 3.14159265358979))
面積分は正しい相対を出力しますが，線積分に関しては出力してくれません．
この問題を解決するにはsqrtの中身を綺麗にする必要があります．sympyのsimplify関数を使いsin,cosやsqrtを整理することで積分をしてくれるようになります．
import sympy as sym
import numpy as np
t = sym.Symbol("t")
a = sym.Symbol("a")
x = a * (1 + sym.cos(t)) * sym.cos(t)
y = a * (1 + sym.cos(t)) * sym.sin(t)
l = sym.integrate(sym.simplify(((sym.diff(x,t) ** 2) + (sym.diff(y,t) **2))**0.5), (t, 0, np.pi))
print("線積分",l)
線積分 8.0*(a**2)**0.5
このように求めることができました．パラメータで積分を行った場合は，解析解と一致しました．
極座標形式面積と線積分
次は極座標形式のまま面積と線積分を求めてみます．
import sympy as sym
import numpy as np
t = sym.Symbol("t")
a = sym.Symbol("a")
r = (1 + sym.cos(t))
S = 2 * 1/2 * sym.integrate(r**2, (t,0,np.pi))
print("面積",S)
l = sym.integrate(sym.simplify(sym.sqrt((r**2) + (sym.diff(r,t) **2))), (t, 0, np.pi))
print("線積分",l)
面積 4.71238898038469*a**2
線積分 sqrt(2)*Integral(sqrt(a**2*cos(t) + a**2), (t, 0, 3.14159265358979))
このように極座標形式では面積を求めることしかできませんでした．Sympyに積分をさせるときはある程度綺麗な形にする必要がありそうです．
Sympyが積分してくれない場合でも，途中経過を出力することで人間が簡単に積分することができる場合があります．
まとめ
- 線形代数と微分積分の計算をPythonを使用する方法を記載しました．
- ベクトル解析や常微分方程式もこれらを応用すれば解けるはずです．
- 大学院入試や講義課題で計算が難しい場合や面倒な場合はPythonを使ってやってみるのはいかがでしょうか？