工学系大学院生のブログ

2人の大学院生による雑記ブログ

第8-2回 一般座標系における方程式[python, julia]




\begin{eqnarray}
\frac{\partial Q}{\partial t}+\frac{\partial E}{\partial x}+\frac{\partial F}{\partial y}+\frac{\partial G}{\partial z}=0
\end{eqnarray}

三次元における一般的な保存形の方程式はこのように書けます.

第8-2回では,一般座標系における上記の方程式を求めていきます.


1 座標の変換


ある座標系((x.y.z))から一般座標系\((\xi,\eta,\zeta)\)に写像を行います.
写像を行うので,次のような関数を考えます.

\begin{eqnarray}
\xi=\xi(x,y,z) \ \ \ (1) \\
\eta=\eta(x,y,z) \ \ \ (2)\\
\zeta=\zeta(x,y,z) \ \ \ (3)
\end{eqnarray}


また,当然ながらこの逆も成り立ちます.

\begin{eqnarray}
x=x(\xi,\eta,\zeta) \ \ \ (4) \\
y=y(\xi,\eta,\zeta) \ \ \ (5) \\
z=z(\xi,\eta,\zeta) \ \ \ (6)
\end{eqnarray}


微小要素の変換として,下記の関係が成り立ちます.

\begin{eqnarray}
\left[
\begin{array}{c}
dx \\
dy \\
dz
\end{array}
\right]
=\left[
\begin{array}{ccc}
\frac{\partial x}{\partial \xi} & \frac{\partial x}{\partial \eta} &\frac{\partial x}{\partial \zeta} \\
\frac{\partial y}{\partial \xi} & \frac{\partial y}{\partial \eta} &\frac{\partial y}{\partial \zeta} \\
\frac{\partial z}{\partial \xi} & \frac{\partial z}{\partial \eta} &\frac{\partial z}{\partial \zeta}
\end{array}
\right]
\left[
\begin{array}{c}
d \xi \\
d \eta \\
d \zeta
\end{array}
\right] \ \ \ (7)
\end{eqnarray}



*** 例 ************

\begin{eqnarray}
x=2\xi+3\eta+5\zeta \ \ \ (8)
\end{eqnarray}

このとき,xの微小要素は\(\xi\)と\(\eta\)と\(\zeta\)の微小要素をそれぞれの傾き分を足したものになるので

\begin{eqnarray}
dx=2d \xi +3d \eta +5d \zeta \ \ \ (9)
\end{eqnarray}

また下記の式が成り立ちます.

\begin{eqnarray}
\frac{\partial x}{\partial \xi}=2 \ \ \ (10) \\
\frac{\partial x}{\partial \eta}=3 \ \ \ (11) \\
\frac{\partial x}{\partial \zeta}=5 \ \ \ (12)
\end{eqnarray}

以上のことから,式(9)~(12)をまとめると

\begin{eqnarray}
dx=\frac{\partial x}{\partial \xi}d \xi + \frac{\partial x}{\partial \eta}d \eta +\frac{\partial x}{\partial \zeta}d \zeta \ \ \ (13)
\end{eqnarray}

よって式(7)の形になります.
**** 例終わり ********




式(7)の逆行列を求め,次のように変換できます.
\begin{eqnarray}
\left[
\begin{array}{c}
d \xi \\
d \eta \\
d \zeta
\end{array}
\right]
=J \left[
\begin{array}{ccc}
y_{\eta}z_{\zeta}-y_{\zeta}z_{\eta} & z_{\eta}x_{\zeta}-z_{\zeta}x_{\eta} & x_{\eta}y_{\zeta}-x_{\zeta}y_{\eta} \\
y_{\zeta}z_{\xi}-y_{\xi}z_{\zeta} & z_{\zeta}x_{\xi}-z_{\xi}x_{\zeta} &x_{\zeta}y_{\xi}-x_{\xi}y_{\zeta} \\
y_{\xi}z_{\eta}-y_{\eta}z_{\xi} & z_{\xi}x_{\eta}-z_{\eta}x_{\xi} &x_{\xi}y_{\eta}-x_{\eta}y_{\xi}
\end{array}
\right]
\left[
\begin{array}{c}
dx \\
dy \\
dz
\end{array}
\right]\ \ \ (14)
\end{eqnarray}


ただし次のように文字を使っています.

\begin{eqnarray}
\frac{1}{J}=x_{\xi}(y_{\eta}z_{\zeta}-y_{\zeta}z_{\eta})+x_{\eta}(y_{\zeta}z_{\xi}-y_{\xi}z_{\zeta})+x_{\zeta}( y_{\xi}z_{\eta}-y_{\eta}z_{\xi} ) \ \ \ (15)
\end{eqnarray}

\begin{eqnarray}
\frac{\partial x}{\partial \xi}=x_{\xi} \ \ \ (16)
\end{eqnarray}


さらに式(7)と同様に,式(1)~(3)から微小区間を考えると次の式が書けます.

\begin{eqnarray}
\left[
\begin{array}{c}
d \xi \\
d \eta \\
d \zeta\end{array}\right]
=\left[
\begin{array}{ccc}
\frac{\partial \xi}{\partial x} & \frac{\partial \eta}{\partial x} & \frac{\partial \zeta}{\partial x} \\
\frac{\partial \xi}{\partial y} & \frac{\partial \eta}{\partial y} & \frac{\partial \zeta}{\partial y} \\
\frac{\partial \xi}{\partial z} & \frac{\partial \eta}{\partial z} & \frac{\partial \zeta}{\partial z}
\end{array}
\right]
\left[
\begin{array}{c}
dx \\
dy \\
dz
\end{array}
\right] \ \ \ (17)
\end{eqnarray}


式(14)と式(17)より行列の要素は等しいので,以下の9つの式が成り立ちます.

\begin{eqnarray}
\xi_{x}=J(y_{\eta}z_{\zeta}-y_{\zeta}z_{\eta}), \ \ \xi_{y}=J(z_{\eta}x_{\zeta}-z_{\zeta}x_{\eta}), \ \ \xi_{z}=J(x_{\eta}y_{\zeta}-x_{\zeta}y_{\eta}) \ \ \ (18)\\
\eta_{x}=J(y_{\zeta}z_{\xi}-y_{\xi}z_{\zeta}), \ \ \eta_{y}=J(z_{\zeta}x_{\xi}-z_{\xi}x_{\zeta}), \ \ \eta_{z}=J(x_{\zeta}y_{\xi}-x_{\xi}y_{\zeta}) \ \ \ (19) \\
\zeta_{x}=J(y_{\xi}z_{\eta}-y_{\eta}z_{\xi}), \ \ \zeta_{y}= J(z_{\xi}x_{\eta}-z_{\eta}x_{\xi}), \ \ \zeta_{z}=J(x_{\xi}y_{\eta}-x_{\eta}y_{\xi}) \ \ \ (20)
\end{eqnarray}


また,ここでは詳しく述べませんが,式(15)で定義される\(J\)はヤコビアンといわれており,一般座標つまり計算空間における体積と物理空間における体積比になっています.

$$ J = \frac{[計算空間のセル体積]}{[物理空間のセル体積]}$$

2 一般座標における方程式


三次元における一般的な保存形の方程式は次のように書けます.

\begin{eqnarray}
\frac{\partial Q}{\partial t}+\frac{\partial E}{\partial x}+\frac{\partial F}{\partial y}+\frac{\partial G}{\partial z}=0 \ \ \ (21)
\end{eqnarray}


チェインルールにより次のように変換できます.

\begin{eqnarray}
\frac{\partial E}{\partial x} &=& \frac{\partial \xi}{\partial x}\frac{\partial E}{\partial \xi}+\frac{\partial \eta}{\partial x}\frac{\partial E}{\partial \eta}+\frac{\partial \zeta}{\partial x}\frac{\partial E}{\partial \zeta} \nonumber \\
&=& \frac{\partial}{\partial \xi} \left( \frac{\partial \xi}{\partial x}E \right)+ \frac{\partial}{\partial \eta} \left( \frac{\partial \eta}{\partial x}E \right)+\frac{\partial}{\partial \zeta} \left( \frac{\partial \zeta}{\partial x}E \right) \nonumber \\
&=& J \left[ \frac{\partial}{\partial \xi} \left( \frac{\xi_x}{J}E \right)+ \frac{\partial}{\partial \eta} \left( \frac{ \eta_x}{J}E \right)+\frac{\partial}{\partial \zeta} \left( \frac{\zeta_x}{J}E \right) \right] \ \ \ (22)
\end{eqnarray}


ここで,\(J\)は計算空間のセル体積と物理空間のセル体積の比であるため定数です.

式(22)の変換を三方向\((x,y,z)\)に対して行い,式(21)に代入,そして\((\xi,\eta,\zeta)\)についてまとめると次のようになります.

\begin{eqnarray}
\frac{\partial}{\partial t} \left( \frac{Q}{J}\right)+\frac{\partial}{\partial \xi} \left[ \frac{1}{J} \left( \xi_x E+ \xi_y F + \xi_z G \right) \right]\nonumber \
+\frac{\partial}{\partial \eta} \left[ \frac{1}{J} \left( \eta_x E+ \eta_y F + \eta_z G \right) \right]+ \\
\frac{\partial}{\partial \zeta} \left[ \frac{1}{J} \left( \zeta_x E+ \zeta_y F + \zeta_z G \right) \right]=0 \ \ \ (23)
\end{eqnarray}

\begin{eqnarray}
\Leftrightarrow\frac{\partial \hat{Q}}{\partial t} + \frac{\partial \hat{E}}{\partial \xi}+\frac{\partial \hat{F}}{\partial \eta}+\frac{\partial \hat{G}}{\partial \zeta}=0 \ \ \ (24)
\end{eqnarray}


それそれの文字は下記の通りです.

\begin{eqnarray}
\hat{Q} &=& \frac{Q}{J} \ \ \ (25)\\
\hat{E} &=& \frac{1}{J} \left( \xi_x E+ \xi_y F + \xi_z G \right) \ \ \ (26)\\
\hat{F} &=& \frac{1}{J} \left( \eta_x E+ \eta_y F + \eta_z G \right) \ \ \ (27)\\
\hat{G} &=& \frac{1}{J} \left( \zeta_x E+ \zeta_y F + \zeta_z G \right) \ \ \ (28)\\
\end{eqnarray}

これで一般座標系における方程式に変換できました.



次回は\(\xi_x\)や\(\eta_x\)等の評価方法について進めていきます.

*コードは第8-5回にあります.

by hide


このエントリーをはてなブックマークに追加

コメントを残す

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です

トラックバックURL: https://teru-hide.com/cfd-generalized-coordinate-system/trackback/
にほんブログ村 英語ブログへ
にほんブログ村