EffectiveMass.py ダウンロード/コピー
EffectiveMass.py
EffectiveMass.py
1"""
2概要:
3 バンド構造から有効質量と運動量空間の波数kの関係を計算します。
4詳細説明:
5 このスクリプトは、CSVファイルからバンド構造データ（運動量空間の波数kとエネルギーE）を読み込みます。
6 読み込んだデータを用いて、エネルギーの運動量空間の波数kに対する二階微分を計算し、そこから有効質量を導出します。
7 有効質量の計算では、異なる差分幅（nskip）を使用して二階微分を評価し、結果をプロットに表示します。
8 計算された有効質量の逆数が特定の閾値を超えた場合、プロット上でデータの連続性を分断する処理が含まれます。
9 最終的なプロットは、有効質量と運動量空間の波数kの関係を示します。
10主な機能:
11 - CSVファイルからのバンド構造データの読み込み
12 - エネルギーの二階微分に基づいた有効質量の計算
13 - 異なる差分幅での有効質量計算の比較
14 - matplotlibによる結果のグラフ表示
15関連リンク:
16 EffectiveMass_usage
17"""
18import sys
19import numpy as np
20from numpy import sqrt, exp, sin, cos, tan, pi
21import numpy.linalg as LA
22import csv
23from matplotlib import pyplot as plt
24
25
26#===================================
27# physical constants
28#===================================
29pi = 3.14159265358979323846
30pi2 = 2.0 * pi
31h = 6.6260755e-34 # Js";
32hbar = 1.05459e-34 # "Js";
33c = 2.99792458e8 # m/s";
34e = 1.60218e-19 # C";
35e0 = 8.854418782e-12; # C<sup>2</sup>N<sup>-1</sup>m<sup>-2</sup>";
36kB = 1.380658e-23 # JK<sup>-1</sup>";
37me = 9.1093897e-31 # kg";
38R = 8.314462618 # J/K/mol
39a0 = 5.29177e-11 # m";
40
41
42#===================================
43# parameters
44#===================================
45a = 4.0e-10 #m
46
47infile = 'band.csv'
48cutline = 1
49
50#===================================
51# figure configuration
52#===================================
53fontsize = 12
54legend_fontsize = 8
55
56
57#=============================
58# other functions
59#=============================
60def savecsv(outfile, header, datalist):
61 """
62 概要:
63 リストのリスト形式のデータをCSVファイルに保存します。
64 詳細説明:
65 指定されたファイル名でCSVファイルを開き、ヘッダー行を書き込んだ後、datalistの各要素を列として、
66 行ごとにデータを書き込みます。ファイル書き込みに失敗した場合はエラーメッセージを表示します。
67 引数:
68 :param outfile: 出力するCSVファイルのパスとファイル名。
69 :type outfile: str
70 :param header: CSVファイルのヘッダー行の文字列のリスト。
71 :type header: list of str
72 :param datalist: 出力するデータを含むリストのリスト。各内部リストは1つの列に対応します。
73 例: [[col1_data1, col1_data2], [col2_data1, col2_data2]]
74 戻り値:
75 :returns: なし
76 :rtype: None
77 """
78 try:
79 print("Write to [{}]".format(outfile))
80 f = open(outfile, 'w')
81 except:
82# except IOError:
83 print("Error: Can not write to [{}]".format(outfile))
84 else:
85 fout = csv.writer(f, lineterminator='\n')
86 fout.writerow(header)
87# fout.writerows(data)
88 for i in range(0, len(datalist[0])):
89 a = []
90 for j in range(len(datalist)):
91 a.append(datalist[j][i])
92 fout.writerow(a)
93 f.close()
94
95
96def read_csv(fname):
97 """
98 概要:
99 指定されたCSVファイルから数値データを読み込みます。
100 詳細説明:
101 CSVファイルをオープンし、最初の行をラベルとして読み込みます。
102 残りの行から2列のデータを浮動小数点数として読み取り、それぞれのリストに格納します。
103 数値変換に失敗した行は警告メッセージを表示してスキップします。
104 引数:
105 :param fname: 読み込むCSVファイルのパスとファイル名。
106 :type fname: str
107 戻り値:
108 :returns: 1列目のラベル文字列、2列目のラベル文字列、1列目のデータリスト、2列目のデータリストのタプル。
109 :rtype: tuple of (str, str, list of float, list of float)
110 """
111 x = []
112 y = []
113 with open(fname) as f:
114 fin = csv.reader(f)
115 labels = next(fin)
116 xlabel = labels[0]
117 ylabel = labels[1]
118 for row in fin:
119 try:
120 x.append(float(row[0]))
121 y.append(float(row[1]))
122 except:
123 print("Warning: Invalid float data [{}] or [{}]".format(row[0], row[1]))
124
125 return xlabel, ylabel, x, y
126
127
128def main():
129 """
130 概要:
131 有効質量と運動量空間の波数kの関係を計算し、結果をプロットします。
132 詳細説明:
133 `infile` で指定されたCSVファイルからバンド構造データ（運動量空間の波数kとエネルギーE）を読み込みます。
134 波数kの差分dkとエネルギーEを用いて、エネルギーの運動量空間の波数kに対する二階微分を計算し、
135 そこから有効質量を導出します。
136 異なる差分幅（nskip = 1とnskip = 4）で二階微分を計算し、それぞれ有効質量を求めます。
137 有効質量の逆数 abs(minv) が 1.0e20 以下の場合は、有効質量のデータとして None を格納し、
138 プロット上でデータの連続性を分断します (cutline が 1 の場合)。
139 また、有効質量の符号が反転した場合もプロットを分断します。
140 最終的に、計算された有効質量を電子の静止質量 me で正規化し、
141 matplotlib を用いて運動量空間の波数kに対する正規化有効質量をプロットし、表示します。
142 引数:
143 :param なし: この関数は引数を取りません。
144 戻り値:
145 :returns: なし
146 :rtype: None
147 """
148 klabel, Elabel, k, E = read_csv(infile)
149 print("k=", k)
150 print("E=", E)
151
152 nk = len(k)
153 dk = k[1] - k[0]
154
155 km = hbar * hbar * (pi2 / a)**2.0
156 nskip = 1
157 xk = []
158 ymc = []
159 signprev = None
160 print("")
161 print("Differentiate for nskip = ", nskip)
162 print("i E[i-1] E[i] E[i+1] 1.0/m*")
163 for i in range(nskip, nk - nskip, nskip):
164 d2Edk2c = (E[i+nskip] + E[i-nskip] - 2 * E[i]) * e / pow(nskip * dk, 2.0)
165 minv = d2Edk2c / km
166 print(i, E[i-1], E[i], E[i+1], minv)
167 if abs(minv) <= 1.0e20: # << 1.0/me ~ 1e30
168 if cutline:
169 xk.append(k[i])
170 ymc.append(None)
171 signprev = -signprev
172 continue
173 else:
174 m = km / d2Edk2c
175
176 if signprev is None:
177 signprev = m
178 elif signprev * m < 0.0:
179 if cutline:
180 xk.append(k[i])
181 ymc.append(None)
182 signprev = m
183
184 xk.append(k[i])
185 ymc.append(m / me)
186
187 nskip = 4
188 xk2 = []
189 ymc2 = []
190 signprev = None
191 print("")
192 print("Differentiate for nskip = ", nskip)
193 print("i E[i-1] E[i] E[i+1] 1.0/m*")
194 for i in range(nskip, nk - nskip, nskip):
195 d2Edk2c = (E[i+nskip] + E[i-nskip] - 2 * E[i]) * e / pow(nskip * dk, 2.0)
196 minv = d2Edk2c / km
197 print(i, E[i-1], E[i], E[i+1], minv)
198 if abs(minv) <= 1.0e20: # << 1.0/me ~ 1e30
199 if cutline:
200 xk2.append(k[i])
201 ymc2.append(None)
202 signprev = -signprev
203 continue
204 else:
205 m = km / d2Edk2c
206
207 if signprev is None:
208 signprev = m
209 elif signprev * m < 0.0:
210 if cutline:
211 xk2.append(k[i])
212 ymc2.append(None)
213 signprev = m
214
215 xk2.append(k[i])
216 ymc2.append(m / me)
217
218 plt.plot(xk, ymc, linewidth = 0.5, marker = 'o', markersize = 1.0, label = 'nskip = 1')
219 plt.plot(xk2, ymc2, linewidth = 0.5, marker = 'o', markersize = 1.0, label = 'nskip = 2')
220 plt.xlabel(klabel)
221 plt.ylabel("m$_e$ / m$_e^0$")
222 plt.xlim([-0.5, 0.5])
223# plt.ylim([-0.5, 0.5])
224 plt.tight_layout()
225
226 plt.pause(0.1)
227 print("Press ENTER to exit>>", end = '')
228 input()
229
230
231if __name__ == "__main__":
232 main()