Add MATLAB-compatible delay support (InputDelay/OutputDelay/IODelay/InternalDelay)

[pganguli]

Brings pkg-control's time-delay handling up to parity with MATLAB's Control System Toolbox: full InputDelay/OutputDelay/IODelay/InternalDelay support across the toolbox, including exact (not approximate) discretization for the common case.

What this adds

• InputDelay/OutputDelay/IODelay properties, correctly composed across series/parallel/append/feedback, and honored in c2d/d2c and time-domain simulation.
• InternalDelay: an LFT representation for delays trapped inside a feedback loop, with discretization, simulation, and frequency response.
• c2d() now matches MATLAB exactly for zoh on delays that aren't a whole number of sample periods (SISO and a well-defined MIMO case), rather than only approximating.
• pade() generalized into a full model-level function: approximates any delay (ordinary or internal, continuous or discrete, on tf/zpk/ss) with a rational filter.
• c2d's DelayModeling="state" option, including for InternalDelay.
• transpose() correctly swaps/transposes delay fields; ctranspose() on a delayed system now errors clearly rather than returning a wrong result (a conjugate transpose would need a negative delay, which doesn't exist).

Test suite grew substantially: 163 new tests were added. 10 pre-existing skips unrelated to delays remain.

What this doesn't cover yet

1. No exact pole/zero/stability analysis for InternalDelay. zero(), pole(), and stepinfo() error out on an InternalDelay system and point the user at pade() as a workaround (approximate, not exact). MATLAB computes these exactly via a transcendental root-finder; matching that is a substantial piece of work on its own and is the biggest gap left after this PR.
2. Descriptor (dss) models combined with delay are untested. Nothing actively forbids it, but no test exercises this combination, so correctness there is unknown rather than confirmed.
3. A rare class of algebraic-loop-plus-delay combination isn't handled. Specifically, topologies where closing a delayed loop would require solving an algebraic loop and the delay is zero-valued on that same path; this combination isn't covered by the current connect/feedback widening.
4. Exact discretization is zoh-only. c2d() matches MATLAB exactly for the default zoh method on fractional-sample delays; the foh, tustin, and matched methods still fall back to the older approximation for such delays. It's not yet confirmed whether MATLAB itself does anything more exact for those methods, so this may not be a real gap. Worth checking before treating it as follow-up work.
5. ctranspose() on a delayed system is unsupported by design, not an oversight. A conjugate transpose substitutes s -> -s (continuous) or z -> 1/z (discrete), which turns a delay into an advance. This is unrepresentable, since delay properties are nonnegative-only. It errors clearly instead of silently returning a wrong answer, and is not planned to change.

Illustrative examples

Each snippet below runs unmodified in Octave (with this branch's pkg-control loaded) and produces the shown output.

pkg load control;

%% Exact fractional-delay discretization (zoh)
h = tf (10, [1 3 10], 'IODelay', 0.25);
hd = c2d (h, 0.1)
    %% => tf with IODelay = 3 (samples), i.e. z^-3 *
    %%    (0.0119 z^2 + 0.0641 z + 0.0097) / (z^2 - 1.6551 z + 0.7408)
    %%    Matches MATLAB's c2d exactly.

%% InternalDelay: a delay trapped inside a feedback loop
G = ss (-1, 1, 1, 0, 'IODelay', 0.3);
L = feedback (G, 1);
freqresp (L, [0.05, 0.5, 1.7, 4, 11])
    %% => magnitudes [0.5000 0.4952 0.4483 0.2979 0.0896]
    %%    L carries the delay internally (not as InputDelay/OutputDelay);
    %%    pole(L) or zero(L) would error here, pointing at pade(L, n)
    %%    as the approximate workaround (see limitation #1 above).

%% pade() approximates any delay, ordinary or internal
L_approx = pade (L, 4);
freqresp (L_approx, [0.05, 0.5, 1.7, 4, 11])
    %% => magnitudes [0.5000 0.4952 0.4483 0.2979 0.0896]
    %%    Matches the exact InternalDelay freqresp above to 4 decimal
    %%    places. A 4th-order Pade approximation is enough here.

%% c2d DelayModeling="state" for InternalDelay
opt = c2dOptions ('DelayModeling', 'state');
Ld = c2d (L, 0.1, opt);
freqresp (Ld, [0.1, 1, 5])
    %% => magnitudes [0.4999 0.4886 0.2471], matching c2d(L, 0.1) with
    %%    the default delay-modeling option at this sample rate.

%% transpose() correctly swaps InputDelay/OutputDelay and transposes IODelay
sys = tf ({1, 1; 1, 1}, {[1 1], [1 2]; [1 3], [1 4]});
sys = set (sys, 'InputDelay', [0.1; 0.2], 'IODelay', [1, 2; 3, 4]);
sys_t = sys.';
get (sys_t, 'OutputDelay')   %% => [0.1; 0.2]  (was InputDelay before transpose)
get (sys_t, 'IODelay')       %% => [1, 3; 2, 4]  (transposed)
sys' %% => error: "lti: ctranspose: delay is not yet supported"

[pganguli]

By the way, I would consider any numerical differences in the output of functions that use the new delay parameters from equivalent MATLAB code to be bugs. Please feel free to let me know if you find any such in-scope differences and I will be happy to take a look into fixing them.

Additional discussion may be found in this thread at Octave Discourse