有限体积方法

有限体积方法（精选十篇）

有限体积方法 篇1

目前分析柱体绕流的数值方法主要有基于有限体积法 ( FVM) 、有限单元法 ( FEM) 、有限差分法 ( FDM) 和格子Boltzmann方法( LBM) 等,国内外学者利用传统方法和新兴方法在这方面有很多的研究。Sohankar[3]探究了低雷诺数下,以不同入射角绕流的方柱涡脱落形态以及阻塞比、边界条件对流场的影响。Geller[4]用格子Boltzmann方法、有限单元法和有限体积法对层流流动进行了基本的计算, 比较了LBM与FEM、FVM所得的结果,讨论了LBM算法的有效性。刘天成等[5]利用LB并行算法对方柱绕流进行了计算,探究了雷诺数对流场参数的影响。Arumuga[6,7]用格子Boltzmann方法计算了二维不可压缩流稳态和非稳态下的钝体绕流,研究了雷诺数、阻塞比和流场长度对方柱和圆柱绕流的影响。 Regulski等[8]用格子Boltzamnn方法和谱元法对比研究了钝体绕流流场中涡街特性。陈静涛[9]用Fluent软件计算了二维湍流模型,通过各项流体参数分析流场中圆柱的尾涡特性。虽然学者们分析了阻塞比、立柱位置及流场长度对流场的影响,以及用CFX和LBM计算了钝体绕流。但是基于FVM讨论流场阻塞比,立柱位置和流场长度对绕流的影响的研究较少,也很少具体的对比分析LBM与FVM的结果; 且在用Fluent进行计算的时候,入流边界的速度大都采用的是恒定的值。本文用格子Boltzmann方法和有限体积法分别对方柱绕流进行了数值模拟计算。根据非平衡态外推格式对方柱角点进行LB拐角边界处理,与以往Fluent设置的恒定速度不同, 计算采用能降低流场影响的抛物线形式的来流速度; 并基于有限体积法探究了不同边界条件对方柱绕流结果的影响,通过阻力系数的对比,可以看出传统方法数值模拟要求的流场边界与LBM计算要求的流场边界条件有一定的区别。最后,将两种方法所得的结果进行了对比分析。

1数值方法

1.1控制模型

有限体积法是基于N-S方程,计算的是流体的宏观行为,其控制方程[10]为:

而基于分子动理论的格子Boltzmann方法属于介观方法,考虑的是大量微观粒子的平均行为,计算时采用具有单松弛时间格式的二维九速度( D2Q9) 模型,通过Boltzmann-BGK方程[11]计算流体质点从而获得整个流场宏观物理量。

式( 2) 中,ri表示某一时刻流体质点的位置,eα表示9个粒子速度; τ 为松弛时间,与黏度 υ 有关[2],τ = ( 6υ + 1) /2; feαq为fα对应的平衡态分布函数,表示为:

式( 3) 中,cs为格子声速,有,c 取 1; ωi为权系数。从而可以知道流体模型的宏观密度、速度

格子Boltzmann-BGK方程可以通过ChapmanEnskog展开导出N-S方程[12]。LBM算法采用的是格子单位,而FVM用到的是物理单位,本文中根据何雅玲教授[10]提到的单位转换方法进行转换。其中,无量纲参数阻力系数Cd、升力系数Cl和斯托劳哈数St的计算均是不需要转换,它们分别[7]为:

式( 6) 中,Umax为来流速度最大值; D为特征值尺寸; T为涡脱落周期,Fx、Fy分别为立柱受到的阻力和升力。

1.2LB拐角边界处理

格子Boltzmann模型的边界采用具有二阶精度的Half-way反弹格式,并结合Machado[13]提出的边界拐角处理方式,这种处理能准确计算出流场分布函数,提高数值精度。由于进出口边界是开放的,只有非曲面的方柱存在拐角,因此,只需在方柱的4个拐角处进行处理,采用内插值和外插值格式。如图1所示,我们知道此时刻拐角处的流体质点的速度uc和分布函数f1、f2和f5。因此通过邻近粒子插值得到拐角点处粒子的密度 ρ ,这与平直边界处粒子密度的计算不同。然后根据fα的对称性以及非平衡态反弹格式处理得到余下几个方向的分布函数。为了得到宏观参数和LB粒子群之间较为准确的结果,设f1= 0 ,并满足

即可算出拐角点处粒子的宏观参数。

2FVM模型及边界影响研究

2.1几何模型

绕流中,钝体与流场之间有一层很薄的边界层, 一些流体参数在该边界层影响下剧烈的变化。因此,在对流场模型进行网格划分的时候,采用分块网格画法,并对方柱周围进行加密,减少边界层影响, 保证结果的可靠性和准确性。方柱几何模型如图2所示,其中区域长度为L,高度为H,方柱边长为D, 阻塞率B = D/H ,立柱中心距上边界为H/2 ,立柱位置为S 。

水流速度定义为来流速度Ux和垂直于来流速度的Uy。流场入口边界采用入流速度形式,出口边界采用自由流形式,流场的上下边界、钝体的边界及内部均为固壁无滑移边界。在自由流中,选择合适的来流速度形式和阻塞比有助于减少流场的不稳性带来的影响[14]。因此,选择抛物线形式U = 4Umax( H y) y / H2的来流速度对雷诺数Re = 100进行FVM数值计算,采用UDF编译码导入流体计算软件Fluent中,并与横定值进行了比较,所得结果如图3,表1所示。

当来流的速度为一恒定值时,方柱的Cd值为1. 554,抛物线时为1. 389,从表1可看出,后者与已有文献值更吻合,说明了抛物线形式的来流速度比恒定的来流速度更具有正确性。

2.2进出口边界模拟结果分析

为了研究基于FVM下流场边界条件对钝体绕流的影响,现借助Fluent软件对单方柱进行数值模拟,分别研究模型的横向规模L、阻塞比B、与入流边界距离S对流场的影响,所得结果如图4、图5所示。

控制流场模型横向长度为L = 50D ,方柱中心与入流边界的距离S = 12. 5D ,探究阻塞比B =1/8,1 / 10,1 /15,1 /20,1 /25下的方柱阻力系数和升力系数。图5( a) 是方柱平均阻力系数随不同阻塞比的变化曲线。从图中可以看出,同一雷诺数下,阻塞比对柱体绕流的影响较大,且随着阻塞比的减小,阻力系数先减小后增大,这一变化规律与文献[17]中实验研究变化规律一致。当阻塞比B = 1 /8时,方柱升力系数( 图4) 震荡得比较厉害,B = 1 /10的时候相对其他的阻塞比很早就开始震荡,因此较高的阻塞比有利于抑制流场的不稳定性[7]。

图5( b) 是控制流场模型横向长度为L = 50D , 阻塞比B = 1 /10,模拟方柱中心与入流边界的距离S = 5. 5D ,10D ,12. 5D ,25D不同工况下的方柱绕流。从图5( b) 中Cd变化曲线可以看出,S = 5. 5D情况下的方柱阻力系数最大,S = 10D、12. 5D、25D时的阻力系数大小十分接近,但能看出S = 12. 5D时的阻力系数最小,即所研究中所得的最小值,并且能有效减小模型尺寸,节省CPU计算时间。

图5( c) 是控制阻塞比B = 1 /10,方柱中心与入流边界的距离S = 12. 5D ,分别取不同区域长度L = 30D,40D,50D,60D,得到的不同工况下方柱绕流的阻力系数,获得了计算域大小对阻力系数的影响。从图中Cd变化曲线可以看出,模型的横向规模L对流场的影响不是很大。由于考虑到模型边界对绕流的影响,在取模型的时候,尽量往大的尺寸取, 但是太大又会增加计算的时间,降低了工作效率。

3LBM模拟及结果对比分析

3.1压力等值线分析

从图5讨论的结果来看,所讨论的几种情况下L = 50D、B = 1 /10、S = 12. 5D这组边界条件更能有效减少FVM模拟下边界对流场的影响。而已有文献[14]表明B = 1 /8的LBM更能提高流场的稳定性。根据相关研究文献[6],本文选择LB模型尺寸40D × 8D ,抛物线式的来流速度,最大速度ULmax= 0. 1,计算所得压力如图6和图7所示。

结果显示整个流场中方柱上游部分的压力最大,且压力等值线近似椭圆形分布,如图6所示,取y1= H /2 ,得某一时刻下、一个周期内和同一时刻不同雷诺数下y1线上压力的变化情况[图7( a) 、 ( b) ]。可以看出y1线上的压力随x增加呈指数增加,在前驻点上达到最大值; 在一个周期内压力也几乎相等,变化趋势一致; 同时,随着雷诺数的增加,压力稍微有所增加,但变化趋势一致,在x = 300 ~ 360之间压力几乎相同。因此可以知道,当流场达到稳定状态之后,立柱上游的压力不会受到涡脱落影响,且雷诺数对压力较大值的影响也不是很大。

3.2速度云图对比分析

为研究流场不同位置下的速度,特选取三个位置的速度变化,比较两种方法所得的结果,如图8所示。图8给出了LBM和FVM两种方法在某一时刻所得的速度云图,该时刻一漩涡正从方柱的上边缘脱落,且两种方法模拟下的涡街形态保持一致。在方柱的中心位置、涡正脱落处以及涡已脱落处,我们分别标定x = x1,x2,x3。在这三个位置处我们得到了速度随y向距离的变化曲线和y1= H /2处速度在x方向上的速度变化,如图9所示。

从图9( a) 中我们可以看出在图8中x1处,此时刻下方柱迎流面的上下角点处有两个涡正在分离。 因此,在方柱的内部速度最小( velocity = 0) ,而在方柱的上下边缘速度达到最大值( velocity > 0. 005 m/ s) ,且速度变化曲线关于流场中轴线( y1= H /2 ) 对称。该现象符合方柱绕流流场变化规律。在x = x2处,关于y1= H /2上下流场的对称性不如x = x1处对称性好,因此该位置处速度变化较为复杂,在流场下半部分漩涡处的速度最大,而在上半部分涡脱落处再次达到第二最大速度。因此我们可以总结出, 在漩涡处的速度最大,而在流场其它位置处速度保持均匀,在x = x3处的速度变化也正好吻合这一现象。从图9( b) 可以知道中轴线上速度随x向距离的增加先减少后增大,最后趋于动态稳定。结合图8,发现以前驻点为中心向前扩展的椭圆形内的速度呈抛物线型减小,在中心处减到最小。在椭圆至来流入口处的速度保持一致,为入口最大速度( velocity = 0. 005 m / s2) 。

3.3无量纲参数的对比分析

图10分别是基于FVM与LBM两种方法模拟所获得的单方柱的升阻力系数以及相应的升力频谱图,从图10( a) 中可以看出,同一雷诺数、相同抛物线型的来流速度下的数据变化曲线不是严格的正余弦变化,而FVM的数据变化就较为均匀,平稳。这主要是因为LBM基于的是分子动理论,是通过计算大量分子的平均行为得到的宏观量,而FVM直接计算的是宏观参数,因此更具有连续性、均匀性。从图10( b) 中发现两种方法所得的方柱升力系数积分值均接近零,但LBM获得的方柱涡脱落升力系数振荡幅值稍小于FVM所得的结果,且前者能更快达到平衡,并在图10( c) 中的频谱图中得到了验证。通过表2的数据可以看出,两种方法计算所获得Cd、St值与已有文献对比,结果是比较吻合的,从而证明了两种方法在模拟钝体绕流是准确的。

4结论

基于微观分子动力论与粘性流体理论,本文分别用格子Boltzmann方法和有限体积法模拟了雷诺数Re = 100下的方柱绕流。用有限体积法研究了不同来流速度形式对流场的影响并探讨了不同边界条件对结果的影响,用格子Boltzmann方法探讨了方柱上游压力变化情况; 并比较了两种方法所得的结果。结果显示:

( 1) 抛物线形式的来流速度下的流场更为精确,有限体积法在阻塞比B = 1 /10时边界对流场影响较小,而格子Boltzmann方法是在B = 1 /8 。

( 2) 流场达到稳定后,方柱上游的压力变化呈指数分布,且不受涡脱落的影响,雷诺数对压力影响不大。

( 3) 方柱附近漩涡上的速度最大,大于入口来流的最大速度,在远离方柱的涡街上的速度越来越小,该变化与低雷诺数下的流场涡脱落变化规律一致。

( 4) 基于有限体积法的软件模拟所得结果的精度要受模型网格的影响,而格子Boltzmann则不需要划分网格,因此在网格处理上占有优势。且在相同的时间步内,后者的计算效率更高。

摘要：基于黏性流体理论,分别采用格子Boltzmann方法(LBM)和有限体积法(FVM)建立了黏性流场中方柱绕流模型。探究LBM在非光滑曲面钝体绕流方面的应用;并结合FVM进行对比分析。在FVM模型中,采用局部加密的方法对钝体边界进行处理,而在LBM模型中,除了传统的Half-way边界处理方法,还结合了拐角边界处理方法。为获得较好的可对比数据,根据已发表文献中的理论及UDF编译码技术,分别对两模型的进出口边界条件进行了讨论和设置。对比分析了两模型下的速度云图以及获得的升、阻力系数,Strouhal数。结果发现方柱上游压力不受涡脱落影响,雷诺数对其影响也较小;两种方法下的速度、无量纲参数吻合较好;但两者最适进出口边界不同,且相同条件下,LBM比FVM数值模拟能更快达到稳定状态。

有限体积方法 篇2

亚式期权定价问题的交替方向迎风有限体积方法

考虑亚式期权定价问题的数值求解.对亚式期权定价问题给出了恰当的`边界条件,并提出了一类加权迎风有限体积格式和相应的交替方向格式.对价格漂移占优问题,采用加权迎风技术以避免数值解的非物理震荡;同时,结合亚式期权定价问题特性提出了相应的交替方向格式.从理论上严格证明了提出的格式满足极大值原理,得到了一致误差估计.

作 者：孙鹏 赵卫东 SUN Peng ZHAO Wei-dong  作者单位：山东大学,数学与系统科学学院,山东,济南,250100 刊 名：山东大学学报（理学版）  ISTIC PKU英文刊名：JOURNAL OF SHANDONG UNIVERSITY(NATURAL SCIENCE) 年，卷(期)：2007 42(6) 分类号：O241.82 F830.9 关键词：亚式期权定价   交替方向   有限体积   最大值原理   迎风方法  

浅谈大体积混凝土施工方法 篇3

关键词：大体积混凝凝土施工

0引言

混凝土板厚大于1.0m，混凝土强度等级大干C30，属于典型的大体积混凝土。这种大体积混凝土底板施工具有水化热高、收缩量大、容易开裂等特点，故底板大体积混凝土浇筑应作为一个施工重点和难点认真对待。大体积混凝土施工重点主要是将温度应力产生的不利影响减少到最小，防止和降低裂缝的产生和发展。

1原材料的选择

水泥优先选用普通硅酸盐水泥，粗骨料选用碎石或卵石，粒径5～30mm，含泥量≤1％，具有良好的形状，质地坚硬，细长和片状颗粒不多于10％，不舍风化颗粒。细骨料选用中砂，含泥量≤1％，色泽均匀、干净，细度模数控制在2.7左右，砂率为30～45％。选用细度小、颜色浅、含碳量低、质量稳定的优质1级粉煤灰以及木钙减水剂，用以降低早期水化热，增强后期强度。

2施工要点及做法

混凝土配合比的选择：大体积混凝土配合比可按参照下表或经试验室试配确定。根据规范要求，每立方米混凝土中水泥用量≤550kg。水化热控制在T3D≤230kJ／kg，T3d≤270kJ／kg。砼中有效含碱总量小于5kg／m³，氯离子含量≤水泥重量的0.06％。施工应进行温度应力计算，混凝土水化热温升最高控制在50℃以下，大体积混凝土入模温度严格控制在18℃以下。为满足泵送和施工操作要求，要求混凝土塌落度为160～180±20mm。为保证基础底板水平流水、立体交叉施工，要求混凝土满足初凝时间3～4h。对混凝土的运输要求：混凝土由集中搅拌站运送到现场，时间不得超过0.5h，期间严禁加水。混凝土到工地后，要取样测定塌落度，塌落度达不到入泵要求时，根据配合比要求添加高效减水剂，严禁加水。混凝土搅拌运输车到工地后必须在1h内泵送完毕。大体积混凝土浇筑：将大体积混凝土分成面积大致相等的作业区，按作业区分层退坡浇筑。每作业区浇筑时，严格控制浇筑顺序。采用分段定点，一个坡度，分层浇筑，顺序推进，一次到顶的浇筑方法。混凝土浇筑按分层退坡前进，振捣棒设前后两排，前排振捣棒振捣浇筑点混凝土，后排振捣棒振捣斜坡处混凝土。在构件边角处，采用振动模板的办法解决构件表面的蜂窝麻面。振捣棒插入点间距不大于40cm，并插入下层10cm，每孔振捣时间不宜少于10～15s，不得超过30s，以混凝土泛浆和不冒气泡为准。振捣棒应快插慢拔，使混凝土充分密实以保证混凝土密实性。为防止大体积混凝土表面出现塑性裂缝，浇筑完成后分三次抹压成型，最后一遍用铁抹子搓平表面，全部过程保证在混凝土终凝前完成。当混凝土结构有后浇带时，在后浇带两侧混凝土强度达到100％后，再用高一个强度等级的混凝土浇筑后浇带，内掺15％UEA-M复合膨胀剂。混凝土试块留置。试块现场取样，现场制作。按每作业区100m³取一组，同时混凝土集中搅拌站亦留置与现场同样组数的试块，以便检测及查找问题。大体积混凝土的养护：混凝土最后一遍抹压完毕，随退随盖塑料布，随铺草袋或棉毡，防止气温过高混凝土表面失水过快而产生收缩干裂，达到保温保湿的目的，养护时间不少于28d。侧壁混凝土浇筑完毕后，直接在模板上淋水养护，拆模后及时覆盖1层塑料薄膜，在薄膜上懸挂麻袋片保水养护，养护时间不少于28d。

3控制温度裂缝措施

3.1合理选择配合比，降低水化热温度。水化热温度高，为此使用大粒径粗骨料，严格控制砂、石级配和含泥量，在混凝土中掺加木钙减水剂和粉煤灰等，优选混凝土配合比，以减少水泥用量，使水泥用量减到320kg／m³，降低水化热温升，同时浇筑底板时在混凝土中加10％左右的块石，既节省混凝土，又起吸热降低水化热温升的作用。

3.2降低混凝土入模温度。为降低浇筑温度，采用地下低温水中加入适量冰屑、石子洒水冷却、砂表面护盖等方法降低搅拌温度，尽量缩短混凝土运输时间，混凝土中掺加木钙缓凝剂，使初凝时间延长到6h以上，减缓浇筑速度，并薄层浇筑，通风机强制通风，以加快浇筑期间热量的散发，推迟水化热峰值出现，延长混凝土升温期。

3.3适当配置温度构造筋。对结构薄层和应力集中部位，加设φ12细密钢筋网片，在洞口四周应力集中部位适当增配构造筋，以分散应力。

3.4控制拆模时间，根据测温结果，若混凝土拆模后的表面温度或大气温度与混凝土内部温度差小于25℃，即可拆模；若降低后的表面温度或大气温度与混凝土内部温度差大于25℃时，不仅不能拆模，还应采取模板上覆盖保温材料的保温措施，减小温差。

混合网格和谐有限体积法 篇4

非结构有限体积法因为具有清晰的物理概念和较强的几何灵活性,得到迅速发展和广泛应用。但由于地形处理(即通量梯度项与底坡源项的平衡)是有限体积模型中的一个难点,直接影响到计算结果的合理性和模型稳定性。简单的方法处理会破坏格式的和谐性并影响稳定性。

对于通量梯度项与底坡源项的平衡,很多学者进行了深入研究,并提出了不同的处理方法,并将能够平衡通量梯度项与底坡源项并保持稳态解的格式统称为和谐格式。比较国内外底坡源项处理方法,效果较优异的主要有以下几种:周建国[1,2]提出“水面梯度法”(SGM,Surface Gradient Method),许为厚和潘存鸿提出“水位方程法”(WLF,Water Level Formation)[3],Audusse等[4]提出了一种静水重构法处理源项,并与其他方法比较[5,6,7],静水重构法因为未在算法中使用雅克比矩阵而使得程序编码比较容易[8]。对于以上3种方法,经过综合比较,本文选用静水重构法进行底坡源项处理。

本文最终将静水重构法与Roe有限体积法进行耦合,建立了和谐Roe有限体积计算模型。并采用了一些文献中经典算例,对模型进行了验证和分析。计算模型构建中,计算网格采用三角形与四边形混合网格。

1 二维浅水控制方程

浅水水流具有水深尺度远小于平面尺度,垂向流速小;假设流速沿水深方向平均,并沿水深方向进行积分简化方程,简化后的方程即为二维浅水方程,守恒形式如下[9,10]:

式中:q=[h,hu,hv]T,为守恒物理向量;F(q)=[hu,hu2+gh2/2,huv]T,为x方向通量向量;G(q)=[hv,huv,hv2+gh2/2]T,为y方向通量向量;h为水深;u为x方向流速;v为y方向流速;g为重力加速度;b(q)=[0,gh(sox-sfx),gh(soy-sfy)]T,为源项向量;sox=-z/x,soy=-z/y,分别为x和y方向的底坡源项;z为地面高程;sfx,sfy分别表示x和y方向底摩擦源项。

sfx,sfy采用曼宁公式计算,计算式为:

式中:n为糙率系数;其它符号意义同上。

方程中忽略了科氏力和风力的影响。

2 有限体积法数值离散

2.1 有限体积法

离散节点采用网格中心式,在网格形心布置计算节点,物理变量在节点上的值有网格平均的含义。将二维浅水控制方程(式(1))在控制体(见图1)进行积分可得:

式中:Ψ为控制平面域(面积为A);dA表示面积分微元;其他符号意义同上。

对应用高斯-格林公式转换为沿边界线积分[11],并将q,b(q)采用控制体内平均值代替,式(4)可转化为:

式中:A为控制体面积;Ψ为控制体边界(逆时针方向);dL为线积分微元;α为边界外法向向量与x轴交角;q为物理量向量在控制体内的平均值;b(q)为源项在控制体内的平均值;其他符号意义同上。

设边界上的法向数值通量Hn(q)=F(q)cosβ+G(q)sinα,可得有限体积的基本方程:

式中:m为控制体的总边数;Li为第i条边的边长;Hn(q)i为第i条边的平均法向数值通量;其他符号意义同上。

2.2 坐标旋转变换

Spekreijse[12]证明了F(q)、G(q)在控制体每条边上具有旋转不变性,并且推导出如下关系:

式中:q,qn分别为全局坐标系x-y下和局部法向坐标系ξ-η下的守恒物理量;Hn(q)为控制体边平均法向数值通量;Fn(qn)为守恒物理量qn在局部法向坐标系ξ-η下沿ξ向通量向量;T(α),T(α)-1为转移矩阵和逆矩阵,表示如下:

其他符号意义同上。

将式(7)代入式(6)可得:

根据以上分析,整个浅水方程的最终离散格式可表示如下:

式中:n为时间序号;Δt为时间步长;其他符号意义同上。

3 法向数值通量计算

Roe平均物理量表示如下[13,14,15]:

式中:hL,hR分别为边界面左右边水深;unL,unR分别为局部法向坐标系ξ-η下ξ方向边界面左右边流速;vnL,vnR分别为边界面局部法向坐标系ξ-η下η方向边界面左右边流速。

Roe格式黎曼通量计算式如下:

4 源项离散

根据文献[4]中静水重构法,局部坐标系下守恒物理量能够被重构如下:

控制体每条边的静水水深能够被重构如下:

式中:(zb)LR为边界计算点处地面高程,为了计算时水深总是正值,边界计算点处地面高程采用逆风方法确定其值:

在当前有限体积计算框架下,底坡源项差分能够非常容易的加入到数值通量方程[式(12)]中,因为浅水方程的旋转不变性,仅仅底坡源项x向分量应该被加入数值通量中,表示如下:

上面方程中qLR-HRM,qLR+HRM表示在方程式(12)替换qnL,qnR去计算控制体边界数值通量。

摩擦源项差分如下:

5 实例计算及分析

5.1 算例1:静水问题

为了验证模型二维和谐性,本文选用文献[16]中算例,验证格式的在一个不规则地形完全正方形区域静水稳定情况。正方形区域为[0,1]×[0,1],底部地形为二维非平底,底部高程计算式如下:

初始水深为:

因此,计算区域内水面水位为常数1,x和y方向动量设为零:

计算区域初始状态时静止的,在任何其他时刻都应该保持静止,如果计算格式不具有和谐性,将产生非物理波。为了检验本文中格式的和谐性,文中分别对区域采用50×50,100×100和200×200计算网格。计算到t=0.1 s时,计算区域水位与初始水位最大偏差见表1。显然,本方法的误差是由于程序编程采用VB.NET,数据类型采用双浮点型,在计算机计算中有截断误差引起的,与文献[16]相似,远小于文献[17]的误差。

5.2 算例2:不规则地形缩放型渠道问题

为了验证模型在不规则地形情况下,计算格式是否和谐稳定,本文选用了文献[18]中的经典算例———缩放型开放渠道。渠道长度取3,沿长度中心对称缩放,缩放段长度为1(见图2),外边缘线用下式定义:

本算例中,渠道最小宽度为0.9 m,计算区域为[0.3]×[-yb(x),yb(x)](见图2)。

计算中河底地形有两个高斯椭圆形岛,河底地形高程计算式如下(见图2):

计算网格采用100×300,两端平行段采用正四边形网格,中间缩放段采用三角形网格(如图3)。计算中初始条件与边界条件采用与文献[18]中的取值一样。具体取值如下:

初始水深为:

初始佛汝得数取2,x方向左端入口边界条件设定佛汝得数等于2,x方向右边出口边界取水位值等于1,y方向的边界条件设定为反射边界(陆地边界)。

根据上述输入数据,采用本文方法进行模拟,计算水位等值线面见图4。从图中可以看出,算法成功的捕获了缩放及不规则地形综合引起的复杂波纹及扰动,计算的水位波纹非常平稳且合理,对比文献[18]中格式计算的水位等值线图(图5),可以清楚的看到,本文的计算结果与其是明显一致的。

6 结语

本文采用静水重构法处理通量梯度项与源项的平衡,并且耦合Roe有限体积法,建立了和谐有限体积法。该方法可以适用不规则地形,并且对于单一网格类型和混合网格类型(如本文中三角形网格与四边形网格混合),都能进行计算,并且保持良好的和谐型和稳定性。

文中采用了一个浅水动力学中经典算例对格式进行了检验,模型计算网格采用三角形网格与四边形网格混合网格。根据模型计算结果分析,该和谐有限体积法较好的模拟了算例,表现出良好的稳定性和和谐性,并且与一些文献中其他算法的模拟效果表现出良好的一致性。

本文中计算格式尚存在以下不足,格式在计算精度上尚需进一步提高,并且在干湿边界处理上处理不足,针对这些缺点,有待下一步对算法进行进一步研究改进。

大体积混凝土施工方法论文 篇5

大体积混凝土结构裂缝主要包括干燥收缩裂缝、塑性收缩裂缝、安定性裂缝、温差裂缝、碳化收缩裂缝等。

1.1收缩裂缝。混凝土在逐渐散热和硬化过程中会导致其体积的收缩，对于大体积混凝土，这种收缩更加明显。如果混凝土的收缩受到外界的约束，就会在混凝土体内产生相应的收缩应力，当产生的收缩应力超过当时的混凝土极限抗拉强度，就会在混凝土中产生收缩裂缝。影响混凝土收缩的主要因素主要是混凝土中的用水量、水泥用量及水泥品种。混凝土中的用水量和水泥用量越高，混凝土收缩就越大。水泥品种对干缩量及收缩量也有很大的影响，一般中低热水泥和粉煤灰水泥的收缩量较小。

1.2塑性收缩也是大体积混凝土收缩一个主要来源。在水泥活性大、混凝土温度较高或者水灰比较低的条件下，混凝土的泌水明显减少，表面蒸发的水分不能及时得到补充，这时混凝土尚处于塑性状态，稍微受到一点拉力，混凝土的表面就会出现分布不规则的裂缝。出现裂缝以后，混凝土体内的水分蒸发进一步加快，于是裂缝迅速扩展。所以在这种情况下混凝土浇筑后需要及早覆盖养生。

1.3温差裂缝 混凝土内部和外部的温差过大会产生裂缝。温差裂缝产生的主要原因是水泥水化热引起的混凝土内部和混凝土表面的温差过大。温差的产生主要有三种情况：

第一种是在混凝土浇筑初期，这一阶段产生大量的水化热，形成内外温差并导致混凝土开裂，这种裂缝一般产生在混凝土浇筑后的第3天。

第二种是在拆模前后，这时混凝土表面温度下降很快，从而导致裂缝产生。

第三种情况是当混凝土内部温度高达峰值后，热量逐渐散发而达到使用温度或最低温度，它们与最高温度的`差值即内部温差。这三种温差都会产生裂缝，但最严重的是水化热引起的内外温差。

1.4安定性裂缝 安定性裂缝表现为龟裂，主要是由于水泥安定性不合格而引起。

2、裂缝的防治措施

2.1设计措施。

2.1.1精心设计混凝土配合比。

2.1.2增配构造筋，提高抗裂性能。

2.1.3避免结构突变产生应力集中。

2.1.4在易裂的边缘部位设置暗梁，提高该部位的配筋率，提高混凝土的极限抗拉强度。

2.1.5在结构设计中应充分考虑施工时的气候特征，合理设置后浇缝。

2.2原材料控制措施。

2.2.1尽量选用低热或中热水泥，在条件许可的情况下，应优先选用收缩性小的或具有微膨胀性的水泥。减少混凝土内的拉应力，提高混凝土的抗裂能力。

2.2.2适当搀加粉煤灰。混凝土中掺用粉煤灰后，可提高混凝土的抗渗性、耐久性，减少收缩，降低胶凝材料体系的水化热，提高混凝土的抗拉强度，抑制碱骨料反应，减少新拌混凝土的泌水等。

2.2.3选择级配良好的骨料。应选择线膨胀系数小、岩石弹模较低、表面清洁无弱包裹层、级配良好的骨料，尽量采用中砂，严格控制砂、石子的含泥量。控制水灰比在0.6以下。还可以在混凝土中掺缓凝剂，减缓浇筑速度，以利于散热。

2.2.4适当选用高效减水剂和引气剂。

2.3施工方法控制措施 大体积混凝土施工时内部应适当预留一些孔道，在内部通循环冷水或冷气冷却，降温速度不应超过0.5℃～1.0℃/h。对大型设备基础可采用分块分层浇筑，以利于水化热散发和减少约束作用。还应加强混凝土的浇灌振捣，提高密实度。尽可能晚拆模，拆模后混凝土表面温度不应下降15℃以上。尽量采用两次振捣技术，改善混凝土强度，提高抗裂性。

2.4温度控制措施 混凝土温度和温度变化对混凝土裂缝是极其敏感的。当混凝土从零应力温度降低到混凝土开裂温度时，混凝土拉应力超过了此时的混凝土极限拉应力。因此，通过应降低混凝土内水化热温度和混凝土初始温度，减少和避免裂缝风险。

3、大体积混凝土施工措施

大体积混凝土施工重点主要是将温度应力产生的不利影响减少到最小，防止和降低裂缝的产生和发展。因此本工程采取如下施工措施：

3.1混凝土配合比。考虑到水泥水化热引起的温度应力和温度变形，在混凝土级配及施工过程中选用32.5普通硅酸盐水泥，山砂，山石。外加剂采用WG-I高效复合防水剂，在混凝土中掺入水泥重量0.8%，初凝时间控制在12～14h。掺入粉煤灰，以替代部分水泥用量，推迟混凝土强度的增长，从而减少水泥水化热的不利影响。采用华能Ⅱ级粉煤灰，细度符合国家现行标准的规定。施工期间，要根据天气及材料等实际情况，及时调整配比，并且应避免在雨天施工。提高混凝土抗拉强度，保证骨料级配良好。控制石子、砂子的含泥量不超过1%和3%，且不得含有其他杂质。

3.2温度控制。为控制好混凝土内部温度与表面温度之差不超过25℃，施工中必须采取控温措施。

3.2.1尽量降低混凝土入模浇筑温度，必要时用湿润草帘遮盖泵管。

3.2.2为防止混凝土表面散热过快，避免内、外温差过大而产生裂缝，混凝土终凝后，立即搭设大棚进行保温养护，大棚保温养护时间根据测温控制，当混凝土表面温度与大气温度基本相同时（约4～5d），撤掉大棚保温养护，改为浇水养护。浇水养护不得少于14d；保湿保温养护措施：先铺一层塑料布，上面铺二层草帘子，根据温差来决定草帘子的增加量。

3.3浇筑方案。本工程地下室底板尺寸较大，为防止冷缝出现，采用泵送商品混凝土，施工时采取斜面分层、依次推进、整体浇筑的方法，使每次叠合层面的浇注间隔时间不得大于8h，小于混凝土的初凝时间。现场采用2台砼输送泵。要求施工队准备四组人，结合现场具体浇筑实际情况调动，要求一定确保每一下料口混凝土能很好地覆盖上层已浇筑的混凝土，避免形成冷缝。采用2台砼泵从两侧开始后退浇筑，保证浇筑的整体性。浇筑完核心筒后，整体推进，直至完成。

3.4方案可行性计算。混凝土浇筑方法为斜面分层布料方法施工，即“一个坡度、分层浇筑、循序渐进”。混凝土泵车采取“z”形行走路线。经现场与搅拌站技术人员研究，泵浇筑速度20m3/h。混凝土初凝时间为10～12h。在整个浇筑过程中，混凝土泵车循环最大距离为9m。

3.5混凝土的振捣。在每一下料口，三个振捣手均匀分布在整个斜面，沿图示中小箭头方向推进，确保不漏振，使新泵出的混凝土与上一斜面混凝土充分密实地结合。振捣应及时、到位，避免混凝土中石子流入坡底，发生离析现象。混凝土采用机械振捣棒振捣。振捣棒的操作，要做到“快插慢拔”，上下抽动，均匀振捣，插点要均匀排列，插点采用并列式和交错式均可；插点间距为300～400mm，插入到下层尚未初凝的混凝土中约50～100mm，振捣时应依次进行，不要跳跃式振捣，以防发生漏振。每一振点的振捣延续时间30s，使混凝土表面水分不再显著下沉、不出现气泡、表面泛出灰浆为止。为使混凝土振捣密实，每台混凝土泵出料口配备4台振捣棒，分三道布置。第一道布置在出料点，使混凝土形成自然流淌坡度，第二道布置在坡脚处，确保混凝土下部密实，第三道布置在斜面中部，在斜面上各点要严格控制振捣时间、移动距离和插入深度。

3.6混凝土表面处理。大体积混凝土的表面水泥浆较厚，且泌水现象严重，应仔细处理。对于表面泌水，当每层混凝土浇筑接近尾声时，应人为将水引向低洼边部，处缩为小水潭，然后用小水泵将水抽至附近排水井。在混凝土浇筑后4～8h内，将部分浮浆清掉，初步用长刮尺刮平，然后用木抹子搓平压实。在初凝以后，混凝土表面会出现龟裂，终凝要前进行二次抹压，以便将龟裂纹消除，注意宜晚不宜早。

大体积方形钢煤斗结构的有限元分析 篇6

煤斗结构是火力发电厂煤仓间的一个重要构筑物。根据工艺要求,煤斗一般采用圆形和方形,其应力和变形分布均匀,受力性能好。但是,有时根据工艺要求也采用方变圆的结构形式,这种煤斗由于几何外形存在突变,在外荷载的作用下,煤斗局部会出现应力集中。截至目前,还没有这类形状煤斗明确的应力和变形的计算公式。传统的计算方法只是将钢煤斗各个构件拆分开,每个受力部位均采取近似的受力模型进行力学计算,竖壁按粮仓贮仓近似计算。这种方法虽然安全度较高,但计算繁琐。而利用有限元分析,能准确地得到各个部位的应力和变形值。本文结合湖南华润鲤鱼江发电厂工程,介绍了利用有限元软件SAP2000对钢煤斗进行模拟分析的过程,并将有限元计算结果与传统手算的结果进行了对比分析。

1 煤斗结构特征

华润鲤鱼江电厂煤斗由长方体形和锥形两部分组成,上部截面为11 288 mm×5 644 mm,均采用钢结构,煤斗外形可以分为两部分,上方长方体几何容积399 m3,下方锥体漏斗每个几何容积67.567 m3。该煤斗筒体属高壁浅仓,煤斗中部竖壁与斜壁连接部分为钢筋混凝土大梁,竖壁为混凝土梁支撑,煤斗为钢吊斗。煤斗外形示意图见图1,煤斗的有限元模型见图2。

2 传统理论计算分析

根据《贮仓结构设计手册》,矩形钢仓斜壁计算

过程如下:

(1)斜壁任一水平截面单位宽度上的斜向拉力设计值Ni按(1)式计算:

undefined: (1)

式中:αi——斜壁与水平面夹角;Nvi——斜壁相应水平截面单位宽度上的竖向拉力。

(2)斜壁弯矩设计值和挠度计算。水平加劲肋和垂直加劲肋可视为斜壁的弹性支点,将斜壁划分为双向板或单向板。双向板可近似按四边固定板计算:

undefined (2)

undefined (3)

式中:undefinedn——该块板平均法向压力标准值;α,β——系数,可按四边固定板查表;a,b——板短边和长边的边长;t——斜壁壁板厚度;E——钢材的弹性模量。

(3)一肢水平布置的加劲肋计算。水平加劲肋承受的相邻斜壁引起的水平拉力设计值为:

undefined (4)

undefined (5)

式中:qa,qb——水平加劲肋承受的法向线荷载标准值。

一肢水平布置的加劲肋可按封闭框架计算。漏斗布置对称时,弯矩设计值为:

支座:undefined (6)

跨中: Ma=γQqa·aundefined/8-M0 (7)

Mb=γQqb·bundefined/8-M0 (8)

式中:C=bh/ah。

加劲肋最大挠度可按(9)式和(10)式计算:

undefined (9)

undefined (10)

式中:ξ——系数,与加劲肋两端弯矩有关,可查表求得。

该工程中煤斗下部斜壁横截面不是规则的矩形,而是由方形渐变为圆形,在计算下部斜壁受力时采取了近似处理,按方形进行计算,偏于保守。

贮料物理特征:重力密度γ=11 kN/m3,内摩擦角φ=30°,k=tg2(45-30/2)=0.333,γ0=1.0,c=1.0。设计时假设贮料装至仓顶,钢仓水平加劲肋间距为:3.63/5=0.726(m)。

贮料压力按(11)式计算:

ξ=cos2α+k·sin2α (11)

ξ=cos270°+0.333sin270°=0.411

斜壁法向压力标准值按(12)式计算:

P=ξ·c·r·s (12)

斜壁顶部:P上=0.411×1.0×11×6.264=28.32(kN/m2)

斜壁底部:P下=0.411×1.0×11×9.894=44.73(kN/m2)

水平加劲肋处Pna及壁板中间平均值undefinedna见图3。

详细的计算过程从略,下面给出计算结果。

(1)斜壁壁板最大内力设计值:

4～5板:Mmax=2.94 kN·m,N=74.20 kN/m;

0～1板:Mmax=2.33 kN·m,N=197.23 kN/m。

(2)水平加劲肋最大内力设计值:

Nh=76.27 kN,M0=62 kN·m,Ma=32.54 kN·m。

煤斗竖直侧壁计算方法较斜壁计算简单,通常取单位宽度的竖直侧板,水平加劲肋视为支座进行计算,这里不再赘述。

3 有限元模型分析

钢煤斗的设计包括两部分:煤斗壁和加劲肋。模型中,梁柱采用frame单元,钢煤斗壁采用能够承受平面外弯矩的厚壁shell单元。在分析时,必须对通过预埋件连接起来的混凝土构件与钢煤斗壁之间进行耦合,在连接位置建立耦合对,即可确保混凝土构件与钢煤斗协调变形,加劲肋与煤斗壁之间通过焊缝连接。

(1)面截面的绘制。煤斗竖壁均为矩形面,煤斗斜壁包括三角形面和圆锥面,矩形面及三角形面可直接通过“绘制面单元”命令绘制,圆锥面则主要通过“拉伸点成线”、“拉伸线成面”和“合并点”命令绘制,见图4。

(2)为了能更准确地模拟混凝土大梁对钢煤斗的约束作用,采用实体单元对混凝土大梁进行模拟。

(3)水平加劲肋的布置。将竖壁截面沿高度方向剖分成5段,在截面剖分处布置加劲肋HM 250 mm×175 mm×7 mm×11 mm(按手算设计);将斜壁沿高度方向剖分为5段,在截面剖分处布置加劲肋160 mm×100 mm×12 mm(按手算设计)。图5为煤斗加劲肋布置示意图。

(4)约束和荷载。按照实际的支承情况,煤斗支撑在混凝土大梁上,大梁与柱刚性连接,故在模型上对柱底施加固定约束。

按照GB 50077—2003《钢筋混凝土筒仓设计规范》,煤斗可根据浅仓进行加载,考虑煤的重力密度为11 kN/m3,内摩擦角为30°,冲击影响系数为1.3。原竖壁高6.264 m,斜壁高3.63 m,煤斗按满载计算,计算结果如下:

面压力:竖壁底部:fb=22.97 kN/m2。

斜壁表面压力:斜壁顶部:fb=28.32 kN/m2;斜壁底部:fu=44.73 kN/m2。

竖壁表面压力为三角形,斜壁表面压力为梯形,在sap2000中分别通过节点样式P1和P2来布置。其中:P1=-3.663z+116.11,P2=-4.521z+143.15。图6和图7为有限元计算结果示意云图。

通过对上述有限元模型进行自动网格剖分并求解分析,得到了计算结果。从应力分布(见图8)来看,斗壁板应力分布较大的区域在中部过渡段,一般为170～230 MPa,在竖壁与斜壁连接的位置应力集中比较严重,最大应力约为232 MPa,其他板的应力均<170 MPa。整个煤斗的最大位移在位于煤斗三角形斜板最上面一道加劲肋的中部,其位移值约为12 mm。

从整体上看,有限元模型与实际情况比较接近,计算所得的应力值与位移值是合理的。在竖壁与斜壁连接部位,应力集中比较严重,但其区域比较集中,由于钢材具有良好的塑性变形能力,可以引起应力重分布,从而缓解应力集中,所以,这种较小区域上的应力集中,不会对煤斗的整个结构安全造成影响。具体施工时,应采取适宜的焊接措施,控制焊缝质量。

4 结果对比分析

将传统理论计算结果(手算结果)与有限元计算结果进行对比,结果见表1。由表1可以看出,理论计算结果均比有限元结果偏大,且竖壁部分的理论计算结果与有限元计算结果的偏差较斜壁部分更大。竖壁部分的最大板应力、弯矩与挠度差别均在50 %以上;斜壁部分最大板应力吻合较好,水平加劲肋弯矩相差比竖壁部分还要大。

竖壁部分理论计算与有限元计算结果偏差较大的最主要原因在于两者计算模型的差异。在理论方法中,计算竖壁部分采用的模型是选取单位宽度的侧板将其简化成受梯形或三角形均布荷载的扁梁来计算其内力。这样建模便于手算,但是,忽略了整个侧板作为一个整体受力与局部单位宽度板受力的差异,实际上由于两侧板对计算模型仍有约束作用,也即计算模型平面外的约束条件在手算中未予考虑,平面外刚度的影响也被忽略,所以,计算结果要较有限元分析结果偏差较大,理论计算结果偏于保守。两种计算结果中有关竖壁部分水平加劲肋的差别集中在跨中弯矩上,支座处的最大弯矩两者相差不多,而跨中的最大弯矩手算结果为48.03,有限元结果为21.2,相差比较大。其原因在于理论计算中水平加劲肋是按照封闭钢架进行计算,长轴方向的两个垂直加劲肋对水平加劲肋的约束作用在这里默认为刚性支座,不产生水平位移,但实际情况是这两个垂直加劲肋是会产生水平位移的,不是作为固定支座,在有限元分析中垂直加劲肋的位移被考虑了进去,导致理论计算的跨中弯矩值要大于有限元计算结果,这也是一种偏于保守的计算简化。

斜壁部分对板应力的计算结果两种方法相差不大,最大的差别在于水平加劲肋的内力计算。因为在手算中无法考虑斜壁部分横截面由矩形渐变为圆形的过程,默认横截面一直为矩形,而在有限元模型中则完全可以考虑这种影响,水平加劲肋在四角部位均为圆弧形,所以,斜壁部分的水平加劲肋弯矩结果相差较大,理论计算结果将方变圆简化为完全方形,导致结果偏大,也偏于保守。有限元计算的水平加劲肋位移比理论计算结果大的原因,在于此处水平加劲肋正好布置在壁板最大受力处;而竖壁的水平加劲肋位移小于理论计算的原因,是由于水平加劲肋布置的位移不是很合理,未布置在壁板最大受力处。

5 结论

5.1

对理论计算与有限元计算的对比可以发现,理论计算均采用了偏于保守的简化模型,其计算结果普遍比有限元分析结果偏大。这种保守处理虽然在实际工程中是偏于安全的,但也在一定程度上造成了浪费,增加了工程费用。而且理论计算模型过于简化,不能准确地计算出某些特定部位的内力,例如本工程中斜壁部分的截面变形问题。

5.2

有限体积方法 篇7

有限元法和有限体积法作为两种不同的数值解法, 有着独特的原理和求解思路, 本文对这两种方法分别进行讨论:

下面以Poisson方程第一边值问题为例:

其中G为一多边形域, ) 83 () 82 (G为G的边界。f2L () 为已知函数。

1 有限元法

有限元法是计算流体力学经常采用的一种主要方法。它是以变分原理或加权余量法为基础, 吸取有限差分法的离散化思想而发展起来的一种数值方法。有限元法是采用分块逼近的思想而形成的一种数值计算方法。在使用有限元法中, 如何选择合适的单元形状及其插值函数是很重要的。

以以下问题为例, 有限元法的定义如下。

把边值问题化成变分形式:求uH01 (7) G (8) , 使得

其中

选取有限维子空间VnH01 (7) G (8) , 求unVn, 使得

取样条函数作为Vn的基函数, 即得有限元法。

2 有限体积法

2.1 积分插值法

假定是多边形区域。将做三角剖分Th, 三角形单元记为KhT, h表示所有三角形的最大边长。剖分Th称为原始剖分。还要构造与hT相应的对偶剖分hT*。常用的方法外心对偶剖分, 如图4所示, 设任一三角形单元的内角不大于90°。Qi是三角单元的外心, Mi是各边中点。分别连接外心点形成的多边形, 即围绕内节点0P的一个对偶剖分。

对于问题 (1) - (2) , 积分插值法的构造如下:

如图4, 在对偶单元上积分 (1) , 并利用Green公式, 得

可得差分格式

有限体积法 (Finite Volume Element Methods) , 又称为广义差分法, 为积分法的推广。

2.2 有限体积法及差分格式

求解Poisson方程的有限体积法定义为:求uhUh, 使得

或等价于

其中

对于均匀剖分 (即轴和轴方向的步长均为) , 则得五点差分格式:

3 小结

有限元法的基本思想是分块逼近, 它将连续的求解域离散后, 在各个子域上构造插值函数逼近微分方程的解。在满足微分方程和相应初边值条件下, 对全部子域进行积分, 总体合成, 并建立有限元方程组, 求解该方程组就得到各节点上待求的函数值。有限元法对椭圆型问题有更好的适应性, 它的计算量较有限体积法大。

有限体积法又称为控制体积法, 该方法同有限差分法相似, 通常将计算区域划分为网格, 寻求函数在网格节点上的近似解。它对每一个网格节点周围取一个控制体积, 将待解的微分方程对每个控制体积积分, 从而得到一组离散方程。其中的未知数是网格节点上的因变量。因此就网格划分和求解而言, 有限体积法可视为离散方法不同的一种特殊的有限差分法。该方法吸取了有限差分法以及有限元发的有点, 是一种对区域剖分具有单元特征、能适应复杂边界条件、离散方法简单灵活等特点。与有限元法相比, 有限体积法具有占用内存少、计算效率高以及数值稳定性好等方面的优点。

摘要：本文以Poisson方程第一边值问题为例, 对常见的流体力学数值解法有限元法与有限体积法进行了分析。

有限体积方法 篇8

在工程燃烧设备的传热过程中, 辐射换热约占90%[1]。随着经济的发展, 对锅炉炉膛传热计算实用性、精确性要求也不断提高。辐射换热是大型锅炉炉膛内的主要换热形式, 因而准确的计算炉膛内的辐射换热量对大型锅炉设计和优化有重要意义。

近几十年来, 对于吸收、散射性的炉内介质的辐射换热, 各国学者做了大量的研究, 已发展了许多的辐射换热的计算方法, 主要有蒙特卡洛法、区域法、离散坐标法、有限体积法等[2,3]。离散坐标法和有限体积法具有很好的适应性较高的求解效率, 在大型锅炉炉膛辐射换热分析中很有前景。其中有限体积法易用于非结构网格, 同时满足局部守恒性, 因而更具优势。

本文将有限体积法推广用于求解和分析大型电站锅炉炉膛内的辐射换热过程。首先给出了有限体积法对辐射传递方程进行离散和求解的基本过程。然后对该方法进行了必要验证, 最后用于分析某电厂600 MW锅炉炉膛内的辐射换热。

1 辐射换热有限体积法

对于大型炉膛, 由于存在大量煤、灰粒子, 因而散射效应较强。对于考虑充满吸收、发射、散射性介质内的辐射传递, 辐射传递离散方程可以写为[2]

undefined

对于不透明、漫射边界, 给定如下边界条件

undefined

使用有限体积法对辐射传递方程进行空间离散, 对每一个空间控制单元内积分并使用散度定理可以得到

undefined

其中Np为体元p的边界数, NV为体元总数。使用阶梯格式处理单元表面的辐射强度

(sm·ei) Iundefined=max{sm·ei, 0}Iundefined-max{-sm·ei, 0}Iundefined (4)

方程 (3) 可以写为如下形式

undefined

上式整理为

KpIundefined=Hp (6)

其中

undefined

具体的求解过程:对于每一个离散方向, 逐单元p求解方程直至收敛, 这样就可以获得炉膛内的辐射强度场。

2 方法验证

首先通过一个简单算例对所给出的有限体积法进行验证。考虑二维方形炉膛内的辐射换热。炉膛的边长为10 m, 壁面为黑壁面, 且不考虑煤、灰散射。炉膛光学厚度为1, 下壁面 (y=0) 温度为1 000 K, 其余壁面为500 K。这里使用限体积法 (FVM) 对求解辐射平衡时炉膛内的温度场。有限体积法空间离散采用400个四边形单元 (如图1) , 角度离散为Nθ×Nφ=20×20个方向。

图2给出了有限体积法求解得到的沿底边中线 (x=5) 的温度分布并与文献中区域法 (Zone) 的计算结果进行了比较。有限体积法的计算结果和文献中的结果吻合的很好, 最大相对误差小于2%。为了更好的了解炉膛内的温度场, 图3中给出了该方形炉膛内的温度场分布。可以看出, 由于辐射换热的影响, 随着远离热壁面, 温度逐渐降低。由于炉膛几何形状的影响, 最高温度出现在热壁面的中心处, 而最低温度出现在热壁面正对壁面的两个拐角处。这些结果表明有限体积法应用于炉膛内辐射换热计算具有很好的精度。

3 炉膛辐射换热的应用

下面将有限体积法用于分析某电厂600 MW锅炉内的辐射换热量。炉膛的尺寸外形见图4, 其沿垂直纸面方向的深度尺寸为19 558 mm。炉膛内的温度分布如图5所示, 在靠近燃烧器区有一个温度峰值。水冷壁的温度为600 K。辐射计算时考虑炉膛内煤灰的吸收和散射作用。图6给出了有限体积法进行辐射换热计算所使用的网格, 整个三维炉膛划分为5 304个六面体网格。

通过有限体积法对辐射传递方程进行求解, 得到三维炉膛内的辐射强度场。为了分析水冷壁上的辐射能照射强度, 图7给出了求解得到的壁面上辐射温度分布, 辐射温度表征了投射辐射能量的强度。从图中可以看出, 燃烧器附近区域辐射温度最高, 因而辐射投射能量强度最大。图8进一步给出了壁面上的辐射换热热流, 可以看出, 壁面吸收的辐射热流分布趋势与投射辐射强度分布类似, 即均在燃烧器附近出现峰值, 最大辐射热流约为200 kW/m2。在燃烧器的上部和下部随着远离燃烧区辐射换热热流逐渐降低。这均与锅炉的实际运行工况相一致。说明有限体积法可以用于实际大型电站锅炉炉内辐射换热的计算和分析。

4 结论

本文将有限体积法推广应用于求解和分析大型电站锅炉炉膛内的辐射换热。给出了有限体积法对辐射传递方程进行离散和求解的基本过程, 将有限体积法用于分析某电厂600 MW锅炉炉膛内的辐射换热, 验证了有限体积法求解大型电站锅炉炉膛辐射换热的可行性。结果表明有限体积法可以有效的求解大型电站锅炉炉膛内的复杂辐射换热过程。

参考文献

(1) 王应时.燃烧过程数值计算 (M) .北京:科学出版社, 1986.

(2) 冯俊凯, 沈幼庭, 杨瑞昌.锅炉原理及计算 (第三版) (M) .北京:科学出版社, 2003.

(3) 沈珞蝉, 夏建军, 姚强, 岑可法, 池作和.区域法在大型电站锅炉炉膛传热数值计算中的应用 (J) .浙江大学学报 (工学版) .1996, 30:7-17.

(4) Modest, M.F., Radiative Heat Transfer, 2nd ed., AcademicPress, San Diego, CA, 2003.

有限体积方法 篇9

本文根据已实测大体积混凝土温度发展曲线, 匹配养护混凝土, 测得各龄期混凝土的力学性能;并根据经实测温度校正的ANSYS数值分析的温度场作为计算温度应力场基础, 最后利用劈裂抗拉强度与计算得到的温度应力第一主应力比值, 即开裂指数, 来评价大体积混凝土由于早期温度应力引起的开裂。

1 试验原材料及试验方案

1.1 原材料及性能

采用甘肃祁连山水泥集团股份有限公司生产的P•O42.5普通硅酸盐水泥;砂的细度模数2.6, Ⅱ级区中砂;粒径为5mm~20mm, 连续级配, 其中5mm~10mm占40%, 10mm~20mm占60%;外加剂为高效减水超缓凝泵送剂;国电兰州热电有限责任公司的粉煤灰, 经0.045mm方孔筛筛余8.1%;酒钢磨细矿渣粉, 经0.08mm方孔筛筛分, 其筛余为3%。

1.2 混凝土的养护制度

选择养护曲线原则是承台最容易开裂表面处 (预模拟) 混凝土的温升曲线作为本匹配养护的温升曲线。匹配养护制度是指实测得到大体积温度曲线, 在可调节温度养护箱内, 根据温度发展规律在实验室内进行跟踪养护的混凝土养护制度。具体养护制度, 即实验室内养护的温升曲线见图1。

混凝土表面的温度在浇筑完成后, 36h内的温度上升速率很快, 之后速率减缓, 77.5h达到峰值, 此后逐渐下降。

1.3 混凝土配合比和试验方法

混凝土轴心抗压强度试验、混凝土劈裂抗拉强度试验和混凝土静力受压弹性模量 (简称弹性模量) 试验用试件形状、尺寸以及试验方法依据《普通混凝土力学性能试验方法标准》GB/T50081-2002。配合比见表1。

2 结果与分析

图2和图3中可见, 1d~5d匹配养护混凝土的早期弹性模量和混凝土的劈裂抗拉强度增长很快, 但静弹性模量较混凝土劈裂抗拉强度增长速度要慢。在5-28d这段龄期内劈裂抗拉强度和静弹性模量均增长很小, 仅劈裂抗拉强度0.32MPa, 静弹性模量为3.6GPa。静弹性模量急剧增长, 相应的混凝土温度应力也急剧增长。因此, 在研究混凝土的早期开裂, 要将早期混凝土静弹性模量的快速增长视为一个重要的控制因素。从优化配合比开始, 如利用静弹性模量的较小的骨料, 适当地增加胶凝材料的量, 利用大掺量矿物掺合料, 掺加超缓凝剂使静弹性模量早期缓慢发展, 防止在早期由于干缩和温度突变产生早期裂缝, 另外将混凝土早期徐变完全发挥出来。一方面, 由于工程上养护不及时或不浇水养护, 干缩、外部温度骤降;使混凝土的内表综合应力增加, 使混凝土表面开裂;另一方面, 混凝土内部由于温度降低产生收缩, 而外部温度变化比内部更大, 以及湿度梯度引起的收缩, 也是表面裂缝的一个主要原因。而大体积混凝土对养护的要求比其他混凝土更敏感, 因此, 控制温度裂缝, 良好的外保温和保湿是非常关键的。

3 大体积混凝土温度应力场的ANSYS分析

混凝土温度应力来源于混凝土受约束时热变化而导致的体积变化。在大体积混凝土中, 水泥水化放出热量, 混凝土导热系数小, 而和空气对流换热的换热系数大, 从内部到表面产生温度梯度, 越靠近表面处梯度越大, 尤其是表面到深度为50cm时, 温差尤为可观。此温度梯度导致的应力足以使大体积混凝土产生裂缝。因此, 本文进行应用ANSYS软件温度梯度导致的温度应力场的分析。本ANSYS有限元模型分析中, 采用的是Solid70单元, 其具有热和结构耦合分析的单元。该模型中全部采用六面体八节点单元, 整个模型单元总数为4004, 节点总数为4830。

根据大体积混凝土结构对称性, 取尺寸6m×14m×2.6m为模拟混凝土块。图2中实测静弹性模量, 进行线性插值将混凝土的静弹性模量非线性的作为载荷加入到应力计算模型中。由于承台原有底部具有灌注桩, 下表面为全约束, 另外前两个表面处于对称结构的中心, 也认为是全约束, 所以其计算的约束条件为:下表面、前两个表面全约束。

混凝土浇筑后1d、3d、5d、7d、14d和28d的S1温度应力云图, 如图4所示。

从图4可见, 混凝土温度应力S1在3d时达到最大, 此后逐渐减小。1d、3d、5d、7d、14d和28d承台表面混凝土温度最大温度应力S1分别为:0.57MPa、4.46MPa、4.23MPa、3.79MPa、1.26MPa和0.86MPa。整个承台的混凝土S1出现在下表面靠近四周, 各个龄期承台最大温度应力S1分别为:5.38 MPa、5.17MPa、4.86 MPa、3.29 MPa和2.55 MPa。

4 混凝土开裂的确定

工程结构物由于各种原因而丧失其正常工作能力的现象, 称为失效。对于脆性材料, 强度失效的形式为脆性断裂。脆性断裂, 指材料失效时未产生明显的塑性变形而突然断裂的现象。因大体积混凝土属于脆性材料, 且温度应力为拉伸型和拉应力占优的混合型应力状态, 适合于用第一强度理论来评价温度应力能否导致大体积混凝土开裂。材料发生断裂的主要因素是最大拉应力, 不论何种应力状态, 只要其最大拉应力σt max达到极限拉应力σu, 材料就会发生断裂。在复杂应力情况下σt max=σ1, 断裂准则为σ1=σb, 其中σb为断裂破坏的强度极限。图4中的S1为第一主应力, 即σ1。因此, 劈裂抗拉强度和S1进行比较, 只要S1大于劈裂抗拉强度, 混凝土一定会开裂。匹配养护实验得到各龄期混凝土劈裂抗拉强度和ANSYS计算得到各龄期的S1比较见图5。

图5中可以看出, 0-3d随着龄期增加, 承台表面最大温度应力增加;3d-28d, 随着龄期增加, 承台表面最大温度应力S1逐渐下降。而混凝土劈裂抗拉强度随着龄期增加逐渐增大, 3d内增加速率很大, 后期缓慢增长。混凝土劈裂抗拉强度和S1对比得到:3d时fts/S1值最小为1.12, 此后逐渐增大, 说明3d是混凝土最易开裂时刻。因此, 3d是控制混凝土温度应力的关键时刻。此时混凝土内部温度较高, 易形成大的温差, 且此时混凝土静弹性模量也几乎到了峰值, 此时温度应力最大。

此外, 在实际工程上, 进行了近一个月的实际观察, 工程上没有发现任何肉眼可见微裂缝。说明劈裂抗拉强度大于S1可保证大体积混凝土表面不开裂, 可以作为评价混凝土表面是否开裂的判据。

最大温度应力集中在承台下表面四周, 并1d和3d温度应力均超过了混凝土劈裂抗拉强度。下面应该出现一些裂缝, 不过我们无法观测到。因此, 对于底面全约束混凝土, 其温度应力主要集中在承台底部的四周, 这些区域容易开裂, 应在以后实际的工程中着重考虑。

3结语

1) 在匹配养护条件下, 混凝土静弹性模量和劈裂抗拉强度均发展很快, 1d~5d匹配养护混凝土弹性模量和劈裂抗拉强度趋势一样, 但比混凝土劈裂抗拉强度增长速度要慢。7d值为36.56GPa, 达到28d值的96.2%, 之后缓慢增长。在5d时, 混凝土劈裂抗拉强度达5.5MPa, 之后强度增长缓慢。

2) 承台表面温度应力S1 3d时达最大, 为4.46MPa, 5d、7d、14d和28d的温度应力分别为:4.23MPa、3.79MPa、1.26 MPa和0.86 MPa。3d时fts/S1值最小, 此后逐渐增大, 说明3d是混凝土最易开裂时刻。

3) 在混凝土表面养护条件良好的情况下, fts/S1>1.12, 大体积混凝土表面不会开裂。

有限体积方法 篇10

深浅腔动静压混合轴承[1,2]具有抗过载能力强、稳定性好、对负载变化敏感性低等优点,因而被广泛地应用于高速精密磨床的静压主轴上。在设计工作中,选取深浅腔动静压轴承结构参数及工作参数的主要依据,是以油腔压力、承载能力、油膜刚度、进油流量及温升为代表的承载特性[3,4]。如果参数选取不当,会造成油膜刚度不足及温升过高,进而导致磨削质量下降,严重时可导致主轴刮擦、抱轴等事故。因此,选取合适的轴承结构参数及工作参数是获得最优承载特性的关键。

已有学者对优化深浅腔动静压轴承的承载特性进行了研究。孙恭寿等[5]采用解析法及有限差分法对环面缝隙节流深浅腔动静压轴承进行了研究,但未能进一步研究动静压轴承的优化设计问题。郭胜安等[6]采用FLUENT软件研究了深浅腔动静压轴承的工作参数和结构参数对轴承刚度、流量及温升的影响规律,并给出了参数的优选区间。Wang等[7]通过求解动态雷诺方程及能量方程得到了轴承平均温升、承载能力、油流率及涡动频率等参数,但仅研究了两个备选方案。上述研究结果表明,寻找动静压轴承优化设计的最优的参数组合,是一个多目标优化问题[8],它存在计算工作量大、用时长的困难。

本文改变了将正交试验法仅用于试验领域的传统习惯,把它移植到理论计算领域。根据平行平板扩散流动计算模型,依据控制体及其边界上流量守恒的物理本质,得到了微元控制体边界上压力的插值函数,提出了分析小孔节流深浅腔动静压轴承的油腔压力、承载能力、静刚度、进油流量及温升等承载特性的计算方法。进而用该方法研究了动静压轴承的供油压力、主轴转速、进油孔径、浅腔深度、初始油膜厚度等参数对轴承承载特性的影响规律,并确定了参数的优化区间。在此基础上,采用四因素三水平的正交试验法,在满足多目标性能最优的前提下,从理论上研究了该轴承结构参数与工作参数的最优组合。以该组参数试制了试验平台,并测量了不同转速及供油压力下油腔的压力值。

1 有限体积法

1. 1 轴承结构及工作原理

图1 为小孔节流深浅腔动静压轴承的结构图。在轴承的对称中心处建立直角坐标系,水平向右为X轴,竖直向上为Y轴,垂直纸面向外为Z轴。圆周均布的4 个油腔从右下起按逆时针依次编号为1、2、3、4。以第1 油腔为例,轴承的结构参数如下:周向封油边包角为 φb,深腔包角为 φc。深腔油腔深度为hs,浅腔油腔深度为hp。

在供油压力ps作用下,润滑油经过小孔节流器进入深腔,在逆时针旋转的轴颈的迁移作用及深浅腔之间二次节流作用下流向浅腔,并在浅腔末端形成压力远大于深腔压力的动压。最后润滑油流过周向封油边进入回油槽或从轴向端面回油。静止时,在静压效应作用下,主轴通过作用在有效承载面积的上下油腔压差承载;转动时,主轴主要靠远大于深腔油压的动压力工作。

1. 2基本假设

做如下假设:① 润滑油流动状态为层流;②在载荷作用下轴心竖直向下偏移(偏位角始终为0) ;③ 基于等温假设;④ 全部发热量均由润滑油吸收并由端泄带走,发热量与端泄散热平衡,即“绝热流动”。

1. 3物理方程

1. 3. 1平行平板扩散流动方程

一维雷诺方程积分后可得平行平板间流量计算公式[9]:

式中,Q为压差 Δp作用下流过长度为 Δl的平行平板间的流量;h为平行平板间润滑油的油膜厚度;μ 为润滑油的黏度;b为平行平板的宽度;v为运动固体壁面对润滑油的剪切速度。

1.3.2流量连续性方程

流量连续性方程为

式(2) 表明,流入某个区域或截面的流量Qin应与流出该区域或截面的流量Qout相等。

1. 4网格划分及边界条件

按油腔形状将油膜区域进行网格划分,如图2 所示。划分好的网格如图3 所示。其中,第( i,j)个节点的油膜厚度hi,j为

式中,θi为Y轴负方向沿逆时针到第(i,j) 个节点的圆周角度;h0为初始油膜厚度;ε 为偏心率。

下面分析边界条件。

(1)进油孔流量边界条件。若第(i,j)个节点处存在进油孔,则在供油压力ps与该节点压力pi,j之差的作用下,经过节流小孔后进入油腔的润滑油流量Qini,j为

式中,α为小孔流量系数,常取0.6~0.7;dc为节流小孔直径;ρ为润滑油密度。

否则Qini,j=0。

(2) 压力边界条件。计算区域边界上的压力为0。

1. 5离散格式推导

以图2 中的第(i,j) 个节点为中心的微元控制体[10,11,12,13]作为研究对象(图4)。

由平行平板扩散流动方程(式(1)) 及流量连续性方程(式(2)) 可知:

(1)第(i,j)个控制体东侧流出的流量Qei,j应等于第(i+1,j)个控制体西侧流入的流量Qwi+1,j,则

式中,pi,j为第(i,j)个节点的压力;pi+1,j为第(i+1,j)个节点的压力;pei,j、hei,j分别为第(i,j)个控制体东侧边界处压力值和油膜厚度;pwi+1,j、hwi+1,j分别为第(i+1,j)个控制体西侧边界处压力值和油膜厚度;ΔXi为X方向增量,ΔXi=RΔθi;R为主轴半径。

同理可知,在第(i,j) 个控制体西侧,有

在第(i,j) 个控制体北侧,有

在第(i,j) 个控制体南侧,有

(2)第(i,j)个控制体流入的流量Qini,j应等于其4个边界流出的流量之和,即

(3) 考虑到压力连续性条件,有

令

求解式(4) ~ 式(10) 可得:

令

则式(12) 可转换为

采用TDMA迭代法[14,15]求解由式(13) 构成的方程组,便可得到整个油腔的压力分布。

1. 6轴承承载能力

轴承X、Y方向承载能力FX、FY分别为

其中,φ 为轴承偏位角。为了仿真及计算的方便,仅研究主轴沿竖直向下运动后的承载能力,本文取 φ = 0。

1. 7轴承油膜刚度

Y方向油膜刚度k0为

式中,eY为Y方向轴心偏移量。

1. 8轴承温升

动静压轴承消耗的功率主要包括两部分:泵功耗及剪切摩擦功耗。

根据牛顿摩擦力公式,轴承的剪切摩擦力Ff为

轴承的剪切摩擦损耗Nf、泵损耗Np分别为

式中,Qi为轴承的总进油量,即4 个油腔进油量之和。

轴承温升 Δt为

式中,cV为润滑油的质量热容,取2120 J/(kg·K)。

2 基于有限体积法的承载特性分析

2. 1 初选轴承结构参数及供油参数

在理论分析时,初选的轴承结构参数及供油参数见表1。研究承载特性时,除所研究参数外,其余参数按表1 选取。

2. 2 求解结果

在图5 中,深腔区域的油腔压力为均匀分布。静止时在圆周方向,油膜压力从深腔到浅腔逐渐减小,在浅腔末端压力达到最小(图5a);工作时在圆周方向,油膜压力从深腔到浅腔逐渐增大,在浅腔末端压力达到最大(图5b);在轴向,从中间向两端压力逐渐减小。

2. 3 承载特性分析

2. 3. 1 供油压力的影响

动静压轴承是依据供油产生的静压效应及一定转速下楔形油膜的动压效应工作的。在使用过程中,轴承的结构参数是确定的,则供油压力成为制约轴承性能的关键工作参数。当供油压力为1. 5 ~ 4. 0 MPa时,承载能力、油膜刚度、进油流量和温升随偏心率 ε 的变化规律如图6 所示。

由图6 可知,当偏心率一定时,承载能力、油膜刚度及进油流量均随供油压力的增大而增大,而温升却随供油压力的增大而减小。这表明增大供油压力对提高轴承的承载特性是有利的,但是供油压力越大,越易导致进油管道的泄漏。因此,在保证足够的刚度、承载能力及温升条件时,供油压力取2. 5 ~ 3. 5 MPa较为合适。

2. 3. 2 主轴转速的影响

主轴转速是动静压轴承另一个重要的工作参数。高速旋转时产生的动压效应可以提高轴承的刚度。当偏心率为0. 1 ~ 0. 4 时,不同主轴转速下,轴承承载能力、油膜刚度、进油流量和温升的变化规律如图7 所示。

由图7 可知,当偏心率一定时,轴承承载能力、油膜刚度、进油流量及温升均随主轴转速的增加而增大。这是因为,随着主轴转速升高,动压效应增强;深腔压力减小,进油流量增加;剪切摩擦损耗增加导致温升增加。

2. 3. 3 进油孔径的影响

不同偏心率下,进油孔径对轴承承载能力、油膜刚度、进油流量和温升的影响规律如图8 所示。

由图8 可知,在某一偏心率下,随着进油孔径增大,承载能力、油膜刚度先增大后减小,进油流量不断增大,温升逐渐降低。这表明随着进油孔径增加,进油流量增大,温升降低。但进油流量越大,回油越困难,越容易造成漏油问题,所以进油孔径不宜过大,取0. 6 ~ 0. 7 mm较为合适。

2. 3. 4 浅腔深度比的影响

为了增强动压效应,提高动静压轴承的承载能力和刚度,工程上常采用深浅腔的油腔结构。浅腔深度与初始油膜厚度的比值成为衡量浅腔深度大小的重要参数。浅腔深度与初始油膜厚度之比hp/ h0分别为0. 5、1、2、4、6、12、25 时,轴承的承载能力、刚度、进油流量和温升随浅腔深度的变化规律如图9 ~ 图11 所示。其中,图9 ~ 图11 的工作转速分别为0、3000、6000 r/min。

在图9 中,当转速为0,偏心率 ε = 0. 2,hp/ h0从0. 5 增大至4 时,承载能力急剧增至2. 78 k N,刚度从365 N/μm急剧增至521 N/μm,流量从1. 96 L / min增至2. 34 L / min;hp/ h0> 5 后,承载能力、刚度及进油流量趋于平稳。温升不随hp/ h0的改变而改变,始终为1. 36 ℃。这是因为,在转速为零时轴承不存在动压效应,浅腔深度越大则静压效应越明显;由于不存在速度剪切摩擦损耗,只有泵损耗转化为热量,故温升保持不变。

图10 中,当转速为3000 r/min,偏心率 ε =0. 2,hp/ h0从0. 5 增大至1 时,承载能力从3. 28k N增大至3. 61 k N,刚度从676 N / μm增大至734N / μm,进油流量从2. 33 L / min增大至2. 55 L /min,温升从6. 47 ℃ 降低至5. 52 ℃ ;当hp/ h0继续增大至6 时,承载能力减小至2. 93 k N,刚度减小至576 N/μm,进油流量减小至2. 43 L/min,温升降至4. 56 ℃。当hp/ h0从10 增大至20 时,承载特性缓慢减弱。这表明浅腔深度比为2 ~ 3 时,承载能力、刚度及进油流量最大而温升最小;浅腔深度比大于20 时,动压效应几乎消失。为了发挥深浅腔轴承在动压效应作用下增强油膜刚度及降低温升的优势,hp/ h0取1 ~ 3 较为合适。

2. 3. 5 初始油膜厚度的影响

初始油膜厚度从0. 018 mm增大至0. 030 mm时(浅腔深度比hp/ h0取2),承载能力、刚度、进油流量和温升随油膜厚度的变化规律如图12 所示。

由图12 可知,当偏心率一定时,随着初始油膜厚度的增大,承载能力、油膜刚度及温升不断减小,而进油流量逐渐增大。这是因为,随着初始油膜厚度的增加,液阻减小,进油流量增加。流量的增加增强了端泄散热能力,从而引起温升降低。但初始油膜厚度过大,轴承承载能力及刚度过小。因此,初始油膜厚度取0. 02 ~ 0. 03 mm较为合适。

3 结合正交试验法的结构优化设计

通过以上分析,得到了轴承工作参数及结构参数的优选范围,但并未找到使承载特性最优的一组参数。若以承载能力最大、刚度最大、流量最小及温升最低这4 个指标为优化目标进行优化计算,则计算量会特别大。众所周知,正交试验法在寻找最优工艺参数时具有减少试验次数的优势[16]。所以,下文借助正交试验法进行理论分析。ε = 0. 2 时,研究初始油膜厚度、浅腔深度比、进油孔径及供油压力4 个因素的影响(表2)。

根据工程经验,动静压轴承的承载特性按照其重要程度依次为温升、刚度、承载能力及供油流量。下面从这4 个指标来研究其最优参数组合。

(1)从温升看,最敏感的影响因素为初始油膜厚度h0,最弱的因素为浅腔深度比hp/h0。为使温升最小,较优的因素组合为:ps=3 MPa,dc=0.7 mm,h0=0.03 mm,hp/h0=2。

(2) 从刚度及承载能力看,最敏感的影响因素为初始油膜厚度h0,最弱的因素为进油孔直径dc。为使刚度最大,较优的因素组合为:ps= 3MPa,dc= 0. 7 mm,h0= 0. 02 mm,hp/ h0= 1。

(3) 从供油流量看,最敏感的影响因素为初始油膜厚度h0,最弱的因素为进油孔直径dc。为使供油流量最小,较优的因素组合为:ps= 2MPa,dc= 0. 6 mm,h0= 0. 02 mm,hp/ h0= 2。

由以上分析可知,4 个承载特性指标的优化结果是相互矛盾的,需要根据其重要程度对相关参数做出取舍。因此,最终确定的参数为:ps=3 MPa,dc= 0. 65 mm,h0= 0. 025 mm,hp/ h0= 2。

4 试验分析

4. 1 试验对象

试验在湖南大学与湖南海捷主轴科技有限公司联合研制的液体悬浮电主轴试验台上进行。该试验台由液体悬浮电主轴( 图13)、液压站及JC-3000型冷却器( 图14) 、电气控制系统( 图15 )组成。根据所确定的最优参数值,试制了该款电主轴的前后轴承。经配磨加工及塞尺检测,轴承油膜的单边间隙为0. 02 ~ 0. 03 mm。液压站用来提供稳定的供油压力;冷却器用于控制电主轴的进油温度。电气控制系统用于控制该款电主轴在额定频率200 Hz以下实行恒转矩调速,而在额定频率200 Hz以上实行恒功率调速。

4. 2 试验原理

试验原理及压力测试原理如图16 及图17 所示。按照图16,将电气控制系统、液压冷却系统、液体悬浮电主轴连接在一起。试验时,首先检查电主轴各个元器件及连接系统是否正常工作;然后启动液压系统,检查主轴是否浮起;最后调节电气控制系统,在主轴自重下记录不同转速时前轴承的压力值。如图17 所示,在前轴承第2 个油腔的浅腔与封油边交汇处钻出一个测压孔,再通过管道把测压孔与压力表相连,试验时直接读取压力表的数值即可。

4. 3 试验及结果讨论

在主轴自重(主轴质量为70 kg)下,当供油压力分别为2. 0 MPa、2. 5 MPa、3. 0 MPa,转速分别为0、2000 r / min、4000 r / min时,油腔压力值见表3。

通过理论压力和实测压力的对比,发现轴承浅腔末端处的压力值随供油压力及转速的变化趋势一致,而局部数据点存在一定误差,其中最大相对误差为10. 7% 。该误差产生的可能原因如下:①由初始油膜厚度的实际值与计算时参数的设定值存在一定差距;②按经验选择的小孔流量系数与实际情况存在一定差距。但就总体趋势而言,理论计算结果与试验测试结果基本吻合。因而该试验在一定程度上验证了优化结果的正确性。

5 结论

(1)基于平行平板扩散流动计算模型及流量守恒原理,推导了微元控制体边界压力的插值函数,建立了分析小孔节流深浅腔动静压轴承的油腔压力、承载力、静刚度、进油流量及温升等承载特性的有限体积计算方法。

(2)定量研究了动静压轴承供油压力、主轴转速、进油孔径、浅腔深度、初始油膜厚度等参数对轴承承载特性的影响规律,进而提出了可增大轴承刚度和降低油膜温升的结构参数及工作参数的优选区域。

(3)采用四因素三水平的正交试验设计方法,在满足多目标性能最优的前提下,理论研究了该轴承结构参数与工作参数的最优组合。最优的一组参数为:ps=3 MPa,dc=0.65 mm,h0=0.025 mm,hp/h0=2。