python - 有没有更清晰的方法来用 numpy 求解线性方程组？

Question

我有一组 6 个方程，我希望 numpy 帮我求解。所以我构建了一个 6x6 的系数矩阵，并用各种值填充它。然而，我最终为此编写的代码非常难以辨认，并且几乎没有向我的代码读者传达我想要求解的方程式。

例如，填充系数矩阵看起来像这样:

# Coefficients matrix
# Order of variables: w, X, Y, Z, s, t
A = np.mat( np.zeros((6,6)) )

A[0:3,0] = cam_inv[...,2]
A[0:3,1:4] = -np.identity(3)
A[3:6,1:4] = np.identity(3)
A[3:,4] = -eigvecs[...,0]
A[3:,5] = -eigvecs[...,1]

# Constants matrix (RHS of equation)
b = np.mat( np.zeros((6,1)) )
b[0:3,0] = -cam_inv[...,0:2] * point
b[3:,] = mean.T

res = np.linalg.solve(A,b)

(其中 cam_inv、eigvecs、mean 和 point 是在别​​处计算的一些其他矩阵。)

显然上面的代码可以有更多的注释，但我觉得即使有一些注释，它仍然无法真正传达正在解决的基础方程式。是否有更好的方法将方程输入求解器，从而生成更清晰的代码？

Answer 1

问题是表示等式的 A 行没有到代码行的一对一映射。我在自己的工作(经济学)中所做的就是对A的每一行都有一个明确的英文名称(或者至少是一行代码，没有功能表示)的函数。必要时，我有一个明确的但与我最终使用的代码执行相同操作的代码版本较慢或可能更长。

例如(来自 Bretscher 的 Linear Algebra with Applications 1997，ex. 37 p. 29 一个简单但不切实际的例子)，考虑一个具有三个行业的经济体，I1、I2、I3，每个行业都将其他两个行业的输出作为输入.他们应该生产什么产品以满足消费者和工业需求？

A =np.zeros((3,3)) 
#Each unit of production by I1 requires 0.1 units of good 2 and .2 of good 3
A[:,0] = [0, 0.1, 0.2] 
#Each unit of production by I2 requires 0.2 units of good 1 and .5 of good 3
A[:,1] = [0.2, 0, 0.5] 
#Each unit of production by I3 requires 0.3 units of good 1 and .4 of good 2    
A[:,2] = [0.3, 0.4, 0] 
#The required production for consumers. 
b = np.array([320,150,90]).reshape(-1,1)
#The optimal production levels of x1, x2, and x3
res = np.linalg.solve(A,b)

按照我的建议做可能会更慢或更不简洁，但阅读起来会更清晰。