计算科学 - 一维不可压缩非定常 Couette 流显式有限差分 CFD - 吾爱随笔录

我目前正在关注 J.Anderson Jr. 的 CFD 基本应用程序，在编写我的第一个 CFD 问题时遇到了一些麻烦。正如标题所示，我正在使用显式有限差分方法求解不可压缩的 Couette 流。如何在不爆炸的情况下增加步数？我试图更改 delta_t 但它不会比我在代码中规定的更大，而不会出现混乱，任何建议都将不胜感激。哦，上板的 Bc 是 2 单位/吨

这是控制方程。无量纲化后

$\frac{\partial u }{\partial t} = \frac{1}{Re} \frac{\partial u^2 }{\partial y^2}$

并使用二阶中心差

$u_{i}^{n+1}=\frac{\Delta t}{Re}\left ( \frac{u_{i+1}^{n}-2u_{i}^{n}+u_{i -1}^{n}}{(\Delta y)^{2}} \right )+u_{i}^{n}$

像这样的抛物线方程的稳定性准则是：

$\Delta t \leq \frac{1}{2} Re \left ( \Delta y \right )^{2}$

所以这是我为中心差分方法编写的代码，参考 Lorena Barb 的 CFD 的 12 个步骤

import numpy as np
from matplotlib import pyplot as py
D=1
N=21
Rey=5000
dy=D/(N-1)
dt=Rey*(dy)**2

u= np.zeros(N)
un=np.zeros(N)
un[20]=2
u[20]=2
for n in range(N-1):
    un = u.copy()
    for i in range(1, N-1):
       u[i] = un[i] + dt/Rey/dy**2 * (un[i+1] - 2 * un[i] + un[i-1])


py.plot(u/2, np.linspace(0, 2, N)/2);

这是我得到的图表

解图是

我意识到我哪里出错了...我的矩阵都搞砸了

import numpy as np
from matplotlib import pyplot as plt

y_max=1
t_max=1
N=101
Rey=5
dy=y_max/(N-1)
dt=0.5*Rey*(dy)**2
u0=2

def Couette_Flow(dt,dy,y_max,t_max,u0,Rey):
    c=dt/Rey/dy**2
    y=np.arange(0, y_max+dy,dy)
    t=np.arange(0, t_max+dt,dt)
    r=len(y)
    s=len(t)
    u=np.zeros([s,r])
    u[:,0]=u0     
    for n in range(0,s-1):
         for j in range(1, r-1):
             u[n+1,j] = u[n,j] + c*(u[n,j+1]-2*u[n,j]+u[n,j-1])
    return y,u,r,s




y,u,r,s = Couette_Flow(dt,dy,y_max,t_max,u0,Rey)

这仍然不是完美的，尽管我认为我对方程进行无量纲化的事实扰乱了我的解决方案，但是对于雷诺数的小值，它似乎生成了我想要的图形。