昨日はカテナリになるのをみたけれど、ほぼ同じようなプログラムでFPUT問題（過去はFPU問題と言われていた）が計算できる。
やっているのはこちらと同じく、
d2xi/ dt2 = (xi+1 + xi-1 - 2xi) + α[(xi+1 - xi)2- (xi - xi-1)2]
の計算。通常はエネルギーのモードごとの図を描きますが、
実際の動きをGIFアニメにした。
こちら。
これだけ見てても再帰現象はわからんな、、、
Fermat Libraryでこういうのを見た。
あれ？カテナリ―（懸垂曲線）になるんじゃ？と思いますが、張力が場所によって違う効果を入れると
— tomo (@tonagai) September 27, 2020
d²y/dx²=ρg/T √(1+(dy/dx)²)
となってカテナリ(coshの形）になります。この話はほぼまっすぐになっていて√の中身が1になる近似を使っている。 https://t.co/gk3vjTDntA
ばねでつながれた質点系でカテナリになるのかどうか試してみようか、とやってみた。使うのはPython＋ScipyのDOP853.
計算する式はこちらで。
d2xi/ dt2 = -γdxi/dt + k(xi+1 + xi-1 - 2xi) -g
結果：黒いのがカテナリ。
おお、なるほどカテナリになってる。
IEEE Microwave magazineの特集はマイクロ波”顕微鏡”の特集。物質の構造だけでなく、生物や、SAWのようなアコースティックデバイスの表面まで見られるとか。これは今後も進む分野だと期待してます。
https://ieeexplore.ieee.org/xpl/mostRecentIssue.jsp?punumber=6668
一方、Microwave Journalはミリ波部品とパッケージ。今後ますますパッケージは重要になりますが、AiPやFOWLPよりPolyStrata押しの記事。
https://www.microwavejournal.com/publications/1/editions/275
Vayyarの4Dイメージングレーダの仕組みとは？
https://www.eetimes.com/under-the-hood-vayyaars-4d-imaging-radar/
MIPIのA-PHY（自動車向けインターフェース）の紹介
Microsoft Surface Duoの分解
https://www.ifixit.com/Teardown/Microsoft+Surface+Duo+Teardown/136576
ドローンで衛星通信のテストするとはなかなかうまいこと考えてるな。
https://www.microwavejournal.com/publications/1/editions/274
そしてアンリツの220GHzのテストシステム。これはいつ見ても高そう（というか私が使ったら壊しそう、、、）。
で、高周波になると様々な導体の粗さモデルを使ったりするが、それと因果律の関係。これはちゃんと考えたことなかったな。
Qualcommとワシントンポストの5Gマスタークラス。
まずはLesson1から。
https://www.washingtonpost.com/brand-studio/wp/2020/06/12/feature/5g-master-class/
続いてiFixitのSamsung Galaxy Watch3分解。Apple Watchと比べると集積度がとても低く見えるが、これで十分な機能という話もある。
https://www.ifixit.com/Teardown/Samsung+Galaxy+Watch3+Teardown/135891
今月のMicrowave Magazine。
https://ieeexplore.ieee.org/xpl/mostRecentIssue.jsp?punumber=6668
シリコンフォトニクスをマイクロ波、ミリ波通信に使う話。これは知らなかった（のと、SiPっていうのか、、、普通はSiPというとSystem in Packageですが、Silicon Photonics)。なかなか面白い。
あとダイレクトRFサンプリングの記事もあるが、コロナと5Gの健康問題についての話も。
電波塔とか燃やされてたな、、、でもこの記事の結論は、何かを判断するようなデータはない、というもの（肯定も否定もしてない）。
次はIMECの発表。Infineon(Google)のProject Soliとどう違うのかな？
Qorvoの5G向けフィルタのホワイトペーパー。
https://www.microwavejournal.com/articles/34267-a-new-generation-of-5g-filter-technology
Mixed Mode Sパラメータの解説。
SiTimeのMEMSクロックICの記事。
すぐiPhoneをパクるSamsungのUWB搭載。
https://www.engadget.com/samsung-galaxy-note-20-ultra-wide-band-connectivity-151223521.html
時間依存のシュレーディンガー方程式も空間微分を差分化すると普通の常微分方程式になるので、ルンゲクッタ8次のDormand-Prince法で計算できる。
具体的には空間の拡散項は2次の差分
(φi+1 - 2*φi + φi-1) / (dx*dx)
で近似。中央に井戸型ポテンシャルを置いて、左から波束が入ってくるときの計算。
ただし、複素数になるのでライブラリがodeからcomplex_odeに変わる。
ずっと前にExcel VBAでやったもののPythonでの書き直し。これは”計算物理”にでてたもので、もともとはシッフの例題にあるものだそうです。
クリックするとGIFアニメが始まります。
ソースはこんな感じで（インデントがないけどココログにコピペすると消えるので、、、）
%matplotlib notebook
import numpy as np
from scipy.integrate import complex_ode
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib.animation import FuncAnimation
from matplotlib.animation import PillowWriter
import numba
INTERVAL = 0.05 # [s] interval time
FRAME_INTERVAL = 1000 * INTERVAL # [msec] interval between frames
FPS = 1000/FRAME_INTERVAL # frames per second
Nt=6000
Nx=2000
t0=0
dx = 0.001
dt = dx * dx / 2
xmax = dx * float(Nx) / 2
xmin = -dx * float(Nx) / 2
tmax=Nt*dt
v_width = 0.064
v0 = 0.6*((70.7 * np.pi) **2)
deltax = 0.035
x0 = -0.3
p0 = np.pi * 50
t=np.linspace(0,tmax,Nt)
x=np.linspace(xmin,xmax,Nx)
v=np.empty(Nx,dtype=np.complex)
u0=np.empty(Nx,dtype=np.complex)
for i in range(Nx):
u0[i]=np.exp(-((x[i]-x0)**2) / (4 * deltax * deltax))*np.exp(p0*x[i]*1j)/((2 * np.pi * deltax * deltax) ** 0.25)
if x[i]>=-v_width/2 and x[i] <=v_width/2:
v[i]=v0
else:
v[i]=0
sol= np.empty((Nt, Nx),dtype=np.complex)
sol[0] = u0
def Schroedinger(t, u): #odeintのときとt,xの並びが逆
f=np.zeros(Nx, dtype=np.complex)
for i in range(Nx):
if i < Nx-1:
p1=i + 1
else:
p1=0
if i > 0:
m1=i - 1
else:
m1=Nx-1
diff2 = (u[m1] + u[p1] -2.0* u[i]) / ( dx ** 2)
f[i] = (-1j)*(-diff2 + v[i]*u[i])
return f
numba_f = numba.jit(Schroedinger,nopython=True)
solver=complex_ode(Schroedinger)
solver.set_integrator('dop853')
solver.set_initial_value(u0,t0) #なぜか関数と並びが逆
k=1
while solver.successful() and solver.t < tmax:
solver.integrate(t[k])
sol[k] = solver.y
k+= 1
def update(i):
plt.cla() #現在描写されているグラフを消去
plt.xlim(-1,1)
plt.ylim(0,4)
plt.plot(x,np.abs(sol[i]),c='Red',lw=3)
plt.plot(x,4*abs(v)/v0,c='blue',lw=3)
plt.xlabel("x")
plt.ylabel("Phi")
plt.title("Schroedinger equation by DOP853")
plt.grid()
# Plot
fig = plt.figure(figsize=(10,6))
ani = FuncAnimation(fig, update, frames=range(0,Nt,15), interval=FRAME_INTERVAL, blit=True)
plt.show()
ani.save("sh01.gif", writer = 'pillow',fps=FPS)
先日は蔵本・シバシンスキー方程式でしたが、今日は蔵本モデル。
dφi/dt = ωi + K/N Σsin (φj - φi ) （和はj=1 ～ N)
位相がどんどんそろっていくのがわかる。
クリックでGIFアニメが始まります。
蔵本先生と言えば、最近は千葉先生の研究で蔵本モデルが有名ですが、私にはこっちのほうがなじみがある蔵本・シバシンスキー方程式（Kuramoto-Sivashinsky equation）。1次元の乱流がシンプルに見られる。これはKdV方程式ができればすぐにできる。
∂u/∂t +u*∂u/∂x = -∂^2u/∂x^2 -∂^4/∂x^4
空間微分をO(h^6)という高次の差分で近似、かつ時間微分をルンゲクッタ8次のスキームであるDormand-Prince法(DOP853)で計算。
結果はこんな感じ。クリックでGIFアニメが始まります。（最初が動きがないので面白くないがだんだん動く）
ソースはKdVと同じなので関数部分だけ。
def kuramoto(t, u): #odeintのときとt,xの並びが逆
f=np.zeros(Nx)
for i in range(Nx):
if i < Nx-1:
p1 = i + 1
else:
p1 = 0
if i < Nx - 2:
p2 = i + 2
elif i == Nx - 2:
p2 = 0
else:
p2 = 1
if i < Nx - 3:
p3 = i + 3
elif i == Nx - 3:
p3 = 0
elif i == Nx - 2:
p3 = 1
else:
p3 = 2
if i > 0:
m1 = i - 1
else:
m1 = Nx-1
if i > 1:
m2 = i - 2
elif i == 1:
m2 = Nx - 1
else:
m2 = Nx - 2
if i > 2:
m3 = i - 3
elif i == 2:
m3 = Nx - 1
elif i == 1:
m3 = Nx - 2
else:
m3 = Nx - 3
d1 = (u[p3] / 9 - u[p2] + 5* u[p1] - 5 * u[m1] + u[m2] - u[m3] / 9) / (20 * dx / 3)
d2 = (2 * u[p3] - 27 * u[p2] + 270 * u[p1] - 490 * u[i] + 270 * u[m1] - 27 * u[m2] + 2 * u[m3])/(180 * dx * dx)
d4 = (-u[p3] +12*u[p2]-39*u[p1]+56*u[i]-39*u[m1]+12*u[m2]-u[m3])/(6*dx*dx*dx*dx)
f[i] = -u[i] * d1 - d2 - d4
return f
まずIEEE Microwave magazineの特集は、Quantum computing for microwave engineers。
https://ieeexplore.ieee.org/xpl/mostRecentIssue.jsp?punumber=6668
なるほどこういうのを見るとRFとも深くかかわるのが分かる。
次はMicrowave Journal。
Artificial Intelligence and Machine Learning Add New Capabilities to Traditional RF EDA Tools
LTEアンテナをAIで設計改善する例。
どちらも流行りものとRFのお話でちょっと面白いというかネタが切れているというか。
あと、TF-SAW vs. BAW: RF Companies Taking Different Strategies
についてはノーコメント。
Samsung Researchが6Gについてのホワイトペーパーだしてる。まあまだ先なので好き勝手に、、、
あとはこれに驚く、、、
KdV方程式
∂u/∂t +u*∂u/∂x+δ^2 ∂^3u/∂x^3 = 0
δ=0.022
ですが、これも空間を差分化すると普通の常微分方程式系と見なせる。そこでルンゲクッタ8次のDOP853を使おう。
空間差分はオリジナルのザブスキーとクルスカルのものと同じとします。
クリックでGIFアニメが始まります。
ソースコードはこんな感じ。
%matplotlib notebook
import numpy as np
from scipy.integrate import ode
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib.animation import FuncAnimation
from matplotlib.animation import PillowWriter
import numba
INTERVAL = 0.05 # [s] interval time
FRAME_INTERVAL = 1000 * INTERVAL # [msec] interval between frames
FPS = 1000/FRAME_INTERVAL # frames per second
Nt=512
Nx=256
xmax=2.0
xmin=0.0
tmin=0.0
tmax=16/np.pi
t0=0
dt = (tmax-tmin)/Nt
dx = (xmax - xmin) / Nx
d=0.022
t=np.linspace(0,tmax,Nt)
x=np.linspace(0,xmax,Nx)
u0=np.cos(np.pi*x)
sol= np.empty((Nt, Nx))
sol[0] = u0
def kdv(t, u): #odeintのときとt,xの並びが逆
f=np.zeros(Nx)
for i in range(Nx):
if i < Nx-1:
p1 = i + 1
else:
p1 = 0
if i > 0:
m1 = i - 1
else:
m1 = Nx-1
if i < Nx - 2:
p2 = i + 2
elif i == Nx - 2:
p2 = 0
else:
p2 = 1
if i > 1:
m2 = i - 2
elif i == 1:
m2 = Nx - 1
else:
m2 = Nx - 2
diff3 = (-u[m2] + 2 * u[m1] - 2 * u[p1] + u[p2]) / (2 * (dx ** 3))
diff1 = (-u[m1] + u[p1]) / (2 * dx)
f[i] = -(u[m1]+u[i]+u[p1]) * diff1/3 - d*d*diff3
return f
numba_f = numba.jit(kdv,nopython=True)
solver=ode(kdv)
solver.set_integrator('dop853')
solver.set_initial_value(u0,t0) #なぜか関数と並びが逆
k=1
while solver.successful() and solver.t < tmax:
solver.integrate(t[k])
sol[k] = solver.y
k+= 1
def update(i):
plt.cla() #現在描写されているグラフを消去
plt.xlim(0,2)
plt.ylim(-1,3)
plt.plot(x,sol[i],c='Red',lw=3)
plt.xlabel("t")
plt.ylabel("u")
plt.title("KdV equation by DOP853")
plt.grid()
# Plot
fig = plt.figure(figsize=(10,6))
ani = FuncAnimation(fig, update, frames=range(0,Nt,1), interval=FRAME_INTERVAL, blit=True)
plt.show()
ani.save("kdv.gif", writer = 'pillow',fps=FPS)
最近のコメント