10行Python代码实现多种数学运算

发表时间: 2024-06-25 18:04

【CSDN 编者按】在如今这个 AI 飞速发展的时代，Python 以其简洁、易读的语法和强大的表现力，深受程序员和科学家的喜爱。本文将展示在仅仅 10 行代码内，Python 如何实现复杂的数学运算，包括向量加法、点积、矩阵乘法、矩阵转置以及线性系统求解器等。

原文链接：https://wordsandbuttons.online/how_much_math_can_you_do_in_10_lines_of_python.html

翻译 | 郑丽媛

出品 | 程序人生（ID：coder_life）

记得我刚开始用 Python 时，还没有 NumPy 或 SymPy。那时，我们常常在 MatLab 中做研究，在 Delphi 中进行快速原型设计……说出来你可能不信，我们还会在 Prolog 中做符号计算。那时的 Python 类似于一个好用版本的 Perl，但研究人员并不爱用。

不过，用它来编写程序还是挺有趣的，所以我选择用 Python 来实现我的论文中的几个实验，其中需要用到一些线性代数。起初，我写了一些 Pascal 风格的循环嵌套程序，感觉不是很有意思，但随着我对 Python 研究得越来越深入，我发现事情开始变得有趣起来。

某一刻，我感觉这门语言不仅仅是写循环嵌套，于是开始深入研究这个语言，而不仅限于我的论文内容——这可能是我做过的最好的错误决定。

Python 以其友好的学习曲线而闻名，但这也可能会带来问题。你可以用 Python 写 Pascal 程序，可以用 Python 写 C++ 程序，还可以用 Python 写 Lisp 程序，但只有当你用 Python 的方式写 Python 程序时，才能真正发挥它的优势。

好在，就算仅掌握 Python 的基本知识可能会错过精彩部分，对用户来说体验依然很好。

为了说明这一点，我选择了线性代数领域。按理说，Python 有许多优秀的库可以处理线性代数，所以不必重新实现下文中提到的任何内容。而我选择线性代数，是想要证明在 Python 中，以如此少量的代码来表达丰富含义的语法是完全可行的。（免责声明：为了保证可读性，部分代码示例并没有严格限制在一行内，示例中的所有缩进完全是为了提升阅读体验。）

列表解析

列表解析是 Python 单行代码的精髓，这是一种描述列表转换的特殊语法。比方说，如果我们要用一个向量乘以一个标量，在类似 Pascal 的 Python 代码中，可能会这样写：

def scaled(A, x): B = list() for i in range(len(A)): B.append( A[i] * x ) return B

怎么说呢，这种写法没啥问题，也有它的好处，例如总能找到一行放置断点，但这种写法有些“冗长”。在 Python 中，你可以这样简单地写：

def scaled(A, x): return [ai*x for ai in A]

然后它便会这样工作：

List comprehension
 [1.0, 4.0, 3.0]* 2.0 [, ...

Step 1
 [1.0, 4.0, 3.0] * 2.0 [2.0, ...

Step 2
 [1.0, 4.0, 3.0] * 2.0 [2.0, 8.0, ...

Step 3
 [1.0, 4.0, 3.0] * 2.0 [2.0, 8.0, 6.0]

不过，列表解析不是 Python 所独有的，Haskell 和 Clojure 也有类似功能，甚至 C# 的 LINQ 也提供了针对范围的特殊语法。而鲜少有人知道的是，Python 还支持字典解析和元组解析。

>>> [2*v for v in [1.0, 2.0, 3.0]][2.0, 4.0, 6.0]
>>> {k:2*v for k, v in {0:1.0, 1:4.0, 2:3.0}.items()}{0: 2.0, 1: 8.0, 2: 6.0}
>>> tuple(2*v for v in (1.0, 4.0, 3.0))(2.0, 8.0, 6.0)

列表压缩

列表压缩可将多个可迭代对象作为一个整体进行迭代，并将所有对象都转化为一个元组列表。虽然标题上写的是“列表”，但它也同样适用于元组、字典、生成器等任何可迭代的对象。

>>> A = [1, 2, 3]>>> B = ('a', 'b', 'c')>>> C = {1:'a', 2:'b', 3:'c'}>>> D = xrange(123)>>> zip(A, B, C, D)
[(1, 'a', 1, 0), (2, 'b', 2, 1), (3, 'c', 3, 2)]

同样地，你可以用一行代码实现向量加法。

def sum_of(A, B): return [ai+bi for (ai, bi) in zip(A, B)]

对于两组列表，这就能像拉链一样将数据压缩在一起。

List zipping
 A = [1.0, 4.0, 3.0] B = [7.0, 3.0, 1.0]zip(A, B) = [...

Step 1
 A = [1.0, 4.0, 3.0] B = [7.0, 3.0, 1.0]zip(A, B) = [(1.0, 7.0), ...

Step 2
 A = [1.0, 4.0, 3.0] B = [7.0, 3.0, 1.0]zip(A, B) = [(1.0, 7.0), (4.0, 3.0), ...

Step 3
 A = [1.0, 4.0, 3.0] B = [7.0, 3.0, 1.0]zip(A, B) = [(1.0, 7.0), (4.0, 3.0), (3.0, 1.0)]

求和函数

求和函数可以进行简单地求和操作。你可以通过累积向量元素的乘积，用一行代码实现向量点积。

def dot_of(A, B): return sum([ai*bi for (ai, bi) in zip(A, B)])

不过，有时简单的求和并不够用。比如在处理浮点数时，你可能会遇到一些恼人的小误差。为了让你更方便，Python 提供了另一个求和函数，它可以对部分求和进行操作，从而获得更精确的输出。

>>> [0.1]*10[0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1]
>>> sum([0.1]*10)0.9999999999999999
>>> import math>>> math.fsum([0.1]*10)1.0

条件表达式

条件表达式是 Python 中的三元运算符。它们是依赖于条件的表达式，通常也是简单而纯粹的表达式。我们可以用条件表达式来制作一个单位矩阵。

def identity(n): return [[1.0 if i==j else 0.0 for j in range(n)] for i in range(n)]

这实际上是一个带有条件表达式的嵌套列表解析，结果如下：

Conditional expression: 0 == 0
 j = 0[[1.0, ... i = 0

Conditional expression: 1 != 0
 j = 1[[1.0, 0.0, ... i = 0

Conditional expression: 2 != 0
 j = 2[[1.0, 0.0, 0.0], ... i = 0

Conditional expression: 0 != 1
 j = 0[[1.0, 0.0, 0.0], [0.0, ... i = 1

Conditional expression: 1 == 1
 j = 1[[1.0, 0.0, 0.0], [0.0, 1.0, ... i = 1

Conditional expression: 2 != 1
 j = 2[[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], ... i = 1

Conditional expression: 0 != 2
 j = 0[[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, ... i = 2

Conditional expression: 1 != 2
 j = 1[[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, ... i = 2

Conditional expression: 2 == 2
 j = 2[[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]] i = 2

条件表达式在简洁性和表达力方面非常出色，但由于它们是表达式，所以不容易添加副作用。但也并不是说不可能，只是需要一些技巧。在 Python 中，元组计算会从左到右依次计算每个元组元素，所以当你需要在表达式中添加更多内容时，只需使用元组即可。

def identity(n): return [[1.0 if i==j else (0.0, sys.stdout.write("i != j\n") )[0] for j in range(n)] for i in range(n)]

你也不能在表达式中使用 print，因为它们不应该有副作用，但可以使用 sys.stdout.write。

def identity(n): return [[1.0 if i==j else (0.0, sys.stdout.write("i != j\n") )[0] for j in range(n)] for i in range(n)]

同样的技巧也适用于 lambda 表达式，你可以通过元组让你的 lambda 表达式变得尽可能复杂。不过我建议如果不是必要情况，请不要这么做。

我猜肯定会有人反对说，“你的例子根本不需要条件表达式！”确实如此。在 Python 中，你可以显式地将布尔类型变量转换为浮点数变量。

def identity(n): return [[float(i==j) for j in range(n)] for i in range(n)]

这只是风格问题。我个人而言，更喜欢将事实和数字分开，不这样做也完全可以。

将容器内容作为参数传递

假设我们在 Python 中有一个矩阵，用列表来表示，其中每个嵌套列表是一行。如果我们将这些行传递给前面提到的 zip，它会生成一个元组列表，其中每个元组表示矩阵的一列。这样我们就得到了一个现成的矩阵转置。

不过，我们还需要解决两个问题。一个很简单：我们希望矩阵是一个列表，而不是元组，所以需要用一个简单的列表解析来解决这个问题。另一个问题需要使用特殊的 Python 语法，它能让我们把列表变成函数参数的元组，这个语法就是在矩阵前面加一个星号。

def transposed(A): return [list(aj) for aj in zip(*A)]

传递列表和传递作为参数的列表区别如下：

Passing a container
>>> A = [[1, 2, 3],...  [4, 5, 6],...  [7, 8, 9]]
>>> def print_A(A):...  print A
>>> print_A(A)
[[1, 2, 3], [4, 5, 6], [7, 8, 9]]

Passing a container's content
>>> A = [[1, 2, 3],... [4, 5, 6],... [7, 8, 9]]
>>> def print_a_b_c(a, b, c):... print a, b, c
>>> print_a_b_c(*A)
[1, 2, 3] [4, 5, 6] [7, 8, 9]

Using zip to transpose a matrix
>>> A = [[1, 2, 3],...  [4, 5, 6],...  [7, 8, 9]]
>>> def print_zip_a_b_c(a, b, c):...  print zip(a, b, c)
>>> print_zip_a_b_c(*A)
[(1, 4, 7), (2, 5, 8), (3, 6, 9)]

同样的语法也可以反向使用。如果你希望把所有参数接收为一个单独的元组，也可以这样做。

>>> def print_arguments_as_tuple(*A):...  print A
>>> print_A_as_tuple(1, 2, 3)
(1, 2, 3)

矩阵乘法

当你的矩阵表示为列表时，你可以通过一个包含点积的列表解析轻松地将矩阵与向量相乘。

def matrix_vector_of(A, X): return [dot_of(ai, X) for ai in A]

这种操作方式在投影几何中尤为重要。通过一次矩阵乘法，你可以同时实现旋转、缩放、平移、仿射和投影变换。

平面上的点投影

在计算几何中，一个平面通常由空间中的一个点和该点的法向量定义。通过这种表示法，你可以非常简单地将任意点投影到平面上。

def projected_on(A, Pn, Pd): return sum_of(A, scaled(Pn, (Pd — dot_of(Pn, A)) / dot_of(Pn, Pn)))

我们不需要为此使用任何特殊语法，但请稍等，我们一会儿回头再谈。

参数中的列表乘法

接下来，来求两点之间的欧氏距离。我们需要做的就是进行向量减法，计算其自身的点积，然后对结果取平方根。

def distance_between(A, B): return pow(dot_of( *(sum_of(A, scaled(B, -1.0)), )*2 ), 0.5)

这里的技巧是做一个自点积，也就是将相同的参数传递给函数而不重复。Python 有一个很好的语法来实现多次复制，那就是... 乘法。你可以使用整数乘法来复制列表、元组等。

>>> (1, 2, 3) * 2(1, 2, 3, 1, 2, 3)
>>> [1, 2, 3] * 2[1, 2, 3, 1, 2, 3]
>>> ([1, 2, 3],) * 2([1, 2, 3], [1, 2, 3])

在这个例子中，我们需要复制一个列表并将其传递给函数。

List concatenation
>>> A = [1, 2, 3]
>>> def print_A_B(A, B):... print A, B
>>> print_A_B(A*2 ,[])
[1, 2, 3, 1, 2, 3] []

List duplication
>>> A = [1, 2, 3]
>>> def print_A_B(A, B):... print A, B
>>> print_A_B((A,)*2, [])
([1, 2, 3], [1, 2, 3]) []

Passing a content of a list duplication
>>> A = [1, 2, 3]
>>> def print_A_B(A, B):... print A, B
>>> print_A_B(*(A,)*2)
[1, 2, 3] [1, 2, 3]

旋转与线性求解器

求解线性方程本身就是一个重要的工程领域。从概念上来说它很简单，但你会遇到各种各样的问题，比如计算误差、不合理的执行时间、内存耗尽等，这就变得不太简单了。

我想在这里想展示的是，你可以利用 Python 的强大功能，用一行代码实现一个线性求解器。不过请记住，我并不推荐这样做，有些库在这方面做得更好。

这是一个迭代求解器，通过从一个超平面跳到另一个超平面，缓慢但稳定地接近解。不过它也可能耗尽你的堆栈，这一切都取决于具体的系统。

def solve(A, B, Xi): return Xi if distance_between(matrix_vector_of(A, Xi), B) < 1e-5 else solve(A[1:] + [A[0]], B[1:] + [B[0]], projected_on(Xi, A[0], B[0]))

这里唯一有趣的是我们在投影队列中旋转超平面时使用的列表语法。一旦你理解了它，就会发现非常简单：你只需将列表的头部切下来，使其成为一个新的列表，然后将尾部与之前的头部连接起来。

Rotation: step 1
>>> B = [1, 2, 3]>>> print B[1:] + [B[0]]
[2, 3, 1]

Rotation: step 2
>>> B = [2, 3, 1]>>> print B[1:] + [B[0]]
[3, 1, 2]

Rotation: step 3
>>> B = [3, 1, 2]>>> print B[1:] + [B[0]]
[1, 2, 3]

矩阵求逆

再强调一次，这绝不是求逆矩阵的最佳方法。但它确实有效，并且只需一行代码。

我的想法是，通过求解线性方程组，将单位矩阵按单位向量逐个除以输入矩阵。线性方程组基本上就是矩阵除以向量，因为它的左边是乘法，而我们通过求解方程组得到第二个操作数。

def inverted(A): return transpose([solve(A, ort, [0.0]*len(A)) for ort in identity(len(A))])

结论

至此，仅用 10 行 Python 代码，我们就实现了向量加法、向量乘法、点积、矩阵乘向量、单位矩阵、矩阵转置、矩阵求逆、点投影、欧式距离和线性系统求解器。

我认为 Python 之所以强大，不仅是因为它拥有丰富的包。当然，现代 Python 提供了几乎所有的包，完全不需要自己编写很多基础代码。但我相信，它之所以有这么多包，是因为用 Python 写代码实在是太有趣了。

由 CSDN 和 Boolan 联合主办的「2024 全球软件研发技术大会（SDCon）」将于 7 月 4 -5 日在北京威斯汀酒店举行。

由世界著名软件架构大师、云原生和微服务领域技术先驱 Chris Richardson 和 MIT 计算机与 AI 实验室(CSAIL)副主任，ACM Fellow Daniel Jackson 领衔，BAT、微软、字节跳动、小米等技术专家将齐聚一堂，共同探讨软件开发的最前沿趋势与技术实践。

10行Python代码实现多种数学运算

热门阅读

推荐阅读