import numpy as np import matplotlib.pyplot as plt # 设置模拟参数 D = 0.12# 热扩散率 L = 20 # 模拟深度 T = 365 * 10 # 模拟时间 h = 0.1 # 栅格大小 dt = 0.01 # 时间步长 N = int(T / dt) # 时间步数 M = int(L / h) + 1 # 深度格点数 K = int(T / 90*dt) # 画图次数 # 初始化温度分布 T0 = np.ones(M) * 10 T0[0] = 10 + 12 # 表面温度 T0[-1] = 11 # 深度为20米处温度 # 初始化温度矩阵 T = np.zeros((M, N)) T[:, 0] = T0 # 进行数值求解 for n in range(1, N): for i in range(1, M - 1): T[i, n] = T[i, n-1] + D * dt / h**2 * (T[i+1, n-1] - 2*T[i, n-1] + T[i-1, n-1]) # 边界条件 if T.all()!=0: T[0, n] =10 + 12 * np.sin(2 * np.pi * n * dt / T) T[-1, n] =11 else: break # 每隔90天画一次图 if n % int(T / 90*dt) == 0: plt.plot(T[:, n], np.linspace(0, L, M), label=f'{n*dt/365:.0f} year') plt.legend() # 显示温度轮廓图 plt.xlabel('Temperature (°C)') plt.ylabel('Depth (m)') plt.title('Temperature Profile') plt.gca().invert_yaxis() plt.show()这个代码为什么画不出图，有画布但是图中没有变量

经过观察代码，我发现在计算边界条件时，使用了`T.all()!=0`作为判断条件，这是错误的。因为`T`是一个`numpy`数组，`T.all()`会返回一个布尔值，表示所有元素是否都非零。而在循环过程中，`T`数组中的元素是会被赋值为零的，因此会导致循环提前结束，后面的时间步无法计算，进而导致无法正确绘制温度轮廓图。

解决这个问题的方法很简单，只需要将`T.all()!=0`改为`not np.all(T[:, n] == 0)`即可，表示判断当前时间步的温度矩阵是否全为零。

修改后的代码如下：

import numpy as np 
import matplotlib.pyplot as plt 

# 设置模拟参数 
D = 0.12 # 热扩散率 
L = 20  # 模拟深度 
T = 365 * 10  # 模拟时间 
h = 0.1  # 栅格大小 
dt = 0.01  # 时间步长 
N = int(T / dt)  # 时间步数 
M = int(L / h) + 1  # 深度格点数 
K = int(T / 90*dt)  # 画图次数 

# 初始化温度分布 
T0 = np.ones(M) * 10 
T0[0] = 10 + 12  # 表面温度 
T0[-1] = 11  # 深度为20米处温度 

# 初始化温度矩阵 
T = np.zeros((M, N)) 
T[:, 0] = T0 

# 进行数值求解 
for n in range(1, N): 
    for i in range(1, M - 1): 
        T[i, n] = T[i, n-1] + D * dt / h**2 * (T[i+1, n-1] - 2*T[i, n-1] + T[i-1, n-1]) 
    # 边界条件 
    if not np.all(T[:, n] == 0): 
        T[0, n] =10 + 12 * np.sin(2 * np.pi * n * dt / T) 
        T[-1, n] =11 
    else: 
        break 

    # 每隔90天画一次图 
    if n % int(T / 90*dt) == 0: 
        plt.plot(T[:, n], np.linspace(0, L, M), label=f'{n*dt/365:.0f} year') 
        plt.legend() 

# 显示温度轮廓图 
plt.xlabel('Temperature (°C)') 
plt.ylabel('Depth (m)') 
plt.title('Temperature Profile') 
plt.gca().invert_yaxis() 
plt.show()

运行修改后的代码，即可正确绘制温度轮廓图。