= aE hE
k +1
+ aW hW
k +1
+ b
( 7b)
材料 ,其物性不同 , k eff 不同 。 式 ( 1) ～ ( 6) 对于求解一维 ,二维 ,三维相变问题 均适用 。然而由于相变问题的一维模型最容易采用 实验结果加以订正和修正 ,因此 ,我们以一维问题的 数值求解及实验为例验证其正确性 。
太 阳 能 学 报 24 卷 498
方程可进一步简化为 [ 9 ] : ρ 9 h = div ( k ・ grad T ) + S h 9t
( 4)
2 ) 固液相变界面自动满足界面条件 , 无其他
2 数值分析
经典 stefan 问题的求解常用显热容法 。该方法 的主要优点是温度作为可直接求解的变量 。然而该 法有其明显的局限性 ,即对相变问题 ,相变只能发生 在一个固定温度上 [ 10 ] 。而且 ,经典 stefan 问题的求 解中 ,常常假设固液两相的密度不变[ 1 ] [ 4 ] [ 10 ] 。在本 文所研究的问题中 ,使用冰作为相变材料 ,冰和水的 物性参数完全不同 。因而 , 我们采用热焓法求解控 制方程 。 引入热焓法 , 相变问题可得到简化并具以下优 点
第 24 卷 第4期
2003 年 8 月
太 阳 能 学 报 Vol124 ,No14
ACTA EN ER GIAE SOLARIS SIN ICA Aug. ,2003
文章编号 : 025420096 ( 2003) 0420497207
h ( T m) + L + h =
T
∫c d T T >
Tm
T
1
Tm ( 1)
m
c d T T ≤ T ∫
0 s
对相变问题 : 固相 : V = 0 , T < T m ; 液相 : T >
T m ; 界面 : T = T m
2) 连续方程 :
收稿日期 :2001208228
: 1) 两相区控制方程与单相控制方程相似 ;
4 期 曾 艳等 : 固液相变蓄能中有效导热系数的数值分析与实验研究 499
对每一个时间层 ,使用已知的上一轮迭代值 <P 求出式 ( 12 ) 中各项系数 。然后用 TDMA 算法直接 k +1 k +1 k +1 k +1 求出式 ( 7b ) 中各点焓值 h P , h E , h W , <P 。根 据 <P 同法可求出各点焓值 h P , h E , h W 。比 较相邻两轮迭代值 ,如果满足式 ( 13) ,认为该时层迭 代已收敛 ,进入下一时层进行计算 。如果不满足式 ( 13) ,需用式 ( 14) 对所得焓值进行亚松弛 ,使用修正 值再进入该时层的下一轮迭代 。直至满足式 ( 13) 为 止 。最终 ,可求得每个节点在每个时层上的焓值和 液相率 。由此可求得所有相变材料熔解所需时间 , 并得到固液相变位置与时间的变化曲线 ( 请见第 4 部分) 。
δ δ t t 式中 : a P = 1 + 2α 2 , a W = α 2 ,αE = α δx δx δ t 0 k 0 ( 12 ) , b = hP - ρ l L ( <P - <P) δx 2 其中 , 下标 W , P , E 代表计算相邻节点 ( 如图 0 0 1) 。 h P 和 <P 分别代表前一时层的焓和液相率 。上 标 k 表示某一时间层上的第 k 次迭代 。
60 ℃,测定不同加热壁温下 ,相变温度场和相变界面
- h P ) < TOL ( 13 ) + ( 1 - ω) h P
k +1
k +1
( TOL 为收敛判据 , 计算中取为 10 - 4 ) = ωh P
k +2
( ω 为松弛因子 , 0 < ω < 1 , 计算中取为 0 . 8) ( 14 )
这样处理后 ,由于方程中没有了速度相 ,液相的数学 描述和固相完全一样 , 使边界条件和数学模型的求 解大为简化 。为了修正模型在液相计算中的误差 , 公式中的导热系数 k 用自然有效导热系数 k eff 代替 , 即可得方程 ( 5) 如下 : ρ 9 h = div ( k eff ・ grad T ) + S h 9t
( 7a)
N u = 0 . 18
if 1 ≤
l ≤2 X ( t)
3 υ β l g l ΔT Pr = ; Ra = ;Δ T = | T w - T m | α α υ l l l
整理式 ( 7) 得 :
aP hP
k +1
( 6) Pr 数取决于物性 , R a 数与物性有关 。不同的相变
固液相变蓄能中有效导热系数的 数值分析与实验研究
曾 艳 ,田怀璋 ,余 鹏 ,陈林辉 ,陈敬良
( 西安交通大学能动学院 ,西安 710049)
摘 要 :引入有效导热系数经验公式 [ 3 ] ,使固液蓄能数学模型得以简化为仅用能量方程加以描述 。运用热焓法求解能量 方程 ,得到计算相变位置随时间变化曲线 。并通过实验测得相变过程的实际温度场 ,得到实验相变位置随时间变化曲 线 。二者比较 ,证明了自然对流对固液相变换热的影响不可忽略 ,验证了该公式在 R a < 1010 , Pr < 105 下的正确性 。 关键词 : 有效导热系数 ; 固液相变蓄能 ; 数学模型 ; 液相率 ; 界面移动速率 中图分类号 : T K512 + . 4 文献标识码 :A
0 引 言
用相变材料储存能量 ( 或蓄冷) 在下述场合下可 以大幅度节能或降低耗能费用 : 用电低谷时与用电 高峰时有较大的用电差价 ; 能量的供应因时间不同 有较大变化的太阳能 、 风能和余热等 。但是由于相 变材料的传热能力通常较差 , 所以需要通过优化蓄 能装置的结构研究来改善其蓄能和释能规律 。如果 完全依靠实验研究 , 需要花费大量时间 , 人力和财 力 。因而希望通过建立数学模型 , 运用数值计算来 研究蓄能结构的能量存储和释放规律 。
对于液相率 ,基于其本身的物理意义 ,任一节点 必须套用以下关系式 : k +1 k +1 if <P < 0 , <P = 0 <P
k +1
> 1 , <P
k +1
= 1
( 15 )
3 实验方法
对一维问题 ,为验证上述解法的正确性 ,我们设 计了一维恒壁温边界条件的实验装置 。实验装置如 图 2 所示 :1 为一 300 × 200 × 200mm 长方体外壁绝 热的容器 ; 在距容器右端面 100mm 处 , 内置导热隔 板 2 将容器分为 A 、 B 两部分 ;A 内装有相变材料 ;B 内装加热液体和内固定有电加热器的铝板 3 , 由于 铝板表面的温度均匀并且其截面与容器截面的大小 基本相同 ,所以 , 在实验过程中 , 加热液体的温度始 终是均匀的 。热量沿导热隔板由 B 内加热液体均
) , c = cl < + cs (1 - <) ,ρ = ρ l < +ρ s (1 - < ρ c ( 8) keff ( 7)
( 5)
由文献 [ 3 ] [ 4 ] ,得到有效导热系数的经验公式如下 :
1
X ( t)
0 . 25
if
l > 2 X ( t) Pr Ra Pr + 0 . 2
div V = 0
( 2)
( 3) 动量方程 : 9V + (V ・ grad) V = 9t 1 ( grad V ) + S V ρgrad P + div V ・
( 3)
1 固液相变蓄能的数学模型
固液相变问题 ( 又称 Stefan 问题 ) 的完备数学 模型如下 [ 1 ] : 1) 能量方程 : ρ9 h + ( V ・ grad h ) = div ( k ・ grad T ) + S h 9t 其中 :
k eff = N u kl N u = 0 . 22 Pr Ra Pr + 0 . 2
0 . 28
限制 ; 3) 热焓方程允许两相区存在一个熔解域 。 考虑为一维问题 , 将方程 ( 5 ) 转换为热焓表达 式 ,如下 : 9h 9 α9 h 9< = - ρ lL 9t 9X 9x 9t 其中 , α=
0 . 29
边界条件 : x = 0 , h ( 0 , t ) = ρ l cl ( T w - T m) ( 9) dh ( 10 ) x = l , = 0 dx 初始 条 件 : t = 0 , < ( x , 0) = 0 , h ( x , 0) = ρ ( 11 ) s c s ( T m - T init ) ( T init 为给定温度初场) 方程 ( 6) 求解可用全隐控制容积有限差分格式 离散 ,得方程如下 : 0 0 hP - hP h E - 2 h P + hW <P - <P ρ =α L l δ δ t t δx 2
求解上述方程组存在着下列问题 : 方程 ( 1) 中导 热系数 k 应用于液相或固液界面的计算并不准确 , 因为该处的传热不但有导热 ,还有对流换热 ,造成实 际的传热量为导热量的数倍之多[ 2 ] ; 实际的蓄能结 构有大量的传热和加强热交换的肋片 , 且方程组的 边界条件十分复杂 ,造成方程组的求解困难 。 考虑到研究上述问题时 , 我们关心的是蓄能结 构能量储存和释放的规律 , 而不是流场速度分布 。 并且 ,实际相变过程中 , 液相流动速度非常小 , 其影 响可忽略 , N u 数转化为用 R a 数和 Pr 数加以描 述 [ 3 ] [ 4 ] 。从分子热运动的观点分析 , 液相的浮力的 存在加速了分子的运动 , 促进了蓄能单元内温度梯 度的扩大 ; 同时液相的粘性和热扩散作用的存在阻 碍液体分子运动 。无量纲数 R a 很好的反应了浮力 跟流动中的粘性和热扩散作用之间的力量对比关 系 ,所以必须予以考虑 。 Pr 数则是反映了两种扩散 作用 — — — 粘性和热扩散之间的力量对比 。因此可将 含自然对流的 Stefan 问题简化为仅用能量方程加 以描述 ,速度对蓄能和释能的影响可以通过对导热 系数的修正来考虑 。 若假定相变物质材料的热物性不变 , 那么能量
A B S ( hP hE
k +2 k +2 k +1 k +2 k +2 kHale Waihona Puke Baidu+2
k
匀传入 A 内相变材料 ,使其慢慢熔解 ,产生相变 。 在 A 内上 、 中、 下三条沿传热方向上 , 各等间距 布有 20 个测点 , 以测定 A 内各点温度随时间的变 化关系 。实验共用两个铂电阻 。一个 ( Pt1 ) 置于环 境 ,测环境温度 ; 另一个 ( Pt2 ) 置于 B 内测量被加热 液体温度 。铂电阻采用四线法测量 。 实验中 , 温度测量所用的数采系统由数字万用 表 ,多路程控开关 ( 内置两块 40 通 道 扫 描 卡 ) 及 IEEE488 总线接口卡组成 。实验中 , 热电偶的冷端 与铂电阻均接入程控开关的端子板上 , 经选通后输 出 ,送到数字万用表的输入端 , 经数字万用表滤波 , 放大及 AD 转换后 , 送入计算机 。由于数字万用表 与程控开关的远程控制采用 IEEE488 通用仪器总 线 ,而 PC 机 采 用 ISA 总 线 , 所 以 需 要 插 入 一 块 IEEE488 总线接口卡 , 以实现计算机对数字万用表 和程控开关的远程控制 , 完成对热电偶和铂电阻测 温的数据采集 。 实验中 , 通过计算机监控系统测量加热液体的 温度并控制电加热器的温控开关 。用 Pt2 对 B 中被 加热液体的温度实时监控 ,当液温超过设定温度时 , 电加热器断开 ,即停止加热 ; 反之 ,电加热器闭合 ,开 始加热 。如此可保证液温恒定 , 实现等壁温导热板 壁面边界条件 。 由式 ( 5) 可知 , k eff 与 Δ T 有关 ,即与 ( T w - T m ) 有关 。实验中 ,我们分别将 Δ T 控制在 40 ℃,50 ℃,
[4 ]
图1 一维控制容积
Fig. 1 One2dimensional domain
求解全隐有限差分格式 ,需满足以下假设[ 4 ] : ( 1) 相变材料热物性与温度变化无关 。且固液 两相有其各自的热物性 ; ( 2) 相变材料最初处于固态 ; ( 3) 相变材料各向同性 ; ( 4) 所有相变子单元等价且彼此独立 。
