python control控制系统仿真再更新：使用内外环系统进行垂直起降飞行器横向控制仿真

control库的安装方法为pip install control
原文链接地址

使用内/外循环设计进行横向控制

注意，无需逐段复制代码，直接跳转到最后复制完整代码即可。
本节演示了矢量推力飞机示例的环路整形控制器的设计。此示例摘自[http://www.cds.caltech.edu/~murray/amwiki Astrom and Murray]的第11章（频域设计）。此处列出的python代码包含在文件pvtol-nested.py中。

为了设计矢量化推力飞机横向动力学的控制器，我们使用了“内/外”环设计方法。我们首先使用框图表示动态

其中，为了简化，我们令以下物理量：

$$H_{\theta u_1}=\frac{r}{J{s}^2}，H_{x{u}_1}=\frac{J{s}^2-mgr}{J{s}^2（m{s}^2+cs）}$$

通过将动力学和控制器分为两个部分来构造控制器：
1、 由侧倾动力学$P_i$和控制$C_i$组成的内循环
2、由侧向位置动力学$P_o$E和控制器$C_o$组成的外循环。

内环$H_i$使用矢量推力控制飞机的侧倾角，而外环控制器$C_o$命令侧倾角调节侧向位置。

以下代码导入所需的库并定义动力学特性：

from matplotlib.pyplot import * # Grab MATLAB plotting functions
from control.matlab import *    # MATLAB-like functions

＃系统参数
m = 4；＃飞机质量
J = 0.0475；＃围绕俯仰轴的惯性
r = 0.25；＃距力中心的距离
g = 9.8；＃重力常数
c = 0.05；＃阻尼系数（估计值）

＃动态传递函数
Pi = tf（[r]，[J，0，0]）; ＃内部循环（滚动）
Po = tf（[1]，[m，c，0]）; ＃外循环（位置）

对于内层反馈循环，使用超前补偿器

k = 200；a = 2；b = 50
Ci = k * tf（[1，a]，[1，b]）＃超前补偿器
Li= Pi * Ci

内环的闭环动力学由$H_i$给出

Hi = parallel(feedback(Ci, Pi), -m*g*feedback(Ci*Pi, 1));

最后，我们使用另一个主补偿器设计横向补偿器

＃现在设计横向控制系统
a = 0.02；b = 5；K = 2；
Co ＝ -K * tf（[1，0.3]，[1，10]）；＃另一个主要补偿器
Lo = -m * g * Po * Co;

可以使用灵敏度函数和互补灵敏度函数来表征系统的性能

L = Co*Hi*Po;
S = feedback(1, L);
T = feedback(L, 1);

使用以下命令计算环路传递函数的频率响应和奈奎斯特图

bode(L）
nyquist(L, (0.0001, 1000));
axis([-700, 5300, -3000, 3000]);

gangof4(Hi*Po, Co);

相应的图如下所示：

伯德图

奈奎斯特图

四个图

标准答案是这样的，但是我未能复现这个图。上两个图我都可以复现。

代码实现

# pvtol-nested.py - inner/outer design for vectored thrust aircraft
# RMM, 5 Sep 09
#
# This file works through a fairly complicated control design and
# analysis, corresponding to the planar vertical takeoff and landing
# (PVTOL) aircraft in Astrom and Murray, Chapter 11.  It is intended
# to demonstrate the basic functionality of the python-control
# package.
#

from __future__ import print_function

import os
import matplotlib.pyplot as plt  # MATLAB plotting functions
from control.matlab import *    # MATLAB-like functions
import numpy as np

# System parameters
m = 4               # 飞行器质量
J = 0.0475          # 俯仰轴惯性
r = 0.25            # 力与质心距离（原文为：distance to center of force
g = 9.8             # 重力加速度
c = 0.05            # 阻尼系数（估计）

# 传递函数
Pi = tf([r], [J, 0, 0])  # inner loop (roll)
Po = tf([1], [m, c, 0])  # outer loop (position)

# Use state space versions
Pi = tf2ss(Pi)
Po = tf2ss(Po)

#
# 内控制环路设计
#
# This is the controller for the pitch dynamics.  Goal is to have
# fast response for the pitch dynamics so that we can use this as a 
# control for the lateral dynamics
#

# Design a simple lead controller for the system
k, a, b = 200, 2, 50
Ci = k*tf([1, a], [1, b])  # lead compensator
Li = Pi*Ci

Si = feedback(1, Li)
Ti = Li*Si

Hi = parallel(feedback(Ci, Pi), -m*g*feedback(Ci*Pi, 1))

a, b, K = 0.02, 5, 2
Co = -K*tf([1, 0.3], [1, 10])  # another lead compensator
Lo = -m*g*Po*Co



# Finally compute the real outer-loop loop gain + responses
L = Co*Hi*Po
S = feedback(1, L)
T = feedback(L, 1)

plt.figure('时域响应')
y,t=step(T,T=np.linspace(0,10,100))
plt.plot(t,y)
plt.xlabel('time/s')
plt.ylabel('y(t)')
plt.grid()

# 对整个系统绘制伯德图
plt.figure('伯德图')
bode(L, np.logspace(-4, 3))

plt.figure('奈奎斯特图')

nyquist(L, (0.0001, 1000))
plt.axis([-4000, 300, -300000, 300000])

plt.figure('四个图')
gangof4(Hi*Po, Co)
plt.show()