第七色在线视频,2021少妇久久久久久久久久,亚洲欧洲精品成人久久av18,亚洲国产精品特色大片观看完整版,孙宇晨将参加特朗普的晚宴

為了賬號(hào)安全，請(qǐng)及時(shí)綁定郵箱和手機(jī)立即綁定
已解決430363個(gè)問(wèn)題，去搜搜看，總會(huì)有你想問(wèn)的

LU分解中除以零警告-Doolittle算法工作

LU分解中除以零警告-Doolittle算法工作

滄海一幻覺(jué) 2022-11-01 14:45:31
我通過(guò)以下鏈接實(shí)現(xiàn)了矩陣的 LU 分解的標(biāo)準(zhǔn)方程/算法：（1）和（2）這完美地返回了如下方矩陣的 LU 分解。但是，我的問(wèn)題是-它也給出了Divide by Zero warning.代碼在這里：import numpy as npdef LUDecomposition (A):    L = np.zeros(np.shape(A),np.float64)    U = np.zeros(np.shape(A),np.float64)    acc = 0    L[0,0]=1    for i in np.arange(len(A)):        for k in range(i,len(A)):            for j in range(0,i):                acc += L[i,j]*U[j,k]            U[i,k] = A[i,k]-acc            for m in range(k+1,len(A)):                if m==k:                    L[m,k]=1                else:                    L[m,k] = (A[m,k]-acc)/U[k,k]            acc=0    return (L,U)A = np.array([[-4, -1, -2],              [-4, 12,  3],              [-4, -2, 18]])L, U = LUDecomposition (A)我哪里錯(cuò)了？
查看完整描述

1 回答

?
ibeautiful

TA貢獻(xiàn)1993條經(jīng)驗(yàn) 獲得超6個(gè)贊

看來(lái)您可能在第一個(gè)內(nèi)層for循環(huán)方面犯了一些縮進(jìn)錯(cuò)誤：U必須在之前進(jìn)行評(píng)估L；您也沒(méi)有正確計(jì)算求和項(xiàng)acc，也沒(méi)有正確地將對(duì)角項(xiàng)設(shè)置L為 1。在進(jìn)行一些其他語(yǔ)法修改后，您可以按如下方式重寫(xiě)您的函數(shù)：


def LUDecomposition(A):


    n = A.shape[0]

    L = np.zeros((n,n), np.float64)

    U = np.zeros((n,n), np.float64)


    for i in range(n):

        # U

        for k in range(i,n):

            s1 = 0  # summation of L(i, j)*U(j, k) 

            for j in range(i):

                s1 += L[i,j]*U[j,k]

            U[i,k] = A[i,k] - s1


        # L

        for k in range(i,n):

            if i==k:

                # diagonal terms of L 

                L[i,i] = 1

            else:

                s2 = 0 # summation of L(k, j)*U(j, i) 

                for j in range(i):

                    s2 += L[k,j]*U[j,i]

                L[k,i] = (A[k,i] - s2)/U[i,i]


    return L, U

與scipy.linalg.lu作為可靠參考A相比，這一次給出了矩陣的正確輸出：


import numpy as np

from scipy.linalg import lu


A = np.array([[-4, -1, -2],

              [-4, 12,  3],

              [-4, -2, 18]])


L, U = LUDecomposition(A)

P, L_sp, U_sp = lu(A, permute_l=False)


P

>>> [[1., 0., 0.],

     [0., 1., 0.],

     [0., 0., 1.]])


L

>>> [[ 1.          0.          0.        ]

     [ 1.          1.          0.        ]

     [ 1.         -0.07692308  1.        ]]


np.allclose(L_sp, L))

>>>  True


U

>>> [[-4.         -1.         -2.        ]

     [ 0.         13.          5.        ]

     [ 0.          0.         20.38461538]]


np.allclose(U_sp, U))

>>>  True

注意：與 scipy lapack getrf 算法不同，此 Doolittle 實(shí)現(xiàn)不包括旋轉(zhuǎn)，這兩個(gè)比較只有在P返回的置換矩陣scipy.linalg.lu是單位矩陣時(shí)才為真，即scipy 沒(méi)有執(zhí)行任何置換，這確實(shí)是您的矩陣的情況A. 在 scipy 算法中確定的置換矩陣旨在優(yōu)化結(jié)果矩陣的條件數(shù)，以減少舍入誤差。最后，您可以簡(jiǎn)單地驗(yàn)證一下A = LU，如果分解正確，情況總是如此：


A = np.random.rand(10,10)

L, U = LUDecomposition(A)


np.allclose(A, np.dot(L, U))

>>>  True

盡管如此，就數(shù)值效率和準(zhǔn)確性而言，我不建議您使用自己的函數(shù)來(lái)計(jì)算 LU 分解。希望這可以幫助。


查看完整回答
反對(duì) 回復(fù) 2022-11-01
  • 1 回答
  • 0 關(guān)注
  • 106 瀏覽
慕課專(zhuān)欄
更多

添加回答

舉報(bào)

0/150
提交
取消
微信客服

購(gòu)課補(bǔ)貼
聯(lián)系客服咨詢(xún)優(yōu)惠詳情

 幫助反饋 APP下載

慕課網(wǎng)APP
您的移動(dòng)學(xué)習(xí)伙伴

 公眾號(hào)

掃描二維碼
關(guān)注慕課網(wǎng)微信公眾號(hào)