平成5年虔 研究報告 大分県工業試験場

1 CA馴こよる熱交換器配管の乱流熱伝達特性の解析
機械部 大 塚 裕 俊

要 旨

CAEによる熱交換器配管の乱流熱伝達特性の把握のため、k−e乱流モデルによるソフトウェアを用いて単管モデル及び千鳥型の熱交換器管列のモデル化と熱伝達特性解析を実施しその評価を行った。必要に応じて基礎方程式や乱流のモデル化の手法、熟的境界条件や離散化手法についても概説した。
1．緒 言

そのためまず時間平均則と乱流の空間等方性をポイン

近年地球環境の保全やエネルギー資源の枯渇化が話題トするk【∈モデルについて概説するとともに、熱伝達となり、資源の有効活用や省エネルギーに関する技術開特性の把握のため重要な熱的境界条件について、そのモ発が諸処で行われている。

デル化を中心に記述する。また同時に利用したソフトウェ

熱エネルギー関連幾器の中でも、熱交換器はシステムアでの離散化手法についても言及する。次に数値解析ののエネルギー効率を決定的に左右する重要な要素続器で前段階として単管モデルでの熱伝達特性の解析を行ったあるため、以前からその性能向上のため、さまざまな研上で、千鳥型の熱交換器管列のモデルイヒとその熱伝達特究が進められてきた。

性解析へと進み、その評価を行う。

しかしながらその多くは実験的・経験的な研究が主流

なお乱流状態の空間非等方性を前提とするレイノルズであり、最も基本的な管列塾熱交換器についても現実的応力モデル（Reynoldsstressmodel）や、近年研究がにはその対象となる流れは大きなレイノルズ数の支配す進められた統計力学による理論的乱流モデルであるる乱流状態となる場合がほとんどであるため、理論的なRNG（RenomalizationGroupTheory）による手法アプローチは困難な場合が多かった。

は有力であり、その適用と評価も大変興味のあることろ

ところが近年の流体の数値解析手法の進展は、乱流状であるが、今回は別稿にゆだねることにする。［1］態の流れについてもその特性把握において一定の条件下ではあるが、成功をおさめつつある。
一般に乱流のモデル化では、複雑な方程式系を時間平均則を利用した仮定に基づいて簡単化し、経験則で補足

2．基礎方程式とk一合モテ■lし
2−1 基礎方程式
流体の挙動を把握するための物理的な方程式は質量保して解ける形に追い込み計算するのが常道となっている。

存式，運動量保存式，エネルギー保存式の3つである。

工学的な計算では例えばこのような時間平均値による特いずれも直交座標上で流体の微小体積中での各々の物理性把握でも十分な場合が多く、CÅEによる理論的アプ量の保存（つりあい）から導出されるものである。これローチがその力量を発揮する場ほ、乱流を取扱うモデルらの方程式はすべて一般的なものであり、その導出につの特性向上などに伴い今後ますます広がっていくと考えいては流体力学の書籍を参考とされたい。［2］［3］［4］られる。
本論文では、現在乱流状態の流体の数値解析手法として一般的な評価を得ているk−−ざモデルを利用したソフ

質量保存式：㌍＋蓋（1）．

トウェア（FLUENT）を用いて熱交換器配管の乱流

β：密度熱伝達特性の解析を行い、主とLて経験値との比較によ

u：速度りその評価を行う。
ー45−

平成5年度 研究報告 大分農工業試験場となる。pu iu J部分はレイノルズ応力と呼ばれ、乱運動量保存式：流によって生じる項となる。この変動項（相関項）を平（pui）＋怠（puiuj）＝】＋pgi

均流（主流）の特性により表現すれば、1つの乱流モデ（2）ルを構成することが出来る。層流での応力テンソル（3）式との相似によりブシネスク近似として知られる次式をp：圧力仮定する。

g：重力加速度ニ†I

p u．u∵＝

・

なお丁＝ほ応力テンソルであり、流体の粘性による応いては、文献［3］を参照のこと。）

ここで1

但し ∂＝＝， k：乱流エネルギー＝uiu j〃t：乱流粘性係数〃：粘性係数エネルギー保存式：∂㌻㌻

豊 霊†）一号〃t‡計∂ij（7）但し ∂＝＝r▲j＝ト（霊†＋豊〕一号〃記†。∂i∫（3）∂㌻γ＋＋〃t（力項を表し、以下で与えられる。（なお丁＝の導出につkは乱流エネルギーとして定義される速度の変動成分による物理量である。また乱流粘性係数〝tは乱流の速（ph）＋蕊丁（puih）＝度スケールと大きさ（長さ）のスケールの積に比例する（ス）＋＋ui

監丁＋S≠（4）という仮定（Prandti−Kolmogoroff仮定）により（速度スケール∝kl／2、長さスケール∝k3／2／￡）以下のようよ：熱伝導率に与えられる。∈はkの消散率として定義される物理量
S≠：熱源である。

h：エソタルピー〃t＝pCu ミk2（8）

2．2 k−∈モデル
C。：実験による比例定数以上の方程式は非線形項を多く含み、レイノルズ数（Re）が増加すると流体の粘性による効果よりも慣性

このようにして定義された新しい物理量、即ち乱流エによる効果が大となり、大小スケールの渦を伴う乱流へ

ネルギーkとその消散率∈についても傲づ＼体積中での保と移行し、限られたメッシュ数やコンピュータ資源では存式を構成でき、これにより全ての方程式系を解ける形

到底コスト的に実用的でない数値計算となってしまう。に完結させることができる。

そこで乱流を取扱うため、乱流変動のある流れでの物（pk）＋丁理量を時間的平均値と変動成分の和として表現するアプ

ローチ（レイノルズ平均操作）を用いる。［5］［6］

蓋丁（莞豊（puik）ニ）＋Gk一戸g（9）今、速度uiを時間平均値uiと変動成分u 主によって（pり＋（p…）＝㍑上＝ul＋u ま（5）tとあらわす。これを（2）式の運動量保存式に代入し、∂Ⅹi＼已㍉一）＋Cl古意Gk珊一C2吉野む㌃（10）

あらためてuiをuiと書き直すとCl‥C㌃∴実験定数

孟：（pulu∫）＝ぺ怒＋富㌃（pu iu j）＋pgj

−−46−

￡2

げk、げ。：乱流拡散プラントル数（6）平成5年度 研究報告 大分県工業試験場

なおGkはkの生成項であり、次式で与えられる。

しかし乱流においては、これと異なり壁の対数法則を拡張することにより壁面熱流束を下記の式から求める。

Gk＝〃t（† ＋†）†（11）［7］わ（△T／△y）＝⊥旦＿む♂k，k y＋ p rq

αり Cl，C2としては経験的に次の値が用いられ1n（Ey＋）

5る。［7］＋幸（禁）す打

グk＝1．0、♂6＝1．3、Cl＝1．魂、C2＝1．92＋毒（書）（託−1）（15）

ここで ユー ＝流体の熱伝導率

3．流体物性値と熱的境界条件
3，1 流体物性値△T ＝Tw−Tf

理想気体については状態方程式から密度は以下のよう△y ＝壁から壁に隣接するにあらわされる。

格子点pまでの臣巨離

q ＝壁面への熱流束

P＋Po（12）P r ＝分子プラソトル数
R こガス定数

P rt＝乱流プラソトル数（壁面でL2）

P：基準圧力A ＝26．0（van Driest定数）

P。：局所的静圧ガ ニ0．42（Ⅴon Kaman定数）
E ＝9．81（滑らかな壁での実験定数）

圧縮性効果はそれほど大きくないとして今回はP。＝0とした。すなわちpはTのみの関数として扱う。
同様にして

y＋ ＝境界層厚さ

上式によって、得られた熱流束q■と△Tから乱流で熱伝導率 よの熱伝達係数hfが求まる。

定圧比熱 C。

また、下図のように流体の外に外部熱媒体（水）の温粘性係数 郎

度T∞と、壁（固体）の熟電導率或wを仮定し、壁の厚さなどの流体物性債も温度Tのみの関数とし、物性表などをtとする。
すると流体から外部への熱流束は次式であらわされる。

からTの多項式近似であらわす。

［5］

3．2 熱的境界条件

q＝hぎ（Twal卜Tぎ）

流体が層流の場合は流体と壁との問の熱伝達は、流体＝hext（Text鵬Twall）の熱伝導率戌によって局所的に一次式で近似される。［5］q＝汗wall ＝汗i、16、Lq＝h（Text一丁り（17）（下図参照）（13）よw

q：壁面熱流速TwallÅT∵（Tw−Tぎ）壁面温度T∇と近傍の流体

hf Tf

セルの温度Tfの差△y：壁から壁に隣接する格子点pまでの距離

く t 〉

よって、この場合勲伝達係数h fは
h f（16）（17）より総合的なhは1（14）（18）＋としてあらわすことができる。

（図1）

有−47−

平成5年慶 研究報告 大分県工業試験場

壁（固体）の熱伝達係数hextは次のようになる。

JvolumeV

hext＝（pu）dV＝JAPu・dA（23）（19）

Textは壁（固体）の外側の表面温度であるが、これを管内を流れる水の温度T∞と等しいとすれば、外部（水）への熱移動の定常状態を求めることが可能となる。

4．0 離散化手法について（図3）

今回利用したソフトウェア（flu en t）では、コ
図3のコソトロールボリュームに対して（20）式を積分

ントロールボリュームに基づく有限差分法を用いている。

すると次式を得る。

本手法では、微分方程式を各セルに対して積分し、それぞれの物理量がコントロールボリューム基準で保存され（puA）e−（puÅ）w＝0るような有限差分表現式を得ている。セル上では図2に（24）

Je＝（puA）e、Jv＝（puA）wと表せば上武は示すような格子に基づいて変数を保存している。この手法では、すべての方程式のすべての変数（圧力、速度成Je−Jw＝0分、および他のすべてのスカラー量）に対して同じコントロールボリュームを用い、全変数をセルの中心で定義（25）運動量の（21）式を積分すると次式を得る。

する。［5］［8］（pu2Å）e一（pu2A）w＝−（pe−pw）Å

数分方程式の積分について、一次元の方程式を例にとuE￣u p＋（〃〔り、直交座標に対して考える。まずはじめに、流体につu p…uwÅⅩA）e−（〃〔いての3つの方程式（連続の式、運動量の式およびエネ同様に変形してルギーなどスカラー量￠の式）の一次元方程式を考える。

Jeue−Jwuw＝−（pe…pw）A

十〔蔓慧（uE−up上蓋：（up−uw）Å（27）

スカラー量の（22）式を積分すると同様に次式を得る。
×CeiiStorage− −Location−−一一−−＋S￠

仁・一ニ：トGて・iしiL．irle

訂妄
∂
㌻妄
∂
㌻妄

㍗GriLin∂h lG eⅠ−1thGridLine

これらの代数方程式系は一般的な数値計算の手段により解くことが出来る。さらにこの原理をBFC（Bo dyFit t e d Co o din a t e）など曲線座漂系を用いた柔軟な座襟系に拡張することが可能である。

（p u）＝05．単管モデルでの解析（puu）＝二託＋孟♭（霊〕（21）（pu￠）＝蕊㌢霊＋S￠単管モデルによる乱流熱伝達解析を行った結果を以下に示す。
（22）

気体流（空気）〜管内流（水）問の伝熱をモデル化し、（20）〜（22）或は発散定理（ガウス定理）を用いて各コソ

図1に示すようなグリッド分割と境界条件で解析を行っトロールボリュームに対して積分できる。

た。結果は以下のような無次元数で表示する。

−48−（28）

平成5年度 研究報告 大分県工業試験場

Zu s a u s k a sの式による円柱まわり流れの平均

また得られた結果を図2〜図12に示す。流れ図等でヌセルト数の経験式［9］はRe＝3000の場合を示す。

Num＝0．26ReO・6prO・37（29）

レイノルズ数 Re＝宝P rO・彷
（P rw

（32））

（1000＜Re＜200000）
ヌセルト数

Nu＝（30）Prw：管壁での流体のプラソトル数
P r ：主流での 〝圧力降下係数 f＝Ap／‡pu喜（31）も同時に示す。（図8）u。：入口での流速

D ：配管の直径∠一OUTLET〃 ：粘性係数p ：密度グリッド＄2×62×1

hf：熱伝達係数△p‥圧力低下（上流入口〜下流出口間）
熟伝達境界

よ こ領域平均温度TAに対する熱伝達率

対称境界

図1単管モデルの計算格子と境界条件．︑■−49−−−

三三毒毒室当表ヨ喜藩誓書 ト

図7 等乱流エネルギ繚図

11

諸表雷撃要叫ル

仙≡叫≡≡叫叫叫叫還還還ヨ当室瑠彗≡≡‰≡≡≡雲叫慧3叫≡盟‡滞
図5 等圧力線図

図6 等温慶練図

H二−こ＝．；．章一二三至こ．二∴≡妻＝≡言責ヨ要撃ごニ‰ 一⁝
こ≡蔓誉ぎこ鮎Ⅹ＝緩逐喜善書署ぎ≡三菱表責室室還貞一一一一−一言rT︐．ぎ．．言︼≡．ゝ．邑宴．ごy−≠irL

鮎は図4 等速度線図

コ≡こ≡∃霊叫叫当∃室叫一室塊

箆‡封

図3 速度ベクトル屈﹁≡三図2 流れ線図J≡≡■芸芸≡■ヽ■−・■・−一■−−ノー、一−、−・−■−7三≡三≡芸

叫1￣■

∴

平成5年度 研究報告 大分県工業試験場

N u4288 ユ0×1040・3

図8 ReによるNuの変化（単管モデル）

110×1044

図9 Reによる圧力降下係数fの変化（単管モデル）

T白U見F−To
TAVヱーで○

ヰβ①
ロ R8＝ 8000＋ Re＝＝10000◇ R8＝40000△ R8＝ユ00000

口 Re＝ 3008＋ Re＝10000◇ R◎＝皇0000△ R8岩100000

角度（8）

角度（8）

180

図10 管壁での圧力係数Cpの変化（単管モデル）

80

図11管整での温度の変化（単管モデル）

流れ方向増 ↓↓ ↓ゆ ↓↓ ↓50

図12 管壁でのNuの変化（単管モデル）

平成5年度 研究報告 大分県工業試験場

れはReの増大とともに配管後流部でのNuも増大して6．千鳥配管列モデルでの解析

千鳥配管列モデルによる乱流熱伝達解析を行った結果いくという経験値データとはかなり異なった結果となっを以下に示す。グリッド分割と境界条件は前のケースとている。［8］今後の乱流伝熱モデルの構成や適用についほぼ同等である。（図13〜図14）て検討を要する部分となろう。また管壁温度Tsurfについての無次元化したものが図11であるが、Nuと同様配管相互の影響による流れや熱伝達性能の変化がポイにβ＝900 付近で最大値があらわれる。

ントであるが、計算の経済性などを考慮した上で3列のまた、圧力降下係数fは図9のようにReとともに低千鳥配管ブロックを1単位としてモデル化する。

そしてこの2列目の配管についてデータを収集する。

下．している。しかし、管壁での圧力変化を無次元化した結果は以下のような無次元数で表示する。得られた結果

Cpの分布（図10）から、Reが増大するほど後流のを図16〜図32に示す。流れ図等ではRe＝3000

剥離位置が後退するとともに、∂＝900 付近でのCpの場合を示す。の最小値が大きくなりCp＝−3に近づいていくのがわ

レイノルズ数 Re＝かる。

umax・D
〃／p

ヌセルト数

Nu＝＝

千鳥配管列モデルでは、Reと管壁での平均ヌセルト

hf・D数Nuの関係は、、図22のように単管モデルで良好なj結果を示したReの範囲においてG rims o n式によ

圧力降下係数 f＝△p／‡pu害る経験値とその増加傾向においてよく合致している。但し、経験式によれば後部列の酉己管のNuは前部列の配管流れ方向のによる乱れの増大によって増加する傾向を示す。

配管ピッチ比 Ⅹ＝Ⅹ／D

また、管壁でのNuの分布では図26のようにReの流れ方向に直角の配管ピッチ比増加によりその最大値をとる位置が♂＝700 付近からY＝y／D♂＝1000 付近へと後退するとともに管壁の先端部（β＝00 付近）を除いて全体的にNuが増大してゆく傾

ひmax：最小流路での平均速度向を示している。これは図25のように管壁温度Tsurf

Ⅹ ：奥行配管ピッチy ：幅 〝の分布でも同様である。
また囲23のように圧力降下係数＝′ま、配管ピッチ比

Ⅹ，Yについては図15を参照のこと。

によって差はあるもののいずれもReの増大とともに減G rims o nの空気流についてのⅩ，Yによるヌセ少する傾向にある。管壁での圧力係数Cpの変化は図24のようにReが」＼さい場合、配管相互の影響が顕著にルト数の変化の経験値［9〕も同時虻示す。（図22）

見られ、流れの上流部管壁と下流部管壁での変化が複雑
Num＝0．522偽ReO・溺2（Ⅹ＝1。25、Y＝3のとき）（33）になる。
配管ピッチ比を変化させた場合を見ると、流れと直交する方向の配管ピッチ比（Y）の変化による管壁でのN

ア．考 察

uの変化ほⅩに比べて非常に複雑な変化傾向を示す。

単管モデルではReと管壁での平均ヌセルト数Nuの（図27〜因28）流れと同じ方向の配管ピッチ比Ⅹを

関係は囲8に示すように5000＜Re＜40000程変化させると、Ⅹの減少とともにNuは全体的に増加す度の範囲ではZu k a u s k a sの経験式とよい一致をることがわかるがYでは一概には判断できない。

示す。しかし40000＜ReではNuは経験値よりも（図27〜図28）大きな値となる。

このことは配管ピッチ比Ⅹ，Yによる平均ヌセルト数N

また、管壁でのNuの分析をみると、図12のようにuの変化をとっても同様である。

Reの増加とともにβ＝900（流れに直角方向）付近

図30のようにⅩの増加によりNuは減少しているが、を中心としたNuの最大値が増加している。しかし、こYを変化させた場合では図29のようにReの大きい範−51−

平成5年虔 研究報告 大分県工業試験場

囲でNuの増加傾向が示されるだけである。

意味することから明らかである。またⅩの増加によって圧力降下係数fはYの増加により顕著に減少する。 もfは減少する傾向を示す。（図32）（図31）これは流れと直交する方向の流路幅の拡大を

グリッド182×42×1 但しこの例は Ⅹ＝Ⅹ／D＝2 のケース
Y＝y／D＝2

図14 千鳥配管列モデルの計算格子（部分拡大）

流れ方向

X＝X／D
Y＝y／D

図15 配管ピッチ比 X Y−52−平成5年虔 研究報告 大分県工業試験場

岳こ一芸＝￣＝≡堅塁
図16 流れ線図（X＝2、Y＝2）

図17 速度ベクトル図（X＝2、Y＝2）

図18 等速度繰回（X＝2、Y＝2）

む王

図柑 等圧力線図（X＝2、Y＝2）

図21等乱流工ネ捗ギ緑園（X＝2、Y＝2）−53−

平成5年虔 研究報告 大分県工業試験場

寵管ピッチ比□ Ⅹ＝1．25 Y＝3＋ Ⅹ；2 Yコ＝3
◇ Ⅹ＝3 Y＝33000

図22 ReによるNuの変化

10000

30000

Ro

図23 Reによる圧力降下係数ずの変化（配管ピッチX比X＝1．25、Y＝3）ll￠1▲Rll台

0 nV n＞
ハリ nU ▲‖＞
nU nV O3 ∩︶ ∩＞R一lR3爪U￠＜U300000

図24 管壁での圧力係数Cpの変化

ロ十◇0∵で0TTO∴∴

1￡0O

e e Q
R R Rロ十◇80

図25 管壁での温度の変化（配管ピッチ比 X＝Y＝1．25）（配管ピッチ比 X＝Y＝1．25）

耽れカ向↓↓一﹁・−・．三・ヰ
↓．↓

ヨ己管ピッチ比↓

図26 管壁でのNuの変化⊂〕

Y三n ユ．25

十Yヨ 2く＞Y＝ き

図27 菅壁でのNuの変化（配管ピtソテ比 X＝Y＝1．25）（配管ピッチ比 ×＝ さ．25 Re＝30000）

ヰ∂∈⊃ロ R（】＝ 3000＋ Reニ10000◇ R①＝30000

配一】紆ピッチ比〔コⅩロE l．25

十Ⅹ出≡ 2◇

Ⅹ壬＝ 31．25

図28 管壁でのNuの変化（配管ピッチ比 Y＝3 Re＝30000）
−54椚

Y23

図29 配管ピッチ比YによるNuの変化（配管ピッチ比X＝L25）

平成5年虔 研究報告 大分県工業試験場□ Re＝ 3000＋ Re；王10000◇ Re＝ 30000□ Re＝ 3000
十 Re＝10000◇ Re＝30000

1。25

X
21．25

3Y
23

図31配管ピッチ比Yによる圧力降下係数fの変化

図30 配管ピッチ比XによるNuの変化（配管ピッチ比×＝1．25）（配管ピ・ソテ比Y＝3）□ Re＝ 3000＋ R8＝10000◇ R針＝30000

1．ヱ5

X
23

図32 配管ピッチ比×による圧力降下係数Fの変化（配管ピッチ比Y＝3）

参考文献

8．まとめ

ぷ−∈乱流モデルを用いた単管モデルでの乱流熱伝達

1）C．G．Speziale and S。Thangam， Anaiysis of an

解析（空気流）を行った結果、Nuの平均値は、およそRNGBased Turbulence ModelforSeparated Flows ，5000＜Re＜40000の範囲で経験値と良好な合Int．J。Engng Sci．Vol．30，No．10，1992

致を示した。しかし管壁のNu分布については経験値と2）日本機械学会，流れの数値シュミレーショソ，コロナ社，異なり、後流側でのNuの原著な増大が認められなかっ

1988

た。

3）原田幸夫，流体力学覆水力学演習，棋書店，1977

同様の千鳥塾管列モデルでは上記のReの範囲でNu

4）H．P．スウ，ベクトル解析，森北出版，1980の平均値は経験値と同様の増加傾向を示した。また配管

5）Fluent。Inc，FLUENT User s Guide，Fltlent．‡nc，

ピッチ比の影響については、ⅩとNu及びfの増減関係

1991

や分布傾向が明らかなのに比べ、YではReの大小など

6）毛利昌康， 汎用熱流体解析プログラムFluentア ，配管によりNuとの増減関係や分布傾向は複雑な変化を示し

技術，2。1988

た。但しYとぎは顧著な相関を示Lた。

7）B．E．Launderarまd D。臥Spaldirlg， The Numerical
Computation ofTurbulent Flows

9．おわりに本論文で利用した流体解析ソフトFLUENてについ，Computermethods

in applied mechanics arid Er唱ineering，3．1974
8）スハス Vノミタンカー，コソビュー一夕による熱移動とては、その利用や技術指導で（株）流体コンサルタントの流れの数値解析，森北出版，1985

毛利昌康氏、小林克雄氏にお世話になりました。お礼申9）日本機械学会，伝熱工学資料（改訂4版），丸善，1991

し上げます。

鵬55−