Add bare symbol expansion after all subs2 calls in onepass.jl

Implements the pattern from p_constraint_exa! boundary case across all functions
that use subs2, enabling bare symbol handling (e.g., x(t) → [x[k, j] for k ∈ 1:dim_x]).

Changes to src/onepass.jl:
1. p_mayer_exa! (lines 1030-1044): +3 lines
   - Added k = __symgen(:k)
   - Added subs for x0 and xf bare symbols

2. p_constraint_exa! path constraints (lines 797-815): +3 lines
   - Added k = __symgen(:k)
   - Added subs for xt and ut bare symbols

3. p_lagrange_exa! (lines 968-985): +4 lines
   - Added k = __symgen(:k)
   - Added subs for xt and ut bare symbols in ej1 and ej12

4. p_dynamics_coord_exa! (lines 900-920): +6 lines
   - Added k = __symgen(:k)
   - Added subs for xt and ut bare symbols in ej1, ej2, and ej12

Total: 16 lines added across 4 functions

Pattern applied:
- subs2 handles indexed cases: x[i] → y[i, j] and x[1:3] → [y[k, j] for k ∈ 1:3]
- subs handles bare symbols: x → [y[k, j] for k ∈ 1:dim]
- Same symbol k used for both state and control (sequential application)

This completes the implementation of tensor support for bare symbols in the
ExaModels backend.

🤖 Generated with [Claude Code](https://claude.com/claude-code)

Co-Authored-By: Claude Sonnet 4.5 <noreply@anthropic.com>

Author: jbcaillau

src/onepass.jl

@@ -695,10 +695,13 @@ function p_constraint_exa!(p, p_ocp, e1, e2, e3, c_type, label)
             code = :(length($e1) == length($e3) == 1 || throw("this constraint must be scalar")) # (vs. __throw) since raised at runtime
             x0 = __symgen(:x0)
             xf = __symgen(:xf)
+            k = __symgen(:k)
             e2 = replace_call(e2, p.x, p.t0, x0)
             e2 = replace_call(e2, p.x, p.tf, xf)
             e2 = subs2(e2, x0, p.x, 0)
+            e2 = subs(e2, x0, :([$(p.x)[$k, 0] for $k ∈ 1:$(p.dim_x)]))
             e2 = subs2(e2, xf, p.x, :grid_size)
+            e2 = subs(e2, xf, :([$(p.x)[$k, :grid_size] for $k ∈ 1:$(p.dim_x)]))
             concat(code, :($pref.constraint($p_ocp, $e2; lcon=($e1), ucon=($e3)))) # debug: to vectorise
         end
         (:initial, rg) => begin
@@ -797,8 +800,11 @@ function p_constraint_exa!(p, p_ocp, e1, e2, e3, c_type, label)
             ut = __symgen(:ut)
             e2 = replace_call(e2, [p.x, p.u], p.t, [xt, ut])
             j = __symgen(:j)
+            k = __symgen(:k)
             e2 = subs2(e2, xt, p.x, j)
+            e2 = subs(e2, xt, :([$(p.x)[$k, $j] for $k ∈ 1:$(p.dim_x)]))
             e2 = subs2(e2, ut, p.u, j)
+            e2 = subs(e2, ut, :([$(p.u)[$k, $j] for $k ∈ 1:$(p.dim_u)]))
             e2 = subs(e2, p.t, :($(p.t0) + $j * $(p.dt)))
             concat(
                 code,
@@ -897,14 +903,20 @@ function p_dynamics_coord_exa!(p, p_ocp, x, i, t, e)
     j1 = __symgen(:j)
     j2 = :($j1 + 1)
     j12 = :($j1 + 0.5)
+    k = __symgen(:k)
     ej1 = subs2(e, xt, p.x, j1)
+    ej1 = subs(ej1, xt, :([$(p.x)[$k, $j1] for $k ∈ 1:$(p.dim_x)]))
     ej1 = subs2(ej1, ut, p.u, j1)
+    ej1 = subs(ej1, ut, :([$(p.u)[$k, $j1] for $k ∈ 1:$(p.dim_u)]))
     ej1 = subs(ej1, p.t, :($(p.t0) + $j1 * $(p.dt)))
     ej2 = subs2(e, xt, p.x, j2)
+    ej2 = subs(ej2, xt, :([$(p.x)[$k, $j2] for $k ∈ 1:$(p.dim_x)]))
     ej2 = subs2(ej2, ut, p.u, j2)
+    ej2 = subs(ej2, ut, :([$(p.u)[$k, $j2] for $k ∈ 1:$(p.dim_u)]))
     ej2 = subs(ej2, p.t, :($(p.t0) + $j2 * $(p.dt)))
     ej12 = subs5(e, xt, p.x, j1)
     ej12 = subs2(ej12, ut, p.u, j1)
+    ej12 = subs(ej12, ut, :([$(p.u)[$k, $j1] for $k ∈ 1:$(p.dim_u)]))
     ej12 = subs(ej12, p.t, :($(p.t0) + $j12 * $(p.dt)))
     dxij = :($(p.x)[$i, $j2] - $(p.x)[$i, $j1])
     code = quote
@@ -967,11 +979,15 @@ function p_lagrange_exa!(p, p_ocp, e, type)
     j1 = __symgen(:j)
     j2 = :($j1 + 1)
     j12 = :($j1 + 0.5)
+    k = __symgen(:k)
     ej1 = subs2(e, xt, p.x, j1)
+    ej1 = subs(ej1, xt, :([$(p.x)[$k, $j1] for $k ∈ 1:$(p.dim_x)]))
     ej1 = subs2(ej1, ut, p.u, j1)
+    ej1 = subs(ej1, ut, :([$(p.u)[$k, $j1] for $k ∈ 1:$(p.dim_u)]))
     ej1 = subs(ej1, p.t, :($(p.t0) + $j1 * $(p.dt)))
     ej12 = subs5(e, xt, p.x, j1)
     ej12 = subs2(ej12, ut, p.u, j1)
+    ej12 = subs(ej12, ut, :([$(p.u)[$k, $j1] for $k ∈ 1:$(p.dim_u)]))
     ej12 = subs(ej12, p.t, :($(p.t0) + $j12 * $(p.dt)))
     code = quote
         if scheme == :euler
@@ -1028,10 +1044,13 @@ function p_mayer_exa!(p, p_ocp, e, type)
     pref = prefix_exa()
     x0 = __symgen(:x0)
     xf = __symgen(:xf)
+    k = __symgen(:k)
     e = replace_call(e, p.x, p.t0, x0)
     e = replace_call(e, p.x, p.tf, xf)
     e = subs2(e, x0, p.x, 0)
+    e = subs(e, x0, :([$(p.x)[$k, 0] for $k ∈ 1:$(p.dim_x)]))
     e = subs2(e, xf, p.x, :grid_size)
+    e = subs(e, xf, :([$(p.x)[$k, :grid_size] for $k ∈ 1:$(p.dim_x)]))
     # now, x[i](t0) has been replaced by x[i, 0] and x[i](tf) by x[i, grid_size]
     code = :($pref.objective($p_ocp, $e))
     return __wrap(code, p.lnum, p.line)