【第17話】波数ベクトルってなに？【Pythonコピペで量子力学完全攻略マニュアル】

#################################################################
## ２次元自由粒子の運動（E = 1.0[eV]）
#################################################################
import numpy as np
import matplotlib as mat
import matplotlib.pyplot as plt
import matplotlib.animation as animation

#図全体
aspect = 0.9
fig = plt.figure(figsize=(10, 10 * aspect))
#全体設定
plt.rcParams['font.family'] = 'Times New Roman' #フォント
plt.rcParams['font.size'] = 16 #フォントサイズ
plt.rcParams["mathtext.fontset"] = 'cm' #数式用フォント
#カラーリストの取得
colors = plt.rcParams['axes.prop_cycle'].by_key()['color']

#################################################################
## 物理定数
#################################################################
#プランク定数
h = 6.62607015 * 1.0E-34
hbar = h / (2.0 * np.pi)
#電子の質量
me = 9.1093837015 * 1.0E-31
#電子ボルト
eV = 1.602176634 * 1.0E-19
#虚数単位
I = 0.0 + 1.0j

#################################################################
## 物理系の設定
#################################################################
#電子のエネルギー
E = 1.0 * eV
#波数
k = np.sqrt(2.0 * me * E / hbar**2)
#角振動数
omega = E/hbar
#１周期
T = 2.0 * np.pi / omega
#時間間隔
#dt = 0.2 * 10**-16
dt = T / 100
#時間刻み数
NT = 100
#空間刻み数
NX = 500
#空間刻み間隔
nm = 10**-9
#空間幅
x_min = -2.0
x_max = 2.0
#入射角度
theta = np.pi/6

#アニメーションフレーム作成用
def update( tn ):
    #グラフ消去
    fig.clear()
    #両軸位置の設定
    axes = fig.add_axes(    [0.07, 0.05/aspect, 0.78, 0.78/aspect])
    #カラーバー位置の設定
    bar_axes = fig.add_axes([0.87, 0.05/aspect, 0.05, 0.78/aspect])

    #実時間
    t = tn * dt
    #描画間隔
    delta_x = (x_max - x_min) / NX
  #x軸・y軸のメッシュ生成
    x = np.arange(x_min, x_max + delta_x, delta_x) * nm
    y = np.arange(x_min, x_max + delta_x, delta_x) * nm
    X, Y = np.meshgrid(x, y)
    #波数の指定
    kx = k * np.cos(theta)
    ky = k * np.sin(theta)
    #平面波の実部
    Z = np.exp( I * ( kx * X + ky * Y ) - I * omega * t ).real
    Z = Z[:-1, :-1]

    #グラフ描画
    img = axes.pcolormesh( X/nm, Y/nm, Z, vmin=-1.0, vmax=1.0, cmap=cmap)
    #カラーバーの追加
    fig.colorbar(img, cax=bar_axes)
    #見出し
    axes.set_title( u"時刻:" + r"$" + str(tn/10) + r"[{\rm fs}]$" , fontname="Meiryo" )
    #軸ラベル
    #axes.set_xlabel("x")
    #axes.set_ylabel("y")

#カラーマップの配色配列
color_list = []
color_list.append( [ 0,   "blue"  ] )
color_list.append( [ 0.5, "black" ] )
color_list.append( [ 1,   "red"   ] )
#カラーマップ用オブジェクトの生成
cmap = mat.colors.LinearSegmentedColormap.from_list('cmap', color_list)

#アニメーションの設定
ani = animation.FuncAnimation( fig, update, range(NT), interval=10)
#アニメーションの保存
#ani.save("output.html", writer=animation.HTMLWriter())
#ani.save('output.gif', writer="imagemagick")
#ani.save("output.mp4", writer="ffmpeg", dpi=300)


#グラフの表示
#plt.show()