python有限差分法求解一维热传导方程

1、方程及其离散

1.1一维热传导方程

$\left\{\begin{matrix} \frac{\partial u\left(x,t \right )}{\partial x}=\frac{\partial u\left(x,t \right )}{\partial x^2}\\ u\left(x,0 \right )=4x\left(3-x \right ) \\ u\left(0,t \right )=0 \\ u\left(3,t \right )=0 \end{matrix}\right.$

$1\leqslant t\leqslant 1000,0\leqslant x\leqslant 1$

1.2离散化

设定步长 $h=0.1,\tau =0.0001$ ，依据上述方程得到递推关系：

2、python求解实现

import numpy as np
import matplotlib.pyplot as plth = 0.1#空间步长
N =30#空间步数
dt = 0.0001#时间步长
M = 10000#时间的步数
A = dt/(h**2) #lambda*tau/h^2
U = np.zeros([N+1,M+1])#建立二维空数组
Space = np.arange(0,(N+1)*h,h)#建立空间等差数列，从0到3，公差是h#边界条件
for k in np.arange(0,M+1):U[0,k] = 0.0U[N,k] = 0.0#初始条件
for i in np.arange(0,N):U[i,0]=4*i*h*(3-i*h)#递推关系
for k in np.arange(0,M):for i in np.arange(1,N):U[i,k+1]=A*U[i+1,k]+(1-2*A)*U[i,k]+A*U[i-1,k]
#不同时刻的温度随空间坐标的变化
plt.plot(Space,U[:,0], 'g-', label='t=0',linewidth=1.0)
plt.plot(Space,U[:,3000], 'b-', label='t=3/10',linewidth=1.0)
plt.plot(Space,U[:,6000], 'k-', label='t=6/10',linewidth=1.0)
plt.plot(Space,U[:,9000], 'r-', label='t=9/10',linewidth=1.0)
plt.plot(Space,U[:,10000], 'y-', label='t=1',linewidth=1.0)
plt.ylabel('u(x,t)', fontsize=20)
plt.xlabel('x', fontsize=20)
plt.xlim(0,3)
plt.ylim(-2,10)
plt.legend(loc='upper right')
plt.show()