想当爷爷的孙子 2020-05-04 12:34 采纳率: 66.7%
浏览 334
已采纳

python 用雅可比方法求解稀疏矩阵,完全不会写,求助a

图片说明

import numpy as np
def Jacobi(A, b, iter_n, initial_guess=0):
n = len(A)

D = np.diag(A)
R = A - np.diag(D)
x_i = initial_guess * np.ones(n)
for i in range(iter_n):
    print('x_',i,'=',x_i)
    x_i = (b - R.dot(x_i)) / D

return x_i

def A_ij(n):
A = np.empty((n, n))
for i in range(n):
A[i, i] = 2
for i in range(n-1):
A[i,i+1]=A[i+1,i]=1
return A
def b_i(n):
b=np.empty(n)
b[0]=1
b[n-1]=-1
return b
def x0(n):
return np.zeros(n)
print(Jacobi(A,b,100,x0))

  • 写回答

1条回答 默认 最新

  • threenewbee 2020-05-04 12:56
    关注
    # -*- coding: utf-8 -*-
    
    #Jacobi迭代法 输入系数矩阵mx、值矩阵mr、迭代次数n、误差c(以list模拟矩阵 行优先)
    
    def Jacobi(mx,mr,n=100,c=0.0001):
        if len(mx) == len(mr):  #若mx和mr长度相等则开始迭代 否则方程无解
            x = [] #迭代初值 初始化为单行全0矩阵
            for i in range(len(mr)):
                x.append([0])
            count = 0 #迭代次数计数
            while count < n:
                nx = [] #保存单次迭代后的值的集合
                for i in range(len(x)):
                    nxi = mr[i][0]
                    for j in range(len(mx[i])):
                        if j!=i:
                            nxi = nxi+(-mx[i][j])*x[j][0]
                    nxi = nxi/mx[i][i]
                    nx.append([nxi]) #迭代计算得到的下一个xi值
                lc = [] #存储两次迭代结果之间的误差的集合
                for i in range(len(x)):
                    lc.append(abs(x[i][0]-nx[i][0]))
                if max(lc) < c:
                    return nx #当误差满足要求时 返回计算结果
                x = nx
                count = count + 1
            return False #若达到设定的迭代结果仍不满足精度要求 则方程无解
        else:
            return False
    
    #调用 Jacobi(mx,mr,n=100,c=0.001) 示例
    mx = [[8,-3,2],[4,11,-1],[6,3,12]]
    
    mr = [[20],[33],[36]]
    print(Jacobi(mx,mr,100,0.00001))
    
    

    https://blog.csdn.net/cswfqxs_/article/details/84067711

    本回答被题主选为最佳回答 , 对您是否有帮助呢?
    评论

报告相同问题?

问题事件

  • 请详细说明问题背景 4月14日

悬赏问题

  • ¥15 数据库获取信息反馈出错,直接查询了ref字段并且还使用了User文档的_id而不是自己的
  • ¥15 将安全信息用到以下对象时发生以下错误:c:dumpstack.log.tmp 另一个程序正在使用此文件,因此无法访问
  • ¥15 速度位置规划实现精确定位的问题
  • ¥15 代码问题:df = pd.read_excel('c:\User\18343\Desktop\wpsdata.xlxs')路径读不到
  • ¥15 为什么视频算法现在全是动作识别?
  • ¥15 编写一段matlab代码
  • ¥15 用Python做岩石类别鉴定软件
  • ¥15 关于调取、提交更新数据库记录的问题
  • ¥15 之前删了盘从下vs2022遇见这个问题 搞了一整天了
  • ¥15 从Freecad中宏下载的DesignSPHysics,出现如下问题是什么原因导致的(语言-python)