本节书摘来自华章社区《机器人与数字人:基于MATLAB的建模与控制》一书中的第2章,第2.4节对偶数、对偶向量及代数形式,作者[美]顾友谅(Edward Y.L.Gu),更多章节内容可以访问云栖社区“华章社区”公众号查看
2.4对偶数、对偶向量及代数形式
在数学运算中,两个实数a和b可以与某一单位量j结合产生一个新的数c=a+jb。通常情况下,仅在下面三种不同的情况下,单位量j才可以单独定义为
j2=+1,j2=0,j2=-1
在过去,组合数的三种形式已经有相关的研究[9, 10, 12],可以通过表2-1来简要描述。
![2dfd675e4c6cc3e1e99a5d72aa7d33cd7c91f08e](https://yqfile.alicdn.com/2dfd675e4c6cc3e1e99a5d72aa7d33cd7c91f08e.png?x-oss-process=image/resize,w_1400/format,webp)
在表21中,每种类型的组合数构成一个加法群,但只有复数满足乘法群,因此它能形成代数域。而双数和对偶数仅能形成乘法半群,因为双数在a=b时以及对偶数在a=0时逆运算不存在。因此,这两种代数只能形成代数环而不是域,并分别称为双环和对偶环。
值得注意的是,第三种类型是著名的复数域,已经被当作最流行的数学工具广泛使用。由于双数代数既复杂又没有意义,所以很少使用。而对偶数在过去研究中被忽视了。事实上,对偶数是非常独特和有吸引力的简单代数式,尤其是在三维空间变换应用中。
上述定义的三种结合代数也可以线性映射到2×2矩阵中,并保留代数域(或者环)运算[4]。也就是说,令
Γ(a+jb)=ab
j2ba=C(214)
这样,对于双数有Γ(a+jb)=ab
ba,对于对偶数有Γ(a+jb)=ab
0a,对于复数有Γ(a+jb)=ab
-ba。 例如,如果给定c1=a1+jb1和c2=a2+jb2是两个复数且j2=-1, 那么,
Γ(c1+c2)=Γ(a1+a2+j(b1+b2))=a1+a2b1+b2
-(b1+b2)a1+a2
并且
Γ(c1·c2)=Γ(a1a2-b1b2+j(a1b2+a2b1))
这等价于
30
a1b1
-b1a1a2b2
-b2a2=a1a2-b1b2a1b2+a2b1
-(a1b2+a2b1)a1a2-b1b2
同样,如果c1和c2都是对偶数且j2 = 0,那么加法显然成立,而乘法为
Γ(c1·c2)=Γ(a1a2+j(a1b2+a2b1))
这也等价于
a1b1
0a1a2b2
0a2=a1a2a1b2+a2b1
0a1a2
很容易验证关于两个双数的映射是有效的。总之,这种映射Γ满足
Γ(c1c2)=Γ(c1)Γ(c2)
其中,“”可以是加法或者乘法。因此,映射Γ是保留代数结构的同态。此映射也适合三个组合数的逆运算求解。
用式(214)对2×2矩阵求逆,得到
C-1=ab
j2ba-1=adjCdetC=a-b
-j2baa2-j2b2
其中,adj C是C的伴随矩阵。对于复数j2=-1,使得
C-1=a-b
baa2+b2
这很明显等价于通过其共轭复数求得的复数的逆:
(a+jb)-1=a-jb(a+jb)(a-jb)=a-jba2+b2
与此相反,对偶数求逆且j2=0变为
C-1=a-b
0aa2
或者
(a+jb)-1=a-jb(a+jb)(a-jb)=a-jba2
这表明如果实数部分a=0,对偶部分b为非零,那么对偶数将不存在逆解。因此,对偶数集合仅能构成乘法半群,并且只能形成一个代数环——对偶环。
同样,对于双数且j2=+1,有
C-1=a-b
-baa2-b2
或者
(a+jb)-1=a-jb(a+jb)(a-jb)=a-jba2-b2
这也表明如果a=b,双数将不存在逆解,双数集合也仅能形成环而不是域。
尽管如此,对偶环具有许多独特和优越的代数性质,这使得它在空间变换应用中具有很大吸引力,尤其在线向量的表达方面。为了在下一节中能够明确区分对偶数和其他代数,令j重新用表示为a+b,且2=0。
2.4.1对偶环微积分
1873年,Clifford首先提出了对偶数[9],Study在1901年进一步研究了在欧几里得空间中如何表示对偶角[10, 1118]。假设在三维空间中两条直线之间的最短距离为d,这两条直线之间的夹角为α,则对偶角定义为
α∧=α+d
前提条件是2=0。这样定义对偶角的几何意义被认为是螺旋运动。也就是说,当用螺丝刀拧螺丝旋转α角度时,同时向下旋进的垂直距离为d。
就像复数域中的通用函数一样,光滑的对偶函数f(a+b)可以在a点根据标准泰勒级数展开,由于2=0,展开式化简为
f(a+b)=f(a)+bf′(a)+2b22!f″(a)+…
=f(a)+bf′(a)(215)
其中,f′(a)和f″(a)分别是在a点的一、二阶导数。值得注意的是,以上对偶环的一般泰勒级数展开式生成了新的对偶函数。基于这样简洁的展开,很容易看到
ea+b=ea·eb=ea+bea(216)
sin(α+d)=sinα+dcosα(217)
cos(α+d)=cosα-dsinα(218)
更重要的是,一般展开式(215)可以直接使用,通过对偶数微分来确定任何光滑函数的一阶导数,这样既没有数值错误也不需要符号推导。正如我们都熟悉的数值计算一样,集成的计算机函数计算可以很容易达到至少四阶精度,如众所周知的龙格库塔算法。然而,在缺少先前采样数据存储的情况下,要计算函数导数数值往往很难达到令人满意的高阶近似值。现在,MATLAB等计算软件可以执行内部子程序来完成对偶数微分运算,就像在MATLAB中已经建立的复数计算一样。这样我们就可以对任何函数进行精确的数值求导计算,在计算机中通过方程(215)对瞬时点求导,且不需要经过符号推导。
作为算例分析,令
f(x)=xexp(-x2)
我们希望计算函数f(x)在x=3处的一阶导数,即f′(3)的数值。令x=3+1=3+,那么
f(x)=(3+)exp(-(3+)2)=(3+)e-9-6=(3+)e-9(1-6)=(3-17)e-9
显然,取以上的对偶部分为
Du{f(3+)}=-17e-9=-2098×10-3
这正是f′(3)的值,因为在x中有意设定b=1。为了验证数值结果,可以用符号表示
f′(x)=(1-2x2)exp(-x2)
将x=3代入以上导数,得到完全相同的答案,即-17e-9=-2098×10-3。
在后面的章节中,将需要计算更多复杂函数的导数,包括高维向量场和矩阵域。在编程进行数值计算时,可能会很困难或者无法人为找到其符号导数。在这种情况下,由于对偶环拥有以上属性,相对其他算法它将以绝对的优势作为数值求导解决方案[1]。
此外,假设用q+替换方程(215)中的a+b,其中q∈瘙綆n可以为n维向量,且n=1,2,…;是关于时间的导数;函数f(·)可能是高维向量或矩阵(注意,在这种情况下,矩阵乘法维数必须一致)。那么,
f(q+)=f(q)+f(q)q=f(q)+(q)
(219)
其中(q)=ddtf(q)=f(q)q是函数f(q)基于链式法则求解的时间导数。方程(219)中新的性质将会扩展方程(215)的特殊优势,用来直接对大规模函数进行时间数值求导,在以后的机器人学和数字人建模应用中将会更多应用。当然,这种导数也可以推广到关于任意时间对任何参数的求导,这也将在以后的研究中更多地展现出来。
现在学习一个实例。假设机器人连杆坐标系相对固定基座的定向旋转矩阵如下:
R=c1c23s1c1c23
s1c23-c1s1s23
s230-c23
其中,s1=sinθ1,c1=cosθ1,s23=sin(θ2+θ3)以及c23=cos(θ2+θ3)为传统的简写方式。在某个时刻,机器人的三个关节角度的弧度θ1=π/6, θ2=π/2和θ3=0,并且给定关节速度1=05rad/s,2=-04rad/s和3=08rad/s。通过MATLAB编程,可以直接得到旋转矩阵R的数值结果:
现在,用θi+i取代每个关节角度θi,i=1, 2, 3,并且从最终数值结果R∧=R+S中取出对偶部分,可以直接从MATLAB工作窗口中看到输出结果:
根据式(219),此矩阵S是矩阵R的瞬时时间导数,即S=。
此外,如果此刻希望确定R关于θ2的数值偏导数,即S=Rθ2,那么仅需要在R中用θ2+1替换θ2,而其他两个关节角度维持不变。显然,这相当于在上述MATLAB程序中设置了关节速度2=1,1=0和3=0。最终显示如下结果:
这只是R关于θ2的瞬时偏导数。使用同样方法,也可以计算出数值偏导数Rθ1和Rθ3。
2.4.2对偶向量和对偶矩阵
对偶数定义可以扩展到两个向量或者两个矩阵组合,以此形成的对偶向量或者对偶矩阵如下:
a∧=a1+b1
a2+b2
a3+b3=a1
a2
a3+∈b1
b2
b3=a+b
或者对于两个同维矩阵A和B,有
A∧=A+B
由于很容易发现对偶数的高维扩展,它相对于复数和双数有特别简洁的性质,因此理应进一步详细讨论。
首先,两个向量的内积类似于两个对偶数乘法。令a∧=a+b,c∧=c+d,且每个a、b、c和d都是3×1向量。那么,a∧和c∧的点积为
a∧Tc∧=aTc+(bTc+aTd)
为了完成叉积a∧×c∧,令a∧的斜对称矩阵定义为
S(a∧)=0-a3-b3a2+b2
a3+b30-a1-b1
-a2-b2a1+b10=S(a)+S(b)
那么,
a∧×c∧=S(a∧)c∧=S(a)c+(S(b)c+S(a)d)=a×c+(b×c+a×d)
同样,若A、B、C和D为n×n方矩阵,两个对偶矩阵A∧=A+B和C∧=C+D点积类似于对偶向量:
A∧C∧=AC+(BC+AD)
并且毫无例外,
A∧a∧=Aa+(Ab+Ba)
此外,为了确定对偶方阵A∧=A+B的逆,必须首先假设实数部分A是非奇异的。那么,A∧的逆可以确定为
A∧-1=(A+B)-1=A-1-A-1BA-1
其很容易通过测试验证,
A∧-1A∧=(A-1-A-1BA-1)(A+B)=I+O=I
其中,I和O分别是n×n单位矩阵和零矩阵,所以它们的变换结果也一样,A∧A∧-1=I。
记得在线性代数中有一个著名的方程结合了实矩阵A的行列式和迹:
det(eA)=etrA(220)
事实上,通过实方阵在线性代数中的某些主要性质,此方程很容易实现变换。首先,实方阵A总可以对角化或者通过相似变换矩阵P的约旦分解如下:
P-1AP=λ1X
0λn
36
式中,λ1到λn是A的特征值,基于约旦分解的右上三角区域X可能含有一些1,但是左下三角区域的元素都由零组成。因此,
P-1eAP=exp(P-1AP)=eλ1X
0eλn
对于右上角X,无论X为何值,左下角所有元素值都始终为零,这保证了行列式
det(eA)=det(P-1eAP)=∏ni=1eλi=exp(λ1+…+λn)=etrA
现在,任何对偶方阵M应用指数函数的泰勒展开式为
exp(M)=I+M
这意味着如果A是非奇异的,且exp(A-1B)=I+A-1B,那么
det(A∧)=det(A+B)=det(A)·det(I+A-1B)=det(A)·det(eA-1B)
利用方程(220),获得对偶方阵A∧行列式新的性质如下:
det(A∧)=det(A+B)=det(A)·etr(A-1B)=det(A)·(1+tr(A-1B))(221)
这意味着对偶矩阵A∧的行列式是对偶标量,其实数部分正好是A∧的行列式的实数部分,而对偶部分是det(A)tr(A-1B)。显然,其运算条件必须是det(A)≠0。
例如,令一个3×3对偶矩阵为
A∧=3+23-4
1-2+0
-1+43-2-4+3
=300
1-20
-13-4+23-4
010
4-23=A+B
由于det(A)=24,则
A-1=0333300
01667-050000
00417-03750-02500
有
37
A-1B=0666710000-13333
033330-06667
-0916702500-09167
因此,tr(A-1B)=-025,det(A∧)=24(1-025)=24-6。
可以通过对偶数的加法和乘法直接计算A∧的行列式:
det(A∧)=-41-2+
-1+43-2+(-4+3)3+23
1-2+=24-6
因此,这两种不同方法计算的结果完全一致。
2.4.3单位旋量和特殊正交对偶矩阵
定义3一个对偶向量a∧0=a0+b0,如果‖a∧0‖2=a∧T0a∧0=1,那么它被认为是单位旋量。
定义4一个对偶矩阵R∧=R+S,如果R∧TR∧=R∧R∧T=I,I为单位实矩阵,并且det(R)=+1,那么它被认为是特殊正交对偶矩阵。在这种情况下,R∧∈SO(3, d),即对偶环属于特殊正交群。
根据定义3,由于a∧T0a∧0=aT0a0+2aT0b0,单位旋量必须满足aT0a0=1且aT0b0=0。实际上,定义表明了如何创建对偶向量的对偶部分,使它能够较好地表示线向量。众所周知,如果一个力f施加到一个旋转刚体上,则扭矩τ一般定义为τ=p×f,其中p是径向向量,其尾部固定在旋转中心上,箭头指向力的作用线。显然,这意味着该力f和扭矩τ总是互相垂直的。扭矩的定义启发我们可以采用类似的概念来定义对偶向量,可以唯一地表达坐标系相对参考坐标系的位置和方向。
在以往机器人学中,经常使用旋转矩阵R来表示坐标系的方向,以及用3×1(点)向量P表示坐标系的原点相对参考坐标系的位置。虽然这个坐标系的三个坐标轴可以通过旋转矩阵R相对应的列ri(i=1, 2, 3)来表示,每一个都是点向量,但仍然需要添加一个位置向量p来唯一完全地表示坐标系的位置和姿态。现在,定义一个“瞬时”向量si,同扭矩τ的定义一样,则有si=p×ri,其中R中的三个单位列向量ri都理解为力f。于是每个坐标轴的完整的对偶向量为
r∧=ri+si=ri+(p×ri)(222)
通过依次增加这样的三个对偶向量,最终获得用来唯一表示坐标系相对于参考坐标系的位置和姿态的3×3对偶旋转矩阵:
R∧=R+S
(223)
38
换句话说,现在把每个坐标轴都当作线向量,而不是点向量,并通过式(222)表示对偶向量。
现在很清楚的是,每个ri都是单位向量,并且瞬时向量si=p×ri与ri垂直,即rTiri=1和rTisi=0(i=1,2,3),那么R∧中的每个对偶列向量r∧i很明显都是单位旋量,而在式(223)中对偶旋转矩阵R∧为特殊正交对偶矩阵,即R∧TR∧=R∧R∧T=I,也意味着RTR=I且RTS+STR=0或者RST+SRT=0。换句话说,RTS必须是一个3×3的斜对称矩阵,且迹tr(RTS)=0。将这个零迹代入到行列式方程(221)中,并注意到式(221)的实数部分A现在是实数旋转矩阵R,且对偶部分B是旋转矩阵S,得到det(R∧)=det(R)=+1。这充分地证明了如下两个引理。
引理1如果一个3×3对偶矩阵R∧=R+S∈SO(3,d),那么R∧中的每一列或每一行都是单位旋量。
引理2一个3×3对偶矩阵R∧=R+S,当且仅当R∈SO(3)和SRT是一个斜对称矩阵时,它是特殊正交对偶矩阵,即R∧∈SO(3,d)。
根据引理2,对应R∧=R+S∈SO(3,d),SRT必须是斜对称矩阵,且使得RST=-SRT或者RST+SRT=O。通过回顾式(222)定义如何增加三个对偶向量,并以此形成式(223)的对偶旋转矩阵,能得到
S=PR或SRT=P(224)
矩阵P=SRT正好是位置向量p的斜对称矩阵,即P=p×。这也表明如果一个坐标系的位置和姿态由式 (223)定义的3×3的特殊正交对偶矩阵R∧表示,则坐标系的姿态完全由实数部分R描述,而坐标系的原点位置通过叉积运算P=SRT确定。
例如,给定一个3×3对偶矩阵
R∧=157
7-2-1
-51-2=R+S
39
根据引理2,此时不能判断它是特殊正交对偶矩阵,除非证明它的实数部分R和SRT满足R∈SO(3)并且SRT是斜对称矩阵。而事实上,
RRT=100
00-1
010100
001
0-10=I并且det(R)=+1
那么
SRT=057
7-20
-50-2100
001
0-10=0-75
70-2
-520
显然是一个斜对称矩阵。因此,这两个条件成立,这个R∧是一个特殊正交对偶矩阵,也就是R∧∈SO(3,d)。此外,根据式(224),所求坐标系原点的位置向量可以直接通过以上的斜对称矩阵SRT=P产生,即p=(257)T。
李群中SO(3,d)的特殊正交对偶矩阵涵盖了对偶环及其变换,将在第4章中重新讨论和进一步研究,并且在机器人运动学建模和路径规划应用中,它是一个很重要而且很有用的理论。