# -*- coding: utf-8 -*-
"""
Created on Tue Dec 23 16:17:33 2014

@author: grivet
Listing 6.3 : Factorisation de Cholesky d'une matrice symétrique

"""
import numpy as np

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

N = 5;
#
# formation d'une matrice symétrique à diagonale dominante
#
A = np.random.rand(N,N)-0.5
A = A + A.transpose()+2.0*np.eye(N,N)
affiche_mat('A = ',A)
print
A0 = np.copy(A)
G = np.zeros((N,N)) 
v = np.zeros(N)
for j in range(N):
    v[j:] = A[j:,j]
#   print (v)
    for k in range(j):
        v[j:] = v[j:] - G[j,k]*G[j:,k]
    G[j:,j] = v[j:]/np.sqrt(v[j])
affiche_mat('matrice G = ',G)
affiche_mat('produit G*G^T', np.dot(G,G.transpose()))