# -*- coding: utf-8 -*-
"""
Created on Mon Jan 19 16:29:21 2015

Listing 10.2 : factorisation QR
@author: grivet

"""
from pylab import *

def affiche_mat(titre,M):
    print (titre)
    for row in M:
        for val in row:
            print (' %10.5f' %val, end = "")
        print()

N = 5
M = rand(N,N) - 0.5*ones((N,N))
V = zeros((N,N)); Q = eye(N,N)
A = diag(arange(N))
A = dot(M.transpose(),dot(A,M))
A0 = copy(A)
affiche_mat('matrice de depart: ',A)
for k in range(N-1):
    x = A[k:N,k]
    e1 = zeros_like(x)
    e1[0] = 1
    v = x + sign(x[0])*norm(x)*e1
    v /= norm(v)
    V[k:N,k] = v
    A[k:N,k:N] = A[k:N,k:N] -2*outer(v,dot(v.transpose(),A[k:N,k:N]))
affiche_mat('matrice R: ',A)
affiche_mat('matrice des vecteurs: ',V)
for k in range(N):
    Q = Q - 2*outer(dot(Q,V[:,k]),V[:,k].transpose())
affiche_mat('matrice Q: ',Q)
affiche_mat('produit QR: ',dot(Q,A))