ベッセル関数とその導関数の零点
モジュールのインポート
import numpy as np
from scipy.special import jv, yv # ベッセル関数
from scipy.special import jn_zeros, jnp_zeros, jnjnp_zeros # ベッセル関数の零点
import matplotlib.pyplot as plt
import scienceplots
#Warning : As of version 2.0.0, you need to add import scienceplots before setting the style (plt.style.use('science')).
plt.style.use(['science', 'notebook'])
ユーザ関数
第1種ベッセル関数の導関数は，
def dbesj(m,x): # 第1種ベッセル関数の導関数
f = (jv(m-1,x)-jv(m+1,x))*0.5
return f
フォーマット付出力のために，
def table2(sx,sy1,sy2,x,y1,y2):
nx = len(x)
print(f"{sx:>8s}",f"{sy1:>8s}",f"{sy2:>8s}")
for i in range(nx):
print(f'{x[i]:8.3f}',f'{y1[i]:8.3f}',f'{y2[i]:8.3f}')
return
ベッセル関数の計算
0次および1次の第1種ベッセル関数を計算すると，
x_min, x_max, n_x = 0.0, 30.0, 11
x = np.linspace(x_min, x_max, n_x)
y0 = jv(0,x) # 0次の第1種ベッセル関数
y1 = jv(1,x) # 1次の第1種ベッセル関数
table2('x','J_0(x)','J_1(x)', x, y0, y1)
x J_0(x) J_1(x)
0.000 1.000 0.000
3.000 -0.260 0.339
6.000 0.151 -0.277
9.000 -0.090 0.245
12.000 0.048 -0.223
15.000 -0.014 0.205
18.000 -0.013 -0.188
21.000 0.037 0.171
24.000 -0.056 -0.154
27.000 0.073 0.137
30.000 -0.086 -0.119
1次の第1種ベッセル関数とその導関数を計算すると，
y2 = dbesj(1,x) # 1次の第1種ベッセル関数の導関数 (TE1m-mode)
table2('x','J_1(x)','dJ_1(x)/dx', x, y1, y2)
x J_1(x) dJ_1(x)/dx
0.000 0.000 0.500
3.000 0.339 -0.373
6.000 -0.277 0.197
9.000 0.245 -0.118
12.000 -0.223 0.066
15.000 0.205 -0.028
18.000 -0.188 -0.003
21.000 0.171 0.028
24.000 -0.154 -0.050
27.000 0.137 0.068
30.000 -0.119 -0.082
プロットして，
fig = plt.figure() # グラフ領域の作成
x = np.linspace(0.0, 30.0, 151)
plt.grid(color = "gray", linestyle="--")
plt.minorticks_on() # #補助目盛りをつける
plt.xlabel(r"$x$") # x軸のラベル設定
plt.ylabel(r"Bessel functions") # y軸のラベル設定
plt.ylim(-1.0, 1.0) # y軸範囲の設定
plt.plot(x, jv(1,x), label=r"$J_1(x)$")
plt.plot(x, dbesj(1,x), label=r"$J_1'(x)$")
plt.legend(ncol=1, loc='upper right', fancybox=False, frameon=True)
fig.savefig('p3_tlt_s14_1.pdf')
plt.show()
ベッセル関数の零点
整数次$n$の第1種ベッセル関数の零点$J_n(x)=0$ は，
https://docs.scipy.org/doc/scipy/reference/generated/scipy.special.jn_zeros.html
nt = 10
n_max = 3
for n in range(n_max+1):
print(n, jn_zeros(n ,nt))
0 [ 2.40482556 5.52007811 8.65372791 11.79153444 14.93091771 18.07106397
21.21163663 24.35247153 27.49347913 30.63460647]
1 [ 3.83170597 7.01558667 10.17346814 13.32369194 16.47063005 19.61585851
22.76008438 25.90367209 29.04682853 32.18967991]
2 [ 5.1356223 8.41724414 11.61984117 14.79595178 17.95981949 21.11699705
24.27011231 27.42057355 30.5692045 33.71651951]
3 [ 6.3801619 9.76102313 13.01520072 16.22346616 19.40941523 22.58272959
25.7481667 28.90835078 32.06485241 35.21867074]
整数次$n$（$n\le 1200$）までの第1種ベッセル関数および導関数の零点
$J_n(x)=0$，$J_n'(x)=0$ は，
https://docs.scipy.org/doc/scipy/reference/generated/scipy.special.jnjnp_zeros.html
nt = 31
zo, n_, m_, t_ = jnjnp_zeros(nt)
n_ # 次数
array([0, 1, 0, 2, 1, 0, 3, 2, 4, 1, 0, 3, 5, 2, 1, 0, 6, 4, 3, 2, 1, 7,
0, 5, 4, 8, 3, 6, 2, 1, 0], dtype=int32)
m_ # m番目の零点
array([0, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 1, 1, 2, 2, 2, 1, 1, 2, 2, 3, 1,
3, 1, 2, 1, 2, 1, 3, 3, 3], dtype=int32)
t_ # Jn(x)=0のとき0, Jn'(x)=0のとき1
array([1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1,
0, 0, 1, 1, 0, 0, 1, 0, 1], dtype=int32)
zo # 零点の値の配列
array([ 0. , 1.84118378, 2.40482556, 3.05423693, 3.83170597,
3.83170597, 4.20118894, 5.1356223 , 5.31755313, 5.33144277,
5.52007811, 6.3801619 , 6.41561638, 6.70613319, 7.01558667,
7.01558667, 7.50126614, 7.58834243, 8.0152366 , 8.41724414,
8.53631637, 8.57783649, 8.65372791, 8.77148382, 9.28239629,
9.64742165, 9.76102313, 9.93610952, 9.96946782, 10.17346814,
10.17346814])
a = 14.3
dpi = 2*np.pi*a
wlc = dpi/zo[1:]
fc = 300/wlc
print("a=",a,'[mm]')
print(f"{'No.':>7s}",f"{'mode':>7s}",f"{'xi':>10s}",\
f"{'wlc[mm]':>10s}",f"{'fc[GHz]':>10s}")
for i in range(1,nt):
print(f'{i:6.0f}',end=' ')
if t_[i]==1:
print(f'{"TE":2s}',end=' ')
elif t_[i]==0:
print(f'{"TM":2s}',end=' ')
print(f'{n_[i]:2.0f}',f'{m_[i]:2.0f}',f'{zo[i]:10.6f}',\
f'{wlc[i-1]:10.6f}',f'{fc[i-1]:10.6f}')
a= 14.3 [mm]
No. mode xi wlc[mm] fc[GHz]
1 TE 1 1 1.841184 48.799881 6.147556
2 TM 0 1 2.404826 37.362190 8.029508
3 TE 2 1 3.054237 29.418003 10.197837
4 TM 1 1 3.831706 23.448968 12.793740
5 TE 0 1 3.831706 23.448968 12.793740
6 TE 3 1 4.201189 21.386696 14.027412
7 TM 2 1 5.135622 17.495358 17.147406
8 TE 4 1 5.317553 16.896785 17.754857
9 TE 1 2 5.331443 16.852765 17.801234
10 TM 0 2 5.520078 16.276862 18.431071
11 TM 3 1 6.380162 14.082644 21.302818
12 TE 5 1 6.415616 14.004820 21.421197
13 TE 2 2 6.706133 13.398116 22.391208
14 TM 1 2 7.015587 12.807133 23.424447
15 TE 0 2 7.015587 12.807133 23.424447
16 TE 6 1 7.501266 11.977918 25.046089
17 TM 4 1 7.588342 11.840471 25.336830
18 TE 3 2 8.015237 11.209844 26.762193
19 TM 2 2 8.417244 10.674462 28.104462
20 TE 1 3 8.536316 10.525565 28.502034
21 TE 7 1 8.577836 10.474617 28.640666
22 TM 0 3 8.653728 10.382757 28.894061
23 TM 5 1 8.771484 10.243370 29.287238
24 TE 4 2 9.282396 9.679564 30.993131
25 TE 8 1 9.647422 9.313323 32.211920
26 TM 3 2 9.761023 9.204932 32.591225
27 TM 6 1 9.936110 9.042729 33.175824
28 TE 2 3 9.969468 9.012472 33.287205
29 TM 1 3 10.173468 8.831752 33.968344
30 TE 0 3 10.173468 8.831752 33.968344
これより，遮断周波数が最も低いTE$_{11}$モードが基本モードであることがわかる．
この計算例では，円形導波管の半径が $a=14.3$ [mm]の場合を求めており，この形状の導波管では，一つのモードのみが伝搬する周波数帯域は，6.15 GHzから8.03 GHzである．
ただし，円形導波管の形状が軸対称の変化のみであれば，基本モードに対して，TM$_{01}$，TE$_{21}$モードは発生しないため，第一次高次モードはTM$_{11}$モードとなり，シングルモードとなる周波数帯域は，6.15 GHz から12.79 GHzとなり，これはX帯（8.2-12.5GHz）に適した標準的な円形導波管といえよう．
円形導波管モードの遮断波数の計算
# imode = 1 (TE), -1 (TM)
def circular_waveguide_xi2(imode,m, n):
if imode == 1:
print('TE mode:',end=' ')
xi = jnp_zeros(m,n)
elif imode == -1:
print('TM mode:',end=' ')
xi = jn_zeros(m,n)
print(xi)
return xi[n-1]
circular_waveguide_xi2(imode=1, m=1, n=6)
TE mode: [ 1.84118378 5.33144277 8.53631637 11.7060049 14.86358863 18.01552786]
18.015527862681804
circular_waveguide_xi2(imode=-1, m=1, n=6)
TM mode: [ 3.83170597 7.01558667 10.17346814 13.32369194 16.47063005 19.61585851]
19.615858510468243
前のページに戻る
「目次」のページに戻る