0 引言

越野车通常应用于抢险救灾、军民运输、资源勘探等方面，考虑到越野车经常遇到桥头、断头险路、窄巷、拥挤停车等转向困境，为此需具备良好的转向灵活性^[1-2]。而四轮毂电机相较于传统的内燃机驱动越野车具有低速大转矩、转矩输出响应极快等优异特性^[3-4]，同时各轮采用一体化的轮毂电机系统，凭借其各轮转矩独立可控的优势，借鉴履带式车辆滑移转向原理^[5]，能够实现车辆原地转向，为此相较于传统车辆可以极大地提升车辆的转向灵活性。

为了实现车辆的原地转向功能，众多学者基于不同工况采用不同的控制算法，控制各轮驱动转矩以满足对期望横摆角速度的快速响应。Hou^[6]提出了一种适用于道路附着的分级转向稳定性控制策略，有效提高了4WID车辆在变附着路面上的鲁棒性。孟祥明^[7]针对轮毂电机驱动车辆在不同附着系数路面上横摆稳定性下降的问题，提出了一种横向稳定性控制策略，基于附着系数合理分配驱动力矩, 有效增强了车辆的行驶稳定性。康亚彪^[8]以轮毂电机驱动的无人车为研究对象，分析其在原地转向运动时轮胎的受力变化，计算出实现原地转向所需克服的阻力矩，进而推导出所需的驱动力矩和输出功率。Chen^[9]针对双电履带车辆通过方向盘转角及车速指令利用BP神经网络算法拟定期望横摆角速度，并采用PID算法对滑转率与横摆角速度误差进行反馈，能够一定程度上提高转向时的稳定性。陈泽宇^[10]为研究履带式车辆原地转向的控制算法，将电机的输出力矩定义为横摆角速度的负反馈函数，负反馈增益则与方向盘转角直接相关，该算法在方向盘阶跃输入及动态调整输入中都具有良好的响应性和优秀的控制效果。张杰^[11]通过将方向盘转角及其变化率、转向半径偏差及其变化率作为模糊系统控制输入，对两侧转矩分别进行前馈补偿与反馈校正，能够有效地提升车辆在不同路面下转向的动态响应性与鲁棒性。综上，关于车辆原地转向控制相关的研究中，多数聚焦在特定用途的履带式车辆，且缺少对车辆原地转向的机理分析及运动过程中的轮胎特性变化。因此为了较好地实现分布式驱动车辆的原地转向控制效果，需要考虑车辆横摆转向时的响应性、稳定性以及鲁棒性等问题，同时结合转向时的轮胎特性变化、启动力矩决策与各轮转矩协调分配等方面设计一种适用性高且稳定可靠的原地转向控制算法。

本研究以四轮毂电机驱动越野车为研究对象，首先建立车辆原地转向动力学模型，分析原地转向机理，其次基于路面条件决策出横摆角加速度和期望横摆角速度，得到原地转向启动力矩，并利用模型预测控制算法跟踪期望横摆角速度并预分配各轮期望力矩，再次基于横摆稳定性控制协调各轮转矩，控制转向中心偏移量，最后完成原地转向控制策略的仿真验证和实车验证。

1 原地转向动力学分析

车辆在原地转向运动过程中，各轮胎处于纵横向运动强耦合状态，且轮胎附着利用率和路面附着力较大，因此为了明确原地转向控制策略的优化控制目标与稳定性约束条件，需要通过车辆动力学与轮胎力学模型探究转向横摆稳定性与轮胎纵横向力之间的耦合关系^[12]。

1.1 四轮毂电机驱动整车参数

试验车辆为课题组与相关合作企业联合研制的四轮毂电机驱动越野车如图 1所示，车辆和轮毂电机系统相关参数如表 1所示。

1.2 车辆原地转向运动轮胎受力分析

当车辆进行原地转向时，轮胎变形如图 1所示，轮胎接地印迹是长为2a，宽为2b的矩形。

当轮胎滑移率与侧偏角同时存在时，可通过地面垂向力以抛物线分布的联合轮胎模型对轮胎受力进行描述，合力大小为：

$ F= \begin{cases}\mu F_{\mathrm{z}}\left[3 \sigma \theta-\frac{(3 \sigma \theta)^2}{3}+\frac{(3 \sigma \theta)^3}{27}\right], & \sigma \leqslant \sigma_{\mathrm{m}} \\ \mu F_{\mathrm{z}}, & \sigma>\sigma_{\mathrm{m}}\end{cases} $

(1)

式中，μ为路面附着系数；F_z为地面法向作用力；σ为轮胎总滑移率; θ为轮胎固有参数与地面垂向力之间的函数；σ_m为轮胎处于完全滑动状态时滑移率能达到的最大值，可通过式(2)~(3)计算：

$ \sigma_{\mathrm{m}}=\frac{1}{\theta}, $

(2)

$ \theta=\frac{4 k a^2 b}{3 \mu F_{\mathrm{z}}}, $

(3)

式中k为轮胎单位面积的侧向刚度。

σ可由纵向滑移σ_x与侧向滑移σ_y合成得到：

$ \sigma=\sqrt{\sigma_{\mathrm{x}}^2+\sigma_{\mathrm{y}}^2}。$

(4)

由式(4)可知定义的总σ为正值，轮胎纵向力与侧向力分别可表示为：

$ \left\{\begin{array}{l} F_{\mathrm{x}}=\frac{\sigma_{\mathrm{x}}}{\sigma} F \\ F_{\mathrm{y}}=\frac{\sigma_{\mathrm{y}}}{\sigma} F \end{array}, \right. $

(5)

式(5)中σ_x，σ_y定义为：

$ \left\{\begin{array}{l} \sigma_{\mathrm{x}}=\frac{\omega \cdot r-v_{\mathrm{x}}}{v_{\mathrm{x}}} \\ \sigma_{\mathrm{y}}=\frac{v_{\mathrm{x}}}{\omega \cdot r} \tan \alpha \end{array}, \right. $

(6)

式中，ω为车轮转速; r为车轮半径; v_x为纵向车速; α为轮胎侧偏角。

1.3 车辆原地转向稳定域确定

车辆进行原地转向运动要求各动力执行器产生的横摆力矩克服转向阻力矩，还需精确控制各轮处于稳定区域内，车辆原地转向各轮受力如图 2所示。

$ F_{\mathrm{x} i}=F_{\mathrm{t} i}-F_{\mathrm{f}i}, $

(7)

$ F_{\mathrm{y} 1}+F_{\mathrm{y} 2}+F_{\mathrm{y} 3}+F_{\mathrm{y} 4}=0, $

(8)

式中，F_xi，F_ti，F_fi，F_yi分别为各车轮纵向力，驱动力，滚动阻力，侧向力，i=1, 2, 3, 4分别为左前、右前、左后以及右后车轮。

原地转向运动过程中车辆纵向速度与加速度皆为零，将车辆的纵向力和侧向力对质心取矩，即产生车辆原地转向所需的横摆力矩M_s以及转向中的摩擦阻力矩M_f：

$ M_{\mathrm{s}}=\left(-F_{\mathrm{x} 1}+F_{\mathrm{x} 2}-F_{\mathrm{x} 3}+F_{\mathrm{x} 4}\right) \cdot B / 2, $

(9)

$ M_{\mathrm{f}}=a\left(F_{\mathrm{y} 1}+F_{\mathrm{y} 2}\right)+b\left(F_{\mathrm{y} 3}+F_{\mathrm{y} 4}\right) 。$

(10)

通过选取车轮纵向滑转率λ作为自变量，以路面附着系数作为影响因素拟合得原地转向中M_s与M_f的动态变化曲线和理论$ \dot{\gamma}$与λ的关系曲线，如图 3，4所示，进而研究车辆在不同路面附着系数下原地转向时稳定性和响应性的控制权重。

将刚开始横摆转向时的λ作为最小滑转率λ_min，$ \dot{\gamma}$峰值处的λ设为最大滑转率λ_max，为提高原地转向的稳定性，要求各轮滑转率保持在[λ_min, λ_max]区间内。同时针对不同附着路面设计不同控制目标，在低附路面以提升车辆稳定性来减小原地转向中心的偏移量，在高附路面以提升车辆的横摆响应性跟踪驾驶员意图。

1.4 原地转向影响因素分析 1.4.1 车辆轴距-轮距比

为了实现车辆原地转向运动，车辆横摆力矩需克服转向阻力矩，且横摆力矩M_s，转向阻力矩M_f分别和轮距B，轴距L呈正相关关系，考虑到车辆L和B的参数变化是相对而言的，分开研究意义不大，为此将车辆轴距-轮距的比例关系(L/B)设为自变量来研究。由图 5可得最低电机转矩与(L/B)呈非线性增加关系，从图 6可看出μ_c与(L/B)呈非线性正相关关系。本次试验车辆(L/B)为1.67，计算得原地转向运动所需的最低路面附着系数μ_c为0.3，为此在设计车辆控制策略时，将0.3设为最低路面附着系数。

1.4.2 路面附着系数

在车辆(L/B)为1.67的前提下，车辆进行原地转向运动所需的最低路面附着系数μ_c为0.3，当路面附着系数低于0.3时，各车轮会出现过度滑转进而导致车辆持续失稳。随着路面附着条件的提升，车辆在原地转向时的轮胎滑转/侧偏现象在一定程度上有所下降，M_f也随之增大，车辆各驱动电机需输出更大的驱动力矩来克服M_f。为此本小节将μ作为自变量，得到图 7中车辆在不同μ下的原地转向中的最低电机转矩以及转向角的仿真试验曲线。

1.4.3 道路坡度

为了增加原地转向功能在复杂工况下的适用性，基于现有履带式车辆坡道转向的动力学分析^[13-14]，考虑坡度变化对越野车辆滑移转向的影响，以判断车辆坡道上原地转向是否能够安稳实现。

车辆在坡度为θ，方位角为φ的路面上，纵倾坡度α和侧倾坡度β为：

$ \tan \alpha=\frac{G \sin \theta \cos \varphi}{G \cos \theta}=\tan \theta \cos \varphi, $

(11)

$ \tan \beta=\frac{G \sin \theta \sin \varphi}{G \cos \theta}=\tan \theta \sin \varphi。$

(12)

设定法向反力作用中心的纵向偏移量为x₀；侧向偏移量为y₀；车辆重心高度为h_g；可得下列车辆力矩平衡方程表达式为：

$ \left\{\begin{array}{l} X h_{\mathrm{g}}=F x_0 \\ Y h_{\mathrm{g}}=F y_0 \end{array}。\right. $

(13)

法向反力作用中心的纵向、侧向偏移量为：

$ x_0=h_{\mathrm{g}} \tan \alpha=h_{\mathrm{g}} \tan \theta \cos \varphi, $

(14)

$ y_0=h_{\mathrm{g}} \tan \beta=h_{\mathrm{g}} \tan \theta \sin \varphi。$

(15)

车辆在坡道接地面受到的法向载荷为纵倾坡与侧倾坡上法向载荷之和，各轮法向载荷为：

$ \left\{\begin{array}{l} F_{\mathrm{z} 1}=\frac{b m g \cos \theta}{2 L}\left(1-h_{\mathrm{g}} \tan \theta \cos \varphi\right)+b m g \frac{h_{\mathrm{g}} \sin \theta \sin \varphi}{B L} \\ F_{\mathrm{z} 2}=\frac{b m g \cos \theta}{2 L}\left(1-h_{\mathrm{g}} \tan \theta \cos \varphi\right)-b m g \frac{h_{\mathrm{g}} \sin \theta \sin \varphi}{B L} \\ F_{\mathrm{z} 3}=\frac{a m g \cos \theta}{2 L}\left(1+h_{\mathrm{g}} \tan \theta \cos \varphi\right)+a m g \frac{h_{\mathrm{g}} \sin \theta \sin \varphi}{B L} \\ F_{\mathrm{z} 4}=\frac{a m g \cos \theta}{2 L}\left(1+h_{\mathrm{g}} \tan \theta \cos \varphi\right)-a m g \frac{h_{\mathrm{g}} \sin \theta \sin \varphi}{B L} \end{array}, \right. $

(16)

式中F_zi (i=1, 2, 3, 4)为各轮的法向载荷。

由于接地面瞬时转向中心会在车辆原有质心位置上产生一定偏移量ρ，对此侧向力Y对瞬时转向中心产生一个额外力矩，如式(17)表达：

$ \begin{gathered} M_{\mathrm{ts}}=M_{\mathrm{s}}-Y \rho= \\ \frac{\mu m g \cos \theta}{L}\left(\frac{L^2}{2}-\rho h_{\mathrm{g}} \tan \theta \cos \varphi+\rho^2+\frac{2 \rho^3 h_{\mathrm{g}} \tan \theta \cos \varphi}{L^2}\right) \\ -m g \rho h_{\mathrm{g}} \sin \theta \sin \varphi。\end{gathered} $

(17)

为验证不同的道路坡度角对车辆原地转向性能的影响，将坡道路面阻力系数设为0.05，路面附着系数设为0.85，得到如图 8所示的在道路坡度为5°，10°，15°，20°以及25°时车辆的转向中心偏移量随车辆的转向方位角变化的动态曲线。从图 8可以看出，当车辆所处坡度角超过20°时，车轮极易发生过度滑转，导致车辆失控侧翻现象发生，为此本研究为了保证车辆的安全性与稳定性，设定车辆原地转向适用坡度低于20°的高附路面。

2 原地转向动力学模型建立 2.1 车辆七自由度动力学模型

搭建车辆七自由度动力学模型来分析车辆原地转向时的稳定性与响应性的耦合关系，并作以下假设：

(1) 路面水平坚实，轮胎无下沉现象。

(2) 忽略整车悬架动态特性，无侧倾和俯仰，转向中心即为车辆质心。

通过图 3可得到如下各自由度的动力学方程：

车辆纵向动力学方程：

$ m a_{\mathrm{x}}=m\left(\dot{v}_{\mathrm{x}}-\gamma v_{\mathrm{y}}\right)=F_{\mathrm{x} 1}+F_{\mathrm{x} 2}+F_{\mathrm{x} 3}+F_{\mathrm{x} 4}。$

(18)

车辆侧向动力学方程：

$ m a_{\mathrm{y}}=m\left(\dot{v}_{\mathrm{y}}+\gamma v_{\mathrm{x}}\right)=F_{\mathrm{y} 1}+F_{\mathrm{y} 2}+F_{\mathrm{y} 3}+F_{\mathrm{y} 4} 。$

(19)

车辆横摆动力学方程：

$ I_{\mathrm{z}} \dot{\gamma}=\left(F_{\mathrm{x} 2}+F_{\mathrm{x} 4}-F_{\mathrm{x} 1}-F_{\mathrm{x} 3}\right) \cdot B / 2+ \\ \;\;\left(F_{\mathrm{y} 1}+F_{\mathrm{y} 2}\right) \cdot a-\left(F_{\mathrm{y} 3}+F_{\mathrm{y} 4}\right) \cdot b。$

(20)

上述式中，v_x为纵向车速；a_x为纵向加速度；v_y为侧向车速；a_y为侧向加速度；γ为车辆横摆角速度；I_z为横摆转动惯量；F_xi为各轮纵向力；F_yi为各轮侧向力；B为车辆轮距；a为车辆质心到前轴的距离；b为质心到后轴的距离。

2.2 轮胎模型

由于车辆在原地转向运动转向过程中车轮滑移与侧偏现象同时存在，为了提高车辆原地转向时的响应性和稳定性，本研究选取形式较为简单，计算效率高的Dugoff轮胎模型，方便原地转向过程中轮胎与路面相关的参数估计以及控制策略的设计。

$ \left\{\begin{array}{l} \bar{F}_{\mathrm{x} i}=\mu F_{\mathrm{z} i} \frac{C_{\mathrm{x} i} \lambda_i}{1-\lambda_i} f\left(L_i\right) \\ \bar{F}_{\mathrm{y} i}=\mu F_{\mathrm{z} i} \frac{C_{\mathrm{y} i} \tan \alpha_i}{1-\lambda_i} f\left(L_i\right) \end{array}, \right. $

(21)

$ f\left(L_i\right)=\left\{\begin{array}{ll} L_i\left(2-L_i\right), & L_i <1 \\ 1, & L_i \geqslant 1 \end{array}, \right. $

(22)

$ L_i=\frac{\left(1-\lambda_i\right)\left[1-\varepsilon v_{\mathrm{x} i} \sqrt{\left(C_{\mathrm{y} i} \tan \alpha_i\right)^2+\left(C_{\mathrm{x} i} \lambda_i\right)^2}\right]}{2 \sqrt{\left(C_{\mathrm{y} i} \tan \alpha_i\right)^2+\left(C_{\mathrm{x} i} \lambda_i\right)^2}}, $

(23)

式中，$\bar{F}_{\mathrm{x} i}, \bar{F}_{\mathrm{y} i} $分别为轮胎纵向力和侧向力；λ_i，α_i分别为各轮的纵向滑转率与轮胎侧偏角；ε为校正轮胎滑移对模型计算值造成影响的系数；L_i为判断是否位于线性区间的系数；C_xi，C_yi分别为轮胎纵向刚度与侧偏刚度。

2.3 轮毂电机模型

为了提高模型的计算效率及相对准确说明轮毂电机的动态响应特性，选用永磁同步水冷电机，并采用带有一阶惯性环节的驱动电机简化模型：

$ T_{\text {out }}=T_{\text {in }} \cdot \frac{1}{1+\tau s}, $

(24)

式中，T_out为电机输出转矩；T_in为输入目标转矩；τ为时间常数。

3 多目标协调优化的原地转向控制策略设计 3.1 原地转向控制策略架构设计

为了实现车辆在不同坡道与附着路面上都有较好的原地转向性能，原地转向控制整体架构如图 9所示，由目标决策层、转矩优化分配层与控制执行层组成。

3.2 峰值横摆角加速度和期望横摆角速度确定

在车辆进行原地转向时，轮胎达到最大可利用附着力，车辆附着极限下的纵向力F_{x_max1}为：

$ F_{\mathrm{x} \_\max 1}=\sqrt{(\mu m g)^2-F_y^2} 。$

(25)

利用车轮旋转动力学方程可得到起动时的车辆纵向力F_{x_max2}：

$ F_{\mathrm{x}\_\max 2}=\frac{T_{\mathrm{x}} i_{\mathrm{n}}-F_{\mathrm{f}} r}{J_{\mathrm{w}} \dot{w} r}, $

(26)

式中，T_x为车轮转矩；i_n为传动比；F_f为车轮滚动阻力；r为车轮半径；J_w为车轮的转动惯量；w_r为车轮转动的加速度。

在车辆横摆力矩克服阻力矩后开始转向运动，此时车轮突破路面附着极限，可将峰值横摆角加速度设定在附着椭圆上，并满足以下条件^[15]：

$ \left\{\begin{array}{l} F_{\mathrm{y}}=G\left(F_{\mathrm{x}}\right) \\ F_{\mathrm{x}}=\min \left(F_{x \_\max 1}, F_{\mathrm{x} \_\max 2}\right) \end{array}。\right. $

(27)

求解得峰值横摆角加速度$ \dot{\gamma}_{\max }$：

$ \begin{aligned} &\;\; \dot{\gamma}_{\max }=\left[\left(F_{\mathrm{x} 2}+F_{\mathrm{x} 4}-F_{\mathrm{x} 1}-F_{\mathrm{x} 3}\right) \cdot \frac{B}{2}+\right. \\ & \left.\left(F_{\mathrm{y} 1}+F_{\mathrm{y} 2}\right) \cdot a-\left(F_{\mathrm{y} 3}+F_{\mathrm{y} 4}\right) \cdot b\right] \frac{1}{I_{\mathrm{z}}} 。\end{aligned} $

(28)

将驾驶员输入的加速踏板开度信号进行归一化处理：

$ \xi_{\mathrm{s}}= \begin{cases}0, & P_{\text {acc }} \leqslant P_{\text {acc } \_0} \\ \frac{P_{\text {acc }}-P_{\text {acc } \_0}}{P_{\text {acc } \_ \text {max }}-P_{\text {acc } \_0}}, & P_{\text {acc }}>P_{\text {acc } \_0}\end{cases}, $

(29)

式中，ξ_s为加速踏板开度的比重系数；P_acc为加速踏板的开度信号量；P_{acc_0}为加速踏板开度的自由行程；P_{acc_ max}为加速踏板开度信号的最大值。

为了实现车辆在不同路面下的多目标协调优化，将路面附着系数与加速踏板开度信号设为系统输入量，基于一阶惯性环节，参考车辆最大横摆角加速度，设计期望的横摆角速度为：

$ t_{\mathrm{s}}=\frac{4 \tau_{\mathrm{s}}}{\sum \hat{\mu}_i}, $

(30)

$ \dot{\gamma}_{\mathrm{d}}=\lambda_{\mathrm{s}} \dot{\gamma}_{\max } \mathrm{e}^{-t / t_{\mathrm{s}}}, $

(31)

$ \gamma_{\mathrm{d}}=\left\{\begin{array}{ll} \int_0^t \dot{\gamma}_{\mathrm{d}} \mathrm{~d} t, & t<3 t_{\mathrm{s}} \\ \frac{\left(\mathrm{e}^3-1\right) t_{\mathrm{s}}}{\mathrm{e}^3} \dot{\gamma}_{\max }, & t \geqslant 3 t_{\mathrm{s}} \end{array}, \right. $

(32)

式中，t_s为不同附着路面极限下单位响应时间；τ_s为控制系统时间常数；$ \dot{\gamma}_\text{d}$为期望横摆角加速度；$ \dot{\gamma}_{\max }$为峰值横摆角加速度。

3.3 原地转向起动力矩决策

由于在起步阶段，轮胎附着仍处于限值以下导致路面附着系数不能收敛，且期望横摆角加速度起始即达到峰值引起电机转矩的快速提升^[16]，为此设计一种多模式判断的起动转矩决策方法，通过设定2 s的起动缓冲时间，并结合附着系数估计标志位slip_flag判定选择哪种力矩决策方案。

(1) t≤2&slip_flag=0

车辆处于静止时，各轮胎力尚未突破路面附着极限，此时路面附着系数估计不收敛。系统默认路面附着系数估计初始值$\hat{\mu}_0 $为1，且$ \dot{\gamma}_\text{d}$=γ_d=0，进而可求解单轮最大需求转矩为：

$ T_{\mathrm{req}}=\hat{\mu} F_{\mathrm{z}}。$

(33)

为了避免驱动轮在起动瞬间达到峰值转矩，车辆出现过度滑转，引入一阶惯性环节确保车辆起动力矩T_s平缓提升：

$ T_{\mathrm{s}}=T_{\mathrm{req}}\left(1-\mathrm{e}^{-t / t_{\mathrm{s}}}\right) 。$

(34)

(2) t≤2&slip_flag=1

车辆在缓冲时间内，车轮附着系数收敛，车辆极有可能发生过度滑转，由于此时的实际$\hat{\mu} $远低于$\hat{\mu}_0 $。为此需更新$\hat{\mu} $与最大T_req来抑制车轮过度滑转：

$ \left\{\begin{array}{l} T_{\mathrm{req}}=\hat{\mu} F_{\mathrm{z}} \\ T_{\mathrm{s}}=\xi_{\mathrm{s}} T_{\mathrm{req}}\left(1-\mathrm{e}^{-t / t_{\mathrm{s}}}\right) \end{array}。\right. $

(35)

(3) t>2&slip_flag=0

车辆虽度过缓冲期，但$\hat{\mu} $还未收敛，此时$ \dot{\gamma}_\text{d}$与γ_d跟随期望轨迹，$\hat{\mu} $保持初始值不变，目标转矩将作为模型预测控制器优化控制约束后求解的结果得出。

$ \left\{\begin{array}{l} \hat{\mu}=\hat{\mu}_0 \\ \dot{\gamma}_{\mathrm{d}}=\frac{\xi_{\mathrm{s}} \dot{\gamma}_{\max }}{\mathrm{e}^{t / t_{\mathrm{s}}}} \\ \gamma_{\mathrm{d}}=\int_0^t \dot{\gamma}_{\mathrm{d}} \mathrm{~d} t \\ T_{\mathrm{d}}=M P C_{\left(\gamma-\gamma_{\mathrm{d}}\right)} \end{array}。\right. $

(36)

(4) t>2&slip_flag=1

车辆达到原地转向稳定运动状态，$\hat{\mu} $收敛并更新初始值$\hat{\mu}_0 $，依据$\hat{\mu} $重新优化$ \dot{\gamma}_\text{d}$与γ_d，并通过模型预测控制算法求解最优转矩T_d。

3.4 基于模型预测控制的期望力矩决策

为了达到车辆原地转向时对$\dot{\gamma} $的准确跟踪，利用车辆动力学模型进行模型预测控制(Model Predictive Control，MPC)求解^[17]，控制器初步采用转矩平均分配的策略进行求解，基于$ \dot{\gamma}_\text{d}$求解得到驱动电机目标输出转矩：

$ T_{\mathrm{d}}=\left(T_{\mathrm{m}}\right)_{\mathrm{r}}=\frac{\left(\dot{\gamma}_{\mathrm{d}} I_z+M_{\mathrm{f}}\right) r i_{\mathrm{m}}}{\boldsymbol{B}} 。$

(37)

采用向前欧拉法对式(36)进行离散化处理，可得到离散空间状态方程表达式：

$ \tilde{x}(k+1)=\boldsymbol{A}_{k, t} \tilde{x}(k)+\boldsymbol{B}_{k, t} \tilde{u}(k), $

(38)

式中，$\boldsymbol{A}_{k, t}=\left[\begin{array}{ll} 1 & 0 \\ 0 & 1 \end{array}\right] $；t为系统采样时间；B_{k, t}=$\left[\begin{array}{llll} \frac{T i_{\mathrm{m}}}{m r} & \frac{T i_{\mathrm{m}}}{m r} & \frac{T i_{\mathrm{m}}}{m r} & \frac{T i_{\mathrm{m}}}{m r} \\ \frac{T B i_{\mathrm{m}}}{2 r I_{\mathrm{z}}} & \frac{T B i_{\mathrm{m}}}{2 r I_{\mathrm{z}}} & \frac{T B i_{\mathrm{m}}}{2 r I_{\mathrm{z}}} & \frac{T B i_{\mathrm{m}}}{2 r I_{\mathrm{z}}} \end{array}\right] $。

将状态误差和控制误差相结合，可得到新的系统状态量和状态空间表达式：

$ \boldsymbol{\xi}=\left[\begin{array}{l} \tilde{x}(k) \\ \tilde{u}(k-1) \end{array}\right], $

(39)

$ \left\{\begin{array}{l} \boldsymbol{\zeta}(k+1)=\tilde{\boldsymbol{A}} \boldsymbol{\zeta}(k)+\tilde{\boldsymbol{B}} \Delta \tilde{u}(k) \\ \boldsymbol{\eta}(k)=\boldsymbol{C} \boldsymbol{\zeta}(k) \end{array}, \right. $

(40)

式中，$\tilde{\boldsymbol{A}}=\left[\begin{array}{cc} \boldsymbol{A} & \boldsymbol{B} \\ 0 & \boldsymbol{I}_{N_{\mathrm{u}}} \end{array}\right] ; \quad \tilde{\boldsymbol{B}}=\left[\begin{array}{c} \boldsymbol{B} \\ \boldsymbol{I}_{N_{\mathrm{u}}} \end{array}\right] $；$ \boldsymbol{C}=\left[\begin{array}{ll} \boldsymbol{I}_{N_{\mathrm{x}}} & 0 \end{array}\right]$；C为观测矩阵；N_u为控制量维数；N_x为状态量维数。

预测时域内的系统输出方程为：

$ Y=\psi \boldsymbol{\xi}(k)+\boldsymbol{\varTheta} \Delta U, $

(41)

式中 ξ(k)为系统当前时刻的状态量。

$\begin{array}{c} \boldsymbol{\xi}(k)=\left[\begin{array}{c} \tilde{x}(k) \\ \tilde{u}(k-1) \end{array}\right], Y=\left[\begin{array}{l} \eta(k+1) \\ \eta(k+2) \\ \;\vdots \\ \eta\left(k+N_{\mathrm{c}}\right) \\ \eta\left(k+N_{\mathrm{p}}\right) \end{array}\right], \\ \psi=\left[\begin{array}{l} C \tilde{A} \\ C \tilde{A}^2 \\ \;\;\vdots \\ C \tilde{A}^{N_{\mathrm{c}}} \\ C \tilde{A}^{N_{\mathrm{p}}} \end{array}\right], \quad \Delta U=\left[\begin{array}{c} \Delta \tilde{u}(k) \\ \Delta \tilde{u}(k+1) \\ \vdots \\ \Delta \tilde{u}(k+c) \\ \tilde{u}(k+p) \end{array}\right], \\ \boldsymbol{\varTheta}=\left[\begin{array}{cccc} C \tilde{B} & 0 & \cdots & 0 \\ C \tilde{A} \tilde{B} & \tilde{C} \tilde{B} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ C \tilde{A}^{N_{\mathrm{c}}-1} \tilde{B} & C \tilde{A}^{N_{\mathrm{c}}-2} \tilde{B} & \cdots & C \tilde{A}^0 \tilde{B} \\ C \tilde{A}^{N_{\mathrm{p}}-1} \tilde{B} & C \tilde{A}^{N_{\mathrm{p}}-2} \tilde{B} & \cdots & C \tilde{A}^{N_{\mathrm{p}}-N_{\mathrm{c}}} \tilde{B} \end{array}\right] 。\end{array} $

设计包括期望状态的跟踪能力、控制量和控制增量极值约束问题两个子项的多约束优化条件下的目标函数：

$ \begin{cases}J(k) & =\sum\limits_{i=1}^{N_{\mathrm{p}}}\left\|\eta(k+i)-\eta_{\mathrm{r}}(k+i)\right\|_\boldsymbol{Q}^2+ \\ & \sum\limits_{i=1}^{N_\mathrm{c}-1}\|\Delta u(k+i)\|_\boldsymbol{R}^2+\rho \varepsilon^2 \\ \text { s. t. } & \left\{\begin{array}{l} \Delta u_{\min } \leqslant \Delta u_t \leqslant \Delta u_{\max } \\ u_{\min } \leqslant u_t+A_{\mathrm{I}} \Delta u_t \leqslant u_{\max } \end{array}\right.\end{cases}, $

(42)

式中，ρ为权重系数；ε为松弛因子；Q与R为权重矩阵；Δu_t为在t时刻的控制增量，Δu_min与Δu_max分别为其上下限值；u_t为在t时刻的控制量；u_min与u_max分别为其大小限值。

通过求解目标函数min(J)得到控制序列中的首个元素作为系统最佳控制增量，可得到当前时刻状态反馈控制律为：

$ u(t)=u(t-1)+\Delta u_t^*, $

(43)

式中u(t)为当前时刻车辆跟踪γ_d所需的驱动电机输出转矩。

3.5 基于转向中心偏移量的各轮转矩优化控制

车辆在不同路面工况下进行原地转向运动时车轮可能出现过度滑转或车身横摆过快等失稳现象，为此需要通过控制各轮输出转矩将各轮滑转率控制在一定范围内。考虑到车辆在坡道上转向时会发生载荷转移，需设计基于坡道轮荷转移的转矩补偿策略，优化各轮胎的附着利用率。

3.5.1 驱动转矩调节的介入与退出机制

在车辆原地转向过程中设定基于自适应最大转矩判断系统(Adaptive Maximum Torque Determination，AMTD)触发机制如图 10所示，当车轮滑转率过大时(λ_i>λ_opt+ε)，此时车轮出现过度滑转，为了增加车辆的稳定性控制，利用当前路面附着条件的轮胎纵向力F_xi计算得最大限制转矩T_{AMTD_ i}，并将其设置在AMTD的触发条件中，计算得各轮的稳定工作区间：

$ T_{\text {AMTD }_{-} i}=F_{\mathrm{x} i} r+J_\omega \overline{\dot{\omega}}_i+\xi, $

(44)

式中ξ为最大驱动转矩的校准参数。

3.5.2 基于滑模控制的驱动转矩调节

为了抑制各轮的过度滑转，确保车辆原地转向时的稳定可控，设计一种基于滑模控制(Sliding Mode Control, SMC)算法的转矩矢量调节策略^[18]，使各轮滑转率保持在[λ_min, λ_max]区间内。当触发转矩调节机制时，基于SMC算法对各轮电机输出转矩进行矢量调节，能够得到单个车轮的转矩控制模型为：

$ \left\{\begin{array}{l} F_{\mathrm{z} i} \dot{v}_{\mathrm{x} i} / g=F_{\mathrm{x} i}-F_{\mathrm{r} i} \\ \mathrm{~J} \dot{\omega}_i=T_i-F_{\mathrm{x} i} r \end{array}, \right. $

(45)

式中，$\dot{v}_{\mathrm{xi}} $为车轮的加速度；F_ri为行驶阻力；T_i为驱动轮输出转矩; ω_i为各车轮转动的加速度。

分别选取λ_i，T_i作为系统的控制变量x，系统的输入量u可得：

$ \dot{x}=\dot{\lambda}_i=\frac{v_{\mathrm{x} i} \dot{\omega}_i}{\omega_i^2 r}-\frac{\dot{v}_{\mathrm{x} i}}{\omega_i r}=\frac{v_{\mathrm{x} i}}{\mathrm{~J} \omega_i^2 r} u-\frac{F_{\mathrm{x} i}-F_{\mathrm{r}}}{\omega_i r m_{1 / 4}}-\frac{v_{\mathrm{x} i} F_{\mathrm{x} i}}{\mathrm{J} \omega_i^2}。$

(46)

定义系统的滑模面：

$ S=e=\lambda-\lambda_{\mathrm{opt}} 。$

(47)

对S求导得：

$ \dot{S}=\dot{\lambda}_i+c_\lambda\left(\lambda-\lambda_{\mathrm{opt}}\right), $

(48)

式中，S为滑模面函数；e为滑移率误差；c_λ为误差常数。

为了削弱系统的抖振，采用饱和函数sat(S/Φ)代替符号函数sgn(S)的指数趋近律：

$ \operatorname{sat}(S / \varPhi)=\left\{\begin{array}{l} S / \varPhi|S|<\varPhi \\ \operatorname{sgn}(s)|S| \geqslant \varPhi \end{array}, \right. $

(49)

$ \dot{S}=-\varepsilon \operatorname{sat}(S / \varPhi)-k S, $

(50)

式中，Φ为边界层的厚度；ε为速度常数；k为滑模面的调节系数。

计算可得单轮的驱动转矩调节量：

$ \begin{aligned} \Delta T_i & =\frac{\mathrm{J} \omega_i}{v_{\mathrm{x} i}}\left(F_{\mathrm{x} i}-F_{\mathrm{r} i}\right)+\frac{\mathrm{J} \omega_i^2 r \dot{\lambda}_i}{v_{\mathrm{x} i}}+ \\ F_{\mathrm{x} i} r & -\frac{\mathrm{J}_{\omega i}^2 r}{v_{\mathrm{x} i}}(\varepsilon \operatorname{sat}(S / \varPhi)+k S)。\end{aligned} $

(51)

3.5.3 基于坡道轮荷转移的转矩补偿控制

为了减小车辆在原地转向过程中由于载荷转移导致的轮荷较小的车轮发生滑转失稳风险的概率，采用基于轮荷转移的转矩补偿策略对各轮转矩进行重构分配，基于道路坡道信息与车身运动状态的各轮载荷计算^[19]如下表示：

$ \begin{gathered} {\left[\begin{array}{c} F_{z 1} \\ F_{z 2} \\ F_{z 3} \\ F_{z 4} \end{array}\right]=\left[\begin{array}{c} b \\ b \\ a \\ a \end{array}\right] \frac{m g \cos \theta}{2 L}+\left[\begin{array}{c} -b \\ -b \\ a \\ a \end{array}\right] \frac{m g h_{\mathrm{g}} \sin \theta}{2 L}+} \\ {\left[\begin{array}{c} -1 \\ -1 \\ 1 \\ 1 \end{array}\right] \frac{m a_{\mathrm{x}} h_{\mathrm{g}}}{2 L}+\left[\begin{array}{c} -b \\ b \\ -a \\ a \end{array}\right] \frac{m g h_{\mathrm{g}} \sin \varphi}{B L}+\left[\begin{array}{c} -b \\ b \\ -a \\ a \end{array}\right] \frac{m a_{\mathrm{y}} h_{\mathrm{g}}}{B L} 。} \end{gathered} $

(52)

此时后轴车轮转矩分配系数ζ如下表示：

$ \begin{gathered} \zeta=\frac{T_{\mathrm{r}}}{T_{\mathrm{req}}}=\frac{F_{\mathrm{zr}}}{F_{\mathrm{zf}}+F_{\mathrm{zr}}}= \\ \frac{a \cdot g \cos \theta+\left(a_x+g \sin \theta\right) h_{\mathrm{g}}}{L \cdot g \cos \theta+(a-b) g h_{\mathrm{g}} \sin \theta}, \end{gathered} $

(53)

式中，F_zf为前轴载荷；F_zr为后轴载荷；T_r为后轴所需力矩；T_req为总需求力矩。

基于ζ进行各轮转矩重分配：

$ \left[\begin{array}{c} T_{\mathrm{r} 1} \\ T_{\mathrm{r} 2} \\ T_{\mathrm{r} 3} \\ T_{\mathrm{r} 4} \end{array}\right]=\left[\begin{array}{c} (1-\zeta) / 2 \\ (1-\zeta) / 2 \\ \zeta / 2 \\ \zeta / 2 \end{array}\right] T_{\mathrm{r}}, $

(54)

式中T_ri (i=1，2, 3, 4)为基于实时轮荷转移补偿的各车驱动转矩分配值。

当车辆在稳态工况下转向时，车辆的横摆力矩ΔM_s为：

$ \begin{gathered} \left.\Delta M_{\mathrm{s}}=\frac{B}{2}\left[\left(F_{x 2}-F_{x 1}\right) \sin \theta+F_{x 4}-F_{x 3}\right)\right]= \\ \frac{B}{2 r}\left(\Delta T_{\mathrm{R}} \sin \theta-\Delta T_{\mathrm{L}}\right), \end{gathered} $

(55)

式中，B为轮矩；ΔT_R，ΔT_L分别为车辆右侧和左侧两车轮的驱动转矩差。

基于各轮荷转移的横摆力矩分配为：

$ \left[\begin{array}{c} \Delta T_{\mathrm{L}} \\ \Delta T_{\mathrm{R}} \end{array}\right]=\left[\begin{array}{c} 1-\zeta \\ \zeta \end{array}\right] \frac{2 r \cdot \Delta M_{\mathrm{s}}}{B} 。$

(56)

联立上式可得：

$ \left[\begin{array}{c} T_{\text{r}1} \\ T_{\text{r}2} \\ T_{\text{r}3} \\ T_{\text{r}4} \end{array}\right]=\left[\begin{array}{c} (1-\zeta) / 2 \\ (1-\zeta) / 2 \\ \zeta / 2 \\ \zeta / 2 \end{array}\right] T_\text{r}+\left[\begin{array}{c} -(1-\zeta) \\ 1-\zeta \\ -\zeta \\ \zeta \end{array}\right] \frac{2 r \cdot \Delta M_\text{s}}{B} 。$

(57)

4 仿真与实车验证

在不同坡道与附着系数路面上进行车辆原地转向控制算法的离线仿真测试与实车道路验证，综合验证所设计的控制策略的有效性。

4.1 原地转向控制策略仿真验证

基于Carsim与Matlab/Simulink软件搭建对应控制模型，在本次仿真试验中，车辆在坚实平坦路面上进行原地转向运动，将高附、低附、坡道路面的路面附着系数分别设置为0.9，0.3，0.85，坡道路面坡度设置为10°，仿真试验中车辆的转向方向由方向盘转角控制，横摆运动速度由加速踏板控制，且与加速踏板开度信号成正相关关系，假设横摆运动方向为逆时针，前轮转角保持为零，试验结果如下。

如图 11所示，随着油门踏板与方向盘转角的逐渐增加，车辆各轮毂电机输出转矩值从0 N·m快速增加，达到峰值后车辆开始由静止状态切换至横摆运动，最后输出转矩值维持稳定。控制系统在识别到路面附着系数与道路坡度后，自适应调节各轮输出转矩，高附路面下横摆角速度快速增加在5.4 s时达到峰值94.5 rad/s，此时车辆能够较好地跟随期望横摆角速度，而在MPC单独控制的结果中出现了较大的超调量，其中在2.8 s时达到最大超调量66.4 rad/s。低附路面下横摆角速度逐渐增加于8.7 s时稳定在65.8 rad/s附近，相比MPC单独控制，转矩协调控制系统能够更快速稳定地跟随目标横摆角速度。坡道路面下横摆角速度逐渐增加于8.23 s达到峰值55.6 rad/s，且能较为稳定地跟踪期望横摆角速度。而在MPC单独控制中出现很大的超调与波动量，虽后续也能跟踪期望横摆角速度，但依旧存在一定波动量。

从图 12可得，车辆在高附路面条件下，由MPC单独控制于2.44 s时车轮出现完全滑转，其后虽有所缓解，由于驱动转矩持续波动，致使车轮转速不断震荡无法收敛，而在联合控制下的车轮转速峰值大幅降低，最终车辆稳态转向时的转速也趋于稳定。车辆在低附路面条件下，由MPC单独控制于4.8 s时车轮出现完全滑转，随后滑转现象有所缓解，但其收敛速度明显更慢。而在联合控制下车轮转速峰值大幅降低，随即车轮发生过大滑转，转矩调节迅速介入控制，最终稳态转向阶段各轮轮速趋于稳定。车辆在坡道路面条件下由MPC单独控制于1.42 s时车轮出现完全滑转，在车辆转过180°后，车轮再次完全滑转且不断震荡，最终各轮速趋于稳定。而联合控制的车轮转速略微降低，转向初期车轮同样会出现过度滑转，此时转矩调节系统及时介入控制，滑转过度现象开始缓解，最终各轮速趋于稳定且系统的收敛性明显更快。

由图 13可以看出，相较于MPC单独控制而言，加入转矩协调优化的联合控制车辆在高附、低附、坡道路面上的纵向车速分别降至0.38，0.19，0.45 km/h，侧向车速分别降至0.32，0.12，0.54 km/h，车辆质心位置沿x轴最大偏移量分别减小至0.11，0.12，0.41 m，沿y轴最大偏移量分别减小至0.25，0.19，0.22 m，为此可知轮胎纵、侧向力方面的稳定性与收敛速度皆优于MPC单独控制，进而验证了联合控制的必要性。

4.2 实车试验验证

通过选取坡道路面利用联合控制进行实车原地转向试验，同时与仿真试验结果对比分析，综合验证所设计的控制策略的真实有效性。

在本试验中，选取坚实的混凝土纵坡路面作为坡道路面进行实车原地转向测试，当前路面附着系数约为0.8，为了确保试验时车辆的安全稳定可控，为此选取10%的路面坡度，车辆初始静止在坡道路面上，其余试验条件与高附路面保持一致，试验结果如下。

图 14可以看出，车辆在坡道路面各驱动转矩随着车辆转向方位角的增加，各轮荷补偿控制开始介入调节，在1.27 s时达到峰值转速为193 rpm，峰值转矩为885 N·m，峰值车速为1.15 km/h，横摆角速度于5.6 s时达到峰值为56.5 rad/s，在系统转矩协调优化控制下，各车轮未出现明显滑转现象，并在7.1 s时完成原地转向，后续车辆趋于稳定转向。虽然试验结果波动率及其范围超过仿真工况，但考虑到实际环境噪声的影响，也符合实际情况。试验结果整体趋势与仿真试验结果基本一致，同时整体误差控制在可接受的范围内，从而说明了所设计的控制策略在各种路面差异下的有效性与适用性。

5 结论

(1) 首先通过建立车辆原地转向动力学模型，对运动车轮受力进行分析，计算得原地转向时的横摆力矩与阻力矩，确定原地转向时的车身稳定域。其次基于理论分析与仿真试验总结车辆自身的结构参数(L/B)、路面附着与地形坡度等条件对原地转向性能的影响规律，明确原地转向的稳定域约束条件及控制目标。

(2) 针对车辆原地转向功能在不同路面上的适用性问题，设计了一种多目标协调优化的原地转向控制策略，首先依据不同坡度与路面附着差异决策出当前路面下的峰值横摆角加速度与期望横摆角速度。其次根据不同模式判断决策出原地转向的起动力矩，并利用MPC算法进行期望横摆角速度跟踪与期望力矩决策。再次为了保证车辆原地转向时的稳定性并减小转向中心的偏移，基于横摆稳定性进行转矩协调优化控制。最后对原地转向控制策略进行仿真与实车试验，结果表明本研究所设计的控制策略使期望横摆角速度的跟随误差控制在7%以内，转向中心最大偏移量控制在470 mm以内，从而实现车辆在高附路面上提升其横摆响应性、低附与坡道路面提升其横摆稳定性的控制目标。