Engineering Sciences Report, Kyushu University
Vol.28, No.1 pp.17–22, JUN. 2006

九州大学大学院総合理工学報告第 28 巻 第 1 号 17-22 頁 平成 18 年 6 月

局所熱非平衡状態を考慮した充填層伝熱モデル
田中 雄一郎*1,†・横峯 健彦*2 ・江原 真司*2 ・清水 昭比古*2(平成18年4月28日 受理)

Model for Heat Transfer in Packed Bed
under Local Thermal Non‐equilibrium Conditions
Yuichiro TANAKA, Takehiko YOKOMINE, Shinji EBARA, Akihiko SHIMIZU
E-mail of corresponding author: tanaka@aees.kyushu-u.ac.jp†The focus of this work is the forced convective heat transfer in a packed bed under local thermalnon-equilibrium conditions. It is considered that local thermal non-equilibrium conditions, which are
made by the temperature difference between pebbles and their surrounding continuum media, cause
the anisotropic effective thermal conductivity of packed bed by means of thermal dispersion. In thisstudy, new model for the thermal dispersion is established with the same procedure to constructturbulence heat flux model. This model represents that the temperature difference between thephases directly influences thermal dispersion under local thermal non-equilibrium conditions.
Key words: Packed bed, Convective heat transfer, Thermal dispersion, Local volume averaging

method, Local thermal non-equilibrium condition

1. 緒

言

熱分散とは、微視的な速度場の変動に起因した熱混

近年、多孔質の応用範囲は広がりを見せており、様々合を表し、従来の研究から、充填層内の有効熱伝導率な分野で研究が行われている。しかし、多孔質内で起は分子熱拡散と熱分散の相乗効果として表されること、こる種々の現象の解析は何れも容易ではなく、特に、流れが速くなると充填層内の伝熱は熱分散に支配され充填層内の熱流動現象はその複雑性からまだ未知な部ること、また、熱分散が有効熱伝導率の非等方性に影分が多く、新たな研究課題として注目されている。

響を及ぼしていることなどが分かっている1)。しかし、充填層内の熱流動解析は、その流れ場の解析に基づ

これらの研究は流体と固体粒子の温度を平衡とみなすいて行われる。充填層内の流体は、粒子の存在による局所熱平衡状態におけるものが主であり、現在、充填影響を受け複雑な流路を持ち、その流路断面積は流れ層の応用が検討されている触媒、熱交換器、燃焼器な方向に変化し、また流路間には流れの干渉も考えられどの開発においては、従来の局所熱平衡状態でのモデる。機器設計等を考えた場合、このような個々の流路ルに替わる新しいモデルが必要となるであろう。何故について複雑現象を逐一追うことは現実的とは言えななら、局所熱非平衡状態における熱分散は流体と粒子い。よって、充填層内の流れ場の理論解析には、ある間の熱伝達による影響を受け、局所熱平衡状態でのそ種の平均操作が有効となる。その主流と言えるのが、れとは異なる挙動を示すと考えられるからである。本

物理量を充填層内の検査体積内で空間平均する局所体研究では、局所熱非平衡状態における充填層伝熱モデ積平均法である。ここ数十年間での理論解析の多くがルの構築の為に必要となる、熱分散流束のモデリング目的としたのは、この局所体積平均操作により見落とを行った。

してしまう空隙スケールの現象が、巨視的な系に及ぼす影響を評価することにあったと言える。特に、微視2. 理 論 背 景的座標系における流体のエネルギー方程式を局所体積2.1 局所体積平均平均した際に現れる熱分散項の評価に関する研究は、局所体積平均とは、固体粒子と空隙部を流れる流体充填層内の熱流動の解明に大きく寄与してきた。

により占められている充填層内の検査体積V（Fig.1）*1 先端エネルギー理工学専攻 修士課程
*2 エネルギー理工学部門座標空間（x空間）の物理量に変換するもので、この操において、微視的座標空間（ξ空間）の物理量を巨視的

18

局所熱非平衡状態を考慮した充填層伝熱モデル

作により多孔質内の複雑な輸送現象の記述を簡易化し、やWhitaker4)等により行われ、次の関係が示されていその解析を効率よく行うことができる。検査体積Vがる。

以下の関係式を満たす時、局所体積平均は有効とな

1
ψ dA
V ∫A fs(7)(8)表す。この条件を満たす検査体積Vにおいて、流体に∂ψ∂ψ
1=− ∫ ψ w idA∂t∂t
V A fs

関するある物理量ψの局所体積平均値は次式で定義さAfsは粒子と流体の界面における表面積で、dAは流体側れる。

から見て外向き方向の面積要素ベクトルである。wは∇ψ = ∇ ψ +る2)。
dl(1)

L

dは粒子径、lは検査体積Vの、Lは充填層の代表寸法を

ψ =

1
ψ dV
V ∫V f(2)

流体と固体粒子の界面の速度を表す。
2.2 巨視的支配方程式

Vf は流体が占める体積である。また、その局所相平均値は

ψf=

微 視 的 座標 系に お け る支 配方 程 式 （連 続の 式 、Navier-Stokes方程式、各相のエネルギー方程式）は、

1
ψ dV
V f ∫V f(3)∂u j∂ξ jで定義され、充填層内の空隙率をε =Vf /Vとおくと、(2)、1 ∂p∂ui∂∂ ⎧⎪ ⎛ ∂ui ∂u j ⎞ ⎫⎪+++u j ui ) = −⎟⎬(⎨ν ⎜ρ ∂ξ i ∂ξ j ⎩⎪ ⎜⎝ ∂ξ j ∂ξ i ⎟⎠ ⎭⎪∂t ∂ξ j(3)式で定義される平均値の間には次の関係式が成り立つ。

ψ =ε ψf⎧ ∂T f⎫pfj⎪⎩

fjj(5)⎪⎭j⎛sp s+ ψˆ∂T f ⎞⎝

fj(11)⎠⎞( ρ c ) ∂∂Tt = ∂∂ξ ⎜⎜ k ∂∂ξT ⎟⎟の和として表す。
f

⎛(10)

( ρ c ) ⎪⎨ ∂t + ∂∂ξ ( u T )⎪⎬ = ∂∂ξ ⎜⎜ k ∂ξ ⎟⎟(4)

更に、ψをその局所相平均値とその平均値からの偏差

ψ= ψ(9)=0⎝ssj(12)⎠ρは密度、νは動粘性係数、cpは定圧比熱、kは熱伝導率を、添え字fとsはそれぞれ流体と固体粒子に関する物

ψˆf

(6)

=0

微視的支配方程式から巨視的支配方程式を導くには、ξ座標におけるψの勾配や時間微分を局所体積平均し、x座標の関数で表す必要がある。この議論はSlattery3)

理量を表す。簡単の為に、全ての物性値は定数とみなし、空隙率ε は一定、更に充填層構造が変形しない（w= 0）と仮定する。これらの仮定の下、(9)から(12)式を検査体積Vに関して局所体積平均し、(5)から(8)式を用いて整理すると、以下の巨視的支配方程式を導くことができる。

∂ uj

ξ3

∂x j
D uix3

Dtf

ξ1s

ξ2

the control volume Vx+

D Tfx2

Fig.1 Microscopic and Macroscopic coordinatesf(13)=0ff⎧ ⎛∂ u j ⎞ ⎫⎪∂ ⎪ ⎜ ∂ ui
⎟⎬++
⎨ν f
∂x j ⎪ ⎜ ∂x j
∂xi ⎟ ⎪⎠⎭⎩ ⎝⎧⎪ p f⎛ ∂ui ∂u j ⎞ ⎪⎫f1∂+uˆiuˆ j (14)⎟⎟ ⎬ n j dA −⎨− δ ij + ν f ⎜⎜ξξ∂∂∂V f ∫Afs ⎩⎪ ρ fxi ⎠⎭j
⎪
⎝ j

Dt

x1

1 ∂ pf=−ρ f ∂xiffff⎧⎫⎫⎪∂ ⎪ ∂ T f ⎪ ∂ ⎪⎧ 1
=α f T f n j dA ⎬⎨α f⎬+⎨∂x j ⎪∂x j ⎪ ∂x j ⎩⎪V f ∫A fs⎭⎪⎩⎭+∂T1∂uˆ jTˆf
α f f n j dA −V f ∫A fs∂ξ j∂x jf(15)平成 18 年∂ Ts∂ts=

第 28 巻九州大学大学院総合理工学報告∂ ⎪⎧ ∂ Ts⎨α s∂x j ⎩⎪∂x j+s⎪⎫
⎬
⎭⎪

第1号

19(13)から(16)式の巨視的方程式を閉じる為には、(14)
f

⎫ 1
∂ ⎧1
∂T
α sTs n j dA⎬ + ∫ α s s n j dA (16)
⎨A fs
ξjV
∂x j ⎩Vs ∫Afs∂s⎭njは流体から粒子を指す方向の単位法線ベクトル、αは式と(15)式の右辺第四項、 − ∂ uˆi uˆ j
∂x j と熱分散の寄fˆˆ与を表す − ∂ u T∂x をモデル化すれば良い。本研
jfj

究では、次節から、乱流熱流束のモデリングと同様のf手順で熱分散流束 uˆ Tˆに関する輸送方程式を導出if

し、熱分散流束のモデルを構築していく。

温度伝導率を表す。上式においてf ∂ψDψ ∂ψ=+ uj
Dt∂t∂x j(17)3. 熱分散流束のモデリング
3.1 熱分散輸送方程式f

熱分散流束 uˆiTˆf

に関する輸送方程式を導く為に、流体の速度偏差 uˆ と温度偏差 Tˆ に関する輸送方程式を用いている。
界面における境界条件は、以下の通りである。

k f njifui = 0(18)を導く。(5)式を(9)、(10)、(11)式に代入し、それぞれT f = Ts(19)(13)、(14)、(15)式との差を取ると、∂T f∂ξ j= ks n j∂uˆ j∂Ts∂ξ j(13)から(16)式の巨視的支配方程式を解く為には、それぞれの面積分項をモデル化し、この方程式群を閉じる∂ uiDuˆi+ uˆ jDt∂x j必要がある。

⎧⎪ p f⎛ ∂ui ∂u j ⎞ ⎪⎫1+⎟⎟ ⎬ n j dA⎨− δ ij + ν f ⎜⎜V f ∫A fs ⎩⎪ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪DTˆfDt⎧⎪ pˆ f⎛ ∂uˆ ∂uˆ ⎞ ⎫⎪1=− δ ij + ν f ⎜ i + j ⎟ ⎬ n j dA (21)⎨∫⎜⎟V f Afs ⎩⎪ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪と表すことができ、この項が固体粒子の存在による流R2 = −∂x j−fuj⎞f(22)
⎟ ui
⎠∂Tˆf
∂ξ j=

∂ ⎛ ∂Tˆf ⎞⎜α f⎟ + R 2 (27)∂ξ j ⎜⎝ ∂ξ j ⎟⎠ujfujff⎞∂fuˆ j uˆi(28)⎟ ui +x∂⎠j∂T1∂uˆ jTˆf
α f f n j dA +V f ∫Afs∂x j∂x jDTˆf ˆ DuˆiDuˆiTˆf = uˆi+ Tf
Dt
Dt
Dt()

足し合わせると、次式を導くことができる。

寄与は他の項の寄与と比べて充分に小さくなると類推∂ ui
D(uˆiTˆf ) + Tˆf uˆ j
∂x jDt

検査体積 V における平均熱伝達率 h を用いると境界条+uˆi uˆ j

∂ Tf

件(20)式より、
f) (23)f1∂T1
α s s n j dA =hAfs ( T fVs ∫A fs ∂ξ jVs ( ρ c p )sfs− Ts ) (24)(30)fると仮定する。(15)、(16)式の右辺第三項の面積分は、∂T
1
1shAfs ( Ts − T f
α f f n j dA =∫AV f fs
∂ξ jV f ( ρcp )(29)

より、(26)式の両辺に Tˆf を、(27)の両辺に uˆi を掛けてるにつれ、(15)、(16)式における右辺第二項の面積分のされる。そこで本研究では、この項の寄与を無視でき

f

ここで、

Kは透過率、FはForchheimer係数である。
局所熱平衡状態における研究1)から、流れが速くな+ uˆ j⎞
∂ ⎛ 1
α f T f n j dS ⎟⎜⎟∂x j ⎜⎝ V f ∫Afs⎠⎧⎪ pˆ f⎛ ∂uˆi ∂uˆ j ⎞ ⎪⎫1+⎟⎟ ⎬ n j dA⎨− δ ij + ν f ⎜⎜V f ∫A fs ⎩⎪ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪ff

∂ Tf

⎛ν εR1 = ⎜ f + F ε 2⎝ K

ル化される5)。

uj+ uˆ j⎛ ∂uˆi ∂uˆ j ⎞ ⎫⎪∂ ⎧⎪ pˆ f+⎟⎟ ⎬ + R1 (26)⎨− δ ij + ν f ⎜⎜∂ξ j ⎩⎪ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪(26)、(27)式において、体抵抗の寄与を表わしていることが分かる。この項は一般にDarcy-Forchheimerの法則を用いて次式でモデ∂uˆi∂ξ j+ uˆ j=(14)式の右辺の面積分は(5)、(6)式より⎛ν f ε= −⎜+ Fε 2K
⎝

f(25)

=0∂ξ j(20)= Tˆf∂x jf+ Tˆf uˆ j∂Tˆ∂uˆi+ uˆi uˆ j f∂ξ j∂ξ j⎛ ∂uˆi ∂uˆ j ⎞ ⎪⎫∂ ⎪⎧ pˆ f+⎟⎟ ⎬⎨− δ ij + ν f ⎜⎜∂ξ j ⎩⎪ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪+uˆi∂ ⎛ ∂Tˆf ⎞ ˆ⎜α f⎟ + T f R1 + uˆi R 2 (31)∂ξ j ⎜⎝ ∂ξ j ⎟⎠20

局所熱非平衡状態を考慮した充填層伝熱モデル(31)式の局所体積平均を取ると、熱分散流束 uˆiTˆffに関する輸送方程式を導出できる。(7)、(8)式と境界条件(18)式より、

DuˆiTˆf
Dtf∂uˆi uˆ jTˆf∂x j+

f=f+ Tˆf uˆ j+ uˆi uˆ j

∂x jf

∂Tˆ
∂
+α f uˆi f∂x j∂ξ j∂ Tff∂ ui⎛ ∂uˆi ∂uˆ j ⎞ ⎫⎪∂ ⎪⎧ pˆ f
= Tˆf+⎟⎟ ⎬⎨− δ ij + ν f ⎜⎜∂ξ j ⎪⎩ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎪⎭ff

⎛ν f ε= −⎜+ Fε 2K⎝

f+ uˆi R 2uˆi R 2= Tˆf

f

= uˆi

fff

(32)

R1 = 0(33)R2 = 0(34)f+ α f uˆi

∂Tˆff

∂ξ jf

∂ ui
∂x j

f

− uˆiuˆ jf∂ Tf⎛ν f ε= −⎜+ Fε 2
⎝ Kφi

⎛ν f ε= −⎜+ Fε 2
⎝ Kf

∂x j

生成項(5)と(6)式より uiTˆf⎛ ∂uˆi ∂uˆ j ⎞ ⎫⎪∂ ⎧⎪ pˆ f+ Tˆf+⎟⎟ ⎬⎨− δ ij + ν f ⎜⎜∂ξ j ⎩⎪ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪ffujfuj⎞⎟ ui + φi (38)⎠f(39)=0⎛ ∂Tˆ ⎞1∂uˆ ∂Tˆfuˆi ⎜ α f f ⎟ n j dA − α f i∫
AV f fs ⎝⎜ ∂ξ j ⎠⎟∂ξ j ∂ξ j(B)項

ujfff

uj= uˆiTˆfff⎞ ˆ⎟ uiT f⎠f+ φiTˆff(40)が成立するので、(A)項を次式で近似し、そのモデルとする。

⎛ ∂uˆi ∂uˆ j ⎞ ⎫⎪∂ ⎧⎪ pˆ f+Tˆf⎨− δ ij + ν f ⎜⎜⎟⎟ ⎬∂ξ j ⎪⎩ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪(A)項+⎛ ∂uˆi ∂uˆ j ⎞ ⎪⎫∂ ⎪⎧ pˆ f+⎨− δ ij + ν f ⎜⎜⎟⎟ ⎬∂ξ j ⎪⎩ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪⎛ ∂uˆi ∂uˆ j ⎞ ⎫⎪∂ ⎧⎪ pˆ f+Tˆf⎨− δ ij + ν f ⎜⎜⎟⎟ ⎬∂ξ j ⎪⎩ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪⎫
⎪⎬
⎪⎭

拡散項

慣性項

= − Tˆf uˆ j

⎞
f(37)
⎟ ui
⎠

(38)式の両辺に Tˆf を掛け、その局所体積平均を取ると⎧∂ ⎪+⎨ uˆi uˆ jTˆf∂x j ⎪⎩

ffuj

φは積分定数であり次式を満たす。

以上より、
DuˆiTˆfDtfuj

よって、f(28)、(29)式より R1 と R2 はξに依存しないので、ff⎛ ∂uˆi ∂uˆ j ⎞ ⎫⎪∂ ⎧⎪ pˆ f+⎨− δ ij + ν f ⎜⎜⎟⎟ ⎬∂ξ j ⎪⎩ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪f+ Tˆf R1⎛ ∂uˆ ∂uˆ ⎞ ⎪⎫1⎪⎧ pˆ f− δ ij + ν f ⎜ i + j ⎟ ⎬ n j dA (36)⎨∫⎜ ∂x j ∂xi ⎟V f Afs ⎩⎪ ρ f⎝⎠ ⎭⎪(22)式より(36)式は次式のように変換できる。

∂x j⎛ ∂Tˆ ⎞1+ ∫ uˆi ⎜ α f f ⎟ n j dAAV f fs ⎜⎝ ∂ξ j ⎟⎠∂uˆ ∂Tˆf
− αf i∂ξ j ∂ξ jTˆf R1f⎛ ∂uˆi ∂uˆ j ⎞ ⎫⎪∂ ⎧⎪ pˆ f+⎨− δ ij + ν f ⎜⎜⎟⎟ ⎬∂ξ j ⎪⎩ ρ f⎝ ∂ξ j ∂ξ i ⎠ ⎭⎪ff(35)

散逸項

各項の物理的意味は乱流熱流束モデルから類推した6)。

⎛ν f ε≈ −co ⎜+ Fε 2K⎝ujfujf⎞ ˆ⎟ uˆiT f
⎠

f(41)

co はモデル定数である。
一方、(B)項は境界条件(18)式より

本研究では、拡散項と散逸項を無視して(A)と(B)項⎛ ∂Tˆ ⎞∂Tˆ
1f 1uˆi ⎜ α f f ⎟ n j dA = − ui
α f f n j dA (42)∫∫AAV f fs ⎜⎝ ∂ξ j ⎟⎠V f fs∂ξ jのモデル化について考える。まず、(A)項のモデル化に(23)式の関係より、(B)項は次式で表すことができる。

3.2 モデリングついて示す。グリーンの定理より

平成 18 年

九州大学大学院総合理工学報告⎛ ∂Tˆ ⎞1uˆi ⎜ α f f ⎟ n j dA∫
V f A fs ⎜⎝ ∂ξ j ⎟⎠第 28 巻第1号

が提案した代数熱分散モデルと一致する。以下、中山らが行った検証手順を追うことにより、局所熱平衡状=−態における本モデルの妥当性を示す。(46b)式において、

fui

V f ( ρcp ) fs

hAfs ( Tsf− Tf) (43)

局所熱平衡状態を仮定し流れと平行な成分に注目すると、上式より、局所熱平衡状態では(B)項は無視できるが、α xx =局所熱非平衡状態では無視できないことが分かる。よって、(B)項が局所熱平衡状態と局所熱非平衡状態の違いを陽的に表していると言える。

f= − Tˆf uˆ jff− uˆiuˆ j∂x j⎛ν f ε−co ⎜+ Fε 2
⎝ K+∂ Tff

∂ uiujuishAfs Ts

f− Tf) (44)(47)

150 (1 − ε )(48a)2ui = ( uD ,0,0 ) の系について行う。この時、∂ uif(48b)f∝(1 − ε )34

εuD(49)

この結果を引用すると、局所熱平衡状態では

α xx∂x jkfの寄与は無視できる程に小さいとみなすことができる。
更に慣性項を無視すると

1− εε 3d

また、中山らが行った数値実験8)より

証

tf

ε 3d 2

F = 1.75

本モデルの検証を、定常状態で巨視的一軸方向流れ− uˆiTˆfff2uˆ x4. 検2uˆ xf⎛ν⎞c0 ⎜ f + Fu D ⎟K⎝⎠K=⎞ ˆ⎟ uˆiT f⎠f(fV f ( ρcp ) fpKとForchheimer係数Fは次式で表される。

f∂x jfuj( ρc )

Ergunの経験則9)より、粒子径dを持つ充填層の透過率

以上の議論より、(35)式はDuˆiTˆf
Dt21cxxε (1 − ε ) Pe d12≈2150 (1 − ε ) Pr + 1.75Pe d(50)cxxはモデル定数を、Pedは粒子径dに基づいたぺクレ数⎪⎧ ⎛ ν f= ⎨c0 ⎜ ε + F ε 2⎩⎪ ⎝ Kfujfuj⎞ ⎪⎫⎟⎬⎠ ⎭⎪−1Ped=uD d/αf を表す。上式より、低ぺクレ数においてαxx/kfがPedの二乗に、高ぺクレ数においてPedに比例する事がわかる。これは局所熱平衡状態における、巨視的一軸ff⎧h ui Afsf ∂ Tfs⎪ˆˆ× ⎨ ui u j−Ts − T f∂xVρc( p )f fj⎪⎩(⎫f ⎪⎬ (45)⎪⎭)上式を巨視的熱物性テンソルαijとベクトルλiを用いて方向流れ系における従来の知見と一致する。

5. 結論

本研究では、局所熱非平衡状態における熱分散モデ
ルの構築を試みた。得られたモデルは、局所熱非平衡整理すると、− ( ρ c p ) f uˆiTˆf = α ij

α ij =f∂ Tf∂x j( ρc )ph ui⎛ν fs− Tff) (46a)状態における熱分散流束が、流体と粒子の温度差によ(46b)する。また、本モデルから求めた局所熱平衡状態でのり影響を受ける可能性を示唆しており、Quintard 7)らが提案した two-equation モデルにおける示唆と一致f⎛ν
⎞c0 ⎜ f + Fu D ⎟⎝K⎠λi =(− λi Tsfuˆi uˆ jf

巨視的一軸方向流れ系における熱分散係数は、従来の知見と一致することが分かった。

Afs⎞Vc0 ⎜ + FuD ⎟ f
⎝K⎠

(46c)

(46a)式はQuintard7) らが提案したtwo-equationモデルと類似した形となっていることが分かる。これは局所熱非平衡状態における熱分散流束が、流体と固体粒子の温度差に影響を受ける可能性を示唆している。
また、局所熱平衡状態を仮定すると、上式は中山ら8)

参 考 文 献1) Kaviany, M., 1995, “Principles of heat transfer inporous media, 2nd ed.,” Springer Verlag, New York.
2) Whitaker, S., 1969, “Advances in theory of fluid
motion in porous media,” Ind. Engng Chem., 61,14-28.
3) Slattery, J. C., 1967, “Flow of viscoelastic fluids
through porous media,” AIChEJ., 13, 1066-1071.
4) Whitaker, S., 1967, “Diffusion and dispersion inporous media,” AIChEJ., 13, 420-427.

22

局所熱非平衡状態を考慮した充填層伝熱モデル

5) Hsu, C. T. and Cheng, P., 1990, “Thermal dispersionin a porous medium,” Int. J. Heat Mass
Transfer, 29, 2002-2006.
6) 数値流体力学編集委員会 編, 数値流体力学シリーズ
第三巻 乱流解析, 東京大学出版会, 1995, 259-269
7) Quintard, M., Kaviany, M. and Whitaker, S., 1997,“Two-medium treatment of heat transfer in porous

media : Numerical resultes for effective properties,”
Adv. Water. Resour., 20, 77-94.
8) Nakayama, A. and Kuwahara, F., 2005, “Algebraic
model for thermal dispersion heat flux within porous
media,” Wiley Inter Science, DOI 10.1002/aic.10503
9) Ergun, S., 1952, “Fluid flow through packed column,”
Chem. Eng. Prog., 48, 89-94.