記事リバイナル②(台風の進路(コリオリの力))

※再掲載第２弾！！2006年8月16日の記事です。

そろそろ台風が頻発する季節になりました。

　今日は,北半球では赤道付近で発生し,海域から多量の水分やエネルギーを吸収しながら発達して北上する台風が,なぜ右（東）の方に進路を変えていくのか？

　そして,なぜ上空から見て左巻き（反時計回り）の風が吹くのか？

　という,ごくありふれた疑問について解説してみます。

　例として,ちょっと古いけど左回転しているLPレコードがあり,その上に"一寸法師"よりも小さい小人が乗っているという仮想的な状況を考えてみます。

(↑左回転は仮定であり,実際のLPレコードは,裏から見ない限り右回転(時計回り)です。)

　レコードの中心は地球の北極に相当し,レコードの最大半径の部分は地球の赤道に相当します。

　

　まず,レコードの回転している"最大半径＝赤道"の上にいる小人が"レコードの中心＝北極"めがけて真っ直ぐに小石を投げたとします。

　本人は真っ直ぐ中心に向かって投げたつもりでも,小石が手を離れた瞬間には慣性によって小石はレコードの回転スピードと同じ速度で右に向かう接線速度を持ちますから,実はそれは中心の方向に向かって真っ直ぐに飛ぶのではなく　て,次第に右の方に逸れていくということになりますね。

　↑ここで右というのは,小石を投げた小人にとっての右です。(わかってるとは思いますが念のため)

　次に逆に"中心＝北極"の上に小人がいて今度は"最大半径＝赤道"めがけてやはり小石を投げたとします。

　今度は北極で小人は自転しているかもしれませんが,スピンの回転半径はゼロなので,その慣性による小石の左右方向への速度はゼロですから確かに"真っ直ぐ"進むはずです。

　ところが,レコードの上,つまり北半球の地球上にいる人は"左から右＝西から東"に回転しています。その人から見ると"上＝北"から真っ直ぐ飛んでくる小石は"下＝南"から見て"左に左に（西に西に）"逸れていくように見えます。

　逆に"小石を投げた方＝北"から見ると,見かけ上はやはり右の方に逸れていくわけですね。

　小石を台風だとみなし地球の自転の角速度をΩとすると　Ωは"360度＝2πラジアン(rad)"を24時間で回転する角速度です。

　地球の半径をRとし,緯度をθとすると,そこでの回転半径はRcosθですから,回転の接線速度はΩRcosθです。

　したがって"赤道"での接線速度は最大速度"ΩR＝時速1667ｋｍ＝秒速463ｍ"ですが,日本付近の緯度:θ＝35度での接線速度は"ΩRcosθ＝秒速379ｍ"です。

　日本付近では回転速度が赤道より"約２割＝秒速80ｍ"程度減っています。

　したがって,赤道付近で発生した台風は地球のまさつにより"ΩR＝時速1667ｋｍ＝秒速463ｍ"の地球に引きずられて慣性による右向きの速度を持っていて,その右向き速度は北上しても全く変わらないものです。

　しかし,地球自身の回転速度は緯度が上がるにつれて次第に小さくなるものですから,日本付近では１秒間に80ｍくらいの割合で,右（東）へ右へと逸れていくことになるわけです。

　先に,ＬＰレコードの例で述べたように仮に北極で台風が発生して南下したとしてもそれは右に逸れていきます。

　実は北半球ではどこから投げた石も見かけ上,右に逸れていくわけです。

　例えば,スナイパー:「ゴルゴ13」が１ｋｍ遠方の標的を狙って狙撃しても,弾丸は僅かに右に逸れていくのでそれを勘定に入れて狙う必要があるわけです。

　もしも南半球なら逆に左に逸れるわけですね。

　こうした北半球で右にそれる現象は結局,遠心力などと同じく"見かけの力＝慣性力"が働いていると考えることができて,それを発見者の名前にちなんでコリオリ(Coｒiolis)の力といいます。

　次に,北半球での台風を考えると,台風ですから"中心＝目の部分"の気圧が最低でまわりの気圧は目の部分のそれより高いわけです。

　風はどのように吹くか,というと水が高いところから低いところへと流れるように,風も気圧の高いところから低いところ目指してその気圧のスロープに沿って吹いていきます。

　もしもコリオリの力がなければ,風は"外周部から中心に向かって一直線に進む＝落下していく"はずなのですが,コリオリの力によって気圧のスロープも右にねじれてしまっています。

　それ故,風は外周部から中心に向かっていくときに,右にフックして逸れていきながら,最後は中心の気圧最低の目に向かっていくことになり,そのために左巻き（反時計回り）になるのです。

　南半球での台風は逆に右巻き（時計回り）ですね。

　どこかの"バカな大学教授"が風呂の水が排水口へと流れていき排水されるときに,北半球では左巻き（反時計回り）だが,赤道を越えて南半球に入ったとたんに右巻き（時計回り）に変わる,"と主張したと聞きましたが,それは誤りですね。

　風呂の水程度の規模では地球自転の影響などは出てきません。

　たまたま排水口付近で左巻きの角運動量を持っていたら左巻きになり,逆なら右巻きになるというだけで,それはカオス現象,偶然の産物にすぎません。

　しかし台風くらいの規模になると地球の自転がもろに効いてきます。

　遠心力の加速度は緯度θでΩ²Ｒcosθですから,最大でも重力の加速度の0.3％程度です。

　北極で体重100kg重の人が赤道で体重を測っても300ｇ重くらいしか軽くはなりません。

　一方,コリオリの力の加速度は台風の北上速度をｖとして2Ωsinθ×ｖです。

　Ω＝7,2×10^-5/ｓですから,緯度θが35度で台風の北上速度が100mを10秒で走る程度の時速36km程度なら,加速度ａ＝8.3×10^-4m/ｓ²であり,重力ｇ＝9.8m/ｓ²の約0.01％程度です。

　そこで,コリオリの力の加速度は,最大で重力加速度の0.3％程度しかない遠心力のさらに1/30程度にすぎません。

　しかし,台風程度の規模だとそれがかなり効いてきます。

　地球自転の証拠であるとされるフーコー(Foucauｌt)の振り子をこのコリオリ力で説明することもできます。

　ニュートン(Isac Newton)は"慣性系の同等性＝ガリレイの相対性原理"は認めても"回転系を含む非慣性系の同等性＝マッハ(Mach)原理 → 一般相対性原理"を認めることをあきらめました。

　そして,彼が"絶対座標系＝絶対空間"に固執せざるを得なかったのも,こうした"遠心力やコリオリ力の絶対性"を解消する道はない,という考えからだったという話もあります。

　こうした"見かけの力＝慣性力"の扱いはとても悩ましいところがあります。 

※追記　(２０１８年10/26)

今年の夏は北東に進む常識的なモノとは異なる台風がいくつかありました。しかし少し考えて見るに,これは,別に上記記事と矛盾するものではなく,非常にゆっくりと北上していれば,右に右に曲がるのですからやがては1周,あるいはそれ以上回転して西方に向かうこともあり得るし,そこから急に移動速度が大きくなったなら,そのまま西進する行路をとっても不思議ではないと思います。

地震に関する過去の科学記事(バックナンバー)

　この際ですから,地震関係の過去の科学記事を紹介しておきます。

(↑ははーん？なるほど手抜きだね。)

　まず,2006年11/４の「結晶内での弾性波(地震波)」です。

　そして,2009年９/８の「定量的地震学１」,９/15の「定量的地震学２」,９/23の「定量的地震学３」,９/29の「定量的地震学４」,11/17の「定量的地震学５」,11/19の「定量的地震学６」です。

　この後,このシリーズをなぜ途中で中断したのかという理由はよく覚えていませんが,４から５の間があいてるのは10月から介護ヘルパー２級をとるために学校に通っていた時期と丁度一致しています。

　さらには,,2009年12/10の「震源の探知(大森公式等)」があります。

　最後の小記事については丸々転載しておきます。

　(転載開始)↓

　2009年12月10日　「震源の探知(大森公式等)」

今日は,ちょっと手抜きのツナギで2006年11/14の記事 「結晶内での弾性波（地震波）」の続きとして,中学か高校の地学で習った簡単な知識を披瀝してお茶をにごします。

　年末,体調はよくも悪くもないですがヒマ人なりに予定がこみ合っていて貧乏なこともありフリ－マンをエンジョイというわけにはいきません。

　ではまず,記事の再掲から始めます。

　(※再掲)

　今日は等方的とは限らない一般の結晶内での弾性波,特に地震波について考察してみます。

弾性体の密度をρ,応力テンソルを{σ_jk},それを構成する部分の局所的速度をｕ＝{ｕ_j}とすると,この弾性体が従うべき運動方程式は流体方程式と同じくρｄ²ｕ_j/ｄｔ²＝∂σ_jk/∂ｘ_kで表与えられます。

そして,歪み速度テンソル{ｕ_jk}をｕ_jk≡(∂ｕ_j/∂ｘ_k＋∂ｕ_k/∂ｘ_j)/2で定義すれば,これは{σ_jk}と同じく反対称の単なる回転を除いた対称テンソルです。これは６個の独立成分を持ちます。

線形弾性体近似では,フック(Hooke)の法則:σ_jk＝Ｃ_jklmｕ_lmが成立するため,先の運動方程式はρｄ²ｕ_j/ｄｔ²＝Ｃ_jklm(∂²ｕ_m/∂ｘ_k∂ｘ_l)となります。

Ｃ_jklmは,一般に81個の成分を持ちますが,{σ_jk}も{ｕ_jk}も対称テンソルなのでｊとｋ,ｌとｍについて対称ですから,独立成分は６×６＝36個となります。

また,歪みエネルギーＷはΔＷ＝σ_jkΔｕ_jkから求まり,対称２次形式Ｗ＝(1/2)Ｃ_jklmｕ_jkｕ_lmになるので,Ｃ_jklmはｊ,ｋとｌ,ｍの交換についても対称であり,結局,その独立成分は(30/2)＋6＝21個となります。

この方程式の平面波の解をｕ_j＝ｕ_0jexp[i(ｋｒ－ωｔ)]として代入すると,ρω²ｕ_j＝Ｃ_jklmｋ_kｋ_lｕ_mとなります。これを書き直すと(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l)ｕ_m＝0 となります。

この１次方程式が自明でない解を持つためには,３行３列の行列係数:(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l) (ｊ,ｍ＝1,2,3)の行列式がゼロ,すなわち,det(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l)＝0 が成立する必要があります。

そして,ω²を未知数としてこの方程式を解けば,ω²の３個の解が得られ,それらはｋの関数となります。

ところで,もしも結晶が等方性弾性体であれば,歪みエネルギーＷ＝(1/2)Ｃ_jklmｕ_jkｕ_lmは座標系の回転に対して不変なスカラーであるはずですが,ｕ_jkの２次形式の形で得られる独立な不変スカラーはｕ_kk²とｕ_jkｕ_jkのみです。

そこで,適当な係数λ,μを選んでＷ＝(1/2)λｕ_kk²＋μｕ_jkｕ_jkと書けるはずです。そこで,応力テンソル{σ_jk}＝{Ｃ_jklmｕ_lm}はσ_jk＝λｕ_llδjk＋2μｕ_jkと書けます。ここに弾性定数λ,μはラメ(Lame)の定数と呼ばれています。

このときには運動方程式もρｄ²ｕ/ｄｔ²＝(λ＋μ)∇(∇ｕ)＋μ∇²ｕと,簡単になります。

先に挙げた１次方程式の係数の行列要素はρω²δ_jm－Ｃ_jklmｋ_kｋ_l＝(ρω²－μｋ²)δ_jm－(λ＋μ)ｋ_jｋ_mとなります。

ｋの向きをｘ軸の正の向きに取ると,ｋ₁＝ｋ＝|ｋ|,ｋ₂＝ｋ₃＝0 ですから,行列[(ρω²－μｋ²)δ_jm＋(λ＋μ)ｋ_jｋ_m]は対角成分が[ρω²－(λ＋2μ)ｋ²]と２つの(ρω²－μｋ²)の対角行列になります。

したがって,det(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l)＝0 の解はＶ_P²≡(ω_P/ｋ)²＝(λ＋2μ)/ρ,Ｖ_S²≡(ω_S/ｋ)²＝μ/ρとなります。

それぞれの速度に対応する平面波は,速度Ｖ_Pで波の伝播方向への振動である"縦波＝Ｐ波"と,速度Ｖ_Sで波の伝播方向に垂直な方向の振動である"横波＝Ｓ波"です。

しかし,もしも異方性の弾性体の場合ならω＞0 の３つの解ω(ｋ)はｋの関数として１次の同次式ではあっても,単なる１次関数になるとは限らず,弾性波は一般に分散性の波でもあると思われます。

そこで,伝播速度は"単一波の速度＝位相速度"ではなく,群速度Ｖ＝∂ω/∂ｋで与えられると考えられます。

ところで,σ_jkやｕ_jkを[σ]＝^t(σ₁₁,σ₂₂,σ₃₃,σ₂₃,σ₃₁,σ₁₂)のように６次元の列ベクトルで表わすボイト(Vogit)の表記という表記法があります。

この表記で表現すると,フックの法則は[σ]＝Ｃ[ｕ]と簡単になります。ここでＣは６行６列のテンソル行列です。

また,このとき,特に結晶が等方性弾性体ならＣの成分のうちで独立なのはやはり２つだけで,先に述べたラメの定数はλ＝Ｃ₁₂,μ＝Ｃ₄₄＝(Ｃ₁₁－Ｃ₁₂)/2と表わされます。

例えば結晶が立方体構造をしている立方晶系では独立成分は３つです。つまり,Ｃ₁₁＝Ｃ₂₂＝Ｃ₃₃,Ｃ₁₂＝Ｃ₂₃＝Ｃ₃₁,かつ４～６行の成分は４～６列しか成分のない対角行列であってＣ₄₄＝Ｃ₅₅＝Ｃ₆₆であり,結局,この結晶構造を示すにはＣ₁₁とＣ₁₂とＣ₄₄の３つだけあれば十分である,ということになります。

この条件ではｘ軸をｋの向きに取り,先のdet(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l)＝ 0 でのＣ_jklmをボイトの表記でＣを表わせば,ξ＝(ω/ｋ)²,ａ＝Ｃ₁₁/ρ,ｂ＝Ｃ₄₄/ρ,c＝Ｃ₁₂/ρとして(ξ－ａ)(ξ－ｂ)²＝0 ですから,解はξ＝ａ,ｂとなります。

つまり速度をＶ_P,Ｖ_Sとすると,Ｖ_P²＝Ｃ₁₁/ρとＶ_S²＝Ｃ₄₄/ρの２つの速度が得られます。これは等方性弾性体でのＣ₄₄＝μ,Ｃ₁₁＝λ＋2μのケースと一致しています。

一方,六方晶系ではＣ₁₁＝Ｃ₂₂で,Ｃ₁₃＝Ｃ₂₃,また４～６は対角行列でＣ₄₄＝Ｃ₅₅,Ｃ₆₆＝(Ｃ₁₁－Ｃ₁₂)/2ですから,結局のところＣ₁₁,Ｃ₁₂,Ｃ₁₃,Ｃ₃₃,Ｃ₄₄の５つだけがあれば十分となります。

結果だけ書くと,Ｖ_P²＝Ｃ₁₁/ρ,およびＶ_S²＝Ｃ₄₄/ρ,(Ｃ₁₁－Ｃ₁₂)/(2ρ)となり,"横波＝Ｓ波"には２つの速度があるので地震の観測ではＳ波のほうに２重の波が観測されると予測されます。

ただし,異方性の結晶では一般に波は完全にＳ波とＰ波になるように行列が対角形になるとは限らないので,これらの計算においては偶々波の進行方向が結晶の対称軸と一致したり直交していたりする特別なケースだけしか扱っていません。

一般的な扱いについては,暇があってその気になればまた記事にするかもしれません。

参考文献；ランダウ＝リフシッツ 著「弾性理論」(東京図書),角谷典彦 著「連続体力学」(共立出版)

(再掲終了※)

　と書きましたが,ここからもっと身近な地震の話題へとつなげます。

等方性媒質を仮定すると,Ｖ_P＝{(λ＋2μ)/ρ}^1/2,Ｖ_S＝(μ/ρ)^1/2ですからＶ_P＞Ｖ_Sです。

そこで,震源Ｏから地震を感知する場所:ＡまでＰ波,およびＳ波が到着するまでの時間を,それぞれＴ_P,およびＴ_Sと書けばＴ_P＝∫_O^A(1/Ｖ_P)ｄｓ,Ｔ_S＝∫_O^A(1/Ｖ_S)ｄｓです。

Ｖ_P＞Ｖ_Sですから,Ｔ_S＞Ｔ_Pとなります。

つまり,まず縦揺れのＰ波だけが到着しその後に横揺れＳ波も到着して両方が重ね合わされた大きい揺れが始まるのですね。

特にＯＡ間の到るところで媒質が同じ均等な弾性体であれば,この区間でＶ_P,Ｖ_Sが一定です。

そこで,このときにはＯＡ間の波の進行距離(直線距離でなくてもよい)をｄ_A≡∫_O^AｄｓとするとＴ_P＝ｄ_A/Ｖ_P,Ｔ_S＝ｄ_A/Ｖ_Sです。

これから,引き算するとＴ_S－Ｔ_P＝(ｄ_A/Ｖ_S－ｄ_A/Ｖ_P)＝ｄ_A(1/Ｖ_S－1/Ｖ_P)です。

それ故,Ｔ_SP≡Ｔ_S－Ｔ_P,ｋ≡(1/Ｖ_S－1/Ｖ_P)^-1＝Ｖ_PＶ_S/(Ｖ_P－Ｖ_S)と定義すればｄ_A＝ｋＴ_SPという比例関係の公式を得ます。

　この式は発見者の大森房吉氏(1918)の名前を取って,大森公式(Ohmori's law)と呼ばれます。

　　　　(下図は中越地震本震の地震計記録からです。)

　　　　　

こうした等方的で均質な媒質を伝わる場合,弱いＰ波だけが来て「あ,ひょっとして地震かな？」と思ってから本格的で大きな揺れがくるまでの初期微動の時間をＴ_SPとするとその地点Ａから震源Ｏまでの距離ｄ_AはＴ_SPに比例します。

初期微動時間が長いほど震源までの距離(地理的な距離だけでなく震源深さまでも含めた距離)が長いということが言えます。もしも震源が直下にあれば縦揺れは上下震動,横揺れは左右震動になります。

さて,標高ｈが違う３地点Ａ,Ｂ,Ｃがありその座標がそれぞれ(ｘ_A,ｙ_A,ｈ_A),(ｘ_B,ｙ_B,ｈ_B) (ｘ_C,ｙ_C,ｈ_C)であるとき,上記のような公式に基づいて初期震動時間の観測から観測点から震源までの距離ｄ_A,ｄ_B,ｄ_Cが全てわかったとします。

このとき,未知の震源Ｏの座標を(ｘ_O,ｙ_O,ｈ_O)と仮定すれば地震波が全て直進ならｄ_A²＝(ｘ_A－ｘ_O)²＋(ｙ_A－ｙ_O)²＋(ｈ_A－ｈ_O)²,ｄ_B²＝(ｘ_B－ｘ_O)²＋(ｙ_B－ｙ_O)²＋(ｈ_B－ｈ_O)²,ｄ_C²＝(ｘ_C－ｘ_O)²＋(ｙ_C－ｙ_O)²＋(ｈ_C－ｈ_O)²が成立します。

これは３つの未知数ｘ_O,ｙ_O,ｈ_Oに対する３つの連立２次方程式ですから解くことが可能です。つまり,正確に震源までの距離を観測できる点が３個あれば原理的には震源の位置がわかります。これは３点法と呼ばれるものですね。

普通,震源の高さｈ_Oは海抜で測って負の数であり震源は地中または海中にあります。

大森公式ｄ_A＝ｋＴ_SPは地震に対する公式で大森係数と呼ばれるｋ＝Ｖ_PＶ_S/(Ｖ_P－Ｖ_S)は８(ｋｍ/ｓ)程度の値ですが,この式は雷における音と光の２つのように速さの異なる２つの信号が届く現象では全く同じ公式が成立します。

ただし,雷の場合は光速がｃ～30万km/ｓ,空気中の音速がｖ～340ｍ/ｓでｃ＞＞ｖなのでｋ＝ｃｖ/(ｃ－ｖ)～ｃｖ/ｃ＝ｖですから,事実上ｄ_A～ｖＴなのでほとんど光速は関係なしです。(Ｔは光が見えてから音が聞こえるまでの時間です。)

例えば,ピカッと光ってから雷の音が聴こえるまでの時間が１秒なら,自分から340ｍ程度離れたところ(空中かもしれません)で"放電＝落雷"があったと推測されます。

PS：＞私のマノン・レスコーたちへ　

　病気とか事故に会っていなければどこで誰と遊んでいてもいいけれど,連絡つかない状況だととにかく心配です。大丈夫かなあ。。。。)

　("マゾ, ロリコン, 準デブ専, メガネフェチ＝変態"のＴＯＳＨＩより。。)

(Ｔ.ウッズも聖人君子ではなく適当にストレスを解消しているらしいので安心しました。あれだけの収入を得る能力があれば大勢の家族を養えますね。

　性欲のある健康な男だし,大有名人で品行方正キャラが売り物なら日本だと風俗通いもままならないので女遊びは囲い込むしか仕方ないとして,もしもプロテスタントならマドンナの養子のように大勢を扶養するのはむしろ義務かな？)

PS2:ドクちゃんって日本が好きなんだ。。？？

　国ではなくて故郷(くに)が好きなんでしょ？

　http://news.aimu-net.com/read.cgi/wildplus/1257045443/

PS3:「陳情政治」からいきなり「非陳情政治」へと革命的に移行するのではなく,途中で１人(小沢さん)に権力が集中するシステムを経る。。

　共産革命でも過渡期ではプロレタリア独裁(現実には中国では毛沢東,ロシアではレーニンの１人独裁)です。それはそれでマヌーバとしてありなのでしょうが独裁者がプラトン的(プラトニック)に正しい哲人でいるうちはいいけど,うまく移行できないと"元の木阿弥"か,もっと悪くなります。

　小沢氏が中国に大挙して出かけていって米国にブラフかけて「ポッポ政権」を裏で後押ししてるのかも知れないけど,うまくいくのかねえ。基地問題。。

　朝っぱらから女子アナが「吉宗は何将軍？」と聞かれて「暴れん坊将軍？」と答えていました。それでもいいんだろうけど,ニュース番組なので一応「米(コメ)将軍だろ？」とツッコミたくなりました。。。アハハ

　↑(転載終了)

水滴の成長と蒸発(2)

　さて水滴の成長と蒸発の続きです。

　色々としあわせな状態にかまけてはいますが,基本的にやることは変わっていません。

　水滴が媒質を含む希薄溶液から成る場合について記述します。 

§3.希薄溶液の水滴の凝結と蒸発

①定常状態(平衡状態)の熱と質量の関係

　蒸発の潜熱をＬ_e(単位質量当たり)とすると,凝結によって単位時間当たり解放される熱量はＬ_e(ｄｍ/ｄｔ)です。

そこで,この潜熱を熱源とするｕ＝0での熱方程式は(ｄｑ/ｄｔ)­₀＝－∫_r=aｊ_hｎｄＳ＋Ｌ_e(ｄｍ/ｄｔ)­₀ (31)となります。ただし,ａは水滴半径,ｑは水滴の熱量です。

この場合,熱平衡(定常)の条件:(ｄｑ/ｄｔ)­₀＝0はＬ_e(ｄｍ/ｄｔ)­₀＝∫_r=aｊ_hｎｄＳ,またはＶを水滴の体積として凝結による潜熱の解放が均質に生じるとすれば,∇ｊ_h＝Ｌ_e(ｄｍ/ｄｔ)­₀/Ｖ＝∇ｊ_hです。

　あるいは,ｊ_h＝－ｋ∇Ｔ(25)より

∫_r=aｊ_hｎ_rｄＳ＝－∫_r=aｋ_a(∂Ｔ/∂ｒ)ｄＳ＝－4πａｋ_a(Ｔ_∞－Ｔ_a)

なので,水滴表面温度Ｔ_aが一定条件下では,

Ｌ_e(ｄｍ/ｄｔ)­₀＝∫_r=aｊ_hｎｄＳ＝4πａｋ_a(Ｔ_a－Ｔ_∞)

と表現できます。

ただし,前述したようにｋ_aは水滴の熱伝導度です。

最後の式は,(蒸発に必要な熱量)＝－Ｌ_e(ｄｍ/ｄｔ)­₀＝4πａｋ_a(Ｔ_∞－Ｔ_a),または(凝結によって解放される熱量)＝Ｌ_e(ｄｍ/ｄｔ)­₀＝＝4πａｋ_a(Ｔ_a－Ｔ_∞)＝(伝導により流出する熱量)を示しています。

このＬ_e(ｄｍ/ｄｔ)­₀＝4πａｋ_a(Ｔ_a－Ｔ_∞)なる式は,(ｄｍ/ｄｔ)­₀＝(4πａｋ_a/Ｌ_e)(Ｔ_a－Ｔ_∞)(32),Ｔ_a＝Ｔ_∞＋{Ｌ_e/(4πａｋ_a)}(33)と書き直すことができます。

故に,蒸発水滴:(ｄｍ/ｄｔ)­₀＜0 の温度はまわりよりも低く,凝結水滴:(ｄｍ/ｄｔ)­₀＞0 の温度はまわりよりも高いといえます。

②Clausius-Clapeyronの公式

　さて２相(液相,気相)平衡状態の飽和水蒸気ｅ_sat(Ｔ)の温度保存式を求めます。

　熱力学によれば,平衡の条件は液相,気相のGibbs自由エネルギーについてＧ_w＝Ｇ_v(34),またはμ_w＝μ_vが成立することです。

　Ｇ_w,Ｇ_vは１モル当たりの水,水蒸気のGibbs自由エネルギーです。

　また,μ_w,μ_vは１モル当たりの水,水蒸気の化学ポテンシャルですが,Ｇ_w＝μ_w,Ｇ_v＝μ_vであり,同じものを別記号で表わしているだけです。

　温度ＴがＴ＋ｄＴ,圧力ｐがｐ＋ｄｐになっても(34)の関係が維持されるためには,その際のＧ_w,Ｇ_vの増加分ｄＧ_w,ｄＧ_vについてもｄＧ_w＝ｄＧ_v(35)が成立することが必要です。

　ところで,熱平衡状態ではｄｕ＝Ｔｄｓ－ｐｄｖでＧ＝ｕ＋ｐｖ－Ｔｓですから,ｄＧ＝－ｓｄＴ＋ｖｄｐです。

　ただし,ｕ,ｓはそれぞれ１モル当たりの内部エネルギー,エントロピー,ｖは１モル当たりの体積です。

　したがって,ｄＧ_w＝ｄＧ_v (35)は－ｓ_wｄＴ＋ｖ_wｄｐ＝－ｓ_vｄＴ＋ｖ_vｄｐを意味します。

　つまりｄｐ/ｄＴ＝(ｓ_v－ｓ_w)/(ｖ_v－ｖ_w)(36)を意味します。

　ところが,Ｔ(ｓ_v－ｓ_w)＝Ｌ_eＭ_w (37)ですから,これは

  ｄｐ/ｄＴ＝Ｌ_eＭ_w/{Ｔ(ｖ_v－ｖ_w)}とも書けます。

同じ１モルでは水の体積は水蒸気の体積よりはるかに小さいので,この式の右辺で水蒸気の体積ｖ_vに比して水の体積ｖ_wを無視します。

左辺のｐにｐ＝ｅ_sat(Ｔ),右辺のｖ_vにｖ_v＝ＲＴ/ｅ_satを代入すると,

(1/ｅ_sat)(ｄｅ_sat/ｄＴ)＝Ｌ_eＭ_w/(ＲＴ²)(38)を得ます。

これが,有名なClausius-Clapeyronの公式(クラウジウス・クラペイロンの公式)です。

念のため,改めて追記するとＬ_eは単位質量当たりの蒸発の潜熱,Ｍ_wは水の分子量です。

　　　　　(↓下図はネット検索で入手した図の転載です。)

　　

※(注):エンタルピー:ｈ≡ｕ＋ｐｖを導入するとｄｕ＝Ｔｄｓ－ｐｄｖよりｄｈ＝Ｔｄｓ＋ｖｄｐです。

　２相平衡にあるときの圧力ｐは,ｄｐ/ｄＴ＝(ｓ_v－ｓ_w)/(ｖ_v－ｖ_w)(36)を満たしますから,ｄｈ＝Ｔｄｓ＋ｖｄｐ＝Ｔｄｓ＋ｖ(ｄｐ/ｄＴ)ｄＴと書けます。

そこで,ｈ_v－ｈ_w＝Ｔ(ｓ_v－ｓ_w)＋∫ｖ(ｄｐ/ｄＴ)ｄＴですが,相変化が定温(ｄＴ＝0)の下で移行する場合の潜熱はエンタルピーの変化として表現されます。

つまり,ｈ_v－ｈ_w＝Ｔ(ｓ_v－ｓ_w)＝Ｌ_eＭ_wです。

したがって,ｐ＝ｅ_satに対するClausius-Clapeyronの公式:

(1/ｅ_sat)(ｄｅ_sat/ｄＴ)＝Ｌ_eＭ_w/(ＲＴ²)(38)は,

(1/ｅ_sat)(ｄｅ_sat/ｄＴ)＝ｄ(lnｅ_sat)＝(ｈ_v－ｈ_w)/(ＲＴ²)とも表現されることがわかります。(注終わり)※

　蒸発の潜熱:Ｌ_eがＴに依らず一定値を取ると仮定して,

　(1/ｅ_sat)(ｄｅ_sat/ｄＴ)＝Ｌ_eＭ_w/(ＲＴ²)(38)を積分すると,

　ｅ_sat(Ｔ_a)＝ｅ_sat(Ｔ_∞)exp[{Ｌ_eＭ_w/(ＲＴ_∞)}(Ｔ_a－Ｔ_∞)/Ｔ_a](39)

　が得られます。

③Kelvinの公式,Raoultの規則

　半径ａの純水滴と湿潤空気の混合系で,気相と液相が平衡にあるためにはμ_w＝μ_v(40),かつｐ_w＝ｐ＋2σ/ａ(41)が成立する必要があります。

　ただし,σは表面張力です。

(このブログのバックナンバーで表面張力について言及しているのは2009年4/28の「水の波(2)(浅水波,深水波,表面張力波)」だけですね。

　一般の球でない水面の場合にはｐ_w＝ｐ＋σ/ａ₁＋σ/ａ₂です。)

　ここで,ｐ,Ｔの変化の下で,μ_w/Ｔ＝μ_v/Ｔの関係が維持されるためにはｄ(μ_w/Ｔ)＝ｄ(μ_v/Ｔ)(42)が成立することが必要です。

　ところが(∂μ_w/∂ｐ_w)_T＝ｖ_w,[∂(μ_w/Ｔ)/∂Ｔ]_pw＝－ｈ_w/Ｔ²より,

ｄ(μ_w/Ｔ)＝(－ｈ_w/Ｔ²)ｄＴ＋(ｖ_w/Ｔ)ｄｐ_wです。

　ただしｈはエンタルピー:ｈ≡ｕ＋ｐｖです。

　一方,水蒸気のみの１成分系では,Ｇ_v＝μ_vであり,

混合系ではＧ_xa＝ｘ_aＧ_a＋ｘ_vＧ_v(ｘ_a＋ｘ_v＝1)です。

ここで,ｘ_aの添字ａは(乾燥)空気(air)のａ,

ｘ_vのｖは水蒸気(vapor)のｖです。

そして,水蒸気の分圧はｅ＝ｘ_vｐですから,Ｇ_v＝μ_v＝μ_v⁺(Ｔ)＋ＲＴlnｅよりμ_v＝μ_v⁺(Ｔ)＋ＲＴlnｐ＋ＲＴlnｘ_vです。

これは,μ_v＝μ_v(Ｔ,ｐ,ｘ_v)＝μ_v⁰(Ｔ,ｐ)＋ＲＴlnｘ_vと書けます。

μ_v⁰(Ｔ,ｐ)は空気が純水蒸気(ｘ_v＝1)で,その気圧がｐのときの化学ポテンシャル(Gibbs自由エネルギー)です。

(∂μ_v/∂ｐ)_T＝ｖ_v,[∂(μ_v/Ｔ)/∂Ｔ]_pw＝－ｈ_v/Ｔ²より,ｄ(μ_v/Ｔ)＝(－ｈ_v/Ｔ²)ｄＴ＋(ｖ_v/Ｔ)ｄｐ＋Ｒｄ(lnｘ_v)ですから,

条件ｄ(μ_w/Ｔ)＝ｄ(μ_v/Ｔ)(42)は,(－ｈ_w/Ｔ²)ｄＴ＋(ｖ_w/Ｔ)ｄｐ_w＝(－ｈ_v/Ｔ²)ｄＴ＋(ｖ_v/Ｔ)ｄｐ＋Ｒｄ(lnｘ_v)と書けます。

また,ｐ_w＝ｐ＋2σ/ａ(41)よりｄｐ_w＝ｄｐ＋2ｄ(σ/ａ)でありｘ_v＝ｅ/ｐですから,{－(ｈ_v－ｈ_w)/Ｔ²}ｄＴ＋{(ｖ_v－ｖ_w)/Ｔ²}ｄｐ－(2ｖ_w/Ｔ)ｄ(σ/ａ)＋Ｒｄ{ln(ｅ/ｐ)}＝0 です。

つまり,(－Ｌ_eＭ_w/Ｔ²)ｄＴ＋{(ｖ_v－ｖ_w)/Ｔ²}ｄｐ－(2ｖ_w/Ｔ)ｄ(σ/ａ)＋Ｒｄ{ln(ｅ/ｐ)}＝0 (43)が成立することが必要です。

ここでのｅはｅ＝ｅ_sat(Ｔ)を意味しますが,今の混合系では定温,定圧下で半径ａの変化に対してｅ＝ｅ_sat(Ｔ)が変化すると予想されるため,ｅ＝ｅ_sat(Ｔ)を,改めてｅ＝ｅ^a _sat(Ｔ)と書きます。

定温,定圧下ではｄＴ＝ｄｐ＝0なので,(43)は－(2ｖ_w/Ｔ)ｄ(σ/ａ)＋Ｒｄ{ln(ｅ^a_sat/ｐ)}＝0(44)となります。

ｐ＝一定ですから,ln{ｅ^a_sat(Ｔ)}＝2ｖ_wσ/(ＲＴａ)＋const.より,

ｅ^a _sat(Ｔ)＝ｅ^∞_sat(Ｔ)exp{2Ｍ_wσ/(ＲＴρ_wａ)}(45) を得ます。

これがKelvinの公式です。

　一方,溶解液に対するRaoult(ラウール)の規則から,ｅ'_satをこの溶液の溶媒(水)の飽和蒸気圧とし,ｍ_w,ｍ_s,ｍを溶媒,溶質,(溶媒＋溶質)の質量,Ｍ_w,Ｍ_sを溶媒,溶質の分子量とすると,

　ｉをvan'tHoff factorとして,　

　ｅ'^∞_sat/ｅ^∞_sat＝ｎ_w/(ｎ_w＋iｎ_s)＝ｍ_wＭ_s/(ｍ_wＭ_s＋iｍ_sＭ_w)

　＝[1＋iｍ_sＭ_w/{(ｍ－ｍ_s)Ｍ_s}]^-1(46)と書けます。

(45),(46)よりｅ^∞_satをｅ_satと書けば,ｅ'^a_sat(Ｔ)＝

ｅ_sat(Ｔ)exp{2ｖ_wσ/(ＲＴρ_wａ)}[1＋iｍ_sＭ_w/{(ｍ－ｍ_s)Ｍ_s}]^-1(47)

です。

④溶解液の蒸発,凝結

　ここまでの公式を以下のように要約します。 

すなわち,ａ(ｄａ/ｄｔ)₀＝{Ｄ_vＭ_w/(ρ_wＲ)}{ｅ_∞/Ｔ_∞－ｅ'^a_sat(Ｔ_a)/Ｔ_a}(7)',ａ(ｄａ/ｄｔ)­₀＝{ｋ_a/(Ｌ_eρ_w)}(Ｔ_a－Ｔ_∞)(32)',

および,

ｅ_sat(Ｔ_a)＝ｅ_sat(Ｔ_∞)exp[{Ｍ_wＬ_e/(ＲＴ_∞)}(Ｔ_a－Ｔ_∞)/Ｔ_a](39)',　ｅ'^a_sat(Ｔ_a)＝ｅ_sat(Ｔ_a)exp{2Ｍ_wσ/(ＲＴ_aρ_wａ)}[1＋ｉｍ_sＭ_w/{(ｍ－ｍ_s)Ｍ_s}]^-1 (47)'です。

そして,(32)'をＴ_a＝Ｔ_∞(1＋δ),

δ≡{Ｌ_eρ_w/(ｋ_aＴ_∞)}ａ(ｄａ/ｄｔ)­₀ (48)と表現すれば,

(7)',(39)',(47)',(48)によって,

ａ(ｄａ/ｄｔ)₀＝{Ｄ_vＭ_w/(ρ_wＲ)}(ｅ_∞/Ｔ_∞－ｅ_sat(Ｔ_∞)exp[Ｍ_wＬ_eδ/{ＲＴ_∞(1＋δ)}＋2Ｍ_wσ/{ＲＴ_∞ρ_wａ(1＋δ)}][1＋iｍ_sＭ_w/{(ｍ－ｍ_s)Ｍ_s}]^-1) (49)が得られます。

δ＝{Ｌ_eρ_w/(ｋ_aＴ_∞)}ａ(ｄａ/ｄｔ)­₀ですが,δ＜＜1なので

exp[Ｍ_wＬ_eδ/{ＲＴ_∞(1＋δ)}～ 1＋Ｍ_wＬ_eδ/(ＲＴ_∞),

かつexp[2Ｍ_wσ/{ＲＴ_∞ρ_wａ(1＋δ)}]～ exp{2Ｍ_wσ/(ＲＴ_∞ρ_wａ)}

と近似できます。

さらに,

exp{2Ｍ_wσ/(ＲＴ_∞ρ_wａ)}[1＋iｍ_sＭ_w/{(ｍ－ｍ_s)Ｍ_s}]^-1～ 1＋ｙ(50),

ｙ≡2Ｍ_wσ/(ＲＴ_∞ρ_wａ)－iｍ/{Ｍ_s(ｍ－ｍ_s)}(51)と置くことで,

近似式:ｅ'^a _sat(Ｔ_a)＝ｅ _sat(Ｔ_∞)(1＋ｙ)(52)

を得ます。

以上から,

ａ(ｄａ/ｄｔ)₀＝{Ｄ_vＭ_w/(ρ_wＲ)}{ｅ_sat(Ｔ_∞)/Ｔ_∞}][ｅ_∞/ｅ_sat(Ｔ_∞)－{1/(1＋δ)}{1＋Ｍ_wＬ_eδ/(ＲＴ_∞)}(1＋ｙ)]です。

さらに,δの２次以上とδｙを無視すると,

ａ(ｄａ/ｄｔ)₀＝{Ｄ_vＭ_wｅ_sat(Ｔ_∞)/(ρ_wＲＴ_∞)}[ｅ_∞/ｅ_sat(Ｔ_∞)－(1＋ｙ)－δ{Ｌ_eＭ_w/(ＲＴ_∞)－1}]

です。

そこで,δ＝{Ｌ_eρ_w/(ｋ_aＴ_∞)}ａ(ｄａ/ｄｔ)­₀より,

ａ(ｄａ/ｄｔ)₀[ρ_wＲＴ_∞/{Ｄ_vＭ_wｅ_sat(Ｔ_∞)}＋{Ｌ_eρ_w/(Ｔ_∞ｋ_a)}{Ｌ_eＭ_w/(ＲＴ_∞)－1}]

～{ｅ_∞－ｅ_sat(Ｔ_∞)(1＋ｙ)}/ｅ_sat(Ｔ_∞),

あるいは,

ａ(ｄａ/ｄｔ)­₀～[{ｅ_∞－ｅ'_sat(Ｔ_∞)}/ｅ_sat(Ｔ_∞)]/[{Ｌ_eρ_w/(Ｔ_∞ｋ_a)}{Ｌ_eＭ_w/(ＲＴ_∞)－1}＋ρ_wＲＴ_∞/{Ｄ_vＭ_wｅ_sat(Ｔ_∞)}](53)

を得ます。

　ただし,ｅ'_sat(Ｔ_∞)＝ｅ_sat(Ｔ_∞)exp{2Ｍ_wσ/(ＲＴ_∞ρ_wａ)}[1＋ｉｍ_sＭ_w/{(ｍ－ｍ_s)Ｍ_s}]^-1(54)です。もちろん,ｍ＝4πａ³ρ/3です。

　ここまでは,もっぱら(ｄｍ/ｄｔ)­₀,ａ(ｄａ/ｄｔ)­₀の静止水滴に関する話でしたが,運動水滴の場合には単にＤ_v→ Ｓ_hＤ_v/2,ｋ_a→ Ｎ_uｋ_aとすればいいだけです。

内容を思い出すのにてこずりましたが,ノートはここで終わっているので,この項目についてはこれで終わりにします。

PS:一昨年の暮れ頃には観音様の御出現,今年は天使の御出現と,ここのところ浮世離れしたことが続いています。

翌日(21日)も天使がそばにいてとてもよかったです。

夜には今年最後の手話講習(初心コース)に出席,終わって21時頃から椎名町駅近くで応用,上級コース合同の忘年会でした。

初対面のろう(あ)の人との手話もはずみましたが,翌日が早いので夜中の零時半頃には帰宅してそのまま就寝しました。

今日(22日)は,仕事をお休みして朝から眼底出血していた右目を含む両目眼底の総合的検査のため,予約していた帝京病院眼科診療に行ってきます。

しかし,今は右目も左目と同じくらい見えるようですから,医者の予想通り幸い自然に出血が止まって治癒しているようで,レーザー手術もしないで済みそうです。

水滴の成長と蒸発(1)

　前の一連の記事のプルーム上昇式について書いていた1983年８月の覚書きノートには続きがあります。それは水滴の成長と蒸発という項目です。

これも1000ページ近い"英文のANLモデル(アルゴンヌ国立研究所模型;後のSACTI)のドキュメントを全部翻訳して内容を完全に把握する"という当時のジョブの一環として勉強したものです。

冷却塔(cooling tower)からの温排気排出の影響には,水蒸気による白煙発生の他に,それが周辺大気で冷却され凝結した水滴の飛散による濃霧,また冷却水が海水のように塩分などを含む場合には飛散水滴による農作物等への塩害etc.も考えられるからです。

そこで,この項目についても覚書きとして記述しておきます。

(当時は,仕事の都合上,タイムリミットもあって,とにかく当座の計算アルゴリズム等に必要な式に到達するのが目的だったので,途中のプラス,マイナスなど細かいところはアバウトだったりテキトーだったりなところもあるので,あとで修正が必要なのです。)

以下は本文です。

§1.拡散流束による水滴質量の増減

　水蒸気の密度をρ_v,その密度流束(局所流束密度)をｊ_vとします。

この局所流束密度ｊ_vは,移流風速(水蒸気を含む空気の平均流速)をｕとすると,そのゆらぎ(fluctuation)に相当する拡散流束を加えて,

ｊ_v＝ρ_vｕ－Ｄ_v∇ρ_vと表現されます。

ただし,Ｄ_vは拡散係数です。

この流束を用いた水蒸気の質量保存を示す連続の方程式は,

∂ρ_v/∂ｔ＋∇ｊ_v＝0 or ∂ρ_v/∂ｔ＝－∇ｊ_vです。

そして,ｊ_v＝ρ_vｕ－Ｄ_v∇ρ_vにより∇ｊ_v＝∇(ρ_vｕ－Ｄ_v∇ρ_v)です。

しかし,この水蒸気を含む空気は非圧縮(または,高々Boussinesq近似)を仮定しているため,∇ｕ＝0 を満たしますから,

これは,∇ｊ_v＝ｕ∇ρ_v－∇(Ｄ_v∇ρ_v) と書けます。

さらに,拡散係数Ｄ_vは空間位置によらない定数であると仮定すると,

∇ｊ_v＝ｕ∇ρ_v－Ｄ_v∇²ρ_vです。

結局,上記連続の方程式:∂ρ_v/∂ｔ＝－∇ｊ_vは,

∂ρ_v/∂ｔ＋ｕ∇ρ_v＝Ｄ_v∇²ρ_v (1)　を意味します。

これは,Fick型の(移流)拡散方程式(diffusion equation)と呼ばれる方程式です。

　一方,ρ_wを液体の水の密度として半径ａの球水滴の質量をｍとするとｍ＝∫ρ_wｄＶ＝4πａ³ρ_w/3です。

水滴の中心を原点:Ｏとするような極座標で考えると,ｒ＞ａの領域では水蒸気に対して,上述の拡散方程式:∂ρ_v/∂ｔ＋ｕ∇ρ_v＝Ｄ_v∇²ρ_v (1)が成立しますが,ｒ＜ａの球内は(1)が成立しない特異領域と考えられます。

　そして,空気は水滴表面から内部に浸入することはないので,表面ｒ＝ａでは速度の表面への垂直成分:ｕ_r＝0 であるべきですが,現実の水滴表面ではより強い粘着条件:ｕ＝0 が満たされていると思われます。

　つまり,表面ｒ＝ａではｕ＝0 なる境界条件が満たされると考えます。

　水滴の質量保存から,明らかにｄｍ/ｄｔ＝－∫_r=aｊ_vｎ_rｄＳですが,

　上記境界条件を要求すると,ｒ＝ａでは

　ｊ_v＝ρ_vｕ－Ｄ_v∇ρ_v＝－Ｄ_v∇ρ_vなので,

　ｄｍ/ｄｔ＝Ｄ_v∫_r=a(∂ρ_v/∂ｒ)_r=aｄＳ (2)です。

　特に,表面だけでなく全空間でｕ＝0 を満たす,水滴に対して静止した空気中でのｄｍ/ｄｔを(ｄｍ/ｄｔ)₀と書けば,場は中心対称ですから,

　(ｄｍ/ｄｔ)₀＝4πａ²Ｄ_v(ｄρ_v/ｄｒ)_r=a (3)と書けます。

他方,ｕ＝0 なら∂ρ_v/∂ｔ＝Ｄ_v∇²ρ_vです。さらにｔ→ ∞ の極限では拡散場は定常(∂ρ_v/∂ｔ＝0)になり,しかも球対称ですから,

∇²ρ_v＝ｒ^-1ｄ²(ｒρ_v(ｒ)/ｄｒ²＝0 です。　

これを,ｒ＝ａでρ_v＝ρ_v(ａ);ｒ＝∞でρ_v＝ρ_∞の境界条件下で解けば

ρ_v(ｒ)＝ρ_∞＋ａ{ρ_v(ａ)－ρ_∞}/ｒ (4)を得ます。

これにより,(ｄρ_v/ｄｒ)_r=a＝{ρ_∞－ρ_v(ａ)}/ａ (5)ですから,

ｕ＝0の静止大気中での水滴の質量保存の式は,　

(ｄｍ/ｄｔ)₀＝4πａ²Ｄ_v(ｄρ_v/ｄｒ)_r=a＝4πａＤ_v{ρ_∞－ρ_v(ａ)}

と書けます。

水蒸気が理想気体なら,それが従う状態方程式:ｅ＝ρ_vＲＴ/Ｍ_wから,

ρ_∞＝ｅ_∞Ｍ_w/(ＲＴ_∞),ρ_v(ａ)＝ｅ_aＭ_w/(ＲＴ_a)を代入して,

(ｄｍ/ｄｔ)₀＝(4πａＤ_vＭ_w/Ｒ)(ｅ_∞/Ｔ_∞－ｅ_a/Ｔ_a) (6)

です。

または,ｍ＝4πａ³ρ_w/3なる式により,

方程式:ａ(ｄａ/ｄｔ)₀＝{Ｄ_vＭ_w/(ρ_wＲ)}(ｅ_∞/Ｔ_∞－ｅ_a/Ｔ_a) (7)

を得ます。

ただし,Ｔ_aは水滴表面(ｒ＝ａ)の温度,ｅは水蒸気圧,Ｒはモル気体定数(Ｒ～ 8.31J/(molＫ)),Ｍ_wは水の分子量(Ｍ_w～ 18ｇ/mol)です。

しかし,水滴半径ａが分子半径と比較してそれほど大きくない場合には,空気や水蒸気を連続体とみなす巨視的流体力学の適用範囲の外に出ると考えられるため,必要な小水滴に対する補正をします。

以下,Fuchs(1959)によって与えられた方法に従います。

　拡散方程式とその解は水滴表面からの水蒸気の反跳長さ(recoiling length):Δ_vを超えないなら正しくないと考えます。

すなわち,ｒ≧ａ＋Δ_vではρ_v(ｒ)＝ρ_∞＋Ａ/ｒ(8)が成立しますが,

ａ≦ｒ≦ａ＋Δ_ｖの境界領域では水蒸気流束を構成する分子はMaxwell-Boltzmannの速度分布に従う粒子として気体運動論的に扱うことにします。

水分子の質量をｍ_wとし,分子の平均速度をｖ＝(ｖ_x,ｖ_y,ｖ_z),その大きさをｖ≡|ｖ|,Boltzmann定数をｋ_B≡Ｒ/Ｎ_A(Ｎ_AはAvogadro数)とすると,温度がＴのときの分子の平均速さは,

＜ｖ＞＝(8ｋ_BＴ)^1/2/(πｍ_w)^1/2＝(8ＲＴ)^1/2/(πＭ_w)^1/2 (9)

で与えられます。

※(注):＜ｖ＞＝∫ｖexp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄ³ｖ/∫exp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄ³ｖ＝∫₀^∞ｖ³exp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄｖ/∫₀^∞ｖ²exp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄｖ　です。

ところが,∫₀^∞exp(－αｖ²)ｄｖ＝(1/2)(π/α)^1/2ですから,両辺をαで微分することで,∫₀^∞ｖ²exp(－αｖ²)ｄｖ＝(1/4)π^1/2/α^3/2を得ます。

さらに,これをαで微分すると,∫₀^∞ｖ⁴exp(－αｖ²)ｄｖ＝(3/8)π^1/2/α^5/2です。

一方,∫₀^∞ｖ³exp(－αｖ²)ｄｖ＝(1/2)∫₀^∞ｕexp(－αｕ)ｄｕ＝(1/2)(1/α²)より,∫₀^∞ｖ³exp(－αｖ²)ｄｖ/∫₀^∞ｖ²exp(－αｖ²)ｄｖ＝2/(πα)^1/2です。

＜ｖ＞＝∫₀^∞ｖ³exp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄｖ/∫₀^∞ｖ²exp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄｖの右辺は,上記のα＝ｍ_w/(2ｋ_BＴ)の場合に当たりますから,

　＜ｖ＞＝(8ｋ_BＴ)^1/2/(πｍ_w)^1/2＝(8ＲＴ)^1/2/(πＭ_w)^1/2です。

一方,二乗平均速度を評価するのなら,

＜ｖ²＞^1/2＝[∫ｖ²exp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄ³ｖ]^1/2＝[∫₀^∞ｖ⁴exp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄｖ/∫₀^∞ｖ²exp{－ｍ_wｖ²/(2ｋ_BＴ)}ｄｖ]^1/2

＝{3/(2α)}^1/2＝(3ｋ_BＴ/ｍ_w)^1/2＝(3ＲＴ/Ｍ_w)^1/2です。

後者はエネルギーの等分配則＜ｍｖ²/2＞＝(3/2)ｋ_BＴから予想される結果に一致しています。

(8/π)^1/2と3^1/2 には大差がないので,＜ｖ＞と＜ｖ²＞^1/2は同じオーダーの量です。(注終わり)※

　さて,ｎを分子濃度,つまり単位体積当たりの分子数とします。

分子の従うMaxwell-Boltzman分布を,

"ｖ＝(ｖ_x,ｖ_y,ｖ_z)とｖ＋ｄｖ＝(ｖ_x＋ｄｖ_x,ｖ_y＋ｄｖ_y,ｖ_z＋ｄｖ_z)を対角線とする速度空間の微小直方体内に,単位体積当たりＦ(ｖ²)ｄｖ_xｄｖ_yｄｖ_z個の分子があること"で表現すれば,

∫Ｆ(ｖ²)ｄｖ_xｄｖ_yｄｖ_z＝ｎ,∫ｖＦ(ｖ²)ｄｖ_xｄｖ_yｄｖ_z/ｎ＝＜ｖ＞

です。

すると,速度空間では速度ｖがｖ～ ｖ＋ｄｖの間にあり,かつ座標空間の任意の向きの平面の微小立体角:ｄΩ＝cosθｄθｄφに衝突する単位時間,単位面積当たりの分子数は,

{ｖsinθｄΩ/(4π)}Ｆ(ｖ²)ｄｖ_xｄｖ_yｄｖ_z

＝ｖＦ(ｖ²)ｄｖ_xｄｖ_yｄｖ_zsinθcosθｄθｄφ/(4π)

で与えられます。

ただし,この平面がｖとなす角がθであるように極座標軸を取ってます。

平面の一方の側の単位面積に衝突する単位時間当たりの分子の総数は,これを角度θについて 0 からπ/2まで,φについて 0 から2πまで,さらにｖについて 0 から∞まで積分した値です。

∫₀^π/2sinθcosθｄθ∫₀^2πｄφ＝πですから,この値:[∫ｖＦ(ｖ²)ｄｖ_xｄｖ_yｄｖ_z∫₀^π/2sinθcosθｄθ∫₀^π/2ｄφ]/(4π)は,(1/4)∫ｖＦ(ｖ²)ｄｖ_xｄｖ_yｄｖ_z＝(1/4)ｎ＜ｖ＞に等しいことがわかります。

したがって,(分子数密度流束)＝(任意の面に衝突する単位時間,単位面積当たりの分子数)は,　

(1/4)ｎ＜ｖ＞ (10)　と書けます。

右辺の因子:＜ｖ＞は先に書いたように,

＜ｖ＞＝(8ｋ_BＴ)^1/2/(πｍ_w)^1/2＝(8ＲＴ)^1/2/(πＭ_w)^1/2 (9)

で与えられます。

また,水蒸気の温度を改めてＴの代わりにＴ_gと書くと,水蒸気はｅ＝ρ_wｋ_BＴ_g/ｍ_w＝ｎｋ_BＴ_gなる状態方程式に従います。

故に,水滴の周囲の水蒸気の分子数密度流束は

ｗ^↓＝(1/4)ｎ＜ｖ＞＝(1/4){ｅ/(ｋ_BＴ_g)}(8ｋ_BＴ_g)^1/2/(πｍ_w)^1/2

＝ｅ/(2π_wｋ_BＴ_g)^1/2です。

他方,水滴表面に付着する水蒸気の分子数密度流束をｗ^↑と書き,平衡状態にあるときのｗ^↑の水滴の周囲の分子数密度流束ｗ^↓に対する比を凝結定数と呼んでα_cとします。

すると,ｗ^↑＝α_cｗ^↓です。

そうすると,水滴表面の分子数密度流束は,

ｗ^↑＝α_cｅ/(2πｍ_wｋ_BＴ_a)^1/2 (13)

で与えられることになります。

実験によれば,α_cの値は,α_c≒0.01～ 0.07です。

ただし,水蒸気分子の平均速さは,＜ｖ＞＝(8ＲＴ_g)^1/2/(πＭ_w)^1/2でも

分子速度ｖの分布は等方的で,そのベクトル平均はゼロと考えられるので

分子数密度流束ｗ^↓の向きはあらゆる方向について対称です。

水滴表面の水蒸気圧:ｅが飽和水蒸気圧:ｅ_satに等しいときに"定常状態＝蒸発平衡"に達するとすれば,

平衡時の水滴表面での飽和分子数流束:ｗ^↑_satは,　

ｗ^↑_sat＝α_cｅ_sat/(2πｍ_wｋ_BＴ_a)^1/2 (14)　です。

ただし,Ｔ_aは水滴表面の絶対温度です。

それ故,流入する分子数流束と流出する分子数流束は互いに独立であるとしてＴ_g～ Ｔ_a ＝Ｔと置けば,水滴に流入する総流束密度の大きさは

ｗ_sat＝ｗ^↑－ｗ^↑_sat＝α_c(ｅ－ｅ_sat)/(2πｍ_wｋ_BＴ)^1/2(15)

と書けます。

この(15)式は,ｅ＜ｅ_satなら水滴が蒸発し,ｅ≧ｅ_satなら蒸発しないことを意味しています。

さて,水滴中心からの距離ｒにおける水蒸気の分子数密度をＮ_v(ｒ),密度をρ_v(ｒ)とすればρ_v(ｒ)＝ｍ_wＮ_v(ｒ)です。

したがって,ｒ＝ａ＋Δ_ｖにおける水蒸気が,流出入する水蒸気を代表するとすれば,その水滴表面へと向かう凝結流束密度は,

α_cｍ_wＮ_v(ａ＋Δ_ｖ)＜ｖ_v＞/4＝α_c＜ｖ_v＞ρ_v(ａ＋Δ_ｖ)/4　です。

他方,水滴表面ｒ＝ａで蒸発に向かう蒸発流束密度は,

α_c＜ｖ_v＞ρ_v(ａ)/4です。

そこで,表面での総蒸発流束Ｊ_v(ａ)は,

Ｊ_v(ａ)＝4πａ²ｊ_v(ａ)＝πａ²α_c＜ｖ_v＞{ρ_v(ａ)－ρ_v(ａ＋Δ_c)} (16)

で与えられることがわかります。

一方,冒頭で書いたように,巨視的流体力学の見方では,

ｒ＝ａでの流速密度はｊ_v＝－Ｄ_v∇ρ_vなので,総流束は

Ｊ_v(ｒ)＝4πａｒ²ｊ_v(ｒ)＝－4πｒ²(ｄρ_v/ｄｒ)Ｄ_v

です。

ところで,今の気体運動論的考察では,

ｒ≧ａ＋Δ_ｖでρ_v(ｒ)＝ρ_∞＋Ａ/ｒ (8)が成立する,

と仮定しているので,

ｒ＝ａ＋Δ_ｖにおいてはρ_v(ａ＋Δ_v)＝ρ_∞＋Ａ/(ａ＋Δ_v)です。

故に,(ｄρ_v/ｄｒ)_r=a+Δv＝－Ａ/(ａ＋Δ_c)²ですから,

Ｊ_v(ａ＋Δ_v)＝－4π(ａ＋Δ_v)²(ｄρ_v/ｄｒ)_r=a+ΔvＤ_v＝4πＤ_vＡ (17)

を得ます。

"定常状態＝蒸発平衡"では,Ｊ_v(ａ＋Δ_v)＝Ｊ_v(ａ) (18)であるべきで,

ρ_v(ａ＋Δ_v)＝ρ_∞＋Ａ/(ａ＋Δ_v)ですから(16),(17)より,

4πＤ_vＡ＝πａ²α_c＜ｖ_v＞{ρ_v(ａ)－ρ_v(ａ＋Δ_c)}

＝πａ²α_c＜ｖ_v＞{ρ_v(ａ)－ρ_∞－Ａ/(ａ＋Δ_v)}

が成立します。

したがって,

[ａ²α_c＜ｖ_v＞/(ａ＋Δ_v)＋4Ｄ_v]Ａ＝ａ²α_c＜ｖ_v＞{ρ_v(ａ)－ρ_∞},

および＜ｖ_v＞＝(8ｋ_aＴ_a)^1/2/(πｍ_w)^1/2＝(8ＲＴ_a)^1/2/(πＭ_w)^1/2 (9)

が成立します。

以上から,

Ａ＝ａ{ρ_v(ａ)－ρ_∞}/[ａ/(ａ＋Δ_v)＋4Ｄ_v/(ａα_c＜ｖ_v＞)]

＝ａ{ρ_v(ａ)－ρ_∞}/[ａ/(ａ＋Δ_v)＋Ｄ_v{2πＭ_w/(ＲＴ_a)}^1/2/(ａα_c)] (19)　を得ます。

これは,冒頭で述べた巨視的流体力学の考察:

(4)ρ_v(ｒ)＝ρ_∞＋ａ{ρ_v(ａ)－ρ_∞}/ｒ,

および(ｄｍ/ｄｔ)₀＝4πａＤ_v{ρ_∞－ρ_v(ａ)}から導かれた式,

(6)(ｄｍ/ｄｔ)₀＝(4πａＤ_vＭ_w /Ｒ)(ｅ_∞/Ｔ_∞－ｅ_a/Ｔ_a),

または(7)ａ(ｄａ/ｄｔ)₀＝{Ｄ_vＭ_w /(ρ_wＲ)}(ｅ_∞/Ｔ_∞－ｅ_a/Ｔ_a)

を次のように補正することに相当します。

すなわち,(4)はａ{ρ_v(ａ)－ρ_∞}→ Ａ＝ａ{ρ_v(ａ)－ρ_∞}/[ａ/(ａ＋Δ_v)＋Ｄ_v{2πＭ_w/(ＲＴ_a)}^1/2/(ａα_c)]として,

ρ_v(ｒ)＝ρ_∞＋ａ{ρ_v(ａ)－ρ_∞}/ｒ→ (4)'ρ_v(ｒ)＝ρ_∞＋Ａ/ｒ

と補正します。

そこで,(6),および(7)は拡散係数をＤ_v→ Ｄ_v'＝Ｄ_v/[ａ/(ａ＋Δ_v)＋Ｄ_v{2πＭ_w/(ＲＴ_a)}^1/2/(ａα_c)](20)として,

それぞれ(ｄｍ/ｄｔ)₀＝(4πａＤ_v'Ｍ_w/Ｒ)(ｅ_∞/Ｔ_∞－ｅ_a/Ｔ_a)(6)',

およびａ(ｄａ/ｄｔ)₀＝{Ｄ_v'Ｍ_w/(ρ_wＲ)}(ｅ_∞/Ｔ_∞－ｅ_a/Ｔ_a)(7)'

に補正します。　

これで小水滴に対して補正された質量保存式を得ました。

　次に,水蒸気を含む空気に対して運動している水滴に対する補正を考察します。

　これは,水滴が運動する代わりに静止水滴に対して空気が運動する,つまりｕ≠0 の等価な問題として設定できます。

この場合も,ｒ＞ａにおける拡散は定常(∂ρ_v/∂ｔ＝0)に達しているとすれば,質量保存の拡散方程式はｕ∇ρ_v＝Ｄ_v∇²ρ_vです。

水滴表面での流束密度がｊ_v＝－Ｄ_v∇ρ_vなので,

やはりｄｍ/ｄｔ＝－∫_r=aｊ_vｎ_rｄＳ＝∫_r=aＤ_v(∂ρ_v/∂ｒ)ｄＳです。

　ここで,平均通風係数と呼ばれる量ｆ_vを,

ｆ_v≡ｄｍ/ｄｔ/(ｄｍ/ｄｔ)₀(21)で,

　平均輸送係数ｋ_vをｋ_v≡－∫_r=aＤ_v(∂ρ_v/∂ｒ)ｄＳ/[4πａ²{ρ_v(ａ)－ρ_∞}] (22)　で定義します。

　すると,(2)'と(21),(22)からｋ_v＝Ｄ_vｆ_v/ａです。

さらに,無次元のSherwood数Ｓ_hを次式で定義します。

すなわち,Ｓ_h≡2ｋ_vａ/Ｄ_v＝2ｆ_v (23)です。

こう定義すると,ｄｍ/ｄｔ＝Ｓ_hＤ_v2πａ(ρ_∞－ρ_sat) (24)です。

ｕ＝0 のときにはｆ_v＝1 or Ｓ_h＝2です。

ただし,これらは単に定義した,または言葉を置き換えただけで何ら新しいことを述べていません。文献を読む際の便宜だけでしょう。

§2.水滴中の熱の移動

　ｋを媒質中の熱伝導度,Ｔを絶対温度とすると,Fourierによれば熱流束ベクトル:ｊ_hはｊ_h＝－ｋ∇Ｔ (25)で与えられます。

そして,一般的な流体中のエネルギー保存の式は,粘性による散逸が無視できる場合には流体単位質量あたりのエントロピーをｓとしてρＴ(∂ｓ/∂ｔ＋ｕ∇ｓ)＝－∇ｊ_hと書けます。ρは流体の密度です。

今の場合のように,媒質流体(水,空気)が非圧縮と見なせるときには

∇ｕ＝0,かつＴｄｓ＝ｃ_pｄＴ(ｃ_pは単位質量当たりの定圧比熱)なので,

これはρｃ_p(∂Ｔ/∂ｔ＋ｕ∇Ｔ)＝κ∇²Ｔ

or ∂Ｔ/∂ｔ＋ｕ∇Ｔ＝κ∇²Ｔ(26)　と変形されます。

ただし,κ≡ｋ/(ρｃ_p)でκは熱拡散係数と呼ばれる量です。

　ｊ_h＝－ｋ∇Ｔ(25)より,水滴の熱量をｑとし,乾燥空気の熱伝導度をｋ_aとすると,

　ｄｑ/ｄｔ＝－∫_r=aｊ_hｎ_rｄＳ＝∫_r=aｋ_a(∂Ｔ/∂ｒ)ｄＳ　(27)

　です。

特に,ｕ＝0(水滴に対して水蒸気静止)のときのｄｑ/ｄｔを(ｄｑ/ｄｔ)₀と書きます。

ｕ＝0 の熱伝導方程式:∂Ｔ/∂ｔ＝κ∇²Ｔのｔ→∞での定常状態(∂Ｔ/∂ｔ＝0)の∇²Ｔ＝ｒ^-1ｄ²(ｒＴ)/ｄｒ²＝0 の解:Ｔ＝Ｔ_∞＋ａ(Ｔ_∞－Ｔ_a)/ｒから,(ｄｑ/ｄｔ)₀＝4πａｋ_a(Ｔ_∞－Ｔ_a)(28)です。

ａが小さい小水滴のときには,先に拡散係数Ｄ_hの代わりにＤ_h'を用いるよう補正したのと同じく,熱伝導度ｋ_aの代わりにｋ_a'≡ｋ_a/[ａ/(ａ＋Δ_T)＋ｋ_a{2πＭ_w/(ＲＴ_a)}^1/2/(ａα_T＋ρｃ_p)](29)を用います。

なお,Δ_T,α_Tは拡散での反跳長さΔ_v,凝結定数α_cに相当する熱伝導(熱拡散)での量です。

また,水滴に対して水蒸気が運動するときには,Nusselt数と呼ばれる数Ｎ_uをＮ_u≡2ｄｑ/ｄｔ/(ｄｑ/ｄｔ)₀で定義します。

(ｄｑ/ｄｔ)₀＝4πａｋ_a(Ｔ_∞－Ｔ_a) (28)からｄｑ/ｄｔ＝Ｎ_uｋ_a 2πａ(Ｔ_∞－Ｔ_a)(30)と書けます。ｕ＝0 のときにはＮ_u＝2です。

今日はここまでですが,まだ続きます。

なお,もちろんオリジナルではなく当時,参考文献を見たはずなのですが,それをノートに書いてないので今は不明です。

確か1983年当時に勤めていた会社の蔵書だったと記憶してます。

PS:新しい試みとしてワードのルビ(ふりがな)を記号としてココログにコピーしてみたのですが取り合えずのところは失敗でした。おかげで文章を校正する時間を浪費しました。

貧乏ひまなしで師走は日曜日でも時間がありません。バスに遅れるので続きはあとで。

PS:民主党だけドタバタしてるように見えるのは与党だから仕方ないとはいえ,こんなとき野党は何のためにあるのか？責任追及だけやるしかないのか？そもそも代議士とは何をするのが仕事なのか？

政界再編はいずれ必要なんだろうけど,今やられてもそれが終わる頃には私自身の寿命などは終わってるだろうなあ。。。

イヤ,またまたマスコミや評論家の多くと同じく,自分のこと抜きで他人を揶揄しちまったか。

なんとか法とかの法律に賛成反対とか,国や公共団体の迷惑行為,大きなお世話はやめて欲しいけど,基本的に政治のお世話には期待してないアナーキーを自認している自分。。

自分の世話さえままならず,社会的にも無力で有象無象の自分がいまさら何の能書きをたれてるのか。お笑いだ。。。

PS2：機械文明や都会のコンビニエンスに依存し切っている身はいまさらどうしようもないという気分もあるけれど,もしもここ10年くらいの間,PCやネットがなかったら,既に廃人になってたかも。。　

今のような状態でも,自分が生き生きと生きてるように感じることが多々あるのは悲しいことかも知れないけど,血の通ってないPCやネット依存症に身を委ねて合理化しているおかげでしょうネ。

プルーム上昇のモデル式(3)

　プルーム上昇式の残りです。このシリーズはこれで終わりです。　

§5.逆転層貫通の条件

プルームが逆転層(inversion layer)を貫通する条件を考察します。

※(補足)：逆転層とは通常大気と逆に高度ｚが上昇するにつれて気温が増加していくような層のことです。

　断熱冷却されながら上昇していく周辺大気より暖かいプルームが上昇しても,こうした層があると自身よりも暖かい周辺大気に遭遇してプルームが浮力を失なうために,反射されて層を貫通できない可能性があります。

こうした事情は,2006年９/25の記事「ベナール対流の安定性とレイリー数」の序論では次のように書いています。

層の底ｚ＝0 では温度がＴ＝Ｔ₀ 天井でＴ＝Ｔ₁ として定常熱伝導の方程式∇²Ｔ＝ 0 の解はＴ(ｚ)＝Ｔ₀－βｚ＝Ｔ₀＋{(Ｔ₁－Ｔ₀)/Ｈ}ｚ です。

ここでβ＞0 ,つまりＴ₀＞Ｔ₁としないと地球上の重力場の中では対流が起きません。そればかりでなく対流が起きるためにはβの値が断熱温度減率(湿度によって変化)よりも大きくないと対流は起きません。(例えばランダウ「流体力学１」を参照)

(補足終わり)※

さて,本題に戻ります。

　ｚ＝ｚ_iに微小厚さΔｚ_iの鋭い逆転層があって,大気はこの高さｚ＝ｚ_i(と厚さΔｚ_iの極く薄い層)を除けば中立,つまり層の上下では共にｓ＝(ｇ/Ｔ_a)(∂θ_a/∂ｚ)＝0 であるとします。

そして逆転層より下:ｚ＜ｚ_iの中立状態の温位はθ_a＝θであり,逆転層より上のｚ＞ｚ_iの同じく中立状態の温位はθ_a＝θ＋Δθ_iであるとします。

(ⅰ)鉛直浮力プルームの場合

ｚ＝ｚ_iの層を通過できたと仮定したときにθ_p≧θ_a,つまりθ_p≧θ＋Δθ_iが満足される,すなわちθ'≡θ_p－θと置いたときにθ'≧Δθ_iが満足されるならプルームは逆転層を貫通できるはずです。

何故なら,プルームが逆転層を貫通する仮想変位を行ったとき,受ける浮力は上向き,つまりｇ(θ'－Δθ_i)/Ｔ_a≧0 なので,それは下には押し戻されないからです。

(ただし,θ'≧Δθ_iは十分条件であって必要条件ではありません。つまりθ'＜Δθ_iで浮力が負の向きであっても上向きの運動量が十分であれば貫通します。)

　微小厚さ:Δｚ_iの逆転層に対して,逆転層の強さ:ｂ_iを,(32)ｂ_i≡ｇΔθ_i/θ_a＝∫_Δzi(ｇ/θ_a)(∂θ_a/∂ｚ)ｄｚ＝∫_Δziｓｄｚによって定義します。

また,改めてθ'≡θ_p－θに対する浮力流束:Ｆ_zと容積流束:Ｖを定義し直しておきます。Ｆ_z≡∫_P(ｇρ_pθ'ｗ'/Ｔ_a)ｄｘｄｙ/(πρ_a),Ｖ＝∫_Pρ_pｗ'ｄｘｄｙ/(πρ_a)です。

　こう定義すると,プルーム断面上でのθ'の平均値は,＜θ'＞＝(Ｔ_a/ｇ)(Ｆ_z/Ｖ)で与えられます。

　すぐ上で述べた逆転層貫通の(十分)条件θ'≧Δθ_iはプルーム断面について,＜θ'＞≧Δθ_i,つまり(Ｔ_a/ｇ)(Ｆ_z/Ｖ)≧Δθ_i⇔ (Ｆ_z/Ｖ)≧(ｇΔθ_i/Ｔ_a)と解釈されます。

一方,ｂ_i＝ｇΔθ_i/θ_aですから,θ_a～Ｔ_aより結局逆転層の貫通条件は(Ｆ_z/Ｖ)≧ｂ_iと表わせることがわかります。

ところで,中立大気ではｓ＝(ｇ/Ｔ_a)(∂θ_a/∂ｚ)＝0 ですから方程式の１つ:(22)ｄＦ_z/ｄｔ＝－ｓ(ｗＶ)からＦ_z＝Ｆ (一定)です。

そこで,ｚ＝ｚ_i～ｚ_i＋Δｚ_iでもＦ_zi＝Ｆ(一定)と考えられます。

すると,(23)ｄ(ｗＶ)/ｄｔ＝Ｆ_zは,ｄ(ｗＶ)/ｄｔ＝Ｆ(一定)となりますから,これもすぐに解けてｗＶ＝Ｆｔ＋Ｆ_mを得ます。

それ故,最後の(24)ｄＶ/ｄｔ＝2α(ｗＶ)^3/2/Ｖはｄ(Ｖ²/2)/ｄｔ＝2α(Ｆｔ＋Ｆ_m)^3/2となります。

ｔ＝0でＶ＝0 の点源なら,解はＶ²/2＝(4α/5)Ｆ^-1{(Ｆｔ＋Ｆ_m)^5/2－Ｆ_m^5/2},すなわちＶ＝(8α/5)^1/2Ｆ^-1/2{(Ｆｔ＋Ｆ_m)^5/2－Ｆ_m^5/2}^1/2と書けます。

　今の場合は浮力プルーム(Ｆ_m＝0)を想定しているので,ｗＶ＝Ｆｔ＋Ｆ_m＝Ｆｔであり,Ｖ＝(8α/5)^1/2Ｆ^3/4ｔ^5/4です。

　したがってｗ＝Ｆｔ/Ｖ＝(8α/5)^-1/2Ｆ^1/4ｔ^-1/4を得ます。

　これが,ｔ＝ｔ_iでｚ＝ｚ_iに到達するとすれば,ｚ_i＝∫₀^tiｗｄｔ＝(8α/5)^-1/2(4/3)Ｆ^1/4ｔ_i^3/4ですから,ｔ_i＝(8α/5)^2/3(3/4)^4/3Ｆ^-1/3ｚ_i^4/3です。

　そこで,Ｖ_i＝(8α/5)^1/2Ｆ^-1/4ｔ_i^5/4＝(8α/5)^4/3(3/4)^5/3Ｆ^1/3ｚ_i^5/3を得ますから,結局Ｖ_i＝(4・3⁵α⁴/5⁴)^1/3Ｆ^1/3ｚ_i^5/3です。

　そして,実験によればＶ_i～ 0.0717Ｆ^1/3ｚ_i^5/3 ⇔ α～ 0.125です。

　以上から,浮力プルームの貫通条件:(Ｆ_zi/Ｖ_i)≧ｂ_iは,{5⁴/(4・3⁵α⁴)}^1/3Ｆ^2/3ｚ_i^-5/3≧ｂ_iと書けることがわかります。

したがって,無風時(calm)に中立状態で浮力プルームが逆転層(ｚ＝ｚ_i)を貫通する条件として,(32)ｚ_i≦(4・3⁵α⁴/5⁴)^-1/5Ｆ^2/5ｂ_i^-3/5が得られました。

(ⅰ)折れ曲がりプルームの場合

　この場合も,浮力プルームについては(Ｆ_zi/Ｖ_i)≧ｂ_iを逆転層貫通の条件であるとします。

　周辺大気は中立(ｓ＝0)なのでＦ_zi＝Ｆ(一定)です。

そして,前記事で書いたように,(21)ｄＶ/ｄｚ＝2β(ｕＶ)^1/2のｚ＝0(ｔ＝0)でＶ＝0 の解は,Ｖ＝ｕ(βｚ)²ですからＶ_i＝ｕ(βｚ_i)²より,Ｆ_zi/Ｖ_i＝{Ｆ/(ｕβ²)}ｚ_i^-2です。

　そこで,折れ曲がり浮力プルームの逆転層貫通条件は,Ｆ/(β²ｕ)}ｚ_i^-2≧ｂ_iと表現できます。

したがって,有風時に中立状態で浮力プルームが逆転層(ｚ＝ｚ_i)を貫通する条件として(33)Ｆ≧β²ｕｂ_iｚ_i²～ 0.25ｕｂ_iｚ_i², or (33)'ｚ_i≦β^-1{Ｆ/(ｕｂ_i)}^1/2～ 2.0{Ｆ/(ｕｂ_i)}^1/2が得られました。

(ⅲ)鉛直プルームで浮力がゼロのジェット・プルーム,および浮力流束は不足だが十分な運動量の鉛直流束を持つ場合

　初期浮力がゼロのジェット(Ｆ＝0)の場合を考えると,ｓ＝0 なら(22)ｄＦ_z/ｄｔ＝－ｓ(ｗＶ)＝0 よりＦ_z＝Ｆ＝0 です。

　そこで,ｚ＜ｚ_iではＦ_z＝0 (一定)ですがｚ＝ｚ_iで浮力流束:－ｂ_iＶ_i＝－(∫_Δziｓｄｚ)Ｖ_iを獲得して,ｚ＞ｚ_iではＦ_z＝－ｂ_iＶ_i(一定)と考えることができます。

　逆転層:ｚ＝ｚ_iより上ではプルーム内部の温位θ_pが外気のそれθ＋Δθ_iよりも小さく,エントレインメントが無視できてＶ＝Ｖ_i(一定)と仮定すれば次のように書けるはずです。

ｚ＜ｚ_i(ｔ＜ｔ_i)ではＦ_z＝0 (一定),ｗＶ＝Ｆ_m(一定),ｄＶ/ｄｔ＝＝2α(ｗＶ)^3/2/Ｖであり,一方,ｚ＞ｚ_i(ｔ＞ｔ_i)ではＦ_z＝－ｂ_iＶ_i(一定),ｗＶ＝－ｂ_iＶ_i(ｔ－ｔ_i)＋Ｆ_m,Ｖ＝Ｖ_i(一定)です。

そこで,ｔ＞ｔ_iではｗ＝－ｂ_i(ｔ－ｔ_i)＋Ｆ_m/Ｖ_iより,ｗ＝0 となるのはｔ_max＝ｔ_i＋Ｆ_m/(ｂ_iＶ_i)のときです。

このときの高さ:ｚ＝ｚ_maxは,ｚ_max＝ｚ_i＋∫_ti^ｔmaxｗｄｔ＝ｚ_i－(ｂ_i/2)(ｔ_max－ｔ_i)²＋(Ｆ_m/Ｖ)(ｔ_max－ｔ_i)で与えられます。よってｚ_max－ｚ_i＝Ｆ_m²/(2ｂ_iＶ_i²)です。

さらに,Ｖ_iをＦ_m,ｚ_iで表わすためｔ＜ｔ_iでのＶを求めます。

従う方程式はｄＶ/ｄｔ＝2α(ｗＶ)^3/2/Ｖ or ｄ(Ｖ²/2)/ｄｔ＝2α(ｗＶ)^3/2＝2αＦ_m^3/2です。

ｔ＝0 でＶ＝Ｖ₀≡Ｆ_m/ｗ₀の解は,Ｖ²/2＝2αＦ_m^3/2ｔ＋(Ｆ_m/ｗ₀)², or Ｖ＝{αＦ_m^3/2ｔ＋(Ｆ_m/ｗ₀)²}^1/2で与えられます。

そして,ｗ＝Ｆ_m/Ｖ＝よりｚ_i＝∫₀^ti(Ｆ_m/Ｖ)ｄｔ＝2Ｆ_m(Ｖ_i－Ｖ₀)/(4αＦ_m^3/2)ですからＶ_i＝2αＦ_m^1/2ｚ_i＋Ｖ₀です。

点源ならＶ₀＝0 なので,Ｖ_i＝2αＦ_m^1/2ｚ_iより,ｚ_max－ｚ_i＝Ｆ_m²/(2ｂ_iＶ_i²)＝(1/2)(Ｆ_m/ｂ_i)(2α)^-2ｚ_i^-2です。

よって,(34)ｚ_max/ｚ_i＝1＋(1/8)α^-2(Ｆ_m/ｂ_i)ｚ_i^-3を得ます。

このケースでは,ｚ_max≧3ｚ_iを貫通条件として採用します。そして,そして,α～ 0.125としてその条件を求めることにします。

無風時(calm)に中立状態でジェット・プルームが逆転層(ｚ＝ｚ_i)を貫通する条件として,(35)Ｆ_m≧16α²ｂ_iｚ_i^-3 ～ 0.25ｂ_iｚ_i^-3, or (35)'ｚ_i≦(4α)^2/3(Ｆ_m/ｂ_i)^1/3～ 1.59(Ｆ_m/ｂ_i)^1/3が得られました。

もしも,わずかの浮力Ｆを持つプルームならｚ＜ｚ_i(ｔ＜ｔ_i)ではＦ_z＝Ｆ(一定),ｗＶ＝Ｆｔ＋Ｆ_m,ｄ(Ｖ²/2)/ｄｔ＝2α(Ｆｔ＋Ｆ_m)^3/2であり,一方,ｚ＞ｚ_i(ｔ＞ｔ_i)ではＦ_z＝Ｆ－ｂ_iＶ_i(一定),ｗＶ＝－ｂ_iＶ_i(ｔ－ｔ_i)＋Ｆｔ＋Ｆ_m,Ｖ＝Ｖ_i(一定)です。

ｔ＞ｔ_iでは,ｗ＝－ｂ_i(ｔ－ｔ_i)＋Ｆｔ/Ｖ_i＋Ｆ_m/Ｖ_iより,ｗ＝0 となるのはｔ＝ｔ_max＝(ｂ_iＶ_iｔ_i＋Ｆ_m)/(ｂ_iＶ_i－Ｆ)のときです。

このときの高さ:ｚ＝ｚ_maxは,ｚ_max－ｚ_i＝∫_ti^ｔmaxｗｄｔ＝(Ｆｔ_i＋Ｆ_m)²/{2Ｖ_i(ｂ_iＶ_i－Ｆ)}で与えられます。

ここで,Ｖ_i＝(8α/5)^1/2Ｆ^-1/2{(Ｆｔ_i＋Ｆ_m)^5/2－Ｆ_m^5/2}^1/2,かつ(Ｆｔ_i＋Ｆ_m)²＝{5Ｖ_i²Ｆ/(8α)＋Ｆ_m^5/2}^4/5なのでｚ_max－ｚ_i＝{5Ｖ_i²Ｆ/(8α)＋Ｆ_m^5/2}^4/3/{2Ｖ_i(ｂ_iＶ_i－Ｆ)}なる式を得ます。

　以下,この場合もｚ_max≧3ｚ_iを貫通条件として中立状態大気中の逆転層(ｚ＝ｚ_i)を貫通する条件式が得られます。

　何故か,ここまででノートが終わってるので,この項はこれでおしまいにします。

参考文献：G.A.Briggs “Plume Rise",U.S.Atomic Energy Commision Division of Technical Information(1969)

PS1:土曜日４日午後に,４～６日に文京区春日の文京シビックセンターで開催されていた「障害者ふれあいの集い」に行ってきました。

　　

　　

　11/30に出席した手話講習会での情報によると,豊島区や北区でも同時に開催されているらしいです。

　実は私もつたない作品＝"牛乳などの紙パックから再生紙を作り折り紙などを切り貼りしてパウチしたコースター"を何枚か出品していて,１枚30円也で販売されてたりしていました。

　　　

　　

　　

　　

PS2：鈴木宗男さん。がんばれ。

　一緒にしたら迷惑かもしれないけど高知のスクールバス運転手の件といい,理不尽だと思いますね。

　ちまたでウワサの市川海老蔵氏とそのケンカ相手も世論に負けないでくださいネ。。

　有名人といえども酒癖はしょうがない面ありますね。

　不特定の乗客を乗せるタクシー運転手と同じく(尤もこちらは遊びじゃなく仕事ですが),常識があろうとなかろうと外での個人的飲酒はキケンと背中合わせなのは.大人なんだから覚悟の上でしょう。

　男でも女でも,普段は物静かだったりリッパに見える方でも,酒の席では基本的に無礼講です。

　酔っ払うと自慢話やホラ話のテンコ盛り,酒乱,泣き上戸,怒り上戸,裸になるなどなど,性格が変わる(実は本性？)ような方々,いっぱい見てきています。本当は羨ましい限りです。

　本気で酒でストレス解消しようとに飲みにきているのなら,その解消法はヒトさまざまです。まあ,そうした愚行のお相手も込みで飲み屋の代金を払ってるのだと思っています。

　有名人とはいえ,かつての草薙つよポンのケースと同じく,普段の姿と重ねて批判されるのはカワイソーですね。

　私だって他人の一方的なストレス解消を迷惑だと感じて,ケンカは弱いけれど思うだけなら「こいつ,どついたろか？」と考えることもしばしばですが。。。

プルーム上昇のモデル式(2)

　今日もプルーム上昇式の続きです。

§2.理論の単純化,エントレインメントの仮定

　これまでの議論では,未知数Ｆ_z,ｗ,Ｖに対して２つの方程式しか与えられておらず,問題が閉じていません。

　これは,プルームの明確な定義が与えられていないため,プルームの大きさを与える式が無いことに起因すると考えられます。

　そこで,連続の方程式∇(ρ_pｖ_p)＝0 を用いてプルームの大きさを規定する,ある方程式を与えます。

　プルーム半径をｒとして,改めてプルームの平均進行速度の大きさをｖ_sとします。さらに,エントレインメント(entrainment)速度,つまり,プルーム境界から乱流により流入する外気の平均速さをｖ_eとします。

　そして,プルームの行路長さ:ｓに対してｄｓ≡ｖ_sｄｔ,またはｄｔ≡ｄｓ/ｖ_sによってｄｔを定義しておきます。

　すると,質量の保存式から明らかにΔ(ρ_pπｒ²ｖ_s)＝ρ_a2πｒｖ_eΔｓ,または(ｄ/ｄｓ)(ρ_pｒ²ｖ_s/ρ_a)＝2ｒｖ_eです。こうした関係があることをエントレインメントの仮定といいます。

　特に,鉛直プルーム:ｓ＝ｚ,ｄｚ＝ｗｄｔの場合にはＶ＝ρ_pπｒ²ｗ/(πρ_a)＝ρ_pｒ²ｗ/ρ_aなのでエントレインメントの仮定はｄＶ/ｄｚ＝2ｒｖ_eを意味します。

　また,(ｄ/ｄｓ)(ρ_pｒ²ｖ_s/ρ_a)＝2ｒｖ_eよりρ_aはｓに依らないとして(ｒ²/ρ_a){ｄ(ρ_pｖ_s)/ｄｓ}＋(2ｒρ_pｖ_s/ρ_a)(ｄｒ/ｄｓ)＝2ｒｖ_eです。

さらに,ｄｔ＝ｄｓ/ｖ_sよりｄ/ｄｓ＝(1/ｖ_s)(ｄ/ｄｔ)を用いればｄｒ/ｄｔ＝ρ_aｖ_e/ρ_p－{ｒ/(2ρ_p)}{ｄ(ρ_pｖ_s)/ｄｓ}を得ます。

ここでTaylorによる自然な仮定:ｖ_e∝ｗを採用し,鉛直プルームに対するエントレインメント係数αを導入してｖ_e＝αｗとします。

また,再びＶ＝ρ_pπｒ²ｗ/(πρ_a)＝ρ_pｒ²ｗ/ρ_aよりｒ²＝ρ_aＶ/(ρ_pｗ)～Ｖ/ｗと近似できるという考察から,逆にプルーム半径:ｒをｒ≡(Ｖ/ｗ)^1/2によって定義します。

この定義によれば,Ｖ＝ｒ²ｗですが,Δ(ρ_pπｒ²ｖ_s)＝ρ_a2πｒｖ_eΔｓ,または(ｄ/ｄｓ)(ρ_pｒ²ｖ_s/ρ_a)＝2ｒｖ_eの関係は変わらず,ρ_p～ρ_aの近似で考えると依然としてエントレインメントの仮定はｄＶ/ｄｚ＝2ｒｖ_eです。

元々,これも仮定ですから以下ではｄＶ/ｄｚ＝2ｒｖ_eをエントレインメントの仮定として採用します。この右辺にｖ_e＝αｗ,ｒ＝(Ｖ/ｗ)^1/2を代入するとｄＶ/ｄｚ＝2α(ｗＶ)^1/2を得ます。

以上から,(抵抗力)を無視すると次のような簡単な基本方程式系が得られました。

(4)ｄＦ_z/ｄｚ＝－ｓＶ,(12)ｄ(ｗＶ)/ｄｚ＝Ｆ_z/ｗ,

および,(13)ｄＶ/ｄｚ＝2α(ｗＶ)^1/2　です。

一応,未知量Ｖ,Ｆ_z,ｗの意味を再び述べておきます。

Ｖは鉛直容積流束,Ｆ_zは鉛直方向の浮力(ｇθ'/Ｔ)の流束,ｗはプルーム面の鉛直速度成分です。

プルームの行路長さと混同するおそれのある同じ変数ｓを用いていますが,(4)式のｓは,もちろんｓ≡(ｇ/Ｔ_a)(∂θ_a/∂ｚ)で定義される安定度パラメータです。

そして,Morton,Taylor,Turner(1956)によれば,エントレインメント係数αの具体的な値の候補はα＝0.093です。一方,Briggs(1968)によればα＝0.075です。また,Ricou,Spalding(1961)によればジェットに対しα＝0.080です。

　§3.折れ曲がりプルーム(bent-over plume)

　進行方向が完全に鉛直ｚ方向でｓ＝(一定)のプルーム断面が水平面内にある鉛直プルームとは異なる状態のプルームを折れ曲がりプルームといいます。

折れ曲がりプルームの場合,積分面をプルームと交わる鉛直面として,この面上ではｓ＝(ｇ/Ｔ)(∂θ_a/∂ｚ),およびｖ_a＝ｕｉは一定であると仮定します。

前と同様,∫_∞∇(ρ_pθ'ｖ_p)ｄｙｄｚ＝－(∂θ_a/∂ｚ)∫_∞ρ_pｗ'ｄｙｄｚですがプルームの外ではθ'＝0,ｗ'＝0 とすれば,これは(ｄ/ｄｘ)∫_Pρ_pθ'(ｕ＋ｕ')ｄｙｄｚ＝－(∂θ_a/∂ｚ)∫_Pρ_pｗ'ｄｙｄｚを意味します。

そこで,浮力流束を(5)Ｆ_z≡∫_P(ｇρ_pθ'ｗ'/Ｔ_a)ｄｘｄｙ/(πρ_a)の代わりに(14)Ｆ_z≡∫_P{ｇρ_pθ'(ｕ＋ｕ')/Ｔ_a}ｄｙｄｚ/(πρ_a)で定義すれば,(15)ｄＦ_z/ｄｘ＝－ｓＭと表わすことができます。

ただし,再びｓはｓ＝(ｇ/Ｔ_a)(∂θ_a/∂ｚ)でありＭはＭ≡∫_Pρ_pｗ'ｄｙｄｚ/(πρ_a)で定義される量です。

この場合も排出源がhot-sourceなら,浮力流束Ｆ_zの初期値Ｆは{(単位時間にプルーム積分面を通過する総排出熱量)×(ｇ/Ｔ_a)/ｃ_p}/(πρ_a)と考えられるので,やはりＦ＝ｇＱ_H/(πρ_aｃ_pＴ_a)で与えられます。

　そして,運動量に対する方程式は,(8)ｄ(ｖＶ)/ｄｚ＝(Ｆ_z/ｗ)ｋ＋ｖ_a(ｄＶ/ｄｚ)＋(抵抗力)の証明の際に遭遇した式:(ｄ/ｄｎ)(∫_Pρ_pｖ_pnｖ_piｄＳ)＝{1＋(ｄｒ/ｄｎ)²}^1/2(∫_Cρ_pｖ_pｖ_piｄｓ)＋δ_i3∫_P(ｇρ_pθ'/Ｔ_a)ｄＳにおいて進行方向ｎをｉに置くと得られます。

　ｄ/ｄｎをｄ/ｄｘに変え,ｖ_p＝(ｕ＋ｕ')ｉ＋ｖ'ｊ＋ｗ'ｋ,ｄｓ＝{1＋(ｄｒ/ｄｘ)²}^-1/2(－ｄｚｊ＋ｄｙｋ)＋{1＋(ｄｒ/ｄｘ)²}^-1/2(ｄｒ/ｄｘ)ｄｌｉとして,必要なｉ＝3の場合のみを書き下します。

　それは,(ｄ/ｄｘ){∫_Pρ_p(ｕ＋ｕ')ｗ'ｄＳ}＝∫_P(ｇρ_pθ'/Ｔ_a)ｄＳ＋{1＋(ｄｒ/ｄｘ)²}^1/2(∫_Cρ_pｖ_pｗ'ｄｓ)＝∫_P(ｇρ_pθ'/Ｔ_a)ｄＳ－∫_Cρ_pｖ'ｗ'ｄｚ＋∫_Cρ_pｗ'²ｄｙ＋∫_Cρ_p(ｕ＋ｕ')ｗ'(ｄｒ/ｄｘ)ｄｌです。

　ｕ'＜＜ｕなのでｕ＋ｕ'～ ｕと近似して有効運動量流束をｕＭ_eff;Ｍ_eff≡{∫_Pρ_pｗ'ｄＳ－∫ｄｘ∫_Cρ_pｗ'(ｄｒ/ｄｘ)ｄｌ}/(πρ_a)と定義しｕ',ｖ',ｗ'の２次以上の項を無視します。

　すると,Ｆ_z＝{ｕ∫_P(ｇρ_pθ'/Ｔ_a)ｄＳ}/(πρ_a)より,結局運動方程式として(16)ｕ(ｄＭ_eff/ｄｘ)＝Ｆ_z/ｕが得られます。

　鉛直運動量の水平方向流束ｕＭ_effの初期値はプルーム面を通過する鉛直運動量ρ_pｗπｒ²/(πρ_a)＝ρ_pｗｒ²/ρ_aの総流束(鉛直方向流束)ですから,その初期値はＦ_m＝(ρ₀ｗ₀ｒ₀²/ρ_a)×ｗ₀＝ρ₀ｗ₀²ｒ₀²/ρ_a＝ρ₀ｗ₀Ｖ₀/ρ_aです。

そして,Briggs(1972),Richards(1963)によれば有効運動量流束ｕＭ_effの運動量流束ｕＭ;Ｍ＝∫_Pρ_pｗ'ｄｙｄｚ/(πρ_a)に対する比率の評価値はＭ_eff/Ｍ ～ 2.3です。

ｚをプルーム軸心の高さ,ｄｔ≡ｄｘ/ｕとしてｗ≡ｄｚ/ｄｔ＝ｕｄｚ/ｄｘと定義します。ｗはプルーム軸心の上昇速度を与え,解(ｘ,ｚ),or ｚ＝ｚ(ｘ)はプルームの軌跡を示します。

ここで,折れ曲がりプルームにおけるエントレインメント係数をβとします。すなわちｖ_e＝βｗです。有効プルーム半径をｒとすると,Ｍ_eff＝πｒ²ρ_pｗ/(πρ_a)＝ｗρ_pｒ²/ρ_aです。

折れ曲がりプルームの場合のエントレインメントの仮定はｄ(ρ_pπｒ²ｕ)/ｄｘ＝2πｒρ_aｖ_e,つまりｄ(ρ_pｒ²ｕ/ρ_a)/ｄｘ＝2ｒｖ_eと表わせます。

そこで,ｕｄ(Ｍ_eff/ｗ)/ｄｘ＝2ｒｖ_e＝2ｒβｗです。

ｒ≡(Ｍ_eff/ｗ)^1/2と定義すれば,(17)ｕｄ(Ｍ_eff/ｗ)/ｄｘ＝2β(ｗＭ_eff)^1/2を得ます。

さらにｗＶ≡ｕＭ_effによってＶを定義すれば,Ｖは水平方向の容積流束でありｒ＝(Ｖ/ｕ)^1/2です。これとｕ/ｄｘ＝ｗ/ｄｚなる関係を用いて基本方程式系:(15),(16),(17)を書き直せば次のようになります。

すなわち,(18)ｄＦ_z/ｄｚ＝－ｓ'Ｖ,(19)ｓ'＝ｓ(Ｍ/Ｍ_eff),(20)ｄ(ｗＶ)/ｄｚ＝Ｆ_z/ｗ,(21)ｄＶ/ｄｚ＝2β(ｕＶ)^1/2です。

式(21):ｄＶ/ｄｚ＝2β(ｕＶ)^1/2は,ｚ＝0 でＶ＝0 という点源(point source)の初期条件に対しては容易に解けて,(Ｖ/ｕ)^1/2＝βｚ,つまりｒ＝βｚです。

Richards(1963),TVA-Report(1967)によれば折れ曲がりプルームのエントレインメント係数はβ～ 0.5です。

§4.プルーム上昇高度の算出

　具体的にプルーム上昇式を求めるに際して,さらに次のような２つの仮定を与えます。

(仮定Ⅰ):鉛直プルームを満足する状態を無風状態(静穏:calm)に対応させ,折れ曲がりプルームを有風時の状態に対応させる。

(仮定Ⅱ):初期時には鉛直運動量か浮力流束のいずれか一方が支配的である。前者はジェットと呼んでＦ＝0 とし,後者を浮力プルームと呼んでＦ_m＝0 とする。

　さて,まずｄｘ＝ｕｄｔ,ｄｚ＝ｗｄｔにより,全ての方程式を独立変数ｔの式に変換します。

鉛直プルームについては,(4),(12),(13)より

(22)ｄＦ_z/ｄｔ＝－ｓ(ｗＶ),(23)ｄ(ｗＶ)/ｄｔ＝Ｆ_z,(24)ｄＶ/ｄｔ＝2α(ｗＶ)^3/2/Ｖです。

また,折れ曲がりプルームについては,(18),(20),(21)より

(25)ｄＦ_z/ｄｔ＝－ｓ'(ｗＶ),(26)ｄ(ｗＶ)/ｄｔ＝Ｆ_z,(27)ｄＶ/ｄｔ＝2βｗＶ(ｕ/Ｖ)^1/2です。

　以下,順に安定大気:ｓ＞0,ｓ'＞0 (ｓ,ｓ'は一定値を仮定)の中でのプルーム上昇高さを求めます。

(ⅰ)鉛直プルーム(無風:calm)

　このケースの基本方程式系は(22)ｄＦ_z/ｄｔ＝－ｓ(ｗＶ),(23)ｄ(ｗＶ)/ｄｔ＝Ｆ_z,(24)ｄＶ/ｄｔ＝2α(ｗＶ)^3/2/Ｖです。

(22),(23)よりｄ²(ｗＶ)/ｄｔ²＝－ｓ(ｗＶ)を得ます。これはよく知られている線形振動の方程式です。

そこで,ω≡ｓ^1/2としてこれをｄ²(ｗＶ)/ｄｔ²＝－ω²(ｗＶ)と書き,初期条件:ｔ＝0 でｗＶ＝Ｆ_m,ｄ(ｗＶ)/ｄｔ＝Ｆ_z＝Ｆの解を求めると,ｗＶ＝Ｆ_mcosωｔ＋(Ｆ/ω)sinωｔを得ます。

鉛直運動量流束ｗＶがゼロとなる最初のｔ＝ｔ₁がプルーム上昇高さを与えると考えられます。

上昇高さはΔｈ＝∫₀^t1ｗｄｔですから,仮定Ⅱの浮力プルーム(Ｆ_m＝0)ではｔ₁＝π/ω,よりΔｈ＝ω^-1∫₀^πｗｄλです。

同様にジェットプルーム(Ｆ＝0)ではｔ₁＝π/(2ω)よりΔｈ＝ω^-1∫₀^π/2ｗｄλ(λ＝ωｔ)です。

①     浮力プルーム(Ｆ_m＝0)の場合

　ｄ(Ｖ²/2)/ｄｔ＝Ｖ(ｄＶ/ｄｔ)＝2α(ｗＶ)^3/2の右辺にｗＶ＝(Ｆ/ω)sinωｔを代入してｄ(Ｖ²/2)/ｄｔ＝2α(Ｆ/ω)^3/2sin^3/2ωｔです。

故に,Ｖ²/2＝2α(Ｆ/ω)^3/2ω^-1{∫^ωtsin^3/2λｄλ＋const.},すなわちＶ＝2α^1/2Ｆ^3/4ω^-5/4{∫^ωtsin^3/2νｄν＋const.}^1/2です。

　ｗ＝ｗＶ/Ｖ＝(Ｆ/ω)sinωｔ/Ｖより,Δｈ＝ω^-1∫₀^πｗｄλ＝(1/2)α^-1/2Ｆ^1/4ω^-3/4[∫₀^πsinλｄλ{∫₀^λsin^3/2νｄν＋const.}^-1/2]を得ます。

したがって,ω＝ｓ^1/2より定数ｃをｃ≡(1/2)α^-1/2[∫₀^πsinλｄλ{∫₀^λsin^3/2νｄν＋const.}^-1/2]で定義すれば,Δｈ＝ｃＦ^1/4ｓ^-3/8なる関係式が成立することになります。

　ｔ＝0 のとき,Ｖ＝0 の点源(const.＝0)ではMorton,Taylor,Turner(1956)の数値計算があり,その示唆するところによればα＝0.132としてｃ～ 5.0ですから,浮力プルーム上昇高さの無風状態での評価式として(28)Δｈ＝5.0Ｆ^1/4ｓ^-3/8を得ます。

②     ジェットプルーム(Ｆ＝0)

　ｗＶ＝Ｆ_mcosωｔよりｄ(Ｖ²/2)/ｄｔ＝Ｖ(ｄＶ/ｄｔ)＝2α(ｗＶ)^3/2の右辺にｗＶ＝(Ｆ/ω)sinωｔを代入すれば,ｄ(Ｖ²/2)/ｄｔ＝2αＦ_m^3/2cos^3/2ωｔです。

故に,Ｖ²/2＝2αＦ_m^3/2ω^-1{∫^ωtcos^3/2λｄλ＋const.},すなわちＶ＝2α^1/2Ｆ_m^3/4ω^-1/2{∫^ωtcos^3/2νｄν＋const.}^1/2です。

　そして,ｗ＝ｗＶ/Ｖ＝Ｆ_mcosωｔ/Ｖより,Δｈ＝ω^-1∫₀^π/2ｗｄλ＝(1/2)α^-1/2Ｆ_m^1/4ω^-1・2[∫₀^π/2cosλｄλ{∫₀^λcos^3/2νｄν＋const.}^-1/2]です。

したがって,ω＝ｓ^1/2より定数ｃ₀をｃ₀≡(1/2)α^-1/2[∫₀^π/2cosλｄλ{∫₀^λcos^3/2νｄν＋const.}^-1/2]と定義すれば,Δｈ＝ｃ₀Ｆ_m^1/4ｓ^-1/4なる関係の成立することがわかります。

　ｔ＝0 のとき,Ｖ＝0 の点源(const.＝0)では,μ＝π/2－νと変数変換すれば,cosν＝sinμより∫₀^λcos^3/2νｄν＝－∫_π/2^π/2－λsin^3/2μｄμ,ｃ₀＝(1/2)α^-1/2[∫₀^π/2cosλｄλ{∫_π/2-λ^π/2sin^3/2μｄμ}^-1/2]＝(1/2)α^-1/2[∫₀^π/2sinσｄσ{∫₀^σsin^3/2μｄμ}^-1/2]となります。

　点源のとき,ｃ＝(1/2)α^-1/2[∫₀^πsinλｄλ{∫₀^λsin^3/2νｄν＋const.}^-1/2]に対してｃ₀＝(1/2)α^-1/2[∫₀^π/2sinλｄλ{∫₀^λsin^3/2νｄν＋const.}^-1/2]であり,ほとんどｃと同じです。

ｃを求めたのと同じ数値計算によれば,α＝0.132に対してｃ₀～ 4.0ですから,無風状態でのジェットプルームの評価式として(29)Δｈ＝4.0(Ｆ_m/ｓ)^1/4を得ます。

(ⅱ)折れ曲がりプルーム(有風時)

　このケースの基本方程式系は(25)ｄＦ_z/ｄｔ＝－ｓ'(ｗＶ),(26)ｄ(ｗＶ)/ｄｔ＝Ｆ_z,(27)ｄＶ/ｄｔ＝2βｗＶ(ｕ/Ｖ)^1/2です。

(25),(26)よりｄ²(ｗＶ)/ｄｔ²＝－ｓ(ｗＶ)なのでω'＝ｓ'^1/2としてｗＶ＝Ｆ_mcosω'ｔ＋(Ｆ/ω')sinω'ｔです。

これを(27)から得られるｄ(2Ｖ^3/2/3)/ｄｔ＝2βｕ^1/2ｗＶの右辺に代入すると,ｄ(2Ｖ^3/2/3)/ｄｔ＝2βｕ^1/2{Ｆ_mcosω'ｔ＋(Ｆ/ω')sinω'ｔ}です。

　故に,Ｖ＝[3βｕ^1/2{(Ｆ_m/ω')sinω'ｔ＋(Ｆ/ω'²)(1－cosω'ｔ)}]^2/3＝[(3βｕ^1/2/ｓ'){(ω'Ｆ_msinω'ｔ＋Ｆ(1－cosω'ｔ))}^2/3です。

　一方,(21)ｄＶ/ｄｚ＝2β(ｕＶ)^1/2のｚ＝0 (ｔ＝0)でＶ＝0 の解はＶ＝ｕ(βｚ)²ですから,上式のＶと等置するとｚ＝(Ｖ/ｕ)^1/2/β＝{3/(β²ｕｓ')}^1/3{(ω'Ｆ_msinω'ｔ＋Ｆ(1－cosω'ｔ))}^1/3を得ます。

　三角関数の合成ω'Ｆ_msinω'ｔ－Ｆcosω'ｔ＝(ω'²Ｆ_m²＋Ｆ²)^1/2sin(ω'ｔ＋δ)を用いると,ｚの最大値がｚ_max＝{3/(β²ｕｓ')}^1/3{Ｆ＋(ω'²Ｆ_m²＋Ｆ²)^1/2}^1/3で与えられることがわかります。

　そして,(19)ｓ'＝ｓ(Ｍ/Ｍ_eff)によりｓ'∝ ｓなので,有風時の浮力プルーム(Ｆ_m＝0)の場合は,Δｈ＝ｃ₁{Ｆ/(ｕｓ)}^1/3,ジェットプルーム(Ｆ＝0)の場合は,Δｈ＝ｃ₂(Ｆ_m/ｕ)^1/3ｓ^-1/6と書けます。

　ただし,ｃ₁＝[6/{β²(ｓ'/ｓ)}]^1/3,ｃ₂＝[3/{β²(ｓ'/ｓ)}]^1/3ですから,例えばβ～ 0.5,(ｓ'/ｓ)～ 2.3を代入して計算するとｃ₁～ 2.19,ｃ₂～ 1.73です。

一方,Briggsによれば,有風時の浮力プルームの場合の上昇式は,(30)Δｈ＝2.4{Ｆ/(ｕｓ)}^1/3,有風時のジェットプルームの場合の上昇式は(31)Δｈ＝1.5(Ｆ_m/ｕ)^1/3ｓ^-1/6で与えられるらしいです。

今日はここまでです。この項目はまだ少し続きます。

参考文献：G.A.Briggs "Plume Rise",U.S.Atomic Energy Commision Division of Technical Information(1969)

プルーム上昇のモデル式(1)

　唐突ですが,過去の資料を整理していたら1983年(33歳当時)の仕事関係の覚書きノートを見つけたので,ブログに覚書きを残す一環としてこれも書いておこうと思います。

内容は,古典的なBriggsによるプルーム(plume)上昇式の紹介です。(←ノートにはプリューム上昇式と書いてありましたが,今はplumeはプルームと読むべきだと思っているので直しました。)

当時は,環境アセス会社の社員で丁度,冷却塔(cooling-tower)からの温排気のシミュレーションについて調べていた頃と思います。

※(この後,最終的にはANLモデル(Argonne国立研究所モデル)と格闘することになったのですが。。この会社をやめて10年以上後にはSACTIというプログラムとも出会いましたね。illinois大学関連かな？)　

ただし,当時は常識だった気象力学の知識は今ではかなり忘れているので,この記事内容だけで私自身がわかるように補足しています。

Ⅰ.プルーム上昇式の導出

§1.基礎理論(鉛直プルーム)

プルーム内の局所密度ρ_pは不変:∂ρ_p/∂ｔ＝0 とします。

すると質量保存の連続の方程式は(1)∇(ρ_pｖ_p)＝0 となります。ただしｖ_p,およびρ_pはプルームガスの局所速度,および密度です。

また,局所プルーム要素の運動は断熱的とします。

断熱条件はプルームの温位をθ_pとして(2)ｄθ_p/ｄｔ＝0 です。

温度がＴ,圧力がｐの気体に対して温位はθ≡Ｔ(ｐ₀₀/ｐ)^(Ｒ/Ｃp)で定義されます。ただしＲは気体定数:Ｒ≒8.31Ｊ/(molＫ),Ｃ_vはモル定積比熱,Ｃ_pはモル定圧比熱です。

Ｃ_vとＣ_pの間にはMayerの関係式:Ｃ_p＝Ｃ_v＋Ｒが成立します。

ｄθ_p/ｄｔ＝0 は中立大気中の運動では浮力が保存すること,移流以外には熱の放射吸収がないことを意味します。

※(補足):ガス(気体)が希薄で理想気体近似ができるとすると,その状態方程式は気体分子の分子量をｍとしてｐ＝ρＲＴ/ｍです。ｐ,ρ,Ｔは,それぞれ気体の圧力,密度,絶対温度です。

断熱変化の場合には,ｃ₁を定数としてPoissonの式:ｐ＝ｃ₁ρ^γが成立します。γは比熱比でγ≡Ｃ_p/Ｃ_vです。

ｐ＝ρＲＴ/ｍの両辺の対数微分を取ると,ｄlnｐ＝ｄinρ＋ｄlnＴですからｄlnＴ＝ｄlnｐ－ｄlnρですが,断熱変化ではｐ＝ｃ₁ρ^γよりlnｐ＝γlnρ＋lnｃ₁が成立するため,ｄlnρ＝ｄlnｐ/γです。

それ故,ｄlnＴ＝(1－1/γ)ｄlnｐ＝(Ｒ/Ｃ_p)ｄlnｐ,すなわちｄln{Ｔ/ｐ^(Ｒ/Ｃp)}＝0 ,またはＴ/ｐ^(Ｒ/Ｃp)＝(一定)です。

そこで,ある定数ｃ₂を与えθ≡ｃ₂Ｔ/ｐ^(Ｒ/Ｃp)によって温位という量θを定義すると,断熱変化はｄθ＝0 を意味するので,断熱変化の場合の温位の時間変化率はｄθ/ｄｔ＝0 です。

特に,圧力ｐの値がｐ₀₀のときθ＝Ｔとなるように定数ｃ₂を定めれば,θ＝Ｔ(ｐ₀₀/ｐ)^(Ｒ/Ｃp)という表現ができます。

さて,プルームの周辺大気について対流が生じない条件はｄｐ/ｄｚ＝－ρｇ(静力学平衡)ですが,これは対数表現ではｄlnｐ/ｄｚ＝－ｍｇ/(ＲＴ)です。ただし,鉛直上向きをｚ軸の正の向きとしています。

さらに,断熱変化:ｄθ＝0 ならｄlnｐ＝(Ｃ_p/Ｒ)ｄlnＴなので,熱の放射吸収のない断熱環境で対流が生じない条件はｄlnＴ/ｄｚ＝－ｍｇ/(Ｃ_pＴ),またはｄＴ/ｄｚ＝－ｍｇ/Ｃ_pと書けます。

そこで,鉛直高さｚと共にこの比率で温度が減衰すれば対流が生じず,しかも高さ方向の変化ｄｚに対する温度変化は断熱変化なので,温位勾配はゼロ:ｄθ/ｄｚ＝0 です。

この勾配:ｄＴ/ｄｚ＝－ｍｇ/Ｃ_pの温度減率(気温減率)を断熱減率といい,この状態にある大気を中立状態といいます。

ところで,通常の地上付近成層圏の大気の体積組成(分子数組成)は分子量28の窒素Ｎ₂が78.08％,分子量32の酸素Ｏ₂が20.95％,その他はアルゴン,二酸化炭素,水蒸気などです。

そこで,地上付近の乾燥大気の分子量:ｍ＝ｍ_dは28.964です。(小倉義光「一般気象学」より)

また,Ｎ₂,Ｏ₂のような２原子分子では理論的にはモル定積比熱は(5/2)Ｒ,モル定圧比熱はＣ_p＝(7/2)Ｒです。

そこで,この大気の断熱減率を計算するとΓ_d≡ｍ_dｇ/Ｃ_p＝28.964(g/mol)×9.8(ｍ/ｓ²)×2/{7×8.31Ｊ/(molＫ)}＝9.8(Ｋ/ｋｍ)です。

これは,100ｍ上昇する毎に気温が約１度下がるような状態です。この気温減率:Γ_d＝ｍ_dｇ/Ｃ_p＝9.8(Ｋ/ｋｍ)を乾燥断熱減率といいます。

一般には,大気は湿度ゼロの乾燥空気ではなく水分を含んでいます。

容積Ｖの中に質量Ｍ_dの乾燥空気と質量Ｍ_vの水蒸気(vapor)が混在しているとき,この空気の密度はρ＝(Ｍ_d＋Ｍ_v)/Ｖ＝ρ_d＋ρ_vです。

そして,ｐ_dを乾燥大気の分圧,ｅを水蒸気圧とすると,乾燥大気,水蒸気の状態方程式はそれぞれｐ_d＝ρ_dＲＴ/ｍ_d,ｅ＝ρ_vＲＴ/ｍ_vです。

全大気圧ｐはｐ＝ｐ_d＋ｅ＝ρ_dＲＴ/ｍ_d＋ρ_vＲＴ/ｍ_vですから,ｐ＝ρＲＴ{1＋(ｍ_v/ｍ_d)(ρ_v/ρ_d)}/{1＋(ρ_v/ρ_d)},またはｐ＝ρＲＴ^*/ｍ_d;Ｔ^*≡Ｔ{1＋(ｍ_v/ｍ_d)(ρ_v/ρ_d)}/{1＋(ρ_v/ρ_d)}と書けます。

ここで,Ｔ^*は仮温度(virtual temperature)と呼ばれる量で具体的にｍ_d＝28.964(ｇ/mol),ｍ_v＝18(ｇ/mol)を代入すると,水蒸気の混合比:ｑ≡ρ_v/ρ_d＝(ｍ_d/ｍ_v){ｅ/(ｐ－ｅ)}に対してＴ^*≒Ｔ(1＋1.61ｑ)/(1＋ｑ)＝Ｔ(1＋0.61ｒ)と表現されます。

ただし,ｒ≡ｑ/(1＋ｑ)＝ρ_v/(ρ_v＋ρ_d)は比湿(specific humidity)と呼ばれる量で,式から明らかにｒ≦1です。

乾燥大気の温度Ｔに対する安定条件:ｄＴ/ｄｚ＝－ｍ_dｇ/Ｃ_p＝－Γ_dは,Ｔを仮温度Ｔ^*に置き換えると近似的にそのまま成立すると考えれば,ｄＴ^*/ｄｚ＝－ｍ_dｇ/Ｃ_p＝－Γ_dなる式が安定条件に対応します。

左辺のＴ^*にＴ^*＝Ｔ(1＋0.61ｒ)を代入すると,(1＋0.61ｒ)(ｄＴ/ｄｚ)＝－Γ_dです。そこで,不飽和で凝結のない大気の断熱減率ΓはΓ＝Γ_d/(1＋0.61ｒ)で与えられ,Γ≦Γ_dです。

さらに,飽和(saturate)して凝結水滴があり潜熱(latent heat)の存在によって温度降下が鈍くなる場合の湿潤断熱減率をΓ_wとします。

飽和水蒸気圧をｅ_sと書き,単位質量当たりの液体の水が水蒸気になるのに必要な潜熱(蒸発熱)をＬとすれば,ｅ_sと温度Ｔの関係は理論的にはClausius-Clapeyronの式:ｄｅ_s/ｄＴ＝Ｌｅ_s/(ＲＴ²)で表わされます。

(上式は,lnｅ_s＝－Ｌ/(ＲＴ)＋const.またはｅ_s＝Ａexp{－Ｌ/(ＲＴ)}を意味しますが,2006年７/30の記事「気液平衡の統計力学」ではこれより精密な修正式:lnｅ_s＝lnＴ－Ｌ/(ＲＴ)－1＋const.またはｅ_s＝ＡＴexp{－Ｌ/(ＲＴ)}を与えています。)

例えばＴ＝100℃ではＬ＝597(cal/ｇ)ですが常温付近でのＬはＴの変動に鈍感なのでＬがＴに無関係と近似してこれを積分すれば,ln(ｅ_s/ｅ_s0)＝(ｍ_vＬ/Ｒ)(1/Ｔ₀－1/Ｔ)です。

これはｅ_sがＴだけに依存して周囲の大気圧ｐには依らない形の関係式です。

さて,乾燥空気と混合比ｑの水蒸気が混合した空気塊を考えます。

まず,単一気体の断熱変化は,0＝ｄlnθ＝(Ｃ_p/Ｒ)ｄlnＴ－ｄlnｐで表わされます。これはエネルギー収支を示す式と見るなら,0＝ＲＴｄθ/θ＝Ｃ_pｄＴ－ＲＴｄｐ/ｐと表現されます。

そこで,乾燥空気と水蒸気が混合しただけの大気において水蒸気が飽和に達していないとき,各成分の上記左辺の温位変化:ＲＴｄθ/θへの寄与は,それぞれＣ_pｄＴ－ＲＴ(ｄｐ_d/ｐ_d)＝0,ｑＣ_pvｄＴ－ｑＲＴ(ｄｅ/ｅ)＝0 です。

ここで,ｑはこの大気の乾燥空気に対する水蒸気の混合比でｑ＝ρ_v/ρ_d＝(ｍ_d/ｍ_v)(ｅ/ｐ_d)＝(ｍ_d/ｍ_v){ｅ/(ｐ－ｅ)}です。また,Ｃ_pvは水蒸気の定圧モル比熱です。このとき,もちろんｐ_d＝ｐ－ｅです。

この空気塊が上昇してｑ＝ｑ_sとなって飽和した後,さらに上昇してｄｑ_sの水蒸気が凝結したとすれば,空気塊はＬｄｑ_sだけの潜熱を受けるので,変化が断熱なら水蒸気が飽和した空気塊のエネルギー収支の式は次のようになります。

すなわち,Ｃ_pｄＴ＋ｑ_sＣ_pvｄＴ－ＲＴ{ｄ(ｐ－ｅ_s)/(ｐ－ｅ_s)＋ｑ_sｄｅ_s/ｅ_s}＋Ｌｄｑ_s＝0 となります。ｑ_sは飽和時の混合比です。

これに,ｄｅ_s/ｄＴ＝Ｌｅ_s/(ＲＴ²)を代入すれば,－Ｌｄ(ｑ_s/Ｔ)＝Ｃ_pｄＴ/Ｔ＋ｑ_sＣ_pvｄＴ/Ｔ－Ｒｄ(ｐ－ｅ_s)/(ｐ－ｅ_s)を得ます。

ｐ_d＝ρ_dＲＴ/ｍ_d,ｅ_s＝ρ_vＲＴ/ｍ_vですから,静力学平衡:ｄｐ/ｄｚ＝－ρｇ, or ｄ(ｐ_d＋ｅ_s)/ｄｚ＝－(ρ_d＋ρ_V)ｇは,ｄ(ｐ_d＋ｅ_s)/ｄｚ＝－ｇ(ｍ_dｐ_d＋ｍ_vｅ_s)/(ＲＴ),またはｄｐ＝－ｇ{ｍ_d(ｐ－ｅ_s)＋ｍ_Vｅ_s}/(ＲＴ)ｄｚと表現できます。

そこで,ｄ(ｐ－ｅ_s)＝－ｇ{ｍ_d(ｐ－ｅ_s)＋ｍ_Vｅ_s}/(ＲＴ)ｄｚ－ｄｅ_sを－Ｌｄ(ｑ_s/Ｔ)＝(Ｃ_p＋ｑ_sＣ_pv)ｄＴ/Ｔ－Ｒｄ(ｐ－ｅ_s)/(ｐ－ｅ_s)に代入すると,－Ｌｄ(ｑ_s/Ｔ)＝(Ｃ_p＋ｑ_sＣ_pv)ｄＴ/Ｔ＋ｇ{ｍ_d＋ｍ_Vｅ_s/(ｐ－ｅ_s)}ｄｚ/Ｔ＋Ｒｄｅ_s/(ｐ－ｅ_s)です。

これから,ｄＴ/ｄｚ＝－ｇ{ｍ_d＋ｍ_Vｅ_s(ｐ－ｅ_s)}/(Ｃ_p＋ｑ_sＣ_pv)－ＲＴ(ｄｅ_s/ｄｚ)/{(ｐ－ｅ_s)(Ｃ_p＋ｑ_sＣ_pv)}＋{ＬＴ/(Ｃ_p＋ｑ_sＣ_pv)}{ｄ(ｑ_s/Ｔ)/ｄｚ}を得ます。

ところが,ｑ_s＝(ｍ_d/ｍ_v){ｅ_s/(ｐ－ｅ_s)}よりｅ_s/(ｐ－ｅ_s)＝ｍ_vｑ_s/ｍ_dです。

そこで,ｄｑ_s/ｄｚ＝(ｍ_d/ｍ_v)[(ｄｅ_s/ｄｚ)/(ｐ－ｅ_s)－{ｄ(ｐ－ｅ_s)/ｄｚ}{ｅ_s/(ｐ－ｅ_s)²}]＝(ｍ_d/ｍ_v)[(ｄｅ_s/ｄｚ)/(ｐ－ｅ_s)－{ｑ_s/(ｐ－ｅ_s)}{ｄ(ｐ－ｅ_s)/ｄｚ]と書けます。

つまり,(ｄｅ_s/ｄｚ)/(ｐ－ｅ_s)＝(ｍ_v/ｍ_d)(ｄｑ_s/ｄｚ)＋{ｑ_s/(ｐ－ｅ_s)}{ｄ(ｐ－ｅ_s)/ｄｚ}です。これらの式を用いて,できるだけｅ_sを消去します。

以上から,湿潤断熱減率は,Γ_w≡－ｄＴ/ｄｚ＝Γ_d/(1＋ｑ_sＣ_pv/Ｃ_p)＋ｍ_V²ｇ/{ｍ_dＣ_p(1＋ｑ_sＣ_pv/Ｃ_p)}＋ｍ_vＲＴ(ｄｑ_s/ｄｚ)/{ｍ_dＣ_p(1＋ｑ_sＣ_pv/Ｃ_p}＋ＲＴ{ｑ_s/(ｐ－ｅ_s)}{ｄ(ｐ－ｅ_s)/ｄｚ}/{Ｃ_p(1＋ｑ_sＣ_pv/Ｃ_p－ＬＴ{ｄ(ｑ_s/Ｔ)/ｄｚ}/{Ｃ_p(1＋ｑ_sＣ_pv/Ｃ_p)}です。

この式は閉じていないので古い正野重方著の岩波全書の「気象力学」などには近似公式もあるようですがここでは主要課題ではないので割愛します。ただ,一般に湿潤断熱減率は一定値ではありません。

(なお,水蒸気の定圧比熱については,2008年１/６の記事「氷,水,水蒸気の比熱」,および2009年２/９の記事「水蒸気の比熱」などを参照してください。)　(補足終わり)　※

さて,浮力以外の気圧傾度力は無視できると考えて外気圧ｐは近似的に静力学平衡:ｄｐ/ｄｚ＝－ρｇを満足しているとします。

また,大気を流体と見るとき,対象とする運動でのReynolds数:Ｒeは非常に大きいため,1/Ｒeに比例する分子粘性項も無視します。

そして,プルームガスの周辺大気(ambient air)の局所速度,圧力,密度,温度,温位を,それぞれｖ_a,ｐ_a,ρ_a,Ｔ_a,θ_aとします。

特に,気圧については,∂ｐ_a/∂ｘ＝∂ｐ_a/∂ｙ＝0,∂ｐ_a/∂ｚ＝－ρ_aｇが成立しており,プルームガス圧については近似的にｐ_p＝ｐ_aが成立するとします。

この場合,周辺大気の運動方程式は,ρ_aｄｖ_a/ｄｔ＝－∇ｐ_a－ρ_aｇｋ＝0 です。ただしｉ,ｊ,ｋは,それぞれｘ,ｙ,ｚ軸の方向単位ベクトルを表わします。

この方程式の解は明らかにｖ_a＝constで与えられます。

一方,プルームガスの満たす運動方程式はρ_pｄｖ_p/ｄｔ＝－∇ｐ_p－ρ_pｇｋ～(－∇ｐ_a－ρ_aｇｋ)－(ρ_p－ρ_a)ｇｋです。

すなわち,この近似での運動方程式はρ_p(ｄｖ_p/ｄｔ)＝(ρ_a－ρ_p)/ｇｋ,またはｄｖ_p/ｄｔ＝{(ρ_a－ρ_p)/ρ_p}ｇｋです。これらの右辺が浮力項(buoyancy term)です。

状態方程式:ｐ_a＝ρ_aＲＴ_a/Ｍ,ｐ_p＝ρ_pＲＴ_p/Ｍから,ｐ_p～ｐ_aの近似ではρ_a/ρ_p ～Ｔ_p/Ｔ_a～θ_p/θ_aなので温位偏差θ'をθ'≡θ_p－θ_aと定義すれば,(ρ_a－ρ_p)/ρ_p～ (θ_p－θ_a)/θ_a＝θ'/θ_aです。

周辺大気の温位θ_aの定義:θ_a≡Ｔ_a(ｐ₀₀/ｐ_a)^(Ｒ/Ｃp)において特にｐ₀₀～ｐ_aなるパラメータを採用すればθ_a～Ｔ_aなので,結局この近似(Boussinesq近似)での運動方程式は,(3)ｄｖ_p/ｄｔ＝(ｇθ'/Ｔ_a)ｋとなります。

連続方程式が質量の保存を意味するのに対して,運動方程式は運動量の保存を意味します。

さらに,ｖ_p≡ｖ_a＋ｖ'と書き,特に速度偏差:ｖ'の鉛直ｚ方向の成分(or第３成分)をｗ'とします。また,周辺大気の風速ｖ_aは鉛直成分を持たないとして,その向きをｘ軸(or１軸)の正の向きに取ります。

すると,ｖ_p≡ｕｉ,ｖ_p3＝ｖ'₃＝ｗ'と書けます。さらにθ_a,ρ_a,ｕはｚだけの関数でありθ_a＝θ_a(ｚ),ρ_a＝ρ_a(ｚ),ｕ＝ｕ(ｚ)のように表わせるとします。

そして水平プルーム断面上ではｗ'と浮力ｇθ'/Ｔ_aが比例する:ｇθ'/Ｔ_a∝ ｗ'と仮定します。比例係数をｃ(ｚ)とするとｇθ'/Ｔ_a＝ｃ(ｚ)ｗ'です。

定常状態が実現されていると仮定し,面積分することによってプルームを水平面で切った断面上での平均量に対する方程式を求めます。

ただし,プルーム断面の外では明らかにθ'＝0,かつｖ'＝0 で特にｗ'＝0 です。

この結果は,まず浮力の効果についての方程式:(4)ｄＦ_z/ｄｚ＝－ｓＶを与えます。ただし,Ｆ_zは浮力流束で,(5)Ｆ_z≡∫_P(ｇρ_pθ'ｗ'/Ｔ_a)ｄｘｄｙ/(πρ_a)で定義されます。

Ｖは鉛直容積流束で,(6)Ｖ≡∫_Pρ_pｗ'ｄｘｄｙ/(πρ_a)です。

また,ｓは安定度パラメータと呼ばれる量で,(7)ｓ≡(ｇ/Ｔ_a)(∂θ_a/∂ｚ)です。

(証明):(2)ｄθ_p/ｄｔ＝0 より,∂θ_p/∂ｔ＋ｖ_p∇θ_p＝0ですから,定常状態:∂θ_p/∂ｔ＝0 ではρ_pｖ_p∇(θ_a＋θ')＝ 0です。

そして,∂θ_a/∂ｘ＝∂θ_a/∂ｙ＝∂θ'/∂ｘ＝∂θ'/∂ｙ＝0 よりρ_pｖ_p∇(θ_a＋θ')＝－ρ_pｗ'(∂θ_a/∂ｚ)です。

これに,(1)∇(ρ_pｖ_p)＝0 にθ'を掛けたものを加えると∇(ρ_pθ'ｖ_p)＝－ρ_pｗ'(∂θ_a/∂ｚ)を得ます。

これを全水平面で積分すると∫_∞∇(ρ_pθ'ｖ_p)ｄｘｄｙ＝－(∂θ_a/∂ｚ)∫_∞ρ_pｗ'ｄｘｄｙですが,無限遠では明らかにθ'＝0 で,プルーム断面を含む全水平面領域では全ての量はｚに依らず一定ですから,左辺の(∂/∂ｚ)も積分記号の外に出せます。

プルーム水平断面の外ではρ_pθ'ｗ'＝ρ_pｗ'＝0 なので,上式は(ｄ/ｄｚ)∫_Pρ_pθ'ｗ'ｄｘｄｙ＝－(∂θ_a/∂ｚ)∫_Pρ_pｗ'ｄｘｄｙとなります。

周辺大気の密度ρ_aは"上昇高度までの層当たり平均密度＝(一定)"と考え,両辺に(ｇ/Ｔ_a)/(πρ_a)を掛けると{∂(1/Ｔ_a)}{(ｄ/ｄｚ)∫_Pρ_pｇθ'ｗ'ｄｘｄｙ}/(πρ_a)＝－(ｇ/Ｔ_a)(∂θ/∂ｚ){∫_∞ρ_pｗ'ｄｘｄｙ}/(πρ_a)です。

よって,ｄＦ_z /ｄｚ＝－ｓＶ－(1/Ｔ_a)(ｄＴ_a/ｄｚ)Ｆ_zを得ました。

ところが,ｄＴ_a/ｄｚ～ｄθ_a/ｄｚ,Ｆ_z ～ｇＶ(＜θ'＞/θ_a）なので(1/Ｔ_a)(ｄＴ_a/ｄｚ)Ｆ_z～ (ｓＶ)・(＜θ'＞/θ_a)＜＜(ｓＶ)ですから,右辺第１項:(ｓＶ)に比べて右辺第２項:(1/Ｔ_a)(ｄＴ_a/ｄｚ)Ｆ_zを無視します。

そこで,結局ｄＦ_z /ｄｚ＝－ｓＶを得ます。(証明終わり)

次に,運動量に対する方程式は,(8)ｄ(ｖＶ)/ｄｚ＝(Ｆ_z/ｗ)ｋ＋ｖ_a(ｄＶ/ｄｚ)＋(抵抗力)です。ただし,ｖＶは(9)ｖＶ≡∫_Pｖ_pρ_pｗ'ｄｘｄｙ/(πρ_a)で定義される運動量流束です。

この式で定義されるベクトルｖをプルームの速度と考えます。ｗはｖのｚ成分です。

(証明):(1)∇(ρ_pｖ_p)＝0 より,ｖ_p・∇(ρ_pｖ_p)＝0 ,つまりΣ_j=1³ｖ_pi∂_j(ρ_pｖ_pj)＝0 です。

　また,(3)ｄｖ_p/ｄｔ＝(ｇθ'/Ｔ_a)ｋより,定常状態ではρ_p(ｖ_p∇)ｖ_p＝(ｇρ_pθ'/Ｔ_a)ｋ,つまりΣ_j=1³ρ_pｖ_pj∂_jｖ_pi＝(ｇρ_pθ'/Ｔ_a)δ_i3です。

したがって,∇(ρ_pｖ_pｖ_pi)＝Σ_j=1³ｖpi∂_j(ρ_pｖ_pj)＋Σ_j=1³ρ_pｖ_pj∂_jｖ_pi＝(ｇρ_pθ'/Ｔ_a)δ_i3です。

　プルームの進行方向をｎとするとき,プルームの断面と側面で囲まれた領域ＶにGaussの定理を適用すると,∫_V∇(ρ_pｖ_pｖ_pi)ｄＶ＝∫_P+ρ_pｖ_pnｖ_piｄＳ－∫_P-ρ_pｖ_pnｖ_piｄＳ－Δｎ{1＋(ｄｒ/ｄｎ)²}^1/2(∫_Cρ_pｖ_pｖ_piｄｓ)です。　

ただし,ここでのｒはプルーム半径(plume-radius)を意味します。プルーム断面は等方的なので,これは断面の形を円と想定したときの半径です。

－Δｎｄｓはプルーム側面の外向き法線の向きを持ちΔｎは側面ｎ方向の厚みです。それ故,ｄｓは周の微小線素の長さを持つ内向き法線方向ベクトルです。

　故に,Δ∫_Pρ_pｖ_pnｖ_piｄＳ≡∫_P+ρ_pｖ_pnｖ_piｄＳ－∫_P-ρ_pｖ_pnｖ_piｄＳ＝Δｎ{1＋(ｄｒ/ｄｎ)²}^1/2(∫_Cρ_pｖ_pｖ_piｄｓ)＋δ_i3Δｎ∫_P(ｇρ_pθ'/Ｔ_a)ｄＳです。

　よって,(ｄ/ｄｎ)(∫_Pρ_pｖ_pnｖ_piｄＳ)＝{1＋(ｄｒ/ｄｎ)²}^1/2(∫_Cρ_pｖ_pｖ_piｄｓ)＋δ_i3∫_P(ｇρ_pθ'/Ｔ_a)ｄＳを得ます。

　一方,連続の方程式から(ｄ/ｄｎ)(∫_Pρ_pｖ_pnｄＳ)＝{1＋(ｄｒ/ｄｎ)²}^1/2(∫_Cρ_pｖ_pｄｓです。

　さらにｖ_pi＝ｖ_ai＋ｖ'_iですから,(ｄ/ｄｎ)(∫_Pρ_pｖ_pnｖ_piｄＳ)＝δ_i3∫_P(ｇρ_pθ'/Ｔ_a)ｄＳ＋ｖ_ai(ｄ/ｄｎ)(∫_Pρ_pｖ_pnｄＳ)＋{1＋(ｄｒ/ｄｎ)²}^1/2(∫_Cρ_pｖ_pｖ'_iｄｓ)を得ます。

　特にｎをｚ方向のｋに選べばｖ_pn＝ｗ'より,(ｄ/ｄｚ)(∫_Pρ_pｗ'ｖ_piｄＳ)＝δ_i3∫_P(ｇρ_pθ'/Ｔ_a)ｄＳ＋ｖ_i(ｄ/ｄｚ)(∫_Pρ_pｗ'ｄＳ)＋(抵抗力)です。

　ここで,(抵抗力)≡{1＋(ｄｒ/ｄｚ)²}^1/2(∫_Cρ_pｖ_pｖ'_iｄｓ)であり,ｄｓ＝{1＋(ｄｒ/ｄｚ)²}^-1/2(－ｄｙｉ＋ｄｘｊ)＋{1＋(ｄｒ/ｄｚ)²}^-1/2(ｄｒ/ｄｚ)ｄｌｋです。(ｄｌ＝(ｄｘ²＋ｄｙ²)^1/2)

　一方,∫_P(ｇρ_pθ'/Ｔ_a)ｄＳ＝(1/ｗ)∫_P(ｇ/Ｔ_a)(ρ_pθ'ｗ)ｄＳ＝(1/ｗ)∫_P(ｇ/Ｔ_a)(ρ_pθ'ｗ')ｄＳ＋(1/ｗ)∫_P(ｇ/Ｔ_a){ρ_pθ'(ｗ－ｗ')}ｄＳです。

ｖＶ≡∫_Pｖ_pρ_pｗ'ｄＳ/(πρ_a),Ｖ≡∫_Pρ_pｗ'ｄＳ/(πρ_a),およびｖ_p3＝ｗ'よりｗ＝∫_Pρ_pｗ'²ｄＳ/∫_Pρ_pｗ'ｄＳですから,ｗ－ｗ'＝{∫_Pρ_pｗ'²ｄＳ－ｗ'∫_Pρ_pｗ'ｄＳ}/{∫_Pρ_pｗ'ｄＳ}です。

　仮定によりｇθ'/Ｔ_a＝ｃ(ｚ)ｗ'を代入すると,(1/ｗ)∫_P(ｇ/Ｔ_a){ρ_pθ'(ｗ－ｗ')}ｄＳ＝{ｃ(ｚ)/ｗ}[{∫_Pρ_pｗ'ｄＳ}{∫_Pρ_pｗ'²ｄＳ}－{∫_Pρ_pｗ'²ｄＳ}{∫_Pρ_pｗ'ｄＳ}]/{∫_Pρ_pｗ'ｄＳ}＝0 です。

　したがって,∫_P(ｇρ_pθ’/Ｔ_a)ｄＳ＝(1/ｗ)∫_P(ｇ/Ｔ_a)(ρ_pθ'ｗ')ｄＳ故,(ｄ/ｄｚ)(∫_Pρ_pｗ'ｖ_piｄＳ)＝(1/ｗ)δ_i3∫_P(ｇ/Ｔ_a)(ρ_pθ'ｗ')ｄＳ＋ｖ_i(ｄ/ｄｚ)(∫_Pρ_pｗ'ｄＳ)＋(抵抗力),

　または,(ｄ/ｄｚ)(∫_Pｖ_pρ_pｗ'ｄＳ)＝(1/ｗ){∫_P(ｇ/Ｔ_a)(ρ_pθ'ｗ')ｄＳ}ｋ＋ｖ_a(ｄ/ｄｚ)(∫_Pρ_pｗ'ｄＳ)＋(抵抗力)なる式を得ます。

　最後の表式は,まさにｄ(ｖＶ)/ｄｚ＝(Ｆ_z/ｗ)ｋ＋ｖ_a(ｄＶ/ｄｚ)＋(抵抗力)そのものです。

※特にプルームが軸対称で完全に鉛直向きなら(抵抗力)＝{1＋(ｄｒ/ｄｚ)²}^1/2(∫_Cρ_pｖ_pｖ'_iｄｓ)＝0でｄ(ｗＶ)/ｄｚ＝Ｆ_z/ｗです。※(証明終わり)

　さて,ｄＦ_z/ｄｚ＝－ｓＶ,ｄ(ｗＶ)/ｄｚ＝Ｆ_z/ｗを解くために次のような初期条件を与えます。

まず,初期速度は吐出速度に等しいとします。すなわち,(10a)ｖ＝ｗ₀ｋです。そして初期運動量流束は,(10b)ｖＶ＝(ｗＶ)ｋ＝Ｆ_mｋ≡(ρ₀/ρ_a)ｗ₀²ｒ₀²ｋとします。ρ₀は排出ガス密度,ｒ₀は煙突口の口径です。

さらに,初期の浮力流束は,(10c)Ｆ_z＝Ｆ≡(1－ρ₀/ρ_a)ｇｗ₀ｒ₀²で与えられるとします。

つまり,初期容積流束をＶ₀≡ｗ₀ｒ₀²とすれば初期の総質量流束はπρ₀Ｖ₀＝πρ₀ｗ₀ｒ₀²ですから,定義により運動量流束はＦ_m＝πρ₀ｗ₀Ｖ₀/(πρ_a)＝(ρ₀/ρ_a)ｗ₀Ｖ₀です。

そして,初期浮力流束はＦ_z＝Ｆ＝(1－ρ₀/ρ_a)ｇＶ₀～{(Ｔ₀－Ｔ_a)/Ｔ₀}ｇＶ₀となるわけです。

ところで,排出源がhot-sourceの場合の初期浮力流束は,(11)Ｆ＝ｇＱ_H/(πρ_aｃ_pＴ_a)で与えられます。

ここで,Ｑ_Hは排出されるガスによって単位時間に排出される熱量(cal/sec)であり,ｃ_pは空気の単位質量当たりの熱容量(定圧比熱)(cal/(kg・Ｋ))です。

なぜなら,Ｆ_z＝Ｆ＝(ｇθ'/Ｔ_a)(πρ₀ｗ₀ｒ₀²)/(πρ_a)＝{ｇ/(πρ_aＴ_a)}(θ'πρ₀Ｖ₀)ですが,πρ₀Ｖ₀は総質量流束なので排出ガスによって単位時間に外気に付加される熱量はＱ_H＝ｃ_pＴ'πρ₀Ｖ₀ ～ｃ_pθ'πρ₀Ｖ₀に等しいからです。　

つまり,θ'πρ₀Ｖ₀＝Ｑ_H/ｃ_pによってＦ＝ｇＱ_H/(πρ_aｃ_pＴ_a)が得られるわけです。

この場合,排熱強度Ｑ_Hと排出ガス密度ρ₀の情報を得ることにより逆にＶ₀＝ｗ₀ｒ₀²が求まります。

そこで,煙突口径ｒ₀のデータから吐出速度ｗ₀が算定されて,結局全ての初期情報を得ることができます。

ただし,後述するように３つの未知数Ｆ_z,ｗ,Ｖに対し微分方程式は２つしかないので解を決めるにはもう１つ方程式が必要です。

今日はここまでですがまだ続きます。

参考文献：G.A.Briggs "Plume Rise",U.S.Atomic Energy Commision Division of Technical Information(1969),小倉義光「一般気象学」(東京大学出版会),小倉義光「気象力学通論」(東京大学出版会)

台風の進路(コリオリの力)(再々掲)

　丁度,台風シーズンです。

　そこで,私自身の息抜き(手抜き？)として,４年前に書いたテキストだけの解説記事:「台風の進路（コリオリの力） 」に図を入れて再掲してみます。

(※再掲記事):2006年８月16日(水)

　そろそろ台風が頻発する季節になりました。

　今日は,北半球では赤道付近で発生し,海域から多量の水分やエネルギーを吸収しながら発達して北上する台風が,なぜ右（東）の方に進路を変えていくのか？そして,なぜ上空から見て左巻き（反時計回り）の風が吹くのかという,ごくありふれた疑問について解説してみます。

　例として,ちょっと古いけど左回転しているLPレコードがあり,その上に"一寸法師"よりも小さい小人が乗っているという仮想的な状況を考えてみます。(左回転は仮定であって実際のLPレコードは裏から見ない限り右回転(時計回り)です。)

　レコードの中心は地球の北極に相当し,レコードの最大半径の部分は地球の赤道に相当します。　

　まず,レコードの回転している"最大半径＝赤道"の上にいる小人が"レコードの中心＝北極"めがけて真っ直ぐに小石を投げたとします。

　本人は真っ直ぐ中心に向かって投げたつもりでも,小石が手を離れた瞬間には慣性によって小石はレコードの回転スピードと同じ速度で右に向かう接線速度を持ちますから,実はそれは中心の方向に向かって真っ直ぐに飛ぶのではなくて,次第に右の方に逸れていくということになりますね。

　↑ここで右というのは,小石を投げた小人にとっての右です。(わかっているとは思いますが念のため))

　次に逆に"中心＝北極"の上に小人がいて今度は"最大半径＝赤道"めがけてやはり小石を投げたとします。

　今度は北極で小人は自転しているかもしれませんが,スピンの回転半径はゼロなので,その慣性による小石の左右方向への速度はゼロですから確かに"真っ直ぐ"進むはずです。

　ところが,レコードの上,つまり北半球の地球上にいる人は"左から右＝西から東"に回転しています。その人から見ると"上＝北"から真っ直ぐ飛んでくる小石は"下＝南"から見て"左に左に（西に西に）"逸れていくように見えます。

　逆に"小石を投げた方＝北"から見ると,見かけ上はやはり右の方に逸れていくわけですね。

　小石を台風だとみなし地球の自転の角速度をΩとすると　Ωは"360度＝2πラジアン(rad)"を24時間で回転する角速度です。

　地球の半径をRとし,緯度をθとすると,そこでの回転半径はRcosθですから,回転の接線速度はΩRcosθです。

　したがって"赤道"での接線速度は最大速度"ΩR＝時速1667ｋｍ＝秒速463ｍ"ですが,日本付近の緯度:θ＝35度での接線速度は"ΩRcosθ＝秒速379ｍ"で,日本付近では回転速度は赤道より"約２割＝秒速80ｍ"くらい減少しています。

　したがって,赤道付近で発生した台風は地球のまさつにより"ΩR＝時速1667ｋｍ＝秒速463ｍ"の地球にひきずられて慣性による右向きの速度を持っていて,その右向き速度は北上しても全く変わらないものです。

　しかし,地球自身の回転速度は緯度が上がるにつれて次第に小さくなるものですから,日本付近では１秒間に80ｍくらいの割合で,右（東）へ右へと逸れていくことになります。

　先にＬＰレコードの例で述べたように仮に北極で台風が発生して南下したとしてもそれは右に逸れていきます。

　実は北半球ではどこから投げた石も見かけ上,右に逸れていくわけです。

　例えばスナイパー:ゴルゴ13が１ｋｍ遠方の標的を狙って狙撃しても,弾丸は僅かに右に逸れていくのでそれを勘定に入れて狙う必要があるわけです。

　もしも南半球なら逆に左に逸れるわけですね。

　こうした北半球で右にそれる現象は結局,遠心力などと同じく"見かけの力＝慣性力"が働いていると考えることができて,それを発見者の名前にちなんでコリオリ(Coｒiolis)の力と言います。

　そして北半球での台風を考えると,台風ですから"中心＝目の部分"の気圧が最低でまわりの気圧は目の部分のそれより高いわけです。

　風はどのように吹くか,というと水が高いところから低いところへと流れるように,風も気圧の高いところから低いところ目指してその気圧のスロープに沿って吹いていくわけです。

　もしもコリオリの力がなければ,風は"外周部から中心に向かって一直線に進む＝落下していく"はずなのですが,コリオリの力によって気圧のスロープも右にねじれてしまっています。

　それ故,風は外周部から中心に向かっていくときに,右にフックして逸れていきながら,最後は中心の気圧最低の目に向かっていくことになり,そのために左巻き（反時計回り）になるのです。

　南半球での台風は逆に右巻き（時計回り）ですね。

　どこかの"バカな大学教授"は,風呂の水が排水口へと流れていき排水されるときに,北半球では左巻き（反時計回り）ですが,赤道を越えて南半球に入ったとたんに右巻き（時計回り）に変わる,などと主張したと聞きましたが,それは誤りですね。

　風呂の水程度の規模では地球自転の影響などは出てきません。

　たまたま排水口付近で左巻きの角運動量を持っていたら左巻きになり,逆なら右巻きになるというだけで,それはカオス現象,つまり偶然の産物でしかありません。

　しかし台風くらいの規模になると地球の自転がもろに効いてきます。

　遠心力の加速度は緯度θでΩ²Ｒcosθですから,最大でも重力の加速度の0.3％程度です。

　北極で体重100kg重の人が赤道で体重を測っても300ｇ重くらいしか軽くはなりません。一方,コリオリの力の加速度は台風の北上の速度をｖとして2Ωsinθ×ｖです。

　Ω＝7,2×10^-5/ｓですから,緯度θが35度で台風の北上速度が100mを10秒で走る程度の時速36km程度なら,加速度ａ＝8.3×10^-4m/ｓ²であり,重力ｇ＝9.8m/ｓ²の約0.01％程度です。

　そこで,コリオリの力の加速度は最大で重力加速度の0.3％程度しかない遠心力のさらに1/30程度にすぎませんが,台風程度の規模だとそれがかなり効いてきます。

　地球自転の証拠であるとされるフーコー(Foucauｌt)の振り子をこのコリオリ力で説明することもできます。

　ニュートン(Isac Newton)は"慣性系の同等性＝ガリレイの相対性原理"は認めても"回転系を含む非慣性系の同等性＝マッハ(Mach)原理 → 一般相対性原理"を認めることをあきらめました。

　そして,彼が"絶対座標系＝絶対空間"に固執せざるを得なかったのも,こうした"遠心力やコリオリ力の絶対性"を解消する道はない,という考えからだったという話もあります。

　こうした"見かけの力＝慣性力"の扱いはとても悩ましいところがあります。

※以上,再掲記事です。

PS:この連休では,2008年11/23の記事「超弦理論(7)(弦の相互作用と頂点演算子)」を再編集するため,図1-5から図1-12まで都合８個も図を入れましたが,これは私には慣れない作業なので疲れました。出来映えはいかに？。。。　

震源の探知(大森公式等)

　今日は,ちょっと手抜きのツナギで2006年11/14の記事 「結晶内での弾性波（地震波）」の続きとして,中学か高校の地学で習った簡単な知識を披瀝してお茶をにごします。

　年末,体調はよくも悪くもないですがヒマ人なりに予定がこみ合っていて貧乏なこともありフリ－マンをエンジョイというわけにはいきません。

　ではまず,記事の再掲から始めます。

　(※再掲)

　今日は等方的とは限らない一般の結晶内での弾性波,特に地震波について考察してみます。 

弾性体の密度をρ,応力テンソルを{σ_jk},それを構成する部分の局所的速度をｕ＝{ｕ_j}とすると,この弾性体が従うべき運動方程式は流体方程式と同じくρｄ²ｕ_j/ｄｔ²＝∂σ_jk/∂ｘ_kで与えられます。

そして,歪み速度テンソル{ｕ_jk}をｕ_jk≡(∂ｕ_j/∂ｘ_k＋∂ｕ_k/∂ｘ_j)/2で定義すれば,これは{σ_jk}と同じく反対称の単なる回転を除いた対称テンソルです。これは６個の独立成分を持ちます。

微小変位に対しては適切な線形弾性体近似では,フック(Hooke)の法則:σ_jk＝Ｃ_jklmｕ_lmが成立するため,先の運動方程式はρｄ²ｕ_j/ｄｔ²＝Ｃ_jklm(∂²ｕ_m/∂ｘ_k∂ｘ_l)となります。

Ｃ_jklmは,一般に81個の成分を持ちますが,{σ_jk}も{ｕ_jk}も対称テンソルなのでｊとｋ,ｌとｍについて対称ですから,独立成分は６×６＝36個となります。

また,歪みエネルギーＷはΔＷ＝σ_jkΔｕ_jkから求まり,対称２次形式Ｗ＝(1/2)Ｃ_jklmｕ_jkｕ_lmになるので,Ｃ_jklmはｊ,ｋとｌ,ｍの交換についても対称であり,結局,その独立成分は(30/2)＋6＝21個となります。

この方程式の平面波の解をｕ_j＝ｕ_0jexp[i(ｋｒ－ωｔ)]として代入すると,ρω²ｕ_j＝Ｃ_jklmｋ_kｋ_lｕ_mとなります。

これは書き直すと(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l)ｕ_m＝0 とも書けます。

この１次方程式が自明でない解を持つためには,３行３列の行列係数:(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l) (ｊ,ｍ＝1,2,3)の行列式がゼロ,

すなわち,det(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l)＝0 が成立する必要があります。

そして,ω²を未知数としてこの方程式を解けば,ω²の３個の解が得られ,それらはｋの関数となります。 

ところで,もしも結晶が等方性弾性体であれば歪みエネルギー:

Ｗ＝(1/2)Ｃ_jklmｕ_jkｕ_lmは座標系の回転に対して不変なスカラーであるはずですが,ｕ_jkの２次形式の形で得られる独立な不変スカラーはｕ_kk²とｕ_jkｕ_jkのみです。

そこで,適当な係数λ,μを選んでＷ＝(1/2)λｕ_kk²＋μｕ_jkｕ_jkと書けるはずです。

そこで,応力テンソル{σ_jk}＝{Ｃ_jklmｕ_lm}はσ_jk＝λｕ_llδjk＋2μｕ_jkと書けます。ここに弾性定数λ,μはラメ(Lame)の定数と呼ばれています。

このときには運動方程式も,

ρｄ²ｕ/ｄｔ²＝(λ＋μ)∇(∇ｕ)＋μ∇²ｕ

と簡単になります。

先に挙げた１次方程式の係数の行列要素は,

ρω²δ_jm－Ｃ_jklmｋ_kｋ_l＝(ρω²－μｋ²)δ_jm－(λ＋μ)ｋ_jｋ_m

となります。

ｋの向きをｘ軸の正の向きに取ると,ｋ₁＝ｋ＝|ｋ|,ｋ₂＝ｋ₃＝0 ですから,行列[(ρω²－μｋ²)δ_jm＋(λ＋μ)ｋ_jｋ_m]は対角成分が[ρω²－(λ＋2μ)ｋ²]と２つの(ρω²－μｋ²)の対角行列になります。

したがって,det(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l)＝0 の解は,

Ｖ_P²≡(ω_P/ｋ)²＝(λ＋2μ)/ρ,Ｖ_S²≡(ω_S/ｋ)²＝μ/ρ

となります。

それぞれの速度に対応する平面波は,速度Ｖ_Pで伝播する波の伝播方向への振動である"縦波＝Ｐ波"と,速度Ｖ_Sで伝播する波の伝播方向に垂直な方向の振動である"横波＝Ｓ波"です。

しかし,もしも異方性の弾性体の場合ならω＞0 の３つの解ω(ｋ)はｋの関数として１次の同次式ではあっても,単なる１次関数になるとは限らず,弾性波は一般に分散性の波でもあると思われます。

そこで,伝播速度は"単一波の速度＝位相速度"ではなく,群速度Ｖ＝∂ω/∂ｋで与えられると考えられます。

ところで,σ_jkやｕ_jkを[σ]＝^t(σ₁₁,σ₂₂,σ₃₃,σ₂₃,σ₃₁,σ₁₂)のように６次元の列ベクトルで表わすボイト(Vogit)の表記という表記法があります。

この表記で表現すると,フックの法則は[σ]＝Ｃ[ｕ]と簡単になります。ここでＣは６行６列のテンソル行列です。

また,このとき,特に結晶が等方性弾性体ならＣの成分のうちで独立なのはやはり２つだけです。

先に述べたラメの定数はλ＝Ｃ₁₂,およびμ＝Ｃ₄₄＝(Ｃ₁₁－Ｃ₁₂)/2と表わされます。

例えば結晶が立方体構造をしている立方晶系では独立成分は３つです。

つまり,Ｃ₁₁＝Ｃ₂₂＝Ｃ₃₃,Ｃ₁₂＝Ｃ₂₃＝Ｃ₃₁,かつ４～６行の成分は４～６列しか成分のない対角行列であってＣ₄₄＝Ｃ₅₅＝Ｃ₆₆です。

結局,この結晶構造を示すにはＣ₁₁とＣ₁₂とＣ₄₄の３つだけあれば十分である,ということになります。

この条件ではｘ軸をｋの向きに取り,先のdet(ρω²δ_jm－Ｃ_jklmｋ_kｋ_l)＝ 0 でのＣ_jklmをボイトの表記でＣを表わせば,ξ＝(ω/ｋ)²,ａ＝Ｃ₁₁/ρ,ｂ＝Ｃ₄₄/ρ,c＝Ｃ₁₂/ρとして(ξ－ａ)(ξ－ｂ)²＝0 です。

この解はξ＝ａ,ｂとなります。

つまり速度をＶ_P,Ｖ_Sとすると,Ｖ_P²＝Ｃ₁₁/ρとＶ_S²＝Ｃ₄₄/ρの２つの速度が得られます。

これは等方性弾性体でのＣ₄₄＝μ,Ｃ₁₁＝λ＋2μのケースと一致します。

一方,六方晶系ではＣ₁₁＝Ｃ₂₂で,Ｃ₁₃＝Ｃ₂₃,また４～６は対角行列でＣ₄₄＝Ｃ₅₅,Ｃ₆₆＝(Ｃ₁₁－Ｃ₁₂)/2ですから,結局のところＣ₁₁,Ｃ₁₂,Ｃ₁₃,Ｃ₃₃,Ｃ₄₄の５つだけがあれば十分となります。

結果だけ書くと,Ｖ_P²＝Ｃ₁₁/ρ,およびＶ_S²＝Ｃ₄₄/ρ,(Ｃ₁₁－Ｃ₁₂)/(2ρ)となり,"横波＝Ｓ波"には２つの速度があるので地震の観測ではＳ波のほうに二重の波が観測されると予測されます。

ただし,異方性の結晶では一般に波は完全にＳ波とＰ波になるように行列が対角形になるとは限らないので,これらの計算においては偶々波の進行方向が結晶の対称軸と一致したり直交していたりする特別なケースだけしか扱っていません。

一般的な扱いについては,暇があってその気になればまた記事にするかもしれません。 

参考文献；ランダウ＝リフシッツ 著「弾性理論」(東京図書),角谷典彦 著「連続体力学」(共立出版)

(再掲終了※)

　と書きましたが,ここからもっと身近な地震の話題へとつなげます。

等方性媒質を仮定すると,Ｖ_P＝{(λ＋2μ)/ρ}^1/2,Ｖ_S＝(μ/ρ)^1/2ですからＶ_P＞Ｖ_Sです。

そこで,震源Ｏから地震を感知する場所:ＡまでＰ波,およびＳ波が到着するまでの時間を,それぞれＴ_P,およびＴ_Sと書けばＴ_P＝∫_O^A(1/Ｖ_P)ｄｓ,Ｔ_S＝∫_O^A(1/Ｖ_S)ｄｓです。

Ｖ_P＞Ｖ_SですからＴ_S＞Ｔ_Pとなります。

つまり,まず縦揺れのＰ波だけが到着しその後に横揺れＳ波も到着して両方が重ね合わされた大きいゆれが始まるのですね。 

特にＯＡ間の到るところで媒質が同じ均等な弾性体であれば,この区間でＶ_P,Ｖ_Sが一定です。

そこで,このときにはＯＡ間の波の進行距離(直線距離でなくてもよい)をｄ_A≡∫_O^AｄｓとするとＴ_P＝ｄ_A/Ｖ_P,Ｔ_S＝ｄ_A/Ｖ_Sです。

これから,辺々を引き算すると,

Ｔ_S－Ｔ_P＝(ｄ_A/Ｖ_S－ｄ_A/Ｖ_P)＝ｄ_A(1/Ｖ_S－1/Ｖ_P) です。

それ故,Ｔ_SP≡Ｔ_S－Ｔ_P,ｋ≡(1/Ｖ_S－1/Ｖ_P)^-1＝Ｖ_PＶ_S/(Ｖ_P－Ｖ_S)と定義すれば,ｄ_A＝ｋＴ_SPという比例関係の公式を得ます。 

　この式は発見者の大森房吉氏(1918)の名前を取って,大森公式(Ohmori's law)と呼ばれます。　

　(下図は中越地震本震の地震計記録からです。)

　　　　　

こうした等方的で均質な媒質を伝わる場合,弱いＰ波だけが来て「あ,ひょっとして地震かな？」と思ってから本格的で大きな揺れがくるまでの初期微動の時間をＴ_SPとするとその地点Ａから震源Ｏまでの距離ｄ_AはＴ_SPに比例します。 

初期微動時間が長いほど震源までの距離(地理的な距離だけでなく震源深さまでも含めた距離)が長いということが言えます。もしも震源が直下にあれば縦揺れは上下震動,横揺れは左右震動になります。 

さて,標高ｈが違う３地点Ａ,Ｂ,Ｃがありその座標がそれぞれ(ｘ_A,ｙ_A,ｈ_A),(ｘ_B,ｙ_B,ｈ_B) (ｘ_C,ｙ_C,ｈ_C)であるとき,上記のような公式に基づいて初期震動時間の観測から観測点から震源までの距離ｄ_A,ｄ_B,ｄ_Cが全てわかったとします。 

このとき,未知の震源Ｏの座標を(ｘ_O,ｙ_O,ｈ_O)と仮定すれば,地震波が全て直進なら,

ｄ_A²＝(ｘ_A－ｘ_O)²＋(ｙ_A－ｙ_O)²＋(ｈ_A－ｈ_O)²,

ｄ_B²＝(ｘ_B－ｘ_O)²＋(ｙ_B－ｙ_O)²＋(ｈ_B－ｈ_O)²,

ｄ_C²＝(ｘ_C－ｘ_O)²＋(ｙ_C－ｙ_O)²＋(ｈ_C－ｈ_O)²

が成立します。

これは３つの未知数ｘ_O,ｙ_O,ｈ_Oに対する３つの連立２次方程式ですから解くことが可能です。

つまり,正確に震源までの距離を観測できる点が３個あれば原理的には震源の位置がわかります。これは３点法と呼ばれるものですね。

普通,震源の高さ:ｈ_Oは海抜で測って負の数であり,震源は地中または海中にあります。

大森公式:ｄ_A＝ｋＴ_SPは地震に対する公式で大森係数と呼ばれる

ｋ＝Ｖ_PＶ_S/(Ｖ_P－Ｖ_S)は８(ｋｍ/ｓ)程度の値です。

しかし,パラメータｋを色々と変えると,雷における音と光のように速さの異なる２つの信号が届く現象では全く同じ公式が成立します。

ただし,雷の場合は光速がｃ～30万km/ｓ,空気中の音速がｖ～340ｍ/ｓでｃ＞＞ｖなのでｋ＝ｃｖ/(ｃ－ｖ)～ｃｖ/ｃ＝ｖですから,事実上ｄ_A～ｖＴなのでほとんど光速は関係無しです。(Ｔは光が見えてから音が聞こえるまでの時間です。)

例えば,ピカッと光ってから雷の音が聴こえるまでの時間が１秒なら,自分から340ｍ程度離れたところ(空中かもしれません)で"放電＝落雷"があったと推測されます。

PS：＞私のマノン・レスコーたちへ　

　病気とか事故に会っていなければどこで誰と遊んでいてもいいけれど,連絡つかない状況だととにかく心配です。大丈夫かなあ。。。。

　("マゾ, ロリコン, 準デブ専, メガネフェチ＝変態"のＴＯＳＨＩより。。)

(Ｔ.ウッズも聖人君子ではなく適当にストレスを解消しているらしいので安心しました。あれだけの収入を得る能力があれば大勢の家族を養えますね。

　性欲のある健康な男だし,大有名人で品行方正キャラが売り物なら日本だと風俗通いもままならないので女遊びは囲い込むしか仕方ないとして,もしもプロテスタントならマドンナの養子のように大勢を扶養するのはむしろ義務かな？)

PS2:ドクちゃんって日本が好きなんだ。。？？

　国ではなくて故郷(くに)が好きなんでしょ？

　http://news.aimu-net.com/read.cgi/wildplus/1257045443/

PS3:「陳情政治」からいきなり「非陳情政治」へと革命的に移行するのではなく,途中で１人(小沢さん)に権力が集中するシステムを経る。。

　共産革命でも過渡期ではプロレタリア独裁(現実には中国では毛沢東,ロシアではレーニンの１人独裁)です。それはそれでマヌーバとしてありなのでしょうが独裁者がプラトン的(プラトニック)に正しい哲人でいるうちはいいけど,うまく移行できないと"元の木阿弥"か,もっと悪くなります。

　小沢氏が中国に大挙して出かけていって米国にブラフかけて「ポッポ政権」を裏で後押ししてるのかも知れないけど,うまくいくのかねえ。基地問題。。

　朝っぱらから女子アナが「吉宗は何将軍？」と聞かれて「暴れん坊将軍？」と答えていました。それでもいいんだろうけど,ニュース番組なので一応「米(コメ)将軍だろ？」とツッコミたくなりました。。。アハハ

← クリックして投票してくだかなりさい。(１クリック＝１投票です。１人１日１投票しかできません。クリックすると「人気blogランキング」に跳びます。）

  ← クリックして投票してください。(ブログ村科学ブログランキング投票です。１クリック＝１投票です。１人１日１投票しかできません。クリックするとブログ村の人気ランキング一覧のホ－ムページームに跳びます。）

ヤーマン　プラチナゲルマローラー　1日3分コロコロエステ！ローラー型プラチナ配合美顔器  

　

お売りください。ブックオフオンラインのインターネット買取 展開へ！ ▼コミック　尾田栄一郎 「ONE PIECE（52）」 ▼コミック 「ONE PIECE」をオトナ買い

三国志特集　▼コミック 横山光輝 「三国志全巻セット」 「三国志（文庫版）全巻セット」   「三国志（ワイド版）全巻セット」 　▼書籍 「三国志」/吉川英治   「三国志」/北方謙三   「三国志」/宮城谷昌光

 

定量的地震学６

　つなぎで地震学の続きを書きます。

前回は内部面上での応力,または変位(歪み)の不連続性で表現される地震源(面源)の体源による等価物の存在を証明しました。

そして,断層面を横切る剪断作用に対する体積力の等価物として次のような表現が得られるという結論で終わりました。 

まず,応力の不連続性[Ｔ]がゼロでない場合には,その寄与－∫_-∞^∞ｄτ∫_Σ[Ｔ_p(ｕ(ξ,τ),ν)]Ｇ_np(ξ,ｔ－τ;ｘ,0)}ｄΣ_ξは,－∫_-∞^∞ｄτ∫_V(∫_Σ{[Ｔ_p(ｕ(ξ,τ),ν)]δ³(η－ξ)ｄΣ_ξ}Ｇ_np(ξ,ｔ－τ;ｘ,0))ｄＶ_ηで表現できます。

それ故,Σ上の応力不連続性の体積力等価物はｆ^[T](η,τ)≡－∫_Σ[Ｔ_p(ｕ(ξ,τ),ν)]δ³(η－ξ)ｄΣ_ξで与えられます。

一方,変位の不連続性による寄与∫_-∞^∞ｄτ∫_Σ[ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_j{∂Ｇ_np(ξ,τ－ｔ;ｘ,0)/∂ξ_q}は,－∫_-∞^∞ｄτ∫_V (∫_Σ[ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_j{∂δ³(η－ξ)/∂η_q}ｄΣ_ξ)Ｇ_np(ξ,ｔ－τ;ｘ,0)ｄＶ_ηで表現されます。

そこで,Σ上の変位不連続性の体積力等価物はｆ_p^[ｕ](η,τ)≡－∫_Σ[ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_j{∂δ³(η－ξ)/∂η_q}ｄΣ_ξで与えられると考えられます。

ただし,鍵括弧[ ]の量はΣ上の各点での不連続性を表わす量です。

すなわち,[Ｔ(ξ,τ)]≡Ｔ(ξ,τ)|_Σ+－Ｔ(ξ,τ)|_Σ-,かつ[ｕ_i(ξ,τ)]≡ｕ_i(ξ,τ)|_Σ+－ｕ_i(ξ,τ)|_Σ-etc.と定義されています。

という内容のことを書きました。

こうすれば応力,および変位の不連続性はそれぞれ－∫_-∞^∞ｄτ∫_Σ[Ｔ_p(ｕ(ξ,τ),ν)]Ｇ_np(ξ,ｔ－τ;ｘ,0)}ｄΣ_ξ＝∫_-∞^∞ｄτ∫_V[ｆ_p^[T](η,τ)Ｇ_np(ｘ,ｔ－τ;η,0)]ｄＶ_η,および∫_-∞^∞ｄτ∫_Σ[ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_j{∂Ｇ_np(ξ,τ－ｔ;ｘ,0)/∂ξ_q}＝∫_-∞^∞ｄτ∫_V[ｆ_p^[ｕ](η,τ)Ｇ_np(ｘ,ｔ－τ;η,0)]ｄＶ_ηと体源で表現できるわけです。

まず,これらについて前回,少し書き残したことを補足することから始めます。

　上記,変位不連続性の表現の被積分関数は,各ｐに対して,添字ｉ,ｊ,ｑのそれぞれが異なる27個の項を含んでいますが,以下では媒質に対称性があって２個か３個を除く全ての項がゼロであるような重要な例を取り上げる予定です。

　しかし,取り合えず上記の恒等式自体は一般的な非均質,非等方媒質に対して成立する式です。

　これらはＶ内の各点における体源応力が断層面の上の弾性媒質の性質のみに依存する表現になっていることは注目すべきことです。

そして,Ｖ内の断層運動は内部プロセス(内力のみの系)なので,全運動量と全角運動量は保存されなければなりません。

すなわち,全てのτに対し∫_Vｆ^[ｕ](η,τ)ｄＶ_η＝0,かつ全てのτとη₀に対し∫_V[(η－η₀)×ｆ^[ｕ](η,τ)ｄＶ_η＝0 であるべきです。

実際,ｆ_p^[ｕ](η,τ)≡－∫_Σ[ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_j{∂δ³(η－ξ)/∂η_q}ｄΣ_ξですから,ガウスの積分定理によって∫_Vｆ_p^[ｕ](η,τ)ｄＶ_η＝－∫_VｄＶ_η∫_ΣｄΣ_ξ{[([ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_j{∂δ³(η－ξ)/∂η_q}＝∫_ΣｄΣ_ξ∫_ＳextｄＳ_η[ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_jδ³(η－ξ)です。

ところがＳ_extとΣは全く共通点がないので,決してη＝ξとは成り得ず,最右辺は確かに消えます。

一方,∫_V[(η－η₀)×ｆ^[ｕ](η,τ)ｄＶ_ηの第ｍ成分は∫_Vε_mnp(η_n－η_0n)ｆ_p^[ｕ](η,τ)ｄＶ_η＝－∫_ΣｄΣ_ξ[ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_j∫_VｄＶ_η(ε_mnp(η_n－η_0n){∂δ³(η－ξ)/∂η_q})＝∫_Σ(ε_mqpＣ_ijpq(ξ)ν_j[ｕ_i(ξ,τ)])ｄΣ_ξで与えられます。

ところがＣ_ijpq(ξ)＝Ｃ_ijqp(ξ),かつε_mqp＝ε_mqpなので,これの右辺もゼロになります。

こうして体積力の表現が内力である条件を満たす合理的なものであることを示すことができました。これが前回書き残した補足事項です。 

さて,野外の不連続性に等価な体積力の簡単な例として,丁度１点と１方向に加えられた体積力のケースを考えます。

以前の記事「定量的地震学１,３」では次のように瞬間体積力を与えました。

“ｘ＝ξに対する１つの個別粒子に時刻ｔ＝τにｎ方向に瞬時的に加えられる１つの体積力ｆ(ｘ,ｔ)があれば,その成分ｆ_i(ｘ,ｔ)は空間位置を与えるのに３次元のディラック(Dirac)のデルタ関数,衝撃の時刻を与えるのに１次元のデルタ関数に比例してｆ_i(ｘ,ｔ)＝Ａδ³(ｘ－ξ)δ(ｔ－τ)δ_inと表現されます。

ただしＡは衝撃の強さを与える定数です。ｆ_i,δ³(ｘ－ξ),δ(ｔ－τ)の次元がそれぞれ[力/体積]＝ＭＬＴ^-2/Ｌ³,1/Ｌ³,1/Ｔであることに着目するとδ_inは無次元なので,衝撃の強さＡは正しく,”衝撃＝力積”の物理的次元を有することがわかります。”

これは,応力成分の不連続性とみなすこともできます。

ここでは,体積力の応力の不連続性との等価性を見るために,ｘ₃を深さ方向とし,ｘ₃＝ｈの１点(0,0,ｈ)に加わる体積力ですがδ(ｔ－τ)に比例する瞬時的な力ではなく,τ＝0から定常的にｘ₃＝ｈの１点(0,0,ｈ)にＦの大きさの体積力とします。

すなわち,ｆ(η,τ)＝Ｆδ(η₁)δ(η₂)δ(η₃－ｈ)θ(τ),Ｆ≡(0,0,Ｆ)とします。ここでθ(τ)はHeaviside関数(階段関数)です。

これは,平面ξ₃＝ｈ上の１点を横切る応力の不連続性:[Ｔ(ξ,τ)]＝Ｔ(ξ,τ)|_{(ξ1,ξ2,ｈ+)}－Ｔ(ξ,τ)|_{(ξ1,ξ2,ｈ-)}＝－Ｆδ(ξ₁)δ(ξ₂)θ(τ)と同一視できることがわかります。

これは,先に得られた表現:ｆ^[T](η,τ)≡－∫_Σ{[Ｔ_p(ｕ(ξ,τ),ν)]δ³(η－ξ)ｄΣ_ξの右辺において,平面ξ₃＝ｈをΣとして,[Ｔ_p(ｕ(ξ,τ),ν)]＝－Ｆδ(ξ₁)δ(ξ₂)θ(τ)を代入すれば,ｆ^[T](η,τ)＝Ｆδ(η₁)δ(η₂)δ(η₃－ｈ)θ(τ)が確かに得られることからわかります。

断層スリップ(地滑り)によって引き起こされる地震波は,モーメントを帳消しにするようなある断層上の力の分布によって引き起こされる地震波と同じです。

問題としている断層スリップに対し,この力の分布は一意的ではありませんが,等方的な媒質内ではそれを常に複数の偶力源の表面分布として選択可能です。

この結論は,地震が単一の偶力源か複数の偶力源のいずれによってモデル化されるかという疑問に関連して長期間論じられた論争の見地からは皮肉な結論です。

単一の偶力源を擁護する人々は地震は断層上のスリップが原因であると強く信じていました。そして彼らは,これが単一の偶力源,すなわち断層の相対する側の運動に対応する２つの力に等価であると直線的に信じました。

しかし,弾性力学では経験上直線的アプローチは危険なことが多いようです。

複数の偶力源を擁護する人々のあるものは,地震が既存の剪断応力下での体積的崩壊であるに違いないと思っていました。 

今や,複数の偶力源に等価であると認識されている近年の断層理論は,莫大な距離で観測された放射パターンによる支持のみならず,震源に非常に近いところで得られたデータによる強い支持を得ています。 

断層Σは平面ξ₃＝0 上にあるものとします。すると,スリップ[ｕ]はξ₃方向には成分を持たず,ベクトル[ｕ]は平面Σに平行と考えられます。 

さらに,平面ξ₃＝0 上でのスリップ[ｕ]の方向をξ₁方向とします。すると,[ｕ₂]＝[ｕ₃]＝0 でありν₁＝ν₂＝0 ,ν₃＝1 です。 

そこで,等価体積力の表現:ｆ_p^[ｕ](η,τ)＝－∫_Σ[ｕ_i(ξ,τ)]Ｃ_ijpq(ξ)ν_j{∂δ³(η－ξ)/∂η_q}ｄΣ_ξは,ｆ_p^[ｕ](η,τ)＝－∫_Σ[ｕ₁(ξ,τ)]Ｃ_13pq(ξ){∂δ³(η－ξ)/∂η_q}ｄξ₁ｄξ₂となります。 　

既述のように,不均質でも等方的な物質では弾性係数Ｃ_ijpqはＣ_ijpq＝λδ_ijδ_pq＋μ(δ_ipδ_jq＋δ_iqδ_jp)なる形をしています。独立定数λ,μはLame(ラメ)の定数と呼ばれます。

よって,今の場合の被積分関数のＣ_13pq(ξ)は,Ｃ₁₃₁₃(ξ)＝Ｃ₁₃₃₁(ξ)＝μ(ξ)を除いて全てゼロです。

それ故,ｆ₁^[ｕ](η,τ)＝－∫_Σμ(ξ)[ｕ₁(ξ,τ)]δ(η₁－ξ₁)δ(η₂－ξ₂){∂δ(η₃)/∂η₃}ｄξ₁ｄξ₂,ｆ₂^[ｕ](η,τ)＝0,ｆ₃^[ｕ](η,τ)＝－∫_Σμ(ξ)[ｕ₁(ξ,τ)]{∂δ(η₁－ξ₁)/∂η₁}δ(η₂－ξ₂)δ(η₃)ｄξ₁ｄξ₂となります。

まず,ｆ₁^[ｕ](η,τ)＝－∫_Σμ(ξ)[ｕ₁(ξ,τ)]δ(η₁－ξ₁)δ(η₂－ξ₂){∂δ(η₃)/∂η₃}ｄξ₁ｄξ₂を見ると,これは±η₁方向内の力でη₂方向に腕がありη₃方向にモーメントを持つΣにわたる単偶力源を示しているとわかります。

実際,積分を実行すると,ｆ₁^[ｕ](η,τ)＝－μ(η)[ｕ₁(η,τ)]{∂δ(η₃)/∂η₃}です。

これは平面η₃＝0⁺上に寄与する点力とη₃＝0^-上に寄与する反対向きの点力のペアであると考えられ,補足事項の最初の式で示したように,内力なのでｆ₁^[ｕ](η,τ)の合力成分は消えます。

一方,力のモーメントはこの力の成分だけでは消えず,この成分によるη₂軸の回りのモーメントは∫_Vη₃ｆ₁^[ｕ]ｄＶ＝－∫_Vη₃μ(η)[ｕ₁(η,τ)]{∂δ(η₃)/∂η₃}ｄη₁ｄη₂ｄη₃＝∫_Σμ(ξ)[ｕ₁(ξ,τ)]ｄΣ_ξとなります。

Σにわたるスリップを平均すると,＜ｕ(τ)＞≡∫_Σ[ｕ₁(ξ,τ)]ｄΣ_ξ/Ｓです。ただし,Ｓ≡∫_ΣｄΣ_ξは断層の全面積です。

そこで,もしも断層域が均質,つまりΣの上でμ(ξ)＝μ(一定)なら,η₂軸の回りのｆ₁^[ｕ](η,τ)による全モーメントは,単にη₂軸の正の向きに∫_Vη₃ｆ₁^[ｕ]ｄＶ＝∫_Σμ(ξ)[ｕ₁(ξ,τ)]ｄΣ_ξ＝μＳ＜ｕ(τ)＞で与えられます。

体積力には,また３軸成分:ｆ₃^[ｕ](η,τ)＝－∫_Σμ(ξ)[ｕ₁(ξ,τ)]{∂δ(η₁－ξ₁)/∂η₁}δ(η₂－ξ₂)δ(η₃)ｄξ₁ｄξ₂＝－∂{μ(η)[ｕ₁(η,τ)]/∂η₁}δ(η₃)があります。

この成分によるη₂軸の回りの力のモーメントは－∫_Vη₁ｆ₃^[ｕ]ｄＶ＝∫_Vη₁∂{μ(η)[ｕ₁(η,τ)]/∂η₁}δ(η₃)ｄη₁ｄη₂ｄη₃＝－∫_Σμ(ξ)[ｕ₁(ξ,τ)]ｄΣ_ξです。

これについても,もしも媒質が均質でΣの上でμ(ξ)＝μ(一定)なら,－μＳ＜ｕ(τ)＞です。

よって,均質,不均質いずれにしても全モーメントはゼロです。これも既に,補足事項として一般的な形で証明した事実に一致しています。

途中ですが今日はここで終わります。 

参考文献:K.Aki,& P.G.Richards 「Quantitative Seismology(Theory and Method)」

　← クリックして投票してください。(１クリック＝１投票です。１人１日１投票しかできません。クリックすると「人気blogランキング」に跳びます。）

← クリックして投票してください。(ブログ村科学ブログランキング投票です。１クリック＝１投票です。１人１日１投票しかできません。クリックするとブログ村の人気ランキング一覧のホ－ムページームに跳びます。）

ヤーマン　プラチナゲルマローラー　1日3分コロコロエステ！ローラー型プラチナ配合美顔器  

　

お売りください。ブックオフオンラインのインターネット買取 展開へ！ ▼コミック　尾田栄一郎 「ONE PIECE（52）」 ▼コミック 「ONE PIECE」をオトナ買い  

三国志特集　▼コミック 横山光輝 「三国志全巻セット」 「三国志（文庫版）全巻セット」   「三国志（ワイド版）全巻セット」 　▼書籍 「三国志」/吉川英治   「三国志」/北方謙三   「三国志」/宮城谷昌光

 

    予約】セックス・アンド・ザ・シティ　ザ・ムービー　コレクターズ・エディション＜限定盤＞