传递矩阵-matlab程序.doc
- 文档编号:79574
- 上传时间:2022-10-02
- 格式:DOC
- 页数:23
- 大小:104KB
传递矩阵-matlab程序.doc
《传递矩阵-matlab程序.doc》由会员分享,可在线阅读,更多相关《传递矩阵-matlab程序.doc(23页珍藏版)》请在冰豆网上搜索。
%main_critical.m
%该程序使用Riccati传递距阵法计算转子系统的临界转速及振型
%本函数中均采用国际单位制
%第一步:
设置初始条件(调用函数shaft_parameters)
%初始值设置包括:
轴段数N,搜索次数M
%输入轴段参数:
内径d,外径D,轴段长度l,支撑刚度K,单元质量mm,极转动惯量Jpp[N,M,d,D,l,K,mm,Jpp]=shaft_parameters;
%第二步:
计算单元的5个特征值(调用函数shaft_pra_cal)
%单元的5个特征值:
%m_k:
:
质量
%Jp_k:
极转动惯量
%Jd_k:
直径转动惯量
%EI:
弹性模量与截面对中性轴的惯性矩的乘积
%rr:
剪切影响系数
[m_k,Jp_k,EI,rr]=shaft_pra_cal(N,D,d,l,Jpp,mm);
%第三步:
计算剩余量(调用函数surplus_calculate),并绘制剩余量图
%剩余量:
D1
fori=1:
1:
M
ptx(i)=0;
pty(i)=0;
end
forii=1:
1:
M
wi=ii/1*2+50;
[D1,SS,Sn]=surplus_calculate(N,wi,K,m_k,Jp_k,JD_k,l,EI,rr);
D1;
pty(ii)=D1;
ptx(ii)=w1
end
ylabel(‘剩余量’);
plot(ptx,pty)
xlabel(‘角速度red/s’);
gridon
%第四步:
用二分法求固有频率及振型图
%固有频率:
Critical_speed
wi=50;
fori=1:
1:
4
order=i
[D1,SS,Sn]=surplus_calculate(N,wi,k,m_k,Jp_k,Jd_k,l,EI,rr);
Step=1;
D2=D1;
kkk=1;
whilekkk<5000
ifD2*D1>0
wi=wi+step;
D2=D1;
[D1,SS,Sn]=surplus_calculate(N,wi,K,m_k,Jp_k,Jd_k,l,EI,rr);
end
ifD1*D2<0
wi=wi-step;
step=step/2;
wi=wi+step;
[D1,SS,Sn]=surplus_calculate(N,wi,K,m_k,Jp_k,Jd_k,l,EI,rr);
End
D1;
Wi;
Ifatep<1/2000
Kkk=5000;
end
end
Critical_speed=wi/2/pi*60
figure;
plot_mode(N,l,SS,Sn)
wi=wi+2;
end
%surplus_calculate,.m
%计算剩余量
%
(1)计算传递矩阵
%
(2)计算剩余量
function[D1,SS,Sn]=surplus_calculate(N,wi,K,m_k,Jp_k,Jd_k,l,EI,rr);
%
(1)计算传递矩阵
%===============
%(a)初值设为0
%===============
fori=1:
1:
N+1
forj=1:
1:
2
fork=1:
1:
2
ud11(j,k.i)=0;
ud12(j,k.i)=0;
ud21(j,k.i)=0;
ud22(j,k.i)=0;
end
end
end
fori=1:
1:
N
forj=1:
1:
2
fork=1:
1:
2
us11(j,k.i)=0;
us12(j,k.i)=0;
us21(j,k.i)=0;
us22(j,k.i)=0;
end
end
end
fori=1:
1:
N
forj=1:
1:
2
fork=1:
1:
2
u11(j,k.i)=0;
u12(j,k.i)=0;
u21(j,k.i)=0;
u22(j,k.i)=0;
end
end
end
%============
%(b)计算质点上传递矩阵―――点矩阵的一部分!
%============
fori=1:
1:
N+1
ud11(1,1,i)=1;ud11(1,2,i)=0;ud11(2,1,i)=0;ud11(2,2,i)=1;
ud21(1,1,i)=0;ud21(1,2,i)=0;ud21(2,1,i)=0;ud21(2,2,i)=0;
ud22(1,1,i)=1;ud22(1,2,i)=0;ud22(2,1,i)=0;ud22(2,2,i)=1;
end
%============
%(c)计算质点上传递矩阵―――点矩阵的一部分!
%============
fori=1:
1:
N+1
ud12(1,1,i)=0;
ud12(1,2,i)=(Jp_k(i)-Jd_k(i))*wi^2;%%%考虑陀螺力矩
ud12(2,1,i)=m_k(i)*wi^2-k(i);
ud12(2,2,i)=0;
end
%============
%(d)以下计算的是无质量梁上的传递矩阵―――场矩阵
%计算的锥轴的us是不对的,是随便令的,在后面计算剩余量时,zhui中会把错误的覆盖掉
%============
fori=1:
1:
N
us11(1,1,i)=1;us11(1,2,i)=1(i);us11(2,1,i)=0;us11(2,2,i)=1;
us12(1,1,i)=0;us12(1,2,i)=0;us12(2,1,i)=0;us12(2,2,i)=0;
us21(1,1,i)=1(i)^2/(2*EI(i));us21(1,2,i)=(1(i)^3*(1-rr(i))/(6*EI(i));us21(2,1,i)=1(i)/EI(i);
us21(2,2,i)=1(i)^2/(2*EI(I));
us22(1,1,i)=1;us22(1,2,i)=1(i);us22(2,1,i)=0;us22(2,2,i)=1;
end
%============
%此处全为计算中间量
%============
fori=1:
1:
N+2
Su(1,1,i)=0;Su(1,2,i)=0;Su(2,1,i)=0;Su(2,2,i)=0;
Sn(1,1,i)=0;Sn(1,2,i)=0;Sn(2,1,i)=0;Sn(2,2,i)=0;
SS(1,1,i)=0;SS(1,2,i)=0;SS(2,1,i)=0;SS(2,2,i)=0;
end
fori=1:
1:
2
forj=1:
1:
2
SS1(i,j)=0;
Ud11(i,j)=0;Ud12(i,j)=0;Ud21(i,j)=0;Ud22(i,j)=0;
Us11(i,j)=0;Us12(i,j)=0;Us21(i,j)=0;Us22(i,j)=0;
end
end
%============
%(e)调用函数cone_modify修改锥轴的传递矩阵
%============
cone_modify(4,wi);
cone_modify(5,wi);
cone_modify(6,wi);
cone_modify(7,wi);
cone_modify(8,wi);
cone_modify(16,wi);
cone_modify(17,wi);
cone_modify(18,wi);
cone_modify(19,wi);
cone_modify(22,wi);
cone_modify(24,wi);
%============
%(f)形成最终传递矩阵
%============
%Ud11Ud12Ud21Ud22为最终参与计算的传递矩阵
fori=1:
1:
N
u11(:
:
i)=us11(:
:
i)*ud11(:
:
i)+us12(:
:
i)*ud21(:
:
i);
u12(:
:
i)=us11(:
:
i)*ud12(:
:
i)+us12(:
:
i)*ud22(:
:
i);
u21(:
:
i)=us21(:
:
i)*ud11(:
:
i)+us22(:
:
i)*ud21(:
:
i);
u22(:
:
i)=us21(:
:
i)*ud12(:
:
i)+us22(:
:
i)*ud22(:
:
i);
end
u11(:
:
N+1)=ud11(:
:
N+1);u12(:
:
N+1)=ud12(:
:
N+1);
u21(:
:
N+1)=ud21(:
:
N+1);u22(:
:
N+1)=ud22(:
:
N+1);
fori=1:
1:
2
forj=1:
1:
2
SS1(i,j)=0;
end
end
fori=1:
1:
N+1
ud11=u11(:
:
i);ud12=u12(:
:
i);ud21=u21(:
:
i);ud22=u22(:
:
i);
SS(:
:
:
i+1)=(ud11*SS1+ud12)*inv(ud21*SS1+ud22);
Su(:
:
i)=ud21*SS1+ud22;
Sn(:
:
i)=inv(ud21*SS1+ud22);%计算振型时用到
SS1=SS(:
:
i+1);
end
%======
(2)计算剩余量======
D1=det(SS(:
:
N+2);
fori=1:
1:
N+1
D1=D1*sign(det(Su(:
:
i));%消奇点
end
%======
(2)不平衡响应值EE======
EE(:
:
n+2)=-inv(SS(:
:
N+2)*PP(:
:
N+2);
fori=N+1:
-1:
1
EE(:
:
I)=Sn(:
:
i)*EE(:
:
i+1)-Sn(:
:
i)*UF(:
:
i);
end
A.2碰摩转子系统计算仿真程序
%main.m
%该程序主要完成完成jeffcott转子圆周碰摩故障仿真
%===========第一步:
设置初始条件
%rub_sign:
碰摩标志,若rub_sign=0,说明系统无碰摩故障;否则rub_sign=1
%loca:
不平衡质量的位置
%loc_rub:
碰摩位置
%Famp:
不平衡质量的大小单位为:
[g]
%wi:
转速单位为:
[rad]
%r:
偏心半径单位为:
[mm]
%Fampl:
离心力的大小单位为:
[kg,m]
%fai:
不平衡量的初始相位[rad]
clc
clear
[rub_signlocaloc_rubFam
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 传递 矩阵 matlab 程序
![提示](https://static.bdocx.com/images/bang_tan.gif)