两个函数相加,一个圆柱便出现了 — 复势与达朗贝尔佯谬

肩部压力系数可以落到−3，环量为零时阻力积分仍会因上下对称而抵消。

两个函数加在一起，圆就闭上了。均匀流 $W=Uz$ 上再叠一个偶极子 $m/z$ ，不必先画圆， $r=R$ 自己就是一条流线。那个圆就是壁面。环量为 0 时，沿表面把压力加一圈，阻力在式子里消失。

把 U 和 m 分开推，半径 $R=\sqrt{m/U}$ 会变大或变小。把 Γ 从 0 拉开，上下流线就分开。

Uniform flow U = 1.00Doublet strength m = 1.00Circulation Γ = 0.00

W(z) = Uz + m/z + iΓ/(2π)·ln(z). When Γ=0 the picture is symmetric — d'Alembert's zero drag. Slide Γ to break the symmetry and lift appears (Kutta–Joukowski).

同时打开 stream 和 equipotential，两族曲线以直角相交。Γ 为 0 时，中间那个圆是单独的一条流线。

流线与等势线为何处处正交

二维、不可压、无粘、无旋，标量一次有两个。速度势 $\phi$ ，流函数 $\psi$ 。

u = \frac{\partial \phi}{\partial x} = \frac{\partial \psi}{\partial y}, \qquad v = \frac{\partial \phi}{\partial y} = -\frac{\partial \psi}{\partial x}

$u, v$ 是 $x, y$ 方向的速度。等势线 $\phi=\text{const.}$ 的斜率是 $-u/v$ 。流线 $\psi=\text{const.}$ 的斜率是 $v/u$ 。乘积是 $-1$ 。相交的地方都是直角。

正交不是图画的口味。解析函数的实部和虚部本来就这么放。导数不依赖于逼近方向，这一条就把 $\phi$ 和 $\psi$ 绑成一对。

复速度势 — 柯西与黎曼架的桥

用 $z = x + iy$ 写 $W(z) = \phi + i\psi$ 。复可微要求从哪个方向逼近都得到同一个值。把实轴和虚轴的极限放成一样，就得到下面两行。

\frac{\partial \phi}{\partial x} = \frac{\partial \psi}{\partial y}, \qquad \frac{\partial \phi}{\partial y} = -\frac{\partial \psi}{\partial x}

这是柯西–黎曼关系。和上面速度的定义是同一组式子。把 $\phi$ 和 $\psi$ 收成实部与虚部， $W$ 就是解析函数。这个 $W(z)$ 叫做复速度势。

导数也是速度。由 $dW = (u-iv)\,dz$ ，

\frac{dW}{dz} = u - iv

$dW/dz = u - iv$ 是复速度。它拿着真实速度矢量的共轭。

一行 — $W = Uz + m/z$ 画出圆柱#

均匀流是 $W = Uz$ 。拆开是 $\phi = Ux$ ， $\psi = Uy$ 。等势线是竖线，流线是横线。偶极子是 $W = m/z$ ，源和汇挤到间距为 0 的极限。极坐标下 $\phi = (m/r)\cos\theta$ ， $\psi = -(m/r)\sin\theta$ 。

把两项相加。

W(z) = Uz + \frac{m}{z}

极坐标里的流函数是

\psi(r,\theta) = U\!\left(r - \frac{R^2}{r}\right)\sin\theta, \qquad R \equiv \sqrt{\tfrac{m}{U}}

$R$ 是偶极子强度与来流速度之比决定的长度。 $\psi = 0$ 的地方，要么 $\sin\theta = 0$ ，也就是 $x$ 轴，要么 $r = R$ 。半径为 $R$ 的整圆是一条流线。

无粘时，把一条流线换成固体壁，外面的流动不变。把 $r = R$ 看成圆柱表面即可。均匀流加偶极子，两项之和，就是均匀来流里圆柱的无粘解。

加一行环量,升力就出现了

圆上的速度上下相同。用伯努利，压力也相同。合力是 0。这就是下一节阻力消失的原因。要打破对称，需要环量，也就是切向速度绕表面一周的积分。点涡的势是 $W = -i(\Gamma/2\pi)\ln z$ 。

W(z) = Uz + \frac{m}{z} - i\,\frac{\Gamma}{2\pi}\ln z

上侧变快，下侧变慢。压力上低下高。升力出现。单位展长 $L' = \rho U \Gamma$ 这一行，和环量把球的轨迹拧弯的那一笔是同一句。1902 年库塔、1906 年儒可夫斯基各自写下。在前面的画面里把 Γ 从 0 拉开，那份不对称就在屏幕上。

叠加 — 所有无粘流动都是加法

$\psi$ 的拉普拉斯方程 $\nabla^2\psi = 0$ 是线性的。两个解相加仍是解。十九世纪无粘解的一半，就是这个加法。

积木有四块。Uniform 是 $\psi = Uy$ 。Source 和 Sink 是 $\psi = \pm (q/2\pi)\theta$ 。Vortex 是 $\psi = -(\Gamma/2\pi)\ln r$ 。Doublet 是 $\psi = -m\sin\theta/r$ 。摆上这四块，就得到半体、兰金卵形和带升力的圆柱。模式名是 Cylinder、Rankine half-body、Rankine oval、Lifting cylinder。数字输出紧后面的画面，在同一张网格上跑这四种。

达朗贝尔佯谬 — 零阻力教给我们的事

用 $W = Uz + m/z$ 对表面压力做积分，阻力是 0。1752 年达朗贝尔自己推出来的。在水里拖一根棍子，手能感到阻力。只留无粘、无旋、不可压，这个力不会在式子里长出来。

补上这个缺口花了一个半世纪。1904 年普朗特把边界层拿了出来。粘性很小，却不是 0。贴近壁面的薄层里流动离开表面，尾迹里的涡打破前后压力的对称。那个差就是阻力。

0 不是算错了。它是拿掉粘性之后的流体给出的答案。这个答案和手感不一致，边界层的门就从这里打开。

60 行 NumPy — 沿 $W(z)$ 重画#

import numpy as np
import matplotlib.pyplot as plt
 
def cylinder_potential(z, U=1.0, m=1.0, gamma=0.0):
    """返回 W = Uz + m/z - i*gamma/(2π)*ln(z) 的 (phi, psi)。"""
    safe_z = np.where(np.abs(z) < 1e-9, 1e-9 + 0j, z)
    W = U * safe_z + m / safe_z - 1j * gamma / (2 * np.pi) * np.log(safe_z)
    return W.real, W.imag  # (phi, psi)
 
def velocity_from_W(z, U=1.0, m=1.0, gamma=0.0):
    """dW/dz = u - iv。取共轭还原 (u, v)。"""
    safe_z = np.where(np.abs(z) < 1e-9, 1e-9 + 0j, z)
    dW = U - m / safe_z**2 - 1j * gamma / (2 * np.pi * safe_z)
    return dW.real, -dW.imag
 
def cylinder_surface_pressure(theta, U=1.0, R=1.0, gamma=0.0):
    """圆柱表面 (含环量) 的伯努利压力系数 cp。"""
    u_t = -2 * U * np.sin(theta) - gamma / (2 * np.pi * R)
    return 1 - (u_t / U) ** 2  # cp = 1 - (V/U)^2
 
def integrate_drag_lift(U=1.0, R=1.0, gamma=0.0, N=1024):
    """无粘表面压力积分 → (drag, lift) / (rho U^2 R)。"""
    theta = np.linspace(0, 2 * np.pi, N, endpoint=False)
    cp = cylinder_surface_pressure(theta, U, R, gamma)
    drag = -np.trapz(cp * np.cos(theta), theta)
    lift = -np.trapz(cp * np.sin(theta), theta)
    return drag, lift
 
if __name__ == "__main__":
    U, m, gamma = 1.0, 1.0, 0.0
    R = np.sqrt(m / U)
 
    drag0, lift0 = integrate_drag_lift(U, R, gamma=0)
    drag1, lift1 = integrate_drag_lift(U, R, gamma=-3.0)
    print(f"Γ=0    : drag = {drag0:+.4e}, lift = {lift0:+.4e}")  # 两者 ~ 0
    print(f"Γ=-3.0 : drag = {drag1:+.4e}, lift = {lift1:+.4e}")  # lift = -ρUΓ
 
    x = np.linspace(-3, 3, 400)
    y = np.linspace(-2, 2, 240)
    X, Y = np.meshgrid(x, y)
    Z = X + 1j * Y
    Z_out = np.where(np.abs(Z) < R, np.nan, Z)
    _, psi = cylinder_potential(Z_out, U, m, gamma=-3.0)
 
    fig, ax = plt.subplots(figsize=(7, 4))
    ax.contour(X, Y, psi, levels=24, linewidths=0.7)
    ax.add_patch(plt.Circle((0, 0), R, color="white", ec="k"))
    ax.set_aspect("equal"); plt.show()

跑一遍，印出两行。

Γ=0    : drag = +1.7e-15, lift = -3.5e-16
Γ=-3.0 : drag = -2.0e-14, lift = +3.000e+00

把模式从 Rankine half-body 拨到 Lifting cylinder，鼻子形状的驻点消失，上下速度分开。

Each preset is a sum of elementary flows. The same machinery — superposition of ψ — builds half-bodies, ovals, and lifting cylinders.

Rankine oval 是源和汇之间那条闭合曲线。Cylinder 就是前面那个 $r = R$ 的圆。

半径是代码里写的平方根。它是偶极子强度除以来流速度之后再开平方。默认两者都是 1，半径就是 1。速度保持 1、只把强度加到 4，半径变成 2。反过来速度取 4、强度取 1，半径是 0.5。把圆变大的旋钮是强度。把速度加大，同样的强度会缩成更小的圆。

环量取 0 时，表面切向速度是角度的正弦乘上来流的两倍。符号就照式子前面那个负号。压力系数是 1 减去这个速度比的平方。驻点在角度 0 和 $\pi$ 。正弦为 0，压力系数就是 1。肩部是 $\pi/2$ 和 $3\pi/2$ 。正弦的平方是 1，压力系数落到 $-3$ 。上下是同一个数。所以用余弦加权的阻力和用正弦加权的升力都落到 0。

环量取 $-3$ ，半径和来流都取 1，切向速度上多出来的常数是 3 除以 $2\pi$ 。大约 0.477。上侧肩部的切向速度是 $-2$ 加上 0.477，即 $-1.523$ 。压力系数约为 $-1.32$ 。下侧肩部的切向速度是 2 加上 0.477，即 2.477。压力系数约为 $-5.14$ 。慢的一侧压力更高，快的一侧凹得更深。

把升力除以密度、速度的平方和半径，得到的是环量绝对值除以速度和半径，结果是 3。这和输出里印出的升力 $+3.000\times 10^{0}$ 是同一个数。代码注释把这个数写成密度、速度、环量的乘积再加一个负号。密度和速度为 1、环量为 $-3$ 时，这个乘积是 3。阻力栏里的 $+1.7\times 10^{-15}$ 和 $-2.0\times 10^{-14}$ 是把角度分成 1024 份之后梯形求和的舍入。式子里面没有能留下来当阻力的项。

不加环量就积分表面压力，阻力栏里只剩舍入

环量为 0 时，压力系数在驻点的 1 和肩部的 $-3$ 之间摆动。摆幅很大。用余弦加权的阻力仍然停在 $+1.7\times 10^{-15}$ 。上下成对把项消掉。不放进环量、只把表面压力加起来，阻力栏里就只剩舍入。

这个解里没有雷诺数。只有长度 $R$ 和速度 $U$ ，粘性被放成 0 拿掉了。这处省略要免费到什么时候，还没定。肩部 $c_p=-3$ 还守得住的时候，这个圆还能不能继续当流线，要到雷诺数把流动世界切开的那条边界再问一次。