高斯 賽德爾迭代法python實現

2021-10-07 01:26:21 字數 1808 閱讀 6525

import scipy

import scipy.linalg

import numpy.matlib

import numpy as np

import time

import warnings

warnings.filterwarnings('ignore')

def gaussseidel_inv(a,b):

l=-np.tril(a,-1)

u=-np.triu(a,1)

d=np.dot(-2,np.identity(n))

invdl=np.linalg.inv(d-l)

b=np.matmul(invdl,u)

g=np.matmul(invdl,b)

x0=g

x=np.matmul(b,x0)+g

while (np.linalg.norm(x-x0)>1e-14):

x0=x

x=np.matmul(b,x0)+g

return x

n=100

a=np.dot(-2,np.identity(n))

for k in range(0,n-1):

a[k,k+1]=1

a[k+1,k]=1

b=np.ones([n,1])

x=np.linalg.solve(a,b)

begin_time=time.time()

x_gaussseidel=gaussseidel_inv(a,b)

end_time=time.time()

print("runtime",end_time-begin_time)

error_gaussseidel=np.linalg.norm(x-x_gaussseidel)

print("error",error_gaussseidel)

遇到的問題

1.runtimewarning: overflow encountered in matmul x=np.matmul(b,x0)+g

我拿到這個錯誤內心是很迷的,這啥錯誤???後來在網上一搜,雖然沒理解為啥出錯,但是你可以簡單粗暴啊

import warnings

warnings.filterwarnings('ignore')

果然加上之後世界都變得美好。

2.我自己感覺我這個**寫的不咋樣是真的,因為是拿老師給的matlab**改的python,雖然原理理解了,但有些東西都是硬湊出來的,而且三種方法我只改了一種,最後附上matlab的**。不說了,滾去複習數字訊號處理,上網課就是這麼難。

function x = gaussseidel_inv( a,b )

%gaussseidel_inv 此處顯示有關此函式的摘要

n=size(a,1);

d=diag(diag(a,0)); %求矩陣a的第1條對角線的元素。

l=-tril(a,-1); %求矩陣a的第-1條對角線以下的元素。

u=-triu(a,1) ; %求矩陣a的第1條對角線以上的元素。

invdl=(d-l)\eye(n);

b=invdl*u;%(d-l)的逆乘以u

g=invdl*b;%(d-l)的逆乘以b

x0=g;

x=b*x0+g;

while(norm(x-x0)>1e-14)

x0=x;

x=b*x0+g;

endend

雅可比迭代法 高斯 賽德爾迭代法

求解方程組 用雅可比迭代法求解方程組ax b 輸入 a為方程組的係數矩陣,b為方程組右端的列向量,輸入 x0為迭代初值構成的列向量,nm為最大迭代次數,eps為誤差精度 輸出 x為求得的方程組的解構成的列向量,k為迭代次數 d diag diag a l tril a,1 u triu a,1 b ...

雅克比迭代法與高斯塞德爾迭代法求解方程組(C語言)

分別用雅可比 迭代法與高斯塞德爾迭代法解下列方程組 雅可比迭代法 include include define eps 1e 6 define max 100 雅可比迭代法 void jacobi float a,int n,float x y i a i n 1 n s a i n 1 i eps...

高斯 賽德爾迭代演算法 C 實現

輸入 係數矩陣a,最大迭代次數n,初始向量,誤差限e 輸出 解向量 定義係數矩陣 double a 50 50 定義x1解的陣列 double rootx1 100 定義x2解的陣列 double rootx2 100 定義x3解的陣列 double rootx3 100 定義x1的迭代公式 dou...