Crank Nicolson方法求解一维非齐次热传导方程的代码,非齐次项为x*exp(t)-6*x,初值条件为x.^3+x,左边界条件为0,右边界条件为1+exp(t)
时间: 2024-04-15 07:30:25 浏览: 20
下面是使用 Crank-Nicolson 方法求解一维非齐次热传导方程的代码示例:
```python
import numpy as np
import matplotlib.pyplot as plt
# 设置参数和网格
L = 1.0 # 空间范围
T = 1.0 # 时间范围
nx = 50 # 空间网格数
nt = 1000 # 时间网格数
dx = L / nx # 空间步长
dt = T / nt # 时间步长
# 初始化网格
x = np.linspace(0, L, nx)
t = np.linspace(0, T, nt)
u = np.zeros((nx, nt))
# 设置初始条件
u[:, 0] = x**3 + x
# 设置边界条件
u[0, :] = 0
u[-1, :] = 1 + np.exp(t)
# 定义非齐次项函数
def f(x, t):
return x * np.exp(t) - 6 * x
# 构建三对角矩阵
A = np.zeros((nx, nx))
A[0, 0] = 1
A[-1, -1] = 1
for i in range(1, nx-1):
A[i, i-1] = -0.5 * dt / (dx**2)
A[i, i] = 1 + dt / (dx**2)
A[i, i+1] = -0.5 * dt / (dx**2)
# 迭代求解
for n in range(1, nt):
b = np.zeros(nx)
b[0] = u[0, n-1]
b[-1] = u[-1, n-1] + np.exp(t[n])
for i in range(1, nx-1):
b[i] = 0.5 * dt / (dx**2) * (u[i-1, n-1] - 2*u[i, n-1] + u[i+1, n-1]) + 0.5 * dt * f(x[i], t[n]) + u[i, n-1]
u[:, n] = np.linalg.solve(A, b)
# 绘制结果
X, T = np.meshgrid(x, t)
fig = plt.figure(figsize=(8, 6))
ax = fig.gca(projection='3d')
surf = ax.plot_surface(X, T, u.T, cmap='viridis')
ax.set_xlabel('x')
ax.set_ylabel('t')
ax.set_zlabel('u')
plt.show()
```
这段代码使用了 numpy 和 matplotlib 库。它首先设置了参数和网格,然后初始化了网格。接着定义了非齐次项函数,并构建了一个三对角矩阵 A。在迭代求解过程中,使用追赶法(`np.linalg.solve`)求解线性系统 Ax = b,并更新解向量 u。
最后,使用 matplotlib 绘制了结果的三维图形。请确保安装了所需的库,然后可以运行此代码以获得结果。希望对你有帮助!如有任何疑问,请随时提问。