我必須在 python 中將許多（大約 700 個）矩陣與一個隨機元素相乘（在下面，我使用的是盒分布）：
#define parameters
μ=2.
σ=2.
L=700
#define random matrix
T=[None]*L
product=np.array([[1,0],[0,1]])
for i in range(L):
m=np.random.uniform(μ-σ*3**(1/2), μ σ*3**(1/2)) #box distribution
T[i]=np.array([[-1,-m/2],[1,0]])
product=product.dot(T[i]) #multiplying matrices
Det=abs(np.linalg.det(product))
print(Det)
對于 μ 和 σ 的這種選擇，我得到了 e^30 量級的量，但是這個量應該收斂到 0。我怎么知道？因為從分析上可以證明它等價于：
Y=[None]*L
product1=np.array([[1,0],[0,1]])
for i in range(L):
m=np.random.uniform(μ-σ*(3**(1/2)), μ σ*(3**(1/2))) #box distribution
Y[i]=np.array([[-m/2,0],[1,0]])
product1=product1.dot(Y[i])
l,v=np.linalg.eig(product1)
print(abs(l[1]))
這確實給出了 e^-60。我認為這里存在溢位問題。我該如何解決？
編輯：
預計這兩個列印數量是等價的，因為第一個是以下行列式的絕對值：
也就是說，根據比內定理（乘積的行列式是行列式的乘積）：
第二個代碼列印最大特征值的絕對值：
很容易看出，一個特征值為 0，另一個等于。
uj5u.com熱心網友回復：
這通常是一個難題。有幾篇關于浮點演算法和精度的好文章。這是一個著名的
一般技巧之一 - 使用比例變數。像這樣：
import numpy as np
#define parameters
μ=2.
σ=2.
L=700
#define random matrix
T=[None]*L
scale = 1
product=np.array([[1,0],[0,1]])
for i in range(L):
m=np.random.uniform(μ-σ*3**(1/2), μ σ*3**(1/2)) #box distribution
T[i]=np.array([[-1,-m/2],[1,0]])
product=np.matmul(product, T[i]) #multiplying matrices
scale *= abs(product[0][0])
product /= abs(product[0][0])
Det=abs(np.linalg.det(product*scale))
print(Det)
它使事情變得更好，但不幸的是沒有幫助。
在這種特殊情況下，您可以做的是將行列式相乘而不是矩陣。像這樣：
import numpy as np
#define parameters
μ=2.
σ=2.
L=700
#define random matrix
T=[None]*L
scale = 1
product=1
for i in range(L):
m=np.random.uniform(μ-σ*3**(1/2), μ σ*3**(1/2)) #box distribution
T[i]=np.array([[-1,-m/2],[1,0]])
product *= np.linalg.det(T[i]) #multiplying matrices determinants
Det=abs(product)
print(Det)
輸出：
1.1081329233309736e-61
所以這可以解決問題。
轉載請註明出處，本文鏈接：https://www.uj5u.com/yidong/519649.html
下一篇：如何確保熊貓日期范圍的間距均勻？