Как реализовать QR-разложение на основе Householder в Python? ⇐ Python

[Anonymous]

В настоящее время я пытаюсь реализовать QR-разложение на основе Хаусхолдера для прямоугольных матриц, как описано в http://eprints.ma.man.ac.uk/1192/1/qrup ... 2nov08.pdf (страницы 3, 4, 5).

Очевидно, я неправильно понял часть псевдокода, поскольку (1) мои результаты отличаются от numpy.qr.linalg() и (2) матрица R, созданная моими подпрограммами, не является верхней треугольной матрицей.

Мой код (также доступен по адресу https://pyfiddle.io/fiddle/afcc2e0e-085 ... d3/?i=true)

Код: Выделить всё

import math
import argparse
import numpy as np
from typing import Union

def householder(alpha: float, x: np.ndarray) -> Union[np.ndarray, int]:
"""
Computes Householder vector for alpha and x.
:param alpha:
:param x:
:return:
"""

s = math.pow(np.linalg.norm(x, ord=2), 2)
v = x

if s == 0:
tau = 0
else:
t = math.sqrt(alpha * alpha + s)
v_one = alpha - t if alpha  Union[np.ndarray, np.ndarray]:
"""
Applies Householder-based QR decomposition on specified matrix A.
:param A:
:param m:
:param n:
:return:
"""
H = []
R = A
Q = A
I = np.eye(m, m)

for j in range(0, n):
# Apply Householder transformation.
x = A[j + 1:m, j]
v_householder, tau = householder(np.linalg.norm(x), x)
v = np.zeros((1, m))
v[0, j] = 1
v[0, j + 1:m] = v_householder

res = I - tau * v * np.transpose(v)
R = np.matmul(res, R)
H.append(res)

return Q, R

m = 10
n = 8

A = np.random.rand(m, n)
q, r = np.linalg.qr(A)
Q, R = qr_decomposition(A, m, n)

print("*****")
print(Q)
print(q)
print("-----")
print(R)
print(r)

Поэтому я не понимаю, как вводить нули в мою матрицу R/какая часть моего кода неверна. Буду рад любым подсказкам!
Большое спасибо за уделенное время.

Подробнее здесь: https://stackoverflow.com/questions/534 ... -in-python