構造力学メモ（その９）２次元極座標の重調和方程式の一般解

　前回は、コンクリートの割裂引張試験の計算式の理論の裏に、エアリーの応力関数や2次元の極座標、重調和方程式などがあることを紹介した。その中で、重調和方程式の解は天下り的に与える形でしか紹介しなかった。今回は、重調和方程式の一般解を紹介する。
　Wikipediaによれば、このような2次元極座標の重調和方程式の一般解はミッシェル解と呼ばれるらしい。
　単純に見えるコンクリートの割裂引張試験の裏にこれほどの理論的な話があるのは驚くべきことだと思う。

１．2次元の極座標の重調和方程式

　2次元の極座標のエアリーの応力関数と重調和方程式は、以下の通りである。

$$
\begin{equation*}
\left(\frac{\partial^2}{\partial r^2}
+\frac{1}{r}\frac{\partial }{\partial r}
+\frac{1}{r^2}\frac{\partial^2}{\partial \theta^2}\right)^2\phi=0
\end{equation*}
$$

　上記の重調和方程式は、参考文献によれば、以下のように分解して、2行目を解いて、その解を1行目に代入して解くことができると紹介されている。

$$
\begin{align*}
\left(\frac{\partial^2}{\partial r^2}
+\frac{1}{r}\frac{\partial }{\partial r}
+\frac{1}{r^2}\frac{\partial^2}{\partial \theta^2}\right)\phi&=\psi\\
\left(\frac{\partial^2}{\partial r^2}
+\frac{1}{r}\frac{\partial }{\partial r}
+\frac{1}{r^2}\frac{\partial^2}{\partial \theta^2}\right)\psi&=0\\
\end{align*}
$$

　その解き方はいまいちよく分からないので、そのまま解くことを考える。上記の重調和方程式は、以下のようにも書くことができる。

$$
\begin{align*}
\left\{\frac{1}{r}\frac{\partial }{\partial r}
\left(r\frac{\partial}{\partial r}\right)
+\frac{1}{r^2}\frac{\partial^2}
{\partial \theta^2}\right\}^2\phi &=0\\
\left[\frac{1}{r^2}\left\{r\frac{\partial }{\partial r}\left(r\frac{\partial}{\partial r}\right)+\frac{\partial^2}{\partial \theta^2}\right\}\right]^2\phi &=0
\end{align*}
$$

以下のように変数変換をするともう少し扱いやすくなる。
（本当は、$${(r,\theta)\rightarrow (t,\theta_1)}$$といったように、ちゃんと2変数の変数変換を行うべきだが、各変数の混じりあった変数変換ではないので、下記のように1変数だけの変数変換を行っている）。

$$
\begin{align*}
r &=e^t \\
t &=\ln r \\
\frac{\partial}{\partial r}
&=\frac{\partial t}{\partial r}\frac{\partial}{\partial t}
=\frac{1}{r}\frac{\partial}{\partial t} \\
r\frac{\partial}{\partial r}&=\frac{\partial}{\partial t}
\end{align*}
$$

上記を重調和方程式に代入すると、以下のように多少簡単になる。

$$
\begin{align*}
\left[\frac{1}{r^2}\left\{\frac{\partial }{\partial t}
\left(\frac{\partial}{\partial t}\right)
+\frac{\partial^2}{\partial \theta^2}\right\}\right]^2\phi &=0\\
\left[\frac{1}{r^2}\left(\frac{\partial^2}{\partial t^2}
+\frac{\partial^2}{\partial \theta^2}\right)\right]^2\phi &=0\\
\left[e^{-2t}\left(\frac{\partial^2}{\partial t^2}
+\frac{\partial^2}{\partial \theta^2}\right)\right]^2\phi &=0\\
e^{-2t}\left(\frac{\partial^2}{\partial t^2}
+\frac{\partial^2}{\partial \theta^2}\right)
e^{-2t}\left(\frac{\partial^2}{\partial t^2}
+\frac{\partial^2}{\partial \theta^2}\right)\phi &=0
\end{align*}
$$

真ん中の$${e^{-2t}}$$が邪魔だが、演算子の交換は以下のように注意が必要である。
感覚的に言えば、$${e^{-2t}}$$を左側に動かす時に、-2が出てくる
（$${e^{-2t}}$$でなくても、$${e^{at}}$$も同様に左側に動かす時に$${a}$$が出てくる）。

$$
\begin{align*}
\frac{\partial}{\partial t}e^{-2t}f
&=e^{-2t}\frac{\partial f}{\partial t}-2e^{-2t}f\\
&=e^{-2t}\left(\frac{\partial}{\partial t}-2\right)f\\
\frac{\partial^2}{\partial t^2}e^{-2t}f
&=e^{-2t}\left(\frac{\partial}{\partial t}-2\right)^2f
\end{align*}
$$

上記を代入すると、以下のようになる。

$$
\begin{align*}
e^{-4t}\left\{\left(\frac{\partial}{\partial t}-2\right)^2
+\frac{\partial^2}{\partial \theta^2}\right\}
\left(\frac{\partial^2}{\partial t^2}
+\frac{\partial^2}{\partial \theta^2}\right)\phi &=0\\
\left\{\left(\frac{\partial}{\partial t}-2\right)^2
+\frac{\partial^2}{\partial \theta^2}\right\}
\left(\frac{\partial^2}{\partial t^2}
+\frac{\partial^2}{\partial \theta^2}\right)\phi &=0
\end{align*}
$$

２．重調和方程式を解く

　重調和方程式の解の候補として、$${\phi=F(t)e^{\mu \theta}}$$を考えて、代入してみる。

$$
\begin{align*}
\left\{\left(\frac{\partial}{\partial t}-2\right)^2
+\frac{\partial^2}{\partial \theta^2}\right\}
\left(\frac{\partial^2}{\partial t^2}
+\frac{\partial^2}{\partial \theta^2}\right)
F(t)e^{\mu \theta} &=0 \\
\left\{\left(\frac{\partial}{\partial t}-2\right)^2
+\mu^2\right\}
\left(\frac{\partial^2}{\partial t^2}+\mu^2\right)
F(t)e^{\mu \theta}&=0
\end{align*}
$$

さらに、$${F(t)=e^{\lambda t}}$$を代入してみる。

$$
\begin{align*}
\left\{\left(\frac{\partial}{\partial t}-2\right)^2
+\mu^2\right\}
\left(\frac{\partial^2}{\partial t^2}
+\mu^2\right)
e^{\lambda t}e^{\mu \theta} &=0 \\
\left\{\left(\lambda-2\right)^2
+\mu^2\right\}
\left(\lambda^2+\mu^2\right)
e^{\lambda t}e^{\mu \theta}&=0
\end{align*}
$$

逆に、$${\phi=e^{\lambda t}G(\theta)}$$から代入すると、

$$
\begin{align*}
\left\{\left(\frac{\partial}{\partial t}-2\right)^2
+\frac{\partial^2}{\partial \theta^2}\right\}
\left(\frac{\partial^2}{\partial t^2}
+\frac{\partial^2}{\partial \theta^2}\right)
e^{\lambda t}G(\theta) &=0 \\
\left\{\left(\lambda-2\right)^2
+\frac{\partial^2}{\partial \theta^2}\right\}
\left(\lambda^2+\frac{\partial^2}{\partial \theta^2}\right)
e^{\lambda t}G(\theta)&=0
\end{align*}
$$

となるが、結局は$${G(\theta)=e^{\mu \theta}}$$を代入すれば、得られる4次方程式は同じになる。

３．4次方程式を解く

　上記の4次方程式はうまく因数分解できないので、下記のように$${\mu}$$は虚数と考える（実は、$${\mu}$$を実数と考えることもできるが、その場合は、$${\lambda}$$を虚数とすることになって$${e^{\lambda t}=r^{\lambda}}$$となって、$${r}$$の虚数乗が出てくるので、物理的な解ではないと考えられる）。

$$
\begin{equation*}
\mu^2=-n^2, \mu=\pm i n
\end{equation*}
$$

　このとき、$${n}$$は$${-\infty \le n \le \infty}$$の実数とする。

$$
\begin{align*}
\left\{\left(\lambda-2\right)^2-n^2\right\}
\left(\lambda^2-n^2\right)
e^{\lambda t}e^{i n \theta}&=0\\
(\lambda-2-n)(\lambda-2+n)(\lambda-n)(\lambda+n) &=0
\end{align*}
$$

　上記を満たす$${e^{\lambda t}e^{in\theta}}$$が、重調和方程式の解となり、その重ね合わせで一般解を求めることができる。
　4次方程式の解は、

$$
\begin{align*}
\lambda_1 &=2+n,& \lambda_2 &=2-n,
& \lambda_3 &=n, & \lambda_4 &=-n
\end{align*}
$$

　よって、解は、これらの$${e^{\lambda t}e^{\pm i n\theta}}$$が解となり、そして一般解は、これらの解の重ね合わせとなる。$${n}$$が整数とは限らないが、実際の境界条件などから$${n}$$が整数になったりする。そもそもこのような解の重ね合わせは、フーリエ変換で解を表現したともいえる。
　二行目への変形は、$${e^{\pm i \theta}}$$だと複素数なので、$${\cos\theta=}$$ $${(e^{i\theta}+e^{-i\theta})/2}$$,$${\sin\theta=}$$ $${(e^{i\theta}-e^{-i\theta})/2}$$を使って変形している。

$$
\begin{align*}
\phi &=\sum_{n=-\infty}^{\infty} (A_n' e^{(n+2)t}e^{in\theta}
+B_n' e^{(-n+2)t}e^{in\theta}+C_n' e^{nt}e^{in\theta}
+D_n' e^{-nt}e^{in\theta})\\
&=\sum_{n=0}^{\infty} (A_n r^{n+2}\cos n\theta
+B_n r^{-n+2}\cos n\theta+C_n r^{n}\cos n\theta
+D_n r^{-n}\cos n\theta)\\
&+\sum_{n=0}^{\infty} (E_n r^{n+2}\sin n\theta
+F_n r^{-n+2}\sin n\theta+G_n r^{n}\sin n\theta
+H_n r^{-n}\sin n\theta)
\end{align*}
$$

４．重解を求める

　求めたかった解$${r\theta\sin\theta}$$は、出てきていない。
　これは、4次方程式を解くときに重解を考えていなかったからである。
　微分方程式を解く際に4次方程式の解が重解になるときは注意が必要である。重解をもつ時、$${e^{\lambda t}}$$だけでなく$${te^{\lambda t}}$$や$${\theta e^{in \theta}}$$も解になる。

　重解となるのは、$${\lambda}$$の4次方程式と考えた時、すなわち、$${\phi=F(t)e^{\pm in\theta}}$$として、$${F(t)}$$を求める場合、

$$
\begin{align*}
(\lambda_1,\lambda_2,\lambda_3,\lambda_4) &=(2,2,0,0),&(n=0)\\
(\lambda_1,\lambda_2,\lambda_3,\lambda_4) &=(3,1,1,-1),&(n=\pm 1)
\end{align*}
$$

$${n=0}$$のとき$${\lambda^2(\lambda-2)^2=0}$$　なので

$$
\begin{equation*}
\frac{\partial^2}{\partial t^2}
\left(\frac{\partial}{\partial t}-2\right)^2 F(t)=0
\end{equation*}
$$

これらの方程式は、以下のように解くことができる。

$$
\begin{align*}
\frac{\partial^2}{\partial t^2}F(t) &=0\\
F(t) &=c_1+c_2t\\
\left(\frac{\partial}{\partial t}-2\right)^2 F(t) &=0\\
F(t) &=c_1e^{2t}+c_2te^{2t}
\end{align*}
$$

　これらの解を重ね合わせた一般解は以下のようになる。

$$
\begin{align*}
\phi &=A_0 e^{0\cdot i \theta}
+B_0 t e^{0\cdot i \theta}
+C_0 e^{2t}e^{0\cdot i \theta}
+D_0 t e^{2t}e^{0\cdot i \theta} \\
&=A_0+B_0 t+C_0 e^{2t}+D_0 t e^{2t}\\
&=A_0+B_0 \ln r+C_0 r^2+D_0 r^2\ln r
\end{align*}
$$

$${n=\pm 1}$$のとき$${(\lambda+1)(\lambda-1)^2(\lambda-3)=0}$$　なので

$$
\begin{align*}
\phi &=A_1' e^{-t}e^{\pm i\theta}
+B_1' e^{t}e^{\pm i \theta}
+C_1' te^{t}e^{\pm i \theta}
+D_1' e^{3t}e^{\pm i \theta} \\
&=A_1'e^{-t}e^{\pm i\theta}+B_1' e^t e^{\pm i\theta}
+C_1' te^te^{\pm i\theta}+D_1' e^{3t}e^{\pm i\theta}\\
&=A_1 r^{-1}\cos\theta+B_1 r\cos\theta
+C_1 (r\ln r)\cos\theta+D_1 r^3\cos\theta\\
&+E_1 r^{-1}\sin\theta+F_1 r\sin\theta
+G_1 (r\ln r)\sin\theta+H_1 r^3\sin\theta
\end{align*}
$$

$${n}$$の4次方程式と考えた時は、

$$
\begin{align*}
(n_1,n_2,n_3,n_4) &=(1,1,-1,-1),&(\lambda=1)\\
(n_1,n_2,n_3,n_4) &=(0,0,2,-2),&(\lambda=0,2)
\end{align*}
$$

$${\lambda=0}$$のとき$${n^2(n+2)(n-2)=0}$$　なので

$$
\begin{equation*}
\frac{\partial^2}{\partial \theta^2}
\left(\frac{\partial}{\partial \theta}+2\right)
\left(\frac{\partial}{\partial \theta}-2\right)G(\theta)
\end{equation*}
$$

$$
\begin{align*}
\phi &=
A_0' e^{0\cdot t} e^{0\cdot i \theta}
+E_0' e^{0\cdot t}\theta e^{0\cdot i \theta}
+D_2' e^{0\cdot t}e^{-2\cdot i \theta}
+H_2' e^{0\cdot t}e^{2\cdot i \theta} \\
&=A_0'+E_0' \theta
+D_2' e^{-2i\theta}+H_2' e^{2i\theta}\\
&=A_0+E_0 \theta +D_2\cos 2\theta+H_2 \sin 2\theta
\end{align*}
$$

$${\lambda=1}$$のとき$${(n+1)^2(n-1)^2=0}$$　なので

$$
\begin{align*}
\phi &=
A_1' e^{t} e^{i \theta}
+I_1' e^{t}\theta e^{i \theta}
+E_1' e^{t}e^{-i \theta}
+K_1' e^{t}\theta e^{-i \theta} \\
&=A_1' e^te^{i \theta}+I_1' e^t \theta e^{i \theta}
+E_1' e^t e^{-i \theta}+K_1' e^t \theta e^{-i \theta}\\
&=A_1 r\cos \theta+I_1 r\theta\cos \theta
+E_1 r \sin \theta+K_1 r\theta \sin \theta
\end{align*}
$$

$${\lambda=2}$$のとき$${n^2(n+2)(n-2)=0}$$　なので

$$
\begin{align*}
\phi &=
C_0' e^{2t} e^{0\cdot i \theta}
+G_0' e^{2t}\theta e^{0\cdot i \theta}
+A_2' e^{2t}e^{2i \theta}
+E_2' e^{2t}e^{-2i \theta} \\
&=C_0' e^{2t}
+G_0' e^{2t} \theta
+A_2' e^{2t} e^{2i \theta}
+E_2' e^{2t}e^{-2i \theta}\\
&=C_0 r^2 +G_0 r^2\theta
+A_2 r^2 \cos 2\theta+E_2 r^2 \sin 2\theta
\end{align*}
$$

　ようやく求めたかった解の$${r\theta\sin\theta}$$が、$${\lambda=1}$$のときに出てきた。
　教科書によっては、$${\sigma_r=f(r)g(\theta)}$$,$${\sigma_{\theta},\tau_{r\theta}=0}$$の条件から求めているものも見たことはある。例えば以下のサイトなど。

極座標系のエアリーの応力関数

５．一般解のまとめ

　さらに、$${te^{\lambda t}\theta e^{\pm in\theta}}$$も解になるらしいが、こちらについてはちゃんと計算していない。

　参考にした教科書などによれば、一般解は以下のようになる。
（$${\theta t=\theta \ln r}$$,$${\theta te^{2t}=\theta r^2\ln r}$$,$${\theta e^{\pm i\theta}te^t=}$$$${\theta \cos\theta r\ln r}$$,$${\theta \sin\theta r\ln r}$$がそれで、それ以外は、上記の説明の中で求めることができている。）

$$
\begin{align*}
\phi &=(A_0+B_0\ln r+C_0 r^2+D_0 r^2\ln r)\\
&+(E_0+F_0\ln r+G_0 r^2+H_0 r^2\ln r)\theta\\
&+(A_1 r+B_1 r^{-1}+C_1 r^3+D_1 r^2\ln r)\cos\theta\\
&+(E_1 r+F_1 r^{-1}+G_1 r^3+H_1 r^2\ln r)\sin\theta\\
&+(I_1 r+J_1 r\ln r)\theta \cos\theta\\
&+(K_1 r+L_1 r\ln r)\theta \sin\theta\\
&+\sum_{n=2}^{\infty} (A_n r^{n+2}
+B_n r^{-n+2}+C_n r^{n}
+D_n r^{-n})\cos n\theta\\
&+\sum_{n=2}^{\infty} (E_n r^{n+2}
+F_n r^{-n+2}+G_n r^{n}
+H_n r^{-n})\sin n\theta
\end{align*}
$$

　参考文献（上記の二次元の極座標の重調和方程式の一般解が記載されているもの）。

線形弾性論の基礎 | コロナ社

材料力学強度設計への応用 - 株式会社 養賢堂

６．おわりに

　コンクリートの割裂引張試験という単純な試験方法の裏に、2次元極座標の重調和方程式などという数学的なものが潜んでいるというのは驚くべきことだと思う。
　ただ、ひとつ気を付けないといけないのが、この重調和方程式が現実の真の姿ではないということである。重調和方程式は、弾性体をモデル化した、いわば近似理論であり、そもそも弾性体力学自体が分子構造を持った物体をどこまでも滑らかな連続体力学として扱う近似理論である。数学的に解いた解が現実と合っているかどうか実験などで確認することが物理では非常に重要だと思う。
　とはいえ、調和方程式や重調和方程式は数学的に面白い話題であるというのも確かである。エアリーの応力関数に関しては、複素関数で扱う方法もあり、複素数が本質的になる量子力学と違って、あくまで解法上の工夫で複素数を使うにすぎないが、難解な複素数や数学的概念が、物体の力や変形などと関係してくる例としては、今回の話題はちょうどいいのではないかと思う。
　教科書などでは、だいたいこのあたりの話にはあまり触れられないが、確かにこのあたりに触れると、本筋から大きくずれてしまうので触れないのはある意味正解だとは思う。