-
摘要: 研究传感器未建模动态对Buck变换器滑模控制系统的性能影响, 提出一种基于奇异摄动理论的稳定性和输出电压谐波分析的新方法.给出滑模控制器的参数整定方法, 选取传感器的上升时间作为摄动时间, 建立其未建模动态的奇异摄动模型, 在多时间尺度框架下, 揭示传感器稳定输出与摄动时间的影响关系.在此基础上, 构造一个类Lyapunov函数分析未建模动态对整个闭环控制系统的稳定性影响, 证明未建模动态诱发谐波的必然性.针对输出电压的谐波, 在频域内利用描述函数法推导出未建模动态摄动时间与其谐波幅值和频率的数学影响关系.仿真结果验证所提方法的正确性和有效性.Abstract: This paper investigates the influence of unmodeled dynamics from the sensor in a sliding mode controlled (SMC) Buck converter, and presents a quantitative analysis approach for the stability and output voltage harmonics on the basis of singular perturbation theory. The design parameter of SMC is first given. Taking the rise time of sensor as a perturbed parameter, the buck converter is modeled as a singularly perturbed system. Then in the frame of multiple-time scale, the influence relationship of stable sensor and the time constant is deduced, and further the stability of the whole closed-loop buck converter system is given by constructing a Lyapunov-like function, proving the existence of harmonics induced by unmodeled dynamics. The describing function method is introduced to analyze the frequency and amplitude of the output voltage harmonics. Simulations validate the proposed method.
-
Key words:
- Sliding mode control /
- unmodeled dynamics /
- sensor /
- singular perturbation /
- describing function method
-
现代防空反导体系已具备全天域、多层次、多维度防御能力; 在此背景下, 传统弹道弹和巡航弹已无法满足复杂对抗态势下的作战需求. 因此人们开始研究具有多空域机动能力的临近空间高超声速滑翔飞行器. 美国于2003年便提出了HTV-2 (Hypersonic technology vehicle 2)[1-2]计划, 并于2017年完成了高超声速滑翔试验[3]. 俄罗斯于2015年和2017年分别完成了高超声速飞行器“Yu-71”和“先锋”的飞行试验[4]. 日本于2019年起加紧研制高速滑翔飞行器[5]; 同年, 印度也完成了高超技术验证机首飞[6].
基于出色的突防能力, 多高超声速飞行器协同作战模式[7]同样备受关注. 但同时规划多条滑翔段轨迹具有一定难度: 1)约束多, 如载荷、热等过程约束以及位置、高度等终端约束; 2)耦合强, 如气动−约束−轨迹−指标; 3)干扰多、速度快、航程远. 综上, 在线规划方法必须具备很强的场景适应性, 同时还要精确控制各飞行器间的相对运动关系.
许多学者针对高超声速飞行器滑翔段轨迹规划问题展开了研究. 孙长银等[8]探讨了高超声速飞行器领域的新方法和新问题, 提出了一种通用的大包络飞行动力学模型; Liao等[9]基于有限时间非线性观测器实现了滑翔段轨迹的在线规划. 在协同规划方面, 王芳等[10]提出了一种时间最优导弹编队方法, Vincent等[11]和Chen等[12]分别用遗传算法和粒子群优化(Particle swarm optimization, PSO)算法实现了无人机编队控制, Xu等[13]提出了一种用于故障条件下多巡航弹协同控制的时变容错算法.
现有方法有效解决了滑翔段轨迹设计问题, 但协同规划方面的研究多以巡航弹/无人机为研究对象. 针对上述问题, 从提高在线规划效率出发, 本文提出一种自动满足终端约束的滑翔段飞行剖面. 推导滑翔段解析解, 实现了飞行剖面的快速重构. 提出一种基于惯性权重智能选择的改进PSO算法, 用于完成在线协同规划. 最后通过数学仿真验证了方法的正确性和有效性.
1. 飞行器协同轨迹规划问题建模
1.1 滑翔段运动数学模型
高超声速飞行器滑翔段运动方程为
$$\left\{ \begin{aligned} & \dot V = - \dfrac{D}{m} - g\sin \gamma \\ & \dot \gamma = \dfrac{{L\cos \sigma }}{{mV}} - \left( {\dfrac{g}{V} - \frac{V}{r}} \right)\cos \gamma \\ & \dot \psi = \dfrac{{L\sin \sigma }}{{mV\cos \gamma }} + \dfrac{V}{r}\cos \gamma \sin \psi \tan \phi \\ &\dot r = V\sin \gamma \\ &\dot \theta = \dfrac{{V\cos \gamma \sin \psi }}{{r\cos \phi }} \\ & \dot \phi = \dfrac{{V\cos \gamma \cos \psi }}{r} \end{aligned} \right.$$ (1) 式中,
$V$ 为速度,$\gamma $ 为飞行路径角,$\psi $ 为航向角,$r$ 为飞行器质心到地心的距离,$\theta $ 为经度,$\phi $ 为纬度;$m$ 为质量,$\sigma $ 为倾侧角,$g$ 为地球引力加速度.$D = {c_D}q{S_c}$ 为气动阻力,$L = {c_L}q{S_c}$ 为气动升力, 其中${c_D}$ 为阻力系数,${c_L}$ 为升力系数,${S_c}$ 为特征面积,$q$ 为飞行动压.$q=\rho {V^2}/2,$ 其中$\rho $ 为大气密度, 采用指数函数计算$$\rho = {\rho _0}{{\rm{e}} ^{ - \beta h}}$$ (2) 式中,
$h = r - {R_e}$ 为飞行高度$({R_e}$ 为地球平均半径),$\beta = 1.4064 \times {10^{ - 4}}\;{{\rm{m}}^{ - 1}}$ 且${\rho _0}=1.225\;{\rm{kg/}}{{\rm{m}}^{\rm{3}}}$ .设阻力系数按如下规律变化:
$${c_D} = {c_{D0}}{V^{{\beta _1}}}{{\rm{e}} ^{ - {\beta _2}h}}$$ (3) 式中,
$\alpha $ 为飞行攻角,${c_{D0}},{\beta _1}$ 和${\beta _2}$ 为给定参数. 取${c_D}_0=0.156,{\beta _1}=2.6 \times {10^{ - 4}},{\beta _2} = 1.5 \times {10^{ - 6}}$ .定义从当前位置到终端位置间的航程为剩余航程, 记为
${R_L}$ , 则有$$\left\{ \begin{aligned} & {R_L} = {\mu _L}{R_e} \\ & {{\dot R}_L} = - V\cos \gamma \cos \zeta \\ & {\mu _L} = \arccos (\sin {\phi _T}\sin \phi + \cos {\phi _T}\cos \phi \cos \Delta \theta ) \end{aligned} \right.$$ (4) 式中,
$\zeta = {\psi _w} - \psi $ 为航向偏差,$\Delta \theta = {\theta _T} - \theta ,{\theta _T}$ 和${\phi _T}$ 为目标经纬度;${\psi _w}$ 为期望航向角$${\psi _w} = \arccos \left( {\dfrac{{\sin {\phi _T} - \sin \phi \cos {\mu _L}}}{{\cos \phi \sin {\mu _L}}}} \right)$$ (5) 设滑翔段约束条件为
$$\left\{ \begin{aligned} & {{\boldsymbol{c}}_{{\rm{end}}}} = {\left[ {{h_f} - {h_{{\rm{ter}}}},{\gamma _f} - {\gamma _{{\rm{ter}}}},{R_{Lf}},{\zeta _f}} \right]^{\rm{T}}} = {\bf{0}} \\ & {{\boldsymbol{c}}_{{\rm{path}}}} = {\left[ {\dfrac{{{N_m}}}{{{N_{\max }}}},\dfrac{q}{{{q_{\max }}}},\dfrac{{\left| \alpha \right|}}{{{\alpha _{\max }}}},\dfrac{{\left| \sigma \right|}}{{{\sigma _{\max }}}}} \right]^{\rm{T}}} < {\bf{1}} \end{aligned} \right.$$ (6) 式中,
${{\boldsymbol{c}}_{{\rm{end}}}}$ 为终端约束,${{\boldsymbol{c}}_{{\rm{path}}}}$ 为过程约束. 下标“$f$ ”表示实际终端状态, 下标“ter”表示期望终端状态. 下标“max”表示该过程变量允许的最大值,${N_m}$ 为法向过载的绝对值.1.2 多飞行器协同轨迹规划模型
设协同飞行器数量为
${N_{{\rm{hyp}}}}.$ 协同规划应考虑各个飞行器的总体性能以及各自飞行约束, 记为${\boldsymbol{c}}_{{\rm{end}}}^{\left( 1 \right)}, \cdots ,{\boldsymbol{c}}_{{\rm{end}}}^{\left( {{N_{{\rm{hyp}}}}} \right)},{\boldsymbol{c}}_{{\rm{path}}}^{\left( 1 \right)}, \cdots ,{\boldsymbol{c}}_{{\rm{path}}}^{\left( {{N_{{\rm{hyp}}}}} \right)}.$ 考虑到各飞行器的初始运动状态各异, 设其剩余飞行时间为$t_{{\rm{go}}}^{(i)},$ 协同规划指标为$${J_{{\rm{syn}}}} = \max \left( {t_{{\rm{go}}}^{\left( i \right)}} \right) - \min \left( {t_{{\rm{go}}}^{\left( i \right)}} \right),\;\;\;i = 1,2, \cdots ,{N_{{\rm{hyp}}}}$$ (7) 2. 滑翔段飞行剖面设计
设飞行高度为剩余航程的函数为
$$h = f( {{R_L}}) = - \frac{1}{{\bar \beta }}\ln F( {{R_L}} ) = - \frac{1}{{\bar \beta }}\sum\limits_{i = 0}^n {\left( {{a_i}R_L^i} \right)} $$ (8) 式中,
${a_i}$ 为未知系数,$\bar \beta = \beta + {\beta _2},F({R_L})$ 为高阶多项式. 结合式(1)和式(8), 可得$$ \frac{{{\rm{d}}h}}{{{\rm{d}}{R_L}}} = - \frac{{F'( {{R_L}} )}}{{\bar \beta F( {{R_L}} )}} = - \frac{{\tan \gamma }}{{\cos \zeta }} $$ (9) $$\frac{{{{\rm{d}}^2}h}}{{{\rm{d}}R_L^2}} = \frac{{{{\left[ {F'( {{R_L}} )} \right]}^2} - F''( {{R_L}} )F( {{R_L}} )}}{{\bar \beta {F^2}( {{R_L}} )}} = - \frac{{{\rm{d}}\left( {\frac{{\tan \gamma }}{{\cos \zeta }}} \right)}}{{{\rm{d}}{R_L}}}$$ (10) 式中,
$$ F'({R_L})=\displaystyle\sum\limits_{i = 1}^5 {i{a_i}R_L^{i - 1}},F''({R_L})=\displaystyle\sum\limits_{i = 2}^5 i\left( {i - 1} \right){a_i}R_L^{i - 2} $$ 定义
${K_h} = \dfrac{{\rm{d}}^2h}{{\rm{d}}R_L^2}$ , 本文认为${\dot \psi _w} = 0,\dot \zeta = - \dot \psi.$ 由式(11)可得升力为$$L = \dfrac{{{K_h}{{\cos }^2}\gamma {{\cos }^2}\zeta + \frac{g}{{{V^2}}} - \frac{1}{r} + \frac{{{K_b}}}{r}}}{\frac{\cos \sigma - \sin \zeta \sin \gamma \sin \sigma }{m{V^2}\cos \gamma }}$$ (11) 式中,
$ {K_b} = 2\sin \gamma \cos \gamma \sin \psi \tan \phi \sin \zeta$ .设
$n = 5,$ 则确定$F({R_L})$ 需建立6个方程. 由于初始和终端高度已知, 由式(8)可得$$\left\{ \begin{aligned} & {{\rm{e}} ^{ - \bar \beta {h_0}}} = \displaystyle\sum\limits_{i = 0}^5 {{a_i}R_{L0}^i} \\ & {{\rm{e}} ^{ - \bar \beta {h_f}}} = \displaystyle\sum\limits_{i = 0}^5 {{a_i}R_{Lf}^i} \end{aligned} \right.$$ (12) 式中, 下标“0”表示滑翔段初始运动状态.
由于初始和终端飞行路径角已知, 认为
$\cos \zeta \approx 1,$ 由式(9)可得$$\left\{ \begin{aligned} & \displaystyle\sum\limits_{i = 1}^5 {{a_i}\left( {R_{L0}^i - \dfrac{{iR_{L0}^{i - 1}}}{{\bar \beta \tan {\gamma _0}}}} \right) + {a_0} = 0} \\ & \displaystyle\sum\limits_{i = 1}^5 {{a_i}\left( {R_{Lf}^i - \dfrac{{iR_{Lf}^{i - 1}}}{{\bar \beta \tan {\gamma _f}}}} \right) + {a_0} = 0} \end{aligned} \right.$$ (13) 设
${h_1},{h_2}$ 分别为剩余航程2/3处(记为${R_L}_1)$ 和1/3处(记为${R_L}_2)$ 的高度值, 由式(8)可得:$$\left\{ \begin{aligned} & {{\rm{e}} ^{ - \bar \beta {h_1}}} = \displaystyle\sum\limits_{i = 0}^5 {\left( {{a_i}R_{L1}^i} \right)} \\ & {{\rm{e}} ^{ - \bar \beta {h_2}}} = \displaystyle\sum\limits_{i = 0}^5 {\left( {{a_i}R_{L2}^i} \right)} \end{aligned} \right.$$ (14) 定义剖面设计变量
${{\boldsymbol{u}}_h} = {[{h_{{\kern 1pt} 1}},{h_2}]^{\rm{T}}}$ , 则式(7)对应的优化变量分别为${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 1 \right)} = [h_1^{\left( 1 \right)}, \cdots ,h_1^{({N_{{\rm{hyp}}}})}]$ 和${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 2 \right)} = [h_2^{\left( 1 \right)}, \cdots ,h_2^{({N_{{\rm{hyp}}}})}]$ , 共$2{N_{{\rm{hyp}}}}$ 个.最后, 设计倾侧角变化规律为
$$\sigma =\frac{ \zeta {\sigma _{\max }}}{{\zeta _{\max }}}$$ (15) 式中,
${\zeta _{\max }}$ 为给定正数.3. 滑翔段高精度解析解推导
滑翔段高度的解析解即式(8). 考虑到滑翔段中航向偏差
$\zeta $ 和飞行路径角$\gamma $ 较小, 由式(1)和式(9)可得飞行路径角的解析解为$$\gamma = - f'\left( {{R_L}} \right) = \dfrac{{\displaystyle\sum\limits_{i = 1}^5 {\left( {{a_i}R_{L0}^{i - 1}} \right)} }}{{\bar \beta \displaystyle\sum\limits_{i = 0}^5 {\left( {{a_i}R_{L0}^i} \right)} }}$$ (16) 考虑到滑翔段中阻力通常远大于引力沿速度方向的分量, 根据式(1)可得
$$\frac{{{\rm{d}}V}}{{{\rm{d}}h}} = - \frac{{{c_{D0}}{S_c}{\rho _0}}}{{2m}}\frac{{{V^{n + 1}}{{\rm{e}} ^{ - \bar \beta h}}}}{\gamma }$$ (17) 定义
$ {\tau _D} = {c_{D0}}{S_c}{\rho _0}/(2m)$ , 因此$$ - {V^{ - \left( {{\beta _1} + 1} \right)}}{\rm{d}}V = {\tau _D}F\left( {{R_L}} \right){\rm{d}}{R_L}$$ (18) 对上式两端积分, 得
$${V^{ - {\beta _1}}} = {\tau _D}\int_0^{{R_L}} {F\left( {{R_L}} \right){\kern 1pt} } {\rm{d}}{R_L} + {C_V}$$ (19) 式中,
${C_V}$ 为积分常数. 定义$$y\left( {{R_L}} \right) = \int_0^{{R_L}} {F\left( {{R_L}} \right){\kern 1pt} } {\rm{d}}{R_L}=\sum\limits_{i = 0}^n {\frac{{{a_i}R_L^{i + 1}}}{{i + 1}}} $$ (20) 故
${C_V} = {V_0}^{ - {\beta _1}} - {\tau _D}y({R_{L0}}),$ 速度解析解为$$V = {\left[ {{\tau _D}y\left( {{R_L}} \right) + {C_V}} \right]^{ - \frac{1}{{{\beta _1}}}}}$$ (21) 根据式(21), 动压的解析表达式为
$$q = 0.5{\rho _0}{{\rm{e}} ^{ - \beta h}}{\left[ {{\tau _D}y\left( {{R_L}} \right) + {C_V}} \right]^{ - \frac{2}{{{\beta _1}}}}}$$ (22) 令
$\cos \zeta = \cos \gamma = 1,$ 则升力的解析解为$$L = m\frac{{{K_h}{V^2} + g}}{{\cos \sigma }}$$ (23) 由式(23), 可以解析估计出法向过载
$${N_m} \le \frac{{g + f''\left( {{R_L}} \right){{\left[ {{\tau _D}y\left( {{R_L}} \right) + {C_V}} \right]}^{ - \frac{2}{{{\beta _1}}}}}}}{{g\cos {\sigma _{\max }}}}$$ (24) 最后, 由式(4)可得
$$ - \frac{{{\rm{d}}t}}{{{\rm{d}}{R_L}}} = \frac{1}{{V\cos \gamma \cos \zeta }} = {\left[ {{\tau _D}y\left( {{R_L}} \right) + {C_V}} \right]^{\frac{1}{{{\beta _1}}}}}$$ (25) 定义积分常数
${C_t},$ 考虑到$t=0$ 时${R_L}={R_L}_0,$ 且$t={t_{{\rm{go}}}}$ 时${R_{Lf}}=0;$ 结合式(20), 可得$${t_{{\rm{go}}}} = {\tau _D}\left[ {\sum\limits_{i = 0}^n {\frac{{{a_i}R_{L0}^{i + 2}}}{{\left( {i + 1} \right)\left( {i + 2} \right)}} + } \frac{{{C_V}}}{{{\tau _D}}}{R_{L0}}} \right]$$ (26) 4. 滑翔段轨迹协同规划算法
4.1 基本PSO算法
PSO算法用
${M_{{\rm{pso}}}}$ 个粒子寻找${D_{{\rm{pso}}}}$ 个参数. 设粒子$i$ 在空间中的位置和速度分别为${{\boldsymbol{x}}_i} = {[{x_{i,1}}, \cdots ,{x_{i,{D_{{\rm{pso}}}}}}]^{\rm{T}}},$ ${{\boldsymbol{v}}_i} = {[{v_{i,1}}, \cdots ,{v_{i,{D_{{\rm{pso}}}}}}]^{\rm{T}}},$ 其历史最佳位置为${{\boldsymbol{p}}_i} = [{p_{i,1}},$ $\cdots ,{p_{i,{D_{{\rm{pso}}}}}}]^{\rm{T}};$ 整个群体的历史最佳位置为${{\boldsymbol{p}}_g} = [{p_{g,1}},$ $\cdots ,{p_{g,{D_{{\rm{pso}}}}}}]^{\rm{T}},$ 相应性能指标为${g_{{\rm{best}}}}.$ PSO算法为[14-15]$$\begin{aligned}[b]v_{i,j}^{\left( {\ell + 1} \right)}=\;&wv_{i,j}^{\left( \ell \right)}+{c_1}{r_1}\left( {{p_{i,j}} - x_{i,j}^{\left( \ell \right)}} \right)+\\ &{c_2}{r_2}\left( {{p_{g,j}} - x_{i,j}^{\left( \ell \right)}} \right)\end{aligned}$$ (27) $$x_{i,j}^{\left( {\ell + 1} \right)} = x_{i,j}^{\left( \ell \right)} + v_{i,j}^{\left( {\ell + 1} \right)}$$ (28) 式中,
$j = 1,2, \cdots ,{D_{{\rm{pso}}}},\ell $ 为进化代数,$w \in [0,1]$ 为惯性权重,${c_1}={c_2}=2,{r_1},{r_2} \in [0,1]$ 由每代随机产生. 令${\phi _1} = {c_1}{r_1} + {c_2}{r_2},{\phi _2} = {c_1}{r_1}{p_{i,j}} + {c_2}{r_2}{p_{g,j}},$ 则粒子进化公式可改写为$$x_{i,j}^{\left( {\ell + 2} \right)} + \left( {{\phi _1} - w - 1} \right)x_{i,j}^{\left( {\ell + 1} \right)} + wx_{i,j}^{\left( \ell \right)} = {\phi _2}$$ (29) 式(29)避免了对粒子速度进行初始化及更新, 算法结构更简洁, 利于提高计算效率.
调整权重是控制粒子寻优过程的常用方法. 借助强化学习方法, 通过快速智能选择惯性权重实现用较少的粒子数和迭代次数找出满足协同规划参数(
${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 1 \right)}$ 和${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 2 \right)}$ ).4.2 最优权重在线生成方法
设惯性权重为进化代数的分段函数如下:
$$w\left( \ell \right) = \left\{ {\begin{array}{*{20}{l}} &{{w_0} + \dfrac{{{w_1} - {w_0}}}{{\frac{\ell _{\max }}{\rm{3}}}}\ell },&{0 \le \ell < \dfrac{{{\ell _{\max }}}}{3}}\\ &{{w_1} + \dfrac{{{w_2} - {w_1}}}{{\frac{\ell _{\max }}{{3}}}}\ell },&{\dfrac{{{\ell _{\max }}}}{3} \le \ell < \dfrac{{2{\ell _{\max }}}}{3}}\\ &{{w_2} + \dfrac{{{w_3} - {w_2}}}{{\frac{\ell _{\max }}{{3}}}}\ell },&{\dfrac{{2{\ell _{\max }}}}{3} \le \ell < {\ell _{\max }}} \end{array}} \right.$$ (30) 式中,
${w_0} \sim {w_3}$ 为未知量,${\ell _{\max }}$ 为最大进化代数.定义自变量
${{\boldsymbol{u}}_w} = {[{w_0},{w_{{\kern 1pt} 1}},{w_2},{w_3}]^{\rm{T}}}$ 以及动作${\boldsymbol{d}} = $ ${[{d_1},{d_2},{d_3},{d_4}]^{\rm{T}}},$ 其中${d_1} \sim {d_4} \in \left[ {0,1} \right]$ 为连续值, 分别对应${w_0} \sim {w_3}$ 的更新:$$w_j^{} = \left( {2{d_j} - 1} \right){n_j}w_j^{}$$ (31) 式中,
$j = 1,2,3,4;{n_1} \sim {n_4}$ 为预置系数.设状态值
${\boldsymbol{s}} = {{\rm{[}}{p_{g,1}}{\rm{,}} \cdots {\rm{,}}{\kern 1pt} {p_{g,{D_{{\rm{pso}}}}}}{\rm{]}}^{\rm{T}}}$ 以及观测值${{\boldsymbol{o}}_t} = $ $ {{\rm{[}}t_{{\rm{go}}}^{{\rm{(a)}}},t_{{\rm{go}}}^{{\rm{(s)}}},t_{{\rm{go}}}^{{\rm{(l)}}}{\rm{]}}^{\rm{T}}},$ 其中$t_{{\rm{go}}}^{{\rm{(a)}}},t_{{\rm{go}}}^{{\rm{(s)}}},t_{{\rm{go}}}^{{\rm{(l)}}}$ 分别为$t_{{\rm{go}}}^{(i)}$ 的平均值、标准差以及最大/小值之差. 定义函数网络${Q_{{\rm{net}}}}$ 来表征状态${\boldsymbol{s}}$ 和动作${\boldsymbol{d}}$ 的优劣, 其中优劣性由奖励值$r = {g_{{\rm{best}}}}$ 表示. 定义策略网络${P_{{\rm{net}}}}$ 用以表征动作${\boldsymbol{d}}$ 的更新; 由式(35)可知,${P_{{\rm{net}}}}$ 即为${{\boldsymbol{u}}_w}$ 的求解器. 根据上述定义, 记${Q_{{\rm{net}}}}$ 和${P_{{\rm{net}}}}$ 的值分别为${{\boldsymbol{\theta}} _Q}$ 和${{\boldsymbol{\theta }}_P},$ 提出最优惯性权重训练方法如下:步骤 1. 随机初始化
${{\boldsymbol{\theta}} _Q}$ 和${{\boldsymbol{\theta }}_P},$ 并令${{\boldsymbol{\theta}} '_Q} = {{\boldsymbol{\theta}} _Q},{{\boldsymbol{\theta}} '_P} = {{\boldsymbol{\theta }}_P};$ 设$k = 1,$ 随机初始化${{\boldsymbol{u}}_w}$ ;步骤 2. 根据
${{\boldsymbol{u}}_w}$ 执行PSO, 得到状态${{\boldsymbol{s}}_k}$ ;步骤 3. 随机初始化扰动量
$\Delta {{\boldsymbol{d}}_k},$ 得到动作${{\boldsymbol{d}}_k}= {P_{{\rm{net}}}}({{\boldsymbol{\theta}} _P},$ ${{\boldsymbol{s}}_k}) + \Delta {{\boldsymbol{d}}_k};$ 执行动作${{\boldsymbol{d}}_k}$ , 并更新${{\boldsymbol{u}}_w},$ 运行PSO得到奖励${r_k}$ 和状态${{\boldsymbol{s}}_k};$ 步骤 4. 将
$({{\boldsymbol{s}}_k},{{\boldsymbol{d}}_k},{r_k},{{\boldsymbol{s}}_{k + 1}})$ 放入经验回放池, 令$ k =$ $ k + 1,$ 返回步骤2, 直至$k = 1\;000$ ;步骤 5. 从经验池中
${N_R}$ 次随机取样; 设定系数${k_Q},$ 计算第$i{\kern 1pt} {\kern 1pt} \;(i = 1,2, \cdots ,{N_R})$ 次取样结果$${y_i} = {r_i} + {k_Q} \times {{\boldsymbol{\theta}} '_Q}\left( {{{\boldsymbol{s}}_i},{{\boldsymbol{d}}_i},{{\boldsymbol{\theta}} _P}} \right)$$ (32) 步骤 6. 以损失函数最小为目标更新
${{\boldsymbol{\theta }}_Q}$ $$\min \frac{1}{{{N_R}}}\sum\limits_{i = 1}^{{N_R}} {{{\left( {{y_i} - {{\boldsymbol{\theta}} _Q}\left( {{{\boldsymbol{s}}_i},{{\boldsymbol{d}}_i},{{\boldsymbol{\theta}} _P}} \right)} \right)}^2}} \Rightarrow {{\boldsymbol{\theta}} _Q}$$ (33) 步骤 7. 通过梯度求解更新
${{\boldsymbol{\theta}} _P} = {{\boldsymbol{\theta}} _P} + \Delta {{\boldsymbol{\theta}} _P}$ $$\Delta {{\boldsymbol{\theta}} _P} = \sum\limits_{i = 1}^{{N_R}} {\left( {\frac{{\partial {\theta _Q}\left( {{{\boldsymbol{d}}_i},{{\boldsymbol{\theta}} _P},{{\boldsymbol{s}}_i}} \right)}}{{\partial {{\boldsymbol{d}}_i}}} \times \frac{{\partial {{\boldsymbol{\theta}} _P}}}{{\partial {{\boldsymbol{s}}_i}}}} \right)} $$ (34) 步骤 8. 给定系数
${\tau _\theta },$ 令${{\boldsymbol{\theta}} '_Q} = {\tau _\theta }{{\boldsymbol{\theta}} _{Q}} + (1 - {\tau _\theta }){{\boldsymbol{\theta}} '_Q}, {{\boldsymbol{\theta}} '_P} =$ $ {\tau _\theta }{{\boldsymbol{\theta}} '_P} + (1 - {\tau _\theta }){{\boldsymbol{\theta}} '_P}$ .4.3 在线协同轨迹规划方法
设从当前位置到起点的航程为
${R_{{\rm{cov}}}},$ 更新周期为$\Delta {R_{{\rm{cov}}}}.$ 如图1所示, 在线协同规划具体步骤如下:步骤 1. 协同飞行剖面重构
计算观测值
${{\boldsymbol{o}}_t}$ 并与上一步获得的${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 1 \right)}$ 和${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 2 \right)}$ (如果是首次执行则随机初始化${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 1 \right)}$ 和${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 2 \right)}$ )共同输入训练好的策略函数网络, 输出惯性权重; 然后用PSO计算${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 1 \right)}$ 和${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 2 \right)}$ 的值, 使指标函数式(7)达到最小.步骤 2. 个体飞行剖面重构
a) 计算
${R_{L1}},{R_{L2}}$ 及其对应高度值${h_1},{h_2};$ b) 求解系数
${a_0} \sim {a_5}$ , 更新飞行剖面.步骤 3. 计算攻角和倾侧角, 更新运动状态.
步骤 4. 计算航程
${R_{{\rm{cov}}}};$ 根据${R_{{\rm{cov}}}},$ 每5个周期$(5\Delta {R_{{\rm{cov}}}})$ 执行一次协同剖面重构, 每1个周期$(\Delta {R_{{\rm{cov}}}})$ 执行一次个体剖面重构; 未到达更新周期时, 各飞行器只执行步骤3.循环执行步骤
$1\sim 4 $ , 直至飞行器到达目标点上空. 步骤2用于控制终端高度和飞行路径角, 但未考虑过程约束; 步骤1既要保证时间协同, 又要通过优化求解来满足过程约束, 同时还负责修正解析解与实际值间的偏差.5. 仿真分析
参考文献[16-17], 设滑翔段过程约束为
${q_{\max }}=$ ${\rm{ 60\;kPa}},{N_{\max }}= 3{\rm{.0}},{\alpha _{\max }}= 25^\circ ,{\sigma _{\max }} = 50^\circ .$ 考虑对4架总体参数不同的高超声速飞行器进行协同规划; 参考X-33[18], 设其特征面积分别为2.0 m2, 1.6 m2, 1.3 m2和1.0 m2, 质量分别为950 kg, 900 kg, 850 kg和800 kg.设PSO算法的最大进化代数为15, 粒子数为20, 约束条件用罚函数法[19]处理, 当
${J_{{\rm{syn}}}}$ < 1或$\ell ={\ell _{\max }}$ 时停止进化.仿真条件见表1, 初始航向角由式(5)求出, 故初始航向偏差为零. 取
$\Delta {R_{{\rm{cov}}}}$ = 30 km, 仿真结果见图2 ~ 5. 仿真结果表明, 各轨迹均能满足终端约束条件和过程约束条件, 且攻角和倾侧角曲线十分平滑. 各飞行器到达目标的时刻分别为735 s, 734 s, 734 s和734 s.表 1 初始状态和终端约束Table 1 The initial states and the terminal constraints初始值 经纬度 (°) 高度 (km) 速度 (m/s) 飞行路径角 (°) 剩余航程 (km) 飞行器 1 E 90, N 45 56 5400 −2.0 2304 飞行器 2 E 65, N 30 54 5300 −2.0 2263 飞行器 3 E 40, N 35 52 5200 −1.0 2324 飞行器 4 E 70, N 70 50 5100 −1.0 2285 终端值 E 60, N 50 25 — −1.0 0 若不进行协同规划, 各飞行剖面设计参数均取
${h_1}$ = 42 km,${h_2}$ = 38 km, 则各飞行器的到达时刻为694 s, 703 s, 733 s和744 s, 说明协同规划方法正确有效.考虑初始状态下第1次计算协同参数
$({\boldsymbol{u}}_{{\rm{syn}}}^{\left( 1 \right)}$ 和${\boldsymbol{u}}_{{\rm{syn}}}^{\left( 2 \right)})$ , 分别用基本PSO和改进PSO求解; 基本PSO中的惯性权重取0.7, 改进PSO中惯性权重由策略网络生成; 其他参数相同. 由于PSO算法受随机数影响, 分别进行30次独立仿真. 用进化代数表征计算速度并进行比较; 表2所示的结果表明, 改进PSO在计算效率上显著高于基本PSO.表 2 基本PSO和改进PSO计算效率对比Table 2 Comparison of the computation efficiency进化代数 最大值 最小值 平均值 标准差 基本 PSO 20 8 14 3.17 改进 PSO 16 6 10 2.69 最后, 为验证协同规划方法的抗干扰性, 设置服从高斯分布的干扰因素(见表3), 并进行500次蒙特卡洛打靶仿真. 500次仿真中, 最晚与最早到达时刻之差的最大值为5.37 s, 平均值为1.57 s, 说明算法能够在干扰因素的影响下完成在线协同规划.
表 3 滑翔段干扰因素设置Table 3 Disturbances in the glide phase序号 干扰因素 $3\sigma $值 1 初始速度 (m/s) 30 2 初始飞行路径角 (°) 0.2 3 初始高度 (m) 500 4 初始航向偏差 (°) 0.2 5 初始剩余航程 (km) 50 6 气动系数 (%) 10 7 大气密度 (%) 10 以飞行器1为例, 给出打靶仿真结果如图6和图7所示. 结果表明, 各种约束条件依然能够得到满足. 在高度和飞行路径角的解析解推导过程中未用到近似假设条件, 因此可忽略它们和真实值间的误差. 对于速度和剩余时间的解析解, 相应误差会随着时间而减小, 并且剩余时间是根据实时运动状态估计的, 因此不存在累计误差.
6. 结束语
本文建立了多高超声速滑翔飞行器协同规划问题数学模型. 提出了一种新的滑翔飞行剖面, 实现了终端约束的自动满足, 并推导了滑翔段解析运动方程. 提出基于强化学习方法的改进PSO算法, 提高了在线协同规划效率. 仿真结果表明, 在初始状态和总体参数各异的情况下, 多弹时间协同偏差不超过2 s, 并且能够克服干扰因素的影响.
-
表 1 Buck变换器的电路参数
Table 1 Circuit parameters of Buck converter
电路参数 数值 电感 $L=50$ mH 电容 $C=100$ ${\rm{ \mathsf{ μ} }}$F 负载电阻 $R=10\, \Omega $ 输入电压 $E=10$ V 给定输出电压 $V_{ref}=5$ V 表 2 $\psi $取不同值时的输出电压$v_{c}$性能对比
Table 2 Different values of $\psi $ and their influence on the output voltage $v_{c}$
时间常数 谐波幅值 谐波频率 稳态误差 相对误差 $\psi\, ({\rm{ \mathsf{ μ} }}$s) (mV) (Hz) (mV) (%) 6.647 0.24 1000 0.12 0.0024 32.09 2.8 4515 1.4 0.028 291.26 224 526.31 112 2.24 623.02 842 263.2 421 8.42 -
[1] 王艳敏, 曹雨晴, 夏红伟. Buck变换器的电压电流双闭环终端滑模控制.电机与控制学报, 2016, 20(08): 92-97 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=djykzxb201608012Wang Yan-Min, Cao Yu-Qing, Xia Hong-Wei. Terminal sliding mode control for Buck converter with structure of voltage and current double closed loop. Electric Machines and Control, 2016, 20(08): 92-97 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=djykzxb201608012 [2] Krishnamurthy P, Khorrami F. A general dynamic scaling based control redesign to handle input unmodeled dynamics in uncertain nonlinear systems. IEEE Transactions on Automatic Control, 2017, 62(9): 4719-4726 [3] 杨天皓, 李健, 贾瑶, 刘腾飞, 柴天佑.虚拟未建模动态补偿驱动的双率自适应控制.自动化学报, 2017, 44(1): 299-310 doi: 10.16383/j.aas.2018.c160623Yang Tian-Hao, Li Jian, Jia Yao, Liu Teng-Fei, Chai Tian-You. Dual-rate adaptive control driven by virtual unmodeled dynamics compensation in industrial heat exchange process. Acta Automatica Sinica, 2017, 44(1): 299-310 doi: 10.16383/j.aas.2018.c160623 [4] 柴天佑, 张亚军.基于未建模动态补偿的非线性自适应切换控制方法, 自动化学报. 2011, 37(7): 773-786 doi: 10.3724/SP.J.1004.2011.00773Chai Tian-You, Zhang Ya-Jun. Nonlinear adaptive switching control method based on unmodeled dynamics compensation. Acta Automatica Sinica, 2011, 37(7): 773-786 doi: 10.3724/SP.J.1004.2011.00773 [5] Levant A. Chattering analysis. IEEE Transactions on Automatic Control, 2010, 55(6): 1380-1389 [6] Boiko I, Fridman L, Iriarte R, Pisano A, Usai E. Parameter tuning of second-order sliding mode controllers for linear plants with dynamic actuators. Automatica. 2006, 42(7): 677-683 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=a4968355d3bc45b64ef9002990e68fa8 [7] Goncalves J M, Megretski A, Dahleh M A. Global stability of relay feedback systems. IEEE Transactions on Automatic Control. 2011, 46(4): 550-562 [8] Boiko I. Oscillations and transfer properties of relay servo systems-the locus of a perturbed relay system approach. Automatica. 2005, 41(7): 677-683 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=b5b9ba81e1a5872c78312b1b261df1f4 [9] 杨晨, 程盈盈, 都海波, 王金平, 何怡刚. Buck型变换器自适应有限时间降压控制算法研究.自动化学报, 2016, 42(1): 315-320 doi: 10.16383/j.aas.2016.c150446Yang Chen, Cheng Ying-Ying, Du Hai-Bo, Wang Jin-Ping, He Yi-Gang. An adaptive finite-time control algorithm for buck converter systems. Acta Automatica Sinica, 2016, 42(1): 315-320 doi: 10.16383/j.aas.2016.c150446 [10] Lee H, Utkin I V. Chattering suppression methods in sliding mode control systems. Annual Reviews in Control, 2007, 31(5): 179-188 [11] Komurcugil H. Non-singular terminal sliding-mode control of DC-DC buck converters. Control Engineering Practice. 2013, 21(5): 321-332 [12] ACS712 Datasheet (PDF)-Allegro MicroSystem[Online], available: https://www.alldatasheet.com/datasheet-pdf/pdf/168326/ALLEGRO/ACS712.htm, August 1, 2018 期刊类型引用(3)
1. 王艳敏,张涵清,吴文谊,董志华. 基于模型解耦的Buck变换器滑模控制. 华中科技大学学报(自然科学版). 2025(01): 42-47+71 . 百度学术
2. 王艳敏,张伟琦,龙云,张涵清,冯勇. 基于过零检测的Buck变换器自适应连续滑模控制. 控制理论与应用. 2023(11): 1931-1939 . 百度学术
3. 陈维乐,都海波. 一种基于齐次系统理论的二阶离散超螺旋控制算法. 控制理论与应用. 2022(04): 761-769 . 百度学术
其他类型引用(5)
-