我的输出T_results和c_results在最后的值为0,这是我不希望的,所以我将这些数组切片,排除它们的最后一列,得到它们各自的数组T_results_cleaned = T_results[:, :-1]
和c_results_cleaned = c_results[:, :-1]
。然而,现在当我用这些切片数组绘图时,我得到以下错误消息TypeError: Dimensions of C (50, 49) should be one smaller than X(50) and Y(50) while using shading='flat' see help(pcolormesh)
下面是我的完整代码,是解决一组方程使用有限差分。
import numpy as np
import matplotlib.pyplot as plt
import time
import pandas as pd
nz = 50
dz = 1 / nz
nt = 50
t = np.linspace(0, 4, nt)
p = np.array([0.01, 0.5, 1.0])
# Initialize concentration and temperature arrays
C = np.zeros(nz)
T = np.zeros(nz)
# Set initial condition for concentration
C[0] = 1.0
# Initialize arrays to store results
c_results = np.zeros((nt, nz))
T_results = np.zeros((nt, nz))
# Define the Arrhenius rate constant function
def k(temperature):
arrh_expr = 4.68 * np.exp(-806 / ((600 - 273) * temperature + 273))
return arrh_expr
# Calculate time derivative at specific spatial points
def rhsC(c0, c1, c2, T1, dz, p):
return (p[0] * (c2 - 2 * c1 + c0) / dz**2 -
p[1] * (c2 - c0) / (2 * dz) -
p[2] * k(T1) * c1**2)
def rhsT(T0, T1, T2, c1, dz, p):
return (0.01 * (T2 - 2 * T1 + T0) / dz**2 -
0.5 * (T2 - T0) / (2 * dz) +
2.45 * k(T1) * c1**2)
# Time-stepping loop
for it in range(nt):
c_results[it, :] = C
T_results[it, :] = T
T_results_cleaned = T_results[:, :-1] #select all rows and all columns except the last column
c_results_cleaned = c_results[:, :-1] #select all rows and all columns except the last column
for iz in range(1, nz - 1):
# Calculate concentration and temperature derivatives
dc_dt = rhsC(C[iz - 1], C[iz], C[iz + 1], T[iz], dz, p)
dT_dt = rhsT(T[iz - 1], T[iz], T[iz + 1], C[iz], dz, p)
# Update concentration and temperature using forward Euler
C[iz] += dz * dc_dt
T[iz] += dz * dT_dt
# Create a mesh for plotting
z = np.linspace(0, 1, nz)
# Plot concentration and temperature profiles
plt.pcolormesh(z, t, c_results_cleaned)
plt.colorbar(label='Concentration')
plt.xlabel('z')
plt.ylabel('t')
plt.title('Concentration Profile')
plt.show()
plt.pcolormesh(z, t, T_results_cleaned)
plt.colorbar(label='Temperature')
plt.xlabel('z')
plt.ylabel('t')
plt.title('Temperature Profile')
plt.show()
我也试着在plt模块中对z和t的数组进行切片,plt.pcolormesh(z[:-1], t[:-1], c_results_cleaned)
来匹配维度,但效果不好。我该如何解决这个问题?
1条答案
按热度按时间7vux5j2d1#
为了扩展评论,如果我正确理解了你的问题,如果使用完整的数据,该图将工作:
但据我所知,你关心的是最后的零栏。如果你必须删除它,你还需要从
z
变量中删除一行,如下所示:得到这样的图: