Gauge Conditions

Understanding and implementing gauge conditions to remove constraint degeneracies in Dedalus

When entering a system of PDEs in Dedalus, the left-hand side (LHS) of the equations is parsed into a sparse linear system. For the solver to succeed, this LHS matrix must be square and nonsingular - it must have a unique solution for any possible values on the right-hand side (RHS).

The Problem

The linear system must constrain all degrees of freedom of the variables, including gauge freedoms. Failure to do so results in singular matrices that cannot be inverted.

Example: Incompressible Flow in a Periodic Domain

Consider incompressible hydrodynamics in a fully periodic domain discretized with Fourier bases:

import dedalus.public as d3
import numpy as np

coords = d3.CartesianCoordinates('x', 'y', 'z')
dist = d3.Distributor(coords, dtype=np.float64)
bases = (d3.RealFourier(coords['x'], size=Nx, bounds=(0, Lx)),
         d3.RealFourier(coords['y'], size=Ny, bounds=(0, Ly)),
         d3.RealFourier(coords['z'], size=Nz, bounds=(0, Lz)))

# Fields
p = dist.Field(name='p', bases=bases)
u = dist.VectorField(coords, name='u', bases=bases)

# Problem - INCORRECT: This will produce singular matrices!
problem = d3.IVP([p, u], namespace=locals())
problem.add_equation("div(u) = 0")
problem.add_equation("dt(u) + grad(p) - nu*lap(u) = - u@grad(u)")

This formulation produces square matrices but they are not nonsingular. The pressure gauge is undetermined - any constant can be added to the pressure without affecting the equations.

Understanding the Degeneracy

Two related issues exist:

1. Underdetermined variable: Any constant can be added to the pressure
2. Degenerate constraint: For the mean $\vec{k} = 0$ Fourier mode, the divergence equation becomes "0 = 0"

While this seems consistent, the system must be solvable for any RHS, not just the one we entered. If we entered "div(u) = 1", the mean mode would give "0 = 1" - clearly impossible.

The Solution: Gauge Conditions with Tau Variables

We expand the system by adding a spatially-constant variable (a tau variable) to the divergence equation to absorb the degeneracy, and impose the pressure gauge as a separate equation:

# Fields
p = dist.Field(name='p', bases=bases)
u = dist.VectorField(coords, name='u', bases=bases)
tau_p = dist.Field(name='tau_p')  # Constant scalar field

# Problem - CORRECT
problem = d3.IVP([p, u, tau_p], namespace=locals())
problem.add_equation("div(u) + tau_p = 0")  # Tau absorbs degeneracy
problem.add_equation("dt(u) + grad(p) - nu*lap(u) = - u@grad(u)")
problem.add_equation("integ(p) = 0")  # Fix pressure gauge

We’ve added one degree of freedom (tau_p) and one constraint (integ(p) = 0), so the system is still square. Now it’s also nonsingular since the mean pressure is fixed.

How It Works

The tau variable mechanism works elegantly:

• The mean pressure is fixed by the gauge condition
• The degeneracy is lifted - the tau variable absorbs/acquires the mean value of any possible RHS
• The mean divergence of velocity is always zero, as required by the periodic discretization

Gauge Conditions in Bounded Domains

Similar modifications work for other types of gauges and geometries. For incompressible hydrodynamics in a bounded domain, we still need this modification so the integral of the divergence equation is compatible with specified inflow boundary conditions.

coords = d3.CartesianCoordinates('x', 'z')
dist = d3.Distributor(coords, dtype=np.float64)
xbasis = d3.RealFourier(coords['x'], size=Nx, bounds=(0, Lx))
zbasis = d3.ChebyshevT(coords['z'], size=Nz, bounds=(0, Lz))

# Fields
p = dist.Field(name='p', bases=(xbasis,zbasis))
u = dist.VectorField(coords, name='u', bases=(xbasis,zbasis))
tau_p = dist.Field(name='tau_p')
tau_u1 = dist.VectorField(coords, name='tau_u1', bases=xbasis)
tau_u2 = dist.VectorField(coords, name='tau_u2', bases=xbasis)

# Substitutions for tau method
ez = coords.unit_vector_fields(dist)[1]
lift_basis = zbasis.derivative_basis(1)
lift = lambda A: d3.Lift(A, lift_basis, -1)
grad_u = d3.grad(u) + ez*lift(tau_u1)

# Problem
problem = d3.IVP([p, u, tau_p, tau_u1, tau_u2], namespace=locals())
problem.add_equation("trace(grad_u) + tau_p = 0")  # Tau for compatibility
problem.add_equation("dt(u) - nu*div(grad_u) + grad(p) + lift(tau_u2) = - u@grad(u)")
problem.add_equation("u(z=0) = 0")
problem.add_equation("u(z=Lz) = 0")
problem.add_equation("integ(p) = 0")  # Pressure gauge

Interpretation with Inflow

From the modified equation, we can see:

• If the prescribed net inflow is nonzero, the tau variable acquires a corresponding nonzero value
• The velocity then has a spatially uniform convergence equal to this tau value
• For properly specified boundary conditions with no net inflow, the tau variable is zero and velocity is divergence-free

General Pattern

The gauge condition pattern applies to various situations:

1. Identify the gauge freedom
2. Determine which variables have underdetermined components (e.g., mean pressure, total angular momentum)
3. Find the degenerate constraint
4. Identify which equation becomes "0 = 0" or is otherwise degenerate
5. Add a tau variable
6. Introduce a constant scalar field to absorb the degeneracy in the affected equation
7. Add the gauge equation
8. Specify an additional constraint that fixes the gauge (e.g., integ(p) = 0, integ(L) = L0)

Common Gauge Conditions

Pressure gauge in incompressible flow

tau_p = dist.Field(name='tau_p')
problem.add_equation("div(u) + tau_p = 0")
problem.add_equation("integ(p) = 0")

Fixes the mean pressure to zero.

Temperature/tracer mean in periodic domains

tau_T = dist.Field(name='tau_T')
problem.add_equation("dt(T) - kappa*lap(T) + tau_T = source")
problem.add_equation("integ(T) = T_mean")

Maintains a specified mean temperature when there’s no natural reference.

Angular momentum in rotating systems

tau_L = dist.Field(name='tau_L')
problem.add_equation("dt(L) - div(stress) + tau_L = 0")
problem.add_equation("integ(L) = L0")

Fixes total angular momentum in a rotating reference frame.

Key Takeaways

Always Check Matrix Properties

The LHS matrix must be square and nonsingular - gauge freedoms break singularity even if dimensions match.

Use Tau Variables

Constant tau fields elegantly absorb constraint degeneracies while maintaining system structure.

Fix All Gauges

Every gauge freedom requires both a tau variable to lift degeneracy and an equation to fix the gauge.

Check Compatibility

Tau variables ensure integral constraints are compatible with boundary conditions and forcing.

See Also

• Tau Method - Understanding tau variables for boundary conditions
• Example scripts in the Dedalus repository demonstrate gauge conditions in various domains