[系列]Python计算CFD问题<7>：2D对流问题
在前面的例子中，我们处理了2D线性对流问题，虽然问题是二维的，但实际上还是一维的，因为只有计算域是2D，物理变量只有X方向速度是变量。在本例中，我们来处理真正的2D对流问题。
2D对流问题可以用此方程组进行表达：
[attach]69974[/attach]
时间项采用向前差分，空间项采用先后差分。离散后的方程为：
[attach]69973[/attach]
换成迭代格式为：
[attach]69972[/attach]
初始条件为：[attach]69975[/attach]
边界条件：
[attach]69977[/attach]
编制计算程序为：
- import matplotlib.pyplot as plt
- import numpy as np
- import matplotlib.cm as cm
-
- nx=101
- ny=101
- nt=80
-
- dx=2.0/(nx-1)
- dy=2.0/(ny-1)
- sigma=0.2
- dt=sigma*dx
-
- x=np.linspace(0,2,nx)
- y=np.linspace(0,2,ny)
-
- u=np.ones((nx,ny))
- v=np.ones((nx,ny))
- un=np.ones((nx,ny))
- vn=np.ones((nx,ny))
-
- u[.5/dx:1/dx+1,.5/dy:1/dy+1]=2
- v[.5/dx:1/dx+1,.5/dy:1/dy+1]=2
-
- for n in range(nt+1):
- un=u.copy()
- vn=v.copy()
- u[1:,1:]=un[1:,1:]-(un[1:,1:]*dt/dx*(un[1:,1:]-un[0:-1,1:]))-vn[1:,1:]*dt/dy*(un[1:,1:]-un[1:,0:-1])
- v[1:,1:]=vn[1:,1:]-(un[1:,1:]*dt/dx*(vn[1:,1:]-vn[0:-1,1:]))-vn[1:,1:]*dt/dy*(vn[1:,1:]-vn[1:,0:-1])
-
- u[0,:]=1
- u[-1,:]=1
- u[:,0]=1
- u[:,-1]=1
-
- v[0,:]=1
- v[-1,:]=1
- v[:,0]=1
- v[:,-1]=1
-
- X,Y=np.meshgrid(x,y)
- plt.figure()
- im=plt.imshow(u,interpolation='bilinear',origin='lower',cmap=cm.rainbow)
- plt.contour(X,Y,v)
- plt.colorbar(im,orientation='vertical')
输出结果：
初始条件（u，v分布）
[attach]69978[/attach]
[attach]69979[/attach]
最终结果（u,v分布）：
[attach]69976[/attach]
[attach]69980[/attach]
【来源：http://nbviewer.ipython.org/gith ... ons/07_Step_5.ipynb】 |
