Multi-Agent Target Tracking via a Time-Expanded Cellular Sheaf

This example uses a time-expanded sheaf to study how coordination constraints determine the space of feasible agent trajectories. Four scenarios are compared; the nullspace dimension of the restricted sheaf Laplacian measures how many free degrees of freedom remain after constraints are imposed.

Vehicle model. Each agent or target is a planar quadrotor. The state is

\[x = [y,\, z,\, \varphi,\, \dot y,\, \dot z,\, \dot\varphi]^\top \in \mathbb{R}^6\]

where $y$ is the lateral position, $z$ is the altitude, and $\varphi$ is the roll angle. The control input is $u = [T_1, T_2]^\top \in \mathbb{R}^2$.

Sheaf structure. Each vertex $(\text{entity}, t)$ carries a stalk $\mathbb{R}^{n_x + n_u} = \mathbb{R}^8$. Three edge families encode:

Dynamics: $x_{t+1} = A_d x_t + B_d u_t$.
Consensus: inter-agent agreement on a chosen coordinate subspace.
Tracking: agent–target alignment in a chosen coordinate subspace.

By varying which coordinates are constrained and at which timesteps, we obtain four qualitatively different solution spaces whose nullspace dimensions tell the story of the engineering constraints.

The helper types (TrackingEdge, TrackingProblem, BobbingTarget, ScenarioResult) and all solver utilities live in CellularSheaves.ControlSheaves.MultiAgentTracking and can be reused by other examples without copying code.

using CellularSheaves
using CellularSheaves.ControlSheaves.MultiAgentTracking
using CellularSheaves.ControlSheaves.TrackingDSL
using CellularSheaves.TrajectorySheaves: continuous_to_discrete_zoh
using LinearAlgebra
using Plots

Planar Quadrotor Dynamics

State: $x = [y, z, \varphi, \dot y, \dot z, \dot\varphi]$. The linearisation is taken around the hover trim condition.

g     = 9.81
m_veh = 0.5
I_quad = 0.01
ell   = 0.25

Ac = [0.0  0.0   0.0   1.0  0.0  0.0;
      0.0  0.0   0.0   0.0  1.0  0.0;
      0.0  0.0   0.0   0.0  0.0  1.0;
      0.0  0.0  -g     0.0  0.0  0.0;
      0.0  0.0   0.0   0.0  0.0  0.0;
      0.0  0.0   0.0   0.0  0.0  0.0]

Bc = [0.0               0.0;
      0.0               0.0;
      0.0               0.0;
      0.0               0.0;
      1.0 / m_veh       1.0 / m_veh;
      ell / (2I_quad)  -ell / (2I_quad)]

h = 0.05
Ad, Bd = continuous_to_discrete_zoh(Ac, Bc, h)
nx = size(Ad, 1)
nu = size(Bd, 2)

State-index constants for the planar quadrotor.

IDX_Y   = 1   # lateral position y
IDX_Z   = 2   # altitude z
IDX_PHI = 3   # roll angle φ
IDX_YDT = 4   # ẏ
IDX_ZDT = 5   # ż
IDX_PHDT = 6  # φ̇

Plot Recipe

A Plots.jl recipe renders a ScenarioResult as a two-panel $y(t)$ / $z(t)$ figure. Agent trajectories use solid/dashed lines; target reference trajectories are shown as dotted lines for comparison.

@recipe function f(sr::ScenarioResult)
    layout := (1, 2)
    size := (800, 380)
    plot_title := "$(sr.label): null dim = $(sr.null_dim),  ||dz|| = $(round(sr.residual; sigdigits=3))"
    agent_colors  = [:steelblue, :darkorange, :green, :crimson]
    agent_styles  = [:solid, :dash, :dashdot, :dot]
    target_colors = [:gray, :black, :darkgreen, :purple]
    for (i, traj) in enumerate(sr.agent_trajs)
        @series begin
            subplot   := 1
            title     := "y(t)"
            xlabel    := "t (s)"
            ylabel    := "y (m)"
            label     := "A$i"
            lw        := 2
            linecolor := agent_colors[i]
            linestyle := agent_styles[i]
            sr.times, traj[:, sr.y_col]
        end
        @series begin
            subplot   := 2
            title     := "z(t)"
            xlabel    := "t (s)"
            ylabel    := "z (m)"
            label     := "A$i"
            lw        := 2
            linecolor := agent_colors[i]
            linestyle := agent_styles[i]
            sr.times, traj[:, sr.z_col]
        end
    end
    for (j, traj_j) in enumerate(sr.target_trajs)
        @series begin
            subplot   := 1
            label     := "T$j"
            lw        := 1
            linecolor := target_colors[j]
            linestyle := :dot
            sr.times, getindex.(traj_j, IDX_Y)
        end
        @series begin
            subplot   := 2
            label     := "T$j"
            lw        := 1
            linecolor := target_colors[j]
            linestyle := :dot
            sr.times, getindex.(traj_j, IDX_Z)
        end
    end
end

Common Setup

All four scenarios use two agents and two targets over a 40-step horizon (2 s). The projection matrices are defined once and reused:

$R_{yz}$: projects the $(n_x + n_u)$-stalk onto $(y, z)$ (Scenario 1).
$R_y$: projects onto $y$ only (consensus in Scenarios 2–4).
$R_z$: projects onto $z$ only (tracking in Scenarios 2–4).

Each scenario is specified as a @tracking_problem DSL program. The DSL declaratively names agents, targets, a time horizon, and the consensus/tracking coordination edges. Calling lower_tracking_program on the resolved program produces the TrackingProblem struct used by the simulation utilities. The boundary data (initial conditions, target trajectories) is still supplied externally as a plain Dict — it varies per scenario but never changes the sheaf topology.

n_agents  = 2
n_targets = 2
k         = 40
times     = h .* collect(0:k)

R_yz = state_projection_matrix([IDX_Y, IDX_Z], nx, nu)
R_y  = state_projection_matrix([IDX_Y],         nx, nu)
R_z  = state_projection_matrix([IDX_Z],         nx, nu)

1×8 Matrix{Float64}:
 0.0  1.0  0.0  0.0  0.0  0.0  0.0  0.0

Bobbing targets used in Scenarios 2–4. Both complete two full vertical cycles over the 2-second horizon.

omega_2periods = 2π * 2 / (k * h)
bt1 = BobbingTarget(0.0, 1.0, 0.3, omega_2periods)
bt2 = BobbingTarget(0.0, 2.0, 0.3, omega_2periods)
traj_bt1 = trajectory(bt1, 0:k, h, nx, nu, IDX_Y, IDX_Z, IDX_ZDT)
traj_bt2 = trajectory(bt2, 0:k, h, nx, nu, IDX_Y, IDX_Z, IDX_ZDT)

41-element Vector{Vector{Float64}}:
 [0.0, 2.0, 0.0, 0.0, 1.8849555921538759, 0.0, 0.0, 0.0]
 [0.0, 2.0927050983124844, 0.0, 0.0, 1.7926992988449335, 0.0, 0.0, 0.0]
 [0.0, 2.176335575687742, 0.0, 0.0, 1.524961107694578, 0.0, 0.0, 0.0]
 [0.0, 2.2427050983124843, 0.0, 0.0, 1.107949098294274, 0.0, 0.0, 0.0]
 [0.0, 2.285316954888546, 0.0, 0.0, 0.58248331161764, 0.0, 0.0, 0.0]
 [0.0, 2.3, 0.0, 0.0, 1.1542024162330739e-16, 0.0, 0.0, 0.0]
 [0.0, 2.285316954888546, 0.0, 0.0, -0.5824833116176402, 0.0, 0.0, 0.0]
 [0.0, 2.2427050983124843, 0.0, 0.0, -1.107949098294274, 0.0, 0.0, 0.0]
 [0.0, 2.176335575687742, 0.0, 0.0, -1.5249611076945777, 0.0, 0.0, 0.0]
 [0.0, 2.0927050983124844, 0.0, 0.0, -1.7926992988449335, 0.0, 0.0, 0.0]
 ⋮
 [0.0, 1.8236644243122582, 0.0, 0.0, -1.5249611076945784, 0.0, 0.0, 0.0]
 [0.0, 1.757294901687516, 0.0, 0.0, -1.1079490982942746, 0.0, 0.0, 0.0]
 [0.0, 1.7146830451114539, 0.0, 0.0, -0.5824833116176374, 0.0, 0.0, 0.0]
 [0.0, 1.7, 0.0, 0.0, -8.079416913631517e-16, 0.0, 0.0, 0.0]
 [0.0, 1.7146830451114539, 0.0, 0.0, 0.5824833116176391, 0.0, 0.0, 0.0]
 [0.0, 1.7572949016875157, 0.0, 0.0, 1.1079490982942735, 0.0, 0.0, 0.0]
 [0.0, 1.823664424312258, 0.0, 0.0, 1.5249611076945773, 0.0, 0.0, 0.0]
 [0.0, 1.9072949016875163, 0.0, 0.0, 1.7926992988449342, 0.0, 0.0, 0.0]
 [0.0, 1.9999999999999998, 0.0, 0.0, 1.8849555921538759, 0.0, 0.0, 0.0]

Scenario 1: Full (y,z) Coordination at Every Timestep

Both consensus and tracking edges use the $(y, z)$ projection, active at every timestep.

The two agents travel to distinct goals while each tracking its own target in full $(y, z)$. But the consensus edge demands that agents agree in $(y, z)$ at every step. Since the two targets travel different $(y, z)$ trajectories, the constraints "A1 tracks T1", "A2 tracks T2", and "A1 agrees with A2 in $(y,z)$" are mutually incompatible unless the targets coincide. The sheaf Laplacian has no kernel (null dim = 0) and harmonic_extension returns the least-squares compromise with a nonzero residual $\|dz\| > 0$.

x0_a1_s1 = [0.0, 1.0, 0.0, 0.0, 0.0, 0.0]
xk_a1_s1 = [2.0, 1.5, 0.0, 0.0, 0.0, 0.0]
x0_a2_s1 = [0.0, 2.0, 0.0, 0.0, 0.0, 0.0]
xk_a2_s1 = [2.0, 2.5, 0.0, 0.0, 0.0, 0.0]

traj_t1_s1 = generate_reference_trajectory(x0_a1_s1, xk_a1_s1, k, Ad, Bd, nx, nu)
traj_t2_s1 = generate_reference_trajectory(x0_a2_s1, xk_a2_s1, k, Ad, Bd, nx, nu)

41-element Vector{Vector{Float64}}:
 [0.0, 2.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
 [-6.395441560233875e-11, 1.999999999998919, 1.842873076363105e-12, -5.558775661995696e-12, -2.1316282072803006e-14, 2.4868995751603507e-14, -192.01431729271462, -182.57382860160047]
 [0.0003014686024541361, 1.0635296352620491, -0.1475076357969141, 0.024117498444256633, -37.4588145894315, -5.900305431946327, 206.19765846938887, 195.84551191323686]
 [0.004191449524761624, 0.1956968317459575, -0.2807706174548289, 0.14237598975755195, 2.745502448831149, 0.5697861656486812, 0.0, -0.6068901773346677]
 [0.014617368413246323, 0.3314547287430971, -0.24279865015081964, 0.27155656001700335, 2.6848134310978162, 0.9490925264828607, 4.703252740928275, 0.0]
 [0.030828352250184182, 0.4774535321492274, -0.1218556997461882, 0.3669957098306263, 3.155138705190831, 3.8886254895629455, 0.0, 5.823405178802746]
 [0.05006361761669083, 0.6497689803546947, -0.01841513117901819, 0.3939586409185184, 3.737479223071347, 0.24899725281086682, 2.4146405462275506, 0.0]
 [0.06985936825464771, 0.8426795428734619, 0.03176348999456626, 0.39376928192684174, 3.9789432776942015, 1.7581475942031817, 0.0, 0.0]
 [0.0887990111388447, 1.041626706758132, 0.11967086970504075, 0.3566300052109233, 3.9789432776942117, 1.758147594203167, 0.0, 4.223147555648007]
 [0.10493858622193444, 1.2511317395319224, 0.1415915688586964, 0.2871609810181629, 4.401258033259025, -0.8813196280768718, 0.0, 0.0]
 ⋮
 [1.9746937542254241, 2.4516789817385254, 0.30200636413369936, 0.3346733858540866, 0.21476008116180534, -1.3422505072610293, 0.0, 0.0]
 [1.987998392925356, 2.4624169857966245, 0.23489383877064862, 0.2029986110917957, 0.21476008116179807, -1.3422505072610293, 0.0, 0.0]
 [1.9955422602294373, 2.4731549898547236, 0.16778131340759814, 0.10424253002008144, 0.2147600811617903, -1.342250507261029, 0.0, 0.0]
 [1.9989712908221973, 2.483892993912821, 0.10066878804454803, 0.03840514263894364, 0.2147600811617822, -1.342250507261029, 0.0, 0.0]
 [1.9999314193881652, 2.4946309979709187, 0.03355626268149836, 0.005486448948382155, 0.21476008116177359, -1.3422505072610285, 0.0, -2.1476008116176453]
 [2.000000000000015, 2.4999999999999725, -2.5028590311393373e-14, -3.2731298106712423e-14, 0.0, 0.0, 0.0, 0.0]
 [2.0000000000000093, 2.499999999999982, -1.7839896226945484e-14, -1.832561880021899e-14, 6.661338147750939e-16, 6.046053414865415e-16, 0.0, 0.0]
 [2.0000000000000044, 2.499999999999991, -9.523631883112671e-15, -7.233221647252286e-15, 8.81239525796218e-16, 8.227874761720932e-16, -3.808103053749275e-15, -2.853235094001666e-15]
 [2.0, 2.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]

The DSL program below declares agents, targets, the time horizon, and the $(y,z)$ consensus and tracking edges active at every timestep. lower_tracking_program builds the TrackingProblem with target dynamics included and a higher tracking weight to steer the least-squares solution toward the tracking objective.

prog1 = @tracking_problem begin
    agent(a1; dynamics=(Ad, Bd), period=h)
    agent(a2; dynamics=(Ad, Bd), period=h)
    target(t1)
    target(t2)
    horizon(K)
    times(Tall = 0:K)
    consensus(c1; agents=(a1,a2), maps=(R_yz,R_yz), at=Tall)
    track(tr1; agent=a1, target=t1, maps=(R_yz,R_yz), at=Tall)
    track(tr2; agent=a2, target=t2, maps=(R_yz,R_yz), at=Tall)
end
ctx1 = Dict{Symbol,Any}(:K => k, :Ad => Ad, :Bd => Bd, :R_yz => R_yz, :h => h)
prob1 = lower_tracking_program(prog1, ctx1;
    include_target_dynamics=true, tracking_weight=5.0,
).problem

bnd1 = Dict{Int,Vector{Float64}}()
bnd1[agent_vertex(prob1, 1, 0)] = vcat(x0_a1_s1, zeros(nu))
bnd1[agent_vertex(prob1, 2, 0)] = vcat(x0_a2_s1, zeros(nu))
bnd1[agent_vertex(prob1, 1, k)] = vcat(xk_a1_s1, zeros(nu))
bnd1[agent_vertex(prob1, 2, k)] = vcat(xk_a2_s1, zeros(nu))
for t in 0:k
    bnd1[target_vertex(prob1, 1, t)] = traj_t1_s1[t + 1]
    bnd1[target_vertex(prob1, 2, t)] = traj_t2_s1[t + 1]
end

result1 = run_scenario("Scenario 1", prob1, bnd1, times;
    target_trajs = [traj_t1_s1, traj_t2_s1], y_col = IDX_Y, z_col = IDX_Z)

plot(result1)

Scenario 2: y-Consensus and z-Tracking, All Timesteps, Aligned Initial Conditions

Targets bob vertically. Consensus edges constrain agents to agree only on lateral position $y$; tracking edges constrain each agent to match its assigned target in altitude $z$. Both edge families are active at every timestep.

Because the two targets share $y = 0$ and the consensus only constrains $y$, the system is consistent: agents can agree in $y$ while independently tracking different altitudes. A non-trivial nullspace (null dim > 0) arises because $z$ (and all other state components) remain unconstrained by the consensus edges. The TrackingEdge API makes it natural to assign different projection matrices to each agent–target pair.

Since all target vertices are pinned as boundary data, target dynamics edges are omitted (include_target_dynamics = false).

The DSL program uses split projection matrices: $R_y$ for consensus (lateral agreement) and $R_z$ for tracking (altitude following). Both edge families span every timestep (at=Tall). Target dynamics are excluded because all target vertices will be pinned as boundary conditions.

prog2 = @tracking_problem begin
    agent(a1; dynamics=(Ad, Bd), period=h)
    agent(a2; dynamics=(Ad, Bd), period=h)
    target(t1)
    target(t2)
    horizon(K)
    times(Tall = 0:K)
    consensus(c1; agents=(a1,a2), maps=(R_y,R_y), at=Tall)
    track(tr1; agent=a1, target=t1, maps=(R_z,R_z), at=Tall)
    track(tr2; agent=a2, target=t2, maps=(R_z,R_z), at=Tall)
end
ctx2 = Dict{Symbol,Any}(:K => k, :Ad => Ad, :Bd => Bd, :R_y => R_y, :R_z => R_z, :h => h)
prob2 = lower_tracking_program(prog2, ctx2;
    tracking_weight=5.0,
).problem

x0_a1_s2 = [0.0, traj_bt1[1][IDX_Z], 0.0, 0.0, traj_bt1[1][IDX_ZDT], 0.0]
x0_a2_s2 = [0.0, traj_bt2[1][IDX_Z], 0.0, 0.0, traj_bt2[1][IDX_ZDT], 0.0]

bnd2 = Dict{Int,Vector{Float64}}()
bnd2[agent_vertex(prob2, 1, 0)] = vcat(x0_a1_s2, zeros(nu))
bnd2[agent_vertex(prob2, 2, 0)] = vcat(x0_a2_s2, zeros(nu))
for t in 0:k
    bnd2[target_vertex(prob2, 1, t)] = traj_bt1[t + 1]
    bnd2[target_vertex(prob2, 2, t)] = traj_bt2[t + 1]
end

result2 = run_scenario("Scenario 2", prob2, bnd2, times;
    target_trajs = [traj_bt1, traj_bt2], y_col = IDX_Y, z_col = IDX_Z)

plot(result2)

Scenario 3: y-Consensus and z-Tracking at the Last Timestep Only

The same projection matrices as Scenario 2, but coordination edges are active only at $t = k$. Agents evolve freely (dynamics only) for $t = 0, \ldots, k-1$ and must converge to the terminal constraint.

Both targets share $y = 0$ (laterally aligned), but agents start misaligned with each other and with the targets in both $y$ and $z$. Removing $k$ timesteps of coordination edges dramatically enlarges the feasible space (null dim >> Scenario 2).

Changing at=Tall to at=t[end] in the DSL restricts both the consensus and tracking edges to the single terminal timestep t = k. The only topological difference from Scenario 2 is this single keyword change.

prog3 = @tracking_problem begin
    agent(a1; dynamics=(Ad, Bd), period=h)
    agent(a2; dynamics=(Ad, Bd), period=h)
    target(t1)
    target(t2)
    horizon(K)
    consensus(c1; agents=(a1,a2), maps=(R_y,R_y), at=t[end])
    track(tr1; agent=a1, target=t1, maps=(R_z,R_z), at=t[end])
    track(tr2; agent=a2, target=t2, maps=(R_z,R_z), at=t[end])
end
ctx3 = Dict{Symbol,Any}(:K => k, :Ad => Ad, :Bd => Bd, :R_y => R_y, :R_z => R_z, :h => h)
prob3 = lower_tracking_program(prog3, ctx3;
    tracking_weight=5.0,
).problem

x0_a1_s3 = [-0.5, 0.5, 0.0, 0.0, 0.0, 0.0]
x0_a2_s3 = [ 1.0, 1.5, 0.0, 0.0, 0.0, 0.0]

bnd3 = Dict{Int,Vector{Float64}}()
bnd3[agent_vertex(prob3, 1, 0)] = vcat(x0_a1_s3, zeros(nu))
bnd3[agent_vertex(prob3, 2, 0)] = vcat(x0_a2_s3, zeros(nu))
for t in 0:k
    bnd3[target_vertex(prob3, 1, t)] = traj_bt1[t + 1]
    bnd3[target_vertex(prob3, 2, t)] = traj_bt2[t + 1]
end

result3 = run_scenario("Scenario 3", prob3, bnd3, times;
    target_trajs = [traj_bt1, traj_bt2], y_col = IDX_Y, z_col = IDX_Z)

plot(result3)

Scenario 4: Last-Timestep Constraints, Targets Not Aligned in y or z

Identical sheaf topology to Scenario 3, but Target 2 is now at $y = 1.5$ m instead of $y = 0$. The targets are not initially aligned with each other in either $y$ or $z$.

Because the sheaf topology (which edges exist at which timesteps) is unchanged, the restricted Laplacian's nullspace has the same dimension as Scenario 3. The two harmonic extensions produce different trajectories, but they live in equally large affine solution spaces, confirming that nullspace dimension is a structural property of the coordination architecture, not of the boundary data.

bt2_s4 = BobbingTarget(1.5, 2.0, 0.3, omega_2periods)
traj_bt2_s4 = trajectory(bt2_s4, 0:k, h, nx, nu, IDX_Y, IDX_Z, IDX_ZDT)

41-element Vector{Vector{Float64}}:
 [1.5, 2.0, 0.0, 0.0, 1.8849555921538759, 0.0, 0.0, 0.0]
 [1.5, 2.0927050983124844, 0.0, 0.0, 1.7926992988449335, 0.0, 0.0, 0.0]
 [1.5, 2.176335575687742, 0.0, 0.0, 1.524961107694578, 0.0, 0.0, 0.0]
 [1.5, 2.2427050983124843, 0.0, 0.0, 1.107949098294274, 0.0, 0.0, 0.0]
 [1.5, 2.285316954888546, 0.0, 0.0, 0.58248331161764, 0.0, 0.0, 0.0]
 [1.5, 2.3, 0.0, 0.0, 1.1542024162330739e-16, 0.0, 0.0, 0.0]
 [1.5, 2.285316954888546, 0.0, 0.0, -0.5824833116176402, 0.0, 0.0, 0.0]
 [1.5, 2.2427050983124843, 0.0, 0.0, -1.107949098294274, 0.0, 0.0, 0.0]
 [1.5, 2.176335575687742, 0.0, 0.0, -1.5249611076945777, 0.0, 0.0, 0.0]
 [1.5, 2.0927050983124844, 0.0, 0.0, -1.7926992988449335, 0.0, 0.0, 0.0]
 ⋮
 [1.5, 1.8236644243122582, 0.0, 0.0, -1.5249611076945784, 0.0, 0.0, 0.0]
 [1.5, 1.757294901687516, 0.0, 0.0, -1.1079490982942746, 0.0, 0.0, 0.0]
 [1.5, 1.7146830451114539, 0.0, 0.0, -0.5824833116176374, 0.0, 0.0, 0.0]
 [1.5, 1.7, 0.0, 0.0, -8.079416913631517e-16, 0.0, 0.0, 0.0]
 [1.5, 1.7146830451114539, 0.0, 0.0, 0.5824833116176391, 0.0, 0.0, 0.0]
 [1.5, 1.7572949016875157, 0.0, 0.0, 1.1079490982942735, 0.0, 0.0, 0.0]
 [1.5, 1.823664424312258, 0.0, 0.0, 1.5249611076945773, 0.0, 0.0, 0.0]
 [1.5, 1.9072949016875163, 0.0, 0.0, 1.7926992988449342, 0.0, 0.0, 0.0]
 [1.5, 1.9999999999999998, 0.0, 0.0, 1.8849555921538759, 0.0, 0.0, 0.0]

The DSL program for Scenario 4 is identical to Scenario 3. The different target trajectory (traj_bt2_s4) is boundary data supplied externally — it never appears in the DSL program and does not change the sheaf topology or nullspace dimension.

prog4 = @tracking_problem begin
    agent(a1; dynamics=(Ad, Bd), period=h)
    agent(a2; dynamics=(Ad, Bd), period=h)
    target(t1)
    target(t2)
    horizon(K)
    consensus(c1; agents=(a1,a2), maps=(R_y,R_y), at=t[end])
    track(tr1; agent=a1, target=t1, maps=(R_z,R_z), at=t[end])
    track(tr2; agent=a2, target=t2, maps=(R_z,R_z), at=t[end])
end
ctx4 = Dict{Symbol,Any}(:K => k, :Ad => Ad, :Bd => Bd, :R_y => R_y, :R_z => R_z, :h => h)
prob4 = lower_tracking_program(prog4, ctx4;
    tracking_weight=5.0,
).problem

x0_a1_s4 = [-0.5, 0.5, 0.0, 0.0, 0.0, 0.0]
x0_a2_s4 = [ 1.0, 1.5, 0.0, 0.0, 0.0, 0.0]

bnd4 = Dict{Int,Vector{Float64}}()
bnd4[agent_vertex(prob4, 1, 0)] = vcat(x0_a1_s4, zeros(nu))
bnd4[agent_vertex(prob4, 2, 0)] = vcat(x0_a2_s4, zeros(nu))
for t in 0:k
    bnd4[target_vertex(prob4, 1, t)] = traj_bt1[t + 1]
    bnd4[target_vertex(prob4, 2, t)] = traj_bt2_s4[t + 1]
end

result4 = run_scenario("Scenario 4", prob4, bnd4, times;
    target_trajs = [traj_bt1, traj_bt2_s4], y_col = IDX_Y, z_col = IDX_Z)

plot(result4)

Comparison: How Constraints Shape the Solution Space

The four scenarios are summarised below. The nullspace dimension and Laplacian residual $\|dz\|$ are displayed via the ScenarioResult show method.

Scenario	Consensus	Tracking	Active timesteps	Initial alignment
1	$(y,z)$	$(y,z)$	all $0:k$	distinct endpoints
2	$y$	$z$	all $0:k$	aligned with targets
3	$y$	$z$	terminal $k$	unaligned; targets share $y$
4	$y$	$z$	terminal $k$	unaligned; targets offset in $y$

Constraint density governs nullspace size. Scenario 1 is infeasible: the $(y,z)$ consensus edge forces agents to share the same trajectory in the $yz$-plane, but tracking edges demand they follow different targets. The result is null dim = 0 with a nonzero least-squares residual — the energy minimiser finds a compromise, not a feasible trajectory. Scenario 2 (only $y$ consensus + $z$ tracking, all timesteps) is consistent and admits a family of solutions (null dim > 0) because fewer coordinates per edge are constrained. Restricting to the terminal timestep only (Scenarios 3–4) removes $k$ sets of constraints and dramatically enlarges the solution space.

Nullspace dimension is a topological invariant. Scenarios 3 and 4 share the same sheaf topology (identical edge structure) and therefore the same nullspace dimension, even though their target trajectories differ. The harmonic extensions produce different specific trajectories, but they live in equally large affine spaces. The nullspace dimension is a property of the coordination architecture, not the boundary data.

Feasibility is separate from trajectory shape. In Scenario 2, agents start consistent with all-time constraints, so the harmonic solution closely tracks the bobbing targets. In Scenarios 3–4, agents start misaligned; the minimum-energy harmonic path converges to the terminal constraint while the large nullspace deforms it freely in the unconstrained directions.

for r in (result1, result2, result3, result4)
    show(stdout, MIME("text/plain"), r); println()
end

bar(["S1\n(y,z) all-t", "S2\ny/z all-t", "S3\ny/z last-t", "S4\ny/z last-t\n(mis-tgt)"],
    [result1.null_dim, result2.null_dim, result3.null_dim, result4.null_dim];
    ylabel = "Nullspace dimension",
    title  = "Nullspace dimension vs. coordination architecture",
    legend = false,
    color  = [:steelblue, :darkorange, :green, :crimson])