基于Python的CFD编程入门（2）一维非线性对流问题和误差来源

上一期基于Python的CFD编程入门（1）— 一维线性对流问题 我们介绍了较为简单的一维线性对流问题，首先离散化了相应偏微分方程，接着使用代码实现了问题的求解以及可视化，同时我们还发现了一个问题，那就是波形随着迭代计算的进行，其形状和振幅都发生了变化。

01

理论基础

在一维线性对流问题的基础上我们可以联想到一维非线性对流问题就是指速度的一阶导系数不为常数，即波形不是以由常系数确定的恒定速度运动。

一维非线性对流（1D Nonlinear convection）方程

可以看到波形以速度u进行传播，因此第二项为非线性项。

使用与线性方程同样的离散化方法：

调整格式以求解未知量：

02

代码实现

代码与一维线性对流问题的代码基本一致，但需要注意的是由于非线性问题的波速是u，因此在初始化时不能为0。

import numpy as npimport matplotlib.pyplot as pltimport sys,time
nx = 41dx = 2 / (nx - 1)nt = 25dt = 0.025
u = np.ones(nx)u[int(0.5/dx):int(1/dx)] = 2grid = np.linspace(0,2,nx)plt.plot(grid,u)
un = np.zeros(nx)for n in range(nt):    un = u.copy()    for i in range(1, nx):    # for i in range(nx):         u[i] = un[i] - un[i] * dt / dx * (un[i] - un[i - 1])    plt.plot(grid, u)

通过绘制各个时间步的波形图像，可以看到非线性对流问题同样存在与线性对流问题相同的波形变化与振幅衰减的现象，造成该现象的原因是一阶迎风方法的耗散问题。

非线性对流

线性对流

03

误差来源

由于一阶对流方程存在精确解，因此理论上来说波形应该保持不变且按照给定的速度移动，前文出现的波形变化与幅值衰减的现象简单来说是由离散方程的截断误差引起的假扩散现象，相当于引入了人工粘性(artificial viscosity)或数值粘性(numerical viscosity)。

稳定的格式会带来两种特性：耗散（dissipation）和色散（dispersion）。一阶迎风格式带来的是耗散，中心格式带来的是色散。

因此扩散更是一种物理现象的描述，而耗散则是造成能量损失的原因。

还是以一阶线性对流方程为例：

离散后为：

我们使用泰勒展开将上式还原为一阶对流方程：

其中：

可以看到还原后等号右侧多出了一项，这多出来的一项其实就是假扩散的来源——耗散项。

低阶格式引起的假扩散会比较明显。实际上，所有的数值格式都会引起假扩散，格式的阶数越高，假扩散现象可以随之减轻；同时，假扩散也会随着网格的细化而减小。

这里我们尝试增加求解线性对流问题时的网格数量（也就是nx值），再次进行求解，可以看到在平移过程中波形变化与幅值衰减的问题明显减小。

—— end ——