Radiation Magnetohydrodynamics 輻射磁流體力學（一）：輻射輸運方程與輻射磁流體控制方程的建立

08-16

來自專欄計算流體力學——從原理到代碼37 人贊了文章

一、動機

磁流體，往往是描述等離子體的宏觀運動，應用廣泛。等離子體的運動，往往具有高速高溫高導的性質，所以近三十年來迫切的希望能講現實中大量存在的熱現象，特別是熱輻射，加入到磁流體的數值模擬模型中去。

我們知道，熱量通過三種方式傳遞：傳導、對流、輻射。那麼，當我們面對數值模擬時，三大熱量輸運會變成怎麼樣的模型加入原來的模型呢：

熱傳導：內能的耗散不是最高效的。但是對於各向異性(anisotropic)的相對論情形的實現可以容易。
熱對流：在流體方程組中自然包含了進去，但是需要很高的精確度
熱輻射：一般由電子、光子及中微子等粒子作為媒介來實現熱能的轉移，但是這些粒子的行為更像是耗散（diffusion）或者流出（streaming）。往往，這些粒子除了運輸能量，還帶走了一部分動量。更為嚴肅的問題是粒子運動的時標遠遠小於流體的運動。——總而言之，輻射是非常複雜且困難的行為。

二、輻射的建模

下面我們以光輻射為例，即考慮光子的行為：

光子之間不發生碰撞collisionless
光子以光速行進

具體闡述collisionless：光子的性質不是局部的，local的，我們需要考慮全局所以光子的性質，不只是光子在空間中的分布，還有光子方向的分布。

由於，能量是守恆量，所以我們還是從能量出發，輻射的能量和什麼有關係？

肯定，第一是與空間位置和方向有關，哪裡光子分布密度大，哪裡光子照射密集，哪裡輻射能高
有空間就有時間，與t有關
其次，是和頻率有關，不同頻率攜帶的能量不一樣，頻率越高，能量越高，定量來說就是 $E = h v = hbar omega$ ，所以與頻率 $u$ 也有關
最後，光子沒有電荷、沒有（rest）質量，所以都不影響。

我們使用光強intensity來描述輻射的大小，那這個量受什麼更加根本的量控制呢：

光的三維空間分布 $f r$
三維方向分布 $f n$
頻率 $u$
時間t

強度是能量的微分，這個大家肯定沒問題

$d E _ { u } = I _ { u } ( mathbf { r } ,{ mathbf { n } } , u,t ) cos heta d u d a d Omega d t$

這裡符號，是取光強首字母作為符號，記為 $mathbf{I}( mathbf{r},mathbf{n}, u,t) ~~[ ext{erg}~~ ext{s}^{-1} ext{cm}^{-2} ext{Hz}^{-1} ext{sr}^{-1} ]$ ，後面是單位，光強由上面的量決定這件事情說明，光強是一個，六個空間自由度+1個時間自由度的量。

三、輻射輸運方程

最簡單的模型，在某個時間微元上，光子行為限制在某個固定方向上，光子進行自由的流動，那對於這個方向的光產生的輻射的變化，就受到光子流量的變化的直接影響：

文字來說：

單位時間裡，輻射強度的變化——（ $frac { partial I _ { u } } { partial t }$ ）等於該方向上，輻射強度流的變化 $- c mathbf { n } abla I _ { u }$ ：

$frac{partial I_ u}{partial t} + c extbf{n} abla{I_ u} = 0$

下一步，模型加入新產生的輻射，和被吸收的輻射，所以輻射強度的變化，除了流入流出的凈變化，還有新產生和被吸收的，我們寫在右端：

$frac{partial I_ u}{partial t} + c extbf{n} abla{I_ u} = c(J_ u-k_ u I_ u)$

右端：新產生的輻射 $J _ { u }$ ，吸收係數 $k_ u$

下面，加入從該方向散射，係數記為 $sigma_ u$ 出去的，散射意味著該方向的能量耗散出去了

$frac { partial I _ { u } } { cpartial t } + mathbf { n } abla I _ { u } = J _ { u } - (k _ { u } +sigma_ u )I _ { u } +sigma_ u ar{I}_ u$

注意，右邊帶一橫的流量 $ar{I}_ u$ ，意思為總的平均的光強，就是旁邊方向散射回來的意思。

四、輻射強度矩 Radiation Intensity Moment

我們考慮單位方向向量 $extbf{n}$ 上的方位角 $Omega$ 的各級平均，即如下定義1、2、3階 of 輻射 $I_ u$

$J_ u =frac{1}{4pi} int I_ u ext{d}Omega\ mathbf{H}_ u =frac{1}{4pi} int I_ u extbf{n} ext{d}Omega\ K_ u =frac{1}{4pi} int I_ u extbf{n} extbf{n} ext{d}Omega$

對於能量、流量、壓力在單位方向向量 $extbf{n}$ 上的方位角 $Omega$ 的平均，可以由單位方向向量 $extbf{n}$ 上的方位角 $Omega$ 的輻射強度 $I_ u$ 的1，2，3級平均得到

$E_ u =frac{1}{c} int I_ u ext{d}Omega\ mathbf{F}_ u =int I_ u extbf{n} ext{d}Omega\ P_ u =frac{1}{c} int I_ u extbf{n} extbf{n} ext{d}Omega$

考慮輻射輸運方程： $frac { partial I _ { u } } { cpartial t } + mathbf { n } abla I _ { u } = J _ { u } - (k _ { u } +sigma_ u )I _ { u } +sigma_ u ar{I}_ u$

如果假設光源光淵在各個方向角都是相同的，於是乎拆成矩的方差，得到：

$frac { partial E _ { u } } { partial t } + abla mathbf{F}_ { u } = c(s_ u - sigma_ u E_ u)\ frac { partial mathbf{F}_ { u } } { cpartial t } + c^2 abla P _ { u } = -c sigma_ u mathbf{F}_ u$

五、（頻率平均）灰度估計 Gray Approximation

因為前面，我們都是固定頻率 $u$ 來分析時空與輻射的關係，現在我們來考慮輻射，我們講輻射平均掉，這樣的估計，叫做灰度估計 Gray Approximation，因為你把顏色混合起來就是灰色= =

(這裡推導我還沒搞懂，希望有人明白的能指點一下)

我們考慮黑體熱輻射， $I_ u = B_ u$ B就是普朗克公式

$frac { partial I } { cpartial t } + mathbf { n } abla I = sigma_alpha [B(T)-I]\ frac { partial E } { partial t } + abla mathbf{F} = c(sigma_alpha aT^4- sigma_E E)\ frac { partial mathbf{F} } { cpartial t } + c^2 abla P = -c sigma_F mathbf{F}$

這裡a=7.5657*10^-15

六、RMHD的控制方程

規定： $ho$ 密度， u 是速度，p 是熱壓力， $σ_R$ 是Rosseland平均不透明度， $F_r$ 是輻射通量，E是流體總能量 $E =ρε+ 1/2ρu^2$ （ε是流體自身特定的內能） $σ_P$ 是普朗克不透明度， $B = B(T)$ 是普朗克函數， $E_r$ 是輻射能量， $mathbb { P } _ { r }$ 是輻射壓力。 $σ_RF_r / c$ 對物質產生輻射力，

質量守恆： $partial _ { t } ho + abla [ ho u ] quad = 0$

牛頓第二定律： $partial _ { t } ho u + abla [ ho u otimes u + P I ] = sigma _ { R } F _ { r } / c$

氣體能量守恆： $partial _ { t } E + abla [ u ( E + P ) ] = sigma _ { R } F _ { r } / c cdot u - sigma _ { P } left( 4pi B - c E _ { r } ight)$

輻射能量守恆： $partial _ { t } E _ { r } + abla left[ u E _ { r } ight] = - abla cdot F _ { r } - mathbb { P } _ { r } : abla u + sigma _ { P } left( 4pi B - C E _ { ext{r} } ight)$

輻射流量守恆： $partial _ { t } F _ { ext{r} } + abla left[ u F _ { r } ight] = - c ^ { 2} abla cdot mathbb { P } _ { ext{r} } - left( F _ { ext{r} } cdot abla ight) u - sigma _ { ext{F} } c F _ { ext{r} }$

（寫得粗略了點，我會不斷修改的）