Numpy矩阵列无法更新——Householder变换实现中的赋值问题

问题分析 #

核心问题是构造的v向量形状与A[:,j]不匹配，导致赋值操作无法正确覆盖原矩阵的列。

你的代码中，v = a - alpha*np.eye(1, A.shape[0], k)生成的是**(1, n)的二维行向量**，而A[:,j]是(n,)的一维数组。尽管Numpy的广播机制允许两者完成计算，但赋值时二维数组无法直接覆盖一维的矩阵列，最终导致原矩阵没有被正确更新。此外，beta = np.inner(v, v)在v为二维时会得到(1,1)的二维数组，后续计算的2*gamma/beta也是二维数组，进一步加剧了形状不匹配问题。

解决方案 #

将标准基向量构造为一维数组，确保v与A[:,j]形状完全一致，同时处理整数矩阵的类型截断问题：

修改后的完整代码 #

import numpy as np
from numpy import linalg

def householder_qr(A):
    """Householder QR Factorization
    
    """
    # 转为浮点型，避免整数截断
    A = A.astype(np.float64)
    for k in range(0, min(A.shape[1], A.shape[0]-1)):
        alpha = -np.sign(A[k, k])*linalg.norm(A[k:, k], 2)
        # 构造v向量
        a = A[:, k].copy()
        a[:k] = 0
        # 用一维数组构造标准基向量e_k
        e_k = np.zeros(A.shape[0])
        e_k[k] = 1
        v = a - alpha * e_k
        beta = np.inner(v, v)
        if beta == 0:
            continue
        for j in range(k, A.shape[1]):
            gamma = np.inner(v, A[:, j])
            # 形状匹配，直接赋值更新
            A[:, j] = A[:, j] - (2*gamma/beta)*v
        print('A: ', A)
        print('--------------------------')

关键修改点 #

1. 标准基向量一维化：用np.zeros构造一维的e_k，确保v是(n,)的一维数组，与A[:,j]形状完全匹配，解决赋值失效问题。
2. 矩阵类型转换：将输入矩阵转为float64类型，避免整数矩阵赋值浮点数结果时的自动截断。
3. 简化冗余操作：去掉不必要的flatten()调用，直接进行计算与赋值。

验证测试 #

使用你提供的输入矩阵测试：

A = np.array([[1, 0, 0],
              [0, 1, 0],
              [0, 0, 1],
              [-1, 1, 0],
              [-1, 0, 1],
              [0, -1, 1]])
householder_qr(A)

修改后，矩阵的列会被正确更新为Householder变换后的结果，赋值操作生效。

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