新型液排渣燃烧器内燃烧特性的数值模拟-(完整版)实用资料

新型液排渣燃烧器内燃烧特性的数值模拟（完整版）实用资料(可以直接使用，可编辑完整版实用资料，欢迎下载）

新型液排渣燃烧器内燃烧特性的数值模拟新型液排渣燃烧器内燃烧特性的数值模拟（完整版）实用资料(可以直接使用，可编辑完整版实用资料，欢迎下载）冉景煜1刘丽娟2黎柴佐3通过数值模拟的方式,研究了新型液排渣燃烧器在不同过量空气系数下的速度、温度以及组分浓度的分布情况.结果表明,在较小的过量空气系数(α=0.7,0.8时,煤的燃尽情况较差;α≥1.0时,煤粉燃烧更完全,但却不利于氮氧化物的控制.采用分级燃烧的方式,控制燃烧器内为欠氧燃烧(取α=0.9以降低局部氧浓度,既能达到液态排渣要求,又可抑制NOx的生成,并在高温烟气进入炉膛降温之后再补充燃尽风,使得可燃成分在炉膛内再次燃烧,提高燃尽率.通过模拟与实验相结合的方式,对燃烧器进行三种不同负荷下的热态实验研究,该燃烧器负荷适应性好,模拟结果与实验结果相吻合.液排渣燃烧器;气固两相流;过量空气系数;燃尽率;NOx排放TQ052.7;TQ534重庆市自然科学基金重点项目(2021BA6067.重庆大学低品位能源利用技术及系统教育部重点实验室,400030重庆2021-03-292021-05-04所采用的方法,总热阻可采用等效热阻的方式代替,NOx排放量急剧增加.热,以液态渣的形式排出燃烧器.@@[1]林伯川,陈恩鉴,吴承康.煤粉旋风燃烧过程流场特性研究[J].工程热物理学报,2001,22(4:523-525.@@[2]纪任山.煤粉工业锅炉燃烧的数值模拟[J].煤炭学报,2021,34(12:1703-1706.@@[3]ChenZC,LiZQ,JingJPetal.TheInfluenceofFuelBiasinthePrimaryAirDuctontheGas/particleFlowCharacteristicsNeartheSwirlBurnerRegion[J].FuelProcessingTechnology,2021,89:958-965.@@[4]秦裕琨,孙锐,李争起等.径向浓淡旋流煤粉燃烧器气流湍流特性的冷态试验研究[J].工程热物理学报,2001,22(1:111-114.@@[5]赵伶玲,周强泰.复杂曲面花瓣燃烧器煤粉燃烧数值分析[J].中国电机工程学报,2007,27(5:39-44.@@[6]EttouatiH,BoutoubA,BentichalHetal.NumericalModellingandSimulationofPulverizedSolid-fuelCombustioninSwirlBurners[J].TheCanadianJournalofChemicalEngineering,2021,87:350-358.@@[7]张力,蔡松,蒲舸.低热值煤层气燃烧器结构设计的数值模拟[J].煤炭转化,2021,33(4:47-51.@@[8]MaL,JonesJM,PourkashanianMetal.ModellingtheCombustionofPulverizedBiomassinanIndustrialCombustionTestFurnace[J].Fuel,2007,86:1959-1965.@@[9]ShangQ,ZhangJ.SimulationofGas-particleTurbulentCombustioninaPulverizedCoal-firedSwirlCombustor[J].Fuel,2021,88(1:31-39.@@[10]唐倩,王晓琴,裴林清等.新型窑炉燃烧器的数值模拟与实验研究[J].中国电机工程学报,2021,30(26:12-16.@@[11]WangXH,ZhaoDQ,JiangLQetal.TheDepositionandBurningCharacteristicsDuringSlaggingCo-firingCoalandWood:ModelingandNumericalSimulation[J].CombustionScienceandTechnology,2021.181(5:710-728.@@[12]RanJY,LiuLJ,LiCZetal.NumericalStudyonOptimumDesigningoftheAirDistributionStructureofaNewCycloneCombustor[C].Shanghai:InternationalConferenceonEnergy,EnvironmentandSustainableDevelopment,2021.@@[13]孙学信.煤粉锅炉燃烧试验技术与方法[M].北京:中国电力出版社,2001:413-449.@@[14]周力行.湍流气力两相流动和燃烧的理论与数值模拟[M].北京:科学出版社,1994:126-195.@@[15]汪小憨,赵黛青,陈恩鉴等.煤粉低尘燃烧器内燃烧特性的数值模拟[J].过程工程学报,2005,5(5:479-484.@@[16]LuQG,ZhuJG,NiuTYetal.PulverizedCoalCombustionandNOxEmissionsinHighTemperatureAirfromCirculatingFluidizedBed[J].FuelProcessingTechnology,2021,89(11:1186-1192.@@[17]张力,张俊广,唐强等.煤层气燃烧器流动及NOx生成特性数值研究[J].煤炭转化,2021,34(3:71-74.NUMERICALSIMULATIONOFBURNINGCHARACTERISTICSOFANOVELLIQUIDSLAG-REMOVALBURNERRanJingyuLiuLijuanLiChaizuo第38卷第2期2021年北京化工大学学报(自然科学版JournalofBeijingUniversityofChemicalTechnology(NaturalScienceVo.l38,No.22021双层半圆管盘式涡轮桨搅拌槽气液分散特性的数值模拟张雪雯李志鹏高正明*(北京化工大学化学工程学院,北京100029摘要:采用基于气泡聚并和破碎机理的群体平衡(PBMMUSIG模型,对双层半圆管盘式涡轮桨搅拌槽内的气液分散特性进行了数值模拟;考察了不同通气量和操作转速下气液搅拌槽内流体流动,局部气含率和气泡尺寸的分布规律。模拟结果表明:通气工况下搅拌槽内的液相流场具有双循环流动形式;采用PBMMUSIG模型预测的局部气含率分布与文献实验数据吻合较好;搅拌槽内气泡尺寸随转速增加而减小,随气量增加而增大;桨叶排出流区域内气泡尺寸较小,近壁区和循环区内气泡尺寸较大。关键词:气液搅拌槽;数值模拟;半圆管盘式涡轮桨;气泡尺寸;群体平衡模型中图分类号:TQ0273收稿日期:2021-05-25基金项目:国家自然科学基金(20576009/20821004第一作者:女,1983年生,硕士生*通讯联系人Emai:lgaozm@mai.lbuct.edu引言具有气含率和相界面积大、容易换热和处理量大等特点的气液搅拌槽被广泛应用于化工,生化,食品等工业生产中。单层桨搅拌槽的气-液两相分散特性也已得到广泛研究,然而随着工业装置越来越大型化,仅仅使用单层桨已经无法满足工业需求。因此工业上常常选择多层桨操作来实现反应器内物料的搅拌混合。近年来多层桨的气液分散特性已有大量实验研究[1-2]。但由于实验研究手段存在投资大、周期长、测量困难等诸多缺点,如何寻找一种能够缩短开发周期、节约成本的研究手段并部分代替实验研究就显得格外令人关注。随着现代科学技术的不断进步,基于计算流体动力学(CFD的理论及方法,借助计算机进行数值模拟的预测技术在搅拌槽设计方面得到广泛应用。相对实验而言,CFD技术能以较少的费用和较短的时间为代价获得大量有价值的研究结果,并能对搅拌槽和搅拌过程进行详细描述,在一定程度上弥补了实验研究手段的缺陷。气液搅拌槽内传质和反应过程都需要通过气泡的表面进行,气泡尺寸及其分布是影响传质和反应速率的重要因素,因此,气泡尺寸分布的数值模拟对搅拌槽的研究、开发和放大都有着非常重要的意义。但由于气液两相流本身很复杂,气液两相间的相互作用、气泡在运动过程中聚并和破碎导致的气泡尺寸分布不均以及气泡对搅拌体系湍动的影响等因素都会影响到气液两相流模拟的准确性。因此,文献中的相关模拟一般简化为单相模拟[3-4],或者将气泡简化为单一尺寸的刚性球体,忽略其聚并和破碎[5-6]。近年来,已有部分研究者开始考虑气泡间破碎和聚并条件下的分布规律,模拟结果较单一气泡假设有明显改善,但主要集中在Rushton涡轮桨[7-12]。本文以文献[2]的实验数据为基础,在考虑气泡尺寸分布和气泡聚并破碎的基础上对双层半圆管盘式涡轮桨气液搅拌槽内的流场、局部气含率和气泡尺寸分布规律进行数值模拟。1搅拌槽结构和模拟工况1.1搅拌槽结构本文搅拌槽模型与文献[2]实验一致,如图1所示。搅拌槽直径T=045m,槽内设有4块挡板,挡板宽度WB=1/12T,挡板距离槽壁为1/60T,选用双层半圆管涡轮搅拌桨(CD6桨,通气装置为环形气体分布器,直径为027T,其上均匀分布12个直径为0003m的小孔,开孔向下,模拟装置各结构具体参数见表1。1搅拌轴;2挡板;3搅拌桨;4搅拌槽;5气体分布器图1搅拌槽(a及搅拌桨(b装置示意图Fig.1Tank(aandCD6impeller(bgeometries表1搅拌槽结构参数Table1Parametersofagitator/vesselconfigurationsT/mH0/TD1/TD2/TCS/TC1/TC2/T045203903902403313312模拟工况物系为自来水-空气体系,具体模拟工况见表2(其中叶端雷诺数Re=ND2/。表2模拟各工况具体操作条件Table2Detailedoperatingconditionsforeachcase工况搅拌转速N/s-1通气流量Q/Lmin-1表观气速VS/ms-1叶端雷诺数Re14213200224748!10625155000164112!10735213200224112!10745276600290112!1072数值模拟21网格划分选取整个搅拌槽内的流体作为模拟计算域,采用网格分块划分技术。为了更好地捕捉桨叶区以及气体分布器小孔处的流体运动,对这两处进行了网格局部加密。桨叶、气体分布器和槽底附近区域采用四面体网格,其他区域采用六面体网格,网格单元数约为110万,如图2所示。22模拟方法本文采用CFX软件对气液搅拌槽内流场、气含图2网格划分示意图Fig.2Illustrationsofthemeshusedinthemodel率和气泡尺寸分布等进行数值模拟研究。液相主流域和气相的流动分别采用标准k~湍流模型和零方程模型。采用多重参考系(MRF方法处理动静区域的相对运动,其中搅拌浆叶位于旋转运动区域,挡板和气体分布器处于静止区域,流动状态为定常态。相间作用力考虑了曳力,重力,表面张力以及湍动扩散力。曳力模型采用适合于稀疏气相流动的Grace曳力模型,较大通气量时,修正的体积分率指数为2。使用LopezdeBertodano模型来考虑相间的湍动扩散力,其作用是湍流涡体携带气泡由高浓度向低浓度流动,湍动扩散系数选择01。利用群体平衡(PBMMUSIG模型处理气泡的聚并和破碎,将气泡尺寸分为1~10mm的五组,每组气泡直径由式(1与(2计算di=dmin+!di-12(1!d=dmax-dminM(2其中M为气泡分组数,d为气泡尺寸。根据上述公式,得到的各组气泡平均尺寸分别为19,37,55,73和91mm。PBMMUSIG是关于数值密度函数的连续形式,其方程式(3如下tn(∀,t+xi(Ui(m,tn(∀,t=BB-DB+BC-DC(3式中n(∀,t表示尺度为∀的颗粒在时间t内的数密度。BB,DB,BC,DC分别表示破碎及聚并后气泡的生成和损失率。2北京化工大学学报(自然科学版2021年3结果与讨论31气液搅拌槽内流场分析工况2下搅拌槽内两挡板中间平面液相速度矢量图如图3所示。CD6为径向流搅拌桨,高速旋转的搅拌桨在桨叶区产生径向射流,该射流携带周围流体向槽壁运动,与槽壁撞击后,沿槽壁分为上下两股,分别在两层桨叶附近区域形成两对循环涡流。从图3可以看出,模拟计算结果成功反映出了该流动结构。32局部气含率对于气液两相流搅拌槽,表征气体在槽内分布的局部气含率是一个非常重要的性能参数,它能够反映槽内局部的分散及传质特性,在反应器的选择和设计中有着极其重要的作用。在搅拌槽内两挡板中间平面上选取一条与轴线平行直线,该直线距离槽壁0035m。将表2中4种工况下该直线处局部气含率的实验和模拟结果进行对比,如图4所示。图3工况2搅拌槽内液相速度矢量图Fig.3LiquidvelocityvectorprofileinthemidplanebetweentwobafflesforCase2∀实验值;模拟值图44种工况下轴向局部气含率分布图Fig.4Axiallocalgasholdupdistributionforthefourcases由图4可知,局部气含率模拟值总体上与实验值[2]较为吻合,即采用PBMMUSIG模型可以较准确模拟气液搅拌槽内局部气含率的分布规律。该分布规律表现为两层桨叶附近区域均存在着峰值,且顶层桨附近峰值始终高于底层桨附近峰值。这一方面是由于桨叶区的气泡受到循环流和浮力的双重作3第2期张雪雯等:双层半圆管盘式涡轮桨搅拌槽气液分散特性的数值模拟用,另一方面则是相对底层桨而言,顶层桨附近区域存在更多的气体循环。由图4可知但模拟结果与实验值仍存在一定偏差,比如底层桨附近区域模拟的气含率峰值高于实验数据,两层桨之间区域模拟的气含率分布下降很快。模拟与实验间存在偏差的可能原因如下:PBMMUSIG模型主要适用于鼓泡塔[13]的气液流动,且使用的聚并和破碎机理比较简单。相对鼓泡塔而言,当通入相同气量时,存在着机械能输入的气液搅拌槽内的气泡平均直径要小得多,因此采用该模型预测搅拌槽内的气液流动会存在偏差。此外,通过指定体积分率指数来修正Grace曳力模型也会导致一定误差。33气泡尺寸分布搅拌槽内的剪切速率与气泡尺寸分布密切相关,图5给出了工况1下搅拌槽内挡板中间平面气泡平均直径分布与剪切速率对比示意图。由图5可知气泡尺寸在剪切速率最大的叶轮排出流区内最小,在远离叶轮的位置,随着剪切速率逐渐减弱,气泡尺寸也随之变大。对于循环区域而言,其剪切速率较之排出区要小的多,所以大气泡主要存在于搅拌槽循环区域内,尤其是循环区域的中心位置。这一趋势与Kerdouss等[11]在双层六直叶涡轮桨气液搅拌槽内模拟所得结果趋势一致。图5搅拌槽内气泡尺寸及剪切速率分布Fig.5Distributionsofbubblesizeandshearstrainsrateinthestirredvessel对比搅拌槽内液相速度矢量图和气泡尺寸图可以看出,低流速区即循环区域的气泡尺寸较大,这是由于此区域内连续相速度较低且与气泡运动方向相反,因此气泡更容易聚并成大气泡。而处于搅拌槽上层区域处的气泡尺寸也较大,这主要是因为此处上升的气泡与下降的液相流相遇,从而导致气泡运动速度变慢,气泡聚并几率变大,从而尺寸变大。为了研究转速和表观气速对气泡平均尺寸的影响,本文计算了工况1~4下气泡平均尺寸d32值,分别为5536,4674,5148和5377mm。由此可知,在相同气速下,随着搅拌转速的增加,气泡尺寸减小,造成这一现象的原因是桨叶转速越高,所带来的能量耗散越高,能够产生剪切力更高,破碎能力也就更大,从而导致气泡破碎现象更加严重,气泡直径变小。当转速不变,随着通气流量的增大,搅拌槽内出现大气泡的几率也随之增大,这是因为随通气流量的增加,气泡之间更容易发生碰撞和聚并,并成为大气泡。沿挡板中间右半平面上各循环回路中液相速度矢量方向各取17点,具体循环回路从下到上依次为循环回路1,循环回路2,循环回路3,循环回路4。如图6所示。图6沿各循环回路取点示意图Fig.6Circulationloopsinthemidplanebetweentwobaffleswithsomepointsshown图7给出了工况1下,沿图6所取点所得到得气泡尺寸变化趋势图。由图7可知,当处于桨叶上方的循环回路2和4时,气泡在循环流和浮力的双重作用下运动,聚并效应比破碎效应更加明显,所以气泡尺寸沿循环方向单调增加;而处在桨叶下方循环回路1和3时,靠近槽壁处的循环区域内存在着气泡的密集地带,从而造成了气泡尺寸先增加,后减小。因为槽底部存在死区,该处的气泡尺寸较小。4北京化工大学学报(自然科学版2021年∀回路1;#回路2;∃回路3;回路4图7循环回路气泡平均尺寸变化趋势图(工况1Fig.7Bubblediameteralongtheliquidcirculationloop(case1全槽内气泡尺寸分布不均匀,桨叶区气泡较小,越靠近槽壁,气泡越大。循区域气泡较大。4结论(1对于双层半圆管盘式涡轮桨气液搅拌槽,其液相流场具有双循环流动型式。(2气液搅拌槽内轴向局部气含率的模拟结果与实验值趋势基本吻合,验证了PBMMUSIG模型和方法的可行性。(3相同气速下,气泡尺寸随搅拌转速的增加而减小;相同转速,气泡尺寸随气量的增加而变大。(4气液搅拌槽叶轮排出流区域内气泡尺寸较小,靠近槽壁处和循环区域内气泡尺寸较大。符号说明C离底高度,mD搅拌桨直径,md32槽内平均气泡直径,mmg重力加速度,m/s2H液面高,mN搅拌转速,s-1Q通气流量,L/minRe叶端雷诺数T搅拌槽直径,mVS表观气速,m/sz轴向坐标,mG局部气含率流体相密度,kg/m3剪切黏度,Pas#体积黏度,Pas参考文献:[1]MouchaT,LinekV,ProkopovaE.Gasholdup,mixingtimeandgasliquidvolumetricmasstransfercoefficientofvariousmultipleimpellerconfigurations:Rushtonturbine,pitchedbladeandtechmiximpellerandtheircombinations[J].ChemicalEngineeringScience,2003,58:1839-1846.[2]GaoZM,SmithJM,MullerSteinhagenH.Voidfractiondistributioninspargedandboilingreactorswithmodernimpellerconfiguration[J].ChemicalEngineeringandProcessing,2001,40:489-497.[3]刘敏珊,张丽娜,董其伍.圆盘涡轮式搅拌槽的数值模拟[J].郑州大学学报:工学版,2007,28(1:122-124.LiuMS,ZhangLN,DongQW.Thesimulationofrushtonturbinestirredtank[J].JournalofZhengzhouUniversity:EngineeringScience,2007,28(1:122-124.(inChinese[4]朱向哲,苗一,谢禹钧.双层涡轮搅拌桨三维流场数值模拟[J].石油化工设备,2005,34(4:26-29.ZhuXZ,MiaoY,XieYJ.Threedimensionalflownumericalsimulationofdoubleturbineimpellers[J].PetroChemicalEquipment,2005,34(4:26-29.(inChinese[5]宋月兰,高正明,李志鹏.多层新型搅拌槽内气液两相流动的实验与数值模拟[J].过程工程学报,2007,7(1:24-28.SongYL,GaoZM,LiZP.Experimentalstudyandnumericalsimulationofgasliquidflowinastirredtankwithanewmultipleimpeller[J].TheChineseJournalofProcessEngineering,2007,7(1:24-28.(inChinese[6]孙会.搅拌罐内气液两相流场的数值模拟[J].上海电机学院学报,2021,11(3:189-192.SunH.Numericalinvestigationofgasliquidflowinstirredvessels[J].JournalofShanghaiDianjiUniversity,2021,11(3:189-192.(inChinese[7]LaakonerM,MoilanenP,AlopaeusV,eta.lModellinglocalgasliquidmasstransferinagitatevessels[J].ChemicalEngineeringResearchandDesign,2007,85(A5:665-675.[8]GimbunJ,RiellyCD,NagyZK.ModellingofmasstransferingasliquidstirredtanksagitatedbyRushtonturbineanCD6impeller:ascaleupstudy[J].ChemicalEngineeringResearchandDesign,2021,87:437-451.[9]LaneGL,SchwarzMP,EvansGM.Numericalmodellingofgasliquidflowinstirredtanks[J].ChemicalEn5第2期张雪雯等:双层半圆管盘式涡轮桨搅拌槽气液分散特性的数值模拟6gineeringScience2005,60:2203-2214,.北京化工大学学报(自然科学版2021年[12]MinJ,BaoYY,ChenL,eta.lNumericalsiulationofmgasdispersioninanaeratedstirredreactorwithmultiplemipellers[J].Industrial&EngineeringCheicalRemsearch,2021,47:7112-7117.[13]SanyalJMarchisioDL,FoxRO,eta.lOnthecom,parisonbetweenpopulationbalancemodelsforCFDsiumlationofbubblecolumns[J].Industrial&EngineeringCheicalResearch,2005,44(14:5063-5072.m[10]LaneGL,SchwarzMPEvansGM.Predictinggasliq,uidflowinamechanicallystirredtank[J].AppliedMathematicalModelling2006,26223-235.,:[11]KerdoussF,BannariA,ProulxPCFDmodelingofgas.dispersionandbubblesizeinadoubleturbinestirredtank[J].ChemicalEngineeringScience2006,61:3313,3322.NumericalsiulationofgasliquidflowinamdualCD6ipellerstirredtankmZHANGXueWenLIZhiPengGAOZhengMing(CollegeofChemicalEngineeringBeijingUniversityofCheicalTechnology,Beijing100029,China,mAbstractThegasdispersioncharacteristicsinadualCD6ipellerstirredtankhavebeennuericallysiulated:mmmusingthepopulationbalancemodel(PBMcombinedwiththemultiplesizegroup(MUSImode,lwhichisGbasedontheprincipleofbubblecoalescenceandbreakupTheeffectsofvaryingthegasflowrateandrotation.speedonthefluidfieldlocalvoidfractiondistributionandbubblesizedistributionhavebeeninvestigatedThere,.sultsshowthattheliquidflowadoptsadualloopflowpatternundergasdispersionoperatingconditionsThelocal.voidfractiondistributionssiulatedusingthePBMMUSImodelareingoodagreementwiththeexperientaldatamGmfromtheliteratureBubblesizedecreasesastherotationspeedincreases,andincreasesasthegasflowrateincrea.sesThebubblesizeissmallerintheipellerdischargeregionandlargerinthenearwallandcirculationregions..m,KeywordsgasliquidstirredtanknumericalsiulationCD6ipellerbubblesizePBMMUSI:;m;m;;G复合材料探索与求实分析模块采用相同的网格划分,如图2所示,以方便它们之间的数据交换。有限元分析中边界条件及加载情况如下:热化学模块:在结构件外表面施加对流边界条件,对流换热系数由实验测定;其他面施加绝热边界条件。树脂流动模块:结构件内表面固定,外表面施加热压罐压力,树脂只能由上表面流向吸胶膜,其他面因为真空带的密封作用,没有树脂流出。1缸c,S图1AS4/3501—6标准固化周期温度与压力曲线图2有限元网格力学分析模块:在固化过程中,结构件只能沿模具表面滑动,与模具接触的结点沿模具表面法线方向位移为零。‘脱膜后,固定一点以限制z方向的刚体位移,在拐角处对称面上所有节点施加对称边界条件。3分析与讨论针对8mm厚,拐角处半径为20mm的L形层合板进行有限元数值模拟,并分别就温度、102复合材料探索与求实ｌ２ｌＯＯ８Ｏ６Ｏ４目善号ｎｐ盘∞器曰１３Ｈ卢一０２ＯＯＡＰｏｓｉｔｉｏｎＢＯ．ＯＯ．２０．４０．６０．８１．０ＲｅｌａｔｉｖｅＴｈｉｃｋｎｅｓｓ（ａ）压缩厚度（ｂ）纤维体积含量（相对厚度：内表面为Ｏ，外表面为１）图４平直和拐角处压缩厚度及纤维体积含量三种方法所得回弹角的对比表ＭｅｔｈｏｄＲａｄｆｏｒｄ’８ｐｒｅｄｉｃｔｉｏｎ【９】ＦＥＡＳｐｒｉｎｇ—ｂａｃｋ锄ｇｌｅ／。一０．５０４一０．６８Ｒｅｌａｔｉｖｅｅｒｒｏｒ／％３４．５１１．７垦堡巫望堡翌！！！４结论＝Ｑ：！！Ｑ在通用有限元分析程序ＡＢＡＱＵＳ平台上，通过添加适当的用户子程序，模拟了热固性树脂基纤维增强复合材料的固化过程。本文的计算结果与相关文献的计算及实验结果对比表明，该方法可以较好地模拟复合材料固化过程中热传导、反应放热、树脂流动和固化变形等现象。该研究工作对于进一步进行大型复杂复合材料结构固化过程的数值模拟奠定了理论基础和技术手段。参考文献［１］杨乃宾．国外复合材料飞机结构应用现状分析［Ｊ］．航空制造技术，２００２，９：２１—２９［２］ＭＳ，Ａｃ＆Ｓｐｒｉｎｇｅｒ，ＧＳ．Ｃｕｒｉｎｇｏｆｅｐｏｘｙｍａｔｒｉｘｃｏｍｐｏｓｉｔｅｓ［Ｊ］．ＪｏｕｒｎａｌｏｆＣｏｍｐｏｓ—ｉｔｅＭａｔｅｒｉａｌｓ，１９８３，１７（２）：１３５—１６９［３］Ｂｏｇｅｔｔｉ，ＴＡ＆Ｇｉｌｌｅｓｐｉｅ，ＪＷ．Ｔｗｏ—ｄｉｍｅｎｓｉｏｎａｌｃｕｒｅｓｉｍｕｌａｔｉｏｎｏｆｔｈｉｃｋｔｈｅｒｍｏｓｅｔｔｉｎｇｃｏｍｐｏｓｉｔｅ［Ｊ］．Ｊｏｕｒｎａｌ［４］Ｎａｊｉ，ＭｅｆｆｅｃｔｏｎｏｆＣｏｍｐｏｓｉｔｅＭａｔｅｒｉａｌｓ，１９９１，２５（３）：２３９～２７３Ｉ＆ＨｏａＳＶ．Ｃｕｒｉｎｇｏｆｔｈｉｃｋａｎｓｌｅ—ｂｅｎｄｔｈｅｒｍｏｓｅｔｃｏｍｐｏｓｉｔｅｐａｒｔ：ｃｕｒｉｎｇｃｙｃｌｅｔｈｉｃｋｎｅｓｓｖａｒｉａｔｉｏｎａｎｄｆｉｂｅｒｖｏｌｕｍｅｆｒａｃｔｉｏｎ［Ｊ］．Ｊｏｕｒｎａｌｒｅｉｎｆｏｒｃｅｄｐｌａｓｔｉｃｓａｎｄｃｏｍｐｏｓ－ｉｔｅｓ，１９９９，１８（８）：７０２—７２２［５］Ｈｕｂｅｒｔ，Ｐ，Ｖａｚｉｒｉ，Ｒ＆Ｐｏｕｒｓａｒｔｉｐ，Ａ．Ａｕｌａｔｉｏｎｏｆｃｏｍｐｌｅｘｓｈａｐｅｔｗｏ－ｄｉｍｅｎｓｉｏｎａｌｆｌｏｗｍｏｄｅｌｆｏｒｔｈｅｐｒｏｃｅｓｓｓｉｍ－ｆｏｒＮｕｍｅｒｉｃａｌＭｅｔｈｏｄｓｉｎｃｏｍｐｏｓｉｔｅｌａｍｉｎａｔｅｓ［Ｊ］．ＩｎｔｅｒｎａｔｉｏｎａｌＪｏｕｒｎａｌＥｎｇｉｎｅｅｒｉｎｇ，１９９９，４４（１）：ｌ～２６［６］Ｈｕｂｅｒｔ，Ｐ＆Ｐｏｕｒｓａｒｔｉｐ，Ａ．Ａｓｐｅｃｔ１０４ｏｆｔｈｅｃｏｍｐａｃｔｉｏｎｏｆｃｏｍｐｏｓｉｔｅａｎｇｌｅｌａｍｉｎａｔｅｓ：ａｎ热固性复合材料层合板固化过程有限元数值模拟ｅｘｐｅｒｉｍｅｎｔａｌｉｎｖｅｓｔｉｇａｔｉｏｎ［Ｊ］，．ＪｏｕｒｎａｌｏｆＣｏｍｐｏｓｉｔｅＭａｔｅｒｉａｌｓ，２００１，２５（３）：２３９—２７３［７］Ｒａｄｆｏｒｄ，ＤＷ．Ｖｏｌｕｍｅｆｒａｃｔｉｏｎｇｒａｄｉｅｎｔｉｎｄｕｃｅｄｗａｒｐａｇｅｉｎｃｕｒｖｅｄｃｏｍｐｏｓｉｔｅｐｌａｔｅｓ［Ｊ］［８］Ｇｕｔｏｗｓｋｉ，ＴＧＡｒｅｓｉｎｆｌｏｗ／ｆｉｂｅｒｄｅｆｏｒｍａｔｉｏｎｍｏｄｅｌｆｏｒ．Ｃｏｍｐｏｓｉｔｅｅｎｇｉｎｅｅｒｉｎｇ，１９９５，５（７）：９２３～９３４ｃｏｍｐｏｓｉｔｅｓ［Ｊ］．ＳＡＭＰＥｑｕａｒ－ｔｅｒｌｙ，１９８５，１６（４）：５８～“［９］Ｒａｄｆｏｒｄ，ＤＷ＆Ｄｉｅｆｃｎｄｏｒｆ，ＲＪ．Ｓｈａｐｅｈａｔｅｉｎｓｔａｂｉｌｉｔｉｅｓｉｎｃｏｍｐｏｓｉｔｅｓｒｅｓｕｌｔｉｎｇｆｒｏｍｌａｍｉ．ａｎｉｓｏｔｒｏｐｙ［Ｊ］．ＪｏｕｒｎａｌｏｆＲｅｉｎｆｏｒｃｅｄＰｌａｓｔｉｃｓａｎｄＣｏｍｐｏｓｉｔｅｓ，１９９３，１２（１）：５８～７５ＦＤⅡＴＥＥＬＥＭ哐ＮＴＮＵＭ凰ＲＩＣＡＬＳＭＩⅡ．ＡＴＩＯＮＯＦＴＨＥＩ泓ｏＳＥＴＴＤ叮ＧＣＯＭ口ＯＳＩＴＥＬＡⅣⅡＮＡＴＥＣＵＲＩＮＧＰＲＯＣＥＳＳＺｈａｎｇＪｉｋｕｉ，ＣｈｅｎｇＸｉａｏｑｕａｎ，ＬｉＺｈｅｎｇｎｅｎｇ，ＧｕａｎＺｈｉｄｏｎｇ（ＳｃｈｏｏｌｏｆＡｅｒｏｎａｕｔｉｃａｌＳｃｉｅｎｃｅａｎｄＥｎｇｉｎｅｅｒｉｎｇ，ＢｅｉｊｉｎｇＵｎｉｖｅｒｓｉｔｙｏｆＡｅｒｏｎａｕｔｉｃｓ＆Ａｓｔｒｏｎａｕｔｉｃｓ，Ｂｅｉｊｉｎｇ１０００８３，Ｃｈｉｎａ）Ａｂｓｔｒａｃｔ：Ｂａｓｅｄｏｎａｏｎｔｈｅｃｕｒｉｎｇｐｒｏｃｅｓｓｃｏｎｔｒｏｌｅｑｕａｔｉｏｎｓ，ｎｕｍｅｒｉｃａｌｓｉｍｕｌａｔｉｏｎｗａｓｃｏｎｄｕｃｔｅｄｃｕｒｅｓｅｒｉｅｓｏｆｃｏｍｐｌｉｃａｔｅｄｐｈｙｓｉｃａｌａｎｄｃｈｅｍｉｃａｌｖａｒｉｅｔｙｐｒｏｃｅｓｓｄｕｒｉｎｇｔｈｅｃｙｃｌｅｏｆｆｉｂｒｅｒｅｉｎ—ｆｏｒｃｅｄｔｈｅｒｍｏｓｅｔｔｉｎｇｒｅｓｉｎｃｏｍｐｏｓｉｔｅｌａｍｉｎａｔｅｓ，ｓｕｃｈａｓａｎｄｔｈｅｈｅａｔｅｘｃｈａｎｇｅ，ｃｕｒｉｎｇｈｅａｔｒｅｌｅａｓｅ，ｍａｔｒｉｘｏｎｆｌｏｗｉｎｇ，ｃｕｒｅｃｈｅｍｉｃａｌｓｈｒｉｎｋａｇｅ，ｃｕｒｅｄｅｆｏｒｍａｔｉｏｎ，ｅｔｅ．ＣｏｍｐｕｔｅｒｃｏｄｅＷａｓｃｏｍｐｉｌｅｄｏｆｔｈｅｂａｓｉｓｔｏＡＢＡＱＵＳｓｏｆｔｗａｒｅ．Ｔｈｅｆｉｎｉｔｅｅｌｅｍｅｎｔｍｅｔｈｏｄｅｆｆｅｃｔｉｖｅｗｉｔｈａｃｏｍｐｕｔｅｒｃｏｄｅｃｏｍｐｉｌｅｄｗｅｒｅｖｅｒｉｆｉｅｄｂｅｆｅａｓｉｂｉｌｉｔｙａｎｄｃｅｒｔａｉｎｅｘａｍｐｌｅ．Ｋｅｙｗｏｒｄｓ：Ｌａｍｉｎａｔｅ；Ｔｅｍｐｅｒａｔｕｒｅ；Ｃｕｒｅｃｈｅｍｉｃａｌｓｈｒｉｎｋａｇｅ；Ｄｅｆｏｒｍａｔｉｏｎ；ＦＥＡ１０５热固性复合材料层合板固化过程有限元数值模拟作者单位：张纪奎，程小全，郦正能，关志东北京航空航天大学航空科学与工程学院,北京,100083相似文献(4条1.期刊论文张纪奎.郦正能.关志东.程小全.ZHANGJikui.LIZhengneng.GUANZhidong.CHENGXiaoquan固化度与固化收缩对非对称复合材料层合板固化变形的影响-复合材料学报2007,24(2采用三维有限元方法研究复合材料非对称层合板在热载荷和固化收缩载荷下的固化变形情况,建立了材料力学特性、固化体积收缩量和温度与固化度之间的函数关系,考察了层合板变形曲率与温度和固化度之间的关系.数值计算结果表明:非对称层合板变形曲率与固化终止时固化度有密切关系;固化变形主要发生在降温阶段;固化收缩对层合板变形曲率影响很小,主要发生在第二个保温平台的前半段.2.期刊论文张振国.张娜.张秀丽.张恒.ZHANGZhen-Guo.ZHANGNa.ZHANGXiu-Li.ZHANGHeng复合材料薄层板常温固化收缩及翘曲变形-河南科技大学学报（自然科学版）2021,30(6研究了常温下树脂基复合材料层合板固化时的收缩以及翘曲变形,结果表明:薄层板各层不同的树脂含量是翘曲变形的主要因素;对称的铺层结构能有效的降低复合材料薄层板收缩与变形.采用有限元加以分析,有限元计算结果与实验结果相近,可以利用有限元的分析计算近似地预估变形.3.期刊论文戴棣.乔新.DaiDi.QiaoXin复合材料层合板的非同步固化翘曲变形分析-南京航空航天大学学报2000,32(1复合材料层合板在实际固化过程中因板内温度梯度的存在,各铺层内树脂的固化反应并非同步进行.这种非同步过程导致层合板内固化后残余应力的分布不同于经典分析方法的计算结果,使按对称铺层设计的层合板亦存在固化后翘曲变形.针对实际固化过程的非同步特征,文中给出一种基于弹性子叠层块的非同步固化过程描述模型,用此模型对非同步固化过程引起的对称层合板翘曲变形进行计算机模拟,分析树脂固化收缩率,固化过程中板内温度梯度,以及层合板铺叠方式对非同步固化翘曲变形的影响作用,为复合材料层合构件的设计和制造工作提供帮助.4.期刊论文张纪奎.郦正能.关志东.程小全.王军.ZHANGJikui.LIZhengneng.GUANZhidong.CHENGXiaoquan.WANGJun热固性复合材料固化过程三维有限元模拟和变形预测-复合材料学报2021,26(1分析了复合材料热固化过程中各种复杂的物理化学变化之间的相互影响,在此基础上建立了复合材料固化过程数值模拟和固化变形预测的三维有限元分析模型.采用整体-子模块方法将固化过程分为热-化学、流动-压实和应力-变形三个相对独立的子模块.热-化学模块的控制方程基于Fourier热传导方程和树脂固化动力学方程建立,解决了温度和固化度之间的强耦合问题.流动-压实模块的控制方程基于Darcy定律和有效应力原理建立,反映了树脂流动和纤维网络紧密压实之间的流固耦合关系.应力-变形模块建立了考虑热载荷和固化收缩载荷时复合材料层合板的有限元方程.各模块之间的相互作用通过它们之间的数据交换来实现,以树脂在固化过程中的凝胶点和玻璃化转化点为判断依据确定是否运行各模块及其子程序.典型结构的计算结果与实验对比验证了本文三维有限元模型的有效性.本文链接：://d.g.wanfangdata/Conference_6232049.aspx授权使用：洛阳工学院（河南科技大学）(wflskd，授权号：37a47ad2-3b07-4189-8947-9e4d010ced2f下载时间：2021年12月15日复合材料学报第22卷第5期10月2005年ActaMateriaeCompositaeSinicaVol.22No.5October2005文章编号:1000-3851(200505-0155-07收稿日期:2004-11-09;收修改稿日期:2005-01-07基金项目:航空科学基金资助项目(04B51045通讯作者:卢子兴,教授,主要研究方向为复合材料力学、泡沫材料力学E-mail:luzixing@buaa.edu三维四向编织复合材料力学性能的有限元分析杨振宇1,卢子兴1,*,刘振国1,李仲平2(1.北京航空航天大学固体力学研究所,北京100083;2.航天材料工艺研究所,北京100076摘要:在已有研究的基础上,提出了一个新的三维编织复合材料单元胞体模型,该模型正确地反映了纤维束的交织方式,十分接近三维编织复合材料的真实结构,可用于三维四向编织复合材料有效模量的有限元数值预报,并合理确定复合材料内部全场应力分布。采用有限元软件对该模型进行了力学分析,得到了相关等效弹性性能参数。结果表明:有限元计算得到的三维编织复合材料的等效弹性性能与实验结果和理论预测值都吻合较好,从而验证了该模型的有效性。此外,基于新的单元胞体模型还确定了三维四向编织复合材料的应力场,为进一步的强度计算奠定了基础。关键词:三维编织;复合材料;有限元;等效弹性性能;应力场中图分类号:TB330.1文献标识码:AFINITEELEMENTANALYSISOFTHEMECHANICALPROPERTIESOF3-DBRAIDEDCOMPOSITESYANGZhenyu1,LUZixing1,*,LIUZhenguo1,LIZhongping2(1.InstituteofSolidMechanics,BeijingUniversityofAeronauticsandAstronautics,Beijing100083,China;2.InstituteofSpaceMaterialTechnology,Beijing100076,ChinaAbstract:Onthebasisofavailablestudies,anewmodelofrepresentativevolumeelementsof3-Dbraidedcom-positesisproposed,whichtrulyreflectsthebraidedmannerandcoincideswiththeactualconfigurationofthebrai-dedcomposites.Itcanbeusedtopredicttheeffectivemodulusof3-DbraidedcompositeswithFiniteElementMeth-od(FEManddetermineitsstressfieldinreason.Thefiniteelementsoftwarewasadoptedtostudythemechanicalpropertiesofthismodel,andtheeffectivepropertiesofelasticitywereobtained.Theresultsshowthattheonesoffiniteelementcalculationfitwellwiththeexperimentalresultsandtheoreticalprediction,whichdemonstratesthevalidityofthismodel.Inaddition,basedonthenewmodel,thestressfieldof3-Dbraidedcompositesintensionorshearingisdetermined,whichprovidesthebasisforstrengthcalculationinthefuture.Keywords:threedimensionalbraiding;composites;FEM;effectiveelasticproperties;stressfield随着计算机应用的普及和各种有限元计算方法及理论的发展,用有限元分析来预报复合材料力学性能的工作不断涌现,但相对而言,对编织复合材料的力学计算还不多见,且由于编织复合材料结构的复杂性,所以目前的分析模型一般均包含有较多的脱离实际的假设。编织复合材料有限元分析的工作主要有:Lei等[1]建立了三维桁架结构的有限元胞体模型,其中用梁单元模拟纤维束,而基体简化为杆单元,该模型虽然在一定程度上反映了复合材料中纤维束的交织情况,但对纤维束和基体的受力、变形方式及它们之间相互作用的处理过于简化。Chen等[2]基于变分原理,提出用有限多相单元法来预测三维编织复合材料的有效性能。Tang[3]等从三维编织复合材料的细观结构入手,采用非线性有限元模拟分析了三维编织复合材料的变形。而孙慧玉等[4]则通过均匀化理论和非协调多变量有限元法研究了三维编织复合材料的微观力学性能,提出采用非协调位移单元和混合应力单元来模拟三维编织复合材料的有效弹性力学性能。最近,庞宝君等[5]根据自己所建立的三维多向编织复合材料的单胞几何模型,建立了有限元计算模型,他们将纤维束和基体分别视为编织复合材料的组分材料,在单胞内采用规则的8节点体单元对其整个单胞结构进行了剖分,所得单元分成三类:一类为基体单元;另一类为纤维束单元;第三类为混合单元(单元内既有纤维束又有基体材料。此外,文献[6,7]中本文作者也对三维四向编织复合材料的参数化建模技术进行了研究,提出了一种“米”字型体胞的有限元计算模型,数值预报了三维四向编织复合材料的有效弹性模量。然而,上述体胞有限元计算模型与编织物的真实细观结构还存在较大的差别,虽然能给出较好的模量预测结果,但给出的全场应力分布却不尽合理,难以用于材料的强度估算。最近,文献[8]中本文作者又提出了一个新的三维编织复合材料的单元胞体模型,该模型正确反映了纤维束的交织方式,十分接近三维编织复合材料的真实结构,可用于三维四向编织复合材料有效模量的有限元数值预报,并合理确定复合材料内部的全场应力分布。1单元胞体的有限元建模三维编织复合材料具有很复杂的细观结构,且计算机的容量仍然是有限的,故不可能以整体结构为对象进行有限元的分析,而只能选取具有代表性的体积单元(RepresentativeVolumeElement作为研究的对象。文献[8]中,讨论了三维四向编织复合材料的单元胞体模型,所选取的单元为具有周期性分布的最小结构单元。通过胞体尺寸与编织参数之间的几何关系,可以在有限元软件中,按照点→线→面→体的顺序进行有限元的建模。纤维束的几何参数(如图1所示与胞体几何尺寸的关系为α=2arcsin11+cos2θ;c=a41+cos2θ;d=a2cosθ(1图1纤维束截面图Fig.1Sectionofthebundle其中:α为纤维束的截面角;c、d为纤维束的截面尺寸;a为单元胞体的几何尺寸;θ为编织角。根据上述结果可以计算得到有限元模型的每一个关键点的坐标(该模型一共具有96个关键点。由这些参数坐标可以发现,每一个关键点的坐标只与胞体的尺寸相关,因此,十分有利于对此模型进行有限元建模,通过Ansys中的APDL(AnsysParame-terDesignLanguages程序设计语言,可完成整个有限元的建模。2有限元计算2.1材料属性设定与有限元网格划分在新的单元胞体模型中,可将三维编织复合材料看作由两种材料组成,一种为树脂基体,另一种为单向纤维树脂复合材料柱体。这里取纤维束中纤维单丝的数目为9000根,纤维单丝的直径为12μm。设基体为各向同性线弹性材料,而按照单向复合材料的理论,纤维复合柱体为横观各向同性的,则其等效的弹性性能参数计算公式可参见文献[6]。计算中所用的材料参数如表1所示。表1所用碳纤维和树脂基体的弹性性能[6]Table1TheelasticpropertiesofcarbonfiberandepoxyresinmatrixE1/GPaE2/GPaG12/GPaG23/GPaμ12μ23T300carbonfiber22013.89.05TDE-85epoxy采用有限元软件对单元胞体模型划分有限元网格后得到图2所示的结果。图2中单元胞体的参数为θ=35°,Vf=55%,一共分出了32773个单元,6574个节点。并且,在单元胞体表面相对面上的节点是一一对应的。图2单元胞体的有限元网格划分Fig.2Thefiniteelementmeshofunitcell156复合材料学报2.2边界条件的确定和弹性性能的计算根据所取单元胞体模型自身在空间上的周期性分布特点,假设编织体足够大,则当其受到外载荷作用时仍可认为胞体保持空间上的周期性。在下面的计算中针对不同的受载情况,对单元胞体边界的处理也有所不同,目的是为了保证胞体在各种载荷下都能保持空间上的周期性和连续性。为了计算三维四向编织复合材料在纵向上的弹性模量,通过对单元胞体(UA×UA×UB(如图3(a的实线部分施加位移载荷,由有限元法计算得到平均应力和应变,最后得到等效弹性性能。由于新的单元胞体模型不满足结构对称性条件,所以经典的平面假设不再适用,边界条件只能按胞体的周期性来确定。当在Z方向施加位移载荷时,根据周期性边界条件的要求,单元胞体模型中X=0与X=UA的面在X方向上应具有相同的变形模式,而在另两个自由度上具有相同的位移;同样,Y=0与Y=UA的面在Y方向上也应具有相同的变形模式,而在其它两个自由度上保持相同的位移;Z=0与Z=UB的面应具有相同的变形趋势,这样才能保证单元胞体在变形以后仍然保持在空间上的周期性。同时,在计算中对胞体施加适当的约束条件,以保证编织复合材料变形后在空间上的连续性,并使单元胞体面的翘曲变形不受限制。图3单元胞体的剪切变形Fig.3Thedeformationoftheunitcellundershearload2.2.1计算弹性模量Ez和泊松比μzx假设δ为已知的位移载荷,通过有限元软件可计算得到面X=UA相对于面X=0的位移u,则有:–εx=uUA,–εz=δUB(2因此有:μzx=-–εx–εz=-u×UBδ×UA(3因为编织复合材料具有横观各向同性的性质,所以有:μxzEx=μzxEz(4此外,通过有限元软件还可以计算得到Z=UB的单元胞体表面上所有节点在Z方向的节点力之和∑Fz,这样可以得到单元胞体在Z方向的平均应力为:–σz=∑FzUA×UA(5故Z方向的弹性模量为:Ez=–σz–εz=∑FzUBUA2δ(6当在X方向施加位移载荷时,对边界条件进行类似的处理,并施加适当的约束条件,用同样的方法也可以计算得到Ex=Ey和μxy。2.2.2计算剪切模量Gxz对于图3(a所示单元胞体模型,当受到剪切作用时,单元胞体除了产生纯剪切变形外,原单元胞体的各个面还可能产生翘曲变形。各面的翘曲使得Gxz与各面为平面假定(经典理论时的结果不同。考虑到整体复合材料在剪切变形过程中,各单元胞体仍会保持空间的周期性与连续性,而且相邻单元胞体在其分界面上应保持位移连续,因此对每个单元胞体而言(如图3(b所示,B1点与B2点应该保持位移连续,A1点与A2点由于剪切的作用,它们的位移在X、Y方向上连续,在Z方向上满足UZA2-UZA1=δz(7在X=0面上的四个角点处施加固定铰约束,这样使图3(a中相对的两个面都产生相同的变形模式,其中虚线表示面可能出现的变形趋势。δz是对单元胞体施加的剪切变形(位移载荷,通过有限元计算可以得到在X=UA平面上所有沿Z方向的节点力的合力∑Fz。由图3(b可知,宏观的剪切应变(平均应变及相应的平均应力为:γxz=δzUA,–τxz=∑FzUA×UB(8因此有:Gxz=–τxzγxz=∑FzUBδz(9若采用相同的方法处理边界条件,并施加剪切位移载荷,进行有限元法计算则也可确定Gxy。3计算结果与讨论根据表1提供的数据,经计算得到了三维四步矩形编织复合材料的等效弹性性能。图4给出了三维四向编织复合材料的等效弹性性能随编织角及纤157杨振宇,等:三维四向编织复合材料力学性能的有限元分析图4等效弹性性能与编织角的变化关系Fig.4Effectiveelasticpropertiesvariationwithbraidedangle维体积含量的变化曲线。由图4(a、图4(b可见,给定纤维体积比,纵向杨氏模量Ez随编织角的增大而减小,而横向杨氏模量Ex随编织角的增大而增大。纵向杨氏模量Ez则在小编织角情况下对纤维体积比的变化比较敏感,在大编织角时随纤维体积比的改变不大,而横向杨氏模量Ex对纤维体积比始终保持相当的敏感性。从图4(c、图4(d可以看出,给定纤维体积比,随编织角的增大,横向158复合材料学报剪切模量Gxy增大,而纵向剪切模量Gxz则先增大后减小,在编织角约等于45°时具有最大值。此外,在大编织角的情况下,横向剪切模量Gxy受纤维体积比的影响较大,而纵向剪切模量Gxz则在编织角约为40°~50°时受纤维体积比的影响最大。在图4(e、图4(f的计算结果中泊松比的变化比较复杂。给定纤维体积比,随编织角的增大,泊松比μxy先减小,然后再增大,并且存在一个特定编织角(约45°,在小于该编织角的条件下,μxy随纤维体积比的增大而增大,而在大于该编织角的条件下,μxy随纤维体积比的增加变化不明显。泊松比μzx随编织角的增大先增加后减小,并且对纤维体积比的变化不敏感,存在一特定的编织角(约30°,在小于该编织角的条件下,μzx随编织角的增大而增加,而在大于该编织角的条件下,μzx随编织角的增大而减小。总体上,泊松比μxy和μzx主要由编织角控制,而纤维体积比对其影响较小;另外,若给定编织角,则随纤维体积含量的增加,弹性模量Ex和Ez、剪切模量Gxy和Gxz、泊松比μxy基本上表现为线性规律的增大[11],而泊松比μzx则具有一个特殊的非线性的变化过程,在编织角等于30°时,它随纤维体积含量的增加而增大,并且趋于平缓;但给定较小编织角时,泊松比μzx表现出先增大后减小的变化趋势[11]。此外,比较图4中的纵向弹性模量和横向模量可以看出,当编织角大于某一特定角度时,横向模量Ex可以大于纵向模量Ez。根据单元胞体的几何模型可知,当编织角趋近于0°时,模型应退化为单向复合材料,故三维四向编织复合材料的纵向弹性模量的上限应该趋近于相同纤维体积含量的单向复合材料的相应模量,计算的结果也表明,三维四向编织复合材料的纵向弹性模量不会超过单向复合材料的相应值;同样地,当编织角趋近于90°时,纤维沿横向的角度趋近于45°,模型应接近[45°/-45°]T的层合板,故横向弹性模量趋近于[45°/-45°]T层合板的纵向弹性模量,并且不会超过该值。同时,不难看到,由新的有限元模型计算得到的三维编织复合材料的等效弹性性能随编织角和纤维体积含量的变化趋势与文献[10]中的理论计算结果十分相似,且在数值上也比较接近,这在某种程度上证明了该有限元计算模型的正确性。需要指出,经过不同密度的网格划分,对同一性能进行计算时,会发现计算结果会随网格的加密不断地变化,呈振荡的趋势,但随网格的加密,计算结果的变化会越来越小。因此,网格的划分会对最后的计算结果有一定的影响。需指出,有限元计算中,同时也得到了单元胞体在拉伸载荷和剪切载荷作用下的应力云图[11]。例如,图5给出了Z向拉伸载荷下单元胞体的全场应力云图和胞体中主要纤维束的应力云图。可见,纤维束是编织复合材料中的主要承载部分,从应力的数值上看,纤维束中应力是基体中应力的数十倍。这主要是由于纤维的纵向弹性模量(E1=220GPa与基体的弹性模量(E=4.5GPa相差比较悬殊,在相同的应变下,纤维束将产生远比基体大的应力。而处于纤维束交叉处的基体也表现出了一定的应力集中现象,越是靠近交叉部位,基体的应力就越大,这主要是由于在受到拉伸载荷时纤维束会发生微小的转动,从而对基体形成了挤压的作用图5单元胞体的应力云图Fig.5Stressnephogramoftheunitcell159杨振宇,等:三维四向编织复合材料力学性能的有限元分析160复合材料学报验值,介于实验值与理论预测之间。由于有限元计算采用了周期性的边界条件,忽略了三维编织复合材料边界的影响,此外,也没有考虑复合材料内部可能存在的缺陷和纤维束弯曲带来的影响,故有限元模型计算的模量一般大于实验值。力。在剪切载荷作用下,纤维束的受力情况更加复杂,纤维束之间有明显的错动,纤维束交叉处的应力集中现象也更加明显,而且纤维束上存在较大的应力梯度[11](限于篇幅,这里略去应力云图。为了进一步验证有限元计算模型的有效性,这里将有限元计算结果与文献[9]的实验和文献[10]的理论预测结果进行了比较。共选取了5种情况:A(Vf=45%,θ=21°、B(Vf=54%,θ=41°、C(Vf=54%,θ=42°、D(Vf=58%,θ=42°和E(Vf=45%,θ=48°,其结果示于图6中。可见,纵向弹性模量Ez的有限元计算结果比实验数据总体上是偏高的,在40°以上编织角的复合材料预测中,有限元的计算值大于实验结果和理论预测结果,但是在小编织角(θ=21°的时候非常接近于实验数据,并优于理论预测的结果。而对泊松比的预测,在大编织角(θ=48°试件的比较中,有限元法的计算结果稍微小于实验值和理论预测值;而在其他几个试件的比较中,有限元法预测结果略大于实4结论综上所述,有限元计算得到的三维编织复合材料的等效弹性性能与实验结果和理论预测都吻合较好,从而证明了新的三维编织复合材料单元胞体模型的有效性。同时,由于新的单元胞体模型正确地反映了纤维束的交织方式,真实地模拟了三维编织复合材料的细观结构,所以基于新的有限元模型得到的三维编织复合材料内部的应力分布,对下一步进行三维编织复合材料的强度预测具有重要的理论意义和实用价值。此外,从有限元分析得到的应力云图不难看出,编织节点处的纤维束,是单向拉伸载荷下增强相中强度的最薄弱点,在较大拉伸载荷作用下,该处有先被局部拉裂的可能性。同样,位于纤维束交织节点处的基体区域存在着明显的应力集中,该处的拉伸应力、剪切应力对应力集中起着主导作用,这是编织复合材料强度最为薄弱、最易损伤的部位。另外,在载荷作用下,由于纤维束在编织复合材料中处于倾斜状态,纤维束与基体、纤维束与纤维束之间存在有相互的作用,所以导致了编织复合材料内部各纤维束总是处于复杂的应力状态中。参考文献:[1]LeiC,CaiYJ,KoFK.Finiteelementanalysisof3Dbraidedcomposites[J].AdvancesinEngineeringSoftware,1992,14:187-194.[2]ChenL,TaoXM,ChoyCL.Mechanicalanalysisof3Dbraidedcompositesbythefinitemultiphaseelementmethod[J].CompositesScienceandTechnology,1999,59(16:2383-2391.[3]TangZX,PostleR.Mechanicsofthreedimensionalbraidedstructuresforcompositematerials———PartⅢ:Nonlinearfiniteelementdeformationanalysis[J].Composite