循环索引错误（索引 100 超出轴 0 的范围，大小为 100）

Question

我试图使用 for 循环通过使用不同的 Peclet 数来计算 theta 变量，但是当我将 Pe 数设置为 100 时，我的循环不起作用。你知道我的代码有什么问题吗？

这是我的代码

import numpy as np
import scipy.sparse as sp
import scipy.sparse.linalg as spl
import matplotlib.pyplot as plt
np.set_printoptions(precision=2)

rho = 1000
c = 4000
Pe = 10000
L = 2
N = 10
deta = 2 / N
k = 0.07

def comolokko(Pe):   
    x = np.linspace(-1,1,N+1)  
    y = np.linspace(-1.25+deta/2,1.25-deta/2,N)
    y2 = np.linspace(deta/2-deta/2,deta/2-deta/2,N)
    ux = -np.cos(np.pi*x/L)
    C = sp.diags([0.5*np.ones(N),0.5*np.ones(N)],[-1,0],(N+1,N)).tocsr() 
    G = sp.diags([-np.ones(N),np.ones(N)],[-1,0],(N+1,N)).tocsr() / deta
    G[0,0] = 0    # No flux at the boundary X direction
    G[-1,-1] = 0 
    C = sp.diags(ux) @ C
    I = sp.eye(N)
    I2 = np.zeros(N)
    Iy = sp.diags(np.sin(np.pi*y/L))
    Iy2 = sp.diags(np.sin(np.pi*y2/L))
    #print(Iy2.todense())
    F = sp.kron(Iy,Pe * C) -sp.kron(I,G) #flux operator for all vertical faces
    F2 = sp.kron(Iy2,Pe * C) -sp.kron(I,G)
    D = sp.diags([-np.ones(N),np.ones(N)],[0,1],(N,N+1)).tocsr()
    Ax = sp.kron(I,D) @ F
    Ax2 = sp.kron(I,D) @ F2
    ####y direc
    G = sp.diags([-np.ones(N),np.ones(N)],[-1,0],(N+1,N)).tocsr() / deta 
    G[0,0] *= 2    # Diriclet at boundary   cuz 2l ?
    G[-1,-1] *= 2  # Diriclet at boundary
    C = sp.diags([0.5*np.ones(N),0.5*np.ones(N)],[-1,0],(N+1,N)).tocsr()
    # Coordinates of the horizontal face centers
    y = np.linspace(-1,1,N+1)
    y2 = np.linspace(0,0,N+1)
    #print(y2)
    x = np.linspace(-1.25+deta/2,1.25-deta/2,N)
    x2 = np.linspace(deta/2-deta/2,deta/2-deta/2,N)
    uy = np.cos(np.pi*y/2) # Velocity part dependent on y only
    #print(uy)

    C = sp.diags(uy) @ C    # The operator u_y * phi on the faces
    C2 = sp.diags(y2) @ C
                            # to 2D
    I = sp.eye(N)
    Ix = sp.diags(np.sin(np.pi*x/2))  
    Ix3=  sp.diags(np.sin(np.pi*x2/2))
    F = sp.kron(Pe * C,Ix) - sp.kron(G,I)
    F3 = sp.kron(Pe * C2,Ix3) - sp.kron(G,I) # Flux operator for all the horizontal faces
    # Conservation operator 
    D = sp.diags([-np.ones(N),np.ones(N)],[0,1],(N,N+1)).tocsr()
    Ay = sp.kron(D,I) @ F
    Ay2 = sp.kron(D,I) @ F3
    A3 = Ay2
    A = Ax + Ay
    A2 = Ax2 + Ay2
    g = np.zeros(N)
    g[0] = G[0,0]
    b = np.kron(g,np.ones(N))
    theta = spl.spsolve(A,b)
    theta2= spl.spsolve(A2,b)
    theta3= spl.spsolve(A3,b)
    return theta, theta2

代码的第一部分可以正常工作，但是当我使用 for 循环时代码不工作。此外，我可以手动更改 comolokko 定义函数的 Pe 数，而不会出现任何索引错误。

thetaq = list()
thetaq2 = list()
Pe = [10, 100, 1000, 10000]
for n in Pe:
    theta, theta2 = comolokko(n)
    thetaq.append(theta[int(n)])
    theta, theta2 = comolokko(n)
    thetaq2.append(theta2[int(n)])
thetaq = np.array(thetaq)
thetaq2 = np.array(thetaq2)
plt.semilogy(Pe, ((thetaq/thetaq2)*2))
plt.xlabel('Peclet Number')
plt.ylabel('q/q_0')
plt.grid()
plt.show()

出现此错误

IndexError                                Traceback (most recent call last)
<ipython-input-38-3b63413cc0ab> in <module>
      4 for n in Pe:
      5     theta, theta2 = comolokko(n)
----> 6     thetaq.append(theta[int(n)])
      7     theta, theta2 = comolokko(n)
      8     thetaq2.append(theta2[int(n/2)])

IndexError: index 100 is out of bounds for axis 0 with size 100

Answer 1

theta 是一个 ndarray 或稀疏矩阵（请参阅此处的文档：https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.spsolve.html#scipy.sparse.linalg.spsolve）。

当 n=100（for 循环的第二次迭代）时，您要求 theta 矩阵中的第 100 个项目，该项目不存在。

我真的不明白你想做什么，但是，如果你只是想遍历矩阵并获得第 i 个项目，你可以做一个枚举 for 循环，如下所示：

for i, n in enumerate (Pe):
   #you code
   #where i is the i-th iteration (like a counter, i runs from 0 to len(Pe) - 1
   #and n is the Pe itself as in the list (10, 100, etc)