nano_exit

基礎的なことこそ、簡単な例が必要だと思うのです。

pythonでペンローズタイリング

ペンローズタイリングとは、ペンローズが開発した非周期的な図形の敷き詰め方である。
以下のサイトで、ペンローズタイリングの方法と、そのpythonスクリプトが公開されている。
Penrose Tiling Explained

図の描画には、pycairoが使われていて、事前にインストールする必要がある。(pip で普通に入ったので、特に問題はないと思う。)
また、スクリプト内には xrange が使われているが、python3では入っていないので、range に置き換える必要がある。

実際に使って、図を作ってみた。
f:id:koideforest:20181113220049p:plain

デフォルトでは色の設定が赤と青だが、自分で色を変えることが出来る。
また図形の細かさも 変数 NUM_SUBDIVISIONS を変更することで調整出来る。

操作として、時計回りに注目する三角形を分割していくのだが、デフォルトでは二つで一組になるように、偶数番目の三角形は奇数番目に鏡映操作を施したものになっている。
その操作を off にして計算したものがこちら。
f:id:koideforest:20181113221427p:plain

自分で実際に試せるところが、pythonの魅力だなぁと改めて感じた。

円運動

原点を中心とした円の軌跡 \vec{r}は以下のように記述出来る。

\displaystyle
\vec{r} = ( x, y ) = ( r \cos\theta, r \sin\theta ) = r ( \cos\theta, \sin\theta ) = r \hat{ r },
\\
\displaystyle
r = \left| \vec{ r } \right| = \sqrt{ x^2 + y^2 } = r \sqrt{ \cos^2 \theta + \sin^2 \theta } = r
\\
\displaystyle
\hat{ r } = \frac{ \vec{ r } }{ \left| \vec{ r } \right| } = ( \cos\theta, \sin\theta )

原点周りの円周上の運動を考える。そのため rは定数で時間変化しないとし、角度 \theta \equiv \theta(t)の時間微分を角速度 \dot{ \theta}として定義する。(円運動を「角度のみが変化する運動」と言い換えても良いだろう。)

運動方程式の用途として、以下の二つが考えらえる。

  1. 力 → 軌跡( \vec{r}
  2. 軌跡( \vec{r})→ 力

通常は 1. のプロセスについて使われることが多い(「力がわかっていて、その場合の軌跡を求めよ」的な)。
ここでは 2. のプロセスを辿り、円運動を起こすために必要な力について考える。

まずは軌跡 \vec{r}を時間微分する。

\displaystyle
\vec{ v } = \frac{ d }{ dt } \vec{ r } = r( - \dot{ \theta } \sin \theta, \dot{ \theta } \cos \theta ) = r \dot{ \theta } ( -\sin \theta, \cos \theta ) = v \hat{ v }
\\
\displaystyle
v = \left| \vec{ r } \right| = r \dot{ \theta } \sqrt{ \sin^2 \theta + \cos^2 \theta } = r \dot{ \theta }
\\
\displaystyle
\hat{ v } = \frac{ \vec{ v } }{ \left| \vec{ v } \right| } = ( -\sin \theta, \cos \theta )

ここで、 \hat{ r } \hat{ v }の間の角度がどうなっているかを調べるために、内積を取ると、

\displaystyle
\hat{ r } \cdot \hat{ v }= -\cos \theta \sin \theta + \sin \theta \cos \theta = 0
したがって、 \vec{ r } \vec{ v }は直交している。
もう少し具体的に言えば、 \hat{ v }は円弧上の点 \vec{ r }における円に対する接線方向を表している。
また、角度変化 \Delta \thetaにおいて r \, \Delta \thetaが円弧の長さを表すから、 v = r \dot{ \theta }は円弧上で軌跡が動く速さを表している。

運動方程式にするために、もう一度時間微分する。

\displaystyle
\frac{ d }{ dt } \vec{ v } = r \ddot{ \theta } \hat{ v } + r \dot{ \theta } ( - \dot{ \theta } \cos \theta, - \dot{ \theta } \sin \theta )
  = r \ddot{ \theta } \hat{ v } + r \left( \dot{ \theta } \right)^2 ( - \cos \theta, - \sin \theta )
\\
\displaystyle
  = r \ddot{ \theta } \hat{ v } + r \left( \dot{ \theta } \right)^2 \left( - \hat{ r } \right)

したがって、運動方程式

\displaystyle
m \frac{ d }{ dt } \vec{ v } = m r \ddot{ \theta } \hat{ v } + m r \left( \dot{ \theta } \right)^2 \left( - \hat{ r } \right)
  \equiv \vec{ f }_v + \vec{ f }_{ -r } 
\\
\displaystyle
\vec{ f }_v = m r \ddot{ \theta } \hat{ v }
\\
\displaystyle
\vec{ f }_{ -r } = m r \left( \dot{ \theta } \right)^2 \left( - \hat{ r } \right)
と表せる。
ここから、円運動をさせるためには、接線方向の力 \vec{ f }_vと、原点方向の力(求心力) \vec{ f }_{ -r } が必要なことがわかる。

特別な場合として、角速度が一定( \dot{ \theta } = \omega  \ddot{ \theta } = 0 )な等速円運動を考える。
この時、 \vec{ f }_v = m r \ddot{ \theta } \hat{ v } = 0であるから、「求心力 \vec{ f }_{ -r } = mr \omega^2 \left( - \hat{ r } \right) のみが働く時、等速円運動が起こる」と言える。

個人的には、求心力の大きさ f_{-r} = m r \omega^2 がいつも覚えられなくて困っていたが、

  •  r極座標にして大きさと方向を分けた名残
  •  \omega^2:二回時間微分した名残

と思えば、いくらか頭に残りそうな気がしている。

頂角が小さい二等辺三角形の底辺について

回転運動の説明を読んでいる時に、角度が小さい時の差ベクトルの近似について気になったので考えて見た。
等速円運動 [物理のかぎしっぽ]

以下の図のような二等辺三角形における底辺(青)、円弧(赤)、そして垂線(緑)を考える。
f:id:koideforest:20181111195750p:plain
図では頂角は30度でプロットしてある。

問題によっては、赤を近似するために青にしたり、その逆もあったり、何を基準にするかややこしい面があるように思う。
とりあえず、それぞれを式で表す。二等辺三角形の頂角および辺の長さを \theta, rとすると

\displaystyle
(red) = r \theta
\\
\displaystyle
(blue) = r \sin( \theta / 2 ) + r \sin( \theta / 2 ) = 2 r \sin( \theta / 2 )
\\
\displaystyle
(green) = r \sin( \theta )

青を表現するのに \sin関数を使うため、比較として、三角形における \sin関数の最も直接的な利用である垂線(緑)を登場させた。
底辺と垂線は、次のようにして差を評価出来る。

\displaystyle
(green) = r \sin( \theta ) = 2 r \sin( \theta / 2 ) \cos( \theta / 2 ) < (blue)
角度が非常に小さい( \theta \ll 1 )時、 \cos( \theta / 2 ) \approx 1 となって、青と緑は(ほぼ)一致する。

赤と青は、 \sin関数のTaylor展開から \theta << 1で両者が一致するのがすぐわかる。

\displaystyle
(blue) = 2 r \sin( \theta / 2 ) \approx 2 r \left( \theta / 2 \right) = r \theta = (red)

よって、最終的には赤青緑が全て近似的に等しくなる。

図を作るのに使ったスクリプト

from matplotlib import pyplot as plt
import numpy as np
import math

fig = plt.figure()
ax = fig.add_subplot( 1, 1, 1 )

d = 30. # degree
radians = np.array([ math.radians(d) for d in np.linspace( 0, d, 100 ) ])
c = np.cos( math.radians(d) )
s = np.sin( math.radians(d) )

# triangle
tri = plt.Polygon( ( ( 0, 0 ), ( 1, 0 ), ( c, s ) ),
                   color = "black",
                   alpha=0.1
      )
ax.add_patch( tri )

# arc
ax.plot( np.cos( radians ), np.sin( radians ), color = "red" )

# base line of triangle
xs = np.linspace( c, 1., 100 )
ys = ( ( s - 0. ) / ( c - 1. ) ) * ( xs - c ) + s 
ax.plot( xs, ys, color = "blue" )

# vertical line in triangle
ax.vlines( c, 0., s, color = "green" )

plt.tick_params( left='off', right='off', bottom='off', top='off',
                 labelleft='off', labelright='off', labelbottom='off', labeltop='off' )

plt.savefig( "triangle.png" )

三角形および補助線の書き方は以下を参照
Python/Matplotlibで三角形の描画
【Python@matplotlib】matplotlib にて横、縦の補助線を描く方法について - Qiita

剛体の運動:ヨビノリの動画の補足

剛体は結局学部ではやらず、研究室のゼミでやるわけでもなく、独学だったのでヨビノリで復習してみた。
多体問題として剛体を見ると、かなり面白いと感じざるを得ない。

ここでは自分の理解のために、そこで出てくる数学的事項等を自分の好みに合わせて補足してみた。

【大学物理】剛体の力学入門①(特徴と魅力)/全6回【力学】 - YouTube

  • 相対位置が不変

相対位置と言ってしまうとと、ベクトルが変わったらアウトなので、並進のみが許され、回転は出来ないと思われる。
なので、ここでは「相対『距離』が不変」の方が個人的にはしっくり来た。

  • 重心について

当たり前のように、我々は「重心」というワードが普段使っているが、よく考えると「重心」は使われ方が曖昧な言葉で、また物理で意味する重心は「剛体の運動方程式」によって意味を持つようになる(定義される)。
例えば、数学で「三角形の重心」と呼ばれる点は、正確には「幾何中心」であり、また三角形の各点がそれぞれ違う重さを持っているとすると、物理的な意味の重心(質量中心)は幾何中心からズレる。
そのため、既に「重心」を知っている体で剛体の話を始めると、後々「重心って何だ?」という罠にハマる気がした。
重心 - Wikipedia

【大学物理】剛体の力学入門②(重心の運動)/全6回【力学】 - YouTube

  • 作用反作用の法則から、相互作用項が落ちる。


\displaystyle
\sum_i \sum_j' \vec{f}_{ij} = \sum_{i \neq j } \vec{f}_{ij} = \sum_{ i > j } \left( \vec{f}_{ij} + \vec{f}_{ji} \right) = \sum_{ i > j } \left( \vec{f}_{ij} - \vec{f}_{ij} \right) = 0

【大学物理】剛体の力学入門③(回転運動)/全6回【力学】 - YouTube

  • 相互作用項が外積の場合でも落ちる。


\displaystyle
\sum_i \sum_j' \left( \vec{r}_i \times \vec{f}_{ij} \right) = \sum_{i \neq j } \left( \vec{r}_i \times \vec{f}_{ij} \right) = \sum_{ i > j } \left( \vec{r}_{i} \times \vec{f}_{ij} + \vec{r}_{j} \times \vec{f}_{ji} \right)
\\
\displaystyle
= \sum_{ i > j } \left( \vec{r}_{i} - \vec{r}_{j} \right) \times \vec{f}_{ij} = \sum_{ i > j } \left| \vec{r}_{i} - \vec{r}_{j} \right| \hat{ f }_{ij} \times \vec{f}_{ij} = 0
\\
\displaystyle
\left( \hat{ f }  = \vec{ f } / | \vec{ f } | \right)
(「位置ベクトルの差」と「力の方向」が並行というのは、式だけ弄っていると気付きにくいので、こういうところは物理やっている感じがして個人的に好き。)

  • 重心からの相対位置 \vec{r}'の満たす性質 \sum_i m_i \vec{r}'_i = 0について

これを、ある種の重心の定義としても良いと思うので、定理として強調した方が個人的には好き。

【大学物理】剛体の力学入門④(慣性モーメントの定義)/全6回【力学】 - YouTube

  • 固定軸を入れることで、自由度が6から1に減る件について

個人的には、最初の相対距離不変のところで出した「自由度6」の話をもう一度出して繋げた方が分かりやすいと思った。

  1. 一点目:自由度3
  2. 二点目:自由度2
  3. 三点目:自由度1
  4. 四点目以降:自由度0

固有軸を入れることで、一点目と二点目が固定されてしまう(自由度0)ので、自由度が1だけ残る。

  • 慣性モーメントを出す時に、角度は積分変数なのに、角速度は積分の外に出せる件について

正直、この部分の説明が無かったために、この記事を書いてる感はある。
角度は各点で異なるが、相対距離が不変に保たれていると、角速度が位置に依らないことを示す。
動画と同様に円筒座標を取り、その時の二点 \vec{r}_1, \vec{r}_2を考える。
z軸周りの回転のみを扱うと、時間に依存するのはxy平面内の角度 \varphiのみだから、

\displaystyle
\left| \vec{r}_1 - \vec{r}_2 \right| = c \rightarrow \left( \vec{r}_1 - \vec{r}_2 \right)^2 = c^2 \rightarrow \frac{ d }{ dt } \left( \vec{r}_1 - \vec{r}_2 \right)^2 = 0
\\
\displaystyle
\left( \vec{r}_1 - \vec{r}_2 \right)^2 = ( x_1 - x_2 )^2 + ( y_1 - y_2 )^2 + ( z_1 - z_2 )^2
\\
\displaystyle
  = r^2 ( \cos \varphi_1 - \cos \varphi_2 )^2 + r^2 ( \sin \varphi_1 - \sin \varphi_2 )^2 + ( z_1 - z_2 )^2
\\
\displaystyle
  = r^2 ( \cos^2 \varphi_1 + \sin^2 \varphi_1 ) + r^2 ( \cos^2 \varphi_2 + \sin^2 \varphi_2 ) -2 r^2 ( \cos \varphi_1 \cos \varphi_2 + \sin \varphi_1 \sin \varphi_2 ) + ( z_1 - z_2 )^2
\\
\displaystyle
  = 2r^2 -2 r^2 ( \cos \varphi_1 \cos \varphi_2 + \sin \varphi_1 \sin \varphi_2 ) + ( z_1 - z_2 )^2
\\
\displaystyle
  = 2r^2 -2 r^2 \cos ( \varphi_1 - \varphi_2 ) + ( z_1 - z_2 )^2
\\
\displaystyle
\therefore
2 r^2 \sin ( \varphi_1 - \varphi_2 ) ( \dot{ \varphi}_1 - \dot{\varphi}_2 ) = 0
特別な二点を選ばない限り \sin( \varphi_1 - \varphi_2 ) \neq 0 は成り立たない。
したがって、 \dot{ \varphi}_1 = \dot{\varphi}_2 、つまり角速度 \dot{ \varphi }は固定軸を通すと位置に依らなくなる。
これにより、角速度を積分の外に出しても問題ない。

【大学物理】剛体の力学入門⑤(慣性モーメントの性質)/全6回【力学】 - YouTube

  • 原点の取り方について、

hを「任意の回転軸」と「それに並行で重心を通る回転軸」との距離とすると、 I = I_G + Mh^2として、慣性モーメントが求まるが、この時に原点は「任意の回転軸」上のどこかに置かないと答えが変わってしまう(or 計算するのが面倒臭い)ため、原点の取り方が地味に重要。

【大学物理】剛体の力学入門⑥(力学的エネルギー)/全6回【力学】 - YouTube

  • 慣性モーメントが質量 m r^2を掛けたような形になる件について

いつも慣性モーメントの形( r^2を掛けなくては行けない理由)を忘れてしまうので、原点の並進運動と比較している時に言及する方が個人的に好み。
回転運動は、速度のうち角速度のみで十分なので、円弧上の速度までは必要ない。そのため運動エネルギーにおいて、円弧から動径を抜き出した「角速度の二乗」で表そうとすれば、抜き出された「動径の二乗」が質量項にくっつき、これが慣性モーメントとなると思えば覚え易いと思った。

sympyで(平面の)Greenの定理を確認

sympyの練習を兼ねて、平面に対するGreenの定理で、問題を解いてみる。
平面のグリーンの定理 [物理のかぎしっぽ]

以下の積分を求めてみる。

\displaystyle
\int \! \int_D ( 1 - ( x^2 + y^2) ) dx dy
\\
\displaystyle
D = \left\{ (x,y) | x^2 + y^2 \le 1 \right\}
原点で最大値を取り、等方向的で、境界でゼロを持つような、何かしらの密度を積分する、と思えば解りやすいだろうか。

普通に極座標に変換して解くと、面積素片

\displaystyle
dx \, dy
=
\left|
\begin{array}{cc}
  \frac{ \partial x }{ \partial r } & \frac{ \partial x }{ \partial \theta } \\ 
  \frac{ \partial y }{ \partial r } & \frac{ \partial y }{ \partial \theta } 
\end{array}
\right|
dr \, d\theta
=
\left|
\begin{array}{cc}
  \cos \theta & - r \sin \theta \\ 
  \sin \theta &  r \cos \theta
\end{array}
\right|
dr \, d\theta
=
r\, dr \, d\theta
であるから(ヤコビアンの中身は転置してても大丈夫)、

\displaystyle
\int \! \int_D ( 1 - ( x^2 + y^2 ) ) dx dy
  = \int^1_0 r dr \int^{2\pi}_0 \, d\theta \, ( 1 - r^2 )
\\
\displaystyle
  = 2 \pi \left[ \frac{ r^2 }{ 2 } - \frac{ r^4 }{ 4 } \right]^1_0
  = \frac{ \pi }{ 2 }

この積分に対する一連の操作をsympyでやると、

from sympy import *
init_session()

r, theta, c, s = symbols( "r θ c s" )
c = cos( theta )
s = sin( theta )

xrt = r * c
yrt = r * s

f0 = 1 - ( x**2 + y**2 )
f1 = f0.subs( [ ( x, xrt ), ( y, yrt) ] )

Jacob = Matrix([
  [ Derivative( xrt, r ), Derivative( xrt, theta ) ],
  [ Derivative( yrt, r ), Derivative( yrt, theta ) ]
])
Jacob = simplify( Jacob.det().doit() )

I = Integral( f1 * Jacob, ( r, 0, 1 ), ( theta, 0, 2 * pi ) )
I.doit()
# pi / 2

積分区間の決定は自分で判断する必要があるが、大方、機械的に処理できるのではないだろうか。
途中式を吐き出すようにすれば、もう少し分かったような気がするようになると思う。

注意点として、*.subs()で関係式を代入して式変形する際、Derivativeにそれで代入しようとしても、「積分した後の変数に代入しようとする」ため、例えば以下のようにしようとするとゼロ行列になる。

Jacob = Matrix([
  [ Derivative( x, r ), Derivative( x, theta ) ],
  [ Derivative( y, r ), Derivative( y, theta ) ]
])
Jacob = Jacob.subs([ ( x, xrt ), ( y, yrt ) ]).doit()
# zero matrix

Integral()も同様に、*.subs()を使っても積分範囲に代入されて、被積分関数を弄ることが出来ない。

この積分を、Greenの定理で線積分に直して計算してみる。
Greenの定理は以下の式で与えられる。

\displaystyle
\int \! \int_D \left( \frac{ \partial Q }{ \partial x } - \frac{ \partial P }{ \partial y } \right) dx dy
  =  \oint_C \left( P \, dx + Q \, dy \right)

したがって、例えば以下のように P, Qを定めると、

\displaystyle
\frac{ \partial Q }{ \partial x } = \frac{ 1 }{ 2 } - x^2
\rightarrow Q = \frac{ x }{ 2 } - \frac{ x^3 }{ 3 }
\\
\displaystyle
  - \frac{ \partial P }{ \partial y } = \frac{ 1 }{ 2 } - y^2
\rightarrow P = - \left( \frac{ y }{ 2 } - \frac{ y^3 }{ 3 } \right)
\\
\displaystyle
\therefore
\oint_C \left( P \, dx + Q \, dy \right) = 
  - \oint_C \left( \frac{ y }{ 2 } - \frac{ y^3 }{ 3 } \right) dx
  + \oint_C \left( \frac{ x }{ 2 } - \frac{ x^3 }{ 3 } \right) dy
\\
\displaystyle
= \int^{2\pi}_0 \left( \frac{ \sin \theta }{ 2 } - \frac{ \sin^3 \theta }{ 3 } \right) \sin \theta d\theta
  + \oint_C \left( \frac{ \cos \theta }{ 2 } - \frac{ \cos^3 \theta }{ 3 } \right) \cos \theta d\theta

これを計算するとちゃんと \pi/2になることを手で示しても良いが、計算ミスだらけで効率が悪いので、sympyにやってもらうことにする。

P = - integrate( 1/2 - y**2, y )
integrand = P.subs([ ( x, xrt ), ( y, yrt ) ]) * diff( xrt, theta )
integrand = integrand.subs([ ( r, 1 ) ])
int_range = ( theta, 0, 2 * pi )
P_int = integrate( integrand, int_range )
# pi / 4

Q = integrate( 1/2 - x**2, x )
integrand = Q.subs([ ( x, xrt ), ( y, yrt ) ]) * diff( yrt, theta )
integrand = integrand.subs([ ( r, 1 ) ])
int_range = ( theta, 0, 2 * pi )
Q_int = integrate( integrand, int_range )
# pi / 4

I = P_int + Q_int
# pi / 2

 O_x, P_yによる被積分関数の分割を、上の例と逆にしても、同じ答えが得られる。

P = - integrate( 1/2 - x**2, y )
integrand = P.subs([ ( x, xrt ), ( y, yrt ) ]) * diff( xrt, theta )
integrand = integrand.subs([ ( r, 1 ) ])
int_range = ( theta, 0, 2 * pi )
P_int = integrate( integrand, int_range )
# pi / 4

Q = integrate( 1/2 - y**2, x )
integrand = Q.subs([ ( x, xrt ), ( y, yrt ) ]) * diff( yrt, theta )
integrand = integrand.subs([ ( r, 1 ) ])
int_range = ( theta, 0, 2 * pi )
Q_int = integrate( integrand, int_range )
# pi / 4

I = P_int + Q_int
# pi / 2

もっと言えば、 Q_x被積分関数を全部押し付けてもちゃんと求まる。

P_int = 0

Q = integrate( 1 - ( x**2 + y**2 ), x )
integrand = Q.subs([ ( x, xrt ), ( y, yrt ) ]) * diff( yrt, theta )
integrand = integrand.subs([ ( r, 1 ) ])
int_range = ( theta, 0, 2 * pi )
Q_int = integrate( integrand, int_range )
# pi / 2

I = P_int + Q_int
# pi / 2

今はsympyで数値計算的なチェックになってしまっているが、代入を駆使して途中式のチェックとかが出来るとかなり便利だとは思う。
(ただ、手で書くのとどっちが速いか、という問題はあるが。。。)

単振り子の問題がパッと受け入れられない件について

単振り子の振動の振舞を求める問題は、単振動の典型的な例題の一つであるが、何かしっくり来ない。
単振り子 ■わかりやすい高校物理の部屋■
単振り子を徹底解説!近似の使用法&運動方程式から周期を導出する方法
単振り子:運動方程式

そのため、しっくり来ない点について考えてみた。

  1. 考えなきゃいけない力が何となく残っている気がする。
  2. 円弧の長さがパッと思い付かない。
  3. 求める軌跡が直交座標 x, yではなく、円弧 s上の運動。
  4. 円弧に対する接線方向と角度の関係がパッとわからない。
  5. 単振動の場合、振動が微小(角度 \thetaが微小)として \sin関数を線形近似するのがズルく感じる。
  6. 当たり前の様に T = 2 \pi / \omegaと説明されるが、スッとそこに行けない。
  7. 単振動を超えた一般の運動方程式を扱うときに、回転運動の話(求心力・遠心力)が急に入って来る。

1. 考えなきゃいけない力が何となく残っている気がする。

例えば、天井と糸の間に働く力などは普通は省略して図が描かれている。これらの力は、作用・反作用で釣り合っているため、運動方程式中では打ち消しあってゼロになり反映されない。
しかし、よくよく考えると、図示出来ない力があることに気が付く。それは、「天井が受ける作用(今の場合、糸に引っ張られる力)を打ち消すための力」である。
例えば、物を平らな床の上に置いたとする。物が動かないのは、重力と床が押し返す反作用が釣り合っているからと説明される。しかし、床は作用を受けっ放しであり、そのまま真に受けると床は沈んでしまう。
実際には、床の下の木材やら基礎やら何やらに全部力が伝わって、最終的には地面に力が伝わることになる。そうすると、次の疑問として、「地面を押したとき、地球を反対側から支える力はないから、(どんなに小さくても)原理的には地球は動くのか?」ということになる。
答えは(あまりに小さ過ぎて観測できないが)"YES"である。
地球上の人間が一斉にジャンプするとどうなるのか : Spookie's
したがって、安心して振り子の運動方程式に重要な力のみを考えれば良い。

2. 円弧の長さがパッと思い付かない。

半径 l の円において、ラジアン単位の角度 \thetaに対する円弧の長さ sは、「 s = l \theta 」で与えられる。
これがしっくり来なかった。
原因は、「円周=円周率 \times 直径 =  \pi \times 2 l 」で覚えていたからだと気が付いた。
そうではなく、「円周 =  2 \pi \times 半径 =  2 \pi \times l 」だから、任意の角度にすぐに一般化して「円弧=角度 \times 半径 =  \theta l 」と理解できる。

3. 求める軌跡が直交座標 x, yではなく円弧 s上の運動。

例えば、それまでバネの単振動とかやっていると、一次元方向の xだったり、鉛直方向の zだったりして分かり易かった。
しかし、振り子になって問題が二次元になり、力の方向が軌跡に沿っていないため、軌跡を何に設定するかという点において、大きな飛躍を感じてしまう。そもそも軌跡を直交座標ではなく何故「基準点から測った円弧の長さ」に設定するのか?という点に対し、「そうすれば簡単に解ける」というようなノリが多い気がする。
個人的には、そこをもう少し掘り下げて、「知りたいのは単振り子の角度 \theta (t)で、その角度と円弧が \theta = s / lで一対一対応するから」とかの説明が欲しかった。

4. 円弧に対する接線方向と角度の関係がパッとわからない。

これが正直最も大きい問題かも知れない。
毎回、図示して確認するが、接線方向の成分が \sinなのか \cosなのか、いまいち自信が持てない。
以下のサイトの図で見れば、同位角の関係だけ使っているので、まだわかる気がする。
単振り子:運動方程式

5. 単振動の場合、振動が微小(角度 \thetaが微小)として \sin関数を線形近似するのがズルく感じる。

角度 \thetaが微小な時、 \sin\theta \approx \theta = s / l として近似することで、単振動の運動方程式が得られる。
多分、「単振り子の問題を解きます!」ってノリで近似が入ると、急に各論っぽくなって脳の負荷が上がる気がする。
そうではなく、先に単振動の運動方程式の一般形を持って来ておいて、「単振り子の中に隠れている単振動を見つける」というような体であれば、目的がハッキリして近似も受け入れ易くなる気がする。

6. 当たり前の様に T = 2 \pi / \omegaと説明されるが、スッとそこに行けない。

個人的には、振動が \sin( \omega t )の形で表されるから、 \omega T = 2 \pi と書かれた方が分かりやすい。
 \omegaが時間の逆数だから、逆数の逆数になるとちょっと脳の負荷が上がってパッとわからない。

7. 単振動を超えた単振り子の運動方程式を扱うときに、回転運動の話(求心力・遠心力)が急に入って来る。

回転運動を学ぶのが後の方だったりするときに、先取りして単振り子の一般的な運動方程式を扱おうとすると、「単振り子ってよく分かんない」というイメージが付いてしまう気がする。
個人的には、回転運動を勉強した後に単振り子の運動方程式を出して、角度の小さい時には単振動をすることを示した方がスッキリして良い気がするが、全体のカリキュラム的には振動の方が先だから、単振り子の単振動としての役割は小さいかもしれない。


結局のところ、こんな風に自分で何が気に食わないかを考えることが、問題を受け入れることの最短経路なんだろうなぁ悟った次第である。

量子力学でのガリレイ変換

特殊相対論(古典力学)の導入でガリレイ変換とローレンツ変換ローレンツブースト)を比較することが多いが、そもそも量子力学でのガリレイ変換ってなんだ?と思い、簡単に考察。
以下のpdfを参考にした。
http://cat.phys.s.u-tokyo.ac.jp/lecture/QM_1_11/quantum1.pdf
https://ocw.kyoto-u.ac.jp/ja/09-faculty-of-engineering-jp/quantum-theory-for-electrical-and-electronic-engineering/pdf/chap11.pdf

一次元を考える。
静止系 Zに対して、速度vで等速度運動している運動系 Z'を定義する。

\displaystyle
x = x' + vt'
\\
t = t'
例えば、静止系 Zにおいて x=vtで運動しているものは、運動系 Z'では x'=0、つまり止まって見える。
このように、運動物体に合わせた運動系 Z'を取ると、問題が簡単になる。
やりたいことは、「(簡単に解ける) Z'  \rightarrow (実際に観測する) Z」の変換である。

時間に依存しないシュレーディンガー方程式において、時間依存性は無視してしまいがちだが、上記の通りガリレイ変換で位置と時間が混ざるので、時間依存性が発生する。
したがって、時間依存性を含めた形で波動関数 \varphi(x,t)を書いておいた方が良い。

ガリレイ変換によって、微分演算は次のように変換される。

\displaystyle
\frac{ d }{ d x' } = \frac{ d }{ d x } \frac{ d x }{ d x' } + \frac{ d }{ d t } \frac{ d t }{ d x' } = \frac{ d }{ d x }
\\
\displaystyle
\frac{ d }{ d t' } = \frac{ d }{ d x } \frac{ d x }{ d t' } + \frac{ d }{ d t } \frac{ d t }{ d t' } = v \frac{ d }{ d x } + \frac{ d }{ d t }

ここで、質点力学と場の理論の間の注意について触れたい。

質点力学では、質点の位置 r, r'を「軌跡」という形で物体の運動を表していた。
つまり、求めたい質点の位置 r, r'時間だけの関数であり、座標の位置 x, x'の関数ではない。
これによって、ニュートン方程式は、力をポテンシャル Uの勾配を使って表せば、

\displaystyle
m \frac{ d^2 }{ d t'^2 } r'( t ) = - \frac{ d U'( x', t' ) }{ d x' }
\\
\displaystyle
m \left( \frac{ d^2 }{ d t^2 } + 2 v \frac{ d }{ d t } \frac{ d }{ d x } + \frac{ d^2 }{ d x^2 } \right) ( r( t ) - vt ) = - \frac{ d U'( x - vt, t ) }{ d x }
\\
\displaystyle
\therefore m \left( \frac{ d^2 }{ d t^2 } \right) r( t ) = - \frac{ d U( x, t ) }{ d x }
\\
\displaystyle
\left( U( x, t ) \equiv U'( x - vt, t ) = U'( x', t' ) \right)
というように、軌跡が時間にしか依存していないことと、等速直線運動が時間に線形であることによって、ガリレイ変換に対して方程式が不変になっている。

一方で、電磁気学量子力学は質点の物理ではなく、各座標の上で何かしらの特徴(密度など)が定義されている。これが場の理論と言われる由縁である。
そのため、座標系が変わっても、それらの間で正しく同じ場所を指し示していれば、得られる値は(直感的には)同じであると考えられる。
上のニュートン方程式の例で言えば、ポテンシャル U,U'の変換がそれに当たる。ややこしいが、ポテンシャルは場の量(位置、時間に依存した量)なので、上記の変換は妥当である。
波動関数で表現すれば、運動系 Z'波動関数 \varphi'が先に求まっているとすると、静止系 Z波動関数 \varphi

\displaystyle
\varphi( x, t ) = \varphi'( x', t' ) = \varphi'( x - vt, t ) \quad (?)
ということが予想される。

しかし、運動量(とエネルギー)に着目すると、これだけでは上手く行かないことがわかる。
それは、Z,Z'の間をただ移す操作を行っても、その座標上の運動量は変わってくれないため、運動座標系 Z'の持つ運動量が合成されない。
平面波に対して確認してみると、規格化定数を Nとして

\displaystyle
\varphi_{k}( x, t ) = \varphi'_k(x',t')
  = N e^{i( k(x - vt) - \omega t )}
  = N e^{i( kx - ( \omega + kv ) t )} \quad (?)
となり、エネルギーのみが変化してしまう全く不満足な結果が得られる。
そのため、座標系における運動量の見え方の違いをどこに求めるか?ということになる。

量子力学では、波動関数は観測量ではなく、その二乗である「密度」が観測量である。
そのため、波動関数そのものの一致を要請するのは条件が厳し過ぎるので、密度が一致するように変更する。
それは、位相因子の付与を許容することと等価である。

\displaystyle
  | \varphi( x, t ) |^2 = | \varphi'( x', t' ) |^2
\\
\displaystyle
\therefore \varphi( x, t ) e^{ i f( x, t ) } = \varphi'( x', t' )
\\
\displaystyle
\left( f( x, t ) = \Re f( x, t ) \right)
仮に、 f複素数の場合、密度の一致を満たさなくなるので、 fは実数であることが必要。

この時に、シュレーディンガー方程式がガリレイ変換を満たすように fを定める。

\displaystyle
  -\frac{ \hbar^2 }{ 2m } \frac{ d^2 }{ d x'^2 } \varphi' + V'( x', t' ) \varphi' = i \hbar \frac{ d }{ d t' } \varphi'
\\
\displaystyle
  -\frac{ \hbar^2 }{ 2m } \frac{ d^2 }{ d x^2 } e^{i f} \varphi + V( x, t ) e^{i f} \varphi = i \hbar \left( v \frac{ d }{ d x } + \frac{ d }{ d t } \right) e^{ i f } \varphi

それぞれの微分係数は、

\displaystyle
\frac{ d^2 }{ d x^2 } e^{i f} \varphi
  = \frac{ d }{ d x } \left( i f_x e^{i f} \varphi + e^{i f} \varphi_x \right)
  = i f_{xx} e^{i f} \varphi - f_x^2 e^{i f} \varphi + 2 i f_x e^{i f} \varphi_x + e^{i f} \varphi_{xx}
\\
\displaystyle
\frac{ d }{ d t } e^{i f} \varphi
  = i f_t e^{i f} \varphi + e^{i f} \varphi_t

元々の静止系 Zシュレーディンガー方程式は、

\displaystyle
  -\frac{ \hbar^2 }{ 2m } \varphi_{xx} + V( x, t ) \varphi = i \hbar \varphi_t

したがって、 e^{-if}を左から掛けて簡単にすると、

\displaystyle
  -\frac{ \hbar^2 }{ 2m } \left( i f_{xx} \varphi - f_x^2 \varphi + 2 i f_x \varphi_x \right)
  = i \hbar \left( v \left( i f_x \varphi + \varphi_x \right) + i f_t \varphi \right)

 \varphi \varphi_xの係数がそれぞれゼロになるように条件を課すと、まず \varphi_xから f_xが求まる。

\displaystyle
  \left( -\frac{ \hbar^2 }{ 2m } \left( 2 i f_x \right) - i \hbar v \right) \varphi_x  = 0
  \rightarrow f_x = - \frac{ m v }{ \hbar }

したがって、 xの二階微分 f_{xx}はゼロになり、簡単に解ける。

\displaystyle
  \left( -\frac{ \hbar^2 }{ 2m } \left( i f_{xx} - f_x^2 \right) - i \hbar \left( i v f_x + i f_t \right) \right) \varphi
  = 0
\\
\displaystyle
  \left( -\frac{ \hbar^2 }{ 2m } \left( - f_x^2 \right) - i \hbar \left( i v f_x + i f_t \right) \right) \varphi
  = 0
\\
\displaystyle
  \left( -\frac{ \hbar^2 }{ 2m } \left( - \frac{ m^2 v^2 }{ \hbar^2 } \right) - i \hbar \left( - i \frac{ m v^2 }{ \hbar } + i f_t \right) \right) \varphi
  = 0
\\
\displaystyle
  \left( \frac{ m v^2 }{ 2 } - m v^2 + \hbar f_t \right) \varphi
  = 0
\\
\displaystyle
\therefore f_t = \frac{ m v^2 }{ 2 \hbar }

よって、位相因子は、積分定数を除いて、

\displaystyle
f = - \frac{ 1 }{ \hbar } \left( m v x - \frac{ m v^2 }{ 2 } t \right)
と求まる。

得られた結果を平面波に対して確認してみると、

\displaystyle
\varphi_{k}( x, t ) = e^{ - i f( x, t ) } \varphi'_k(x',t')
  = N e^{ \frac{i}{\hbar} \left( m v x - \frac{ m v^2 }{ 2 } t \right) } e^{i( k(x - vt) - \omega t )}
\\
\displaystyle
  = N \, {\rm exp} \left( i \frac{ \hbar k + mv }{ \hbar } x - i \left( \omega + kv + \frac{ m v^2 }{ 2 \hbar } \right) t \right)
\\
\displaystyle
  = N \, {\rm exp} \left( \frac{ i }{ \hbar } \left( \hbar k + mv \right) x - \frac{ i }{ \hbar } \frac{ (\hbar k + mv )^2 }{ 2 m } t \right)
きちんと運動量が合成されているのがわかる。

ではこの位相変化はどこから来るのか?
参考にしたpdfに明快な考察があった。
https://ocw.kyoto-u.ac.jp/ja/09-faculty-of-engineering-jp/quantum-theory-for-electrical-and-electronic-engineering/pdf/chap11.pdf
元々、シュレーディンガー方程式は波動関数に対する方程式であり、それを密度の一致に条件を緩めた時点で、ガリレイ変換に対して不変な方程式ではないことがわかる。
観測者間の変換問題を満たすためには、特殊相対論に拡張してローレンツ変換を考える必要がある。
電子の場合には、運動系に対するディラック方程式の解を、ローレンツ変換を使って静止系に繋げば良い。これによって、位相も含めてダイレクトに波動関数同士を繋げられる。

\displaystyle
\Phi_D( x, t ) = \Phi'_D( x', t' )
運動エネルギーが静止エネルギーよりも十分小さい時、静止エネルギーによる位相を括り出すことができる。

\displaystyle
\varphi_D( x, t ) e^{ - i \omega_c t } = \varphi'_D( x', t' ) e^{ - i \omega_c t' }
\\
\displaystyle
\therefore \varphi_D( x, t )  e^{ - i \omega_c ( t- t' ) } = \varphi'_D( x', t' )
\\
\displaystyle
\left( \omega_c = m c^2 / \hbar \right)

ローレンツ変換は以下で与えられる。

\displaystyle
\beta = v / c
\\
\displaystyle
x' = \gamma( x - \beta c t )
\\
\displaystyle
ct' = \gamma( ct - \beta x )
\\
\displaystyle
\gamma = \frac{ 1 }{ \sqrt{ 1 - \beta^2 } } \sim ( 1 + \beta^2 / 2 + \cdots )

したがって、 \betaの二次まで拾えば、

\displaystyle
  -\omega_c ( t- t' )
  \sim -\omega_c \left( t -
  \left( 1 + \frac{ \beta^2 }{ 2 } + \cdots \right) 
  \left( t - \frac{ \beta }{ c } x \right)
  \right)
\\
\displaystyle
  \sim - \omega_c \left(
  \frac{ \beta }{ c } x - \frac{ \beta^2 }{ 2 } t
  \right) 
\\
\displaystyle
  = - \frac{ 1 }{ \hbar } \left( mvx - \frac{ m v^2 }{ 2 } t \right) = f( x, t )
となり、ガリレイ変換で出てきた位相変化を再現する。
つまり、非相対論であるシュレーディンガー方程式でも、静止エネルギーによる影響は無視出来ないレベルで入って来ていてることがわかる。
相対論による僅かな時間のズレが、静止エネルギーを通して、運動量および運動エネルギーのズレに変換されるのは、なかなかに衝撃だと思う。