SharedWork6.ipynbOpen in CoCalc
Author: Vlad Dfgdfg
Views : 15
Description: Jupyter notebook Work6.ipynb
In [59]:
# -*- coding: utf-8 -*-
from __future__ import division, absolute_import, print_function
import numpy as np

def form_a(n, m):
a = np.zeros((n, n))
for i in range(1, n+1):
for j in range(1, n+1):
if i == j:
a[i-1][j-1] = n + m**2 + j/m + i/n
else:
a[i-1][j-1] = (i + j)/(m + n)
return a

def form_b(n):
b = np.zeros(n)
for i in range(1,n+1):
b[i-1] = 200 + 50*i
return b

def p_print(x):
return str(x.tolist()).replace(']','').replace('[', '').replace(',','')

def cholesky(A):
n = len(A)

L = np.zeros((n, n))

for i in range(n):
for j in range(i + 1):
t_sum = 0
for k in range(j):
t_sum += L[i][k] * L[j][k]

if i == j:
L[i][j] = np.sqrt(A[i][i] - t_sum)
else:
L[i][j] = 1.0 / L[j][j] * (A[i][j] - t_sum)

return L

def solve_ls(L, b):
n = len(L)

y = np.matrix([0.0 for i in range(n)]).T
for i in range(n):
y[i, 0] = (b[i] - L[i, :i] * y[:i, 0]) / L[i, i]

L = L.T
x = np.matrix([0.0 for i in range(n)]).T
for i in range(n - 1, -1, -1):
x[i, 0] = (y[i] - L[i, i:] * x[i:, 0]) / L[i, i]
return p_print(x)

if __name__ == "__main__":
n = 20
m = 8

print(solve_ls(cholesky(form_a(n, m)), form_b(n)))

2.132160949534369 2.663065103035645 3.1926648789835226 3.720965078564629 4.247970479431361 4.773685835845899 5.298115878823181 5.821265316272813 6.343138833139963 6.863741091545193 7.3830767309233005 7.901150368161144 8.417966597734434 8.93352999184357 9.447845100548482 9.960916451902452 10.47274855208504 10.983345885533966 11.492712915076117 12.000854082057543