2021年5月4日火曜日

侵徹中に標的が作る抵抗 2 -動的Cavity expansion analysis-

 前回力の釣り合いの式は以下に示すように、右辺が0となっていました。

$$ \frac{d\sigma_{rr}}{dr}+\frac{2}{r}\left(\sigma_{rr}-\sigma_{\theta\theta}\right) = 0$$

上式は力の釣り合いの式なわけですが、標的の各部分についての運動方程式とも言えます。

力が釣り合っていて、右辺が0ということは、侵徹中にも関わらず標的の各部分は静止をしていることを意味しており、

その意味で、前回の結果は静的cavity expansion analysisとも呼ばれます。


では、力が釣り合わず、標的の各部分が運動する状況を考えると、その$\sigma_{rr}$はどうなるでしょうか?というのが今回の話題です。

計算ノートは、githubに置いています。

まず、力の釣り合いの式を、見慣れた運動方程式の形にすると、

$$ \frac{d\sigma_{rr}}{dr}+\frac{2}{r}\left(\sigma_{rr}-\sigma_{\theta\theta}\right)  =- \rho_t \left( \frac{\partial v}{\partial t} +v\frac{\partial v}{\partial r} \right)$$

となります。ここで、右辺が妙な形になっているのは、オイラー座標系での時間微分になっているためです。$v$は各部分の粒子速度であり、$r, u$と同じ向きで外側が正です。また、$v$は$u$の時間微分ですので、同じような手続きを取ると、


$$ v = \left( \frac{\partial u}{\partial t} +v\frac{\partial u}{\partial r} \right) $$

$$ \frac{\partial u}{\partial t} = v \left(1 - \frac{\partial u}{\partial r} \right)$$


として$u$と$v$の関係が得られます。


今、前回の記事から$u$の形は以下の通り

$$u=r\left[1-\left(1-\frac{a^3}{r^3}\right)^{1/3}\right]$$

として与えられているので、$a$を

$$a = a(t)$$

なる関数だと思って微分して

$$ \frac{\partial u}{\partial t} = v \left(1 - \frac{\partial u}{\partial r} \right)$$

に代入すれば、

$$v = \frac{\dot{a} a^{2}{\left(t \right)}}{r^{2}}$$

を得ます。ここで、$\dot{a}$は$a$の時間微分です。

この結果を用いて、前回と同じくまずは降伏領域($a\geq r \geq c$)の$\sigma_{rr}$を求めましょう。

本項第2式に$\sigma_{rr}-\sigma_{\theta\theta} = Y$と$v$を代入して、

$$\frac{d}{d r} \sigma_{rr}{\left(r \right)} = - \frac{2 Y}{r} - \rho_{t} \left(\frac{\ddot{a} a^{2}{\left(t \right)}}{r^{2}} + \frac{2 \dot{a}^{2} a{\left(t \right)}}{r^{2}} - \frac{2 \dot{a}^{2} a^{4}{\left(t \right)}}{r^{5}}\right)$$


という形で得られます。ここで$\ddot{a}$は$a$の時間に対する2階微分です。

前回の記事と同じように積分してやれば、

$$\sigma_{rr}(a) = \sigma_{rr}(c)-2Y\ln{\frac{a}{c}} +\rho_t\left[ -(\ddot{a}a^2+2a\dot{a}^2)\left(\frac{1}{c} - \frac{1}{a}\right) + \frac{a^4\dot{a}^2}{2}\left(\frac{1}{c^4}-\frac{1}{a^4}\right) \right]$$


として$\sigma_{rr}(c)$を得ます。それでは、$\sigma_{rr}(c)$を求めるために、前記事と同じく弾性変形領域($c\leq r$)を解きましょう。

やることは大して変わらなくて、結局、

$$ \frac{d\sigma_{rr}}{dr}+\frac{2}{r}\left(\sigma_{rr}-\sigma_{\theta\theta}\right)  =- \rho_t \left( \frac{\partial v}{\partial t} +v\frac{\partial v}{\partial r} \right)$$

$$\sigma_{rr}-\sigma_{\theta\theta} = \frac{4Ea^3}{3r^4}$$

を代入して、

$$  \sigma_{rr}(a) = \frac{4Ea^3}{9c^3} -2Y\ln{a/c}+ \rho_t\left(\ddot{a}a+\frac{3\dot{a}^2}{2}\right)$$

となります。前回示したように、

$$ \left(\frac{a}{c}\right)^3 =  \frac{3Y}{2E}$$

ですから、

$$ \sigma_{rr}(a) =\frac{2Y}{3}\left[ 1+\ln{\frac{2E}{3Y}}\right] + \rho_t\left(\ddot{a}a+\frac{3\dot{a}^2}{2}\right) $$


となります。また、$\ddot{a}$は極めて小さいこと、$\dot{a}$、すなわちCavityの拡大速度は侵徹速度$V_z$のCavityの$r$方向成分なので$\dot{a}=V_z\cos{\theta}$となることを用いて


$$ \sigma_{rr}(a) =\frac{2Y}{3}\left[ 1+\ln{\frac{2E}{3Y}}\right] + \rho_t\left(\frac{3(V_z\cos{\theta})^2}{2}\right) $$


として、動的なcavity expansion analysisでの$\sigma_{rr}$を得ます。


侵徹中に標的が作る抵抗 1 -静的Cavity expansion analysis-

 侵徹過程における侵徹体の運動は、侵徹体が標的から受ける抵抗によって決まることに疑いの余地はありません。

とはいえ、その抵抗をどうやって求めるか、はなかなか難しい問題です。

本記事と次回、次々回の記事で、標的が作る抵抗を求め、侵徹深さを直接求めてみようと思います。

まず、今考えている状況を図に示します。下の図は前回示した侵徹体の模式図に、侵徹体に加わる応力を加えています。


さて、侵徹体は標的を押し広げるように変形させながら侵徹をするので、侵徹体周辺で標的は以下の図に示すように、侵徹体の周囲に塑性変形領域があり、その外側には弾性変形領域があります。ここで$a$は侵徹体半径で、$c$は塑性変形領域の半径ですがまだわかっていません。
今、この$\sigma_{rr}$を求めてみましょうというのが今回の話題です。
図に示したように侵徹体先端付近の標的は塑性/弾性変形しています。$a \leq r \leq c$では標的は塑性変形しており、$c\gt r$では弾性変形のみが生じています。侵徹体の運動を考えるためには、変形する標的の侵徹体/標的界面で作る侵徹体を押し返す力($F_z$)を求める必要があります。
本項と次回では、単位面積あたりに動径方向に働く力、$\sigma_{rr}$を求めてみます。次々回では、$\sigma_{rr}$から、$z$方向に働く力$F_z$を求めます。

さて、本項では、球対称な極座標系(f(r))での標的の変形を考える必要がありますが、出発点となるのは以下の2式です。

$$\frac{d}{dr}\left[\left(r-u\right)^3\right] = 3r^2$$
$$ \frac{d\sigma_{rr}}{dr}+\frac{2}{r}\left(\sigma_{rr}-\sigma_{\theta\theta}\right) = 0$$

ここで、第1式は標的の各部分の変形前後での質量(体積)保存則、第2式は標的の各部分に成り立つ力の釣り合いの式です。
第1式では、圧力に対して標的の体積は不変であるという、標的を非圧縮性の材料と仮定しています。
また、$r,u,\sigma_{rr},\sigma_{\theta\theta}$はそれぞれ、侵徹体中心からの距離、$r$における変位、動径方向の応力、円周方向の応力です。

また、動径方向のひずみ、円周のひずみをそれぞれ、

$$ \varepsilon_{rr} = \ln{\left(1-\frac{du}{dr}\right)} \approx -\frac{du}{dr}$$
$$ \varepsilon_{\theta\theta} = \ln{\left(1-\frac{u}{r}\right)} \approx -\frac{u}{r}$$

と定義します。
$\sigma_{rr}$と$\sigma_{\theta\theta}$、標的の降伏強度の間にはトレスカの降伏条件から、
$$\sigma_{rr}-\sigma_{\theta\theta} = Y$$
の関係があります。(加工硬化を仮定した形も解くことが出来るが、簡単のため、完全弾塑性体を仮定します。)

さて、侵徹体/標的界面($r=a$)では、標的は$a$まで押しのけられているので$u(r=a)=a$です。
また、$r<a$には標的はないので考える必要はありません。
この初期条件を用いて、第1式を$r:a \rightarrow r (r\lt c)$まで解くと、
$$ (r-u)^3 = r^3-a^3$$
となるから、$u$は$r$について
$$u=r\left[1-\left(1-\frac{a^3}{r^3}\right)^{1/3}\right]$$
となります。(今回は使わないのですが、一応導出しておきます。完全弾塑性でない仮定を置くと使います。)

次に、$\sigma_{rr}$を求めたいわけですが、最初に述べたように、標的は塑性変形領域と弾性変形領域に分かれています。塑性変形領域については、トレスカの条件から第2式中の$\sigma_{rr}-\sigma_{\theta\theta}$は$ Y$であることがわかっているので、直ちに解くことが出来ます。

第2式に$\sigma_{rr}-\sigma_{\theta\theta}=Y$を代入して$r:a\rightarrow c$まで積分すれば、
$$\frac{d\sigma_{rr}}{dr} = \frac 2 r Y$$
$$\sigma_{rr}(a) = \sigma_{rr}(c) -2Y\ln{\frac{a}{c}}$$
を得ます。
$\sigma_{rr}(c)$については、弾性変形領域について考えることで明らかになります。

ということで、$r\gt c$の弾性変形領域について考えてみます。
弾性変形領域ではひずみが小さい(<0.2%以下)なので、常に$r \gg u$です。ということで、第6式中の$u$の2次以上の項を無視すれば、
 
 $$ r^3-3r^2u+3ru^2-u^3 = r^3-a^3$$
 $$ -3r^2u = -a^3$$
 $$ u = \frac{a^3}{3r^2}$$

が得られます。ここで得られた$u$を第3,4式のひずみの定義(第3項)に代入すれば、

$$ \varepsilon_{rr} = -\frac{du}{dr} = \frac{2a^3}{3r^3}$$
$$ \varepsilon_{\theta\theta} = -\frac u r = - \frac{a^3}{3r^3}$$

となります。さて、$\varepsilon_{rr},\varepsilon_{\theta\theta}$間のモールの応力円を考えると、最大せん断ひずみは$\varepsilon_{rr}-\varepsilon_{\theta\theta}$ですから、
$$ \sigma_{rr} - \sigma_{\theta\theta} = 2k = 2G\left(\varepsilon_{rr}-\varepsilon_{\theta\theta}\right)$$
となり、上式に$\varepsilon_{rr},\varepsilon_{\theta\theta}$を代入すると、
$$\sigma_{rr}-\sigma_{\theta\theta} = \frac{4Ga^3}{9r^4} = \frac{4Ea^3}{3r^4}$$
となります。ここで、今、標的は非圧縮性であることを仮定しているのでポアソン比$\nu=0.5$であること、$G=\frac{E}{2(1+\nu)}$であることを用いています。
これで、$\sigma_{rr}-\sigma_{\theta\theta}$を第2式を代入すると、
$$\frac{d\sigma_{rr}}{dr} = -\frac{4Ea^3}{3r^4}$$
となり、これを$r:c\rightarrow \infty$について積分すると、
$$\sigma_{rr}(c)=\frac{4Ea^3}{9c^3}$$
を得ます。
今、明らかに$c$において$\sigma_{rr} - \sigma_{\theta \theta}= Y$であることから、
$$ \sigma_{rr} - \sigma_{\theta \theta} =\frac{3}{2} \sigma_{rr} = Y$$
$$ \sigma_{rr}(c) = \frac{2}{3} Y = \frac{4Ea^3}{9c^3}$$
$$ \left(\frac{a}{c}\right)^3 =  \frac{3Y}{2E}$$
となります。これらを、第9式に代入すると
$$\sigma_{rr}(a) = \frac{4Ea^3}{9c^3}+2Y\ln{\frac{a}{c}} = \frac{2Y}{3}+\frac{2Y}{3}\ln{\frac{3Y}{2E}}$$
$$ \sigma_{rr}(a)= \frac{2Y}{3}\left[1+\ln{\frac{2E}{3Y}}\right]$$
が得られます。
これが、準静的な仮定から得られる動径方向の応力です。
この一連の導出をcavity expansion analysisと呼びます。

この結果と、侵徹深さの間の関係は、また次回(次々回)





2021年5月1日土曜日

侵徹体の形状から侵徹体の重量を求める

 本記事では典型的な侵徹体の質量を、その形状因子のみから計算することを目的とします。質量は体積が明らかであれば密度を掛ければ算出できるので、ここでは特に体積の導出方法について検討します。(逆に形状と質量があれば平均的な密度を求めることが出来ます。


今回はsympyで展開したので、興味がある方はこちらをご覧ください。


侵徹体の形状は様々なものがありますが、ここでは以下のように、半径$a$、長さ$L$の円柱部に半径$s$の円を侵徹体軸上で回転させた先端部を貼り付けた侵徹体を考えましょう。この侵徹体形状は低速度侵徹を解析的に解くことが可能な形状であり、侵徹の試験において多く用いられています。




このような侵徹体の先端形状を定義する際に、CRH(Calibre-Head radius)と呼ばれる以下の指標を用います。

$$ \psi = \frac{s}{2a}$$

この指標を用いることで、侵徹を簡潔に表すことができます。(その辺の詳細はまた後日

本項の目標は、上述の特徴を有する侵徹体の体積$V$を

$$ V = \pi a^2 (L+ka)$$

のように、円柱部の断面積を用いて簡潔に表現することです。この$k$を求めてみましょう。


今、上の図には、$\theta_0, l$が定義されていますが、その具体的な形はわかっていません。まず、図から明らかに、

$$\sin{\theta_0} = \frac{s-a}{s} = \frac{2\psi-1}{2\psi}$$

なので、$\theta_0$は$\sin^{-1}{\frac{2\psi-1}{2\psi}}$であることがわかります。また、$l$は明らかに$s\cos{\theta_0}$ですが、

$$\cos{\theta_0} = \sqrt{1-\sin^2{\theta_0}}=\sqrt{\frac{4\psi^2-(4\psi^2-4\psi+1)}{4\psi^2}}=\frac{\sqrt{4\psi-1}}{2\psi}$$

であること、また、$s=2\psi a$であることから、

$$l=a\sqrt{4\psi-1}$$

であることがわかります。


では、侵徹体先端部の重さ(=体積)を求めましょう。

侵徹体軸上に、$\theta_0 \rightarrow \pi/2$の範囲で存在する半径$s$の円の体積を求めるためには、

$$y = s \sin{(\theta)} - (s-a) $$

$$ x = s \cos{(\theta)}$$

$$ V_o = \int_{0}^{l} \pi y^2 dx$$

を解けばいいわけですが、$x$は$\theta$の関数なので、$\theta$についての導関数を求めて置換すると、

$$ dx=-s\sin(\theta) d\theta$$

$$l:\ 0 \rightarrow l$$

$$\theta:\ \pi/2 \rightarrow \theta_0$$

から、

$$ V_o = \int_{\pi/2}^{\theta_0} -\pi s\sin(\theta)  y^2 d\theta$$

となり、これを積分して気合いで整理すると、

$$ V_o = \pi a^{3} \left(- 8 \psi^{3} \operatorname{asin}{\left(\frac{\sqrt{4 \psi - 1}}{2 \psi} \right)} + 4 \psi^{2} \operatorname{asin}{\left(\frac{\sqrt{4 \psi - 1}}{2 \psi} \right)} + \sqrt{4 \psi - 1} \left(- 4 \psi^{2} + \frac{4 \psi}{3} - \frac{1}{3}\right)\right) $$


となります。ここまでくれば、

$$ V = \pi a^2 (L+ka)$$

と式の形を合わせるように、$k$を

$$k = \left(- 4 \psi^{2} \left(2 \psi - 1\right) \operatorname{asin}{\left(\frac{\sqrt{4 \psi - 1}}{2 \psi} \right)} + \sqrt{4 \psi - 1} \left(- 4 \psi^{2} + \frac{4 \psi}{3} - \frac{1}{3}\right)\right)$$

と定めれば、当初欲しかった$k$が得られます。

後は$V$に侵徹体の材料の密度$\rho$を掛ければ、侵徹体の重量がわかりますし、逆に重量がわかっていれば、侵徹体の平均的な密度がわかります。

2021年1月9日土曜日

高速度侵徹の侵徹効率

前回の記事では、高速度侵徹におけるエネルギー効率を考えました。この例では装甲の抵抗が衝突速度と同じオーダーで増加するので、その極限はエネルギー効率は密度から定まる非0な値となりました。

よく知られているように、高速度侵徹の高速度側では侵徹深さは侵徹体の速度に依存しなくなります。ということで、運動エネルギーに対する侵徹深さの効率は高速度側で0になります。ここで、侵徹効率$\eta_p$をどのように考えるかは議論の余地がありますが、

\[ \eta_p(V_0) = P(V_0)\times R_t(V_0)/KE(V_0) \]

として定義してみます。ここで、$V_0$は衝突速度であり、$R_t$は標的強度項です。つまり、投入された運動エネルギーのうち、どれだけ標的強度項が消耗したエネルギーを占めるかという定義になります。

以下に、Anderson-Walkerモデルを用いて計算した$\eta$、$\eta_p$と$V_0$の関係を示します。ここで、標的は鉄($\rho$:7800g/cm3, Y:1GPa)、侵徹体はWHA($\rho$:17800, Y:1GPa)です。


結論としては見てのとおりですね。

2021年1月2日土曜日

高速度侵徹のエネルギー効率

侵徹体の運動エネルギーが標的の侵徹にすべて利用される低速度侵徹と異なり、高速度侵徹では侵徹体が消耗するために初期の運動エネルギー$KE0$は標的が受け止めるエネルギー$W$とは一致しません。
そこで、高速度側の極限では初期の運動エネルギーの何割くらいが侵徹に利用されているのかふと気になりました。

高速度侵徹の極限では侵徹深さ$P$は、侵徹体長さ$L$、侵徹体密度$\rho_p$、標的密度$\rho_t$を用いて

$$P = \sqrt{\frac{\rho_p}{\rho_t}}L $$

と表されます。一方、衝突時の侵徹体の運動エネルギー$KE$は、衝突速度$V_0$を用いて

$$KE = \frac 1 2 L\rho_p V_0^2$$

と表されます。標的の侵徹に使用されるエネルギーは、少し考えどころです。エネルギーは力×距離ですから、侵徹深さがわかってる現在、標的が作る力がわかれば直ちに求められます。標的が作る力はベルヌーイの式を参考にすれば、標的の侵徹に使用されるエネルギー$W$は

$$ W= \frac{1}{2}\rho_t u^2\times P = \frac{1}{2}\rho_t u^2\times \sqrt{\frac{\rho_p}{\rho_t}}L $$

として簡単に得られます。$W$と$KE$の比がエネルギー効率ですから、比を取って適当に解くと、エネルギー効率$\eta$は

$$ \eta = \frac{W}{KE} = \frac{\rho_{t} \sqrt{\frac{\rho_{p}}{\rho_{t}}}}{\rho_{p} \left(\sqrt{\frac{\rho_{t}}{\rho_{p}}} + 1\right)^{2}} $$

として簡単に得られます。まあ適当に、侵徹体としてタングステンの密度(19.3g/cm3)、
標的の密度として鉄の密度(7.8g/cm3)を代入してみれば、$\eta=0.24$程度となり、
高速度侵徹のエネルギー効率は24%程度になることがわかります。

以上、ふと思って計算してみましたが、たぶん昔の偉い人が似たような計算をしている気がします。

2020年2月1日土曜日

AWモデルにおけるDynamic cavitation analysisの解析解

完全にメモ書きです。

ちょっと頑張ってL/D効果も含めた侵徹深さの評価をしたい話

\[ \left(1+\frac{\rho_tu^2}{Y_t}\right)\sqrt{K_t-\rho_tu^2\alpha^2}= \left(1+\frac{\rho_tu^2\alpha^2}{2G_t}\right)\sqrt{K_t-\rho_tu^2} \]

に対してニュートンラフソン法、あるいは二分法を用いれば$\alpha$を求められると書いていました。これは、この式が平方根を含むために直接解析解を求めることができないためで、AWモデルを用いた計算の大部分を$\alpha$を求めるために消費していました。

一方で、$\alpha$の値域は$\alpha \gt 1$の実数というのは明らかですので、

\[ \left(1 + \frac{\rho u^{2}}{Y}\right)^{2} \left(K_{t} - \alpha^{2} \rho u^{2}\right) = \left(1 + \frac{\alpha^{2} \rho u^{2}}{2 G}\right)^{2} \left(K_{t} - \rho u^{2}\right) \]

とした解を求めても基本的には等価です。ここで$\alpha^2$については明らかに解析解が存在していて、sympyで適当に求めると、

…ここに示しておきたかったんですが長すぎて数式に表示されませんでした。

import sympy 
rho,u2, Y, K_t, alpha, G = sympy.symbols("rho, u^2, Y, K_t, alpha^2, G")
lh = (1 + rho*u2 /Y) * (1 + rho*u2 / Y)* (K_t - rho*u2 * alpha)
rh =  (1 + rho*u2 * alpha * 0.5 / G) * (1 + rho*u2 * alpha /2 / G) * (K_t - rho*u2)
eq = sympy.Eq(lh,rh)
s1,s2 = sympy.solve(eq,alpha)
s2

と書くことが出来ます、とだけ示しておきます。これの平方根を取れば、$\alpha$になります。数値解との対応は以下のとおりです。

個人的には計算時間が4msから1msまで減らせたので満足です。

2019年9月1日日曜日

ノイマン効果ってなんだ?な話

ノイマン効果ってなんだ


いつの間にか、ブログ記事の多くが終末弾道の話になっていて、自分の出自と無関係な方向にどんどんと進んでいるな、という感じがあります。
まぁ、それはともかくとして、終末弾道を考えるとき簡単の場合には侵徹体はバルクと考えてよく、侵徹体は砲から発射されたままの形状を初期条件として侵徹を考えることができます。このような例として、低速度のAPからAPFSDSまでが該当します。
一方で、侵徹体をバルクとして考えられない例として、侵徹体をその場で形成するHEATが考えられます。

HEATの侵徹体の形成について説明を行う際には、モンロー効果、ノイマン効果といった効果が引き合いに出されるかと思います。
こういう効果の説明を簡潔に行っているのは、やはりWikipediaですので、ノイマン効果に関する記述を引っ張ってきますと、

ノイマン効果(Neumann effect)とは、ドイツの科学者、エゴン・ノイマン(Egon Neumann)が1910年に発見した、モンローの円錐形のくぼみに金属板で内張り(くぼみと同じ形の金属の円錐をはめ込むこと)をすると穿孔力がさらに強くなる現象。
爆轟が進行して金属の内張りに達すると、爆轟波によりライナーは動的超高圧に晒されユゴニオ弾性限界を超える圧力に達すると固体の金属でも可塑流動性を持つようになり、液体に近似した挙動を示すようになる。これにより、融着体と呼ばれる金属塊となって前方へ超音速で飛び出していく。
モンロー/ノイマン効果-Wikipedia

のように説明されるようです。結局の所、ノイマン効果というのは「円錐形のくぼみに金属板で内張り(くぼみと同じ形の金属の円錐をはめ込むこと)をすると穿孔力がさらに強くなる現象。」として説明される現象と考えていいようです。
さて、では実際にどのような現象が起こっているかという話ですがWikipediaの記事はなかなか言っていることが難しくてどういうことが起きているのかなかなかわかりません。
そこで、今回の記事では1948年にBirkhoff, MacDougall, Pugh, Taylorにより提案された最も単純なモデルに則ってノイマン効果がどのようにして生じているかを概観したいと思います。
今回の記事はBallisticsのShaped chargeの項を元にしたものになります。

ライナーの中で成り立つ幾何学的関係

今回の記事では以下の図のような状況を考えます。
 Fig. 1 BMPTモデルで考える前提

爆轟速度$D$[m/s]の炸薬中に角度$2\alpha$[rad]の開き角で挿入された円錐状のライナーが挿入されており、爆轟は炸薬中を平面波的に右へ進んでいます。
爆轟が通過し終わったライナーは円錐中央部に叩き出されます。爆轟はライナー形状に関係なくまっすぐ右に進んでいますので、ライナーもまた円錐中央部から爆轟面まで一定の角度をなして分布していると考えられます。そこで、爆轟によって叩き出されたライナーの角度を$2\beta$[rad]とします。

Fig. 2 BMPTモデルで成り立つ各部分の速度と角度の関係

BMPTモデルでは、Fig. 1に示した初期ライナー角と、叩き出されたライナー角の関係から、Fig. 2のような速度の関係を幾何学的に示しました。

ここで、点Aは叩き出されたライナーが円錐中央部で集まっている点(ジャンクション)であり、$V_1$[m/s]はジャンクションが前方に進む速度です。また、点Pは今まさに爆轟が到達した点であり、点Bに位置するライナーは速度$V_0$[m/s]で点Bに向かって叩き出されます。この時、点Pに立った読者は点Aに向かって$V_2$[m/s]で接近します。
また、$D$[m/s]は炸薬の爆轟速度です。

$V_0$は炸薬によって叩き出されたライナーの速度であり、この図形中で唯一不明なパラメータです。$V_0$はGurney methodという炸薬のエネルギーと、ライナー質量/成形炸薬質量比の関係と成形炸薬の特性で決まるとされていますが、今回の記事ではこれは魔法の力によって所与のものとしましょう。

さて、そうすると未だ明らかでないパラメータは$\theta$[rad]と$V_1$、$V_2$、$\beta$[rad]となりますが、$\theta$についてはFig. 2から明らかに

$$\theta = \frac{\pi}{2} -\frac{\beta -\alpha}{2} $$

となります。そうすると、$V_1$は正弦定理から

$$\frac{V_0}{\sin \beta} = \frac{V_1}{\sin\left(\frac{\pi}{2} -\frac{\beta -\alpha}{2}\right) } = \frac{V_1}{\cos \left(\frac{\beta -\alpha}{2}\right)} $$

が直ちに得られます。
さらに$V_2$はオレンジ色の補助線に従うと、

$$V_2 = V_1\cos \beta + V_0 \cos \left(\frac{\pi}{2} -\frac{\beta -\alpha}{2}\right) =V_1\cos \beta + V_0 \sin \left(\frac{\beta -\alpha}{2}\right)  $$

が同様に得られます。
最後に、$\beta$についてですが、これは右側の三角形について考えると、

$$\frac{V_0}{{\sin \left(\pi-2\theta\right)}}=\frac{\frac{D}{\cos \alpha}}{\sin \theta}$$
$$\frac{V_0\cos\left(\left(\beta-\alpha\right)/2\right)}{\sin\left(\beta-\alpha\right)} = \frac{D}{\cos \alpha}$$

という$\beta$のみが未知の式が得られるので、これを数値的に解くことで得られます。
以上でFig. 2の項目はすべて明らかになりました。

ジャンクションでの速度関係

次に、ジャンクションで成り立つ関係に移りましょう。以下にジャンクションで成り立つ速度の関係について示します。

Fig. 3 ジャンクションで成り立つ速度の関係

点Pから速度$V_2$でライナーが点A(ジャンクション)に流れこむ時、ある部分は真後ろ(AC方向)に、ある部分は真正面(AB方向)に打ち出されるでしょう。この時、理想的な円錐であればその他の速度成分は持ちません。これは円錐系のライナーが軸対称に流れ込むために、軸方向以外の速度成分がすべてキャンセルされるためです。この、ジャンクションから放出されるライナーは軸方向以外の速度成分を持たず、真正面に打ち出されるという現象こそがノイマン効果をもたらす現象と言えるでしょう。

さて、一般に用いられるライナー材料は柔らかく、材料の強度項は無視していいでしょう。一方で、ジャンクションに流れ込む速度$V_2$はそれほど高くないので*、ここでは非圧縮性の流体として扱ってみましょう。

そうすると、ジャンクション(点A)に入る圧力と出る圧力はベルヌーイの式から釣り合っており、AC方向にでる速度$V_{AC}$と、AB方向に出る速度$V_{AB}$と$V_2$の間に成り立つ関係は、

$$\frac 1 2 \rho V_2^2 = \frac 1 2 \rho V_{AC}^2 = \frac 1 2 \rho V_{AB}^2$$

となり、Fig. 3の関係が得られます。

メタルジェット速度、スラグ速度

さて、Fig. 3の関係から、ジャンクションから出ていく物質の速度はAC方向、AB方向、いずれの方向にも$V_2$であることがわかりました。
一方で、ジャンクションは$V_1$で動いているので、AB方向に出る物質と、AC方向に出る物質を外から見ると、

$$V_{AB} = V_1 + V_2$$
$$V_{AC}= V_1 - V_2 $$

だけの速度で動いていることがわかります。後ろ向きの方向(AC方向)の速度をスラグ速度$V_s$と呼び、前向きの方向(AB方向)の速度をメタルジェット速度$V_j$と呼んでみると

$$V_{AB} =V_j = V_1 + V_2$$
$$V_{AC}=V_s = V_1 - V_2 $$

と書けます。
以上がBMPTモデルに基づいたメタルジェット速度とスラグ速度の導出です。

BMPTモデルは円錐中央部に流れ込む速度関係を考えるという単純な仮定で見通しよく成形炸薬で生じるスラグとメタルジェットの速度を求められるのでいいですね。

メタルジェット速度と侵徹深さ

口径6 in.のTNTにライナー径5 in.の鉄ライナーを$\alpha$=45度の角度で埋め込んだときのメタルジェット速度を以下に示します。計算に必要なデータはBallisticsより抜粋しました。

Fig. 4 口径6 in.のTNTにライナー径 5 in.の鉄ライナー($\alpha$=45度)のメタルジェット速度とスラグ速度

メタルジェットの速度は期待したとおり高速であることがわかります。このような速度域の侵徹深さを計算してみると、ベルヌーイの式で期待される、侵徹深さが密度のみで決定される侵徹様式に到達していることが以下の図からわかります**。



Fig. 5 軟鋼に鉄侵徹体を衝突させたときの侵徹深さ(修正ベルヌーイの式により計算)

結局ノイマン効果ってなんだ

と言われれば、現象だけを見ればやはり、
炸薬に叩き出された円錐型のライナーが円錐中心部で真正面に高速で打ち出されることで高速なメタルジェットを形成し、そのメタルジェットが深い侵徹孔を形成するという現象がノイマン効果と言えるでしょう。
そしてなぜそのような現象が生じるか、というのは、まさに本記事で議論した内容が最も基礎にあると考えていいと思います***。

本記事の内容は以上です。

*炸薬をTNT、ライナー材料を鉄とした場合、$V_2$はたかだか2km/s程度です。
**ここでは形成されるメタルジェット長さや実際の侵徹深さの議論は一度置いておいて、メタルジェットくらいの速度域での侵徹様式にのみ注目しています。
***本モデルや個々から発展したモデルはスラグとメタルジェットの質量比やメタルジェットの長さについても議論できるんですが筆者がそこまで興味がないのでこのへんで。