（应用数学专业论文）基于边界元的高炉侵蚀线反演.pdf

复旦大学硕士学位论文 摘要 在高炉炼铁过程中由于化学反应以及热一机械压力，炉壁会受到一定程度的 腐蚀，我们需要监测高炉内壁侵蚀情况以保障生产安全。因而在实际生产中， 通常是在炉内的某些位置埋放热电偶，通过测量出的温度值来推测腐蚀情况。 在这篇论文中我们考虑了轴对称稳态热传导反问题高炉侵蚀线监测问 题。本文利用误差平方和最小的原则把问题归结为一个优化问题。并采用边界 元方法计算正问题。边界元方法不需要进行网格剖分，只需对边界进行离散， 有效的减少了计算量。另外为避免大量重复计算我们引进了摄动方法。 由于反问题本身的不适定性，因此在测量数据可能有误差的情形，我们采 用了吉洪诺夫正则化方法。 数值结果表明本方发是稳定有效的。 关键字：边界元，基本解，摄动方法，正则化方法。 复旦大学硕士学位论文 a b s t r a c t i nt h i s p a p e rw ec o n s i d e ra ni n v e r s eb o u n d a r yp r o b l e mf o rt h ea x i s y m m e t f i c s t e a d y s t a t eh e a te q u a t i o n ，w h i c ha r i s e si nm o n i t o r i n gt h eb o u n d a r yc o r r o s i o nf o rt h e b l a s t f u r n a c e m e a s u r e t e m p e r a t u r ea t s o m el o c a t i o n sa r eu s e dt o i d e n t i f yt h e i n t e r f a c eb e t w e e nt h el i q u i di r o na n dt h ef u r n a c el i n i n g t h e p a p e rg i v e sam e t h o dt os o l v et h ei n v e r s ep r o b l e mb ym i n i m i z i n gt h es u mo f a b s o l u t er e s i d u a la tt h em e a s u r i n gl o c a t i o n s ，i nw h i c ht h ed i r e c tp r o b l e m sa r es o l v e d b yb o u n d a r ye l e m e n tm e t h o d ( b e m ) f u n d a m e n t a ls o l u t i o nf o rt h ea x i s y m m e t r i c p r o b l e mi su s e di nb e m w h i l et h eb o u n d a r yi nt h eh e a tp r o b l e mi sn o tf i x e d ，t h en u m e r i c a li n v e r s i o ni s c o m p l i c a t e da n dc o n s u m i n g p e r t u r b a f i o nm e t h o di s s u g g e s t e d t ot r a n s f e rt h e d y n a m i cb o u n d a r yi n t ot h ef i x e db o u n d a r yp l u sas m a l lp e r t u r b a t i o n b e c a u s et h ei n v e r s ep r o b l e mi si l l p o s e d ，t h et i k h o n o vm e t h o di ss u g g e s t e d t h ea c t u a lc a l c u l a t i o ns h o w st h a tag o o da g r e e m e n tb e t w e e nt h ec o m p u t e da n d m e a s u r e dr e s u l t si so b t a i n e d k e y w o r d s ：b e m ，f u n d a m e n t a ls o l u t o n ，p e r t u r b a t i o n ，t i k h o n o vm e t h o d 复旦大学硕士学位论文 第 1 章引言 1 1 高炉炉缸炉底的侵蚀监测问题 冶金高炉在使用的过程中，其内壁，尤其是炉缸和炉底，长期受到高温、 渣铁冲刷和化学作用等等的侵蚀，容易发生破损。当炉缸、炉底所受的侵蚀达 到一定的程度时，高炉就可能需要停产大修，以保障生产安全。因此，为了延 长高炉寿命，保障生产安全，提高经济效益，监测高炉侵蚀线就成了一项很重 要的课题。 1 2 问题的研究现状和主要结果 早在上世纪七十年代，国外就已经开始高炉炉底侵蚀线监测方面的研究和 开发，并陆续应用于生产实践中，取得了一定的效果，进一步的研究正在持续 进行中。例如，奥地利林茨大学建立了这方面的模型。日本新日铁2 0 0 0 年还专 门立项，与东京大学合作开展研究。 我国在这方面的研究起步较晚，但最近也有不少这方面的研究成果。例如： 1 9 9 6 年， 9 刘薇等通过求解控制高炉炉底传热过程的热传导方程，解得在假想 侵蚀边界条件下的有限个边界点的温度值，再通过正交试验确定满足实测边界温 度的侵蚀边界。 1 0 庄玉如等提出一种人工调节侵蚀边界的假设位置、使用有限 元方法计算炉衬温度分布，推断炉衬侵蚀状态的方法； 7 】在1 9 9 9 年利用有限 元法建立了二维数学模型，采用网格变形拟合的方法处理移动边界问题： 1 9 9 6 年的 9 1 5 n2 0 0 4 年的 6 都使用边界元方法计算温度分布，用正交实验进行反演 推测侵蚀边界，但在边界条件的离散处理上，前者使用常数元离散，后者采用基 尔霍夫变换把非线性问题转化为线性问题。 目前国内大部分使用的模型，都将问题简化成了二维的平面传热问题，很 少有使用三维模型的。计算温度分布时使用的数值算法有差分法、边界元法、 有限元法等；而在反演时，一般使用正交实验、神经网络等方法。 本文对炉衬温度分布的计算使用的是三维传热模型。我们认为，炉衬热传 导问题是三维问题，虽然它的方程可以利用炉衬的旋转对称性，转化成柱坐标 复旦大学硕士学位论文 系从而实现简化，但是如果将整个传热问题本身简化为平面热传导问题，则是 一种过度简化，对计算结果的精确性会造成较大的影响。 本文使用边界元方法计算正问题，推导了轴对称拉普拉斯方程的基本解。 另外为了有效的避免大量重复计算，提高运算效率，我们引进了摄动方法。并 运用拟牛顿优化算法进行优化问题的求解。值得指出摄动方法的使用是有条件 的，只有在扰动极小的情况下才能保证计算精度。 1 3 本文的主要内容 高炉在建造的过程中，为了监测炉衬的侵蚀情况，会在炉衬中埋设热电偶， 探测炉衬中的温度。本文将确定炉衬侵蚀线的问题归结为一个稳态热传导方程 的确定边界的反问题，即，求该方程求解区域的一条边界，使得稳态热传导方程 在此边界所界定的区域上的解在热电偶所在位置，与温度的实际测量值一致。 这一问题又可化为求稳态热传导方程的定解问题的一段边界，使该定解问题的 解在测温点的值与温度测量值的误差极小。离散化后，这一问题成为一个优化 问题。 另外，由于利用优化方法反演移动边界，需要进行大量稳态热传导方程边 值问题( 正问题) 的计算，这些正问题的边界尤其是对应于侵蚀线的部分 是不同的。这就需要重复的进行节点配置，矩阵组装，求解代数方程组，这是 相当耗时繁琐的。而且其中大量的计算工作是重复工作，为了简化问题的计算， 我们引进了摄动方法，摄动方法巧妙的把变边界问题转化为两个固定边界的边 值问题的求解，减少了重复计算。但是摄动方法也是有条件的，它只有在扰动 较小或者线性性较强的情况下才是精确高效的。 考虑到反问题通常是不稳定的，我们采用了吉洪诺夫正则化的方法，选取 了适当的正则化项、确定了正则化参数，增加了算法的稳定性。 本文分三个章节，第一章引言，介绍本文选题的意义，和在这个问题中的 研究现状；第二章介绍高炉炉衬热传导问题的模型，和本文的解决途径，包括 方程的推导，和我们使用边界元方法求解的过程，并较为详细的介绍了边界元 复旦大学硕士学位论文 方法和整个计算过程；第三章描述本文解决炉衬侵蚀状态监测的方法，包括炉 衬侵蚀问题的描述和模型的建立、摄动方法的介绍、以及本文的数值结果。 复旦大学硕士学位论文 第2 章高炉炉衬热传导问题的模型和正问题计算方法 宝钢2 号高炉呈旋转对称结构，为简便计，我们采用柱坐标系画出宝钢2 号高炉副丁午面。如下图所示。 6 5 d 3 2 0 图2 1 宝钢2 号高炉半截面图 整个区域分成两块，是由不同传热材料组成。其中区域q 的热传导系数 为3 3 w m k ，区域q 2 的热传导系数为l o w m k 。 2 1 正问题的描述和模型的建立 高炉炉衬的热传导问题是三维区域中的问题，高炉及外围冷却设施都呈中 心旋转对称。下面为表述方便，采用柱坐标系画出求解区域的简单示意图，如 下图所示。 4 复旦大学硕士学位论文 图2 2 计算区域示意图 考虑到高炉内部的侵蚀过程，相对于一代高炉的使用寿命而言非常缓慢 热传导方程中的时间项也可以简化略去。 所以，高炉炉衬的热传导问题的方程可以写成如下形式： 导( 北詈 + 未( t 老 = 。m q c z - ， 其中，和z 分别是柱坐标系中的半径和竖直方向的坐标，七是炉衬材料的热 传导系数，“表示温度。 在对称边界f 处，热流为0 ，故“满足绝热边界条件： 娑= o d n r i ； ( 2 2 ) d n 在炉衬的底部边界r 和外侧边界r 1 处，分别有风冷却和水冷却措施，所以 “满足热交换边界条件： 七芸= 嚏( 一“) 册r 2 ， d n k 竺= h w ( u - u ) 。nf 3 ； d r ( 2 3 ) ( 2 4 ) 其中吃是边界r ：与冷空气之间的热交换系数，丸是边界r 3 与冷水之间的 热交换系数。“。，“。是相应冷却物质的温度。 在顶部边界l 处，热交换可以忽略，因此“满足绝热条件，即 复旦大学硕士学位论文 挈：0 o n f 。； ( 2 5 ) 最后，u 在侵蚀边界r 5 上，实际上是一条等温线，等于熔化铁水的温度丁。， 因此有 u = t o o n f 5 ；u = t o o i l f 5 ； ( 2 6 ) 本文中r o = 1 1 5 0 。c 。 至此，我们建立了高炉炉底炉衬热传导的数学模型，这是一个轴对称的稳 态热传导模型。在此模型中，若已知方程的系数、区域的几何形状和边界条件， 就可以求得炉衬的温度分布，这就是本文的正问题。 实际上，在对特定的高炉进行计算时，方程中的系数和边界条件总是确定 的，而区域的边界除侵蚀边界外也是己知的，所以，它们都可以看做正问题的 常量参数。这样一来，只要设定了侵蚀线的位置，就可以对正问题进行计算。 2 2 数值计算方法简介 关于稳态热传导问题的数值求解，目前已经发展了相当多的数值求解方法， 比如有限差分法( f d m ) ，有限元法( f e m ) ，边界元法( b e m ) 。 有限差分法是一种很古老的数值计算方法，理论比较完整。它直接从微分 方程出发，将求解区域划分成网格，近似的用差分、差商代替微分、微商，于 是无限自由度的问题化成了有限自由度问题。可以通过代数方程组进行求解。 这种方法在解决规则边界问题时极为方便，但是对非规则边界的问题适用性较 差。 有限元方法是从微分方程所对应的泛函出发，用变分原理结合区域剖分得 到离散算式代数方程组。它克服了有限差分法对区域形状的限制，对各种 形状的边界都能灵活处理，有限元法具有适应性强，应用范围广的特点，而且 在工程应用中己经开发了相当多的有限元软件，如s a p ，n a s t r a n ，a s k a ， a n s y s 等等，m a t l a b 中也已经有了有限元工具箱，使用方便，是目前工程计 算的主要手段。有限元方法计算的主要困难在于：由于区域的剖分随着网格的 复旦大学硕士学位论文 加细而使方程组的维数急速增大，计算负荷较大。尤其是对于高维问题，计算 量相当庞大。 边界元法是从上世纪七十年代发展起来的数值计算方法。相对有限差分法、 有限元法，是种较新的数值计算方法。这种方法是用控制微分方程的基本解 建立相应的边界积分方程，再对它结合边界的剖分而得到离散算式。由于只在 边界上剖分，因此实际上是将问题降维处理，降维的结果必然减少代数方程组 的未知数。从而减少存储空间的占用，节省计算时间。另外，与有限元法稍有 不同，边界元法只对边界离散，因此计算误差只来源于边界，这就使我们将减 少误差的全部注意力放在边界上。求域内变量时，只需改变其数量和坐标位置 即可。对于只关心边界上的函数值和梯度值的问题有其优越性。 总的来说，三种数值计算方法各有优劣。考虑到本文研究的热传导问题， 不需要了解整个区域的温度场，只需要知道边界上的温度和外法相导数即可。 这正好能够发挥边界元方法的优点。我们选用边界元作为稳态热传导问题正问 题的数值求解方法。 2 3 边界元法介绍 边界元法是在经典的积分方程理论的基础上，吸收了有限元法的离散技术 而发展起来的计算方法。边界元法的数学基础是积分方程理论，早在1 0 0 多年 前，阿贝尔( a b e l ) 、亥姆霍兹( h e l m h o l t z ) 等就对积分方程理论进行了深入 的研究。弗雷德霍姆( f r e d h o l m ) 在对积分方程研究的基础上，首先将这种理 论用于弹性力学问题的求解。2 0 世纪5 0 年代，前苏联的米赫林( m i k h l i n ) 、 穆什海里什维里( m u s k h e l i s v i l i ) 的工作为积分方程在工程上的应用开辟了道路。 从边界元法计算格式形成的全过程来看，关键问题有两个。个是问题的 边界化，即将给定区域上的定解问题化为可以只考虑边界取值的问题。这一步 的关键是格林( g r e e n ) 公式，这是边界元法的基石。边界化的结果使问题降维。 边界元法的封闭方程组中只有边界节点上的未知量，解完方程组再根据需要有 目的地求内部值，因此减少了计算的盲目性。形成计算格式的第二个关键步骤 是边界的离散化。单就离散技术而言，与有限元法没有很大的不同。 复旦大学硕士学位论文 1 边界积分方程 2 3 1 边界积分关系式 用边界元方法来求偏微分方程的数值解，首先需要把偏微分方程转化为边 界积分方程，称之为边界归化。同一边值问题可以用不同的方式得到几个不同 形式的边界积分方程。 考虑到本文研究的问题是一个稳态热传导问题。区域是轴对称的。温度函 数满足一个轴对称拉普拉斯方程。本质上是一个三维拉普拉斯方程。我们首先 考虑如下个三维边值问题。 a u ( x 、= 0 石q “( x ) i r o = ，( z ) 石f d i b u ( x ) l ，。= 删x r 。 ( 2 7 ) 其中，和g 都是已知函数。 利用第二格林公式( v a u - u a v ) d y = ( ”芸一“妾) 搬和三维拉普拉斯方程 nr 的基本解可以把问题( 2 7 ) 表述成边界积分方程。 其结论可归结为如下两个定理。 定理2 1 设“c 2 ( c z ) l q c l ( q ) ，则有 s ( u 芸n n 担- i n ) d s 一，* a u d v = 茹晏 s ， 其中* ( z ，y ) 是三维拉普拉斯方程的基本解。三维拉普拉斯方程的基本解将 在下文介绍。五。表示三维区域孬的补集。 定理2 2 若u c 2 ( q ) n c ( q ) 是调和函数，则当y r 复旦大学硕士学位论文 c ( 咖( y ) ：弦( 训) 旦掣一吣) 旦蜜型) 峨 ( 2 9 ) ： o no n 其中c ( y ) ：墨盟，o ( y ) 是过y 点的1 1 的切平面围成的三维角( 包含在q 的 4 石 部分) 。 上述两个定理的证明过程可利用调和函数的性质证得。 2 三维拉普拉斯方程的基本解 把物理现象所对应的微分方程转化为积分方程，基本解起着重要作用。下 面我们不加阐述的引进基本解的概念以及三维拉普拉斯方程的基本解。 首先简要的介绍与基本解密切相关的万函数。 所谓j 函数是定义在r 。( d 表示维数) 中并且具有以下性质的函数： ( f ) 8 ( x - y ) ：j o z y o o z=v ( i f ) j j ( x y ) d x = l 占函数的个重要特征是对任意一个连续函数妒( x ) ，都有 p ( x ) 万 一y ) 出= 妒( y ) ，y q 当研究某个线性微分算子l 或相应的线性微分方程时，称满足方程 l u = 8 ( x - y ) 的解h + ( 石，y ) 为算子l 或相应的微分方程的基本解。 可以证明三维拉普拉斯方程的基本解是： “4 ( 五y ) 2 碲1 习 2 1 。) 在下文中如无特别说明，“( z ，y ) 表示三维拉普拉斯方程的基本解。 2 3 2 轴对称拉普拉斯方程的基本解 在本文研究的稳态热传导问题中，温度函数u ( r ，z ) 满足轴对称拉普拉 9 复旦大学硕士学位论文 斯方程h ：0 ，其中l ：一i _ 0 ( r 拿) + 芝以下如无特别说明，l 表示轴对称拉 r0 r0 rd z 普拉斯算子。只要找出了轴对称热传导方程的基本解“+ ( 工) ，) 满足 l “( z ，y ) = 万 一y ) ，就可以把区域内的任意一点的函数值表示成边界积分的形 式。我们从三维拉普拉斯万干呈出发推导轴对称挝晋拉劢万程明基本前年。 由前面定理2 2 ，我们得到三维拉普拉斯方程的积分关系式( 2 9 ) 。 对轴对称拉普拉斯方程使用式( 2 9 ) 。 c ( ) ，m 炉肼“伽) 等叫神塑字肌脚( 圳f = f 南删掣肌，一聘硐1 脚m 剐瑚m r 。南棚( 曲 = f 南坝曲 = i 2 ；= = = = = = = = = = = = = = = = = = = 二= = = = = = = = = = = = = = ；= = = 一d o ( x ) 由石r 2 ( 工) + 尺2 ( y ) - 2 r ( x ) r ( y ) c o s ( o ( x ) 一目( y ) ) + ( z ( x ) 一z ( y ) ) 2 ：1 2 ；：；：；：；：一d 妒 ( 2 1 2 ) 如万( r ( 工) + r ( y ) ) 2 + ( z ( z ) 一z ( y ) ) 2 - 4 r ( x ) r ( y ) s i n 2 妒 其中d l 是沿轴截面的边界l t t 线积分。 式( 2 1 2 ) 就是轴对称拉普拉斯方程的基本解。 对轴对称区域的光滑边界有： c = 等= 等 其中口( y ) 是立体角，当y 在q 内部时目( y ) = 4 z c ， 当y 处于q 的光滑边界时p ( y ) = 2 石。 占( ) ) 是轴截面上的平面角，当y 在q 内部时谷( y ) = 2 万， 当v 处于q 的光滑边界时占( v ) = 7 1 。 1 0 复旦大学硕士学位论文 这样就可以把积分式( 2 1 1 ) 看成是沿轴截面的边界的积分。 2 3 3 边界元离散、求解 由( 2 1 2 ) 式我们得到轴对称拉普拉斯方程的基本解。为了表述习惯 我们用点s 代替点x ，点f 代替点ya 把轴对称拉晋拉斯方程的基本解重新表述 如下： n 雠卜顽丽菰挚篙乒丽 他脚 牙+(s，。=攀3丽丽丽10n丽，f ( 尺( f ) + r ( s ) ) 2 + ( z ( 舌) 一z ( s ) ) 2 t 志c 鬻黑群襄学黯跏m c 帆 眨 2 r ( s ) 。( r ( 孝) 一r ( s ) ) 2 + ( z ( 善) 一z ( s ) ) 2 、“、 + 而订黑害黑i 两跏) ( s ) ) ( r ( 孝) 一r ( s ) ) 2 + ( z ( 亭) 一z ( s ) ) 2 l v “以、。川 其中m = 石瓦可五专算等笺赫，k ( m ) ，e ( m ) 分别表示第一类、第 二类完全椭圆积分。 由于高炉炉底炉衬在各个部位的材质是不完全相同的，例如，在不需要有 良好的导热性能，不易渣化的部位可以使用粘土砖或高铝砖。在需要有良好的 导热性能的部位可以使用碳砖。所以要把炉底区域划分成若干个单一材质的部 分区域，对各个部分区域应用上述的方法建立积分方程组。再用边界元法把它 们离散化为若干个代数方程组联立起来求解。 这样需要在各子区域的连接处补充边界条件。我们根据热传导的物理特性 补充如下边界条件。 复旦大学硕士学位论文 图2 3 分片子区域示意图 f 。：是子区域q ，、q ：的共同边界，我们在f ：上补充边界条件： i“f 吨2 “l r 隅驴如氛 其中k l , k 2 分别是区域q 。，q ：的热传导系数。n 表示区域q 。在r 。：上的外 法线方向。 为了进行数值计算，必须对边界进行离散化。为了叙述方便和简单明了， 我们假定炉底区域是同一材质的情况说明离散化方法。对于分片予区域的情形 可以同样推出。 爿成等) + 未 要= 。 ko 、 = 吃( 一“) d n o n k d 、u = ( u - u ) o n i d u = 0 o i l f 4 u = t o o n f ； 根据不同的要求及具体情况，将边界分割成一些小单元。可以使用常单元、 线单元、二次单元及更高阶单元。在单元上所取的插值函数应分别是常数、线 c ；no 、厂r r 蔷一破 。卜d )b2( 巳 l 复旦大学硕士学位论文 性函数、二次函数及更高阶函数。本文选用常单元，常插值函数。 把边界r ，f ：，f 3 ，f 。，f 5 分别划分成。，n 2 一n l ，n 3 一n 2 ，n 4 一n 3 ，n 5 一n 4 ( - 单- 元。这些小单元称为边界元素，记为e l 。在边界元素上把“和热通量日看作常 数，记为“，q ，并用它们来代表边界元素上的节点( 元素中点) 上的值。此外， 把全部边界元素的节点取作毒。这样积分方程就可以写成如下形式的方程组 c ，“。= q “7 r ( s ) d f 。一p 西尺( s ) d r 。 = 芝目“净( s ，) d f e s 一艺“f q ；r ( s ，) d f e s i 。1 t e j 。t r e i 铮( c + h ) u = g q 式中 “? = h ( s ，毒) ， g = g 1 。g l ： g 2 。g 2 ： 吒z 吒： u = l “lh 2 c = 。d “ 吼2 蔷 c e g l 地 g 2 n g n 一： c n , h = 日1 2 一 日2 2 “儿 7 ，q = q q ： - 日虬 7 两 = f q ；d f e j i = 1 ，2 ，n 5 ，j = l ，2 ，n 5 f e i g = u ：d f e j i = l ，2 ，n 5 ，j = l ，2 ，m f e j ( 2 1 6 ) 另外根据问题中给出的边界条件有下述等式 q = a o + ( 2 1 7 ) 式中 a i = a = o o 如= 厶= 0 复旦大学硕士学位论文 。 a d = 吃 k o 丸 k 0 1 一l ，心一 虬一n 2 ，n 3 - n 2 a s = e 肾n s 吨是ns n 4 除单位晦。 u =8 = i 4 n i + i u 划kr 1 4 0 n 一h 1 n 一h 卜仁 吃一女 一 吃一女 一 k 一0 一 k 一 屈屈屈厦屈玑玑玑叽易 ! | 饥 ，。l | i g 复旦大学硕士学位论文 f 0 属= i l ，屈= l oj m 吃“。 k h a u 。 k 把( 2 1 7 ) 式代入( 2 1 6 ) 式得 ( c + f 1 ) u = g ( a u + ) 数。 屈= h w u 。 k h w u 。 k 0 f 0 肛h a , 4 - 肛h 地 ( 2 1 8 ) 求解线性代数方程组( 2 1 8 ) ，即可得到0 ，边界上的函数值和外法相导 解决奇性的两种尝试 在用边界元法求边值问题时往往会遇到奇性的困难。从式( 2 1 3 ) 可以看 出。弘m 积分p 肘，删) 硼中的积分鲥( m ，2 东黯r “v 、“ 在m ：a 是- - + n a 。本文尝试了两种方法解决这个奇性问题。 其中n = r 2 ( 肘) + 尺2 ( 肘。) + ( z ( m ) 一z ( m 1 ) ) 2 , b = 2 r ( 肘) ( m a m = 羔 足( m ) 2 皓是第一类完全椭圆积分。 一种方法是：由于该奇点是可积奇点，所以可以用数值逼近的方法来做积 分的近似值。即挖去m 。的某个微小邻域仃( m ，占) 。用积分fh + ( m ，m ，) d s r i ，口抽 来做为积分“( 肼，m ，) d s 的近似值，从数值实验来看，是可行的。 r 。 另外一种方法，考虑到温度场的连续性，我们用m 。某个邻域内的外点巧来 代替点m 。这样我们计算得到的是在外点的温度值和法相导数值 “( 瓦) ，芸“( 瓦) 作为h ( m 1 ) o ( m ，) 的近似。 d h 复旦大学硕士学位论文 “+ c m ，磁，塑笋船= ，c 肘，未“+ c m ，巧，豳 本文选用了前一种方法来进行数值计算。实验结果表明效果是理想的。 2 4 边界元计算结果和精确解的比较 为了验证边界元方法的有效性，我们用边界元法试算了一个有精确解的边 值问题，计算结果显示精确度较高。 幽方程组： 挈 ：o 伽 o , 1 1 o ，1 d 己 o n 0 x 0 ，1 】 ( 2 1 9 ) o n l x 0 ，1 _ d u ：0d n o ，1 x 0 d z u ( r ，z ) = r 2 2o n 0 ，1 x l 它的个精确解是“= r 2 2 z2 ，其中，r = 4 x 2 + y 2 。 对其使用边界元求解，得到解的计算值，与精确解在节点上的值做比较。 这个算例中，我们运用边界元方法求出了边界上的值和外法向导数。在边 界0 x 0 ，1 ，1 x 0 ，1 】， o ，1 x o 上由于给出了n e u m a n n 条件，我们计算出了函数值， 另外在边界【0 ，l l x l 上由于给出了d i f i c h l e t 条件，我们计算出了外法向导数。计 算结果平均相对误差为0 0 0 3 7 。这里为节约篇幅，这里不把全部节点处的比较 数据列出，仅以部分节点为例，列出如下表 坐标理论值 计算值 误差相对误差 ( o ，0 9 5 ) 1 8 0 51 8 0 5 40 0 0 0 40 0 0 0 2 ( o ，0 7 5 ) 1 1 2 51 1 2 8o 0 0 3o 0 0 2 7 ( o ，0 5 5 ) 0 6 0 50 6 1 0 20 0 0 5 20 0 0 8 6 ，l飞0 2 a 一砖 = = 卜孔一扣孔一打 孔一甜-0 井1 一r 复旦大学硕士学位论文 ( o 0 3 5 ) 0 2 4 5一o 2 5 1 60 0 0 6 60 0 2 7 1 ( 0 3 5 ，0 ) o 1 2 2 50 1 1 4 4o 0 0 1 9o 0 1 5 l ( 0 7 5 ，0 ) 0 5 6 2 50 5 5 1 90 0 1 0 60 0 1 8 9 ( 1 ，0 2 5 ) 0 8 7 50 8 6 4 10 0 1 0 9o 0 1 2 5 ( 1 ，0 6 5 ) 0 1 5 50 1 5 4 10 0 0 0 90 0 0 5 8 ( 1 ，0 8 5 ) 一0 4 4 5一0 4 4 60 0 0 10 0 0 2 2 ( 0 8 5 ，1 ) 一4 3 9 0 7 70 0 9 2 30 0 2 3 1 ( 0 6 5 ，1 ) 一43 9 7 2 20 0 2 7 80 0 0 7 ( o 4 5 ，1 ) 一4 3 9 8 l0 0 1 9一o 0 0 4 8 ( o 2 5 ，1 ) 一43 9 8 4 4o 0 1 5 60 0 0 3 9 ( 0 0 5 ，1 ) 4 3 9 8 5 40 0 1 4 60 0 0 3 7 表2 1 边界元算例，计算结果比较 第3 章高炉炉缸侵蚀线监测问题的求解 在本文中，高炉炉衬侵蚀线的监测问题是作为高炉炉衬热传导问题的反问 题来求解的。本文用摄动方法结合优化方法来解决这个问题。 3 1 反问题的描述和模型的建立 炉衬侵蚀监测问题的目标是确定侵蚀边界的位置，以便指导高炉冶金生产。 为此，我们需要借助炉衬中预设的热电耦的信息。如果在炉衬中所考察的截面 中，沿半径方向埋设了个热电耦，那么可以读取这些位置上的温度测量值 五，7 2 ，一五，将他们与我们用正问题计算得到的相应位置上的温度计算值 “，“：“。作比较，将使得两者的参数取到最小值时的侵蚀边界形状作为实际侵 蚀线的推测位置。 采用图2 2 中的记号，我们的侵蚀监测问题就是：求l ，使得对应的炉衬 温度函数在测量点的值恰等于热电偶测得的温度值。由于这个温度函数与r 有 1 7 复旦大学硕士学位论文 关，记为“= u ( r ，z ；r 5 ) ，于是确定r s 的问题归结为：求( “( r ，z ；f 5 ) ，r 5 ) 使得 u ( r ，z ；r 5 ) 满足偏微分方程的边值问题( 2 1 ) 一( 2 6 ) ，且成立 简记 “( o ，z ，；r 5 ) = u ，( j = 1 2 ，) “( r 5 ) = “( o ，z ，；r 5 l ( j = 1 2 ，) 反问题就转化为一个优化问题 呼n n0 ，一uj ) 2 5 j = l 其中“i 是边值问题( 2 1 ) 一( 2 6 ) 的解在测温点的值。 由于反问题的不稳定性，温度测量的微小误差可能引起r 5 的巨大误差，因 此，我们采用吉洪诺夫正则化的方法，引入正则化项r ，优化问题成为： 1 1 1 i n ，( r 5 ) = 哩n 善o ，( r 5 ) 一u ，) 2 + 尺( r 5 ) ( 3 1 ) 为了求解优化问题( 3 1 ) ，我们还需将侵蚀线r 5 离散化。采用本文正问题部 分引入的离散方法，给定r 5 上的s 个点只，p 2 ，只，把由它们作分段插值得到 的曲线作为r 5 的近似，这样优化问题就成为 m i n ，( 只，只，只) ：。i l i 。艺0 ，( 鼻，尸2 ，只) 一u ，) 2 + 尺 ( 3 2 ) j = l 由于测量点的数目一般比较少，而侵蚀线的形状又比较复杂，一次必需恰 当地给r 5 分段，并适当地选取只，p 2 ，只的允许范围。我们设定各个只的允许 集合为条适当的线段。对只，只，只进行插值后我们可得到r 5 上的m 个特征 点用来进行边界元离散和构造正则化项。我们采用如下的正则化项 r = 口( 口州- 0 其中口，是侵蚀线上第j 到j + 1 个控制点间的有向线段的倾角，而口是正则 复旦大学硕士学位论文 化系数。 至此，我们建立了高炉炉衬热侵蚀监测的数学模型，这是优化问题。在此 模型中，对每一个给定的高炉，若已知特定的点集处的温度测量值，就可以反 演推测炉衬的侵蚀位置，这就是本文的反问题。 3 2 摄动方法的引入 由于侵蚀效应，高炉内边界是不断变化的，现在通常使用的计算方法每进 行一次正问题的计算，就需要对区域边界重新进行划分，配置节点，从本文正 问题引入的边界元离散过程来看，这是相当繁琐的。另外也不能很好的应用前 面已经得到的计算结果，增加计算量的同时也对存储空间提出了要求。事实上， 大量的机时开销就是花费在重复计算上面。本文尝试使用摄动方法来解决这个 困难。 摄动方法的基本原理是把函数值在初始边界上作渐近展开并把移动边界 表示成初始边界加上微小摄动。我们以二维移动边界问题为例说明这个方法。 假设移动边界y = j ( x ；占) 可表示成如下渐进形式： y = s ( x ；e ) 一s o ( x ) + ( 工) + - s 是一个微小参数，本文中取0 0 2 。 相应的，我们可以把移动边界上的函数u ( x ，s ( 置s ) ) 作如下展开 “( x ，s ( z ；f ) ) “( 工，s 。( 工) ) + s _ ( x ) 宴兰( x ，( x ) ) + 因此当我们把区域上的解写成如下形式时 u ( x ，y ；占) u o ( 石，y ) + 占“i ( 石，y ) + 我们可以得到： “( z ，s ( z ；) ；s ) “。( x ，s o ( x ) ) + f ( u i ( x , s 0 ( x ) ) + s 1 ( x ) ! ：2 ( x , s 0 ( z ) ) ) + d ( 占2 ) d v 这样我们可以通过分别求解固定边界问题，“。来求解移动边界问题“。 回到本文研究的宝钢2 号高炉侵蚀线监测问题。 复旦大学硕士学位论文 由于宝钢2 号高炉的移动边界是一段拐角形边界。无法简单的表示成 ( r ，z ( r ) ) 的形式。 我们用参数方法表示移动边界亍，。用f ，记初始边界。 r s 慝曷 寻j f = r ( s ) + g l ( s ) ie ： 3 【- i = z ( s ) + f 9 2 ( s ) 实际上在冶炼过程中移动边界是一条等温线，温度等于铁水的熔化温 瓦= 1 1 5 0 。c 。 u ( r ( s ) + g i ( s ) ，z ( s ) + s9 2 ( s ) ) = r o 运用上文所述的摄动方法： h l l = “( 。) i l + 占“( 1 ) i l + _ 纠一( 等卜球卅警卜硝州” m 7 o 代入原方程后把原问题求解转化成如下两个边值问题求解： 爿北等h t 警 n n 挲：o 。r 。 t 等= 吃( ) 。n r ： ( 3 ，) 警地( u w - - u r 0 1 、1 彻r 3 娑：o 。凡r 。 “( o ) = l o n f 。 复旦大学硕士学位论文 未( 心等 + 扑等 ：o 觑q 娑：o 。r i 一挲：h a ，。r ： 一k 娑： m 。r 3 挲：o 。r 。 “( 1 ) + f 旦掣o rs 。( s ) + 掣o z - s ：( s ) i = 瓦。n r ， 解上述两个边值问题即可求得h “，u “，利用关系式u u + 占“1 即能得到 u 的近似解。 3 3 计算区域的表示 要反演侵蚀线的形状，需要对侵蚀线进行离散。我们用首尾相接的折线段 来作为侵蚀线的近似。假设侵蚀线被离散成了m 一1 段的折线，那么这些折线段 共有膨个端点，这膨个端点的位置决定了折线段的形状，从而也决定了计算 区域的形状。我们把这些点称为“控制点”。 更进步，假设这些端点的位置总是在某些预先设定的线段上移动，把这 些线段称为控制点的“轨道”。这些轨道设定以后，侵蚀线的形状就可以用这m 个参数来表示了。 为了很好的反演出侵蚀线的形状，一般需要在移动边界上配置尽可能多的 控制点。但由于在高炉炉衬上埋设的热电耦有限，提供的附加条件是有限制的， 我们无法在移动边界上配置很多的移动控制点。一般来说为了保证优化问题能 够顺利进行，需要保证移动控制点的数目少于边界上测温点的数目。另外过多 的移动控制点也对优化提出了很高的要求，计算比较困难。 假设这m 个控制点的初始位置为( f ，乙)i = l ，2 ，m 2 1 复旦大学硕士学位论文 这些轨道的参数方程为( 其中z 为参数) f ，= i + f c o s 舅 f _ 1 ，2 ，m z = z + z s i n 谚 一 各个轨道参数z 的取值f i ，t 2 ，t 。就是正问题的自变量。设t = ( f ，t ：，0 ) ， 旦t 给定，就能运用前面正问题中引入的边界元方法计算出区域内各点的温 度。 假设测温点的实测温度为：( v ( f ) 一r ) ( v ( f ) 一r ) 7 ，t = ( 互，正，t ) 这些测温点的计算值为：v ( t ) = ( v ( f ) ，v 2 ( f ) ，v l ( t ) ) 而反问题可以表示为：求解t ，使得( v ( f ) 一丁) ( v ( r ) 一丁) 7 最小 3 4 数值结果 为了考察方法的有效性，我们做了大量的数值实验。计算参数为：高炉底 部外壁与空气的热交换系数吃= 3 0 w m k ，高炉侧壁与冷却水的热交换系数 h = 7 0 w m - k ，空气温度u a = 3 5 。c ，冷却水温度u w = 3 3 。c ，铁的熔化温度 矗= 11 5 0 。c ，微小参量e 取0 0 2 。在下文所述的数值算例中我们都采用这组 计算参数。 必 在给定侵蚀线的情况下，运用边界元计算高炉炉衬在测温点的温度。以这 个计算得到的温度作为已知条件反演侵蚀线。 我们把反演得到的控制点坐标列表如下： 反演得到控制点坐标控制点坐标 控制点编号径向坐标轴坐标径向坐标轴坐标 16 8 5 7 45 0 16 8 2 25 0 1 26 5 4 0 13 7 5 56 5 7 23 7 5 5 36 3 9 52 5 16 3 2 22 5 1 46 1 92 3 7 4 3 6 1 9 2 3 7 54 6 52 3 5 8 4 6 5 2 3 7 3 复旦大学硕士学位论文 63 12 ，4 1 23 12 4 0 9 71 5 52 4 0 3 51 5 52 4 1 1 802 3 9o2 3 7 2 反演得到的侵蚀线如下图虚线所示： 图3 1 无噪声情况下反演结果 另外我们把反演得到的区域作为计算区域重新计算测温点温度如下表 测温点编测温点温反演得到 县 度温度 11 7 0 1 31 6 9 1 3 21 6 4 0 3 1 6 3 6 1 31 4 2 5 51 4 1 6 49 2 1 7 69 2 0 3 9 54 5 3 6 14 5 0 8 7 67 7 2 2 97 7 1 8 79 9 5 7 59 8 7 2 9 81 1 7 1 31 1 6 9 91 3 2 9 61 3 2 3 4 根据反演结果重新计算得到的温度值与初始温度值作误差分析，结果如下 图所示： 复旦大学硕士学位论文 图3 2 相对误差 地 反问题一般是不适定的，为了检验算法的稳定性，我们在给定侵蚀边界的 条件下，给测温点的温度值加最大为千分之十的随机噪声，考察算法的数值结 果。 反演得到的侵蚀线如下图虚线所示： 图3 3 有噪声情况下反演结果 复旦大学硕士学位论文 图3 4 有噪声情况下相对误差 倒里 另外我们对反问题的计算作了正则化处理。计算结果如下图 图3 5 带正则化处理的反演结果 复旦大学硕士学位论文 图3 6 相对误差 我们利用工业实测数据对高炉侵蚀线进行了反演，得到如下图所示的侵蚀 线形状。跟经验推测比较吻合。 图3 7 工业实测算例 结论与展望：