たか松

理系大学4年です

【解析力学】ラグランジュ方程式の導出

takaphysics22.hatenablog.jp

↑前回導出した、ハミルトンの原理
\displaystyle \delta S=\int_{t_{1}}^{t_{2}}\delta L dt=0 \tag{1}
からオイラーラグランジュ方程式を導きます。

ラグランジアンLは運動エネルギーK\left(=\dfrac{1}{2}m\dot{\boldsymbol{x}}^2\right)とポテンシャルエネルギーUを用いて、
\begin{eqnarray}
L=K-U=\dfrac{1}{2}m\dot{\boldsymbol{x}}^2-U(\boldsymbol{x})
\end{eqnarray}
と表され、位置ベクトル\boldsymbol{x}=(x,y,z)と速度ベクトル\dot{\boldsymbol{x}}=(\dot{x},\dot{y},\dot{z})に依ります。

これより、仮想変位\delta\boldsymbol{x}とそのときの速度の変分\delta\dot{\boldsymbol{x}}(1)式で考えると
\begin{equation}
\displaystyle\delta\int_{t_{1}}^{t_{2}}L dt = \\
\int_{t_{1}}^{t_{2}}\left(\dfrac{\partial L}{\partial x}\delta x+\dfrac{\partial L}{\partial y}\delta y+\dfrac{\partial L}{\partial z}\delta z+\dfrac{\partial L}{\partial\dot{x}}\delta\dot{x}+\dfrac{\partial L}{\partial\dot{y}}\delta\dot{y}+ \dfrac{\partial L}{\partial \dot{z}}\delta\dot{z}\right)dt \\ \tag{2}
\end{equation}
となり、変分\delta\dot{\boldsymbol{x}}の部分を取り出すと
\displaystyle \int_{t_{1}}^{t_{2}}\left(\dfrac{\partial L}{\partial\dot{x}}\delta\dot{x}+\dfrac{\partial L}{\partial\dot{y}}\delta\dot{y}+ \dfrac{\partial L}{\partial \dot{z}}\delta\dot{z}\right)dt \\ =\displaystyle \int_{t_{1}}^{t_{2}}\left(\dfrac{\partial L}{\partial\dot{x}}\dfrac{d}{dt}\delta x+\dfrac{\partial L}{\partial\dot{y}}\dfrac{d}{dt}\delta y+\dfrac{\partial L}{\partial\dot{z}}\dfrac{d}{dt}\delta z\right)dt \\=\displaystyle\bigg\lbrack \dfrac{\partial L}{\partial\dot{x}}\delta x+\dfrac{\partial L}{\partial\dot{y}}\delta y+\dfrac{\partial L}{\partial\dot{z}}\delta z \bigg\rbrack_{t_{1}}^{t_{2}}-\int_{t_{1}}^{t_{2}}\left(\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{x}}\delta x+\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{y}}\delta y+\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{z}}\delta z\right)dt \\ =-\displaystyle\int_{t_{1}}^{t_{2}}\left(\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{x}}\delta x+\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{y}}\delta y+\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{z}}\delta z\right)dt

途中で部分積分を用いました。この結果を(2)式に戻し、整理すると
\begin{eqnarray}
\int_{t_{1}}^{t_{2}}\delta Ldt=\displaystyle\int_{t_{1}}^{t_{2}}\bigg\lbrace\left(\dfrac{\partial L}{\partial x}-\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{x}}\right)\delta x+\left(\dfrac{\partial L}{\partial y}-\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{y}}\right)\delta y+\left(\dfrac{\partial L}{\partial z}-\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{z}}\right)\delta z\bigg\rbrace dt=0
\end{eqnarray}

仮想変位\delta\boldsymbol{x}=(\delta x,\delta y,\delta z)は任意であるから、\displaystyle\int_{t_{1}}^{t_{2}}\delta Ldt=0が成り立つためには

x成分: \dfrac{\partial L}{\partial x}-\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{x}}=0

y成分: \dfrac{\partial L}{\partial y}-\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{y}}=0

z成分: \dfrac{\partial L}{\partial z}-\dfrac{d}{dt}\dfrac{\partial L}{\partial\dot{z}}=0

である必要があります。
これよりオイラーラグランジュ方程式が導出できました。

【解析力学】ハミルトンの原理



解析力学のハミルトンの原理

\delta S \displaystyle=\int_{t1}^{t2}L dt=\int_{t1}^{t2}\left(K-U\right) dt=0

を導出します。

ここでは、時刻t=t_{A}で点Aに存在した質点が、ニュートン運動方程式m\ddot{\boldsymbol{x}}=\boldsymbol{F}に従って運動し、時刻t=t_{B}に点Bに到達したとき、質点がたどる経路として実際に実現される運動を決める原理を導出することを試みます。

実際の経路は、仮想変位\delta \boldsymbol{x}にたいして
\begin{eqnarray}
\int_{t_{1}}^{t_{2}}\left(\boldsymbol{F}-m\ddot{\boldsymbol{x}}\right)\cdot \delta\boldsymbol{x}dt=\int_{t_{1}}^{t_{2}}0dt=0
\end{eqnarray}
が成り立つと考えられます。
左辺を
\begin{eqnarray}
\int_{t_{1}}^{t_{2}}\boldsymbol{F}\cdot \delta\boldsymbol{x}dt-\int_{t_{1}}^{t_{2}}m\ddot{\boldsymbol{x}}\cdot \delta\boldsymbol{x}dt=0
\end{eqnarray}
と書き直して、それぞれ別々に考えましょう。

左辺第一項

ここでは、\boldsymbol{F}が保存力である場合\left(\boldsymbol{F}=-\nabla U\right)であるものと仮定していますので、左辺第一項は
\begin{eqnarray}
\int_{t_{1}}^{t_{2}}\boldsymbol{F}\cdot \delta\boldsymbol{x}dt=-\int_{t_{1}}^{t_{2}}\nabla U\cdot \delta\boldsymbol{x}dt=-\int_{t_{1}}^{t_{2}}\delta U dt
\end{eqnarray}
ただし、\delta U=\displaystyle \nabla U \cdot \delta\boldsymbol{x}=\dfrac{\partial U}{\partial x}\delta x+\dfrac{\partial U}{\partial y}\delta y+\dfrac{\partial U}{\partial z}\delta z です。

左辺第二項

\begin{eqnarray}
\displaystyle \int_{t_{1}}^{t_{2}}m\dfrac{d^2\boldsymbol{x}}{d t^2} \cdot \delta\boldsymbol{x}dt &=& m\dfrac{d\boldsymbol{x}}{dt}\cdot\delta\boldsymbol{x}\Bigg|_{t_{1}}^{t_{2}}- \int_{t_{1}}^{t_{2}}m\dfrac{d\boldsymbol{x}}{dt}\cdot\dfrac{d(\delta\boldsymbol{x})}{dt}dt \\ &=& -\int_{t_{1}}^{t_{2}}m\dfrac{d\boldsymbol{x}}{dt}\cdot\dfrac{d(\delta\boldsymbol{x})}{dt}dt
\end{eqnarray}

ここで、ベクトルの関係式(\boldsymbol{A}+\boldsymbol{B})^2=A^2+B^2+2\boldsymbol{A}\cdot\boldsymbol{B}より、
\begin{eqnarray}
\boldsymbol{A}\cdot\boldsymbol{B}=\dfrac{(\boldsymbol{A}+\boldsymbol{B})^2-A^2-B^2}{2}
\end{eqnarray}
\boldsymbol{A}=\dfrac{d\boldsymbol{x}}{dt},\boldsymbol{B}=\dfrac{d(\delta\boldsymbol{x})}{dt}とそれぞれ代入すれば、
\begin{eqnarray}
-\int_{t_{1}}^{t_{2}}m\dfrac{d\boldsymbol{x}}{dt}\cdot\dfrac{d(\delta\boldsymbol{x})}{dt}dt=-\dfrac{1}{2}m\int_{t_{1}}^{t_{2}}\bigg\lbrace\dfrac{d}{dt}\left(\boldsymbol{x}+\delta\boldsymbol{x}\right)\bigg\rbrace^2-\left(\dfrac{d\boldsymbol{x}}{dt}\right)^2-\bigg\lbrace\dfrac{d(\delta\boldsymbol{x})}{dt}\bigg\rbrace^2 dt
\end{eqnarray}
ここで、2次の微小量\bigg\lbrace\dfrac{d(\delta\boldsymbol{x})}{dt}\bigg\rbrace^2を無視しました。
すると、
\begin{eqnarray}
-\int_{t_{1}}^{t_{2}}m\dfrac{d\boldsymbol{x}}{dt}\cdot\dfrac{d(\delta\boldsymbol{x})}{dt}dt=-\int_{t_{1}}^{t_{2}}\dfrac{1}{2}m\Bigg\lbrack\bigg\lbrace\dfrac{d}{dt}(\boldsymbol{x}+\delta\boldsymbol{x})\bigg\rbrace^2-\bigg(\dfrac{d}{dt}\boldsymbol{x}\bigg)^2\Bigg\rbrack dt=-\int_{t_{1}}^{t_{2}}\delta K dt
\end{eqnarray}
よって、左辺第二項\displaystyle \int_{t_{1}}^{t_{2}}m\dfrac{d^2\boldsymbol{x}}{d t^2} \cdot \delta\boldsymbol{x}dt
\begin{eqnarray}
\displaystyle \int_{t_{1}}^{t_{2}}m\dfrac{d^2\boldsymbol{x}}{d t^2} \cdot \delta\boldsymbol{x}dt=-\int_{t_{1}}^{t_{2}}\delta K dt
\end{eqnarray}
と変形することができました。

以上をまとめると、
\begin{eqnarray}
-\int_{t_{1}}^{t_{2}}\delta U dt+\int_{t_{1}}^{t_{2}}\delta K dt=0
\end{eqnarray}
となり、これを変形すると
\begin{eqnarray}
\int_{t_{1}}^{t_{2}}(\delta K-\delta U) dt=\int_{t_{1}}^{t_{2}}\delta(K-U)dt=\delta\int_{t_{1}}^{t_{2}}(K-U)dt
\end{eqnarray}
を表すことがわかります。ここで、L=K-Uと定義すれば(Lラグランジアンと呼ばれる量です。)
\begin{eqnarray}
S &=& \int_{t_{1}}^{t_{2}}L dt \\
\therefore \delta S &=& \delta \int_{t_{1}}^{t_{2}}L dt=0
\end{eqnarray}
これより、積分S極値をとる経路が実現されることを示すハミルトンの原理を導出することができました。

【統計力学】1次元 ゴム弾性統計力学モデル


温度Tにおいて1本の長さがaである棒がN(N\gg1)個つながったゴム鎖が、右向きまたは左向きで水平に連なっている。
このときゴム鎖の両端に働く張力fを求めます。

ボルツマン定数k_{B}および
公式 \ln N!=N\ln N-Nを用います。
構成は以下の通りです。

  1. 場合の数 W
  2. 右向きの本数n_{1}、左向きの本数n_{2}
  3. エントロピーS
  4. 自由エネルギーF
  5. 張力f

1.場合の数 W

ゴム鎖の場合の数はN本の中から任意のn_{1}本を選び出すことに等しいので、
\begin{eqnarray}
W= {}_N C_{n_{1}}=\dfrac{N!}{n_{1}!n_{2}!}
\end{eqnarray}
と求められます。

2.右向きの本数n_{1}、左向きの本数n_{2}

全体の本数NN=n_{1}+n_{2}であり、ゴム鎖の長さll=(n_{1}-n_{2})dなのでこの2式から右向きの本数n_{1}、左向きの本数n_{2}はそれぞれ、
\begin{eqnarray}
n_{1} &=& \dfrac{N}{2}+\dfrac{l}{2a}=\dfrac{N}{2}\left(1+x\right)\\
n_{2} &=& \dfrac{N}{2}-\dfrac{l}{2a}=\dfrac{N}{2}\left(1-x\right) \\
\end{eqnarray}
となります。
以下x=\dfrac{l}{Na}を用いることにします。

3.エントロピーS

1で求めた、場合の数Wに対数をとり、
\begin{eqnarray}
\ln W &=& \ln N!-\ln n_{1}!-\ln n_{2}! \\&=& (N\ln N-N)-(n_{1}\ln n_{1}-n_{1})-(n_{2}\ln n_{2}-n_{2}) \\ &=& n_{1}\ln \dfrac{N}{n_{1}}+n_{2}\ln \dfrac{N}{n_{2}} \\ &=& -n_{1}\ln \dfrac{n_{1}}{N}-n_{2}\ln \dfrac{n_{2}}{N}
\end{eqnarray}
2のn_{1},n_{2}を用い、さらにボルツマン定数k_{B}をかけると
エントロピーS
\begin{eqnarray}
S &=& k_{B}\ln W=-k_{B}\left(n_{1}\ln \dfrac{n_{1}}{N}+n_{2}\ln \dfrac{n_{2}}{N}\right) \\
&=& -k_{B}\biggl\{\dfrac{N}{2}(1+x)\ln\left(\dfrac{1+x}{2}\right)+\dfrac{N}{2}(1-x)\ln\left(\dfrac{1-x}{2}\right)\biggr\} \\ &=& Nk_{B}\biggl\{\ln2-\dfrac{1}{2}(1+x)\ln(1+x)-\dfrac{1}{2}(1-x)\ln(1-x)\biggr\}
\end{eqnarray}

4.自由エネルギーF

熱力学関係式F=U-TSを用いてHelmholtzの自由エネルギーFを求めます。
ここで、内部エネルギーUは棒が自由に折れ曲がることにより、ゴム鎖全体の長さlには依らないことがわかります。よって内部エネルギーは一定値とすることができるので、その値をU=0とすると自由エネルギーF
\begin{eqnarray}
F(x)=-TS=-Nk_{B}T\biggl\{\ln2-\dfrac{1}{2}(1+x)\ln(1+x)-\dfrac{1}{2}(1-x)\ln(1-x)\biggr\}
\end{eqnarray}
と求められます。

5.張力f

張力fは等温環境下において、ゴム鎖の全長lに対する自由エネルギーFの変化率で与えられるので、
\begin{eqnarray}
f=\left(\dfrac{\partial F}{\partial l}\right)_{T}=\dfrac{1}{Na}\dfrac{dF(x)}{dx}\Bigg|_{T} 
\end{eqnarray}
4で求めた自由エネルギーF(x)の表式を代入し、整理すると
張力f
\begin{eqnarray}
f &=&\dfrac{1}{Na}\Bigg\lbrack-Nk_{B}T\biggl\{-\dfrac{1}{2}\ln(1+x)-\dfrac{1}{2}+\dfrac{1}{2}\ln(1-x)+\dfrac{1}{2}\biggr\}\Bigg\rbrack \\
&=& \dfrac{k_{B}T}{2a}\ln\left(\dfrac{1+x}{1-x}\right) \\
&=& \dfrac{k_{B}T}{2a}\ln\left(\dfrac{1+l/Na}{1-l/Na}\right)
\end{eqnarray}
と求めることができました。

【統計力学】3次元 状態密度

統計力学で登場する状態密度D(\epsilon)を導出します。
今回は3次元です。

3次元において波数ベクトル\bf k(大きさ:k=\sqrt{\dfrac{2m\epsilon}{\hbar^2}})で張られる空間内の要素数\left(2\pi/L\right)^3であり、またその面積は\dfrac{4}{3}\pi k^3であるので、\bf k空間内の状態数W
\begin{eqnarray}
W=2\dfrac{\dfrac{4}{3}\pi k^3}{\left(2\pi/L\right)^3}=\dfrac{L^3}{4\pi^2}\dfrac{4}{3}\left(\dfrac{2m\epsilon}{\hbar^2}\right)^{3/2}
\end{eqnarray}
求めた状態数Wをエネルギー\epsilon微分することにより状態密度D(\epsilon)が求められるので、
\begin{eqnarray}
D(\epsilon)  &=& \dfrac{d}{d\epsilon}W  \\&=&  \dfrac{d}{d\epsilon}\dfrac{L^3}{4\pi^2}\dfrac{4}{3}\left(\dfrac{2m\epsilon}{\hbar^2}\right)^{3/2} \\ &=& \dfrac{V}{2\pi^2}\left(\dfrac{2m}{\hbar^2}\right)^{3/2}\epsilon^{1/2}
\end{eqnarray}
3次元の状態密度はエネルギー\epsilon1/2乗に比例することがわかります。

【統計力学】2次元 状態密度

統計力学で登場する状態密度D(\epsilon)を導出します。
今回は2次元です。(前回の一次元の続きです)

2次元において波数ベクトル\bf k(大きさ:k=\sqrt{\dfrac{2m\epsilon}{\hbar^2}})で張られる空間内の要素数\left(2\pi/L\right)^2であり、またその面積は\pi k^2であるので、\bf k空間内の状態数W
\begin{eqnarray}
W=2\dfrac{\pi k^2}{\left(2\pi/L\right)^2}=\dfrac{L^2}{2\pi}\dfrac{2m}{\hbar^2}\epsilon
\end{eqnarray}
求めた状態数Wをエネルギー\epsilon微分することにより状態密度D(\epsilon)が求められるので、
\begin{eqnarray}
D(\epsilon)=\dfrac{d}{d\epsilon}W=\dfrac{d}{d\epsilon}\dfrac{L^2}{2\pi}\dfrac{2m}{\hbar^2}\epsilon=\dfrac{L^2}{2\pi}\dfrac{2m}{\hbar^2}
\end{eqnarray}
2次元の状態密度はエネルギーの大きさに依らないことがわかります。

【統計力学】1次元 状態密度

統計力学で登場する状態密度D(\epsilon)を導出します。
まずは1次元から考えます。

 

長さLの1次元の空間にある自由粒子シュレディンガー方程式は、エネルギーを\epsilon波動関数\psi(x)として
\begin{eqnarray}
-\dfrac{{\hbar}^2}{2m}\psi(x) &=& \epsilon\psi(x) \\
\therefore\dfrac{\partial^2}{\partial x^2}\psi(x) &=& -k^2\psi(x)  
\end{eqnarray}
ただし、k=\sqrt{\dfrac{2m\epsilon}{{\hbar}^2}}で、このkは波数を表します。

波数ベクトル\bf{k}(大きさ:k)で張られる体積要素は2\pi/Lであることから、\bf{k}空間での長さkにある状態数W
\begin{eqnarray}
W=2\dfrac{k}{(2\pi/L)}=\dfrac{L}{\pi}\sqrt{\dfrac{2m\epsilon}{{\hbar}^2}}
\end{eqnarray}
ただし、2はスピン自由度に由来する因子です。
これより、状態密度D(\epsilon)は状態数Wをエネルギー\epsilon微分することにより求められるので、
\begin{eqnarray}
D(\epsilon)=\dfrac{d}{d\epsilon}W=\dfrac{d}{d\epsilon}\dfrac{L}{\pi}\sqrt{\dfrac{2m\epsilon}{{\hbar}^2}}=\dfrac{L}{2\pi}\sqrt{\dfrac{2m}{{\hbar}^2}}{\epsilon}^{-\frac{1}{2}}
\end{eqnarray}
1次元の状態密度はエネルギー\epsilonの-1/2乗に比例することがわかります。

【統計力学】1次元の量子力学的調和振動子

振動数\nuをもつ1次元の量子力学調和振動子のエネルギー準位は
\begin{eqnarray}
E_{n}=h\nu\left(n+\dfrac{1}{2}\right) \tag{1}
\end{eqnarray}
となります。この調和振動子が温度Tで熱平衡状態にあるとき、
平均エネルギー\langle E \rangleと比熱Cを求めます。

まず、分配関数Z(1)を用いることにより、
\begin{eqnarray}
Z &=& \sum_{n=0}^{\infty}e^{-\beta E_{n}} \\&=& 
\sum_{n=0}^{\infty}\exp(-\beta\left(n+\dfrac{1}{2}\right)h \nu) \\&=&
e^{-\beta h\nu/2}\sum_{n=0}^{\infty}e^{-n\beta h\nu}
\end{eqnarray}
と表すことができますが(\beta=\frac{1}{k_{B}T}、k_{B}ボルツマン定数)、
ここでe^{-\beta h\nu} \ll 1であることから、
無限等比級数の公式\sum_{n=0}^{\infty}x^n=\frac{1}{1-x}を利用すると
\begin{eqnarray}
Z=\dfrac{e^{-\beta h\nu/2}}{1-e^{-\beta h\nu}}
\end{eqnarray}
と求めることができます。

これより、平均エネルギー\langle E \rangleは、
\begin{eqnarray}
\langle E \rangle &=& -\dfrac{\partial}{\partial\beta}\ln Z
\\&=& -\dfrac{\partial}{\partial\beta}\left(-\beta h\nu/2-\ln (1-e^{-\beta h\nu})\right) \\&=& \dfrac{h\nu}{2}+\dfrac{h\nu e^{-\beta h\nu}}{1-e^{-\beta h\nu}}
\\&=& \dfrac{h\nu}{2}+\dfrac{h\nu}{e^{\beta h\nu}-1}
\\&=& \dfrac{h\nu}{2}\dfrac{e^{\beta h\nu}+1}{e^{\beta h\nu}-1}
\end{eqnarray}
となります。
比熱CC=\dfrac{\partial E}{\partial T}で求められますが、\beta=\dfrac{1}{k_{B}T}の変数変換に注意し、
\begin{eqnarray}
C &=& \dfrac{\partial E}{\partial T} \\&=& -k_{B}\beta^2\dfrac{\partial E}{\partial\beta}
\\&=& -k_{B}\beta^2\dfrac{\partial}{\partial\beta}\left(\dfrac{h\nu}{2}+\dfrac{h\nu}{e^{\beta h\nu}-1}\right)
\\&=& -k_{B}\beta^2\dfrac{-(h\nu)^2 e^{\beta h\nu}}{(e^{\beta h\nu}-1)^2}
\\&=& k_{B}(\beta h\nu)^2\dfrac{e^{\beta h\nu}}{(e^{\beta h\nu}-1)^2}
\end{eqnarray}
とわかります。