デルタ列関数によるステファン問題の解析

NAOSITE: Nagasaki University's Academic Output SITE
Title
デルタ列関数によるステファン問題の解析
Author(s)
福地, 信義
Citation
長崎大学工学部研究報告, 15(24), pp.59-64; 1985
Issue Date
1985-01
URL
http://hdl.handle.net/10069/24177
Right
This document is downloaded at: 2014-11-01T13:20:11Z
http://naosite.lb.nagasaki-u.ac.jp
長崎大学工学部研究報告 第15巻 第24号.昭和60年1月
59
デルタ列関数によるステファン問題の解析
福地 信義*
Analysis of Stefan Problem by Using the
Function
of Delta−Sequence
by
Nobuyoshi FUKUCHI*
The oil tankers¢arrying high pour point crude oil almost are equiped with the tank heating system
consisted of meandering steam pipe to reduce the.viscosity of crude oil. However, the accident of
oil solidification occasionally occurs during pumping up cargo o量l from an oil tank of tanker to land for
the sake of.unsuitable oil heating system. For a settlement of this case, the heat energy of steam
is transfered to crude oil in solid state to melt it.
In.ord6r to量nvestigate the heat transfer with moving.solid phase boundary, so・called Stefan problem,・
two dimensional numerical analysis is carried out in which the enthalpy is expressed continuously by us童ng
the function of delta_sequence at interface of both phases. Furthermore, the spreaded melting zone and
the oil temperature by using paraffin are measured in quasi−two dimensional tank with two heat3rs.
The results of theory and experiment are compared for judging the lustice of analytical method.
1.はじ.めに
解析し,固液両相内の温度分布と液相内の対流速を計
高流動点原油を輸送するタンカー.において,その装
算する。また凝似2次元タンクによるパラフィンの溶
備した加熱システムが適当でない場合には船底部に原
融実験を行い,溶融域の進展状態,温度分布などにつ
油を凝固させ,揚油不能となる事故が起ることがある.
いて数値計算の結果との比較を行い,解析法の妥当性
この場合には,船底付近に設置された蒸気管に通気
を調べる.
することにより凝固状態の高流動点油を溶融し,さら
に低粘度化を行なう.
2.ステファン問題
このような,固液相境界が移動する場合(わ熱拡散,
2..1 基礎式
いわゆるステファン問題1)を解析するためには相境界
高流動点油を加熱により溶融する場合のように,相
における融解潜熱の扱いが問題となる..本論文では相
変化があり,固液相境界が移動する場合の熱拡散の状
変化時におけるエンタルピーをデルタ列関数を用いて
態は以下の諸方程式により表わされる.
連続的に扱い,固定メッシュの有限要素法により数値
(1)液相内
昭和59年10月1日受理
*構造工学科(Department of Structural Eng量neering)
60
デルタ列関数によるステファン問題の解析
液相内では自然対流による熱伝達が起こる.油を近
2.2 デルタ列関数を用いた熱拡散方程式
似的にニュートン流体とし,流れを層流とすると,速
(1)エンタルピー関数
度場と温度場の状態を表わす式は次の通りである.
(5)式で表わされるステファン条件を相境界で厳密
なお式はtensor表示とし,添字は総和規約に従うもの
に満足させ,熱拡散方程式を解くことは困難なので,
とする.
次のエンタルピー関数H(θ)を用いて近似的に解析を
1)流れの連続式
行う2
H(θ)イ・(θ)μθ)dθ+μθ・)仇Y(θ一θ・)(6)
器一・ (δ=1,2,3)(1)
ここに,π‘は流速である.oσ妥は座標を表わし,」σ、,コσ2
ここに,Y(コσ)はHeavisideのステップ関数であり,
を水平方向,」C、を鉛直方向(上向きを正)とする.
θmは油の融点である。
iD運動方程式
高流動点原油のような混合体では融点にある程度の
重力場にお・けるNavier−Stokesの運動方程式は次式
幅があり,また比熱,比重量の温度依存性が小さいこと
で表わされる.
を考慮するとエンタルピー関数は次のように表わされ
讐+毒(・幽)一一去藷
る。
+素{・(∂%‘十∂賜∂:r」 ∂:r‘)}一β9・θ(2)
H(θ)一・ρθ+呵コω・(θ一θ・)dθ(7)
ここに,Pは圧力,孟は時間,θは油温である.
ここに,ωk@)はデルタ列関数であり,Diracのデ
またg、=g2=0, g3=一g(重力の加速度)であり,レ
ルタ関数δ(勾との間には次の関係がある.
は動粘性係数,βは体積膨脹率,ρは比重量である.
δ(oσ)=limωk(oじ) (8)
(2)式の境界条件としでは固液相境界およびタンク周
k→oσ
壁においてno−slipの条件をとり,自由表面では鉛直
デルタ列関数の例を3種類Fig.1に示す。ここでは
方向の流速を0とする・
最も扱い易い,
liD熱拡散方程式
ωk(oo)=κ/π(1十κ『コσ2) (9)
器+素(ω,θ)「妾(畷)(3)
を用い,κとωk(00)との関係およびこれに対応するエ
ここに,κ、(=λ、/cρ)は液相内の熱拡散係数,λ・
ンタルピー関数正1(のをFig.2に示す。デルタ列関数を
は熱伝導率,cは比熱である.
用いることにより,融点において潜熱のため不連続と
(2)固相内
なるエンタルピーを連続的に扱うことができる.
凝固相の熱伝導を表わす式は(3)式において移流項
がない次の方程式である.
誓一素( ∂θκ8 ∂:r」) (4)
Wk
亜
30
L一冷㎏蝉
2.5
ここに,κ・(=λ、/cρ)は凝固相内の熱拡散係数で
ある.
20
(3)固液相境界
相境界では,流出入する一流束の差が融解潜熱に等
t5
しくなる,いわゆるステファン条件と呼ばれる次の熱
1.0
’
収支方程式を満足しなければならない・
擁論、
レ響
∼
’
為(∂θ∂η)。一為(誓)、一ρσ・誓(5)
1
Q5
、1
1’
ここに,ηは固相境界に対し外向き法線方向を表わ
0
積当りの融解量を意味する.またσ、は融解潜熱である.
一3 、・乏 ∀1
)
し,添字8,Zをそれぞれ固相側,液相側の値を示す.
呪は・η方向の境界移動量,すなわち境界上の単位面
へ /へ
一〇,5
Fig.1
‘’
、 ’・ 〆・、
1 》冨・3
0
Three examples of function of delta−
sequence。
福 地信 義
35
60
Wk
誰Wk=禰
30
(x104)
r票∫/
詮50
ll//
H(θ)
ii
ii
1i
HeQV
54ρ
/
oil
@ イ
口
O0
(200
む
ii
20
RユrQ箔
O0
妻45
lレk=10
2.5
00
.鋼艦:翌一
55
61
緬ノ
1.う
ii
35
i:k.5
ao /H・c・θ+q・践
t5
25
20
1.0
k=2
Q5
ゆ ツ
15
k=1
/1
./1
;100
O0
3
ヤム へ
暮
/
亀
.逃
〉
ぐ
T0
’ぬ。
lo目e oil
勤 亀
\
Mi @ ノ\
pS
// ソ1
11
_7豪ニーメ’
@ ,ノ ノ
Q0
モ
\
[email protected]{
、℃
@ で)[ず一_
tO
一一 一
0
受心、
,’窪,ノ
melting point
O5
一4 −2 0 2 4
x
3D
10
ヘt轟,.,罰θおec)
Fig.2 A fUnction of delta・sequence and enthalpy
噂
ShQi
・くy
,Rユ十N十†十
手
「\ 辱
ffin
B貢
モ窒浮
@5
ゆ十
●
十十幸
「
一、
fUnCt孟On USing it
210
2Q 30 40 50 60 70 80
temperQture of oil(。C)
〔fil㎏d up by pαrqffin〕
ceUing:(1cryb乏e
…m/ †,・。、d認聰)
item
uhit
symbo1 丁ach拍g
[email protected]
Minas
[email protected]
1 (8柵kne$)
暢〕 層ゴ 漏
rqd iQtive 卜edねr セGnslbrmer
thermoイbw meter チ
bottom:「x)1ystyro【
fOQm
m・Qsu…g・㎝・・・…. 門ヒ…一r
T
initial using
sur口
heGter
秩Bunding 盾奄s
・・≠狽・
code
(34.2W)
狽・撃o.
狽・香D
陛議
P。ur p。int
o C
Melting point
3Kg/m
Paraffin
[email protected]
[email protected]
o C
θP
30
35
鼻
一
一
臼5
30.86×10
30.85×10
一
42’》44
30。81・X10
Speciffc weight
3Kg/m
ρ
ρ・9,(1一・.7・1・一3日9
= 15
(1−0.9Xτ09
Cubical
[email protected]
Specific heat
1/K
」/Kg・K
Latent heat
[email protected]
Themal
モ盾獅р浮モ狽lty
・」/K
.W/m・K
β
C
qL
λ
0.57/ρ
C司ア40//繭δ+3.8(θ一15)
2.。9∼2.26×105
λ・117(1一.。.54・10−3θ)/ρ
・・t・・θ・.・i1・・mper・t・re(.C)・θ15・θ一15
(right)35.9W
M2−1
28。7。C
21.3。C
i1eft)34。2W
M2−2
36.1
26.0
M2噂3
42.5
19.6
Fig.4 Physical properties of high・pour piont
o三1
(unit:mm)
Fig.3 Quasi−two dimensional tank model and
ここに,D/D孟は流体力学的微分であり,1固相では
measuring condition.
D/1冗=∂/∂ち液相ではD/D置=∂/餌十篤∂/∂筋となる.
さらに(11)式は次のように書き改められる.
なお,(9)式におけるκの値}よH(θm)=ρθm(C+αL/
[1+σLω。(θ一θ。 C)P器一蹴(・署)(12)
θm)とωhの最大値ωk(o)ニκ/πとの比較からκ=πσL/
θmとした.
(2)熱拡散方程式
3.解析結果
(7)式のエンタルピー関数を顕熱容量で割ると次の
3.1 パラフィンによる溶融実験2}
ようになる.
(1)実験方法 . 一
ん(θ)一
∴齠サπω・(θ一㊧dθ(1①
高流動点油の加熱による溶融範囲の拡がり方と温度
分布を調べ,.さらに2..2で述べたデルタ列関数を用
このん(θ)を用いることにより,、(3)式および(4)式で表わ
される熱拡散方程式は,相境界における潜熱吸収の取
いた解析法による解とを比較するために,ガラス製の
し
凝似2次元タ.ンクに融点42∼44℃の固形パラフィンを
扱いを考慮した,すなわち弱い意味においてステファ
入れて,2本の円柱状電気ヒーター(加熱管と称する)
ン条件を満足する次の方程式に変換される.
により加熱を行い計測を行った.この実験装置と実験
号}一素( ∂θκ ∂lr」)r.(11)
要目をFig.3に示す.
溶融域の拡がり方はガラス面に貼ったフィルム上に
デルタ列関数によるステファン問題の解析
62
相境界を直接記録した。温度分布は非接触型熱放射計
考えられる.
を用いて,縦横11×13の格子点上の熱放射量を計測し,
速度と温度に関する方程式(1),(2),(12に対しGaler・
経験的に求めた熱放射量一温度の関係式により温度へ
kin法を適用して有限要素法のための定式化を行い,
換算した.なお熱放射量は周囲の熱的環境に影響され
2次元問題として扱った.有限要素としては8節点の
るので,これを考慮して温度への換算式を作った.
2次要素を用い,定式化後の矯,P;θに関する非線形
実験に用いたパラフィンの粘度一温度特性およびそ
連立方程式はNewton−Raphson法により数値的に解き,
の他の物性値を高流動点原油の代表例であるミナス原
時間微分項は差分で近似した.
油(インドネシア産)および大慶原油(中国産)と比
計算対象はミナス原油とし,Fig.7に示す計算条件
較して3)・4)Fig.4に示す。
により解析を行った.また加熱開始時には加熱管まわ
(2)溶融状態と温度分布
りの3×3の要素は液相要素とし,他を固相要素とし
加熱時の溶融域の進展の様子をFig.5(a)に示す.
た.固相要素はその平均温度が流動点以下であるか否
なお,Fig.5(b)は,3.2で述べる数値解析法により
かを判定し,超えている要素は液相要素に転換した.
求めた溶融域の時間的変化である.溶融は最初加熱管
(2) 溶融状態 ・
上方に向って起こり,上面に到達すると液相が固相よ
液相内対流について計算した流速分布の時間的変化
り軽いことと膨脹により上面横方向に液相が拡がり,
の様子をFig.7に示す.対流は加熱管上方への上昇流
次に上面から下方に向って溶ける.加熱管より下部は
と側方の下降流から成り,溶融域の広さに応じて循環
対流が起らず移流の効果がないため,伝熱量が少ない
流の数が変化するのが見られる.
ので,加熱管より上部が完全に溶解するまで融解しな
Fig.8はF量g.7に対する温度等高線の時間的変化で
ある.パラフィンによる計測でも見られるように,固
い.
Fig,6は弾弓分布の経時変化を示す.いずれも固液
相境界では融解潜熱が大きいために,両相とも境界に
対し法線方向の温度勾配が大きくなるのが見られる,
3.2 数値解析
(1)計算法
相境界が移動する場合の熱拡散状態を解析する方法
として,縦横8×9の固定メッシュによる有限要素法
を用いる.ただし,固定メッシュのため移動する境界
②竜一120mln
を完全には追従できないが,相境界付近では非常な高
⑤t・330min,
40
粘度のため流れはほとんど起らず,流速計算において
弓。
境界位置のわずかな食い違いあまり問題にならないと
60 tOP 60
50
6
4
む
15018・
45
1:。
30
4
;:1
③t−180mln
90、
15
15
1i。
⑥t・420min.
:1:
8
270
1:よ
27
@8
0
30010
ll。
330
00
ii§
P2
iii
330
il崔
420
heαter
meGsurement M2_t
㎞甘
12
bottom
(unit of time :minute)
(G) π灘3urernent
h・・tq・antity・Q1(θH一θ1階12x2
(・nit・f ti…λ・/・戸2・。.、・・1・一3)
(b) CGIcutGtion
(unit:。C )
Fig.5 Growth of melting zone in tank model
Fig.6 Distribut量ons of oil temperature in
with two heaters.
measurement M2−1
福 地 信 義
63
ρ。0.ヨx1σ3 ④入tlc
①撞ノ。
28
9
【2。2鴉x1σ3
0
3
3
3
8
4
37
0
38
4
一
36
入t1・ρ匹147・1σ3
⑤入ゼcρ匹331x1σ3
7
2
2
7
8
② 入t1、12.11。.1σ3 ⑤入t1、12。3.31.1♂
4
52
2
ロロロロ
54
44
4
34
・L… 亀■←隔。79 噛し■ … ..・
[email protected]:;1獄でな駕ニモ㌧二=;::
i癩鱗髪磁;ii:・・
46
56
6
ロののうちもロキロヲロも メの ロほコ 40
3
3
58
0
’曹召’ @ゆゼ民」rガ. ○、\庵’メ」・や・・
’・ @・腎¶ . ,o・9摂 r r も煽○貝ρ,
コ や ロ
3
5
ゆサロやひコそセゆ ロサリヘヨぱサ ロロげヒ
’:;=や、こ=:::::1二二三ナこ:::
コ コ ロ ものぼ ロ コ ロ ヒ コロ リ ロ ロロサ ヒ ヒ コ
::::∵::::::::::∵}::,:
8
③入tlcρ【≧220x1σ
Calculating condition of ta6k mode1
2
3獄,c l2.1.脳.1(戸 (unit・・C)
1) Assumed length :1旨 2 閉
heGt quqntity
2) Heating oi1 =Minas crude oi1
(pour point : 340C )
@ Q/(θH一θ1)KA=6x2
3)Oil temperature near
heati・g・・i1・θH・6・.C
8
4)011t・mper珈re at b・gimi・g・θ1・3・.C
40
日ote….
j:Heat・transm1ゴsi。n c。efficient
5)S・一di・g t・mper・ture・醜・25.C
36
‘)日eat B辮翻llgC3:{1’鷲.8、,.2K
@A = Heat transmission ared
i…τ・mper・t・re・e・r heati・g・・i1
ニ1・Oil t・mper・加爬・t b・gimi・g
下
6
0.6875 275 Q275
掴45奇
Q3125 . heαter
⊥
Fig.7 Calculated flow pattgrn in model C2−2
Fig.8 Distributions of calculated oil temper・
and calculating condition.
ature in model C2_2
1.0
qL冨0
qL=218(J/㎏)
弩
お04
一L.
latent heat
M2〆
total
窪1・0
9・6
M2−2(36.10C)
9Q8
M2−3(42.50C)
9
M2−1(28.70C)
歪oβ
te「nperoture
、
surrounding
ち04
220.2
ε
he(1t qLkユntity
2Q2
Q/(θH一θ【)KA=6x2
0
(G)meQSUrement
qL=1・09 U!㎏)
轟・・
0 .Q5
10 C、㎜t
ッ。轟、225
01234567
0
3&σ・)蕊
8 9 10
セime(hour)
Fig.9
(b) CGヒu[otion
Variations of calculated melting zone
1.0
岩、,
corresponding to three kinds of latent
heat.
9
婁Q6
\
も04
液相境界において両相の温度勾配は融解潜熱吸収のた
め大きい.これに反し,境界から離れた固相部では伝
熱量が少ないため,また液相部では対流による掩乱の
ため温度が均一化し,温度勾配が小さくなっている.
\1;二1羅}iW
£
2Q2
0
C2−1(q=3x2)
ロ
note二q=Q!織一θp KA
0
α51 唐謔T..遣2,、225
30
35
(x1σ3)
2.2で述べた,デルタ列関数を用いた熱拡散方程
Fig.10 Variations of melting zone:correspond・
式により相境界の熱収支において潜熱の扱いが妥当で
ing to supplied heat quantity.
デルタ列関数によるステファン問題の解析
64
あるかを調べるために,高流動点原油の融解潜熱σL=
4.結 言
2.18×105(J/kg)に対しその%および0の場合について解
高流動点油を対象に;凝固状態から加熱により溶融
析を行い,この3種類の場合についてタンク全体を1.0
する際の相変化を含む熱拡散の状態を解析するために,
とした溶融量比の経時変化をFig.9に示す.溶融量比
熱拡散方程式においてエンタルピーをデルタ列関数を
の曲線は固定メッシュのために波うっているが,潜熱
用いて温度に関する連続関数として表わし,固定メッ
の大きさの差による溶融速度の違いが現われている.
シュの有限要素法による数値解析を行った.またζれ
Fig.10(a)はパラフィン溶融実験にお・ける外周温度
に対応するパラフィンの溶融実験を行って解析解との
と溶融量比との関係であり,加熱開始2∼3時間は初
比較を行った.これにより完全ではないが本解析法に
期油温に影響されるが,その後は周囲に逃げる熱量の
よるステファン問題の解析が可能であることが分った.
多少により溶融速度が決っている.これに対して,外
しかし,溶融域の進展が階段状になるなど,移動境
部への貫流でなく供給熱量を変化させた場合について
界の取扱いやデルタ列関数におけるκ値の選択法など
計算した溶融量比の時間的変化をFigユ0(b)に示す.
に問題が残されている.
固相と液相およびタンク全体の各平均温度の経時変
本研究は文部省科学研究補助金を得て行ったことを
化を,計測値と計算値を比較して,Fig.11に示す.タ
付記する.
ンク全体の平均温度は溶融域が小さい初期には固相の
温度に近く,溶融が進むに従って液相の温度に近づく
が,加熱管より上部が全て溶け溶融速度が鈍ると潜熱
参考文献
に奪われていた熱が液相温度の上昇に使われ,油温が
1)山口昌哉,野木達夫:ステファン問題(1977),
急激に上昇する.固相温度の時間的変化は小さく,こ
産業図書
の図からも伝導のみによる伝熱は少ないことが分かる.
2)福地信義:高流動点油の加熱による溶融と熱拡散
について,日本造船学会論文集,第156号(1984)
(1) meqSUrement
3)渡辺,他:貨油タンクの高粘度原油洗浄,昭58春
60
70
季船舶技術研究所発表会講演集(1983)
鯉P蟹ノ
55
65
ノト虚』’ルぜ
0
劉
量45
4)日本造船研究協会:タンカーのタンクヒーティン
月目。
/
_ノ 説・鏡
謹
議/
ヨ
量55
塁、。
警
〆
崔go
Aノ 簿
蓉35
45
一一鰻鰹∼’
30
グに関する研究,協会資料,No.79(1968),No.99(1969)
ゲ
/ び
40
25
35
20
30
15
25
99[id phqse
/ 阜一}響,・・ひ●一一
01234567
0 1 2 3 4 5
time(採)url
(b)meosuremenい42−2
tr爬(㎞r)
(Q)meqsurement M 2−1
(2)co【culotion
48
50
48
Heat quantity
Gy(像一θ,)漁=6x2
ε46
46
44
睾44
42
§42
540
ε
1iquid phase
豊暴2
譜
ェび
36
’豊
Heat quantity
Q/(一一a)KA=12x2
評5e
、、Φ、も
40
二
538
,〆
38
説漣
36
巨
9
34
34
=
32
SOlid
phase
o 32
28
solid phase
30
30
0
0.5 1。0 1.5
@ time
.二又,,、誕5・・1δ9・・
28
む ら しむ の 段1げ・}2’‘
、b、、。。,u、謝謝
(o) ooヒu1αtion l C2−1
Fig.11 Variations of averaged oil temperatures
in solid phase and liquid phase.