如何用Python求解含尽可能多零解的线性方程组？

1. 让NumPy生成含尽可能多零的解 #

NumPy自带的最小二乘是基于L2范数最小化，本来就是找L2范数最小的解，这种解的非零元素通常分布比较均匀，没法直接设置让它优先生成零元素。你要的是最稀疏的解，得换用专门的稀疏优化方法：

• 用L1正则化诱导稀疏性：可以用scipy.optimize.lasso_path或者sklearn.linear_model.Lasso，通过给损失函数加L1惩罚项，迫使尽可能多的系数变为零，同时满足方程组的拟合要求，适合你这种方程数远少于变量数的欠定场景。
• 基追踪（Basis Pursuit）：如果要严格满足方程组（而非近似拟合），可以用cvxpy构造L1范数最小化的优化问题，直接找解空间里最稀疏的解，示例代码如下：

  import cvxpy as cp
  import scipy.sparse as sp
  
  # 假设A是10x10000的稀疏矩阵，b是10维向量
  A = sp.csr_matrix(A)  # 用稀疏格式存储，大幅节省内存
  x = cp.Variable(10000)
  objective = cp.Minimize(cp.norm(x, 1))
  constraints = [A @ x == b]
  prob = cp.Problem(objective, constraints)
  prob.solve(solver='ECOS')  # 选择适合稀疏问题的求解器
  sparse_solution = x.value
  

2. 用SymPy直接求解大型方程组 #

完全不现实。SymPy是符号计算工具，擅长处理小规模的符号推导问题，10000个变量的规模远超它的能力范围：

• 符号计算的复杂度随变量数指数增长，这么多变量会直接耗尽内存，计算时间更是无法忍受。
• SymPy的设计目标不是大规模数值稀疏问题，所以这条路走不通。

内容的提问来源于stack exchange，提问作者Gauss