import numpy as np
n=3 #dimenze, tj. počet rovnic a neznámých
A=np.array([[1,2,2],[1,1,1],[1,-1,1]]) #matice soustavy
b=np.array([1,1,5]) #vektor pravých stran
Jacobiho iterace:
x=np.array([2,0,0])
for k in range(4):
y=x
x=np.array([(b[i]-sum([A[i][j]*y[j] for j in range(i)])-sum([A[i][j]*y[j] for j in range(i+1,n)]))/A[i][i] for i in range(n)])
print(x)
# if np.linalg.norm(y-x,np.infty)<0.01: #ukončovací podmínka: rozdíl po sobě jdoucích iterací
# break
# if np.linalg.norm(np.matmul(A,x)-b,np.infty)<0.01: #ukončovací podmínka: reziduum
# break
print(k+1)
[ 1. -1. 3.] [-3. -3. 3.] [1. 1. 5.] [-11. -5. 5.] 4
Gauss–Seidel:
x=[2,0,0]
for k in range(4):
for i in range(n):
x[i]=(b[i]-sum([A[i][j]*x[j] for j in range(i)])-sum([A[i][j]*x[j] for j in range(i+1,n)]))/A[i][i]
print(x)
# if np.linalg.norm(y-x,np.infty)<0.01: #ukončovací podmínka: rozdíl po sobě jdoucích iterací
# break
# if np.linalg.norm(np.matmul(A,x)-b,np.infty)<0.01: #ukončovací podmínka: reziduum
# break
print(k+1)
[1.0, 0.0, 4.0] [-7.0, 4.0, 16.0] [-39.0, 24.0, 68.0] [-183.0, 116.0, 304.0] 4
Zkusíme zanalyzovat iterační matice
B=np.array([[-A[i][j]/A[i][i] for j in range(n)] for i in range(n)])
for i in range(n):
B[i][i]=0
#omega=1
#B=omega*B+(1-omega)*np.identity(n)
print("B=\n",B)
c=np.array(b/np.diag(A))
print("c=\n",c)
B= [[ 0. -2. -2.] [-1. 0. -1.] [-1. 1. 0.]] c= [1. 1. 5.]
x=[2,0,0]
x=np.matmul(B,x)+c
print(x)
[ 1. -1. 3.]
print("Spektrální poloměr:",max(np.abs(np.linalg.eigvals(B))))
print("Řádková norma: ",np.linalg.norm(B,np.infty))
print("Sloupcová norma: ",np.linalg.norm(B,1))
Spektrální poloměr: 1.732050807568877 Řádková norma: 4.0 Sloupcová norma: 3.0
B=np.linalg.inv(np.tril(A))
c=np.matmul(B,b)
B=-np.matmul(B,np.triu(A,1))
#omega=1
#B=omega*B+(1-omega)*np.identity(n)
print("B=\n",B)
print("c=\n",c)
B= [[-0. -2. -2.] [-0. 2. 1.] [-0. 4. 3.]] c= [1. 0. 4.]
x=[2,0,0]
x=np.matmul(B,x)+c
print(x)
[1. 0. 4.]
print("Spektrální poloměr:",max(np.abs(np.linalg.eigvals(B))))
print("Řádková norma: ",np.linalg.norm(B,np.infty))
print("Sloupcová norma: ",np.linalg.norm(B,1))
Spektrální poloměr: 4.561552812808831 Řádková norma: 7.0 Sloupcová norma: 8.0