The vortex filament model

The standard vortex filament model (VFM) describes the motion of thin vortex lines in three-dimensional space.

Vortex lines are assumed to be very thin with respect to the scales of interest, such that they can be effectively described as spatial curves.

The Biot–Savart law

In the VFM, each vortex line induces a velocity field around it given by the Biot–Savart law:

$$\mathit{v}(\mathit{x})=\frac{\Gamma }{4\pi }{\oint }_{\mathcal{C}}\frac{(\mathit{s}-\mathit{x})\times \mathrm{d}\mathit{s}}{|\mathit{s}-\mathit{x}{|}^3}$$

where $\Gamma$ is the vortex circulation (equal to $\kappa$ for quantum vortices), and $\mathit{s}$ is a point along the vortex line. Here $\mathcal{C}$ denotes the whole set of vortex lines in the system.

Mathematically, the above equation derives from the vorticity field:

$$\mathit{\omega }(\mathit{x})\equiv \mathbf{\nabla }\times \mathit{v}(\mathit{x})=\Gamma {\oint }_{\mathcal{C}}\delta (\mathit{s}-\mathit{x})\mathrm{d}\mathit{s}$$

where $\delta$ is Dirac delta function. That is, the vorticity field is singular and localised at the locations of quantum vortices.

The Biot–Savart law describes in particular the motion induced by vortex filaments on themselves and on surrounding vortex lines. The VFM thus describes the collective motion of a set of mutually-interacting vortex filaments which obey the Biot–Savart law. Note that the Biot–Savart integral is singular when evaluated at a vortex location ${\mathit{s}}^{\prime }\in \mathcal{C}$, and the integral must be desingularised by taking into account the finite thickness of the vortex core. The VFM also accounts for vortex reconnections, which occur when two vortex segments are sufficiently close to each other and which change the topology of the vortex system.

Desingularisation

In the VFM, one usually wants to evaluate the Biot–Savart law at locations $\mathit{x}={\mathit{s}}_0$ on the vortex. It is clear that the Biot–Savart integral, as written above, is singular when evaluated at a point ${\mathit{s}}_0$ on the curve.

The divergence of the Biot–Savart integral is of course unphysical, and is related to the fact that the actual thickness of vortex lines not really infinitesimal but finite. The standard way of accounting for the radius $a$ of the vortex core is to split the integral into local and non-local parts:

$$\mathit{v}({\mathit{s}}_0)=\frac{\Gamma }{4\pi }{\int }_{{\mathcal{C}}_0}\frac{(\mathit{s}-{\mathit{s}}_0)\times \mathrm{d}\mathit{s}}{|\mathit{s}-{\mathit{s}}_0{|}^3}+\frac{\Gamma }{4\pi }{\int }_{\mathcal{C}\setminus {\mathcal{C}}_0}\frac{(\mathit{s}-{\mathit{s}}_0)\times \mathrm{d}\mathit{s}}{|\mathit{s}-{\mathit{s}}_0{|}^3}={\mathit{v}}_{\text{local}}({\mathit{s}}_0)+{\mathit{v}}_{\text{non-local}}({\mathit{s}}_0)$$

Here ${\mathcal{C}}_0$ denotes a portion of the set of curves $\mathcal{C}$ which is in the neighbourhood of the point of interest ${\mathit{s}}_0$.

To illustrate this, the figure below shows a trefoil knot curve, or rather its projection on the XY plane. Note that here we represent a discretised version of the curve, where the number of degrees of freedom is finite and controlled by the positions of the markers.

Code for this figure

using CairoMakie
CairoMakie.activate!(pt_per_unit = 1.0)
Makie.set_theme!()

using VortexPasta.Filaments
using VortexPasta.Filaments: Vec3
using VortexPasta.PredefinedCurves: define_curve, TrefoilKnot

trefoil = define_curve(TrefoilKnot())
N = 32
i = (N ÷ 2) + 3
refinement = 8
colours = Makie.wong_colors()
f = Filaments.init(trefoil, ClosedFilament, N, CubicSplineMethod())

fig = Figure(; Text = (fontsize = 24,), size = (600, 500))
ax = Axis(fig[1, 1]; aspect = DataAspect(), xlabel = L"x", ylabel = L"y")
hidexdecorations!(ax; label = false, ticklabels = false, ticks = false)
hideydecorations!(ax; label = false, ticklabels = false, ticks = false)

# Plot complete filament in grey
let color = (:grey, 0.6)
    plot!(ax, f; refinement, color, linewidth = 1.5, arrows3d = (;))
    text!(
        ax, f(i ÷ 2 + 1, 0.6);
        text = L"\bar{\mathcal{C}}_{\!i}", align = (:right, :bottom), color,
        fontsize = 28,
    )
end

# Plot point of interest and local segment around it
let color = colours[2], arrowcol = colours[3]
    # Plot nodes (i - 1):(i + 1)
    scatter!(ax, getindex.(Ref(f), (i - 1):(i + 1)); color)
    text!(
        ax, f[i - 1] + Vec3(0.08, 0.0, 0.0);
        text = L"\mathbf{s}_{i - 1}", align = (:left, :center), color,
    )
    text!(
        ax, f[i] + Vec3(0.08, 0.0, 0.0);
        text = L"\mathbf{s}_i", align = (:left, :center), color,
    )
    text!(
        ax, f[i + 1] + Vec3(-0.08, 0.08, 0.0);
        text = L"\mathbf{s}_{i + 1}", align = (:right, :center), color,
    )

    # Plot local segments
    let ζs = range(0, 1; length = refinement + 1)
        lines!(ax, f.(i - 1, ζs); color, linewidth = 3)
        lines!(ax, f.(i, ζs); color, linewidth = 3)
        text!(
            f(i - 1, 0.5);
            text = L"ℓ^{-}", align = (:right, :bottom), color,
        )
        text!(
            f(i, 0.6) + Vec3(-0.08, 0, 0);
            text = L"ℓ^{+}", align = (:right, :center), color,
        )
        text!(
            ax, f(i - 1, 0.7) + Vec3(0.2, 0.0, 0.0);
            text = L"\mathcal{C}_{\!i}", align = (:left, :top), color,
            fontsize = 28,
        )
    end

    # Plot tangent and curvature vectors (assuming this is a 2D plot...)
    s⃗ = f[i]
    t̂ = f[i, UnitTangent()]
    ρ⃗ = f[i, CurvatureVector()]
    lengthscale = 1.2
    arrows2d!(ax, [s⃗[1]], [s⃗[2]], [t̂[1]], [t̂[2]]; color = arrowcol, shaftwidth = 1.5, tipwidth = 10, lengthscale)
    arrows2d!(ax, [s⃗[1]], [s⃗[2]], [ρ⃗[1]], [ρ⃗[2]]; color = arrowcol, shaftwidth = 1.5, tipwidth = 10, lengthscale)
    text!(
        ax, s⃗ + lengthscale * t̂;
        text = L"\mathbf{s}′", align = (:left, :center), color = arrowcol,
        offset = (2, -4),
    )
    text!(
        ax, s⃗ + lengthscale * ρ⃗;
        text = L"\mathbf{s}″", align = (:right, :center), color = arrowcol,
        offset = (4, 6),
    )
end

save("trefoil_local.png", fig)

Here, to evaluate the velocity induced by the trefoil vortex on its point ${\mathit{s}}_i$, we split the curve into a local part ${\mathcal{C}}_i$ (orange) and a non-local part ${\overline{\mathcal{C}}}_i=\mathcal{C}\setminus {\mathcal{C}}_i$ (grey). The non-local part is far from the singularity, so there is no need to modify the Biot–Savart integral as written above. As for the local part, we can approximate it using a Taylor expansion of the Biot–Savart integral about ${\mathit{s}}_i$ and truncating the integral at a small distance $ϵ\propto a$ from the singularity. See for instance Arms and Hama (1965).

More explicitly, from a Taylor expansion of $\mathit{s}(\xi )$ close to ${\mathit{s}}_i=\mathit{s}({\xi }_i)$, one can show that the Biot–Savart integrand is approximately

$$\frac{[\mathit{s}(\xi )-{\mathit{s}}_i]\times {\mathit{s}}^{\prime }(\xi )}{|\mathit{s}(\xi )-{\mathit{s}}_i{|}^3}\approx \frac{{\mathit{s}}_i^{\prime }\times {\mathit{s}}_i^{″}}{2|\xi -{\xi }_i|},$$

where $\xi$ is the curve arc length and derivatives (primes) are with respect to $\xi$. Note that ${\mathit{s}}_i^{\prime }$ and ${\mathit{s}}_i^{″}$ are respectively the unit tangent and curvature vectors at ${\mathit{s}}_i$ (see figure).

Now, if one integrates e.g. from $\mathit{s}({\xi }_i+ϵ)$ to ${\mathit{s}}_{i+1}=\mathit{s}({\xi }_{i+1})$, one gets:

$${\mathit{v}}_i^{+}\approx \frac{\Gamma }{4\pi }\frac{{\mathit{s}}_i^{\prime }\times {\mathit{s}}_i^{″}}{2}{\int }_{{\xi }_i+ϵ}^{{\xi }_{i+1}}\frac{\mathrm{d}\xi }{\xi -{\xi }_i}=\frac{\Gamma }{4\pi }\frac{{\mathit{s}}_i^{\prime }\times {\mathit{s}}_i^{″}}{2}\ln \left(\frac{{\ell }^{+}}{ϵ}\right),$$

where ${\ell }^{+}={\xi }_{i+1}-{\xi }_i$ is the length of the curve segment ${\mathit{s}}_i\to {\mathit{s}}_{i+1}$ (see figure). Doing something similar for the segment ${\mathit{s}}_{i-1}\to \mathit{s}({\xi }_i-ϵ)$ and adding both contributions, we obtain the local velocity:

$${\mathit{v}}_{\text{local}}({\mathit{s}}_i)=\frac{\Gamma }{4\pi }{\mathit{s}}_i^{\prime }\times {\mathit{s}}_i^{″}\ln \left(\frac{\sqrt{{\ell }^{-}{\ell }^{+}}}{ϵ}\right)=\frac{\Gamma }{4\pi }{\mathit{s}}_i^{\prime }\times {\mathit{s}}_i^{″}[\ln \left(\frac{2\sqrt{{\ell }^{-}{\ell }^{+}}}{a}\right)-\Delta ],$$

where we have taken the cut-off distance to be $ϵ=\frac{e^{\Delta }a}{2}$ (Saffman, 1993; §11.2). Here $\Delta$ is a coefficient which depends on the form of the vorticity profile within the vortex core (see the vortex ring tutorial).

In the code, the local velocity is sometimes referred to as the LIA (local induction approximation) term, as it has been historically used as a fast (and incomplete) approximation to the full Biot–Savart integral.

Streamfunction and energy

Everything that has been discussed until now applies to the velocity derived from the Biot–Savart law. Sometimes we may also be interested in the streamfunction vector $\mathit{\psi }$. In particular, the streamfunction values on vortex lines can be used to estimate the total energy of the vortex filament system.

Streamfunction from vortex lines

In three dimensions, the streamfunction is a vector field which is directly related to the velocity and vorticity by

$$\mathit{v}=\nabla \times \mathit{\psi }\Rightarrow \mathit{\omega }=\nabla \times \mathit{v}=-\nabla ²\mathit{\psi }$$

That is, if one knows the vorticity field $\mathit{\omega }$, then the streamfunction is the solution of a Poisson equation with the vorticity as source term. The solution can be written in terms of the Green's function $G(\mathit{r})=1/(4\pi r)$ associated to the 3D Poisson equation:

$$\mathit{\psi }(\mathit{x})=(G\ast \mathit{\omega })(\mathit{x})=\int G(\mathit{x}-\mathit{y})\mathit{\omega }(\mathit{y}){\mathrm{d}}^3\mathit{y}=\frac{1}{4\pi }\int \frac{\mathit{\omega }(\mathit{y})}{|\mathit{x}-\mathit{y}|}{\mathrm{d}}^3\mathit{y}$$

Using the fact that the vorticity has the form $\mathit{\omega }(\mathit{x})=\Gamma {\oint }_{\mathcal{C}}\delta (\mathit{s}-\mathit{x})\mathrm{d}\mathit{s}$ leads to an expression relating the streamfunction to the geometry of vortex lines:

$$\mathit{\psi }(\mathit{x})=\Gamma {\oint }_{\mathcal{C}}G(\mathit{x}-\mathit{s})\mathrm{d}\mathit{s}=\frac{\Gamma }{4\pi }{\oint }_{\mathcal{C}}\frac{\mathrm{d}\mathit{s}}{|\mathit{x}-\mathit{s}|}$$

Note that the usual Biot–Savart law for the velocity can be recovered by taking the curl of the above expression:

$$\mathit{v}(\mathit{x})=\nabla \times \mathit{\psi }(\mathit{x})=\Gamma {\oint }_{\mathcal{C}}\nabla G(\mathit{x}-\mathit{s})\times \mathrm{d}\mathit{s}=-\frac{\Gamma }{4\pi }{\oint }_{\mathcal{C}}\frac{(\mathit{x}-\mathit{s})\times \mathrm{d}\mathit{s}}{|\mathit{x}-\mathit{s}{|}^3}$$

where we have used $\nabla G(\mathit{r})=-\mathit{r}/(4\pi r^3)$.

Connection with kinetic energy

The kinetic energy of the system is given by

$$E=\frac{\rho }{2}\int \mathit{v}(\mathit{x})\cdot \mathit{v}(\mathit{x}){\mathrm{d}}^3\mathit{x}=\frac{\rho }{2}\int \mathit{\psi }(\mathit{x})\cdot \mathit{\omega }(\mathit{x}){\mathrm{d}}^3\mathit{x}$$

where the second equality is obtained using integration by parts and assuming we're in a periodic domain (or, alternatively, that the velocity goes to zero far from the vortices), and $\rho$ is the fluid density. As before, using $\mathit{\omega }(\mathit{x})=\Gamma {\oint }_{\mathcal{C}}\delta (\mathit{s}-\mathit{x})\mathrm{d}\mathit{s}$ leads to

$$E=\frac{\rho \Gamma }{2}{\oint }_{\mathcal{C}}\mathit{\psi }(\mathit{s})\cdot \mathrm{d}\mathit{s}$$

That is, the energy is simply the line integral of the tangential component of the streamfunction along vortex filaments.

Note that, if one uses the expression for $\mathit{\psi }(\mathit{x})$ from the previous section, one can directly write the total energy from the configuration of the vortex lines:

$$E=\frac{\rho {\Gamma }^2}{8\pi }{\oint }_{\mathcal{C}}{\oint }_{\mathcal{C}}^{\prime }\frac{\mathrm{d}\mathit{s}\cdot \mathrm{d}{\mathit{s}}_1}{|\mathit{s}-{\mathit{s}}_1|}$$

where the prime over the integral is to account for the desingularisation of the streamfunction integral (see section below).

Desingularisation of the streamfunction integral

As for the velocity, the Biot–Savart integral for the streamfunction is singular when evaluated on a vortex line. We can use the same approach to desingularise it, by splitting the integral onto a local and a non-local part.

Using the same notation as for the velocity, a leading-order Taylor expansion of the integrand leads to:

$${\mathit{\psi }}_i^{+}\approx \frac{\Gamma }{2\pi }\frac{{\mathit{s}}_i^{\prime }}{2}{\int }_{{\xi }_i+{ϵ}^{\prime }}^{{\xi }_{i+1}}\frac{\mathrm{d}\xi }{\xi -{\xi }_i}=\frac{\Gamma }{2\pi }\frac{{\mathit{s}}_i^{\prime }}{2}\ln \left(\frac{{\ell }^{+}}{{ϵ}^{\prime }}\right)$$

and then to:

$${\mathit{\psi }}_{\text{local}}({\mathit{s}}_i)=\frac{\Gamma {\mathit{s}}_i^{\prime }}{2\pi }\ln \left(\frac{\sqrt{{\ell }^{-}{\ell }^{+}}}{{ϵ}^{\prime }}\right)$$

Note that, here, the cut-off distance ${ϵ}^{\prime }$ is not necessarily the same as for the velocity. In fact, if one considers a system of non-colliding vortex rings, one can show that, for the system to be Hamiltonian, one needs to take ${ϵ}^{\prime }=\frac{e^{\Delta -1}a}{2}$ (as verified in the example below). This finally leads to:

$${\mathit{\psi }}_{\text{local}}({\mathit{s}}_i)=\frac{\Gamma {\mathit{s}}_i^{\prime }}{2\pi }[\ln \left(\frac{2\sqrt{{\ell }^{-}{\ell }^{+}}}{a}\right)-(\Delta -1)]$$

Example: vortex ring

For a circular vortex ring of radius $R$, one can analytically compute the non-local Biot–Savart integrals. Under the approximation that the local segments are much smaller than the vortex radius (${\ell }^{-}\sim {\ell }^{+}\ll R$), one obtains:

$$\begin{array}{rl}{\mathit{v}}_{\text{non-local}}({\mathit{s}}_i)& =\frac{\Gamma }{4\pi }{\mathit{s}}_i^{\prime }\times {\mathit{s}}_i^{″}\ln \left(\frac{8R}{{\ell }^{\ast }}\right)\\ {\mathit{\psi }}_{\text{non-local}}({\mathit{s}}_i)& =\frac{\Gamma }{2\pi }{\mathit{s}}_i^{\prime }[\ln \left(\frac{8R}{{\ell }^{\ast }}\right)-2]\end{array}$$

where we have defined ${\ell }^{\ast }=2\sqrt{{\ell }^{-}{\ell }^{+}}$. In the expression for the velocity, ${\mathit{s}}_i^{\prime }\times {\mathit{s}}_i^{″}=\hat{\mathit{z}}/R$, where $\hat{\mathit{z}}$ is the direction orthogonal to the plane where the ring lives. As a reminder, ${\mathit{s}}_i^{\prime }$ is the local unit tangent vector, while the second derivative ${\mathit{s}}_i^{″}$ corresponds to the curvature vector (oriented towards the "interior" of the curve) and whose magnitude is $1/R$.

Adding the local integrals obtained by desingularisation leads to:

$$\begin{array}{rl}\mathit{v}({\mathit{s}}_i)& =\frac{\Gamma \hat{\mathit{z}}}{4\pi R}[\ln \left(\frac{8R}{a}\right)-\Delta ]\\ \mathit{\psi }({\mathit{s}}_i)& =\frac{\Gamma {\mathit{s}}_i^{\prime }}{2\pi }[\ln \left(\frac{8R}{a}\right)-(\Delta +1)]\end{array}$$

The first expression obtained is the well-known self-induced velocity of a thin vortex ring.

From the second expression, we can easily obtain the energy associated to the vortex ring:

$$E=\frac{\rho {\Gamma }^{2}R}{2}[\ln \left(\frac{8R}{a}\right)-(\Delta +1)]$$

We can finally check that, with the desingularisation procedure described in the previous sections, the vortex ring obeys Hamilton's equation $v=\partial E/\partial p$ (Roberts and Donnelly, 1970). Here $\mathit{p}$ is the vortex impulse (Saffman, 1993):

$$\mathit{p}=\frac{\rho \Gamma }{2}{\oint }_{\mathcal{C}}\mathit{s}\times \mathrm{d}\mathit{s}$$

For a vortex ring, $\mathit{p}=\rho \Gamma \pi R^2\hat{\mathit{z}}$, that is, the impulse is directly related to the area enclosed by the vortex.

Noting that

$$\frac{\partial E}{\partial R}=\frac{\rho {\Gamma }^2}{2}[\ln \left(\frac{8R}{a}\right)-\Delta ]\text{and}\frac{\partial \mathit{p}}{\partial R}=2\rho \Gamma \pi R\hat{\mathit{z}},$$

we finally obtain that

$$\frac{\partial E}{\partial p}=\frac{\Gamma }{4\pi R}[\ln \left(\frac{8R}{a}\right)-\Delta ]=v_{\text{ring}}$$

Energy in open domains

In open (non-periodic and unbounded) domains, assuming that the velocity field goes to zero sufficiently far from the vortices, a commonly used alternative expression for the energy is (Saffman, 1993):

$$E₀=\rho \Gamma {\oint }_{\mathcal{C}}\mathit{v}\cdot (\mathit{s}\times \mathrm{d}\mathit{s})$$

which only requires knowing the velocity $\mathit{v}$ on the vortex filaments. Besides the above assumptions, this assumes that the fluid is incompressible everywhere.

However, this definition does not satisfy Hamilton's equation, meaning that this "energy" may present fluctuations in time in cases where the previous definition will be conserved.

This can be easily shown for the case of a vortex ring, in which case applying this definition to the expression for the vortex ring velocity leads to

$$E₀=\frac{\rho {\Gamma }^{2}R}{2}[\ln \left(\frac{8R}{a}\right)-\Delta ]$$

which overestimates the actual vortex ring energy by a $\rho {\Gamma }^{2}R/2$ term (and thus doesn't satisfy Hamilton's equation).

These differences might be explained by the fact that this second definition doesn't take into account the internal structure of the vortices, unlike the first definition which does this via the local contribution to the streamfunction.

References

A classical reference for the VFM is Schwarz (1985), while a more modern review is given by Hänninen and Baggaley (2014). Details on the streamfunction vector, its desingularisation and the choice of cut-off to ensure energy conservation can be found in Polanco (2025).