热辐射问题的边界元算法.pdf

2 0 1 1 年第l 期导弹与航天运载技术 N o 12 0 1 1 总第3 1 1 期M I S S I I E SA N DS P A C EV E t t l C L E SS u mN o 3 1 l 文章编号：1 0 0 4 7 1 8 2 ( 2 0 11 ) 0 1 0 0 4 6 - 0 8 热辐射问题的边界元算法 王静u ，高效伟2 ( 1 东南大学工程力学系，南京，2 1 0 0 9 62 大连理工大学航空航天学院，大连，1 1 6 0 2 4 ) 摘要：介绍一种处理介质热辐射问题的纯边界元算法该方法在漫射表面以及不考虑介质散射的条件下，只需要在 边界划分单元就可求解介质辐射换热问题利用立体角与微元表面的转换关系，可将沿立体角的积分转化到边界上，形 成边界积分方程在非等温介质问题中，边界积分内部还含有线积分，采用径向积分法计算此线积分二维和三维算例 证明，边界元法是计算辐射换热问题的一种有效数值方法 关键词：热辐射；参与性介质；漫射表面；边界元法 中图分类号：V 4 1文献标识码：A A B o u n d a r yE l e m e n tM e t h o df o rS o l v i n gH e a tR a d i a t i o nP r o b l e m s W a n gJ i n 9 1 一G a oX i a o w e i 2 ( 1 D e p a r t m e n to f E n g i n e e r i n gM e c h a n i c s ，S o u t h e a s tU n i v e r s i t y ，N a n j i n g ，2 1 0 0 9 6 ； 2 S c h o o lo fA e r o n a u t i c sa n dA s t r o n a u t i c s ，D a l i a nU n i v e r s i t yo fT e c h n o l o g y ，D a l i a n , 116 0 2 4 ) A b s t r a c t ：Ab o u n d a r ye l e m e n tm e t h o dw h i c hc a nb eu s e dt os o l v eh e a tr a d i a t i o nt r a n s f e rp r o b l e m sw i t h o u ti n t e r n a lc e l l si n t h ea s s u m p t i o n so fd i f f u s i v es u r f a c e sa n dn os c a t t e r i n gm e d i u m si sd e s c r i b e d T h er e l a t i o n s h i pb e t w e e nt h ei n f i n i t e s i m a ls o l i d a n g l ea n dt h ei n f i n i t e s i m a ls u r f a c ec a nb eu s e dt oc o n v e r tt h ei n t e g r a l so v e rt h es o l i da n g l ei n t oi n t e g r a l so v e rt h es u r f a c e b o u n d i n gt h em e d i u m ，a n db o u n d a r yi n t e g r a le q u a t i o n sc a nb eo b t a i n e d F o rt h en o n i s o t h e r m a lg r a ym e d i u m ，l i n ei n t e g r a l se x i s t i nt h eb o u n d a r yi n t e g r a l s ，w h i c hc a r lb ee v a l u a t e du s i n gr a d i a lb a s i sf u n c t i o n s 2 Da n d3 Dn u m e r i c a le x a m p l e sa r eg i v e nt o d e m o n s t r a t eab o u n d a r ye l e m e n tm e t h o di sv e r ye f f e c t i v et oc o m p u t er a d i a t i o nt r a n s f e rp r o b l e m K e yW o r d s ：H e a tr a d i a t i o n ；P a r t i c i p a t i n gm e d i u m ；D i f f u s i v es u r f a c e ；B o u n d a r ye l e m e n tm e t h o d O 引言 辐射传热计算是工业领域和科学研究领域中有关 温度场计算方面的重要组成部分。表面之间的辐射换 热研究在2 0 世纪6 0 年代已经相对成熟，之后偏重研 究具有发射、吸收、散射介质的换热计算。发射、吸 收和散射介质的辐射与固体的表面辐射相比要复杂得 多。一方面由于介质内部参与了辐射，在研究辐射换 热问题时需要了解介质中每个局部位置上介质的温 度、辐射强度和物理性质：另方面，介质辐射中波 长的影响也十分显著，需要详细地研究其光谱辐射的 性质。对介质辐射的研究最早主要是针对地球大气层 中的辐射、吸收和散射。在2 0 世纪2 0 年代，工业上 气体辐射的研究主要是与燃烧过程中的传热问题相联 系的，发现燃烧物中的二氧化碳、水蒸汽、一氧化碳 等气体对辐射都具有一定的吸收与发射作用。核动力、 火箭技术、等离子发生器及高温高压燃烧室等科学技 术的发展更推动了对介质辐射性质及传热的研刭。 热辐射的产生与传输机理与导热、对流换热的不 同，导致描述它们的控制方程有很大的差异。由于方 程不同，因而导热、对流换热数值计算中一系列行之 有效的方法都不适用于辐射传热的数值计算1 2 J 。目前求 解辐射传热问题的数值方法有：扩散近似法、热流法、 蒙特卡洛法、有限体积法等。a ) 扩散近似法：将光子 经历看成一种类似于分子输运中的扩散过程，从而把 积分微分形式辐射传递方程转变成形如导热问题的 扩散方程，这样就可以用求解热传导方程的常用方法 来求解：b ) 热流法：将微元体界面上复杂的半球空间 热辐射简化成垂直于此界面的均匀强度或热流，使积 分微分形式的辐射传递方程简化为一组有关辐射强 度或热流密度的线性微分方程，然后用通用的输运方 程求解方法求解p J ；c ) 蒙特卡洛法：将传输过程分解 为发射、投射、反射、吸收和散射等一系列的予过程， 收稿口期：2 0 0 9 - 0 9 1 4 ：修刨日期：2 0 1 0 - 0 1 - 0 8 作者简介：王静( 1 9 8 3 ) 女。博士研究生。主要研究方向为飞行器热防护系统热分折 万方数据 第1 期 王静锋热辐射题的边界兀筇法 4 7 并把它们化成随机问题，即建立每个子过程的概率模 型，令每个单元( 面元或体元) 发射一定量的光束。跟踪、 统计每柬光束的归宿( 被介质和界面吸收，或从系统中 投射出或溢出) ，从而得到该单元辐射能量分配的统计 结果。2 0 世纪4 0 年代中期，蒙特卡洛法首先在核武器 的研制中得到应用，1 9 6 4 年，H o w e l l 将其引入辐射换 热计算领域l M 】；d ) 有限体积法：是R a i t h b y 和c h u i 首先提出的1 3 , 7 ，此后，不少研究学者用该方法处理了 各种热辐射传输与耦合换热等问题p ，有限体积法 的基本思想是在控制体积和控制立体角内对辐射传递 方程进行积分，可得控制体积和控制立体角内辐射能 量守恒方程的有限体积表达式。扩散近似法和热流法 由于对控制方程进行近似，因此精度较低：蒙特卡洛 法由于需要进行繁重的光束传播轨迹计算，因此计算 速度极慢，而有限体积法虽然较蒙特卡洛法计算速度 要快，但是它们既要划分体单元又要划分面单元，因 此在建模和计算方面工作量较大。 边界元法( B E M ) 是一种比较新颖的数值分析方 法，所形成的积分方程通常只在区域的边界上，因此 只需要将边界离散，相当于把问题的维数降低了一维， 这样一来在数据准备和计算方面都大大节约了时间。 然而，目前所发表的文献中，利用B E M 求解参与性介 质辐射换热问题时，最终形成的B E M 边界积分方程由 于在面积分中含有沿辐射射线方向的线积分，因此并 没有将积分完全转化到边界上，仍然需要将计算区域 划分成体单元【1 1 1 。本文利用径向基函数计算边界积分 方程中的线积分，从而形成不需要内部体单元的热辐 射边界元法。在透明介质和等温灰体介质条件下，仅 仅通过微元立体角与微元表面的转换关系式就可以将 积分转化到边界上，形成边界积分。对于非等温灰体 介质，所得边界积分中存在体现介质吸收对辐射能减 弱的线积分，这个线积分在发表的边界元文献中是通 过划分体单元来求解的，而本文采用径向基函数来近 似，从而不需要划分体单元。对于非灰体介质，采用 谱带近似的方法，对每一个谱带列方程求解即可。 l B E M 形式的辐射积分方程 辐射表面和介质都为灰体，且辐射表面为漫射表 面，忽略介质散射假设下，关于辐射强度的积分形式 辐射传递方程为【1 1 , 1 2 1 州力篇篇) 小却刚) ( 1 )L 屯( r ) 厶 P ( ，) f ( 一p ) ( r ) 式中 妒为辐射传播路径与源点P 处外法线方向的夹 角；，：( p ) 为源点P 处妒方向上投射辐射强度：r 为沿妒 方向投射到源点P 的边界表面上的点；办为点，处外 法线方向与辐射传播路径之间的夹角；C ( ，) 为点，处 沿办的有效辐射强度；k 为介质的吸收系数：k 为点 r 和点P 之间射线所经过的路径：，为该路径上的点： 厶 r ( r ) 为点，。处所在温度下黑体辐射强度；f ( ，P ) 和r ( r ；p ) 为定向穿透率，表示经过介质吸收后的辐射 强度与初始辐射强度的比例系数。 r ( r 。P ) 和f ( ，：p ) 表达式为 。 厂1 f 【，p ) = c x p 卜Ik ( r ) ( r ) l ( 2 ) 广1 f ( ，- ：p ) = e x p l 一，k a ( r ) 吆( ，。) I 3 ) 点P ，以及，- 之间的关系如图1 所示。 圈l辐射方程中积分方向上各点及相互关系 式( 1 ) 的物理意义为：从点，发射的辐射能，经 过传播路径上介质吸收和介质发射辐射所到达点P 的 能量。 对于漫射灰体表面，有式( 4 ) 所示的热平衡关系 式，物理意义为：半球空间内的有效辐射等于发射辐 射与反射辐射之和。 L ，；c o s 妒Q = L 厶c o s 阳Q + p L 0 c o s a t o ( 4 ) 式中d Q 为微元立体角：为表面发射率；P 为表面 反射率，并且根据基尔霍夫定律有p = 1 一。 灰体表面的净辐射热通量可以表示为被吸收的投 射辐射与发射辐射之差( 半球空间内) ： q 7 = I 。( 以一l 。) c o s 舳 ( 5 ) 利用式( 4 ) 和式( 5 ) 消去，二，可得如下关系式： I 。，；咖脚= L c o s 妒d j O + 半矿 ( 6 ) 根据漫射灰体表面性质有： 万方数据 4 8 导弹与航天运载技术2 0 1 1 年 L 枷蝴= 彬 茎苫耋霎翥嚣，积分表达式利用传统的边界元 并且，黑体辐射力毛和黑体辐射强度厶之间有如 下关系式【1 2 ，1 3 】： 厶：三毛 ( 8 ) 丌 于是，式( 6 ) 可转化为 巧= 三( 毛+ 字