与点类似,欧几里得平面上的向量也可以表示成两个数字坐标,它们可以转化成向量的大小和方向。例如,向量(3,5)可以理解为沿横轴正方向移动3个单位,沿纵轴正方向移动5个单位后的总位移。图4-3描绘了欧几里得平面上的向量 。
许多物理量都是向量:需要用大小和方向来完全定义。例如,速度、加速度和力都是向量。由于向量很常见,所以让我们创建一个类来表示它们。
图4-3 平面上的一个向量
右击geom2d软件包并选择New(新建)→Python File(Python文件)。将其命名为Vector,并单击OK(确定)按钮。然后输入清单4-6中的代码。
清单4-6 Vector类的代码
Vector类的代码与Point类相似。只不过坐标被命名为u和v,而不是x和y。这是一种约定,避免在无意中混淆点和向量。
在继续学习之前,让我们重构Point类中的__sub__方法,使它返回一个向量。回想一下,两点相减 P-Q 会得到一个从 Q 到 P 的向量。将point.py文件修改成清单4-7中的代码。
清单4-7 重构Point类的__sub__方法
我们将在4.4.3节仔细研究这个操作方法,在那里我们会使用这个操作来创建向量。
现在让我们为Vector类创建一些有用的方法。
与点一样,向量加减是常见操作。例如,对代表两个力的向量求和可以得到两个力的和(也是向量)。
在__init__方法之后,输入清单4-8的代码。
清单4-8 向量的加法和减法
在__add__和__sub__方法中,我们都创建了新的向量实例,来保存投影的加或减的结果。
图4-4描绘了两个向量 和 的加减法运算。注意, 为何可以理解成 和- 的和。
你可能会好奇,我们是否会为其他操作符构建类似的代码。点和向量的加减非常直观,但像__mul__(用于重载乘法操作)这样的操作符要复杂一些。向量的乘法分为点乘、叉乘和标量乘法(对应向量的缩放)。因此,我们不使用操作符,而是直接将这些操作创建为带有描述性名称的方法:scaled_by、dot和cross。
我们将从向量缩放开始。
图4-4 两个向量相加和两个向量相减
要缩放(scale)一个向量 ,你可以将它与标量 k 相乘,从而放大或缩小这个向量。数学上,标量乘法用公式(4.2)表示:
让我们在Vector类中创建一个缩放方法。在__sub__方法下面输入清单4-9中的代码。
清单4-9 向量的缩放
在上述代码中,只返回一个新Vector,其 u 和 v 等于原向量的 u 和 v 与参数factor——缩放因子——相乘的结果。
使用缩放方法,我们可以实现另一个操作,即用给定向量 对 P 点进行 k 次移动,见公式(4.3):
图4-5 用向量 移动 P 点 k 次(图中 k =2)
示意图见图4-5。
因为移动的对象是点,因此让我们在Point类中编写对应的代码(见清单4-10)。
清单4-10 用向量 移 动 P 点 k 次
该方法接收两个实参:向量vector和标量times。向量缩放的次数等于times,最终得到总位移。例如,向量(3,5)按times=2缩放将得到(6,10)。注意,参数times的默认值为1,因为通常输入的vector长度就是实际长度,不需要进行缩放。函数返回的点坐标等于原来的点坐标加总位移向量的坐标。
让我们尝试在Python的shell中移动点。重新启动控制台,以清除之前导入的Point类和Vector类,并输入以下内容:
你可以用计算器确认该运算结果是否正确。
向量的范数(norm)是指它的长度。单位范数的长度为一个单位。拥有单位范数的向量在确认向量方向时非常有用,因此,我们经常会想知道一个向量的范数是否为单位范数(它是否是单位向量)。我们也经常需要归一化(normalize)一个向量:方向不变,长度变为1。二维向量的范数由公式(4.4)给出:
让我们创建一个返回向量范数的特性,然后创建另一个特性来检验某向量是否为单位向量。两者的代码都包含在清单4-11中。
清单4-11 计算向量的范数
norm特性的返回值完全符合公式(4.4)的定义。为了确认向量的范数是否为1,我们使用数值比较函数is_close_to_one并将向量的范数传递给它。
我们将创建另外两个重要的方法:第一个将向量 归一化,生成一个方向与 相同但长度为1的向量 ;第二个将向量缩放到给定长度。对向量进行归一化后得到的向量,我们称之为单位向量(unit vector)或单位方向向量(versor),可以使用公式(4.5)得到:
计算结果是长度为1的向量。将该向量乘以标量 k ,得到向量 ,它的方向和原向量相同,长度等于标量的值,如公式(4.6)所示:
清单4-12是上述公式对应的代码。
清单4-12 计算单位长度或指定长度的向量
为了归一化向量,我们将它乘以它范数的倒数(相当于用向量的长度除以范数)。当需要将向量缩放到指定长度时,我们先将其归一化,然后乘以指定的标量。
你可能已经意识到,我们从不改变任何对象的属性,而是创建并返回一个新的Point或Vector实例。例如,为了归一化向量,我们也可以使用清单4-13中的代码。
清单4-13 原地归一化向量
调用该方法将导致原地归一化(normalization in place),即直接更改当前对象的属性。原地归一化更快,需要的内存更少,但也更容易出错。你的程序更改的对象很可能被其他程序调用,而后者不希望对象变更。这类错误非常难以发现,需要进行大量调试。此外,使用不可变数据的程序更易于理解和解释,因为你不需要考虑对象的状态随时间的变化。
看看下面的代码。它以与前面类似的方式创建了normalize方法,但它包含一个小错误。在这种情况下,归一化将会产生错误的结果。你能找出原因吗?
这个问题很棘手。第一行修改的self.x属性,会在第二次调用self.norm属性时用到。两次调用self.norm会产生不同的结果。这就是必须将self.norm的值存储在变量中的原因。
当对象的数据量很小时,最好避免可变对象。这样,程序在并发执行时也能够正确运行,并且代码也更易于理解。将可变性降到最低将使你的代码更加稳健,在本书中,我们会尽可能地坚持这个原则。
请留意方法的命名约定。在调用时改变对象状态的方法命名如下:
normalize 将向量原地归一化
scale_by 在原地缩放向量
创建新对象作为结果的方法命名如下:
normalized 返回一个新的归一化向量
scaled_by 返回一个新的缩放后的向量
接下来,我们将在Vector类中实现点乘和叉乘。这些简单的乘法将为一些有用的操作(如计算两个向量的夹角或检验是否垂直等)打下基础。
两个向量 和 的点乘(dot product)会得到一个标量,它可以反映两个向量方向的差异。在二维空间中,点乘如公式(4.7)所示,其中 θ 为向量的夹角:
根据两个向量的方向不同,点乘的值也不同,见图4-6。图上有一个参考向量 和另外三个向量: 、 和 。一条垂直于 的直线将空间分成两个半平面。向量 在直线上,因此 和 的夹角 θ 等于90°。而cos(90°)=0,因此 。垂直向量的点乘为零。向量 所在的半平面和 相同,因此, 。最后, 在与 相对的半平面上,因此, 。
图4-6 与不同方向的向量的点乘结果
公式(4.7)可以直接转化为点乘的代码。在Vector类中,输入清单4-14的代码。
清单4-14 向量的点乘
在继续讨论叉乘之前,让我们停一停,分析点乘的一个应用:求向量在指定方向上的投影。
当点乘的其中一个向量是单位向量时,结果是一个向量在另一个向量上投影的长度。为了明白原因,需要参考公式(4.7)。给定向量 和单位向量 ,它们的点乘如下:
其中, 正好是 在 方向上的投影。这对于计算指定方向的投影非常方便,我们可以用它来获得桁架构件上的力的轴向分量,如图4-7所示。在这种情况下,我们只需要计算 来得到轴向分量 。
让我们将这个操作创建为一个新的方法。将清单4-15中的代码输入Vector类中。
图4-7 力 在桁架构件轴向方向 上的投影
清单4-15 一个向量在另一个向量上的投影
请注意,实参direction可能不是单位向量。为了确保算法有效,需要先将其归一化。
两个三维向量的叉乘(cross product)会得到一个垂直于这两个向量所在平面的新向量。向量的顺序很重要,它决定了结果向量的方向。可以使用右手法则得到叉乘的方向。注意,叉乘不满足交换律: ,见图4-8。
在三维空间中,叉乘可以用公式(4.8)来计算:
图4-8 叉乘不满足交换律
在二维平面上,由于所有向量都位于同一平面,因此,叉乘会得到一个垂直于该平面的向量。将 u z = v z =0代入公式(4.8)中,即可得到这个结论:
因此,在二维应用中,叉乘的结果可以输出一个标量,也就是结果向量的 z 坐标。这个坐标也可看作结果向量的长度。由于 x 和 y 坐标为零,因此只需要计算 z 坐标。设向量 与 的夹角为 θ ,应用公式(4.9)可以得到二维向量的叉乘结果:
让我们来创建叉乘算法。输入清单4-16中的代码。
清单4-16 向量的叉乘
二维向量叉乘的一个重要应用是确定角度的旋转方向。从图4-8中可以看到, ,因为从 到 的角度为正(逆时针)。相反,从 到 的角度为负,因此叉乘 。最后,平行向量的叉乘为0,这很显然,因为sin 0=0。让我们仔细研究一下,并在Vector类中创建检验两个向量是否平行或垂直的方法。
使用点乘和叉乘,很容易检验两个向量是平行还是垂直。清单4-17包含了具体的代码。
清单4-17 检查向量是否平行或垂直
检验两个向量是否相互平行,只需要确认它们的叉乘是否为零。同样,检验两个向量是否垂直,只需要检验点乘是否为零。注意,我们使用了函数is_close_to_zero来说明计算中浮点数比较的困难。
计算两个向量之间的夹角可以借助点乘公式:
将左边的点乘除以另一边的范数积,再取反余弦,即得到公式(4.10):
该公式只能计算出角度的大小。如果想知道方向,则需要使用叉乘。角度的符号可以用如下函数:
其中,符号函数sgn定义如下:
要理解为什么公式(4.10)只能得到角度的大小,我们需要记住余弦函数的一个重要性质。回忆一下基本的几何理论,单位向量夹角的余弦正是其水平投影的值。通过图4-9中的单位圆可以看到,两个角度相反的向量(角度的和等于零),其对应的余弦值相同。换句话说,cos α =cos( -α ),这意味着余弦函数会消除掉角度的符号。因此根据点乘的值确定角度的符号也就不可能了。
图4-9 相反角度的余弦值相等
我们创建的许多应用程序,都同时需要角度的大小和符号,在叉乘的帮助下,我们可以恢复这些信息。让我们创建两个新方法,一个方法计算角度的绝对值(针对那些只需要角度大小的情况),另一个方法包含符号。在Vector类中输入清单4-18中的代码。
清单4-18 计算两个向量的夹角
第一个方法(angle_value_to),使用公式(4.10)计算self与other的夹角。我们首先得到点乘值,然后除以范数的乘积。取结果的反余弦,即得到角度值。第二个方法(angle_to),返回带叉乘符号的角度值。Python中的math.copysign(x,y)函数返回x的大小和y的符号。
让我们在控制台中尝试这两个方法。重新加载控制台,然后输入以下内容:
仅供参考,角度值0.78539…对应π/4 rad(45°)。
假设我们有一个向量,想要通过对其旋转一定角度来创建一个新的向量。
如图4-7所示,假设一根杆受到外力作用,我们希望知道力在垂直于杆的方向上的投影,也就是力的剪切分量。为了计算这个投影值,我们首先需要找到与杆 u ˆ的方向垂直的向量,这可以通过旋转该向量π/2弧度得到,如图4-10所示。
旋转变换不影响长度,因此该操作会保留原向量的长度。假设向量旋转的角度为 α ,使用公式(4.11)可以计算旋转后的向量:
图4-10 将杆的方向向量旋转π/2弧度
对应的Python代码如清单4-19所示。
清单4-19 旋转一个向量
rotated_radians函数返回一个新向量,这是将原向量旋转给定弧度后的结果。根据我们的不可变编程原则,我们不会改变原向量,而是返回一个应用旋转后的新向量。
角度π/2 rad(90°)对于旋转向量非常有用。使用π/2 rad,我们可以得到垂直于原向量的新向量。
为了避免重复编写v.rotated_radians(math.pi/2),我们可以在Vector类中定义一个新方法。因为cos(π/2)=0,sin(π/2)=1,公式(4.11)可以简化成如下形式:
我们将这个方法命名为perpendicular。Python代码见清单4-20。
清单4-20 计算垂直向量
另一个我们经常用于旋转的角度是π(180°)。向量旋转π可以得到一个与之共线但反向的向量。因为cos(π)=-1,sin(π)=0。公式(4.11)可以简化成如下形式:
将对应的方法命名为opposite。Python代码见清单4-21。
清单4-21 计算相反向量
perpendicular和opposite这两个方法涉及的知识我们都学过,我们当然也可以继续使用rotated_radians。不过,这两个方法很方便,我们会经常使用它们。
为了在 x 轴和 y 轴上投影一个向量,我们需要使用向量角度的正弦或余弦,如图4-11所示。
图4-11 向量在坐标轴上的投影
这些可以用来计算本书第五部分介绍的桁架结构在全局坐标中的刚度矩阵。杆的刚度矩阵计算是在相对坐标系中进行,其 x 轴是杆的轴向,但我们需要将这个矩阵的方向投影到全局坐标的 x 和 y 轴上,从而构建结构在全局坐标系下的方程组。
如果Vector类没有这两个属性,类的实例可以获取杆的角度,然后计算正弦或余弦。即使这是完全可以接受的,它也需要一些操作:首先计算角度,然后计算正弦或余弦。但你也知道,利用它们的数学定义,我们可以更有效地计算正弦值和余弦值。
假设向量 的范数为 ,投影分别为 u 和 v 。正弦和余弦的计算方法如下:
让我们将其创建为Vector类的属性。输入清单4-22中的代码。
清单4-22 向量方向的正弦和余弦
有了前面的表达式,这个代码非常简单。让我们添加最后几步,来完成我们的Point类和Vector类。