如何在Python项目中实现一个Jacobi迭代算法-创新互联

这篇文章给大家介绍如何在Python项目中实现一个Jacobi迭代算法,内容非常详细,感兴趣的小伙伴们可以参考借鉴,希望对大家能有所帮助。

成都创新互联公司是专业的东台网站建设公司,东台接单;提供做网站、成都网站建设,网页设计,网站设计,建网站,PHP网站建设等专业做网站服务;采用PHP框架,可快速的进行东台网站开发网页制作和功能扩展;专业做搜索引擎喜爱的网站,专业的做网站团队,希望更多企业前来合作!
import numpy as np
import time

1.1 Jacobi迭代算法

def Jacobi_tensor_V2(A,b,Delta,m,n,M):
start=time.perf_counter()#开始计时
find=0#用于标记是否在规定步数内收敛
X=np.ones(n)#迭代起始点
x=np.ones(n)#用于存储迭代的中间结果
d=np.ones(n)#用于存储Ax**(m-2)的对角线部分
m1=m-1
m2=2-m
for i in range(M):
print('X',X)
a=np.copy(A)
#得Ax**(m-2)
for j in range(m-2):
a=np.dot(a,X)
#得d 和 (2-m)Dx**(m-2)+(L'+U')x**(m-2)
for j in range(n):
d[j]=a[j,j]
a[j,j]=m2*a[j,j]
#迭代更新
for j in range(n):
x[j]=(b[j]-np.dot(a[j],X))/(m1*d[j])
#判断是否满足精度要求
if np.max(np.fabs(X-x))

1.2 张量A的生成函数和向量b的生成函数:

def Creat_A(m,n):#生成张量A
size=np.full(m, n)
X=np.ones(n)
while 1:
#随机生成给定形状的张量A
A=np.random.randint(-49,50,size=size)
#判断Dx**(m-2)是否非奇异,如果是,则满足要求,跳出循环
D=np.copy(A)
for i1 in range(n):
for i2 in range(n):
if i1!=i2:
D[i1,i2]=0
for i in range(m-2):
D=np.dot(D,X)
det=np.linalg.det(D)
if det!=0:
break
#将A的对角面张量扩大十倍,使对角面占优
for i1 in range(n):
for i2 in range(n):
if i1==i2:
A[i1,i2]=A[i1,i2]*10
print('A:')
print(A)
return A
#由A和给定的X根据Ax**(m-1)=b生成向量b
def Creat_b(A,X,m):
a=np.copy(A)
for i in range(m-1):
a=np.dot(a,X)
print('b:')
print(a)
return a

1.3 对称张量S的生成函数:

def Creat_S(m,n):#生成对称张量B
size=np.full(m, n)
S=np.zeros(size)
print('S',S)
for i in range(4):
#生成n为向量a
a=np.random.random(n)*np.random.randint(-5,6)
b=np.copy(a)
#对a进行m-1次外积,得到秩1对称张量b
for j in range(m-1):
b=outer(b,a)
#将不同的b叠加得到低秩对称张量S
S=S+b
print('S:')
print(S)
return S
def outer(a,b):
c=[]
for i in b:
c.append(i*a)
return np.array(c)
return a

1.4 实验一

def test_1():
Delta=0.01#精度
m=3#A的阶数
n=3#A的维数
M=200#大迭代步数
X_real=np.array( [2,3,4])
A=Creat_A(m,n) 
b=Creat_b(A,X_real,m)
Jacobi_tensor_V2(A,b,Delta,m,n)

关于如何在Python项目中实现一个Jacobi迭代算法就分享到这里了,希望以上内容可以对大家有一定的帮助,可以学到更多知识。如果觉得文章不错,可以把它分享出去让更多的人看到。


分享题目:如何在Python项目中实现一个Jacobi迭代算法-创新互联
浏览路径:http://scyanting.com/article/hhjii.html