Python 线性方程组求解之:Jacobi迭代算法
程序员文章站
2022-08-03 08:27:36
原理: 注: 实例数据:jacobi_data.txt 。组织方式:[A,b] 2,-1,-1,-51,1,10,111,5,-1,8 ......
- import numpy as np
- delta=0.0001#精度要求
- #数据读取
- data = []
- f = open("h:/notepad/数学实验/jacobi_data.txt")
- for line in f:
- line = line.replace("\n","")
- data.append(list(map(eval, line.split(","))))
- f.close()
- data=np.array(data)
- #对data行初等变换,主元变为每列最大值
- row,column=data.shape
- for i in range(row):
- max_value_index=np.argmax(np.fabs(data[i:row,i]))
- temp=np.copy(data[i,:])
- data[i,:]=data[max_value_index+i,:]
- data[max_value_index+i,:]=temp
- #lu: -(l+u)
- #d:系数矩阵的对角线元素
- #b:ax=b中的b
- lu=np.negative(data[:,0:column-1])
- d=np.zeros(row)
- b=data[:,column-1]
- for i in range(row):
- d[i]=data[i,i]
- lu[i,i]=0
- #迭代求解
- x=np.ones(row)#用于存储迭代过程中x的值
- y=np.ones(row)#用于存储中间结果
- for iteration in range(100):
- print('x:',x)
- #迭代计算
- for i in range(row):
- y[i]=np.vdot(lu[i,:],x)+b[i]
- y[i]=y[i]/d[i]
- #判断是否达到精度要求
- if np.max(np.fabs(x-y))<delta:
- print('iteration:',iteration)
- break
- #将y幅值到x,开始下一轮迭代
- x=np.copy(y)
原理:
注:
实例数据:jacobi_data.txt 。组织方式:[a,b]
2,-1,-1,-5
1,1,10,11
1,5,-1,8
下一篇: 睡前练瑜伽瘦身减压 教你九式睡前瑜伽