一维热方程时间-温度绘图错误排查（有限差分法）

尝试使用一维热方程，以48小时时长为X轴、温度为Y轴进行绘图，采用有限差分法编写Python代码时，运行出现ValueError，提示X与Y维度不匹配。

代码如下：

import numpy as np
import matplotlib.pyplot as plt

L = 0.5    
Ta = 20     
T0 = 40      
D = 1E-4    
tau = 48*3600  
dt = 10*60   
dx = 0.001   
Nt = int(tau/dt)+1  
Nx = int(L/dx)+1    
total_time=48*3600
t=np.arange(0, total_time, 10*60)
x = [j*dt for j in range(Nx)] 
T = np.zeros((Nt,Nx))
          
for j in range(Nx):
    T[0,j] = Ta     # T(x,t=0) = Ta

T[0,0] = T0         # T(x=0,t=0) = T0

for i in range(1,Nt): 
    T[i,0] = T0     # T(x=0,t) = T0
    T[i,-1] = Ta    # T(x=L,t) = Ta

for i in range(Nt-1):
    for j in range(1, Nx-1):
        T[i+1, j] = T[i, j] + D * dt / (dx*dx) * (T[i, j+1] + T[i, j-1] - 2*T[i, j])

    plt.plot(t, T[i,:], label='temperature fluctuation')
    plt.xlabel('t (s)')
    plt.ylabel('T (°C)')
    plt.grid()
    plt.legend(loc='best')
    plt.show()

报错信息：

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[57], line 31
     28 for j in range(1, Nx-1):
     29     T[i+1, j] = T[i, j] + D * dt / (dx*dx) * (T[i, j+1] + T[i, j-1] - 2*T[i, j])
---> 31 plt.plot(t, T[i,:], label='temperature fluctuation')
     32 plt.xlabel('t (s)')
     33 plt.ylabel('T (°C)')

File ~\anaconda3\envs\py11\Lib\site-packages\matplotlib\pyplot.py:2812, in plot(scalex, scaley, data, *args, **kwargs)
   2810 @_copy_docstring_and_deprecators(Axes.plot)
   2811 def plot(*args, scalex=True, scaley=True, data=None, **kwargs):
-> 2812     return gca().plot(
   2813         *args, scalex=scalex, scaley=scaley,
   2814         **({"data": data} if data is not None else {}), **kwargs)

File ~\anaconda3\envs\py11\Lib\site-packages\matplotlib\axes\_axes.py:1688, in Axes.plot(self, scalex, scaley, data, *args, **kwargs)
   1445 """
   1446 Plot y versus x as lines and/or markers.
   1447 
   (...)
   1685 (``'green'``) or hex strings (``'#008000'``).
   1686 """
   1687 kwargs = cbook.normalize_kwargs(kwargs, mlines.Line2D)
-> 1688 lines = [*self._get_lines(*args, data=data, **kwargs)]
   1689 for line in lines:
   1690     self.add_line(line)

File ~\anaconda3\envs\py11\Lib\site-packages\matplotlib\axes\_base.py:311, in _process_plot_var_args.__call__(self, data, *args, **kwargs)
    309     this += args[0],
    310     args = args[1:]
--> 311 yield from self._plot_args(
    312     this, kwargs, ambiguous_fmt_datakey=ambiguous_fmt_datakey)

File ~\anaconda3\envs\py11\Lib\site-packages\matplotlib\axes\_base.py:504, in _process_plot_var_args._plot_args(self, tup, kwargs, return_kwargs, ambiguous_fmt_datakey)
    501     self.axes.yaxis.update_units(y)
    503 if x.shape[0] != y.shape[0]:
--> 504     raise ValueError(f"x and y must have same first dimension, but "
    505                      f"have shapes {x.shape} and {y.shape}")
    506 if x.ndim > 2 or y.ndim > 2:
    507     raise ValueError(f"x and y can be no greater than 2D, but have "
    508                      f"shapes {x.shape} and {y.shape}")

ValueError: x and y must have same first dimension, but have shapes (288,) and (501,)

核心错误点 #

1. 绘图维度完全不匹配：T[i,:]是某一时刻下空间维度的温度分布，长度为501；而t是时间数组，长度为288。两者维度不对应，不能直接用来绘图。如果要做"时间-温度"的图，应该选取固定空间位置的温度序列，即T[:,j]（j为固定空间索引）。
2. 空间坐标计算错误：代码中x = [j*dt for j in range(Nx)]用时间步长计算空间坐标，完全错误，应该用空间步长dx。
3. 绘图逻辑冗余：原代码在时间循环内每次迭代都调用plt.show()，会弹出数百个独立窗口，不符合可视化需求。

修正后的代码示例 #

示例1：绘制固定位置（x=0处）温度随时间变化 #

import numpy as np
import matplotlib.pyplot as plt

L = 0.5    
Ta = 20     
T0 = 40      
D = 1E-4    
tau = 48*3600  
dt = 10*60   
dx = 0.001   
Nt = int(tau/dt)+1  
Nx = int(L/dx)+1    
total_time=48*3600

# 修正时间数组：确保长度和Nt一致
t = np.arange(0, total_time+dt, dt)
# 修正空间坐标计算
x = np.arange(0, L+dx, dx)

T = np.zeros((Nt,Nx))
# 简化初始条件赋值
T[0,:] = Ta
T[0,0] = T0        

# 边界条件赋值
for i in range(1,Nt): 
    T[i,0] = T0     
    T[i,-1] = Ta    

# 有限差分迭代计算
for i in range(Nt-1):
    for j in range(1, Nx-1):
        T[i+1, j] = T[i, j] + D * dt / (dx**2) * (T[i, j+1] + T[i, j-1] - 2*T[i, j])

# 绘制x=0处温度随时间变化
plt.plot(t, T[:,0], label='x=0处温度')
plt.xlabel('时间 (s)')
plt.ylabel('温度 (°C)')
plt.grid()
plt.legend(loc='best')
plt.show()

示例2：绘制多个时刻的空间温度分布 #

# 接上述计算代码，添加以下绘图部分
# 选取关键时刻
time_indices = [0, int(Nt/4), int(Nt/2), int(3*Nt/4), Nt-1]
labels = ['t=0h', 't=12h', 't=24h', 't=36h', 't=48h']

for idx, label in zip(time_indices, labels):
    plt.plot(x, T[idx,:], label=label)

plt.xlabel('位置 (m)')
plt.ylabel('温度 (°C)')
plt.grid()
plt.legend(loc='best')
plt.show()

内容的提问来源于stack exchange，提问作者Elisa909