自行车动力学建模及稳定性分析研究综述1)

王囡囡 熊佳铭 刘才山

(北京大学工学院航空航天工程系，湍流与复杂系统国家重点实验室,北京 100871)

引言

自德国人Karl von Drais 于1817 年发明了世界上第一辆“自行车”以来，自行车已有200 多年的发展史.不同学者从不同角度给出了较为全面的描述.其中，Herlihy[1]详细介绍了自行车本身结构的发展和演化过程；Schwab 等[2]概括性地总结了不同自行车动力学模型及其相关的试验验证工作，以及当考虑骑行人对自行车的控制作用时,不同控制行为方式及其相关的试验工作；而本文主要从约束角度详细总结了不同的自行车动力学模型，以及对自行车特有的自稳定性进行了系统性地阐述.

自行车具有结构简单，使用方便等特点，同时它的动力学行为与一些复杂的力学和数学问题密切相关.静态来看，自行车很难通过前后轮与路面的接触实现稳定直立.因此，其物理本质与“倒立摆”结构类似，都属于一类静态不稳定系统[2].但当自行车在适当的运动速度条件下，却能够实现动态自稳定的直立“自行”.这种自行稳定性现象在历史上引起了众多学者的研究兴趣，并出现了多种多样的理论解释.Rankine[3]将自行车简化为只包含车轮偏转角的单自由度动力学系统，粗略地认为自行车的平衡源于骑行者对车把转角的控制；Bourlet[4]建立了考虑车轮陀螺效应的自行车简单动力学模型，并认为自行车的自稳定性与车轮的陀螺效应有关；文献[5-8]也建立了各自的自行车等效动力学模型，但这些模型在反映自行车的自稳定特征时均存在一定程度上的不足.1899 年法国数学家Carvallo[9]和英国剑桥大学Whipple[10]应用刚体动力学理论，并引入一些线性化的假设，分别独立地建立了自行车做匀速直线运动时的线性化摄动动力学模型.特别是，两位研究者首次将自行车简化为由4 个刚体组成的简单构型，并给定了标准的几何和物理参数.该简化的自行车构型被后来的研究者广泛采用[11-15]，称之为“Carvallo-Whipple 自行车构型”或简称为“Whipple 构型”.

随着对自行车研究的不断深入，以及动力学理论和计算机技术的发展，自行车的研究范畴逐步扩大.当前，关于自行车的研究内容不仅涉及经典力学[16-19]的多个分支，而且与应用数学[20-21]，控制理论[22-25]，机器人[26]等多个学科领域密切相关.

在动力学方面，自行车涉及的问题主要包括：(1)由于两个接触体的运动约束涉及运动物体轮廓面(线)上的一对动点的运动学关系，建立自行车在复杂曲面上的完整和非完整约束方程需要较为复杂的数学表述；(2)经典力学提供了丰富的动力学建模方法和手段，但是不同建模方法将会导致自行车动力学模型具有不同的表现方式.因此，应理解不同动力学模型之间的内在关联关系，揭示自行车内在的动力学规律；(3)自行车的稳定性分析；(4)动力学实验和自行车相关的控制问题.本文将围绕以上4 个方面，期望对自行车动力学的研究进展给出较为全面的综述.

论文的结构如下：第2 节简要描述了Whipple 自行车构型，并对自行车动力学研究中涉及到的完整约束和非完整约束进行了概述；第3 节系统总结了自行车动力学建模中的各类方法；第4 节对自行车3 种动态平衡状态的稳定性问题的研究进行了总结；最后，展望了自行车动力学的研究前景.

1 自行车构型描述及约束方程

图1 为通常采用的Carvallo-Whipple 自行车构型系统.该系统主要由前后框架H 和B(骑行人作为刚体固结于后框架)、前车轮F 和后车轮R 四个刚体组成.忽略铰接处的间隙、摩擦、阻尼等因素，并在建模过程中考虑以下几个方面的假设：

(1)前、后轮边缘为刀刃型结构，与地面的接触关系为点接触；

(2)前、后轮是圆周对称结构，其质心与各自圆心重合；

(3)假设前轮和后轮与地面之间的最大静摩擦系数足够大，以保证前后轮运动时处于纯滚动状态，即与地面无任何滑移；

(4)前框架和后框架分别关于前轮和后轮所在的平面对称，一是几何对称，二是质量分布对称.

图1 Carvallo-Whipple 自行车构型图Fig.1 The benchmark Carvallo-Whipple bicycle in an upright,straight reference configuration

为了描述自行车的运动状态，定义空间惯性坐标系FI={O−ijk}，并将坐标原点设置在后车轮与路面的接触点P处，i轴的正向指向前车轮，k轴正向背离地面，j轴正向指向骑行者的左侧.分别在4 个刚体质心处建立固连坐标系

自行车的4 个刚体通过3 个理想铰链连接，当不考虑路面对车轮的约束时，描述自行车自由运动需要9 个(4 × 6 −3 × 5=9)广义坐标.Xiong等[27-28]将9 个广义坐标定义如下：选择前叉转轴与前框架坐标轴的交点D为参考点，其在惯性坐标系FI中的3 个坐标分量(x,y,z)作为系统的3 个广义坐标.遵循3-1-2 的旋转顺序，选择后框架固连坐标系Fb相对于惯性坐标系FI的3 个欧拉角(ψ,θ,)为系统的另外3 个广义坐标，其中ψ 为绕着z轴的偏航转角，θ 为绕着x轴的倾斜角，为绕着y轴的俯仰角.另外3 个广义坐标分别为车把转角δ(前框架相对于后框架的转角)，后、前轮相对于后、前框架的转角φr和φf.这9 个广义坐标可用向量表示为：q=[x,y,z,ψ,θ,,δ,φr,φf]T.

另外，描述Whipple 自行车构型需9 个几何参数和16 个质量参数.这9 个几何参数包括：自行车的轴距w，前轮尾迹c(即前叉转轴延长线与地面的交点和地面接触点Q之间的距离，当该交点在Q点右侧时，c为正值，否则为负值)，前叉转轴相对于竖直线的倾斜角λ，前后轮的半径(rf,rr)，前框架质心位置(xh,yh)和后框架质心位置(xb,yb).16 个质量惯性参数包括：4 个刚体的质量(mb,mh,mr,mf)以及它们相对于各自质心的惯量张量(Ib,Ih,Ir,If)(共有4×2+2×2=12 个非零分量).Whipple 自行车构型的经典参数可参考文献[17].根据车轮纯滚动的假设条件，容易理解自行车是一个同时受完整和非完整约束的多刚体系统.然而，与Chaplygin 雪橇[29]和独轮机器人[30]等传统的非完整系统相比，自行车中的约束方程不易用显式方程来描述.由于通过铰链约束连接的前后车轮的轮心和前后框架之间形成了复杂的链式传递关系，接触点处的几何约束往往是隐式和非线性的.为建立Whipple 自行车构型下的线性化动力学模型，Meijaard[17]对非线性约束方程进行了线性化近似处理；Psiaki[13]利用解析几何的方法首次推导了完整约束方程的解析形式，其为四次代数方程，且具有高度非线性；Basu-Mandal 等[18]基于寻找两个车轮的最低点为参考条件推导了约束方程；类似地，Peterson 和Hubbard[31-32]考虑轮胎的形状，采用直观的解析几何方法建立了完整和非完整约束方程；Wang 等[33]对两个车轮圆周进行参数化，利用极值条件(车轮与地面的接触点为轮缘的最低点)推导出车轮的完整约束方程.

理论上，接触点处的几何约束方程与两个接触体的几何外形密切相关.因此，在建立由点接触引起的几何约束方程时，需要对两个接触体的几何外形进行参数化处理，进而依赖接触条件，确定接触点处的外形参数与系统广义坐标之间的函数关系，在此基础上，建立由系统广义坐标表示的几何约束方程[30,34].基于Zhao 和Liu[34]关于接触约束的理论描述，Xiong 等[27-28]分别建立了自行车在水平面和旋转曲面上的几何约束方程和速度约束方程.当自行车在水平面上运动时[28]，通过对前后轮轮廓和水平面参数化处理，最后可以得到由广义坐标表示的两个几何约束方程

以上两个方程分别是广义坐标z和关于广义坐标(θ,δ)的隐式函数，详细推导过程参见文献[28].显然，根据以上两个几何约束关系，z和不再是独立的广义坐标，从而可以将自行车构型空间的维数由9 维降为7 维，且定义降维后的独立广义坐标向量为=[x,y,ψ,θ,δ,φr,φf]T.

假设自行车在水平面上作纯滚动，即每个接触点处速度为零.这一速度条件使得每个接触点处对应两个独立的速度约束：该接触点的速度在水平面上的两个分量为零.考虑自行车有两个接触点，故自行车系统共有4 个纯滚动约束.

为恰当表示以上4 个滚动约束方程，需要在接触点处建立局部的接触坐标系.该坐标系不一定要求是正交的，但要求局部接触坐标系中的一个坐标轴对应接触点处共切面的公法线方向[34].利用接触条件得到的接触体的轮廓参数关于广义坐标的函数关系，然后利用两个接触体在接触点处的相对速度的运动学关系，并将该相对速度为零的矢量条件映射到所建立的局部接触坐标系上.这样，在每个接触点处可以建立3 个速度约束方程.容易证明，沿两个接触点处公法线方向的两个速度约束方程一定是两个几何约束方程(1)的全微分形式.自行车两个接触点引起的独立的速度约束方程仅有4 个，它们对应4 个独立的非完整约束方程[30].同时，这也意味着自行车系统的自由度为3，(这里的自由度采取分析力学中的传统定义，即系统的独立虚位移数[35]).

对受约束的多体系统来说，可以将几何约束和微分约束在速度空间中进行统一表示[36].对Whipple 自行车系统来说，可以直接利用两个接触点相对速度为零的条件，并将以上条件映射到接触点的局部坐标系中，进而可以将相应的完整和非完整约束统一表示为

当自行车在曲面上运动时，由于曲面接触点处的几何参数与描述自行车构型之间的广义坐标存在非常复杂的函数关系，通常需要采用计算机符号运算得到系统的约束方程.这一过程可以按照如下的流程进行:在两个接触点P和Q处分别定义三个曲面曲线参数，记为pi=[pi,1,pi,2,pi,3]T,i=r,f，其中，前两个曲面参数表征接触点在路面上的位置信息，最后一个曲线参数表征接触点在轮缘上的位置信息.根据文献[34]给出的接触几何条件，可以建立以上6 个参数p=与系统广义坐标q之间的函数关系以及两个几何约束方程如下

结合以上关系，并利用接触点处相对速度等于零的条件，将接触点的相对速度表示在局部接触坐标系中，可以得到自行车在任意曲面上运动时的速度约束方程[28]

2 动力学建模

在经典动力学框架内，存在多种不同的力学原理和方法建立系统的动力学模型.然而，不同方法可能导致动力学方程的表现形式存在很大不同[37].在缺乏计算机符号运算等便利分析工具的年代，如何建立自行车系统完备的动力学模型是一个具有挑战性的问题[12,17,19,33].Getz 和Marsden[21]采用几何力学的方法建立了一个极为简化的自行车动力学模型，并研究了自行车控制问题.Neimark 和Fufaev[14]由于对自行车势能分析时的错误，得到了一组错误的动力学方程组；Papadopoulos 等[16]为建立摄动后的系统线性化动力学方程，手工推演时不得不引入一些直观的假设；Meajaard 等[17]系统地综述了历史上各种形式的自行车动力学建模方法，并针对经典的Whipple 自行车构型，不仅采用计算机符号推导了自行车在平面上运动的非线性动力学方程组，而且与手工推导得到的线性化动力学模型进行对比，验证了模型的正确性.论文中还提供了在Whipple 自行车构型经典参数下，所得到的线性化方程的系数矩阵的具体表达式以及线性系统的具体特征值，并将其作为其他模型校验基准.Basu-Mandal 等[18]采用矢量力学中的Newton-Euler 方法和分析力学中的Euler-Lagrange 方法，并结合计算机符号运算，得到Whipple 自行车完备的非线性微分代数方程组，并通过数值结果验证了两种模型的等效性.作者还通过研究自行车做匀速直线运动时的线性化方程验证模型的正确性和精度.在此基础上给出了自行车在水平面上做匀速圆周运动的动态平衡解的完备描述.Wang 等[33]同样采用欧拉--拉格朗日方法建立了Whipple 自行车的非线性模型，并通过符号运算工具消除了完整约束和非完整约束，由此得到约化的动力学模型；Huang 等[38]利用Kane 方法推导自行车的控制微分方程.Sharp[12]考虑车轮滑动作用建立了自行车线性化方程；Meijaard 和Schwab[39]在建立自行车动力学模型时考虑了轮胎几何、路况以及系统阻尼等因素；Boyer 和Porez[40]从几何力学角度研究了自行车动力学问题，并通过主纤维丛上的一般约化方法建立动力学模型.由于自行车中所包含的约束方程以及动力学方程组的表现形式异常复杂，人们往往需要通过数值计算[41-42]或实验[32,43-44]验证模型的正确性.

对非完整系统来说，采用Gibbs-Appell 方法可以得到具有最小维数的动力学方程[19,32,45].Xiong等[28]基于约束方程选取了3 个独立的准速度 σ=[σ1,σ2,σ3]T=[θ,δ,φr]T，然后利用Gibbs-Appell 方法推导得到了自行车的约化动力学方程

(1)各项系数只与σ1和σ2有关，即mij=mij(σ1,σ2),,Pi=Pi(σ1,σ2)；

(2)mij(σ1,σ2)=mji(σ1,σ2)；

(4)P3(σ1,σ2)≡0.

然而，这些性质(除了性质(2)是平凡的)很难根据Gibbs-Appell 方程或其他经典力学方程如Lagrange 方程的结构去严格地证明.事实上，这些性质来源于自行车系统内在的对称性，而要阐述对称性，则需要用到几何力学中的一些基本概念[46].对于平面自行车系统，在考虑两个几何约束方程(1)之后，其构型空间是一个七维微分流形，同时受到4 个非完整约束作用，这就形成了该构型空间上的一个约束分布[20].自行车运动在构型空间切丛上形成的相轨迹始终限制在该约束分布之上.取李群G=SE(2)×S1×S1为对称群，它是3 个李子群SE(2)，S1和S1的直积，分别代表自行车后框架的平面刚体运动、后轮的转动以及前轮的转动.在G作用下，自行车系统具有Lagrangian 不变性以及约束分布不变性[20].可以看到，在构型空间每一点处，约束允许的速度子空间与群轨道切空间之和张成了该点的切空间，同时二者交集为一维子空间.在这基础之上，可以应用所谓的非完整Lagrange-d'Alembert 方程以得到约化约束空间上的约化动力学方程[20,47-48].该方程本质上与Gibbs-Appell 方程是等价的，但其考虑了对称性，结构较Gibbs-Appell 方程更为清晰.利用该方程，前述性质很容易得到证明，其中性质(1)是Lagrangian 不变性和约束分布不变性的直接结果，而性质(3)和(4)除了这两个不变性之外，主要来源于前述两个子空间的交集是一维的.

当自行车在空间曲面上运动时，接触点处约束方程的数学描述更为复杂.这对动力学方程组的约化带来了很大的困难.Xiong 等[27]利用第一类拉格朗日方程[49]得到了自行车在曲面路面的运动方程，对应为一组微分代数方程

其中，F(q,q)与自行车的重力和惯性力有关，A=是完整约束和非完整约束的系数矩阵，是前后车轮在接触点处与完整约束和非完整约束相对应的6 个拉格朗日乘子.

3 自行车的相对平衡点及其稳定性

目前关于自行车的相对平衡点及其稳定性的研究主要有3 种情况：

(1)自行车在水平面上的匀速直线运动；

(2)自行车在水平面上的匀速圆周运动；

(3)自行车在旋转曲面上的匀速圆周运动.

为研究以上动态运动的稳定性，首先需要通过系统的动力学方程，得到其对应的平衡解，然后在平衡位置附近建立线性化动力学模型，通过分析该线性系统的Lyapunov 稳定性问题，确定系统相对平衡点的稳定性.历史上，人们对Whipple 自行车做匀速直线运动情况的稳定性问题研究的较多[17,50-51].Carvallo[9]和Whipple[10]于1899 年分别独立推导出自行车的线性化方程组，并针对具体的模型参数计算了自行车自稳定的速度范围.文献[11，15，17，42]相继针对Carvallo-Whipple 自行车构型，分别采用Lagrange 方法和矢量力学方法推导出了线性化方程.Schwab[42]与Meijaard[17]还利用多体系统动力学计算软件分别求出了线性化方程的系数，并与手算结果进行了对比验证.Lennartsson[52]和Aström[53]利用多体动力学软件Sophia 得到的非线性动力学方程在平衡位置附近进行线性化处理，也得到了相应的线性化方程组.Basu-Mandal[18]利用对带乘子的Lagrange 方程线性化，Moore[25]对Kane 方法得到的非线性动力学方程线性化.以上方法均得到了自行车做匀速直线运动时正确的线性化方程组，它具有如下的标准格式

其中，η=[θ,δ]T，f=[P1,P2]T，ω0是自行车匀速直线运动时后轮的角速度.其他系数矩阵M,K0,C和D是与自行车几何和惯量参数相关的常数矩阵.具体表达式见文献[28]的附录.

近些年来，学者们发现自行车的自稳定性与自行车的所有参数相关[26,54]，在经典Whipple 自行车模型基础上还应进一步考虑非理想因素，如可以考虑车轮的滚动摩阻效应[55]；前车架系统对稳定性的影响[36,49,51,56-58]，其中前车架系统包括车头前立柱、前叉和前轮，在该模型中考虑了前车架系统的相对结构运动导致的水平速度自由度；自行车轮胎特征、路面梯度以及阻尼对稳定性也有影响[36,40,59-61]；自行车车轮的侧面滑动影响[11]；还有学者考虑了侧面风对自行车自稳定性影响[62].

自行车匀速直线运动的稳定性可以根据线性化动力学方程的特征根进行判断，Meijaard等[17]得到了Whipple 自行车在运动速度0～10 m/s 范围内动力学方程(7)的特征根分布情况.根据特征根的分布，可以将自行车匀速直线运动受到小扰动时的状态分为weave 运动，capsize 运动和caster 运动.其中，自行车稳定直线运动的速度限定在4.3 m/s

除此之外，也有学者对自行车在水平面上的匀速圆周运动进行了研究.Psiaki[13]和Franke 等[63]分别得到了一族匀速圆周运动解，而Lennartsson[52]和Aström 等[53]则得到了三族匀速圆周运动解.Basu-Mandal 等[18]首次给出了Whipple 自行车在水平路面上的全部匀速圆周运动解族，他们一共得到了4 个解族，而其中的两个解族在某种意义上可以合并为一个解族.该结果得到了Wang 等[33]和Xiong 等[28]的验证.

为了得到自行车在水平面上的匀速直线运动和匀速圆周运动的统一描述，并进行更加严格的稳定性分析，Xiong 等[28]基于Gibbs-Appell 方程对这两种运动进行了研究.这两种运动在某种意义上都可以看作是自行车动力系统的平衡点.注意到σ3=φr是循环坐标，可以将式(5)改写成一阶常微分方程组的形式

其中x=.设x0是方程(8)的平衡点，即其满足f(x0)=0，则必有x0=(θ0,δ0,0,0,ω0)，其中θ0,δ0,ω0是3 个常数.根据前文提到的方程(5)系数所具有的性质，这3 个常数只需要满足两个代数方程.因此，在这两个代数方程独立的情况下，这样的平衡点不是孤立的，而是以单参数解族的形式存在的.进一步，平衡点可以分为两类.第一类平衡点为x0=(0,0,0,0,ω0)，其中ω0为任意常数，这样的平衡点代表自行车的匀速直线运动.第二类平衡点的θ0,δ0不同时为零，这样的平衡点代表自行车的匀速圆周运动.利用Whipple 自行车的经典参数[17]，Xiong 等[28]得到了四族单参数圆周解，与Basu-Mandal[18]及Wang 等[33]的结果完全一致.

对于一些特殊的动力学系统[64]，其平衡点的Lyapunov 稳定性可以通过构造Lyapunov 函数，即采用Lyapunov 直接方法[65]进行判断.但是，对于自行车系统，其动力学方程组系数的复杂性，导致相应的Lyapunov 函数不易被求解.因此，Xiong 等[28]基于Lyapunov 间接方法，即根据系统(8)在平衡点处Jacobi 矩阵的特征根分布来判断平衡点的Lyapunov 稳定性，并且首次给出了自行车匀速直线运动和匀速圆周运动的非线性稳定性分析.事实上，由于x0是非孤立的，可以证明该点的Jacobi 矩阵必有一个零特征根，这意味着x0是非双曲的平衡点.因此，x0的非线性稳定性与线性稳定性并不直接等价[66].Xiong 等[28]结合中心流形定理[67]，第一次严格地证明了如果平衡点x0处的Jacobi 矩阵其余4 个特征根均具有负实部，则x0是稳定的，但不是渐近稳定的.

对于自行车在旋转曲面上匀速圆周运动的稳定性，Xiong 等[27]亦开展了相关研究.首先，他们给出了自行车在旋转曲面上圆周解的定义.该定义给定了广义坐标q，曲面曲线参数p以及Lagrange 乘子A随时间变化的规律，其中一共有21 个待定常数，由微分代数方程组(6)确定.为了消除时变变量，他们引入了一系列旋转变化，将圆周解等价转化为了一个新的微分代数方程组的平衡解，其中一共有21 个待定常数.同时他们发现，这21 个常数只需要满足20 个代数方程.因此与在水平路面上运动的情形一样，这样的平衡点一般来说不是孤立的，而是以单参数解族的形式存在的.他们以旋转抛物面ζ=α(ξ2+η2)，α ≥0 为对象，在不同的曲面系数下都找到了4 个单参数圆周解族.这些圆周解族的结构随着旋转曲面参数的变化而演化，其中包含了不同解族之间的相互碰撞、分裂等现象.

旋转曲面上匀速圆周运动的稳定性指的是相应平衡解的Lyapunov 稳定性.在平衡点处将新的微分代数方程组线性化，并利用线性约束方程消去非独立的变量，可以得到一个7 维的线性系统.同样，平衡点的非孤立性意味着该线性系统必有一个零特征根.结合中心流形定理可以证明，如果线性系统的其余6 个特征根都具有负实部，则该平衡点是稳定的，但不是渐近稳定的.数值计算结果发现，随着旋转抛物面系数α 增大，稳定圆周解的区间越来越小.当α 增大到某个临界值时，所有的圆周解都是不稳定的.

自行车的自稳定性研究除了基于解析模型的定量分析外，还可以从物理角度得到更直观的认识[26,68].历史上，对自行车稳定性的观点主要有3 种：一是以Appell[69]、Timoshenko 和Young[8]为代表的离心力观点，即前轮朝车体倾斜方向转动使车体转为曲线运动，所产生的离心惯性力平衡了重力的倾覆作用；二是以Klein 和Sommerfeld[70]为代表的陀螺效应观点(如图2 所示)，即车轮旋转产生的陀螺力矩在自稳定现象中起到了重要作用；三是以Grammel[71]为代表的“脚轮效应”观点(如图3 所示)，即前轮尾迹c是影响自行车稳定性的关键因素，当c>0 时，前轮的法向支撑力产生对前叉转轴的力矩，推动车把朝倾斜方向转动，使离心力效应的稳定作用自动实现.Jones[72]安装了一个与前轮并列反向旋转的副车轮，以消除陀螺力矩作用，用实验否定了陀螺效应，验证了脚轮效应.刘延柱[73]通过简化的力学分析对这三种效应都作出了解释.在另一篇文献[74]中，他经过估算发现在受控情况下，车轮的陀螺效应仅占稳定性因素的3%左右.

图2 前轮的陀螺效应Fig.2 Gyroscopic effect of the front wheel

图3 车轮脚轮效应Fig.3 Caster effect of the front wheel

但是，Kooijman 等[26]则否认了陀螺效应和脚轮效应是保证自行车实现稳定性的必要条件.他们设计了一款“TBS 自行车”(如图4 所示)，该“自行车”同样包括车身、前叉、前轮和后轮4 个刚体，但是车身和前叉简化为各自带有集中质量的杆件.在该实验自行车中，一方面通过在前轮附近利用反向旋转的副轮彻底消除陀螺效应，另一方面前叉转轴延长线与地面的交点位于前轮触地点的后方，与Jones[72]的自行车正好相反，以此消除脚轮效应.通过实验发现该“自行车”仍然可以实现自稳定，因此，他们认为陀螺效应和脚轮效应不是自行车实现自稳定的必要条件，但也不能否认这两个条件对自稳定的影响效应.他们同时验证了自行车的其他参数，如前车架质量分布，各个几何参数的相互耦合作用等等，均对自行车自稳定性有重要影响.刘延柱[75]认为Kooijman 等[26]的自行车并非无脚轮效应，他通过理论分析发现，由于其重心高度远大于普通自行车，前轮的侧向摩擦力的正面作用超过法向支撑力的负面作用，由此推动前叉朝脚轮效应正确方向转动.

图4 TBS 自行车试验模型[26]Fig.4 The experimental TBS bicycle[26]

为了明确自行车各个参数对稳定性的影响程度，Xiong 等[28]对自行车匀速直线运动状态下各物理参数对稳定性的敏感性进行了理论分析.他们发现轴距w越大，自行车在水平面上作匀速直线运动时的自稳定区域越小；前轮尾迹c越大，自稳定区域越大，这与Jones[72]的观点基本一致；前叉转轴倾角λ 影响相对较小，意味着该参数对自稳定具有很好的鲁棒性；前框架和后框架质心在惯性系中的水平位置对自稳定具有显著影响，这与Kooijman[26]的观点一致；而对于自行车的车轮参数而言，增大前车轮半径，减小后车轮半径则有助于提高自行车的自稳定性；与陀螺效应有关联的前轮惯量参数对自稳定起相反作用.另外，Doria[35]一方面基于方差的敏感度分析方法对自行车不同几何参数对稳定指数的影响进行了分析，另一方面利用相关分析方法，根据不同的相关系数判断几何参数对自行车系统整体稳定指标的影响规律.Tak 等[76]基于自行车线性动力学模型，研究了不同设计参数对自行车偏转速度和车把摇摆速度的灵敏度，结果发现在自行车的25 个参数中有7 个参数对稳定性影响明显.

4 自行车实验研究及控制

为了验证相关理论分析的正确性，许多学者开展了关于自行车的动力学实验研究.Roland 和Lynch[77]实验分析了自行车的横向动力学以及轮胎特征参数，将实验结果与他们所提出的理论模型[78-79]进行对比验证；Miyagishi[80]设计了自行车机器人，用于实现无人骑行自行车的稳定性运动.

为研究自行车的自稳定特性，Kooijman 等[81]在自行车上安装了陀螺仪装置及速度传感器(如图5 所示)，测量自行车运动过程中的车轮角速度、车身偏航速度、车把转角以及后轮转角等相关运动学量，并反演自行车线性化方程的特征根，验证匀速直线运动时的稳定性理论结果.Escalona 等[82]设计了类似的自行车动力学测试系统.通过测量自行车运动过程中关键点的位移、速度和加速度量，验证自行车多体动力学系统模型和的正确性和有效性.

图5 自行车的实验测试系统[85]Fig.5 The instrumented bicycle system[85]

近年来，人们也开展了一系列关于自行车控制的研究.Lowell 和McKell[83]综合考虑了前文所述的自行车的3个稳定因素，将驾车人对前叉的控制用简单的比例规律替代，建立数学模型计算，却得出不稳定的结果.刘延柱[84]对控制规律作了一些修改，即考虑前叉转轴的阻尼因素，且控制规律中除车体倾斜的角度信息以外，还计入了角速度信息，这样计算结果即转变为渐近稳定.有的学者通过在自行车上增加外置平衡装置，如飞轮[85-86]和平衡块[87]实现自行车的稳定持续前行.Tanaka[88]以车把转角为控制变量，采用PD 控制器研究自行车的稳定行驶.Keo[89]提出了outputzeroing 控制器的方法实现自行车在速度为零时的直立状态的稳定性；Owczarkowski[90]和Anjumol[91]采用LQR 控制方法实现自行车的直线跟踪控制；Toshiyuki 和Zhang 等[92]采用带扰动观测器的滑模控制方法研究运动稳定性，Chen 等[93]采用滑模控制器研究自行车的车轮转动角跟踪控制问题.除以之外，有些学者[23,94-98]以电动摩托车为对象，开展了相关动力学的实验研究工作，为未来新型自行车的研制和实验提供了必要的借鉴资料.

5 结论与展望

本文较详细地总结了自行车自发明以来，不同学者对其动力学建模、稳定性分析以及实验验证的研究工作.作为典型的多刚体动力系统，自行车简单的机械结构蕴藏着丰富的动力学现象，在不同路面行驶过程中受到的复杂完整约束和非完整约束导致了其动力学模型的复杂性.为了明晰自行车的自稳定性本质，不同学者采用不同方法分析了不同几何参数和质量参数对稳定性的影响机理，而且通过实验验证了所提出的不同模型的正确性和有效性.研究自行车的稳定性不仅可以探索这个百年来至今都未有明确解释的“自行”问题，而且还可以更深入研究人作为主动控制输入时，人、机与环境之间的复杂耦合关系，通过定量分析和物理直观解释对自行车的自稳定有进一步的认识，对未来真正无人自行车的实现提供一定的参考价值.