【刚度计算】 Weber算法齿轮啮合刚度计算【含Matlab源码 3459期】.zip

%注意求刚度的时候用k=F/delta，而不是用k=F/(B*delta)。 %数据来源：多篇参考文献 %采用方法：Webster方法（材料力学方法） %两个都是外齿轮啮合 %原理参考《机械原理》P184和P185 clear all; clc; format long; F1=170000; %齿轮作用力（摆设，后面求刚度的时候会约掉） F2=170000; pi=3.14; %%%%%%%%%%%%%%齿轮副基本参数及尺寸%%%%%%%%%%%%%%%%% z1=25;z2=39; %《直齿轮轮齿变形计算的数值积分法》颜海燕 B1=20/1000;B2=20/1000; mm=3.0/1000; % z1=37;z2=96; %《齿轮时变啮合刚度改进算法及刚度激励研究》李亚鹏 % mm=11/1000; % z1=45;z2=90; % mm=3/1000; % B1=mm*8;B2=mm*8; ha=1;c=0.25; alpha=20*pi/180; invalpha=0.01490438; d1=mm*z1;r1=d1/2;%分度圆半径 d2=mm*z2;r2=d2/2; % B1=d1*0.4; % B2=d2*0.4; da1=(z1+2*ha)*mm;ra1=da1/2;%齿顶圆半径 da2=(z2+2*ha)*mm;ra2=da2/2; df1=(z1-2*ha-2*c)*mm;rf1=df1/2;%齿根圆半径 df2=(z2-2*ha-2*c)*mm;rf2=df2/2; db1=d1*cos(alpha);rb1=db1/2;%基圆半径 db2=d2*cos(alpha);rb2=db2/2; s1=pi*mm/2;s2=pi*mm/2; %齿厚 alphaa1=acos(rb1/ra1); %齿顶圆压力角 alphaa2=acos(rb2/ra2); a=mm*(z1+z2)/2; %标准中心距 ap=a;%实际中心距 alphap=acos(a*cos(alpha)/ap); invalphap=tan(alphap)-alphap; dp1=d1*cos(alpha)/cos(alphap);%节圆直径1 rp1=dp1/2; sp1=s1*rp1/r1-2*rp1*(invalphap-invalpha);%节圆齿厚1 Hp1=sp1;%节圆齿厚1 dp2=d2*cos(alpha)/cos(alphap);%节圆直径2 rp2=dp2/2; sp2=s2*rp2/r2-2*rp2*(invalphap-invalpha);%节圆齿厚2 Hp2=sp2;%节圆齿厚2 p=pi*mm; Pb=p*cos(alpha); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%齿轮的材料参数%%%%%%%%%%%%%% E1=2.09e11; %单位为Pa E2=2.09e11; miu1=0.26; miu2=0.26; if B1/Hp1>5 Ee1=E1/(1-miu1^2); else Ee1=E1;%轮齿材料的等效弹性模量 end if B2/Hp2>5 Ee2=E2/(1-miu2^2); else Ee2=E2;%轮齿材料的等效弹性模量 end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%传动基本参数区%%%%%%%%%%%%%%%%% PB1=rb1*(tan(alphaa1)-tan(alphap)); PB2=rb2*(tan(alphaa2)-tan(alphap)); B1B2=PB1+PB2; B2C=B1B2-Pb; N1B1=sqrt(ra1^2-rb1^2); N1B2=N1B1-B1B2; N2B2=sqrt(ra2^2-rb2^2); N2B1=N2B2-B1B2; CD=Pb-B2C; N1C=N1B2+B2C; N2C=N2B1+Pb; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% m=9;n=40; N1Bi=zeros(1,m); O1Bi=zeros(1,m); rp1i=zeros(1,m); alphap1i=zeros(1,m); invalphap1i=zeros(1,m); faip1i=zeros(1,m); xp1i=zeros(1,m); yp1i=zeros(1,m); beta1i=zeros(1,m); Lf1i=zeros(1,m); rou1i=zeros(1,m); rK1i=zeros(m,n+1); x1i=zeros(m,n+1); y1i=zeros(m,n+1); Tk1i=zeros(m,n+1); alphaK1i=zeros(m,n+1); invalphaK1i=zeros(m,n+1); A1i=zeros(m,n+1); I1i=zeros(m,n+1); equivalentA1i=zeros(m,n); equivalentI1i=zeros(m,n); N2Bi=zeros(1,m); O2Bi=zeros(1,m); rp2i=zeros(1,m); alphap2i=zeros(1,m); invalphap2i=zeros(1,m); faip2i=zeros(1,m); xp2i=zeros(1,m); yp2i=zeros(1,m); beta2i=zeros(1,m); Lf2i=zeros(1,m); rou2i=zeros(1,m); rK2i=zeros(m,n+1); x2i=zeros(m,n+1); y2i=zeros(m,n+1); Tk2i=zeros(m,n+1); alphaK2i=zeros(m,n+1); invalphaK2i=zeros(m,n+1); A2i=zeros(m,n+1); I2i=zeros(m,n+1); equivalentA2i=zeros(m,n); equivalentI2i=zeros(m,n); N1Bj=zeros(1,m); O1Bj=zeros(1,m); rp1j=zeros(1,m); alphap1j=zeros(1,m); invalphap1j=zeros(1,m); faip1j=zeros(1,m); xp1j=zeros(1,m); yp1j=zeros(1,m); beta1j=zeros(1,m); Lf1j=zeros(1,m); rou1j=zeros(1,m); rK1j=zeros(m,n+1); x1j=zeros(m,n+1); y1j=zeros(m,n+1); Tk1j=zeros(m,n+1); alphaK1j=zeros(m,n+1); invalphaK1j=zeros(m,n+1); A1j=zeros(m,n+1); I1j=zeros(m,n+1); equivalentA1j=zeros(m,n); equivalentI1j=zeros(m,n); N2Bj=zeros(1,m); O2Bj=zeros(1,m); rp2j=zeros(1,m); alphap2j=zeros(1,m); invalphap2j=zeros(1,m); faip2j=zeros(1,m); xp2j=zeros(1,m); yp2j=zeros(1,m); beta2j=zeros(1,m); Lf2j=zeros(1,m); rou2j=zeros(1,m); rK2j=zeros(m,n+1); x2j=zeros(m,n+1); y2j=zeros(m,n+1); Tk2j=zeros(m,n+1); alphaK2j=zeros(m,n+1); invalphaK2j=zeros(m,n+1); A2j=zeros(m,n+1); I2j=zeros(m,n+1); equivalentA2j=zeros(m,n); equivalentI2j=zeros(m,n); %%%%%%%%%%%%%%%%%%%%%%%%%%%%% xM1=rf1; xM2=rf2; yM1=rf1*tan(360/(2*z1)*pi/180); yM2=rf2*tan(360/(2*z2)*pi/180); Hf1=2*yM1; Hf2=2*yM2; %%%%%%%%%%%%%%%%%%%第一个双齿啮合区B2C%%%%%%%%%%%%%%%%%%%%%%% for p=1:m %i啮合点 %第1个外齿轮啮合点 N1Bi(1,p)=N1B2+(p-1)*(B2C/m); O1Bi(1,p)=sqrt(N1Bi(1,p)^2+rb1^2); rp1i(1,p)=O1Bi(1,p); alphap1i(1,p)=acos(rb1/rp1i(1,p));%弧度 invalphap1i(1,p)=tan(alphap1i(1,p))-alphap1i(1,p); faip1i(1,p)=s1/r1-2*(invalphap1i(1,p)-invalpha); xp1i(1,p)=rp1i(1,p)*cos(faip1i(1,p)/2); %草图中的xp yp1i(1,p)=0.5*(s1*rp1i(1,p)/r1-2*rp1i(1,p)*(invalphap1i(1,p)-invalpha)); %草图中的yp beta1i(1,p)=alphap1i(1,p)-faip1i(1,p)/2; %草图中的betap Lf1i(1,p)=xp1i(1,p)-xM1-yp1i(1,p)*tan(beta1i(1,p)); rou1i(1,p)=rb1*tan(alphap1i(1,p)); %每对应一个啮合点齿轮1分成n段 for k=1:n+1 rK1i(p,k)=rf1+(rp1i(1,p)-rf1)*(k-1)/n; %分基圆小于齿根圆和基圆大于齿根圆两种情况讨论 if rb1<rf1 %基圆小于齿根圆时整个齿都视为渐开线齿廓 alphaK1i(p,k)=acos(rb1/rK1i(p,k)); %渐开线上任意点压力角%弧度 invalphaK1i(p,k)=tan(alphaK1i(p,k))-alphaK1i(p,k); x1i(p,k)=rK1i(p,k); %对应草图中的xk y1i(p,k)=0.5*(s1*rK1i(p,k)/r1-2*rK1i(p,k)*(invalphaK1i(p,k)-invalpha)); %对应草图中的yk Tk1i(p,k)=(rp1i(1,p)-rf1)/n; else %基圆大于齿根圆时基圆到齿根圆部分视为矩形，矩形宽为sb alphab1i=acos(rb1/rb1); invalphab1i=tan(alphab1i)-alphab1i; sb1=s1*rb1/r1-2*rb1*(invalphab1i-invalpha); Tk1i(p,k)=(rp1i(1,p)-rf1)/n; x1i(p,k)=rK1i(p,k); if (rK1i(p,k)<rb1) y1i(p,k)=0.5*sb1; else alphaK1i(p,k)=acos(rb1/rK1i(p,k)); %渐开线上任意点压力角%弧度 invalphaK1i(p,k)=tan(alphaK1i(p,k))-alphaK1i(p,k); y1i(p,k)=0.5*(s1*rK1i(p,k)/r1-2*rK1i(p,k).*(invalphaK1i(p,k)-invalpha)); end end A1i(p,k)=2*y1i(p,k)*B1; I1i(p,k)=B1*(2*y1i(p,k)).^3/12; end for k=1:n equivalentA1i(p,k)=2*(A1i(p,k)*A1i(p,k+1)/(A1i(p,k)+A1i(p,k+1)));%小段k的当量截面积 equivalentI1i(p,k)=2*(I1i(p,k)*I1i(p,k+1)/(I1i(p,k)+I1i(p,k+1)));%小段k的当量惯性矩 end %第2个外齿轮啮合点 N2Bi(1,p)=N2B2-(p-1)*(B2C/m); O2Bi(1,p)=sqrt(N2Bi(1,p)^2+rb2^2); rp2i(1,p)=O2Bi(1,p); alphap2i(1,p)=acos(rb2/rp2i(1,p));%弧度 invalphap2i(1,p)=tan(alphap2i(1,p))-alphap2i(1,p); faip2i(1,p)=s2/r2-2*(invalphap2i(1,p)-invalpha); xp2i(1,p)=rp2i(1,p)*cos(faip2i(1,p)/2); %草图中的xp yp2i(1,p)=0.5*(s2*rp2i(1,p)/r2-2*rp2i(1,p)*(invalphap2i(1,p)-invalpha)); %草图中的yp beta2i(1,p)=alphap2i(1,p)-faip2i(1,p)/2; %草图中的betap Lf2i(1,p)=xp2i(1,p)-xM2-yp2i(1,p)*tan(beta2i(1,p)); rou2i(1,p)=rb2*tan(alphap2i(1,p)); %每对应一个啮合点齿轮1分成n段 for k=1:n+1 rK2i(p,k)=rf2+(rp2i(1,p)-rf2)*(k-1)/n; %分基圆小于齿根圆和基圆大于齿根圆两种情况讨论 if rb2<rf2 %基圆小于齿根圆时整个齿都视为渐开线齿廓 alphaK2i(p,k)=acos(rb2/rK2i(p,k)); %渐开线上任意点压力角%弧度 invalphaK2i(p,k)=tan(alphaK2i(p,k))-alphaK2i(p,k); x2i(p,k)=rK2i(p,k); %对应草图中的xk y2i(p,k)=0.5*(s2*rK2i(p,k)/r2-2*rK2i(p,k)*(invalphaK2i(p,k)-invalpha)); %对应草图中的yk Tk2i(p,k)=(rp2i(1,p)-rf2)/n; else %基圆大于齿根圆时基圆到齿根圆部分视为矩形，矩形宽为sb alphab2i=acos(rb2/rb2); invalphab2i=tan(alphab2i)-alphab2i; sb2=s2*rb2/r2-2*rb2*(invalphab2i-invalpha); Tk2i(p,k)=(rp2i(1,p)-rf2)/n; x2i(p,k)=rK2i(p,k); if (rK2i(p,k)<rb2) y2i(p,k)=0.5*sb2; else alphaK2i(p,k)=acos(rb2/rK2i(p,k)); %渐开线上任意点压力角%弧度 invalphaK2i(p,k)=tan(alphaK2i(p,k))-alphaK2i(p,k); y2i(p,k)=0.5*(s2*rK2i(p,k)/r2-2*rK2i(p,k).*(invalphaK2i(p,k)-invalpha)); e