《机器人与数字人：基于MATLAB的建模与控制》——2.4节对偶数、对偶向量及代数形式

本节书摘来自华章社区《机器人与数字人：基于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来简要描述。

在表21中，每种类型的组合数构成一个加法群，但只有复数满足乘法群，因此它能形成代数域。而双数和对偶数仅能形成乘法半群，因为双数在a=b时以及对偶数在a=0时逆运算不存在。因此，这两种代数只能形成代数环而不是域，并分别称为双环和对偶环。
值得注意的是，第三种类型是著名的复数域，已经被当作最流行的数学工具广泛使用。由于双数代数既复杂又没有意义，所以很少使用。而对偶数在过去研究中被忽视了。事实上，对偶数是非常独特和有吸引力的简单代数式，尤其是在三维空间变换应用中。
上述定义的三种结合代数也可以线性映射到2×2矩阵中，并保留代数域(或者环)运算［4］。也就是说，令

Γ(a+jb)=ab
j2ba=C（214)

这样，对于双数有Γ(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
很容易验证关于两个双数的映射是有效的。总之，这种映射Γ满足
Γ(c1c2)=Γ(c1）Γ(c2）
其中，“”可以是加法或者乘法。因此，映射Γ是保留代数结构的同态。此映射也适合三个组合数的逆运算求解。
用式(214)对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， 1118］。假设在三维空间中两条直线之间的最短距离为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)(215)

其中，f′(a)和f″(a)分别是在a点的一、二阶导数。值得注意的是，以上对偶环的一般泰勒级数展开式生成了新的对偶函数。基于这样简洁的展开，很容易看到

ea+b=ea·eb=ea+bea(216)

sin(α+d)=sinα+dcosα(217)

cos(α+d)=cosα-dsinα(218)

更重要的是，一般展开式(215)可以直接使用，通过对偶数微分来确定任何光滑函数的一阶导数，这样既没有数值错误也不需要符号推导。正如我们都熟悉的数值计算一样，集成的计算机函数计算可以很容易达到至少四阶精度，如众所周知的龙格库塔算法。然而，在缺少先前采样数据存储的情况下，要计算函数导数数值往往很难达到令人满意的高阶近似值。现在，MATLAB等计算软件可以执行内部子程序来完成对偶数微分运算，就像在MATLAB中已经建立的复数计算一样。这样我们就可以对任何函数进行精确的数值求导计算，在计算机中通过方程(215)对瞬时点求导，且不需要经过符号推导。
作为算例分析，令

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=－2098×10－3
这正是f′(3)的值，因为在x中有意设定b=1。为了验证数值结果，可以用符号表示
f′(x)=(1－2x2)exp(－x2)

将x=3代入以上导数，得到完全相同的答案，即-17e-9=-2098×10-3。
在后面的章节中，将需要计算更多复杂函数的导数，包括高维向量场和矩阵域。在编程进行数值计算时，可能会很困难或者无法人为找到其符号导数。在这种情况下，由于对偶环拥有以上属性，相对其他算法它将以绝对的优势作为数值求导解决方案［1］。

此外，假设用q+替换方程(215)中的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)基于链式法则求解的时间导数。方程(219)中新的性质将会扩展方程(215)的特殊优势，用来直接对大规模函数进行时间数值求导，在以后的机器人学和数字人建模应用中将会更多应用。当然，这种导数也可以推广到关于任意时间对任何参数的求导，这也将在以后的研究中更多地展现出来。
现在学习一个实例。假设机器人连杆坐标系相对固定基座的定向旋转矩阵如下：

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=05rad/s，2=-04rad/s和3=08rad/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(220)

事实上，通过实方阵在线性代数中的某些主要性质，此方程很容易实现变换。首先，实方阵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)
利用方程(220)，获得对偶方阵A∧行列式新的性质如下：


det(A∧)=det(A+B)=det(A)·etr(A-1B)=det(A)·(1+tr(A-1B))(221)

这意味着对偶矩阵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=0333300
01667－050000
00417－03750－02500
有
37
A－1B=0666710000－13333
033330－06667
－0916702500－09167
因此，tr(A-1B)=-025，det(A∧)=24(1-025)=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)(222)

通过依次增加这样的三个对偶向量，最终获得用来唯一表示坐标系相对于参考坐标系的位置和姿态的3×3对偶旋转矩阵：

R∧=R+S
(223)
38

换句话说，现在把每个坐标轴都当作线向量，而不是点向量，并通过式(222)表示对偶向量。
现在很清楚的是，每个ri都是单位向量，并且瞬时向量si=p×ri与ri垂直，即rTiri=1和rTisi=0(i=1,2,3)，那么R∧中的每个对偶列向量r∧i很明显都是单位旋量，而在式(223)中对偶旋转矩阵R∧为特殊正交对偶矩阵，即R∧TR∧=R∧R∧T=I，也意味着RTR=I且RTS+STR=0或者RST+SRT=0。换句话说，RTS必须是一个3×3的斜对称矩阵，且迹tr(RTS)=0。将这个零迹代入到行列式方程(221)中，并注意到式(221)的实数部分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。通过回顾式(222)定义如何增加三个对偶向量，并以此形成式(223)的对偶旋转矩阵，能得到

S=PR或SRT=P(224)
矩阵P=SRT正好是位置向量p的斜对称矩阵，即P=p×。这也表明如果一个坐标系的位置和姿态由式 (223)定义的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)。此外，根据式(224)，所求坐标系原点的位置向量可以直接通过以上的斜对称矩阵SRT=P产生，即p=(257)T。
李群中SO(3,d)的特殊正交对偶矩阵涵盖了对偶环及其变换，将在第4章中重新讨论和进一步研究，并且在机器人运动学建模和路径规划应用中，它是一个很重要而且很有用的理论。