基于MATLAB的两种通风网络解算方法的编程及实现.pdf

I S S Nl 0 0 5 2 7 6 3C N4 3 一1 2 1 5 T D矿业研究与开发第2 8 卷第2 期M I N I N GR D V 0 1 2 8，N o 22 0 0 7 年4 月A p r 2 0 0 8基于M A T L A B 的两种通风网络解算方法的编程及实现陈梅芳，丁德馨，莫勇刚，张传飞(南华大学核资源与安全工程学院，湖南衡阳市4 2 1 0 0 1)摘要：介绍了通风网络解算的数学模型及拟牛顿算法，针对计算通风网络回路分支风量的拟牛顿法和直接非线性方程组求解法。应用M A 7 r L A B 软件编程实现，并给出了实例解算过程及源程序。这种方法能够方便的进行通风网络解算，为通风网络优化计算提供了一种有效的方法。关键词：通风网络；拟牛顿法；非线性方程法；M A T I A B中图分类号：T D 7 2 2文献标识码：A文章编号：1 0 0 5 一”6 3(2 0 0 8)0 2 0 0 6 5 0 3P r o g r a m m ea n dI m p l e m e n t a t i o no fT w oV e n t i l a t i o nN e t w o r kS o l u t i o n sB a s e dO i lM A T L A BC h e nM e 扣n#，D i n gD e x i n，M oY o n g g a n g，Z h a n gC h u a n f e i(S c h o o lo fU r a n i u mR e s o u r c ea n dS e c u r i t yE n g i n e e r i n g，N a n h u aU n i v e r s i t y，H e n g y a n g4 2 1 0 0 1，H u n a n，C h i n a)A b s t r a c t：I nt h i sp a p e r，am a t h e m a t i c sm o d e la n dq u a s i N e w-t o nm e t h o do fc a l c u l a t i n gv e n t i l a t i o nn e t w o r ki si n t r o d u c e da tf i r s t M A T L A Bs o f t w a r ei su s e dt op r o g r a mf o rq u a s i N e w t o nm e t h o da n dn o n l i n e a re q u a t i o ns e tt oc a l c u l a t et h ea i rv o u l l l eo fb r a n c h e sc o n t a i n e di nv e n t i l a t i o nc i r c u i t，t h ed e t a i l e dc a l c u l a t i o np r o c e s so fae x a m p l ea n dt h es o l l l c ec e d ea r ep r e s e n t e do u t T h em e t h o dc a nc o n v e n i e n t l yi m p l e m e n tt h es o l u t i o no fv e n t i l a t i o nn e t w o r k，s op r o v i d e sa ne f f e c t i v ew a yf o ro p t i m i z e dc a l c u l a t i o no fv e n t i l a t i o nn e t w o r k K e yW o r d s：V e n t i l a t i o nn e t w o r k，Q u s i N e w t o n m e t h o d，N o n l i n e a re q u a t i o ns e tm e t h o d，M A T L A B风网解算是在只有扇风机作用而不加任何调节设施的矿井通风网络中求解各条分支分配的风量。这种网络也叫自然分风网络，网络中所有分支的风阻已知，扇风机的工作风压已知。风网解算方法有图解法、模拟法、数值法等，其中数值解法容易在计算机上实现。本文主要讨论基于M A T L A B 的拟牛顿法和解非线性方程组法的编程求解及实现。l数学模型回路风量法是根据风量平衡定律、风压平衡定律及分支特性方程，得到分析风网所需的方程组，对于一个有m 个节点，几条分支的通风网络，有节点风量平衡方程m 一1 个：弘研2o，i=1，2，m 一1(1)回路风压平衡方程，l m+1 个：蓍c 口一=O,i=1，2，n m+1(2)风压降方程乃个：_=0lc bl 田+h 彤一h 舾一h 巧(3)这样针对n 个未知风量和n 个未知风压的风网解算问题列出了2 n 个方程，其中h 厨，h 所，巧为常量。选择任意一棵生成树，按余支在前树支在后的顺序排列各分支，将式(1)代入：Q=C k 7 Q c，则有：n-m+l劬=i 善c#q i=n m+2，乃一m+3，n(4)把式(3)代入式(2)，去掉h 尉，h 肿，有：Z(Q c)=i c g(oI 毋I 毋一 巧)=o，i=1，2，n m+1(5)把式(4)代人式(5)，得到凡一m+1 个非线性方程组：Z(Q。)=Z(g。，q 2，q。一川)，i=1，2，1,一m+1(6)式中：g，第j 条分支的风量，m 3 s；6。基本关联矩阵曰后中的元素；c。基本回路矩阵傩中的元素；r f 第条分支的风阻，k g m 7；i 第条分支的风压降，P a；彤调节风窗阻力，P a；危朋自然风压，P a；n 风机风压，P a；Q 分支风量向量，m 3 s；Q。余支风量向量，m 3 s。收稿日期：2 0 0 7 0 6 0 r 7作者简介：陈梅芳(1 9 8 3 一)，女，湖南浏阳人，硕士研究生，主要从事矿山通风网络优化研究，E m a i l：z i q i o n 9 0 7 3 1 1 2 6 万方数据矿业研究与开发2拟牛顿算法及M A T L A B 编程计算实例2 1 拟牛顿算法原理拟牛顿法的主要思想是拟定一组余支的初始风量Q。，使得初拟风量不满足风压平衡方程，即Z(Q。)0，存在一组风量增量Q。，使得：(Q。+Q。)=0，按台劳级数展开得：z(Q c+Q c)=z(Q c)+羔她+差她+g。l+0o Q co2=ou Y t l n l+l忽略二阶无穷小，则有：矾矾矾O q lo,1 2却。一。+1矾。a g l矾。向2矾。a g。+lA q l9 2一。+I工。五正一。+。(7)即：，Q c=一，(8)所以有：Q。=一(厂1)F(9)由式(9)解出的第J|次风量校正值Q：，对Q：进行校正：Z X O c “=畦+Q：(1 0)反复迭代，直到m a x IZ(Q：)I 0 0 1)&(k 1 0 0)j=s u b s(J，q l，q 5，q 6 ，Q c(1)，Q c(2)，Q c(3)；一 万方数据陈梅芳，等：基于I I A T L A B 的两种通风网络解算方法的编程及实现6 7f=s u b s(F，q 1，q 5，q 6 ，Q c(1)，Q e(2)，Q e(3)；d e t a Q e=一i n v(j)宰f：Q c=Q e+d e t a Q e；k=k+1：e n dQm-C k 牛Q e；Q c=v p a(Q e，6)Q=v p a(Q，6)在命令窗口运行此M 文件后，显示结果为：Q=7 5 6 0 8 5，1 2 1 5 0 8，1 0 6 3 4 7，2 2 7 8 5 5，1 2 9 0 6 9，9 8 7 8 5 7，2 2 7 8 5 5。余支选择，基本关联矩阵、基本回路矩阵的生成参照文献 3 o如果风网简单，不采用计算软件求解的话，拟牛顿法当然也不失为一个好的迭代计算方法，因为余支及回路的选择对计算结果没有明显影响，收敛速度快。如果用I A T L A B 编程实现，也比较简单且算法稳定，但是由于用到M A T L A B 的符号运算功能，计算速度较慢。3基于M A T L A B 的直接非线性方程组解法对自然分风网络，采用非线性方程组求解命令求解风网各分支风量，以分支风量(个)为未知数，按式(1)列m 一1 个节点风量平衡方程，按式(2)，lm+1 个回路风压平衡方程即可求解出风网的n个未知风量值。针对上述算例，列非线性方法求解的源程序如下：f u n c t i o nY2 r e(x、m=5，n=7；B k=00 0100 1；000 11 10；1100 10O；一10100 10 ；x=x(1)x(2)x(3)x(4)x(5)x(6)x(7)jR=0 50 0 0 50 0 0 40 0 0 50 0 10 0 20 ；I-I f=00 0 00 05 0 0 ；C k=1 1100 00；0101 101；O01 1011 ；列个风量平衡方程f o r i=1：(m 一1)y(i)=B k(i，：)x；e n d列个风压平衡方程f o ri=m-nj=1：nf o r t=2：np(1)=C k(i 一4)，J(1)(R(j(1)x(j(1)2)一H f(j(1)；p(t)=C k(i-4)，j(t)(R(j(t)x(j(t)2)一H f(j(t)；p(1)=p(1)；p(t)=p(t)+p(t 一1)；y(i)=p(t)；e n de n dY=y(1)y(2)y(3)y(4)y(5)y(6)y(7)；然后在命令窗口执行：】【0=1 02 01 0 01 01 09 01 0 0 ；任意给定分支的初始风量f s o l v e(7 f c，1 02 01 0 01 01 09 01 0 0 )；分支1 7 的自然分配风量计算结果：a n s=7 5 6 0 9，1 2 1 5 0 8 5，1 0 6 3 4 6 5，2 2 7 8 5 5 0，1 2 9 0 6 9 3，9 8 7 8 5 7，2 2 7 8 5 5 0。上述程序具有通用性，通风网络分支、节点越多，就越能体现其优越性。对不同的风网，只需在M文件中修改反映风网信息的参数，9 1 和m、基本关联矩阵、基本回路矩阵、风阻向量、风机风压向量。如果通风网络的节点和分支较少，也可以直接列出7个非线性方程，然后用劬l v e 命令求解即可。4结束语(1)用M A I I A B 语言编程解算风网比其他语言简单；拟牛顿法是风网解算诸多方法中收敛速度较快、易于用计算机编程实现的一种方法；基于I I A T I A B 的非线性方程组法，在实现复杂风网非线性方程组的计算机生成后，有可能是最快捷，最直接的风网解算方法。两种方法的计算结果一致。(2)在实际计算中，考虑到自然风压的影响，应把自然风压当成风机风压迭加到风机风压向量中。(3)本文的方法可用于有调节设施的风网的风量计算，对通风网络优化计算提供了一有效的方法。(4)本文提供的源程序均在M A T I A B 6 5 版本下运行通过，可供参考。参考文献：1 徐竹云矿井通风系统优化原理与设计计算方法 M 北京：冶金工业出版社，1 9 9 6 2 王沫然M A I I A B 与科学计算(第二版)M 北京：电子工业出版社。2 0 0 5 3 陈梅芳，丁德馨基于N A r I A B 的矿井通风网络图的矩阵表示及电算方法 J 矿业研究与开发，2 0 1 7 7，2 7(6)：6 5 6 7 万方数据基于MATLAB的两种通风网络解算方法的编程及实现基于MATLAB的两种通风网络解算方法的编程及实现作者：陈梅芳，丁德馨，莫勇刚，张传飞，Chen Meifang，Ding Dexin，Mo Yonggang，Zhang Chuanfei作者单位：南华大学,核资源与安全工程学院,湖南,衡阳市,421001刊名：矿业研究与开发英文刊名：MINING RESEARCH AND DEVELOPMENT年，卷(期)：2008,28(2)被引用次数：2次 参考文献(3条)参考文献(3条)1.陈梅芳;丁德馨 基于MATLAB的矿井通风网络图的矩阵表示及电算方法期刊论文-矿业研究与开发 2007(06)2.王沫然 MATLAB与科学计算 20053.徐竹云 矿井通风系统优化原理与设计计算方法 1996 引证文献(2条)引证文献(2条)1.舒才.郁钟铭 COSMOSFloWorks在矿井通风系统可视化仿真中的应用期刊论文-采矿技术 2010(5)2.李宗翔.王德民.温永宇 矿井3D风网图及基于MATLAB仿真编程实现期刊论文-安全与环境学报 2010(6)本文链接：http:/