[论文解读]二阶非线性不确定系统的PID控制器设计

PID controller design for second order nonlinear uncertain systems
摘要 复现一下论文的推导过程，补充一些推导细节，证明PID可以控制满足一定条件的二阶非线性系统并全局稳定，最后跑个仿真。

  对二阶非线性系统（相对阶为2）
$$\begin{array}{rl}{\dot{x}}_1=& x_2\\ {\dot{x}}_2=& f(x_1,x_2,t)+u(t)\end{array}$$
其中 $f(x_1,x_2,t)$ 满足
$$\left|\frac{\partial f}{\partial x_1}\right|\le L_1,\left|\frac{\partial f}{\partial x_2}\right|\le L_2$$
取控制信号为
$$u(t)=k_{p}e(t)+k_i{\int }_0^{t}e(\tau )\text{d}\tau +k_d\dot{e}(t)$$
其中 $e(t)=y^{\ast }-x_1(t)$ 为待跟踪信号 $y^{\ast }$ 与系统输出 $x_1(t)$ 的差。定义
$$\begin{array}{rl}& x(t)={\int }_0^{t}e(\tau )\text{d}\tau +\frac{f(y^{\ast },0,0)}{k_i}\\ & y(t)=e(t)\\ & z(t)=\dot{e}(t)\\ & g(y,z,t)=-f(y^{\ast }-y,-z,t)+f(y^{\ast },0,0)\end{array}$$
则
$$\begin{array}{rl}\dot{z}(t)=& \ddot{e}=-{\ddot{x}}_1\\ =& -f(x_1,x_2,t)-\left(k_{p}e(t)+k_i{\int }_0^{t}e(\tau )\text{d}\tau +k_d\dot{e}(t)\right)\\ =& -f(y^{\ast }-e,-\dot{e},t)-k_{p}y-k_{i}x+f(y^{\ast },0,0)-k_{d}z\\ =& f(y^{\ast },0,0)-f(y^{\ast }-y,-z,t)-k_{i}x-k_{p}y-k_{d}z\end{array}$$
即
$$\{\begin{array}{rl}& \dot{x}=y\\ & \dot{y}=z\\ & \dot{z}=g(y,z,t)-k_{i}x-k_{p}y-k_{d}z\end{array}$$
其中 $g(y,z,t)$ 为非线性函数，可以使用局部线性化的方法将其中的 $y$ 和 $z$ 分离，即
$$g(y,z,t)=b(y,t)y+a(y,z,t)z$$
其中
$$b(y,t)=\{\begin{array}{rl}\frac{g(y,0,t)}{y}& ,y\ne 0\\ {\frac{\partial g}{\partial y}|}_{(0,0,t)}& ,y=0\end{array},a(y,z,t)=\{\begin{array}{rlr}& \frac{g(y,z,t)-g(y,0,t)}{z}& ,z\ne 0\\ & {\frac{\partial g}{\partial z}|}_{(y,0,t)}& ,z=0\end{array},$$
于是
$$\left[\begin{array}{c}\dot{x}\\ \dot{y}\\ \dot{z}\end{array}\right]=\left[\begin{array}{ccc}0& 1& 0\\ 0& 0& 1\\ -k_i& b(y,t)-k_p& a(y,z,t)-k_d\end{array}\right]\left[\begin{array}{c}x\\ y\\ z\end{array}\right]$$
下面用李雅普诺夫函数证明稳定性。
取PID参数 $k_p$，$k_i$，$k_d$ 满足下面几个条件
$$\begin{array}{rl}& k_p>L_1,k_d>L_2,k_i>0,\\ & (k_p-L_1)(k_d-L_2)-k_i>L_2\sqrt{k_i(k_d+L_2)}\end{array}\text{\hspace{1em}(1)}$$
定义
ψ 0 = inf ⁡ y , z , t { k d − a ( y , z , t ) } ψ 1 = sup ⁡ y , z , t { k d − a ( y , z , t ) } ψ = ψ 0 + ψ 1 2 ϕ ( y ) = k p − b ( y ) ϕ 0 = inf ⁡ y ϕ ( y ) μ = 2 ( ϕ 0 ψ 0 + k i ) 4 ϕ 0 + L 2 2 (2) \begin{aligned} \psi_0 =& \inf_{y,z,t}\{k_d-a(y,z,t)\} \\ \psi_1 =& \sup_{y,z,t}\{k_d-a(y,z,t)\} \\ \psi =& \frac{\psi_0+\psi_1}2 \\ \phi(y) =& k_p-b(y) \\ \phi_0 =& \inf_y\phi(y) \\ \mu =& \frac{2(\phi_0\psi_0+k_i)}{4\phi_0+L_2^2} \end{aligned} \tag{2} ψ0​=ψ1​=ψ=ϕ(y)=ϕ0​=μ=​y,z,tinf​{kd​−a(y,z,t)}y,z,tsup​{kd​−a(y,z,t)}2ψ0​+ψ1​​kp​−b(y)yinf​ϕ(y)4ϕ0​+L22​2(ϕ0​ψ0​+ki​)​​(2)
其中 $\mu$ 是个大于0的常数。取李雅普诺夫函数
$$V(x,y,z)=\left[\begin{array}{ccc}x& y& z\end{array}\right]P\left[\begin{array}{c}x\\ y\\ z\end{array}\right]+{\int }_0^y(\varphi (s)-{\varphi }_0)s\text{d}s$$
其中
$$P=\frac{1}{2}\left[\begin{array}{ccc}\mu k_i& k_i& 0\\ k_i& {\varphi }_0+\mu \psi & \mu \\ 0& \mu & 1\end{array}\right]$$
为半正定矩阵。下面证明 $P$ 为半正定矩阵。根据式(2)的定义和式(1)的约束条件可得
$$\begin{array}{rl}& {\varphi }_0\ge k_p-L_1>0,\\ & k_d+L_2\ge \psi \ge {\psi }_0\ge k_d-L_2>0,\\ & \mu >0,\\ & {\varphi }_0{\psi }_0-k_i>(k_p-L_1)(k_d-L_2)-k_i\\ & >L_2\sqrt{k_i(k_d+L_2)}>L_2\sqrt{k_i{\psi }_0}>0\end{array}$$
$$\begin{array}{rl}{\varphi }_0{\psi }_0-k_i>& L_2\sqrt{k_i{\psi }_0}\\ ({\varphi }_0{\psi }_0-k_i{)}^2>& L_2^2{k}_i{\psi }_0\\ ({\varphi }_0{\psi }_0+k_i{)}^2>& L_2^2{k}_i{\psi }_0+4k_i{\varphi }_0{\psi }_0\end{array}$$
下面用这些条件和结果证明如下3个式子
$$\begin{array}{rl}& \mu <{\psi }_0\\ & 4(\mu {\varphi }_0-k_i)({\psi }_0-\mu )>{\mu }^2{L}_2^2\\ & \mu {\varphi }_0>k_i\end{array}\text{\hspace{1em}(3)}$$
第3个式子可以由前两个得到。
$$\begin{array}{rl}\mu -{\psi }_0=& \frac{2({\varphi }_0{\psi }_0+k_i)-{\psi }_0(4{\varphi }_0+L_2^2)}{4{\varphi }_0+L_2^2}\\ =& \frac{2k_i-2{\varphi }_0{\psi }_0-{\psi }_0{L}_2^2}{4{\varphi }_0+L_2^2}\\ <& \frac{-{\psi }_0{L}_2^2}{4{\varphi }_0+L_2^2}<0\end{array}$$
$$\begin{array}{rl}& 4(\mu {\varphi }_0-k_i)({\psi }_0-\mu )-{\mu }^2{L}_2^2\\ =& -(4{\varphi }_0+L_2^2){\mu }^2+4({\varphi }_0{\psi }_0+k_i)\mu -4k_i{\psi }_0\\ =& \frac{1}{4{\varphi }_0+L_2^2}[-4({\varphi }_0{\psi }_0+k_i{)}^2+8({\varphi }_0{\psi }_0+k_i{)}^2-4k_i{\psi }_0(4{\varphi }_0+L_2^2)]\\ =& \frac{4}{4{\varphi }_0+L_2^2}[({\varphi }_0{\psi }_0+k_i{)}^2-4k_i{\varphi }_0{\psi }_0-k_i{\psi }_0{L}_2^2]>0\end{array}$$
于是矩阵 $P$ 的三阶顺序主子式的行列式
$$\begin{array}{rl}& \mu k_i>0\\ & \det \left[\begin{array}{cc}\mu k_i& k_i\\ k_i& {\varphi }_0+\mu \psi \end{array}\right]=\mu k_i({\varphi }_0+\mu \psi )-k_i^2\ge k_i(\mu {\varphi }_0+{\mu }^2{\psi }_0-k_i)>k_i{\mu }^2{\psi }_0>0\\ & \det \left[\begin{array}{ccc}\mu k_i& k_i& 0\\ k_i& {\varphi }_0+\mu \psi & \mu \\ 0& \mu & 1\end{array}\right]=k_i(\mu {\varphi }_0+{\mu }^2\psi -k_i-{\mu }^3)>k_i({\mu }^2{\psi }_0-{\mu }^3)=k_i{\mu }^2({\psi }_0-\mu )>0\end{array}$$
所以 $P$ 正定。
$$\begin{array}{rl}& A^{\text{T}}P+PA\\ =& \frac{1}{2}\left[\begin{array}{ccc}0& 0& -k_i\\ 1& 0& b-k_p\\ 0& 1& a-k_d\end{array}\right]\left[\begin{array}{ccc}\mu k_i& k_i& 0\\ k_i& {\varphi }_0+\mu \psi & \mu \\ 0& \mu & 1\end{array}\right]+\frac{1}{2}\left[\begin{array}{ccc}\mu k_i& k_i& 0\\ k_i& {\varphi }_0+\mu \psi & \mu \\ 0& \mu & 1\end{array}\right]\left[\begin{array}{ccc}0& 1& 0\\ 0& 0& 1\\ -k_i& b-k_p& a-k_d\end{array}\right]\\ =& \frac{1}{2}\left[\begin{array}{ccc}0& -\mu k_i& -k_i\\ \mu k_i& k_i+\mu (b-k_p)& b-k_p\\ k_i& {\varphi }_0+\mu (\psi +a-k_d)& \mu +a-k_d\end{array}\right]+\frac{1}{2}\left[\begin{array}{ccc}0& \mu k_i& k_i\\ -\mu k_i& k_i+\mu (b-k_p)& {\varphi }_0+\mu (\psi +a-k_d)\\ -k_i& b-k_p& \mu +a-k_d\end{array}\right]\\ =& \left[\begin{array}{ccc}0& 0& 0\\ 0& k_i+\mu (b-k_p)& \frac{1}{2}(b-k_p+{\varphi }_0+\mu (\psi +a-k_d))\\ 0& \frac{1}{2}(b-k_p+{\varphi }_0+\mu (\psi +a-k_d))& \mu +a-k_d\end{array}\right]\end{array}$$
$$\begin{array}{rl}\dot{V}=& {\vec{x}}^{\text{T}}(A^{\text{T}}P+PA)\vec{x}+(\varphi (y)-{\varphi }_0)y\dot{y}\\ =& \left[\begin{array}{ccc}x& y& z\end{array}\right]\left[\begin{array}{ccc}0& 0& 0\\ 0& k_i-\mu \varphi (y)& \frac{1}{2}({\varphi }_0-\varphi (y)+\mu (\psi +a-k_d))\\ 0& \frac{1}{2}({\varphi }_0-\varphi (y)+\mu (\psi +a-k_d))& \mu +a-k_d\end{array}\right]\left[\begin{array}{c}x\\ y\\ z\end{array}\right]+(\varphi (y)-{\varphi }_0)yz\\ =& \left[\begin{array}{cc}y& z\end{array}\right]\left[\begin{array}{cc}{k}_i-\mu \varphi (y)& \frac{1}{2}\mu (\psi +a-k_d)\\ \frac{1}{2}\mu (\psi +a-k_d)& \mu +a-k_d\end{array}\right]\left[\begin{array}{c}y\\ z\end{array}\right]+\left[\begin{array}{cc}y& z\end{array}\right]\left[\begin{array}{cc}0& \frac{1}{2}({\varphi }_0-\varphi (y))\\ \frac{1}{2}({\varphi }_0-\varphi (y))& 0\end{array}\right]\left[\begin{array}{c}y\\ z\end{array}\right]+(\varphi (y)-{\varphi }_0)yz\\ =& -\left[\begin{array}{cc}y& z\end{array}\right]Q\left[\begin{array}{c}y\\ z\end{array}\right]\end{array}$$
其中
$$Q=\left[\begin{array}{cc}\mu \varphi (y)-k_i& -\frac{1}{2}\mu (\psi +a-k_d)\\ -\frac{1}{2}\mu (\psi +a-k_d)& k_d-\mu -a(y,z,t)\end{array}\right]$$
定义
$$\begin{array}{rl}\alpha =& \frac{1}{2}(\mu (\psi +a(y,z,t)-k_d))\\ \beta =& k_d-a(y,z,t)-\mu \end{array}$$
根据式(2)的定义和式(3)的结果得到
$$\begin{array}{rl}& \mu \varphi (y)-k_i\ge \mu {\varphi }_0-k_i>0\\ & {\psi }_1-\mu \ge \beta \ge {\psi }_0-\mu >0\\ & |\psi +a(y,z,t)-k_d|\le L_2\\ & (\mu \varphi (y)-k_i)\beta \ge (\mu {\varphi }_0-k_i)({\psi }_0-\mu )\ge \frac{1}{4}{\mu }^2{L}_2^2\ge {\alpha }^2\end{array}$$
于是 $Q$ 的所有顺序主子式大于0，$Q$ 正定。

补充说明与仿真

取
$$f(x_1,x_2)=2\ln (|x_1|+1)+3\sin x_2$$
则 $L_1=2,L_2=3$，取 $k_p=k_d=4,k_i=0.01$，PID参数需要满足式(1)。但不严格满足也大部分情况下能稳定。另外当待跟踪信号 $y^{\ast }$ 不是常数时不满足 $\dot{x}=y$，证明过程不成立。

$$f(x_1,x_2)=2\ln (|x_1|+1)+3\sin x_2$$
  simucpp 代码如下。

#include <iostream>
#include <cmath>
#include "simucpp.hpp"
using namespace simucpp;
using namespace zhnmat;
using namespace std;

class PIDcontroller: public PackModule
{
public:
    PIDcontroller() {};
    PIDcontroller(Simulator *sim){
        Initialize(sim);
    };
    virtual PUnitModule Get_InputPort(uint32_t n) const {
        if (n==0) return in1;
        if (n==1) return in2;
        if (n==2) return in3;
        if (n==3) return in4;
        return nullptr;
    };
    virtual PUnitModule Get_OutputPort(uint32_t n) const {
        if (n==0) return sum3;
        return nullptr;
    };
    void Initialize(Simulator *sim){
        SUGain(in1, sim);
        SUGain(in2, sim);
        SUGain(in3, sim);
        SUGain(in4, sim);
        SUIntegrator(interr, sim);
        SUSum(sum1, sim);
        SUSum(sum2, sim);
        SUSum(sum3, sim);
        SUGain(gainKp, sim);
        SUGain(gainKi, sim);
        SUGain(gainKd, sim);
        sim->connectU(in3, sum1);
        sim->connectU(in1, sum1);
        sim->connectU(in4, sum2);
        sim->connectU(in2, sum2);
        sim->connectU(sum1, interr);
        sim->connectU(sum1, gainKp);
        sim->connectU(sum2, gainKd);
        sim->connectU(interr, gainKi);
        sim->connectU(gainKp, sum3);
        sim->connectU(gainKi, sum3);
        sim->connectU(gainKd, sum3);
        sum1->Set_InputGain(-1);
        sum2->Set_InputGain(-1);
        gainKp->Set_Gain(4);
        gainKi->Set_Gain(0.1);
        gainKd->Set_Gain(4);
    }

    UGain *in1=nullptr;  // x1
    UGain *in2=nullptr;  // x2
    UGain *in3=nullptr;  // u
    UGain *in4=nullptr;  // u'
    UIntegrator *interr=nullptr;
    USum *sum1=nullptr;
    USum *sum2=nullptr;
    USum *sum3=nullptr;
    UGain *gainKp=nullptr;
    UGain *gainKi=nullptr;
    UGain *gainKd=nullptr;
};

int main() {
    Simulator sim1;
    FUIntegrator(intx1, &sim1);
    FUIntegrator(intx2, &sim1);
    FUFcnMISO(misof, &sim1);
    FUInput(inv, &sim1);
    FUInput(indotv, &sim1);
    FUOutput(outu, &sim1);
    FUOutput(outv, &sim1);
    FUOutput(outy, &sim1);
    auto *controller = new PIDcontroller(&sim1);
    sim1.connectU(intx2, intx1);
    sim1.connectU(misof, intx2);
    sim1.connectU(intx1, misof);
    sim1.connectU(intx2, misof);
    sim1.connectU(controller, 0, misof);
    sim1.connectU(intx1, controller, 0);
    sim1.connectU(intx2, controller, 1);
    sim1.connectU(inv, controller, 2);
    sim1.connectU(indotv, controller, 3);
    sim1.connectU(controller, 0, outu);
    sim1.connectU(inv, outv);
    sim1.connectU(intx1, outy);
    inv->Set_Function([](double t){
        return 1;
    });
    indotv->Set_Function([](double t){
        return 0;
    });
    misof->Set_Function([](double *u){
        double x1 = u[0];
        double x2 = u[1];
        double uin = u[2];
        return 2*log(abs(x1)+1)+3*sin(x2)+uin;
    });
    outu->Set_EnableStore(false);
    outv->Set_EnableStore(false);
    // sim1.Set_SampleTime(0.1);
    sim1.Set_SimStep(0.01);
    sim1.Initialize();
    sim1.Simulate(10);
    sim1.Plot();
    return 0;
}