《飞行控制系统》课程实验报告班级 0314102学号 ********* 姓名孙旭东成绩南京航空航天大学2017年4月（一）飞机纵向飞行控制系统的设计与仿真1、分析飞机纵向动力学模态，求飞机的长周期与短周期阻尼与自然频率。

在MATLAB环境下导入数据文件，输入damp(alon),得出结果：Eigenvalue Damping Freq. (rad/s)-2.29e+000 + 4.10e+000i 4.88e-001 4.69e+000-2.29e+000 - 4.10e+000i 4.88e-001 4.69e+000-3.16e-002 1.00e+000 3.16e-002-7.30e-003 + 3.35e-002i 2.13e-001 3.42e-002-7.30e-003 - 3.35e-002i 2.13e-001 3.42e-002长周期的根为 -7.30e-003 + 3.35e-002i 和 -7.30e-003 - 3.35e-002i阻尼为 2.13e-001自然频率为 3.42e-002(rad/s)短周期的根为 -2.29e+000 + 4.10e+000i 和 -2.29e+000 - 4.10e+000i阻尼为 4.88e-001自然频率为 4.69e+000(rad/s)2、对升降舵及油门单位阶跃输入下的飞机自然特性进行仿真，画出相应的状态曲线。

sys=ss(alon,blon,clon,dlon)[y,t]=step(sys,500)subplot(221)plot(t,y(:,1,1))xlabel('t(s)')ylabel('\Deltau(m/s)')subplot(222)plot(t,y(:,1,2))xlabel('t(s)')ylabel('\Deltau(m/s)')subplot(223)plot(t,y(:,2,1))xlabel('t(s)')ylabel('\Delta\alpha(deg)')subplot(224)plot(t,y(:,2,2))xlabel('t(s)')ylabel('\Delta\alpha(deg)')200400600-10-505t(s)∆q (d e g /s )200400600-4-2024t(s)∆q (d e g /s )200400600-150-100-50050t(s)∆θ(d e g )0200400600-50050100t(s)∆θ(d e g )200400600-2000200400t(s)∆u (m /s )0200400600-6-4-2t(s)∆α(d e g )200400600-2000200400t(s)∆u (m /s )0200400600-2024t(s)∆α(d e g )subplot(221) plot(t,y(:,3,1)) xlabel('t(s)')ylabel('\Deltaq(deg/s)') subplot(222) plot(t,y(:,3,2)) xlabel('t(s)')ylabel('\Deltaq(deg/s)') subplot(223) plot(t,y(:,4,1)) xlabel('t(s)')ylabel('\Delta\theta(deg)') subplot(224) plot(t,y(:,4,2)) xlabel('t(s)')ylabel('\Delta\theta(deg)')subplot(121) plot(t,y(:,5,1)) xlabel('t(s)')ylabel('\Deltah(m)') subplot(122) plot(t,y(:,5,2)) xlabel('t(s)')ylabel('\Deltah(m)')2004006004t(s)∆h (m )200400600-2.5-2-1.5-1-0.54t(s)∆h (m )以上各图为升降舵及油门单位阶跃输入下的飞机自然特性行仿真，左边一列为升降舵的阶跃输入，右边一列为油门的阶跃输入。

3、采用短周期简化方法，求出传递函数()e qG s δ∆∆。

采用根轨迹方法设计飞机的俯仰角控制系统，并进行仿真。

输入命令：a1=alon((2:3),(2:3)) b1=blon((2:3),:)c1=clon((2:3),(2:3)) d1=dlon((2:3),:)[n,d]=ss2tf(a1,b1,c1,d1,1) g1=tf(n(2,:),d)得到传递函数()e qG s δ∆∆为 ：-34.17 s - 82.55----------------------- s^2 + 4.579 s + 22.01根轨迹设计： 输入命令：g1=tf(n(2,:),d)g2=tf([-10],[1 10]) g3=series(g1,g2) sisotool(g3)1010101010Frequency (rad/sec)P h a s e (d e g)Open-Loop Bode E ditor for Open Loop 1 (OL1)M a g n i t u d e (d B )Root Locus E ditor for Open Loop 1 (OL1)1010101010Frequency (rad/sec)P h a s e (d e g )Open-Loop Bode E ditor for Open Loop 1 (OL1)M a g n i t u d e (d B )1010101010Frequency (rad/sec)Open-Loop Bode E ditor for Open Loop 1 (OL1)Root Locus E ditor for Open Loop 1 (OL1)Real AxisI m a g A x i s同样，可得Kth=1在Simulink 中搭建系统仿真模型：进行仿真：0.20.40.60.81t(s)∆θ(d e g )00.20.40.60.811.21.4∆ut(s)4、基于长周期简化方法，求出传递函数()T uG s δ∆∆，设计飞机的速度控制系统，并进行仿真。

输入命令：a1=alon([1,4],[1,4]) b1=blon([1,4],:) c1=clon([1,4],[1,4]) d1=dlon([1,4],:)[n,d]=ss2tf(a1,b1,c1,d1,2); g1=tf(n(1,:),d)得到传递函数为：7.971 s --------------- s^2 + 0.04847 s在Simulink 中搭建系统模型：使用经验试凑法得到PID 控制器参数：Kp=0.9 Ki=0.2 Kd=0仿真结果如下：5、基于纵向线性模型（状态方程），分别对速度控制与俯仰角控制进行仿真。

在Simulink中搭建仿真模型：先在速度通道加阶跃信号，输入命令：subplot(221)plot(t,x1)xlabel('t(s)')ylabel('\Deltau(m/s)')subplot(222)plot(t,x2)xlabel('t(s)')ylabel('\Delta\alpha(deg)')subplot(223)plot(t,x3)xlabel('t(s)')ylabel('\Deltaq(deg/s)')subplot(224)plot(t,x4)xlabel('t(s)')ylabel('\Delta\theta(deg)')和plot(t,x5)xlabel('t(s)')ylabel('\Deltah(m)')得到以下曲线：t(s)∆h (m)0510-0.015-0.01-0.0050t(s)∆u (m /s)0510-0.500.51t(s)∆α(d e g)-10123t(s)∆q (d e g /s)0.511.5t(s)∆θ(d e g)5100.511.5t(s)∆u (m /s)510-0.4-0.200.20.4t(s)∆α(d e g)510-4-2024t(s)∆q (d e g /s)510-0.200.20.40.6t(s)∆θ(d e g )再在俯仰角通道加阶跃信号，重复以上命令，得到如下曲线：0246810-30-20-10010t(s)h (m )（二）飞机侧向滚转角控制系统设计1、求出侧向运动方程的特征根，及对应的模态，求出荷兰滚模态的阻尼及自然频率。

在MATLAB环境下导入数据文件，输入damp(alon),得出结果：Eigenvalue Damping Freq. (rad/s)0.00e+000 -1.00e+000 0.00e+000-6.89e+000 1.00e+000 6.89e+000-1.55e-002 1.00e+000 1.55e-002-1.02e+000 + 5.08e+000i 1.97e-001 5.19e+000-1.02e+000 - 5.08e+000i 1.97e-001 5.19e+000侧向运动方程的特征根为：0.00e+000（航向随遇平衡模态）-1.55e-002（螺旋模态）-1.02e-001 +5.08e+000i，-1.02e-001 – 5.08e+000i（荷兰滚模态）-6.89e+000（侧向滚转收敛模态）荷兰滚模态的阻尼为：1.97e-001自然频率为：5.19e+000(rad/s)2、对副翼与方向舵单位阶跃输入下的自然特性进行仿真sys=ss(alat,blat,clat,dlat)[y,t]=step(sys,400)subplot(221)plot(t,y(:,1,1))xlabel('t(s)')ylabel('\Delta\beta(deg)')subplot(222)plot(t,y(:,1,2))xlabel('t(s)')ylabel('\Delta\beta(deg)')subplot(223)plot(t,y(:,2,1))xlabel('t(s)')ylabel('\Deltap(deg/s)')subplot(224)plot(t,y(:,2,2))xlabel('t(s)')ylabel('\Deltap(deg/s)')得到以下曲线：100200300400-3-2-101t(s)∆β(d e g )0100200300400-0.500.511.5t(s)∆β(d e g )100200300400-30-20-100t(s)∆p (d e g /s )0100200300400-15-10-5t(s)∆p (d e g /s )subplot(221) plot(t,y(:,3,1)) xlabel('t(s)')ylabel('\Deltar(deg/s)') subplot(222) plot(t,y(:,3,2)) xlabel('t(s)')ylabel('\Deltar(deg/s)') subplot(223) plot(t,y(:,4,1)) xlabel('t(s)')ylabel('\Delta\phi(deg)') subplot(224) plot(t,y(:,4,2)) xlabel('t(s)')ylabel('\Delta\phi(deg)')得到以下曲线：100200300400-2.5-2-1.5-1-0.504t(s)∆ψ(d e g )200400t(s)∆ψ(d e g )100200300400-80-60-40-200t(s)∆r (d e g /s )100200300400-40-2020t(s)∆r (d e g /s )100200300400-1500-1000-5000t(s)∆φ(d e g )0100200300400-800-600-400-2000t(s)∆φ(d e g )subplot(121) plot(t,y(:,5,1)) xlabel('t(s)')ylabel('\Delta\psi(deg)') subplot(122) plot(t,y(:,5,2)) xlabel('t(s)')ylabel('\Delta\psi(deg)'得到以下曲线：以上各图中左边为副翼输入单位阶跃响应的曲线，右边为方向舵输入单位阶跃响应的曲线。