（应用数学专业论文）煤自燃过程中对流扩散反应方程组的数值模拟.pdf

论文题目： 专业： 硕士生： 指导教师： 煤自燃过程中对流扩散反应方程组的数值模拟 应用数学 郑宁 褚维盘 ( 签名 ( 签名 赵高长 ( 签名) & 兰兰坠 摘要 本文结合陕西省自然科学基金项目“采空区煤炭自燃耦合场的有限元求解及动态模 拟”，主要结合我校x k 一i 型煤自然发火实验台的实际情况，在能量与质量守衡的原则 下，通过传质传热过程分析，建立了煤在低温氧化放热过程中温度与氧浓度变化的控制 方程，并引入了相应的初始条件与边界条件。可以看出温度场与氧浓度场控制方程格式 相同，二者通过氧化放热反应项发生耦合，组成一个二维的对流扩散反应偏微分方程 组。针对建立的对流一扩散一反应耦合偏微分方程组，分别用交替隐式预估校验方法和特 征有限差分格式进行离散和编程求解，对义马和兖卅煤的自然发火过程进行了数值模 拟。模拟计算得到的义马煤和兖州煤在与实验条件相同的条件下取得结果相符。从而印 证了本论文中所建模型与模拟方法的正确性，可以在预测煤的自然发火过程中发挥作 用。 本文的创新点有： l 、首次在耦合场中动态模拟煤自燃过程中渗流场、温度场及氧浓度场的变化规律。 这比其它单场模拟更科学，更接近现场情况，取得更合理的结果。 2 、本文构造了适合二维半线性耦合对流扩散反应方程组的预估校正格式。 3 、本文首次把秦新强提出的新型特征有限差分法用于解决二维耦合半线性对流扩 散方程组的求解。 关键词：对流扩散反应方程；特征差分；预估校正；煤炭自燃；预测预报 研究类型：应用研究 s u b j e c t ：n u m e r i ca p p r o a c h m e n to fd i f f u s i o n - c o n v e c t i o n r e a c t i o n e q u a t i o n sa n ds i m u l a t i o no fc o a ls p o n t a n e o u sc o m b u s t i o n s p e c i a l t y ：a p p l i e dm a t h e m a t i c s n a m e ：z h e n gn i n g i n s t r u c t o r ：c h uw e i p a n z h a og a o c h a n g a b s t r a c t t h i st h e s i si sj u s to n ep a r to ft h ep r o j e c to f “n u m e r i c a la p p r o a c h i n gm e t h o da n d s i m u l a t i o no ft h ec o u p l ef i e l d si nc o a ls p o n t a n e o u sc o m b u s t i o n ”w h i c hi sf u n d e db yn s fo f s h a a n x ip r o v i n c e b yt h er o l eo fm a s so rh e a tc o n v e r s a t i o n ，t h ec o n t r o le q u a t i o no ft h e r m a l f i e l da n d o x y g e nc o n c e n t r a t i o nm o d e lo fc o a ls p o n t a n e o u sc o m b u s t i o ni nt h ex k - i i id e t e c t i o n o v e nw e r ec o n s t r u c t e d ，a l s ow e r et h eb o u n d a r yc o n d i t i o n si n t r o d u c e di n t o f r o mt h e m a t h e m a t i c a lm o d e l s ，i tc a nb es e e nt h a tt h et w of i e l d se q u a t i o n sh a v es e e mf o r ma n dc o u p l e d e a c ho t h e rb yt h er e a c t i o ni t e m s ot h es i m u l a t i o nb e c o m et h ea p p r o a c h m e n to fuc o u p l e d d i f f u s i o n - c o n v e c t i o n r e a c t i o ne q u a t i o n s i nt h i sp a p e r , t w om e t h o d sw e r eu s e d ，a l t e m a t i v e i m p l i c i tp r e d i c t i o n v e r i f y m e n tf i n i t ed i f f e r e n c ea n dc h a r a c t e r i s t i cf i n i t ed i f f e r e n c em e t h o d s ，t o s i m u l a t i o nt h es e l f - c o m b u s t i o np r o c e s so fy i m ac o a la n dy a n z h o uc o a l t h er e s u l th a v ea g o o dc o h e r e n c yw i t ht h er e a le x p e r i m e n t w h i c hp r o v e dt h a tt h em o d e lc a l lb eu s e df o rt h e p r e d i c t i o no f c o a ls p o n t a n e o u sc o m b u s t i o np r o c e s s t h ei n n o v a t i o n so f t h i sp a p e ra r el i s t e da sf o l l o w i n g ： 1 ) t w of i n i t ed i f f e r e n c em e t h o d sw e r ee s t a b l i s h e df o rt h ec o u p l e df i e l d so fd i f f u s i o n c o n v e c t i o n - r e a c t i o ne q u a t i o n st os i m u l a t et h ec o a ls p o n t a n e o u sc o m b u s t i o n t h er e s u l th a s h i 曲e ra c c u r a c yc o m p a r e d 、v i mt h a tr e s o l v e db yo t h e r s s i n g l ef i e l dm o d e l s 2 ) a l t e r n a t i v ei m p l i c i tp r e d i c t i o n - v e r i f y m e n tf i n i t ed i f f e r e n c em e t h o d sw e r ec o n s t r u c t e d i nt h i sp a p e rt os o l v et h ec o u p l e dd i f f u s i o n - c o n v e c t i o n r e a c t i o ne q u a t i o n s 3 ) e x t e n dq i n x i n q i a n g sc h a r a c t e r i s t i cd i f f e r e n c em e t h o df o ro n ed i m e n s i o ne q u a t i o n i n t ot h er e s o l v e m e n to f t w od i m e n s i o nc o u p l e dd i f f u s i o n - c o n v e c t i o n - r e a c t i o ne q u a t i o n s k e y w o r d s ：c o n v e c t i o n - d i f f u s i o n r e a c t i o ne q u a t i o n p r e d i c t i o n v e r i f y m e n td i f f e r e n c em e t h o d c h a r a c t e r i s t i cd i f f e r e n c em e t h o d c o a ls e l f - i g n i t i o n f o r e c a s t i n gs i m u l a t i o n t h e s i s ：a p p l i c a t i o nr e s e a r c h 要料技大学 学位论文独创性说明 9 2 3 1 1 2 本人郑重声明：所呈交的学位论文是我个人在导师指导下进行的研究工作及 其取得研究成果。尽我所知，除了文中加以标注和致谢的地方外，论文中不包含 其他人或集体已经公开发表或撰写过的研究成果，也不包含为获得西安科技大学 或其他教育机构的学位或证书所使用过的材料。与我一同工作的同志对本研究所 做的任何贡献均已在论文中做了明确的说明并表示了谢意。 学位论文作者签名：鄯予 日期：a 一6 加 学位论文知识产权声明书 本人完全了解学校有关保护知识产权的规定，即：研究生在校攻读学位期间 论文工作的知识产权单位属于西安科技大学。学校有权保留并向国家有关部门或 机构送交论文的复印件和电子版。本人允许论文被查阅和借阅。学校可以将本学 位论文的全部或部分内容编入有关数据库进行检索，可以采用影印、缩印或扫描 等复制手段保存和汇编本学位论文。同时本人保证，毕业后结合学位论文研究课 题再撰写的文章一律注明作者单位为西安科技大学。 保密论文待解密后适用本声明。 学位论文作者签名：指导 1 绪论 1 1 问题的提出 1 绪论 煤炭开采后长时间暴露在空气中与其中的氧发生物理化学作用，释放出热量在煤体 中积聚起来，使煤体温度上升，最终导致煤炭自燃。据统计，我国煤层中一半以上届易 自燃煤层，国有煤矿和地方重点矿井约有4 6 4 9 发生过火灾事故，其中属自燃火灾 的占9 0 以上【l 引。井下煤炭自燃由于高温火源的存在经常引发瓦靳爆炸、煤尘爆炸等 一系列重大事故，严重威胁着矿工的生命安全，给国家带来了惊人的资源损失和设备破 坏。煤炭自燃灾害造成的直接和间接损失每年在5 0 亿元以上，已经成为制约中国煤炭 工业发展的重要因素1 3 j 。为此，中国政府在2 l 世纪议程”中己将煤炭自燃列为重大自然 灾害类型之一，要求加强煤炭自燃机理研究及自燃发生过程中的跟踪临测，建立合理预 测、预报机制，杜绝和预防矿井火灾隐患的发生1 4 j 。 对煤炭自燃过程的研究与预测预报早在煤炭的大规模工业开采以后就已引起了世 界范围内的关注。起初主要通过对少量煤样的氧化自热反应过程进行研究，定性地比较 不同煤种之间自燃倾向性的差异口7 】，但由于现场影响煤自燃过程的因素多种多样，单 一的自燃倾向性只能比较煤发生自燃的难易程度，不能计算煤的自燃发火期为火灾的预 防提供定量的信息，因而实用性受到很大限制。之后，人们根据反应热对煤自燃发火期 进行计算来指导矿井采区的布置和工作面的推进速率【3 】，但因为自燃发火期是通过对实 验室绝热状态下测得的反应热计算得到，忽视了采区实际生产状况下由漏风引起的对流 散热，因而与实际情况误差较大，对火灾预防的指导作用也不明显。 进入2 0 世纪9 0 年代后，世界各主要产煤国先后建立了模拟煤层自燃过程的大型自 然发火实验台，通过相似模拟实验以预测实际条件下松散煤体和自燃危险程度和自然发 火期，实验得到的煤炭自燃发火期与现场实际情况吻合程度较高，因而在煤炭开采的火 灾预防方面得到了应用。但由于实验周期较长( 2 - 5 个月) ，一方面在时间上已不能满足 各大矿区送来的大批煤样自燃倾向性测试，另方面，由于现场生产中风量、工作面推 进速度随时改变，而在煤自然发火台的模拟实验中这些参数的调节不易灵活改变，因此 大多矿区只是将实验得到的发火期作为采区设计时的一个重要参考，还未能针对具体某 一采区的实际生产过程作到及时、准确的预测预报。 由于煤炭自燃过程是煤氧之间发生的氧化放热作用与外部环境对流散热的综合结 果，因而在现场多根据采空区漏风情况对自燃发火过程做出经验性的判断。目前我国还 普遍采用氧浓度和风速的大小将采空区划分为“二区”( 易自燃区、不易自燃区) 和“三 带”( 散热带、氧化自热带及窒息带) ，用以指导井下煤炭自燃火灾可能发生区域的预防 西安科技大学硕士学位论文 性处理【9 。1 2 】。由于不同煤种的氧化反应性不同，对应的i 临界风速也有较大差异，而且随 着开采的进行，“二区、三带”是动态移动的，所以这种按统一的漏风速率参数指标划分 的“二区、三带”只是一个大致范围，无法准确、动念对采空区煤炭自燃的实际过程进行 评判。 因而现行常用的煤炭自燃预测预报主要是利用井下开采中的一些安监信息( 如特征 气体浓度，漏风量，温度，吸氧量) 来形成煤炭自燃倾向性判断标准。这种预测大部分 靠主观经验，并且存在着严重的滞后型，其中的一些特征指标都是在出现暗火后才测得 的，还未能做到真正的提前预测与防止，有很多情况是当安监人员根据工作面释放气体 浓度超标做出火灾的预报时，采空区火情己大面积蔓延到无法控制的地步，只好封闭正 在开采的采区，部分矿区还发生过整套综采设备被烧毁的情况，损失惨重。 随着近年来大型矿难事故频繁发生，国家对矿井生产安全问题的重视己上升到前所 未有的高度。如何根据矿井下的实际条件建立与实际情况相符合的数学模型，实现煤自 燃过程的提前准确预报，防忠于未然，是煤矿安全生产一个迫切需要解决的问题。 1 2 论文背景及国内外研究进展 1 2 1 煤炭自燃过程的数学背景 具有一定氧化活性的煤、充分的空气供给以及适宜的温度环境是决定煤炭自燃所必 不可少的三大因素。其中煤与空气是物质基础，煤在有氧的情况下发生的氧化自热反应 是自燃发生的根本原因，散热环境是决定煤发生自燃还是进行风化的外部条件，当散热 速率大于放热速率时，氧化反应产生的热量将通过传导、对流的方式散失在周围环境中， 这时煤体温度不会上升，因此只能发生缓慢的风化。当散热条件较差时，放出的热量导 致煤体温度上升，加速氧化反应的进行，对氧的消耗加剧，这时能否引起自燃则主要取 决于供风量的大小，若风速较小，将会消耗掉空气中的氧而使反应停止，而风速过大又 会使反应放出的大部分热量被风流带走，进入风化状态，只有当风速适量时，才会发生 持续的氧化自热反应，最终引发自燃。 因此可以看出，煤炭自燃过程实际上是温度场、氧浓度场与空气渗流场三场耦合叠 加的综合作用结果，根据各场的物理意义，可建立其控制方程。 1 ) 空气在煤体多孔介质中的渗流场 由于空气粘度较小，流动比较稳定，由于在研究范围内静压差很小，故空气在煤等 多孔介质中的流动可视为平面流。由于在采空区及自然发火台中煤的形状、大小及分布 的规律都是随机的，因而可以不考虑这些因素对空气在其中流动的影响，即视为各向同 性。 根据多孔介质中空气渗流动力学和d a c c y 定律 1 3 - 1 4 】，空气在煤体中的渗流场可用如 2 下渗流连续性方程( 1 - 1 ) 表示，根据( 1 - 2 ) 可计算不同位置的渗流风速，它决定着体 系的供氧量。 d i v ( k g r a d ( h ) 、= 0 ( 1 - 1 ) v kg r a d ( h )( 1 - 2 ) 式中：丑- _ 各点处的风压；k 为渗透系数，当空气穿过自然发火台气体缓冲区或工 作面的自由空间时，后= 1 ；v 为漏风渗流速度矢量。 由式l 。1 可以看出，渗流风速仅与坐标有关，不随时间变化，是一个典型的椭圆 型偏微分方程( 源项为0 ，拉普拉斯方程) 。 2 ) 温度场和氧浓度分布的控制方程 介质中热量的积聚和散失主要通过热源补偿和对流、传导、辐射三种热传递方式进 行，由于煤自燃研究的温度范围一般在4 0 0 c 以下，且为松散颗粒堆积体，故辐射的影 响很小，可忽略不计。根据能量守衡定律，单位体积内导致煤体温度变化的积聚热量等 于单元内反应放出的热量与通过对流、传导向周围传递散失热量的代数和。 由于传质与热传递具有类似的性质，根据物质守衡定律可知，单位体积内某一组分 物质的积累等于由对流、扩散进入单元内物质的量与反应消耗的量的代数和。 由此可写下温度场和氧平衡的控制方程如下： n i o t ：6 l 掣+ 咖( 以g r a d ( t ) + q 咿刃 ( 1 3 ) 吒i o c ：趁旦嘤+ d i v ( 见寥位f ) + c 2 似c 刃 ( 1 - 4 ) 式1 3 和式1 4 中，方程左边的一项分别对应的是热量和氧在单元中的变化速度， 积分后对应的是一段时间内积聚的热量或氧浓度增量，等式右边的一、二、三项分别对 应的是热量和氧的对流项、扩散项和反应项，其中反应项是与温度t 及氧浓度c 有关的 函数。其中a ，、口2 、b 。、b 2 、q 和c 2 是在平衡计算中各项对应的系数。 式1 3 和1 4 格式非常相似，是含有时变的非稳态扩散对流方程，二者通过反应项 相互耦合，因此有的学者将这一类方程归结为对流一扩散一反应方程。 可见煤炭自燃过程是一个包含有化学反应的传质、传热的动态耦合场求解问题，由 于渗流场是不含时变的椭圆方程，仅与某点的坐标有关，目前求解方法己比较成熟，对 流扩散反应方程是一类非稳态的非线性偏微分方程，是目前数学偏微分方程求解研究 中的一个热点，而由浓度场与温度场构成的对流功“散一反应方程相互耦合，使这一问题 的求解更加困难。本论文的重点是针对对流扩散反应祸合偏微分方程求解方法进行研 究，主要是结合工程背景提出可实现煤炭自燃过程的动态模拟与预测预报的数学方法。 西安科技大学硕士学位论文 1 2 2 煤炭自燃过程数值模拟的进展 由于空气渗流场是稳态的椭圆型偏微分方程，结合现场测得的风压边界条件进行数 值求解的方法也较为成熟，李宗翔等i l5 j 利用该方程对瓦斯及空气在煤中的渗流规律进行 了模拟研究。蒋署光【l 、杜礼明1 1 7 j 等也根据计算得到的漏风速率分布模拟了采空区三 带的划分。正是由于非稳态偏微分方程的求解困难，在研究煤炭自燃过程时，大量文献 虽然建立了数学模型，却没能找到求解办法实现自燃的预测预报，部分文献在求解这类 模型时，往往通过大量简化，将各个场分割开来单独求解，或采用经验公式的办法避免 对祸合偏微分方程组的求解，造成预测的结果与实际情况偏差较大，限制了预测模型的 应用。单独对温度场【18 j 或氧浓度分布i l 训进行求解，由于经验公式比较相糙，无法区分各 煤种自燃倾向性方面的差异，没有考虑各点温度及氧浓度对反应放热的影响，因而结果 与实际情况有较大的差别。 2 0 0 1 年f a k g u n l 2 0 1 和v f i e r r o l 2 l j 首先利用有限差分方法对温度、氧浓度耦合场进行 了求解，对煤堆的氧化自燃过程进行了动态模拟，在低温时与煤堆中温度场的分布变化 规律较为符合，但在快发生自燃时的较高温度下，收敛困难，误差较大，并常出现飞温 现象，在温度超过2 0 0 后的结果不太可靠。 1 2 3 对流扩散偏微分方程求解的数学方法发展现状 对流和扩散现象大量的出现在自然界各个工程领域中，从多孔渗水介质的散布追踪 到可溶物在河口和近海的扩散，从污染物在浅湖的蔓延到河床对化学药品的吸收和各种 生产电力的方法等等都是以对流扩散为基本过程。所有这些对流和扩散过程的数学模型 可归纳为对流扩散方程，一般形式如下： ，j c ( x l ，x 2 ，x n ，) + 6 ( x l ，z 2 ，一，x n “) v u v ( a ( x l ，x 2 ，一，x n , u ) v u ) = f ( x l ，x 2 ，一，x n ，“) 式中 u i u 为通用变量，可以代表流体的速度矢量在空问坐标系上的分量和温度等求解变量； b ( x l ，x 2 ，x n “) = 【6 l ( x i ，x 2 ，x n , “) ，b2 ( x l ，工2 ，一，x 。，甜) ，b 。( x 1 ，x 2 ，一，x 。，“) ) ； 日( x 。，x 2 ，x 。，u ) 为广义扩散系数；厂 ，x 2 ，x 。，“) 为广义源项。这里引入的：“广义” 二字，表示处在a ( x l ，x 2 c x 。，u ) 与f ( x i ，z 2 ，x 。，u ) 上的项不必是原来物理意义上的量， 而是数值计算模型中的一种定义。对于上述对流扩散问题，数学界已经发展了不少获得 其精确解的数学方法。这些精确解是在整个求解区域内连续变化的函数。但是直到目前， 这些分析解还只能对少量的简单的情形得出，对于大量具有工程实际意义的对流扩散问 题，数值计算的方法越来越广泛地得到应用。在过去的几十年内已经发展出了多种数值 方法，其间的主要区别在于区域的离散方式、方程的离散方式及代数方程求解的方法这 三个环节上。在对流扩散问题的计算中应用较广泛的是有限差分法、有限元法、有限容 4 1 绪论 积法等【2 2 郴】。 有限差分法导源于牛顿、欧拉以差商代替导数的思想，于1 9 2 8 年，由 c o u r a n t ，f r i e d r i c h s ，l e w y 首次对偏微分方程的差分方法作了完整的论述，成为现代理论 分析的先导。有限差分方法在2 0 世纪5 0 年代得以蓬勃发展。差分法的优点是容易适应 各种类型的微分方程，形成差分方程的计算量也比较小，且易于编程；缺点是不容易适 应几何形状复杂的求解区域和复杂的边界条件。偏微分方程的差分方法的大致过程为： 1 ) 剖分求解区域，得到离散问题的网格。为便于作导函数的离散，通常采用比较 规则的剖分。 2 ) 通过各种途径，如用差商代替方程中的导数，或者对原方程积分后应用数值积 分或者数值微分，将求解区域上的微分方程离散为网格点上的代数方程，称为差分方程。 3 ) 连同原问题的初始条件、边界条件的离散，得到一个封闭的代数方程组( 差分 方程组) 。该代数方程组的系数矩阵是稀疏的，有利于节省存储量和计算量。如果原来 的微分方程有某些特殊性质，如对称正定性，系数矩阵还可能保持。 4 1 根据系数矩阵性质，选择适当的方法求解差分方程组，得到在网格点上近似解 的数值。 对流扩散方程常用的差分格式有中心显式格式，修正中心显式格式，迎风差分格式， s a m a r s k i i 格式，指数型差分格式，隐式格式等。 有限元法实质上是古典的变分近似法( 如r i t z 法和g a l e r k i n 法) 的推广，是将变分 原理和剖分插值相结合的一类数值方法。利用一般有限元法求解数学物理问题的过程大 致如下： 1 1 利用某种变分原理，将原问题化为变分问题，然后用变分近似法求问题的广义 解。由于同一个数学物理问题可以转化为不同的变分形式，故采取何种变分原理要取决 于原问题的性质和使用者的需要。例如，求解p o i s s o n 方程边值问题时，将未知函数u 看作结构力学中的位移。如果采用最小总位能原理，则原问题化为求能量泛函的极小， 再用r i t z 法求解得到未知量u ；如果采用虚功原理，则原问题化为求解虚功方程，用 g a l e r k i n 法求解也得到未知量i 1 ；如果采用h e l l i n g e r - r e i s s n e r 变分原理。则原问题化为 求h e l l i n g e r - r e i s s n e r 泛函的驻点，用混合有限元法求解，可以同时求出u ( 即位移) 和 它的导数u x ，m ( 即应力) 。这些变分问题中未知函数导数的阶数会降低：且问题以积分形 式出现，便于求它的广义解，从而扩大了解的范围；一些边界条件变为自然边界条件， 对数值求解有利。 2 ) 对求解区域q 作剖分r h ，只要遵守剖分的一般规则，所得的有限个闭集e ，( 称 为单元) 不一定都属于同一类型，更不必是均匀剖分，以便于逼近复杂边界。但是剖分 的方式要和下面的分片插值相协调。 3 ) 采用分片插值的方法构造满足本质边界条件的试探函数类和检验函数空问，即 西安科技大学硕士学位论文 分别在每个单元上插值，并且插值函数能在不同单元间保持某种光滑性连接。这种分片 插值的办法不仅能克服古典变分近似法在构造基函数时遇到的困难，而且构造的基函数 都是小支集。从而使将要求解的代数方程组的系数矩阵为稀疏矩阵，大大节省了计算量。 4 ) 作单元分析和总体合成。即将变分泛函或变分方程( 它们都是在求解区域q 上 的积分) 的计算化为在每个单元e 上积分之和。由于上述分片插值的结果，将变分问题 的近似解表示为试探函数空间基函数如。 的线性组合u 。( p 。后，在e 。上的积分值都可 以用阶数不高的矩阵k i 和向量u ，f i 表示；而总体合成就是将这些单元e 上的积分累加， 采用简单的技巧，其过程可以化为矩阵k i 和向量u f 的简单叠加。 5 1 处理本质边界条件，使在边界上的网格点近似解被强制的满足本质边界条件。 最终得到一个线性代数方程组。 6 1 求解所得到的代数方程组k u = f 。其系数矩阵k 除了稀疏性外，通常还能保持 原问题的一些特性( 如对称性，正定性等) ，有利于求解。 7 1 从k u = f 的解，通过试探函数空间的基函数，得到变分问题的近似解，称为原 问题的有限元解。 混合有限元法是将问题中所求的未知函数，除了原来的u 以外，还将u 的导数( 一 阶或高阶) 作为补充的独立变量一起求解，在解n a v i e r - s t o k e s 方程和重调和方程以及结 构力学等问题中，经常使用这一方法。 有限体积法( f i n i t ev o l u m em e t h o d ，简称f v m 方法) ，是集有限差分法和有限元法 之优点而发展起来的一种新的求解偏微分方程数值方法，尤其适用于满足物理守恒律的 方程。在对求解区域作有限剖分后，将原方程在某个子域上积分：应用g a u s s 定理，将 散度的积分转化为子域边界上的积分；然后选取适当的有限维试探函数空间，在该子空 间上离散含子域边界积分的方程。由于此方法将散度的积分化为子域边界积分后再离 散，数值解满足离散守恒律；而且可以采用非结构网格，所以在计算物理，特别是计算 流体力学上有广泛的应用。 用标准的差分方法或有限元方法来处理对流扩散问题常常失效，会出现较大的数值 扩散或数值振荡等现象，其根本原因就在于对流项的存在。针对对流占优扩散问题这一 数值求解困难，许多作者研究了不同的方法，1 9 8 2 年，j d o u g l a s ，j r 和t fr u s s e l l 等提 出了特征线修正技术求解对流占优的对流扩散问题。这一方法的主要思想是引入特征方 向的导数，把对流扩散方程化为形式上不含有对流项的扩敖方程求解。当以对流为主时， 解沿特征方向的变化比沿时间方向的变化慢的多，可以有效的克服数值震荡，保证数值 解的稳定。与传统的数值计算方法( 差分法、有限元法等) 相结合，提出了特征有限差 分法、特征有限元法( 在不损失精度的情况下，取较大的时间步长，提高计算效率) 、 特征混合有限元法，并给出了理论分析。几十年来，该技术得到了很好的研究和应用。 6 1 绪论 而特征差分格式的成功与否，关键在于插值方式的选择。线性插值可以避免数值振 荡，但却有较大的数值扩散。普通的二次l a n g r a n g e 插值可以减少数值扩散，但却有严 重的数值振荡。高阶单调插值虽然可以避免数值振荡，减少数值扩散，但一般计算量偏 大。为了提高计算的精度，陆金甫、由同顺、秦新强、田振夫、程爱杰等人对特征差分 法进行了修正和改进，提出了精度更高的特征差分算法。 求解对流占优问题的另一有效数值解法就是基于修正等价偏微分方程的差分方法。 修正等价偏微分方程是1 9 7 4 年w a r i m i n g 和h y e t t 提出的，利用这一方法可以简单直接 地判断差分格式的误差，同时可以提出精度更高的有限差分格式。m e h d id e h g h a n 以修 正等价偏微分方程为基础，构造了若干求解常系数对流扩散方程的无振荡无数值扩散高 精度两层加权差分格式。1 ”，j 0 j 另外，待定系数法也是构造高精度差分格式的一个很有效的方法。利用这一方法， 文献p ”构造了一系列求解对流扩散方程的高精度差分格式。 近三十年来，在上述基本方法的基础上，针对不同类型的发展方程问题( 尤其是各 种非线性和耦合问题) ，探寻可靠的高效、高精度的数值计算方法的努力始终没有间断 过，不断地涌现出新的数值方法，如有限体积法与广义插分法、特征和迎风有限元法、 间断有限元法等。值得重视的还有，由我国学者冯康院士倡导的从几何角度出发寻求发 展方程的保结构算法的研究，这是对于构造数值方法的依据和观念上的一个重大革新。 另外，近2 0 年来谱与拟谱方法的研究也取得了得大的进展。发展方程的数值求解问题 在科学与工程计算中处于十分重要的地位，已被广泛的应用于气象预报、地震预测、油 田的勘测与开发技术、机翼与汽轮机叶片等工业产品设计、生态与环境的动态模扭等领 域。适应现代科学技术的突飞猛进，关于线性与非线性发展方程数值计算方法的研究必 将在理论与应用方面得到更加迅速的发展。但迄今为止，还不能提出一种在稳定性、准 确性与经济性各方面对各类对流扩散问题都较理想的格式。目前，人们提出的数值解法 还只能解决带有简单初值条件和边界条件的对流扩散问题，而大量具有工程实际意义的 对流扩散问题，往往都带有复杂初值条件和边界条件，如何求解带有复杂初值条件和边 界条件的对流扩散问题还需要进一步研究。 对于一个实际问题中产生的微分方程，无论是利用有限差分法还是有限元法去寻找 它的数值解，最后总是被归结为一个大型的( 一般是稀疏的) 线性代数方程组的求解问 题。为了确保计算满足一定的精度要求，一般应将网格的剖分尺寸取得比较小，但这样 做有可能使相应的代数方程组系数矩阵的阶数变得很大，更为严重的是导致了它的条件 数也随之变得很坏。这时若采用经典的迭代方法，如r i c h a r d s o n 迭代、j a c o b i 迭代或 g a u s s s e i d e l 迭代去求解，收敛的速度将随步长的缩小而变得缓慢，计算的效率也就很 低了。因此前苏联的计算数学专家r p f e d o r e n k o 在6 0 年代中期提出了多重网格法 ( m u l t i g r i dm e t h o d ) ，它考虑一系列( 从粗到细) 的网格，对应系列的方程组，它们未 7 西安科技大学硕士学位论文 知数的数目各不相同( 从少到多) ，一般用迭代法求解，而且在各方程组之间进行适当 的转换。7 0 年代，多重网格法受到了普遍的重视，被认为是一种行之有效的数值方法， 尤其在8 0 年代之后，这个方法的理论研究取得了重要成果，并在实际应用中显示了它 强大的生命力。正是由于近年来，各种偏微分方程的求解方法研究得到了重视和发展， 尤其是有限元分析中新方法的出现，如混合有限元法、边界元法、多重网络划分、自适 应网格划分技术解决r 方程收敛困难、迭代计算量人的问题，使一些动态的非线性问题 求解成为可能。 1 ，3 论文的研究内容与思路 2 0 0 4 年课题组申请了陕西省自然科学基金项目“采空区煤炭白燃耦合场的有限元求 解及动态模拟”，提出了一种煤炭自燃快速预测预报机制。在前期煤炭氧化反应动力学 研究和宏观传质传热参数测试的基础上，建立煤炭自燃传质、传热数学模型，并根据偏 微分方程求解方法的最新进展对这些耦合方程组进行求解，动态模拟煤炭自燃过程中渗 流场、温度场及氧浓度场的变化规律。对任何种煤，在几天时间内对少量煤样通过化 学动力学研究的物性参数实验，赋予现场边界条件等参数，就能对其从低温氧化到自然 发火过程中的渗流场、温度场、氧浓度分布进行动态的可视化仿真模拟，计算其发火期 和相参数。在此基础上通过边界条件的改变可模拟煤层厚度及工作面供风量和几何参数 等对煤自然发火过程的影响，不仅为寻找有效的火灾预防提供帮助，而目还可对矿井、 采区的布置、设计提供服务。 本论文将根据煤炭自燃发火过程中的氧化反应及传质、传热规律，建立合理的数学 模型，动态模拟煤从低温氧化发展到自燃着火这一过程，利用偏微分方程的理论和方法， 确定温煤炭自燃过程中温度场、氧浓度场等的变化规律，同时还可以通过调节边界条件， 1 绪论 模拟不同参数对煤白燃发火过程的影响，为建立及时可行的自燃预测预报机制奠定基 础。 由于矿井采空区环境条件较为复杂，给煤自燃过程的研究带来了困难，因此利用煤 自然发火实验台先对白燃过程中传质、传热及氧化放热规律进行模拟和研究，再根据这 些规律对井下煤炭自燃情况进行评判和预测是一种行之有效的方法，本论文主要针对我 院建立的煤自然发火实验台进行数值模拟，在总结前人工作的基础，针对煤在自然发火 实验台中的传质传热规律，建立数学模型，并利用已测得的传质、传热参数对耦合场偏 微分方程进行求解，动态地模拟煤的自然发火过程，实现量化的预测预报机制，希望经 过不断的验证和改进，能为井下煤炭自燃的仿真模拟奠定基础。 9 西安科技大学硕士学位论文 2 煤自然发火台数学模型的建立 2 1 煤自然发火实验台的发展与趋势 煤炭自然发火实验台是2 0 世纪9 0 年代确立起来的种煤炭自然发火过程相似模拟 实验技术。一般为一个较大的柱状容器，内部充填松散煤体，并按需要预先布置温度监 控与气体成分的测点，由风机压缩的空气自下而上通过煤体来模拟煤堆或采空区的实际 供风条件，当空气与煤接触时将发生缓慢的氧化放热反应而导致煤体温度上升，最终引 发煤炭自燃，由监控装置连续记录各处的温度、气体成分的变化规律，为贮煤场或井下 采空区煤炭自燃的预测预报提供服务。由于装煤量大，通过风速及环境温度的改变基本 能够反映煤体在实际情况下的氧化反应以及传质、传热过程，实验得到的自然发火期准 确度较高，因而在国内外得到了广泛的应用。 最初由s t o t t ，j b 【3 2 】于1 9 8 0 年在美国矿业局建立了长5 m ，直径o 。6 m 的柱体实验 台；1 9 9 1 年在s t o t t ，j b 的指导下x d c h e n 9 1 3 3 等在新西兰设计建造了长2 m 、直径 o 3 m 、装煤量1 1 0 k g 的一维自燃实验装置；1 9 9 5 年a r e i fa 等【3 4 】根据x d c h e n g 的 实验装置，在澳大利亚昆士兰大学建立了长2 m 、直径0 2 m ，装煤量6 0 k g 的煤自燃 实验台。1 9 9 1 年s m i t h a c ，m i r o ny 和l a z z a r a p 副在美国矿业局建立了装煤量1 3 t 的 实验台；1 9 9 8 年c l i f f d rb e n n e t 和a g a l v i n p6 j 在澳大利亚昆士兰采矿安全测试与研 究中一i l , ( s i m t a r s ) 建立了装煤量1 6 t 实验台。 同时，国内的一些专家也开始从事这方面的研究，其中湘潭工学院李仁发等人 ( 2 0 0 1 ) ( 37 建造了可以装9 0 蜉碎煤的实验装置，此装置装煤量少，人为调节温度、湿 度变化大，不能很好的模拟现场的实际条件；淮南工业学院张国枢等人( 1 9 9 9 ) 【3 别也设 计研制了实验室内模拟煤炭自燃的实验装置及其参数测定系统，模拟和研究煤炭在常温 条件下啻燃的发生、发展过程及其影响因素，研究自燃火源形成及其分布规律等，测试 煤的自燃性；中国矿业大学张瑞新1 3 9 1 ，【4 0 】等建立了模拟煤堆自燃的大型实验台，研究露 天煤体和煤堆的自燃过程。 西安科技大学的徐精彩教授等在上世纪9 0 年代陆续建立了一系列大型的煤低温氧 化自然发火实验台，实验煤量从5 0 0 k g 到1 5 0 0 k g ，通过风量调节模拟煤在接近现场条件 下的自然发火过程1 4 “。课题组对我国各大矿区的主要煤种进行了自燃倾向性与发火期测 试，并在煤低温自燃发火实验模拟和综放面采空区大量现场跟踪观测的基础上，根据对 采空区三带划分和综放面采空区自燃主要影响因素等理论分析结果，提出了综放面及采 区煤自燃危险区域的评判预测模型。 由于自然发火模拟实验得到的煤炭自燃发火期与现场实际情况吻合程度较高，在矿 1 0 2 煤自然发火台数学模型的建立 井开采生产与设计中有重要的参考价值。为了提高精度，近年来自然发火实验台装煤量 越来越大，在2 0 0 4 年，由西安科技大学与兖矿集团合作建立了国内最大的自然发火实 验台，装煤量可达1 5 吨。 目前国内外煤自然发火实验台一般均是通过长时间跟踪煤低温氧化放热过程的模 拟实验估算出煤的放热速率与耗氧速率，实验周期较长( 2 5 个月) ，一方面在时间上己 不能满足各大矿区送来的大批煤样自燃倾向性测试，另一方面，由于现场生产中，风量、 工作面推进速度随时改变，而在煤自然发火台中的模拟实验中这些参数的调节不易灵活 改变，因此大多矿区只是将实验得到的发火期作为采区设计时的一个重要参考，还未能 针对具体条件采区的实际生产过程作到及时、准确的预测预报。 因此如何根据煤自然发火模拟实验中煤体温度及氧浓度的变化规律建立对煤种及 供风条件的改变有普遍适应性的预测模型，可大大减少同一煤种在不同条件下的自然发 火台模拟实验次数，不但节约了成本，而且可有效地缩短实验时间，通过模拟参数的调 节灵活地适应现场开采条件的变化，更加及时地为井下自燃火灾进行预测、预报。 2 2x k - i i i 煤自然发火实验台的实体简介 西安科技大学于2 0 世纪8 0 年代未建立了我国第一个大型煤自然发火实验台( 装 煤量1 吨左右) ，随后相继建立了装煤量o 5 吨和1 5 吨等多个x k 系列煤自然发火 实验台。x k 系列煤自然发火实验台，主要模拟现场实际的蓄散热情况、漏风状况及浮 煤堆积厚度，以井下一般温度( 1 5 。c 3 00 c ) 作为实验起始温度，利用煤自身的氧化放 热引起自燃升温，最终导致自燃，通过连续检测实验炉内各点煤样的温度、气体的变化 情况，计算出不同温度下的煤的放热强度、耗氧速率、各种气体的产生率等参数，以及 测定煤的白燃倾向性及自然发火期。x k 型系列自然发火实验台结构基本相同，只是在 装煤量、控温方式等方面略有不同。目前对煤自燃倾向性分析预测模拟实验主要以 x k - i i i 型煤自然发火台为主。 x k i i i 型煤自然发火实验台结构原理如图2 1 ，主要由炉体、气路及控制检测系统 三部分组成。炉体呈圆柱形，外径1 8 m ，内径1 2 m ，炉内盛煤空间高1 7 5 m ，有效装 煤体积1 6 3 m 3 ，最大装煤量1 5 0 0 k g 。炉体内衬为5 0 n u n 厚的玻璃钢，外层是5 m m 钢筒， 中间填充有1 0 0 m m 厚的聚胺酯泡沫材料保温层。炉体项、底部各留有2 0 c m 自由空间以 保证气流均匀，空气由空压机提供，经三通流量控制阀、从炉体底部通过煤样，从顶盖 出口排出。炉内按一定方式排布温度测点和采气点，温度巡检由计算机自动完成，采集 的气样由气相色谱仪分析。炉体外部结构与检测点的布置如图2 2 所示。当煤体局部温 度急剧上升出现异常时，即认为煤炭发生自燃，得到煤炭自燃发火期等实验数据。 西安科技大学硕士学位论文 图2 1x k i l i 型自燃发火炉结构示意图 f i g2 1e x p e r i m e n t a ls e t u pa n ds k e t c hm o d e lf o rc o a ls e l f - i g n i t i o n t e s t 图2 2 实验台外观与内部温度、气体分析测点的布置 f i g2 2l a y o u to f t h ee x p e r i m e n to v e na n dt h ed e t e c t i o ne l e m e n t sa r r a n g e m e n t x k i i i 型自燃发火炉是我校为各大矿区测试自然发火期用得最多的实验台，对于低 变质程度的烟煤来说，由于反应性较好，从低温氧化到发火需要1 - 2 个月，对于高变质 程度烟煤或无烟煤，由于氧化反应速率较慢，实验周期更长，需要3 - 5 个月。如果现场 要求改变风速，还需重新测试，因而还存在耗时、耗财方面的明显不足。本论文在建立 煤炭自燃过程数学模型的基础上，可以模拟不同风速及环境温度对煤炭自然发火过程的 影响，有效地为现场服务。 1 2 2 煤自然发火台数学模型的建立 2 3 煤自然发火过程数学模型的建立 2 3 1 体系的简化 由于影响煤自燃过程的因素复杂，为了便于求解，在不影响精度条件下对部分影响 因素进行简化，作如下假设： 1 1 松散煤体在所研究的范围内是均匀的，粒度分布、孔隙率、密度和扩散系数和 反应性等物化性质相同，并且在整个反应过程中，煤体的密度和孔隙率不随时间变化； 2 ) 实验中炉体内风速很小，可近似认为通过填充煤粒的风流温度与煤温相同，风 流沿纵向z 轴方向均匀流动，在煤体局部微小单元内，气流温度与煤温相等； 3 ) 由于环境温度随室风气温变化复杂，而为了避免环境温度的影响，炉体外壳中 采用了保温水层和绝热聚酯发泡材料夹套，因此在模拟过程中，可以认为炉体外壁温度 保持不变： 4 1 由于在实验考察的范围内温度不超过3 0 0 ，热辐射效应较小，炉体内的热平衡 主要由煤的氧化放热与传导、对流两种散热作用为