Potential Flow & Flow Visualization

Before computers could solve the full Navier-Stokes equations, engineers designed wings, ship hulls, and nozzles using potential flow theory. By assuming the fluid is ideal — inviscid and incompressible — the governing equations collapse into the elegant Laplace equation, which is linear. Linearity is the gift that keeps giving: we can add simple flows together to build complex ones.

In this lesson we'll learn the velocity potential and stream function, study elementary flows, and superpose them to recreate the classic flow over a cylinder. Along the way we'll use NumPy and Matplotlib to draw streamlines.

Ideal (Inviscid, Incompressible) Flow

A flow is called ideal when:

• Inviscid — viscosity is negligible ($\mu = 0$), so there are no shear stresses.
• Incompressible — density is constant ($\rho = \text{const}$).
• Irrotational — fluid particles do not spin; the vorticity is zero.

These assumptions are surprisingly good outside boundary layers and wakes, which is exactly where lift on a wing is generated.

The vorticity $\vec{\omega} = \nabla \times \vec{V}$ measures local rotation. Irrotational flow means:

$$\vec{\omega} = \nabla \times \vec{V} = 0$$

Velocity Potential $\phi$

If a flow is irrotational, the velocity can be written as the gradient of a scalar velocity potential $\phi$:

$$\vec{V} = \nabla \phi \quad\Longrightarrow\quad u = \frac{\partial \phi}{\partial x}, \qquad v = \frac{\partial \phi}{\partial y}$$

This automatically satisfies irrotationality because the curl of a gradient is always zero.

Stream Function $\psi$

For 2D incompressible flow, the continuity equation $\partial u/\partial x + \partial v/\partial y = 0$ is satisfied automatically if we define a stream function $\psi$:

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

Lines of constant $\psi$ are streamlines — the fluid never crosses them. The difference $\psi_2 - \psi_1$ equals the volume flow rate (per unit depth) between two streamlines.

Equipotential lines and streamlines are orthogonal everywhere.

The Laplace Equation

For flow that is both incompressible and irrotational, both functions satisfy Laplace's equation:

$$\nabla^2 \phi = \frac{\partial^2 \phi}{\partial x^2} + \frac{\partial^2 \phi}{\partial y^2} = 0, \qquad \nabla^2 \psi = 0$$

Because Laplace's equation is linear, the sum of two solutions is also a solution. This is the principle of superposition — the engine of potential flow.

Elementary Flows

A handful of building-block solutions, written most cleanly in polar coordinates $(r, \theta)$:

A source ($m>0$) emits fluid radially outward; a sink ($m<0$) absorbs it. A vortex swirls fluid around a point with circulation $\Gamma$. A doublet is the limit of a source and sink brought infinitely close together.

Superposition Examples

• Source + Uniform stream $\rightarrow$ a Rankine half-body (like a blunt nose).
• Source + Sink + Uniform stream $\rightarrow$ a Rankine oval.
• Doublet + Uniform stream $\rightarrow$ flow over a cylinder.
• Doublet + Uniform stream + Vortex $\rightarrow$ a lifting cylinder (generates lift!).

Visualizing Source + Uniform Flow (Rankine Half-Body)

Let's superpose a source at the origin with a uniform stream and plot the streamlines. The stagnation point sits upstream where the source outflow exactly cancels the free stream.

import numpy as np
import matplotlib.pyplot as plt

# Parameters
U = 1.0          # free-stream velocity (m/s)
m = 5.0          # source strength (m^2/s)

# Cartesian grid
x = np.linspace(-4, 4, 400)
y = np.linspace(-3, 3, 300)
X, Y = np.meshgrid(x, y)

# Polar coordinates
r = np.sqrt(X**2 + Y**2)
theta = np.arctan2(Y, X)

# Stream function: uniform + source
psi = U * Y + (m / (2 * np.pi)) * theta

# Velocity field (analytic)
u = U + (m / (2 * np.pi)) * X / (r**2)
v = (m / (2 * np.pi)) * Y / (r**2)

# Stagnation point: x_s = -m/(2*pi*U) on the x-axis
x_stag = -m / (2 * np.pi * U)
print(f"Stagnation point at x = {x_stag:.3f} m")
print(f"Half-body max half-width (far field) = {m/(2*U):.3f} m")

plt.figure(figsize=(9, 6))
plt.streamplot(X, Y, u, v, density=1.4, color='steelblue', linewidth=0.8)
# Highlight the dividing streamline (the body surface): psi = m/2
plt.contour(X, Y, psi, levels=[m / 2], colors='crimson', linewidths=2)
plt.plot(x_stag, 0, 'ko', markersize=8, label='Stagnation point')
plt.title('Source + Uniform Stream: Rankine Half-Body')
plt.xlabel('x'); plt.ylabel('y'); plt.legend()
plt.axis('equal'); plt.grid(alpha=0.2)
plt.show()

Stagnation point at x = -0.796 m
Half-body max half-width (far field) = 2.500 m

The red curve ($\psi = m/2$) is the body surface — the dividing streamline that separates fluid coming from the source from fluid coming from the free stream.

Flow Over a Cylinder (Doublet + Uniform Stream)

Adding a doublet of strength $\kappa = 2\pi U a^2$ to a uniform stream produces flow over a cylinder of radius $a$. The stream function is:

$$\psi = U\,r\sin\theta\left(1 - \frac{a^2}{r^2}\right)$$

On the cylinder surface ($r=a$) the tangential velocity is $V_\theta = -2U\sin\theta$, peaking at twice the free-stream speed at the top and bottom, and the surface pressure coefficient is:

$$C_p = 1 - 4\sin^2\theta$$

import numpy as np
import matplotlib.pyplot as plt

U = 1.0          # free-stream velocity
a = 1.0          # cylinder radius

x = np.linspace(-3, 3, 500)
y = np.linspace(-2, 2, 350)
X, Y = np.meshgrid(x, y)
r = np.sqrt(X**2 + Y**2)
theta = np.arctan2(Y, X)

# Mask the inside of the cylinder
mask = r < a

# Velocity components (uniform + doublet)
u = U * (1 - (a**2 / r**2) * np.cos(2 * theta))
v = -U * (a**2 / r**2) * np.sin(2 * theta)
u[mask] = np.nan
v[mask] = np.nan

plt.figure(figsize=(9, 6))
plt.streamplot(X, Y, u, v, density=1.6, color='steelblue', linewidth=0.8)
circle = plt.Circle((0, 0), a, color='lightgray', ec='black', zorder=5)
plt.gca().add_patch(circle)
plt.title('Potential Flow Over a Cylinder (Doublet + Uniform)')
plt.xlabel('x'); plt.ylabel('y')
plt.axis('equal'); plt.xlim(-3, 3); plt.ylim(-2, 2)
plt.show()

# Surface pressure coefficient
theta_s = np.linspace(0, 2 * np.pi, 361)
Cp = 1 - 4 * np.sin(theta_s)**2
print(f"Cp at front/rear stagnation (theta=0,180): {1 - 4*np.sin(0)**2:.1f}")
print(f"Cp at top/bottom (theta=90):              {1 - 4*np.sin(np.pi/2)**2:.1f}")

Cp at front/rear stagnation (theta=0,180): 1.0
Cp at top/bottom (theta=90):              -3.0

Notice the front-back symmetry of the streamlines. Integrating this symmetric pressure gives zero drag — the famous d'Alembert's paradox. Real cylinders have drag because viscosity creates a wake, which potential theory ignores.

Kutta-Joukowski Lift

To get lift, add a vortex (circulation $\Gamma$) to the cylinder flow. The flow loses top-bottom symmetry: faster on top, slower on the bottom, so by Bernoulli the pressure is lower on top and a net upward force appears. The Kutta-Joukowski theorem gives the lift per unit span:

$$L' = \rho\, U\, \Gamma$$

This single, beautiful result underlies all of classical aerodynamics — a spinning ball curves (the Magnus effect) for exactly the same reason.

import numpy as np

rho = 1.225      # air density (kg/m^3)
U   = 50.0       # free-stream speed (m/s)
a   = 0.5        # cylinder radius (m)

# Circulation from a rotating cylinder spinning at omega (rad/s)
omega = 60.0
Gamma = 2 * np.pi * a**2 * omega     # bound circulation for a spinning cylinder

L_prime = rho * U * Gamma            # Kutta-Joukowski lift per unit span

print(f"Circulation  Gamma = {Gamma:7.2f} m^2/s")
print(f"Lift per span  L'  = {L_prime:7.1f} N/m")

# Stagnation points move down to: sin(theta) = -Gamma / (4*pi*U*a)
arg = -Gamma / (4 * np.pi * U * a)
if abs(arg) <= 1:
    theta_stag = np.degrees(np.arcsin(arg))
    print(f"Stagnation points at theta = {theta_stag:.1f} deg and {180 - theta_stag:.1f} deg")
else:
    print("Stagnation points have lifted off the cylinder surface.")

Circulation  Gamma =   94.25 m^2/s
Lift per span  L'  =  5772.7 N/m
Stagnation points at theta = -17.5 deg and 197.5 deg

Common Pitfalls

• Confusing $\phi$ and $\psi$. Streamlines are constant-$\psi$ lines and follow the flow; equipotentials are constant-$\phi$ lines and are perpendicular to them.
• Forgetting d'Alembert's paradox. Potential flow predicts zero drag on any closed body. Real drag needs viscosity (boundary layers, separation) — covered in the boundary-layer lesson.
• Sign of source vs sink. A positive $m$ is a source (outflow); a negative $m$ is a sink. Mixing the sign flips the whole picture.
• Dividing by $r=0$. Elementary flows are singular at their center. Always mask or avoid the origin when plotting.
• Using potential flow inside boundary layers. It is only valid in the inviscid outer region.

Key Takeaways

• Ideal flow is inviscid, incompressible, and irrotational — and outside boundary layers it is a great model.
• The velocity potential $\phi$ exists for irrotational flow; the stream function $\psi$ exists for 2D incompressible flow. Both satisfy Laplace's equation.
• Because Laplace's equation is linear, complex flows are built by superposing uniform streams, sources/sinks, doublets, and vortices.
• Doublet + uniform stream = flow over a cylinder, with $C_p = 1 - 4\sin^2\theta$ and zero drag (d'Alembert's paradox).
• Adding circulation produces lift via the Kutta-Joukowski theorem: $L' = \rho U \Gamma$.

Next, we'll leave incompressible flow behind and explore what happens when fluid speeds approach the speed of sound.