4. Python调用MATLAB引擎求解线性方程组

任务拆解：Python 端发起、MATLAB 端求解

本例要完成一条最小但完整的跨语言计算链：Python 负责组织系数矩阵 A 与右端向量 b，通过 MATLAB Engine API for Python 启动 MATLAB 引擎，在 MATLAB 工作区中执行求解，再把解向量 x 取回 Python 并做残差校验。关键点不在“解方程”本身，而在于数据如何跨语言边界传递、类型如何转换、结果如何验证。

环境准备与引擎启动

MATLAB Engine API for Python 需要与已安装的 MATLAB 版本匹配。安装方式是在 MATLAB 安装目录下找到 extern/engines/python，执行 python setup.py install（或使用 pip install matlabengine 指定版本）。安装后在 Python 中 import matlab.engine 不报错即表示可用。

import numpy as np
import matlab.engine

# 启动引擎（首次启动约数秒）
eng = matlab.engine.start_matlab()
print("MATLAB 版本:", eng.version())

注意：start_matlab() 会启动一个独立的 MATLAB 进程，Python 与它通过共享内存/管道通信。引擎对象 eng 是后续所有调用的入口，任务结束后应调用 eng.exit() 释放，否则进程会残留。

构造方程组并传入工作区

选一个 3 元线性方程组作为教学示例，系数矩阵取对角占优的可逆矩阵，保证解唯一且数值稳定：

A = np.array([[4.0, -1.0,  0.0],
              [-1.0, 4.0, -1.0],
              [ 0.0, -1.0, 3.0]])
b = np.array([15.0, 10.0, 10.0])

# 传入 MATLAB 工作区：numpy → matlab.double
eng.workspace['A'] = matlab.double(A.tolist())
eng.workspace['b'] = matlab.double(b.reshape(-1, 1).tolist())

这里有两个容易出错的细节：其一，matlab.double 接受嵌套列表，numpy 数组需先 tolist()；其二，MATLAB 中向量默认按列理解，把 b 转成 (3,1) 的列向量可避免后续维度歧义。若直接传一维 b，A \ b 仍能求解，但返回形状会随版本略有差异，显式列向量更稳妥。

在 MATLAB 端求解并取回结果

# 在 MATLAB 工作区执行左除求解
eng.eval("x = A \\ b;", nargout=0)

# 取回结果：matlab.double → numpy
x_mat = eng.workspace['x']
x = np.array(x_mat).reshape(-1)
print("解向量 x =", x)

eng.eval 的 nargout=0 表示不向 Python 返回表达式结果，只让变量留在 MATLAB 工作区；随后通过 eng.workspace['x'] 读取。也可以直接写 x_mat = eng.eval("A \\ b", nargout=1) 一步取回，两种方式等价，前者更适合需要复用中间变量的场景。

结果校验：残差与回代

跨语言调用最容易出的问题是“看起来跑通了但数据错了”，因此必须做数值校验。用 Python 侧独立计算残差 ‖Ax−b‖₂，与 MATLAB 端结果交叉验证：

residual = np.linalg.norm(A @ x - b)
print("残差 =", residual)
assert residual < 1e-10, "求解结果未通过残差校验"

# 可选：在 MATLAB 端也计算一次残差，做双向确认
eng.workspace['x'] = matlab.double(x.reshape(-1, 1).tolist())
r_mat = eng.eval("norm(A*x - b)", nargout=1)
print("MATLAB 端残差 =", float(r_mat))

对本例给定数据，理论解为 x = [5, 5, 5]，残差应在 1e-12 量级（浮点舍入）。若残差明显偏大，优先检查 b 的维度是否被误当成行向量、A 是否在传入时被转置，以及是否误用了 * 而非 \。

完整可运行脚本

import numpy as np
import matlab.engine

def solve_with_matlab(A, b):
    eng = matlab.engine.start_matlab()
    try:
        eng.workspace['A'] = matlab.double(A.tolist())
        eng.workspace['b'] = matlab.double(b.reshape(-1, 1).tolist())
        eng.eval("x = A \\ b;", nargout=0)
        x = np.array(eng.workspace['x']).reshape(-1)
        return x
    finally:
        eng.exit()

if __name__ == "__main__":
    A = np.array([[4.0, -1.0, 0.0],
                  [-1.0, 4.0, -1.0],
                  [0.0, -1.0, 3.0]])
    b = np.array([15.0, 10.0, 10.0])
    x = solve_with_matlab(A, b)
    print("x =", x)
    print("残差 =", np.linalg.norm(A @ x - b))

脚本用 try/finally 保证引擎一定被释放，避免异常路径下 MATLAB 进程残留。运行后预期输出解向量接近 [5. 5. 5.]，残差在 1e-12 量级。至此，一条“Python 组织数据 → MATLAB 求解 → Python 校验”的跨语言计算链闭环完成。

常见问题速查