【Python】モジュールとパッケージ
モジュールとパッケージの違いが分かったかもしれないのでメモります。
モジュール
- コードをまとめた一つのファイル。
import ファイル名（拡張子なし）でファイルに記述した関数や変数が使えるようになる。モジュールをインポートするという。- Ex. test1 (ファイル) の func1 (関数) を呼び出す
import test1
test1.func1()
- インポートするモジュール（ファイル）をどのように探すかは、sysモジュールのpath変数に記載されいる。基本的にメインプログラムと同一フォルダが最優先に検索される。
パッケージ
- モジュール（ファイル）をまとめ階層構造に組織したもの（フォルダをイメージ）。
from フォルダ名 imort モジュール（ファイル名）でモジュールの関数や変数を使うことができる。- パッケージの中には、
__init__.pyというファイルが必ず含まれる。このファイルがあることによりフォルダがパッケージとして認識される。
参考文献
虚数iのルートは複素数
複素関数の勉強を基礎からしています。岩波の本を使って勉強しています。始めの章に複素数の基礎について書かれていましたが、どのように自然数から複素数まで数を拡張していったか簡潔に書かれていて、説明の上手さに感動しました。
さて、その複素数の基礎部分にちらっと複素数のn乗根は複素数であると書かれていました。特に証明は載っていなくどうやって導出するんだとうと思い、ネットで調べたらでてきました。
サイトには、虚数iを例に、複素数のルートの導出方法について書かれています。
虚数iのルートは複素数 ですね。何だか不思議な感じがしますが、複素平面で考えるとすっきりしますね。複素平面だとn乗根も簡単に求めることができますね。複素平面すごい。
複素数についてもっと勉強を進めていきたいと思います。
参考文献
Pythonで波動方程式の数値計算②
以前作成した一次元波動方程式の数値計算の初期条件を変えてどのような挙動になるか見てみました。
例えば上のアニメーションのように弦の真ん中でインパルス上の波を起こすと、時間が経つにつれ波が左右に移動し（振幅は半分になる）、弦の両端で位相がひっくりかえり逆方向に進みます。二つの波がぶつかると合成し振幅が大きくなります。波の挙動がちゃんと表現できていました。面白いですね。
コード
# matplotlibのアニメーション表示用 %matplotlib nbagg import numpy as np import matplotlib.pyplot as plt import matplotlib.animation as animation import math # 各種定数の設定 t0 = 0 tM = 2 M = 40 dt = (tM-t0)/M # 時間間隔 x0 = 0 xN = 1 N = 20 dx = (xN-x0)/N # 空間間隔 a = (dt/dx)**2 # 差分方程式に出てくる係数 ux0 = 0 # 境界条件 uxN = 0 # 境界条件 A = 2 # 初期条件の式の振幅 # 変数の初期化 u = np.zeros([M+1, N+1]) # 行：時間、列：ｘ # 初期条件1 (t=0でのuの分布) u(x)=A*exp(-5*(x-0.5)^2) for xi in range(N+1): u[0,xi] = A*math.exp(-((xi-N/2)**2)*5) # 初期条件2 (t=dtでのuの分布) u[1,0] = 0 u[1,N] = 0 for xi in range(1,N): u[1,xi] = u[0,xi] + dt*0 + a/2*(u[0,xi+1] -2*u[0,xi] + u[0,xi-1]) # 差分解 (t=2*dt~M*dt) for ti in range(1,M): # 両端の弦の変異は常に0 u[ti+1,0] = ux0 u[ti+1,N] = uxN # 両端以外を差分方程式により更新 for xj in range(1,N): u[ti+1, xj] = 2*u[ti,xj] -u[ti-1,xj] + a*(u[ti,xj+1] -2*u[ti,xj] + u[ti,xj-1]) # matplotlibで表示 x = np.linspace(x0, xN, N+1) # maplotlibで表示するためのx軸の設定 fig = plt.figure(figsize=(6,4)) ax = fig.add_subplot(1,1,1) # アニメ更新用の関数 def update_func(i): # 前のフレームで描画されたグラフを消去 ax.clear() ax.plot(x, u[i,:], color='blue') ax.scatter(x, u[i,:], color='blue') # 軸の設定 ax.set_ylim(-2.22, 2.22) # 軸ラベルの設定 ax.set_xlabel('x', fontsize=12) ax.set_ylabel('u', fontsize=12) # サブプロットタイトルの設定 ax.set_title('Time: ' + '{:.2f}'.format(tM*i/M)) ani = animation.FuncAnimation(fig, update_func, frames=M, interval=100,repeat=True) # アニメーションの保存 #ani.save('wave_eq3.gif', writer='imagemagick') # 表示 plt.show()
参考文献
【Python】 タプルのアンパック
タプルのアンパックを知らなかったのでメモを残します。
タプルのアンパックとは、タプルの各要素を一度に複数の変数に代入できるテクニックです。
people_tuple = ('Taro', 'Jiro', 'Saburo') staff1, staff2, staff3 = people_tuple print('staff1: ' + staff1) print('staff2: ' + staff2) print('staff3: ' + staff3) # staff1: Taro # staff2: Jiro # staff3: Saburo
people_tupleの各要素を変数staff1、staff2、staff3にそれぞれ代入したいとき、一行でシンプルに書けます。使える場面あったら使っていきたいと思います。
参考文献
【Python】 内包表記を知る
入門 Python 3を読んでいたら内包表記という書き方を知りました。内包表記を使うとfor文（+条件付きfor文）を一行でシンプルに書けるようになります。本によると内包表記が使えるかがpython初心者卒業レベル卒業の目安になるらしいです。卒業に向け内包表記が使えるときは使っていきたいと思います。
内包表記とは
内包表記はFor文をシンプルに一行で書く書き方です。さらに、if文が入ったfor文も一行で書くことができます。以下ではリストを対象としたリスト内包表記について書きます。
基本形
基本的な書き方は以下の通り。
- [式 for 変数 in イテラブルオブジェクト]
式は「f(x)」、変数は「x」とすると分かりやすいかもしれません。イテラブルオブジェクトがN個の要素x1~xNから成るとすると、リスト内包表記の出力は[f(x1), f(x2), ・・・, f(xN)]というN個の要素を持ったリストになります。
例：初項1、交差2の等差数列の第10項までのリスト
num_list = [1+2*(i-1) for i in range(1, 11)] print(num_list) # [1, 3, 5, 7, 9, 11, 13, 15, 17, 19]
発展形
for文にさらに条件式（if文）を入れることもできます。書き方は以下の通り。
- [式 for 変数 in イテラブルオブジェクト if 条件式]
これは、「条件式 == True」を満たすイテラブルオブジェクトのみ出力されることを意味します。例を見た方が分かりやすいかもしれません。
例：初項1、交差2の等差数列の第10項まで数列の内、偶数項を抽出したリスト
num_list2 = [1+2*(i-1) for i in range(1, 11) if i % 2 == 0] print(num_list2) # [3, 7, 11, 15, 19]
参考文献
ガウスの消去法の実装
数値解析の本（数値計算 (理工系の基礎数学 8)）に記載されていたガウスの消去法をpythonで実装してみました。部分ピボットはまだ実装していません。
ガウスの消去法とは
（N元）連立一次方程式の解を求める方法のひとつで、大きく分けて前進消去と後退代入の二つの手続きによって解を導出していきます。
一つ目の操作である前進消去は、方程式を上から下に並べたとして、上の方の一次方程式を用いて他の下の方の一次方程式の変数を一つ消去します。これを繰り返していくと、最終的に変数が一つの方程式が導出されます。
二つ目の後退代入は、下の方にある変数が一つの方程式を解き変数の値を確定し、その値を一つ上の方程式に代入しまた別の変数の値を確定していきます。繰り返すことで連立一次方程式の解を導出することができます。
ガウスの消去法の実装
上の3元連立1次方程式を例にガウスの消去法によって解を求めてみました。コードは下記の通りです。
import numpy as np # 各係数の設定 A = np.array([[6, 3, 2], [5, 1, 8], [9, 4, 7]], dtype='float') b = np.array([5, 1, 4], dtype='float') N = A.shape[0] x = np.zeros(3) # 前進消去 for k in range(N-1): for i in range(k+1, N): t = A[i, k] / A[k, k] for j in range(k+1, N): A[i, j] = A[i, j] - t*A[k, j] b[i] = b[i] - t*b[k] # 後退代入 x[N-1] = b[N-1]/A[N-1, N-1] for i in reversed(range(N-1)): u = 0 for k in range(i+1, N): u = u + A[i, k]*x[k] x[i] = (b[i] - u)/A[i,i] print(x)
これを実行すると、 となります。
この解が正しいかnumpyのlinalg.solve関数を用いて確認してみました。
import numpy as np A = np.array([[6, 3, 2], [5, 1, 8], [9, 4, 7]]) b = np.array([5, 1, 4]) x2 = np.linalg.solve(A, b) print(x2)
と、同じ値になっていました。実装に大きな間違いはなさそうでよかったです。
今度は部分ピボット選択をコードに取り入れてみたいと思います。
