航空发动机强度计算

航空发动机强度计算Matlab编程作业姓名：王英杰学号：143403040512班级：1434200101专业：飞行器动力工程学院：创新学院2017年12月27日一数值积分法涡轮叶片强度计算某涡轮的转数为4700r/min，共有68片转子叶片，气流参数沿叶高均匀分布，平均半径处叶栅进出口的气流参数，各个截面的重心位置，截面面积A和ξ轴与y轴的夹角а均列于表1-3中，叶片材料的GH33的密度ρ为8.2x103kg/𝑚3，求叶片各个截面上的离心拉伸应力，气体力弯矩，离心力弯矩与合成弯矩。表1.1叶片各个截面数据012345X（cm）0.530.410.410.400.240.12Y（cm）0.410.380.300.190.110.02Z（cm）62.859.156.053.049.445.8A（cm2）1.802.323.124.105.487.05а31º40′27º49′25º19′22º5′16º57′12º43′C1amC1umρ1mP1mC2amC2umρ2mP2m297m/s410m/s0.894kg/m30.222Mpa313m/s-38m/s0.75kg/m30.178Mpa利用数值积分法进行计算，参考《航空发动机强度计算》书上求涡轮叶片各个截面上离心拉伸应力，气体力弯矩，离心力弯矩与合成弯矩的数值积分法公式，利用MATLAB编程计算结果如下：1．1离心拉伸应力计算MATLAB离心拉伸应力计算计算结果：1.2气体力弯矩计算1.3离心力弯矩计算1.4离心力弯矩和气体力弯矩在x,y轴的合成弯矩计算1.5离心力弯矩和气体力弯矩在主惯性轴合弯矩计算1.6主惯性轴合弯矩计算由上图可以看出适当地设计气离心力弯矩，改变离心力弯矩的大小和方向使之与气体力弯矩方向相反，相互抵消。气体力弯矩和离心力弯矩合成后，合成弯矩在x,y方向上有了明显的降低，降低了叶片的弯曲应力，起到了弯矩补偿的作用。1.6叶片截面主惯性轴（𝜂，𝜉轴）和非主惯性轴（x,y轴）合成弯矩比较。二等厚圆环法涡轮盘强度计算以实心涡轮盘为例，用等厚圆环法计算其应力分布。可以利用。原始数据有：转速：n=12500r/min盘材料密度：ρ=8.0x103kg/𝑚2盘缘径向应力：𝜎𝑎=140𝑀𝑝𝑎泊松比：ν=0.3盘外缘温度:𝑡𝑎=500℃,盘心温度:𝑡0=150℃轮盘厚度ℎ𝑖，弹性模量𝐸𝑖，线膨胀系数𝛼𝑖，以及温度𝑡𝑖沿轮盘半径的分布规律见表1所示。表1：已知轮盘数据序号截面半径R（cm）盘厚H（cm）弹性模量E（Mpa）温度T（℃）线膨胀系数（1/℃）004.861.62x10516516.5x10−615.03.901.62x10516516.5x10−6210.02.971.57x10525017.1x10−6314.02.241.48x10536018.2x10−6415.01.861.40x10540019.0x10−6515.81.601.37x10543019.4x10−6616.61.801.34x10546019.7x10−6717.42.301.30x10550020.3x10−6表2：已经计算的轮盘各平均半径处应力分布值段平均直径（cm）02.57.512.014.515.416.217.017.4段的厚度（cm）4.864.383.432.612.051.731.702.052.30径向应力𝜎𝑟（Mpa）290.3283.3320.5314.8290.0290.2248.3163.2140.0周向应力290.3286.3300.6108.9-167.3-270.3-344.8-428.2-417.1𝜎𝜃（Mpa）利用等厚圆环法进行计算，参考《航空发动机强度计算》书上等厚圆环法的计算公式，利用MATLAB编程进行计算，计算结果如下所示：2.1根据已经给出的圆盘应力画出圆盘沿半径的应力分布图如下所示：2.2通过给出的已知条件进行计算得到的圆盘应力分布如下图所示：2.3将计算得到的圆盘应力分布和已经给出的圆盘应力分布进行比较如下图所示：三传递矩阵法转子临界转速计算简支转子如图所示，总长l=72cm，等截面实心轴直径2.5cm，盘离左支点为24cm，盘质量为m=8.9kg，直径转动惯量𝐼𝑑=319.6𝑘𝑔𝑐𝑚2，极转动惯量𝐼0=639.2𝑘𝑔𝑐𝑚2，轴的弹性模量E为2.06x107𝑁/𝑐𝑚2，两支点均为刚性支撑，试用传递矩阵法计算其临界转速。图3.1单盘转子系统将轴分成9段，每段轴长度均为𝑙𝑖=8cm，各段轴质量按照重心不变的原则离散到各站，每段和站的原始数据如下表：表3.1转子系统原始数据i123456789𝑙𝑖（cm）888888888𝑑𝑖（cm）2.52.52.52.52.52.52.52.52.5𝑚𝑖（kg））0.3080.3089.210.3080.3080.3080.3080.3080𝐼𝑑𝑖（𝑘𝑔𝑐𝑚2）00319.6000000𝐼𝑂𝑖（𝑘𝑔𝑐𝑚2）00639.2000000𝐽𝑖（）1.9171.9171.9171.9171.9171.9171.9171.9171.917𝐽𝑖=𝜋𝑑𝑖464为各截面的惯性矩。3456789021利用传递矩阵法求解该转子系统的临界转速，参考《航空发动机强度计算》一书上传递矩阵法计算公式，用MATLAB编程进行计算，得到的结果如下所示：3.1根据传递矩阵法的原理，画出Δ随Ω（轴的转速）变化的曲线，曲线与Δ=0时的水平轴交点频率就为临界转速。𝜔1cor由上图可以看出当𝜔近似为249rad/s时，Δ曲线第一次与Δ=0轴相交，此时的𝜔值即为转子系统的一阶临界转速𝜔1cor。与教材所给计算临界转速值248.81rad/s基本吻合。由上图可以看出转子在到达115286rpm还未出现第二阶临界转速。在工程实际中往往着重关注转子的第一阶临界转速，工作在第一阶临界转速转子附近的转子振动强烈，危险性最大，而高频振动时转子的振幅却很小，对转子系统的危害较小。3.2也可以采用循环试凑法，设定一个初试的转速值𝜔=250rad/s，代入Δ表达式中进行计算求出Δ值看是否为零，不为零时增加转速值，继续进行上述操作，当Δ值和0接近满足一定的精度时，即可认为此时的转速为临界转速，试凑停止。利用MATLAB进行试凑法计算，选定初始的转速为250rad/s,每次增加1，精度设为0.01，循环计算满足精度条件时，求出一阶临界转速也为249rad/s。附录：一数值积分法涡轮叶片强度计算%AircraftTurbinebladestrengthcalculation%*****************************************************************************%发动机数据n=4700;%r/minw=(4700*2*pi)/60;N=68;%叶片个数RH=8.2e3;%叶片材料GH33密度kg/m3%涡轮叶片各个截面数据：重心坐标，截面面积，坐标轴夹角，进出口速度和压强%x=[0.53,0.41,0.41,0.40,0.24,0.12];%cm%y=[0.41,0.38,0.30,0.19,0.11,0.02];%cm%z=[62.8,59.1,56.0,53.0,49.4,45.8];%cm%A=[1.8,2.32,3.12,4.1,5.48,7.05];;%m2x=[0.0053,0.0041,0.0041,0.0040,0.0024,0.0012];%my=[0.0041,0.0038,0.0030,0.0019,0.0011,0.0002];%mz=[0.628,0.591,0.56,0.53,0.494,0.458];%mA=[0.00018,0.000232,0.000312,0.00041,0.000548,0.000705];%m2afa0=[31.67,27.82,25.32,22.092,16.95,12.72];afa=degtorad(afa0);%角度转弧度C1am=297;C1um=410;RO1m=0.894;P1m=0.222e6;%y叶片进口气动参数C2am=313;C2um=-38;RO2m=0.750;P2m=0.178e6;%y叶片进口气动参数%*****************************************************************************%数值积分法叶片各截面离心拉伸引力计算clcAm=[];zm=[];dz=[];dP=[];xgm=[];fori=1:5Am(i)=(A(i)+A(i+1))/2;zm(i)=(z(i)+z(i+1))/2;dz(i)=z(i)-z(i+1);dP(i)=RH*w^2*Am(i)*zm(i)*dz(i);endxgm(1)=0;%叶片顶部应力为零forj=1:5sdP=0;forjj=1:jsdP=dP(jj)+sdP;endxgm(j+1)=sdP/A(j+1);end%plot(z,xgm);%gridon;%xlabel('叶片高度/m','FontSize',10.5);ylabel('应力/Pa','FontSize',10.5);%legend('离心拉伸应力');%view(90,-90);%*****************************************************************************%数值积分法叶片各截面气体力弯矩计算，气动参数沿叶片径向时均匀分布的px=[];py=[];fori=1:5px(i)=(2*pi*zm(i)/N)*((RO1m*C1am^2-RO2m*C2am^2)+(P1m-P2m))*(dz(i));py(i)=(2*pi*zm(i)/N)*((RO1m*C1am*C1um-RO2m*C2am*C2um))*(dz(i));endMx=[];My=[];Mx(1)=0;My(1)=0;forj=1:5sMx=0;forjj=1:jsMx=sMx+py(jj)*(zm(jj)-z(j+1));endMx(j+1)=sMx;endforj=1:5sMy=0;forjj=1:jsMy=sMy-px(jj)*(zm(jj)-z(j+1));endMy(j+1)=sMy;end%plot(z,Mx);%holdon;gridon;%xlabel('叶片高度/m','FontSize',10.5);ylabel('弯矩/Nm','FontSize',10.5);%plot(z,My);%legend('x方向气体力弯矩','y方向气体力弯矩');%view(90,-90);%*****************************************************************************%数值积分法叶片各截面离心力弯矩计算ym=[];xm=[];%Am(i)=(A(i)+A(i+1))/2;zm(i)=(z(i)+z(i+1))/2;dz(i)=z(i)-z(i+1);fori=1:5ym(i)=(y(i)+y(i+1))/2;xm(i)=(x(i)+x(i+1))/2;endPLy=[];PLz=[];fori=1:5PLy(i)=RH*w^2*dz(i)*Am(i)*ym(i);PLz(i)=RH*w^2*dz(i)*Am(i)*zm(i);endMLx=[];MLy=[];MLx(1)=0;MLy(1)=0;forj=1:5sMLx=0;sMLy=0;forjj=1:jsMLx=sMLx-PLz(jj)*(ym(jj)-y(j+1))+PLy(jj)*(zm(jj)-z(j+1));sMLy=sMLy+PLz(jj)*(xm(jj)-x(j+1));endMLx(j+1)=sMLx;M