提问人:Anushka Tilekar 提问时间:11/5/2023 最后编辑:Anushka Tilekar 更新时间:11/6/2023 访问量:58
在 Python 的输出中“nan”而不是矩阵 elemets?
"nan" instead of matrix elemets in output in Python?
问:
我正在尝试计算这个稀疏矩阵问题 Au=f 的近似解,但我在矩阵的输出中得到“nan”而不是矩阵元素,尽管输入矩阵都包含 np.float64。 我不明白我哪里出了问题。
这是我的代码:-
import numpy as np
import scipy
from scipy.sparse import coo_matrix
from scipy.sparse import random
#Defining a function that takes N as an input
def our_func (N) :
#Taking the given values of parameters h and k
h = 1/N
k = (29*(np.pi))/2
#Defining initial input matrix with domensions N+1
I = np.empty([N+1,N+1], dtype=np.float64)
#[REF :- https://www.geeksforgeeks.org/how-to-create-an-empty-matrix-with-numpy-in-python/]
print("datatype of I is :-")
print(I.dtype)
#Generating the matrix A in coo sparse format
A_coo = coo_matrix(I).reshape((N+1,N+1))
print("datatype of A_coo is :-")
print(A_coo.dtype)
#Converting the matrix A from coo sparse format to CSR sparse format
A = A_coo.tocsr().reshape((N+1,N+1))
print("datatype of A is :-")
print(A.dtype)
# Creating an empty storage array for the right side of our equation Au=f
f = np.empty((N+1), dtype=np.float64)
#The rows of matrix A and f are defined as :-
A[0,0] = 1 #for i==j==0
A[N,N] = 1 #for i==j==N
f[N] = 1 #for i==N
f[0] = 0 #for i==0
for i in range (int(N)) :
for j in range (int(N)) :
A[i,1] = 2 - ((h**2)*(k**2))
A[i,i+1] = -1
A[i,i-1] = -1
A[i,j] = 0 #for i!==j
#Our function Returns matrix A and f, with matrix A in sparse storage format
return A, f
# Defining a new function using the function scipy.sparse.linalg.spsolve
# to solve the sparse matrix-vector problem: Au=f
def our_solution(N):
A,f = our_func(N)
our_sparse_matrix = scipy.sparse.linalg.spsolve(A,f)
return our_sparse_matrix
# Using this defined function to compute the approximate solution for our
# problem, for N=10, N=100, and N=1000.
approx_sol_1 = our_solution(10)
print("approx_sol_1 is :- ")
print(approx_sol_1)
答:
0赞
hpaulj
11/5/2023
#1
为较小的函数运行函数:N
In [4]: A,f = our_func(10)
datatype of I is :-
float64
datatype of A_coo is :-
float64
datatype of A is :-
float64
C:\Users\14256\miniconda3\lib\site-packages\scipy\sparse\_index.py:103: SparseEfficiencyWarning: Changing the sparsity structure of a csr_matrix is expensive. lil_matrix is more efficient.
self._set_intXint(row, col, x.flat[0])
它警告我们,单独设置稀疏数组的元素是没有效率的。csr
In [5]: A
Out[5]:
<11x11 sparse matrix of type '<class 'numpy.float64'>'
with 103 stored elements in Compressed Sparse Row format>
注意的是,矩阵不是特别稀疏,103 个值中有 121 个值是非零的。
In [6]: A.A
Out[6]:
array([[ 0. , -1. , 0. , 0. ,
0. , 0. , 0. , 0. ,
0. , 0. , -1. ],
[ -1. , -18.75084325, -1. , 0. ,
0. , 0. , 0. , 0. ,
0. , 0. , 0. ],...
它可以被修剪,摆脱许多 0。
In [10]: A.eliminate_zeros()
In [11]: A
Out[11]:
<11x11 sparse matrix of type '<class 'numpy.float64'>'
with 28 stored elements in Compressed Sparse Row format>
In [17]: def my_func (N) :
...:
...: #Taking the given values of parameters h and k
...: h = 1/N
...: k = (29*(np.pi))/2
...:
...: #Defining initial input matrix with domensions N+1
...: A = np.zeros([N+1,N+1], dtype=np.float64)
...:
...: # Creating an empty storage array for the right side of our equation Au=f
...: f = np.zeros((N+1), dtype=np.float64)
...:
...: #The rows of matrix A and f are defined as :-
...: A[0,0] = 1 #for i==j==0
...: A[N,N] = 1 #for i==j==N
...: f[N] = 1 #for i==N
...: #f[0] = 0 #for i==0
...: for i in range (int(N)) :
...: for j in range (int(N)) :
...: A[i,1] = 2 - ((h**2)*(k**2))
...: A[i,i+1] = -1
...: A[i,i-1] = -1
...: A[i,j] = 0 #for i!==j
...:
...: #Our function Returns matrix A and f, with matrix A in sparse storage format
...: return A, f
...:
In [18]: A,f = my_func(10)
In [19]: A.shape
Out[19]: (11, 11)
这可以在以下情况下转换为稀疏:
In [21]: M = coo_matrix(A).tocsr()
In [22]: M
Out[22]:
<11x11 sparse matrix of type '<class 'numpy.float64'>'
with 30 stored elements in Compressed Sparse Row format>
除非太大了,否则你会遇到内存错误,我认为这种方法会更快。N
这并不能解决您的价值观问题。我对你的问题了解不够,不知道是否有合理的值。nan
A
编辑
我想知道是不是这样
A[i,j] = 0 #for i!==j
太宽泛了。如果我在 中省略它,推断它已经是 0,则有一个不同的值 (-1)。my_func
A[8,9]
In [51]: A,f = my_func(10)
In [52]: scipy.sparse.linalg.spsolve(A,f)
C:\Users\14256\miniconda3\lib\site-packages\scipy\sparse\linalg\_dsolve\linsolve.py:214: SparseEfficiencyWarning: spsolve requires A be CSC or CSR matrix format
warn('spsolve requires A be CSC or CSR matrix format',
Out[52]:
array([ 0.8987385, -0.1012615, 1. , 0.1012615, 0.8987385,
1.797477 , 1. , 0.1012615, 0.8987385, 1.797477 ,
1. ])
使用自:A1
our_func
In [53]: scipy.sparse.linalg.spsolve(A1,f)
C:\Users\14256\miniconda3\lib\site-packages\scipy\sparse\linalg\_dsolve\linsolve.py:276: MatrixRankWarning: Matrix is exactly singular
warn("Matrix is exactly singular", MatrixRankWarning)
Out[53]: array([nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan])
更正稀疏警告会给出相同的结果:
In [54]: scipy.sparse.linalg.spsolve(coo_matrix(A).tocsr(),f)
Out[54]:
array([ 0.8987385, -0.1012615, 1. , 0.1012615, 0.8987385,
1.797477 , 1. , 0.1012615, 0.8987385, 1.797477 ,
1. ])
我不知道这是否正确,至少它不是单一的。看起来应该有 3 个对角线,“-1,?, -1”。A
A
0赞
Matt Haberland
11/6/2023
#2
我不明白我哪里出了问题。
矩阵是奇异的,因此您尝试求解的线性系统没有唯一的解。您可以看到这一点,因为倒数第二列的所有元素都为 和 为零(或基本上为零)。例如:A
N=10
N=100
A, f = our_func(10)
A.todense()[:, -2].T
# matrix([[0.00000000e+000, 0.00000000e+000, 0.00000000e+000,
# 0.00000000e+000, 0.00000000e+000, 0.00000000e+000,
# 0.00000000e+000, 0.00000000e+000, 0.00000000e+000,
# 0.00000000e+000, 4.69807001e-294]])
# 4.69807001e-294 is essentially 0
当遇到(数值)奇异矩阵时,它会发出警告并返回 NaN 数组,而不是引发错误。sparse.linalg.spsolve
from scipy import sparse
A = sparse.csc_matrix((3, 3))
b = np.zeros(3)
sparse.linalg.spsolve(A, b)
# array([nan, nan, nan])
要解决此问题,请检查用于生成矩阵的代码,并确保其按预期工作。如果您确定这是您尝试求解的线性系统,则需要确定系统没有唯一解意味着什么。
评论
nan