Chapter10.3：描述函数法

该系列博客主要讲述Matlab软件在自动控制方面的应用，如无自动控制理论基础，请先学习自动控制系列博文，该系列博客不再详细讲解自动控制理论知识。
自动控制理论基础相关链接：https://blog.csdn.net/qq_39032096/category_10287468.html?spm=1001.2014.3001.5482
博客参考书籍：《MATLAB/Simulink与控制系统仿真》。

3.描述函数法

3.1 描述函数的基本概念

3.1.1 描述函数的定义

设非线性环节输入输出描述为：
$$y=f(x)$$
当非线性环节的输入信号为正弦信号：
$$x(t)=A\sin \omega t$$
对非线性环节的稳态输出$y(t)$进行谐波分析；

一般情况下，$y(t)$为非正弦的周期信号，可以展开成傅里叶级数：
y(t)=A0+∑n=1∞(Ancos⁡nωt+Bnsin⁡nωt)=A0+∑n=1∞Ynsin⁡(nωt+φn)y(t)=A_0+\\sum_{n=1}^{\\infty}(A_n\\cos{n}\\omega{t}+B_n\\sin{n}\\omega{t})=A_0+\\sum_{n=1}^{\\infty}Y_n\\sin(n\\omega{t}+\\varphi_n) y(t)=A0​+n=1∑∞​(An​cosnωt+Bn​sinnωt)=A0​+n=1∑∞​Yn​sin(nωt+φn​)
其中：A0A_0A0​为直流分量；Ynsin⁡(nωt+φn)Y_n\\sin(n\\omega{t}+\\varphi_n)Yn​sin(nωt+φn​)为第$n$次谐波分量，且有：
Yn=An2+Bn2,φn=arctan⁡AnBnY_n=\\sqrt{A_n^2+B_n^2},\\varphi_n=\\arctan\\frac{A_n}{B_n} Yn​=An2​+Bn2​​,φn​=arctanBn​An​​
其中：An、BnA_n、B_nAn​、Bn​为傅里叶系数；
An=1π∫02πy(t)cos⁡nωtd(ωt)，Bn=1π∫02πy(t)sin⁡nωtd(ωt)(n=1,2,…,)A_n=\\frac{1}{\\pi}\\int_{0}^{2\\pi}y(t)\\cos{n}\\omega{t}{\\rm d}(\\omega{t})，B_n=\\frac{1}{\\pi}\\int_{0}^{2\\pi}y(t)\\sin{n}\\omega{t}{\\rm d}(\\omega{t})(n=1,2,\\dots,) An​=π1​∫02π​y(t)cosnωtd(ωt)，Bn​=π1​∫02π​y(t)sinnωtd(ωt)(n=1,2,…,)
直流分量：
A0=12π∫02πy(t)d(ωt)A_0=\\frac{1}{2\\pi}\\int_{0}^{2\\pi}y(t){\\rm d}(\\omega{t}) A0​=2π1​∫02π​y(t)d(ωt)
若A0=0A_0=0A0​=0且当$n>1$时，YnY_nYn​均很小，则可近似认为非线性环节的正弦响应应仅有基波分量：
y(t)≈A1cos⁡ωt+B1sin⁡ωt=Y1sin⁡(ωt+φ1)y(t)≈A_1\\cos\\omega{t}+B_1\\sin\\omega{t}=Y_1\\sin(\\omega{t}+\\varphi_1) y(t)≈A1​cosωt+B1​sinωt=Y1​sin(ωt+φ1​)
定义正弦输入信号作用下，非线性环节的稳态输出中基波分量和输入信号的复数比为非线性环节的描述函数，用$N(A)$表示，即：
N(A)=∣N(A)∣ej∠N(A)=Y1Aejφ1=B1+jA1AN(A)=|N(A)|{\\rm e}^{j\\angle{N(A)}}=\\frac{Y_1}{A}{\\rm e}^{j\\varphi_1}=\\frac{B_1+{\\rm j}A_1}{A} N(A)=∣N(A)∣ej∠N(A)=AY1​​ejφ1​=AB1​+jA1​​

3.1.2 非线性系统描述函数法分析的应用条件

1. 非线性系统应简化成一个非线性环节和一个线性部分闭环连接的典型结构图，如下图所示：

   

2. 非线性环节的输入输出特性$y(x)$应是$x$的奇函数，即$f(x)=-f(-x)$，或正弦输入下的输出为$t$的奇对称函数，即y(t+πω)=−y(t)y\\left(t+\\displaystyle\\frac{\\pi}{\\omega}\\right)=-y(t)y(t+ωπ​)=−y(t)，以保证非线性环节的正弦响应不含有常值分量，即A0=0A_0=0A0​=0；

3. 系统的线性部分应具有较好的低通滤波性能；线性部分的阶次越高，低通滤波性能越好；欲具有低通滤波性能，线性部分的极点应位于复平面的左半平面；

3.2 典型非线性特性的描述函数

3.2.1 死区饱和非线性环节

输出$y(t)$的数学表达式为：
y(t)={0,0≤ωt≤Ψ1K(Asin⁡ωt−Δ),Ψ1<ωt≤Ψ2K(a−Δ),Ψ2<ωt≤π2y(t)= \\begin{cases} 0,&0≤\\omega{t}≤\\Psi_1\\\\\\\\ K(A\\sin\\omega{t}-\\Delta),&\\Psi_1<\\omega{t}≤\\Psi_2\\\\\\\\ K(a-\\Delta),&\\Psi_2<\\omega{t}≤\\displaystyle\\frac{\\pi}{2} \\end{cases} y(t)=⎩⎨⎧​0,K(Asinωt−Δ),K(a−Δ),​0≤ωt≤Ψ1​Ψ1​<ωt≤Ψ2​Ψ2​<ωt≤2π​​
确定$y(t)$产生不同线性变换的区间端点为：
Ψ1=arcsin⁡ΔA，Ψ2=arcsin⁡aA\\Psi_1=\\arcsin\\frac{\\Delta}{A}， \\Psi_2=\\arcsin\\frac{a}{A} Ψ1​=arcsinAΔ​，Ψ2​=arcsinAa​
由于$y(t)$为奇函数，因此：A0=0,A1=0A_0=0,A_1=0A0​=0,A1​=0，又因为$y(t)$为半周期内对称，因此：
B1=1π∫02πy(t)sin⁡ωtd(ωt)=4π∫0π2y(t)sin⁡ωtd(ωt)=4π[∫Ψ1Ψ2[K(Asin⁡ωt−Δ)]sin⁡ωtd(ωt)+∫Ψ2π2K(a−Δ)sin⁡ωtd(ωt)]=2KAπ[arcsin⁡aA−arcsin⁡ΔA+aA1−(aA)2−ΔA1−(ΔA)2]\\begin{aligned} B_1&=\\frac{1}{\\pi}\\int_{0}^{2\\pi}y(t)\\sin\\omega{t}{\\rm d}(\\omega{t})=\\frac{4}{\\pi}\\int_0^{\\frac{\\pi}{2}}y(t)\\sin\\omega{t}{\\rm d}(\\omega{t})\\\\\\\\ &=\\frac{4}{\\pi}\\left[\\int_{\\Psi_1}^{\\Psi_2}[K(A\\sin\\omega{t}-\\Delta)]\\sin\\omega{t}{\\rm d}(\\omega{t})+\\int_{\\Psi_2}^{\\frac{\\pi}{2}}K(a-\\Delta)\\sin\\omega{t}{\\rm d}(\\omega{t})\\right]\\\\\\\\ &=\\frac{2KA}{\\pi}\\left[\\arcsin\\frac{a}{A}-\\arcsin\\frac{\\Delta}{A}+\\frac{a}{A}\\sqrt{1-\\left(\\frac{a}{A}\\right)^2}-\\frac{\\Delta}{A}\\sqrt{1-\\left(\\frac{\\Delta}{A}\\right)^2}\\right] \\end{aligned} B1​​=π1​∫02π​y(t)sinωtd(ωt)=π4​∫02π​​y(t)sinωtd(ωt)=π4​[∫Ψ1​Ψ2​​[K(Asinωt−Δ)]sinωtd(ωt)+∫Ψ2​2π​​K(a−Δ)sinωtd(ωt)]=π2KA​​arcsinAa​−arcsinAΔ​+Aa​1−(Aa​)2​−AΔ​1−(AΔ​)2​​​
死区饱和特性的描述函数为：
N(A)=2Kπ[arcsin⁡aA−arcsin⁡ΔA+aA1−(aA)2−ΔA1−(ΔA)2],A≥aN(A)=\\frac{2K}{\\pi}\\left[\\arcsin\\frac{a}{A}-\\arcsin\\frac{\\Delta}{A}+\\frac{a}{A}\\sqrt{1-\\left(\\frac{a}{A}\\right)^2}-\\frac{\\Delta}{A}\\sqrt{1-\\left(\\frac{\\Delta}{A}\\right)^2}\\right],A≥a N(A)=π2K​​arcsinAa​−arcsinAΔ​+Aa​1−(Aa​)2​−AΔ​1−(AΔ​)2​​,A≥a
取$\Delta =0$时，饱和特性的描述函数为：
N(A)=2Kπ[arcsin⁡aA+aA1−(aA)2],A≥aN(A)=\\frac{2K}{\\pi}\\left[\\arcsin\\frac{a}{A}+\\frac{a}{A}\\sqrt{1-\\left(\\frac{a}{A}\\right)^2}\\right],A≥a N(A)=π2K​[arcsinAa​+Aa​1−(Aa​)2​],A≥a
死区特性描述函数：
N(A)=2Kπ[π2−arcsin⁡ΔA−ΔA1−(ΔA)2],A≥ΔN(A)=\\frac{2K}{\\pi}\\left[\\frac{\\pi}{2}-\\arcsin\\frac{\\Delta}{A}-\\frac{\\Delta}{A}\\sqrt{1-\\left(\\frac{\\Delta}{A}\\right)^2}\\right],A≥\\Delta N(A)=π2K​​2π​−arcsinAΔ​−AΔ​1−(AΔ​)2​​,A≥Δ

3.2.2 死区与滞环继电特性非线性环节

输出$y(t)$的数学表达式：
y(t)={0,0≤ωt<Ψ1M,Ψ1≤ωt≤Ψ20,Ψ2<ωt≤πy(t)= \\begin{cases} 0,&0≤\\omega{t}<\\Psi_1\\\\\\\\ M,&\\Psi_1≤\\omega{t}≤\\Psi_2\\\\\\\\ 0,&\\Psi_2<\\omega{t}≤\\pi \\end{cases} y(t)=⎩⎨⎧​0,M,0,​0≤ωt<Ψ1​Ψ1​≤ωt≤Ψ2​Ψ2​<ωt≤π​
区间端点：
Ψ1=arcsin⁡hA,Ψ2=π−arcsin⁡mhA\\Psi_1=\\arcsin\\frac{h}{A},\\Psi_2=\\pi-\\arcsin\\frac{mh}{A} Ψ1​=arcsinAh​,Ψ2​=π−arcsinAmh​
$y(t)$为奇对称函数，而非奇函数，则有：
A1=2π∫0πy(t)cos⁡ωtd(ωt)=2π∫Ψ1Ψ2Mcos⁡ωtd(ωt)=2MhπA(m−1)B1=2π∫0πy(t)sin⁡ωtd(ωt)=2π∫Ψ1Ψ2Msin⁡ωtd(ωt)=2Mπ[1−(mhA)2+1−(hA)2]\\begin{aligned} &A_1=\\frac{2}{\\pi}\\int_0^{\\pi}y(t)\\cos\\omega{t}{\\rm d}(\\omega{t})=\\frac{2}{\\pi}\\int_{\\Psi_1}^{\\Psi_2}M\\cos\\omega{t}{\\rm d}(\\omega{t})=\\frac{2Mh}{\\pi{A}}(m-1)\\\\\\\\ &B_1=\\frac{2}{\\pi}\\int_0^{\\pi}y(t)\\sin\\omega{t}{\\rm d}(\\omega{t})=\\frac{2}{\\pi}\\int_{\\Psi_1}^{\\Psi_2}M\\sin\\omega{t}{\\rm d}(\\omega{t})=\\frac{2M}{\\pi}\\left[\\sqrt{1-\\left(\\frac{mh}{A}\\right)^2}+\\sqrt{1-\\left(\\frac{h}{A}\\right)^2}\\right] \\end{aligned} ​A1​=π2​∫0π​y(t)cosωtd(ωt)=π2​∫Ψ1​Ψ2​​Mcosωtd(ωt)=πA2Mh​(m−1)B1​=π2​∫0π​y(t)sinωtd(ωt)=π2​∫Ψ1​Ψ2​​Msinωtd(ωt)=π2M​​1−(Amh​)2​+1−(Ah​)2​​​
因此，死区滞环继电特性描述函数为：
N(A)=2MπA[1−(mhA)2+1−(hA)2]+j2MhπA2(m−1),A≥hN(A)=\\frac{2M}{\\pi{A}}\\left[\\sqrt{1-\\left(\\frac{mh}{A}\\right)^2}+\\sqrt{1-\\left(\\frac{h}{A}\\right)^2}\\right]+{\\rm j}\\frac{2Mh}{\\pi{A^2}}(m-1),A≥h N(A)=πA2M​​1−(Amh​)2​+1−(Ah​)2​​+jπA22Mh​(m−1),A≥h
取$h=0$，理想继电特性描述函数：
N(A)=4MπAN(A)=\\frac{4M}{\\pi{A}} N(A)=πA4M​
取$m=1$，死区继电特性描述函数：
N(A)=4MπA1−(hA)2,A≥hN(A)=\\frac{4M}{\\pi{A}}\\sqrt{1-\\left(\\frac{h}{A}\\right)^2},A≥h N(A)=πA4M​1−(Ah​)2​,A≥h
取$m=-1$，滞环继电特性描述函数：
N(A)=4MπA1−(hA)2−j4MhπA2,A≥hN(A)=\\frac{4M}{\\pi{A}}\\sqrt{1-\\left(\\frac{h}{A}\\right)^2}-{\\rm j}\\frac{4Mh}{\\pi{A^2}},A≥h N(A)=πA4M​1−(Ah​)2​−jπA24Mh​,A≥h

3.2.3 典型非线性特性的描述函数

1. 理想继电特性(库仑摩擦)和有死区的继电特性

   

   • 理想继电特性：
     N(A)=4MπAN(A)=\\frac{4M}{\\pi{A}} N(A)=πA4M​

   • 有死区的继电特性：
     N(A)=4MπA1−(hA)2，A≥hN(A)=\\frac{4M}{\\pi{A}}\\sqrt{1-\\left(\\frac{h}{A}\\right)^2}，A≥h N(A)=πA4M​1−(Ah​)2​，A≥h

2. 有滞环的继电特性和有死区与滞环的继电特性

   

   • 有滞环的继电特性：
     N(A)=4MπA1−(hA)2−j4MhπA2N(A)=\\frac{4M}{\\pi{A}}\\sqrt{1-\\left(\\frac{h}{A}\\right)^2}-{\\rm j}\\frac{4Mh}{\\pi{A^2}} N(A)=πA4M​1−(Ah​)2​−jπA24Mh​

   • 有死区与滞环的继电特性：
     N(A)=2MπA[1−(mhA)2+1−(hA)2]+j2MhπA2(m−1),A≥hN(A)=\\frac{2M}{\\pi{A}}\\left[\\sqrt{1-\\left(\\frac{mh}{A}\\right)^2}+\\sqrt{1-\\left(\\frac{h}{A}\\right)^2}\\right]+{\\rm j}\\frac{2Mh}{\\pi{A^2}}(m-1),A≥h N(A)=πA2M​​1−(Amh​)2​+1−(Ah​)2​​+jπA22Mh​(m−1),A≥h

3. 饱和特性(幅值限制)和有死区的饱和特性

   

   • 饱和特性：
     N(A)=2Kπ[arcsin⁡aA+aA1−(aA)2],A≥aN(A)=\\frac{2K}{\\pi}\\left[\\arcsin\\frac{a}{A}+\\frac{a}{A}\\sqrt{1-\\left(\\frac{a}{A}\\right)^2}\\right],A≥a N(A)=π2K​[arcsinAa​+Aa​1−(Aa​)2​],A≥a

   • 有死区的饱和特性：
     N(A)=2Kπ[arcsin⁡aA−arcsin⁡ΔA+aA1−(aA)2−ΔA1−(ΔA)2],A≥aN(A)=\\frac{2K}{\\pi}\\left[\\arcsin\\frac{a}{A}-\\arcsin\\frac{\\Delta}{A}+\\frac{a}{A}\\sqrt{1-\\left(\\frac{a}{A}\\right)^2}-\\frac{\\Delta}{A}\\sqrt{1-\\left(\\frac{\\Delta}{A}\\right)^2}\\right],A≥a N(A)=π2K​​arcsinAa​−arcsinAΔ​+Aa​1−(Aa​)2​−AΔ​1−(AΔ​)2​​,A≥a

4. 死区特性和间隙特性

   

   • 死区特性：
     N(A)=2Kπ[π2−arcsin⁡ΔA−ΔA1−(ΔA)2],A≥ΔN(A)=\\frac{2K}{\\pi}\\left[\\frac{\\pi}{2}-\\arcsin\\frac{\\Delta}{A}-\\frac{\\Delta}{A}\\sqrt{1-\\left(\\frac{\\Delta}{A}\\right)^2}\\right],A≥\\Delta N(A)=π2K​​2π​−arcsinAΔ​−AΔ​1−(AΔ​)2​​,A≥Δ

   • 间隙特性：
     N(A)=Kπ[π2+arcsin⁡(1−2bA)+2(1−2bA)bA(1−bA)]+j4KbπA(bA−1),A≥bN(A)=\\frac{K}{\\pi}\\left[\\frac{\\pi}{2}+\\arcsin\\left(1-\\frac{2b}{A}\\right)+2\\left(1-\\frac{2b}{A}\\right)\\sqrt{\\frac{b}{A}\\left(1-\\frac{b}{A}\\right)}\\right]+{\\rm j}\\frac{4Kb}{\\pi{A}}\\left(\\frac{b}{A}-1\\right),A≥b N(A)=πK​[2π​+arcsin(1−A2b​)+2(1−A2b​)Ab​(1−Ab​)​]+jπA4Kb​(Ab​−1),A≥b

5. 变增益特性、有死区的线性特性、库仑摩擦加黏性摩擦

   

   • 变增益特性：
     N(A)=K2+2(K1−K2)π[arcsin⁡sA+sA1−(sA)2],A≥sN(A)=K_2+\\frac{2(K_1-K_2)}{\\pi}\\left[\\arcsin\\frac{s}{A}+\\frac{s}{A}\\sqrt{1-\\left(\\frac{s}{A}\\right)^2}\\right],A≥s N(A)=K2​+π2(K1​−K2​)​[arcsinAs​+As​1−(As​)2​],A≥s

   • 有死区的线性特性：
     N(A)=K−2Kπarcsin⁡ΔA+4M−2KΔπA1−(ΔA)2,A≥ΔN(A)=K-\\frac{2K}{\\pi}\\arcsin\\frac{\\Delta}{A}+\\frac{4M-2K\\Delta}{\\pi{A}}\\sqrt{1-\\left(\\frac{\\Delta}{A}\\right)^2},A≥\\Delta N(A)=K−π2K​arcsinAΔ​+πA4M−2KΔ​1−(AΔ​)2​,A≥Δ

   • 库仑摩擦加黏性摩擦：
     N(A)=K+4MπAN(A)=K+\\frac{4M}{\\pi{A}} N(A)=K+πA4M​

3.3 非线性系统简化

3.3.1 非线性特性的并联

若两个非线性特性输入相同，输出相加、减，则等效非线性为两个非线性特性的叠加，并联等效非线性特性的描述函数为各非线性描述函数的代数和；

3.3.2 非线性特性的串联

当两个非线性环节串联时，采用图解法简化。

参数确定：
Δ2=K1(Δ−Δ1)⇒Δ=Δ1+Δ2K1\\Delta_2=K_1(\\Delta-\\Delta_1)\\Rightarrow\\Delta=\\Delta_1+\\frac{\\Delta_2}{K_1} Δ2​=K1​(Δ−Δ1​)⇒Δ=Δ1​+K1​Δ2​​

a2=K1(a−Δ1)⇒a=a2K1+Δ1a_2=K_1(a-\\Delta_1)\\Rightarrow{a}=\\frac{a_2}{K_1}+\\Delta_1 a2​=K1​(a−Δ1​)⇒a=K1​a2​​+Δ1​

当$|x|\le \Delta$时，y(x1)y(x_1)y(x1​)特性可知，$y(x)=0$；

当$|x|\ge a$时，y(x1)y(x_1)y(x1​)特性可知，y(x)=K2(a2−Δ2)y(x)=K_2(a_2-\\Delta_2)y(x)=K2​(a2​−Δ2​)；

当$\Delta <|x|<a$时，y(x1)y(x_1)y(x1​)位于线性区，$y(x)$呈线性，设斜率为$K$，有：
y(x)=K(x−Δ)=K2(x1−Δ2)y(x)=K(x-\\Delta)=K_2(x_1-\\Delta_2) y(x)=K(x−Δ)=K2​(x1​−Δ2​)
当$x=a$时，x1=a2x_1=a_2x1​=a2​，由于：x1=Δ2+K1(a−Δ)x_1=\\Delta_2+K_1(a-\\Delta)x1​=Δ2​+K1​(a−Δ)，因此：a−Δ=a2−Δ2K1a-\\Delta=\\displaystyle\\frac{a_2-\\Delta_2}{K_1}a−Δ=K1​a2​−Δ2​​，因此有：K=K1K2K=K_1K_2K=K1​K2​；

两个非线性环节串联，等效特性还取决于前后次序；

3.3.3 线性部分的等效变换

3.4 非线性系统稳定性分析的描述函数法

3.4.1 变增益线性系统的稳定性分析

设$G(s)$的极点均位于$s$的左半平面，即$P=0,G(j\omega )$的奈奎斯特曲线ΓG\\Gamma_GΓG​如上图所示，闭环系统特征方程为：
1+KG(jω)=0或G(jω)=−1K+j01+KG({\\rm j}\\omega)=0或G({\\rm j}\\omega)=-\\frac{1}{K}+{\\rm j}0 1+KG(jω)=0或G(jω)=−K1​+j0
当ΓG\\Gamma_GΓG​曲线不包围点(−1K,j0)(-\\displaystyle\\frac{1}{K},{\\rm j}0)(−K1​,j0)时，即：$Z=P-2N=-2N=0$，系统闭环稳定；

当ΓG\\Gamma_GΓG​曲线包围点(−1K,j0)(-\\displaystyle\\frac{1}{K},{\\rm j}0)(−K1​,j0)时，系统不稳定；

当ΓG\\Gamma_GΓG​曲线穿越点(−1K,j0)(-\\displaystyle\\frac{1}{K},{\\rm j}0)(−K1​,j0)时，系统临界稳定，将产生等幅振荡；

若设$K$在一定范围内可变，即有：K1≤K≤K2K_1≤K≤K_2K1​≤K≤K2​，则(−1K,j0)(-\\displaystyle\\frac{1}{K},{\\rm j}0)(−K1​,j0)为复平面实轴上的一段直线，若ΓG\\Gamma_GΓG​曲线不包围该直线，则系统闭环稳定，当ΓG\\Gamma_GΓG​包围该直线时，则系统闭环不稳定；

3.4.2 应用描述函数分析非线性系统的稳定性

当非线性特性采用描述函数近似等效时，闭环系统的特征方程为：
1+N(A)G(jω)=0⇒G(jω)=−1N(A)1+N(A)G({\\rm j}\\omega)=0\\Rightarrow{G({\\rm j}\\omega)}=-\\frac{1}{N(A)} 1+N(A)G(jω)=0⇒G(jω)=−N(A)1​
称−1N(A)-\\displaystyle\\frac{1}{N(A)}−N(A)1​为非线性环节的负倒描述函数；

在复平面上绘制ΓG\\Gamma_GΓG​曲线和−1N(A)-\\displaystyle\\frac{1}{N(A)}−N(A)1​曲线时，−1N(A)-\\displaystyle\\frac{1}{N(A)}−N(A)1​曲线上箭头表示随$A$增大，−1N(A)-\\displaystyle\\frac{1}{N(A)}−N(A)1​的变化方向；若ΓG\\Gamma_GΓG​曲线和−1N(A)-\\displaystyle\\frac{1}{N(A)}−N(A)1​曲线无交点，表明上式无$\omega$的正实数解；

非线性系统的稳定性判据：若ΓG\\Gamma_GΓG​曲线不包围−1N(A)-\\displaystyle\\frac{1}{N(A)}−N(A)1​曲线，则非线性系统稳定；若ΓG\\Gamma_GΓG​曲线包围−1N(A)-\\displaystyle\\frac{1}{N(A)}−N(A)1​曲线，则非线性系统不稳定；