用代码求解瞬态扩散方程

瞬态扩散方程是一个偏微分方程，其一般形式为：

$$\frac{\partial u}{\partial t}=D\frac{\partial^2 u}{\partial x^2}$$

其中，$u$ 是随时间 $t$ 和空间 $x$ 变化的函数，$D$ 是扩散系数。

为了求解这个方程，我们需要一个初始条件和边界条件。在这里，我们设定边界条件为 $u(0,t)=0$ 和 $u(L,t)=0$，初始条件为 $u(x,0)=f(x)$。

下面是一个 Python 实现的例子，使用有限差分法进行数值求解：

import numpy as np
import matplotlib.pyplot as plt

# 设置初始条件和模拟参数
L = 1.0  # 空间长度
T = 1.0  # 时间长度
N = 100  # 离散空间点数
M = 1000  # 离散时间点数
D = 1.0  # 扩散系数
dx = L / N  # 空间步长
dt = T / M  # 时间步长

# 初始化网格和初始条件
x = np.linspace(0, L, N+1)
u = np.zeros((N+1, M+1))
u[:, 0] = np.sin(np.pi*x/L)

# 设置边界条件
u[0, :] = 0
u[N, :] = 0

# 使用有限差分法求解
for k in range(M):
    for i in range(1, N):
        u[i, k+1] = u[i, k] + D*dt/dx**2*(u[i+1, k] - 2*u[i, k] + u[i-1, k])

# 绘制结果
fig, ax = plt.subplots()
ax.plot(x, u[:, 0], label='t=0')
ax.plot(x, u[:, int(M/4)], label='t=T/4')
ax.plot(x, u[:, int(M/2)], label='t=T/2')
ax.plot(x, u[:, -1], label='t=T')
ax.legend()
ax.set_xlabel('x')
ax.set_ylabel('u')
plt.show()

该代码使用有限差分法进行数值求解，并绘制出了不同时间点的结果。你可以根据需要，修改模拟参数和边界条件来进行计算。