黎曼问题

黎曼问题（Riemann problem）是守恒律的一种基本初值问题：在某一界面左、右分别给定常状态，求撤去初始间断后出现的波形。它把一个局部接口问题与整个区域的演化连接起来，也是 Godunov 有限体积方法确定界面通量的基础。它与复分析中的黎曼映射、积分中的黎曼和是不同概念。

本条主要讨论一维标量守恒律，先区分凸、凹通量下的波型，再在交通流模型中从一个接口算到三个网格单元的一步更新。方程组的黎曼问题可以含多个特征家族，不能用同一套标量公式逐分量求解。

两个常状态构成的初值

设通量 $f$ 二阶连续可微，一维方程及初值为 $$\begin{array}{cc}{u}_t+f(u{)}_x& =0,\\ u(x,0)& =\{\begin{array}{cc}a,& x<0,\\ b,& x>0.\end{array}\end{array}$$ 常数 $a,b$ 分别为左、右状态。界面也可位于别的 $x_{\ast }$，只需把坐标换成 $x-x_{\ast }$。对于含 $m$ 个守恒量的系统，$a,b$ 是 $m$ 维向量；本条计算中的 $a,b$ 均为实数。

黎曼问题得名于 Bernhard Riemann 在1860年对有限振幅平面波的研究。后来以分段常数重构数值解时，每个相邻网格的交界处都自然出现一个这样的局部问题。原始问题的定义与数值用途相关，但“一个接口的精确解”和“多网格演化的数值近似”是两个层次。

当 $a=b$ 时，常函数就是解。两态不同时，需使用弱解与熵解的含义允许并筛选间断；仅在间断以外代回微分方程，不足以确定唯一波形。

自相似性与波型选择

方程及这组初值在同时缩放 $x,t$ 后保持不变：若 $u(x,t)$ 是解，$u(\lambda x,\lambda t)$ 在 $\lambda >0$ 时具有相同的台阶初值。在有界标量熵解的唯一性框架下，两者相同。因此正时间解可写为 $$u(x,t)=R(x/t).$$ 比值 $\xi =x/t$ 标记传播速度。原点发出的每条射线对应自相似剖面的一个位置。

若 $f^{″}>0$ 在两态之间成立，通量严格凸，特征速度 $f^{\prime }(u)$ 随状态递增：

• $a<b$ 时，两侧特征展开，解为稀疏波，扇区内 $R(\xi )=(f^{\prime }{)}^{-1}(\xi )$。
• $a>b$ 时，特征汇聚，解为压缩激波，速度由兰金—雨贡纽条件给出。

若 $f^{″}<0$，则特征速度随状态递减，以上的状态顺序反转：$a<b$ 产生激波，$a>b$ 产生稀疏波。两种曲率下，激波均满足 $$\begin{array}{cc}s& =\frac{f(b)-f(a)}{b-a},\\ f^{\prime }(a)& >s>f^{\prime }(b),\end{array}$$ 而稀疏扇的边界按空间顺序为 $x=f^{\prime }(a)t$、$x=f^{\prime }(b)t$，必须有 $f^{\prime }(a)<f^{\prime }(b)$。这比记忆“左态大就出现某种波”可靠，因为同一个状态顺序在凸、凹模型中的结果不同。若通量在两态之间改变凸性，解可能含多个波段，上述二分规则不再完整。

交通流：相同密度对的两种接口

在无量纲的交通流模型中，令 $\rho \in [0,1]$ 为相对于拥挤密度归一化的车辆密度，并假设车速为 $1-\rho$，则流量为 $$\begin{array}{cc}f(\rho )& =\rho (1-\rho ),\\ f^{\prime }(\rho )& =1-2\rho .\end{array}$$ 通量严格凹。车速 $1-\rho$ 与特征速度 $1-2\rho$ 不同；前者描述车辆，后者描述密度信息的传播。这里是教学用的理想化模型，下面的数值不对应道路观测。

先取左态 $a=1/5$、右态 $b=4/5$。通量两端相等： $$f(1/5)=f(4/5)=4/25.$$ 因此激波速度为零，低密度与高密度的界面保持在原处。两侧特征速度分别为 $3/5,-3/5$，都指向界面，符合熵方向。界面传递的流量为 $4/25=0.16$。

若先取两态平均值 $(a+b)/2=1/2$，再代入通量，则得到 $f(1/2)=1/4$。这个流量与黎曼解不符。驻留的跳跃在左右有两个不同的迹值，但二者流量相同；界面处单点怎样赋值不影响守恒，实际通量取两侧共同的 $4/25$。

金色虚线连接两端相等流量，斜率零对应驻激波。曲线顶点的流量较大，却不能替代这组初值的接口流量。

再交换左右状态，取 $a=4/5,b=1/5$。两侧特征速度变为 $-3/5,3/5$，形成跨过零速度的稀疏扇。在扇区内解 $f^{\prime }(\rho )=x/t$： $$\rho (x,t)=\{\begin{array}{cc}4/5,& x/t\le -3/5,\\ (1-x/t)/2,& -3/5<x/t<3/5,\\ 1/5,& x/t\ge 3/5.\end{array}$$ 在原界面 $x=0$，正时间状态始终为 $1/2$，流量为 $1/4$。这一回的最大流量由展开的波形选出，并非因为所有接口都应取两态平均。

本图把上一组初值的左右状态交换。零速度射线位于稀疏扇内部，因而接口通量为曲线容量 1/4。

Godunov 通量的来源

令 $G(a,b)$ 表示把左态 $a$、右态 $b$ 放在一个接口后，由其熵黎曼解决定的流量。自相似性使 $x=0$ 处的流量不随正时间改变。当该处没有驻激波时，直接取 $G(a,b)=f(R(0))$；若有驻激波，则取两侧相等的通量。

对本条讨论的严格凸或严格凹标量通量，可统一写成 $$G(a,b)=\{\begin{array}{cc}\underset{a\le v\le b}{min}f(v),& a\le b,\\ \underset{b\le v\le a}{max}f(v),& a>b.\end{array}$$ 这一公式不是先平均左右状态。以凸通量的上升初值为例：若整个稀疏扇向右走，零速度处是左态，最小值也在左端；若扇全向左走，则取右态与右端最小值；若扇跨过零速度，则取 $f^{\prime }(v)=0$ 的内部最小点。凸通量下降初值产生激波；激波速度的符号等于 $f(a)-f(b)$ 的符号，所以零速度射线选择的是两端流量中较大的一端。对凹通量，激波与稀疏波的位置互换，同理得到上述最小、最大规则。

交通通量的凹性使计算更简单：当 $a\le b$ 时，区间最小值在端点；当 $a>b$ 时，检查区间是否含容量点 $1/2$。记 $f_a=a(1-a)$、$f_b=b(1-b)$，则 $$G(a,b)=\{\begin{array}{cc}min\{f_a,f_b\},& a\le b,\\ 1/4,& b\le \frac{1}{2}\le a,a>b,\\ max\{f_a,f_b\},& a>b,\frac{1}{2}\notin [b,a].\end{array}$$ 前两组例子分别得到 $G(1/5,4/5)=4/25$ 和 $G(4/5,1/5)=1/4$。若两态相同，$G(a,a)=f(a)$，这是界面公式对常状态的相容性检查。

从三个局部问题更新三个单元

将一维区域分成等宽单元 $C_j$，宽度为 $\Delta x$。用 $U_j^n$ 表示第 $n$ 个时刻的单元平均密度，在每个单元内重构成常数。相邻单元交界处便有一组黎曼初值 $U_j^n,U_{j+1}^n$。界面通量为 $$F_{j+1/2}^n=G(U_j^n,U_{j+1}^n).$$ 对单元积分守恒律，右边流出、左边流入，得到 Godunov 更新 $$\begin{array}{cc}{U}_j^{n+1}& =U_j^n-\lambda (F_{j+1/2}^n-F_{j-1/2}^n),\\ \lambda & =\frac{\Delta t}{\Delta x}.\end{array}$$ 两个相邻单元使用同一份界面通量，求和时内部流出、流入逐项相消。这是离散守恒的直接原因。黎曼解在有限体积方法中提供通量，而不是把一张波形图当作新的单元平均值。

现在采用周期边界，三个单元依次为 $$(U_1^0,U_2^0,U_3^0)=(0.8,0.2,0.2),$$ 并取 $\Delta x=1$、$\Delta t=1/2$。周期条件使第三格的右邻居回到第一格。

接口方向	左、右密度	Godunov 通量
$C_1\to C_2$	$0.8,0.2$	$0.25$
$C_2\to C_3$	$0.2,0.2$	$0.16$
$C_3\to C_1$	$0.2,0.8$	$0.16$

把每格的流出减流入分别代入，得到 $$\begin{array}{cc}{U}_1^1& =0.8-\frac{1}{2}(0.25-0.16)=0.755,\\ U_2^1& =0.2-\frac{1}{2}(0.16-0.25)=0.245,\\ U_3^1& =0.2-\frac{1}{2}(0.16-0.16)=\mathrm{0.2.}\end{array}$$ 初始和更新后的密度和均为 $1.2$；乘以共同的格宽才是区域内归一化的总车辆量。第一格流出较多，平均密度降低；第二格接收这份差额；第三格本步净收支为零。

上下两行采用相同的密度高度比例。外框高度表示密度一，填充高度表示单元平均密度；上方箭头标记流量，首尾重复的 0.16 是同一周期接口。

时间步与单元平均的含义

交通模型在整个允许区间上有 $|f^{\prime }(\rho )|\le 1$。为让相邻接口发出的波在这一步内不交叠，本条采用便于直接证明的充分条件 $$\frac{\Delta t}{\Delta x}\le \frac{1}{2}.$$ 每个接口的波传播距离不超过 $\Delta t$；两边最多合占 $2\Delta t\le \Delta x$，所以本步各局部黎曼波形能在每格内部相接。界面流量在这段时间保持相应的 $G$，时间积分便是 $\Delta tG$。因此更新式正好给出本步分段常数初值演化后的单元平均。该条件是这里使用的保守充分条件，不是所有 Godunov 实现的必要步长上限。

本例中实际波速的最大绝对值为 $0.6$，实际 CFL 数为 $0.6\times 0.5=0.3$，也符合这个非交叠要求。局部熵解取值位于两态之间，所以本步均值仍在 $[0,1]$ 内。下一步把这些均值重新当作单元内常数，会丢掉上一时刻单元内部的波形；多步数值结果因而一般是近似，而非原始连续初值问题的逐点精确解。

Python 实现与边界

以下代码只实现通量 $f(\rho )=\rho (1-\rho )$、周期边界和上述保守时间步条件。输入为非空的密度序列、正格宽与正时间步；返回一个新列表，不修改原序列。每个接口只求一次通量，时间与额外存储均为 $O(n)$。有限性检查同时拒绝 NaN、无穷和超出 $[0,1]$ 的密度。

Python 3

from math import isfinite


def _density(value):
    value = float(value)
    if not isfinite(value) or not 0.0 <= value <= 1.0:
        raise ValueError("density must be finite and between 0 and 1")
    return value


def traffic_flux(left, right):
    a, b = _density(left), _density(right)
    fa, fb = a * (1.0 - a), b * (1.0 - b)
    if a <= b:
        return min(fa, fb)
    if b <= 0.5 <= a:
        return 0.25
    return max(fa, fb)


def godunov_step(values, dx, dt):
    q = [_density(value) for value in values]
    if not q:
        raise ValueError("at least one cell is required")
    dx, dt = float(dx), float(dt)
    if not (isfinite(dx) and isfinite(dt) and dx > 0 and dt > 0):
        raise ValueError("dx and dt must be finite and positive")
    ratio = dt / dx
    if not isfinite(ratio) or not 0.0 < ratio <= 0.5:
        raise ValueError("this implementation requires 0 < dt/dx <= 0.5")
    n = len(q)
    flows = [traffic_flux(q[j], q[(j + 1) % n]) for j in range(n)]
    return [q[j] - ratio * (flows[j] - flows[(j - 1) % n])
            for j in range(n)]


if __name__ == "__main__":
    before = [0.8, 0.2, 0.2]
    after = godunov_step(before, dx=1.0, dt=0.5)
    print([round(value, 3) for value in after])
    print(round(sum(before), 12), round(sum(after), 12))

运行输出为 [0.755, 0.245, 0.2]，两次总和均为 1.2。浮点表示可能在末位产生微小差异，不能用强行截断密度或事后改总量来掩盖通量错误。

变式。 若周期网格所有密度都相同，每个接口的通量相同，流出与流入相消，更新应保持常数。只有一个周期单元时也有这一性质。空序列没有接口定义，本实现报错；左右状态分别为零和一时通量为零，而从一到零的释放接口通量为 $1/4$，两种输入不能交换处理。

非周期道路还需给出与入流、出流相配的边界数据；把首尾相连不能代表实际道路端点。一般方程组的精确黎曼求解可能需要迭代求中间状态，近似黎曼求解器可降低代价，但必须另外核对物理可容许性、稳定性和熵相容。上述交通代码不处理变道、路口、车辆惯性或空间变化通量。

相关条目

• 稀疏波、激波、兰金—雨贡纽条件：黎曼解中的基本波形与守恒关系。
• 熵解、双曲型方程：解的选择与传播结构。
• 交通流模型、微分方程数值解：模型与离散计算。