Fix lqi/lqi_controller: explicit augmentation, discrete integral gain, index order, LQGProblem weights - #155
Merged
Merged
Conversation
…, index order, LQGProblem weights
`add_output_integrator` built the augmentation as `tf(M)*sys` followed by a
`similarity_transform` with a hard-coded state permutation. That construction was
correct only by virtue of properties the transfer-function conversion does not
guarantee, and it produced two incorrect results in `lqi`/`lqi_controller`. It is
replaced by an explicit state-space augmentation
A_aug = [A 0; h*C[ind,:] λI], B_aug = [B; h*D[ind,:]], C_aug = [C 0; 0 ±I]
with `h = 1, λ = -ϵ` in continuous time and `h = Ts, λ = 1-ϵ` in discrete time, so
that the integrator state is a time integral of the selected outputs in both time
domains and the added outputs and states follow the order of `ind`.
Fixes resulting from the rewrite:
* Discrete-time integral gain was a factor `Ts` too small. The augmented plant's
integrator state was the unscaled running sum `Σy` (the conversion placed the
factor `Ts` in the added output row), whereas the integrator assembled by
`lqi_controller` emitted `Ts*Σe`, so `L_i` multiplied a signal `Ts` times too
large. The realized loop was the designed loop with `L_i` scaled by `Ts`; on the
example in the test suite the closed-loop pole at 0.647 was realized at 0.972.
The controller's integrator is now built from an explicit realization rather than
from `ss(tf(...))`, whose state depends on the internal scaling chosen by the
conversion.
* A non-ascending `integrator_outputs` silently produced the wrong controller. The
augmentation ordered the integrator states by output index while `lqi_controller`
ordered the reference and error channels as given, so the two disagreed whenever
the indices were not sorted. `integrator_outputs=[2,1]` on an asymmetric 2x2 plant
yielded a closed-loop pole at +1.04. The order given is now honoured throughout,
for the integrator states, the reference channels and the integrator entries of
`Q1` alike.
* `add_output_integrator` added a single output equal to the sum of the selected
integrals for a vector `ind`, rather than one output per index as documented,
while still adding one state per index.
* `lqi_controller(::LQGProblem, Qi)` spliced `prob.Q1` directly into a state penalty.
`LQGProblem.Q1` penalizes the performance output `C1*x`, so the method threw a
dimension error whenever `size(C1,1) != nx` and would have used the wrong weight
whenever `C1 != I`. It also discarded `qQ` and `SQ` silently. The augmented weights
are now formed as in `lqr(::LQGProblem)`, i.e. `blkdiag(C1'Q1*C1 + qQ*C2'C2, Qi)`
with `SQ` extended by zero rows for the integrator states.
Validation added to `lqi`: `integrator_outputs` must be non-empty, in range and free
of duplicates; `Q1` and `Q2` must be square and of the right size (only the first
dimension was checked); `length(integrator_outputs) <= nu`, since integrating more
outputs than there are control inputs leaves the augmentation unstabilizable and
previously returned a non-stabilizing gain without complaint. A transmission zero at
the integrator pole is reported as a warning, replacing the downstream error "The
Hamiltonian matrix is not dichotomic". `lqi_controller` now checks that `obs` has the
shape produced by `observer_predictor(G, K; output_state=true)`.
The docstrings stated the augmented dynamics with `-C` in the lower-left block and
`[B; 0]` as the augmented input matrix. Neither was correct: `neg=true` never affected
the state dynamics, only the sign of the added output row, so it was a no-op in `lqi`
and is no longer passed there, and a nonzero `D` feeds the control signal into the
integrator. The sign bookkeeping in the closed-loop assembly is now documented where
it happens.
Tests: `test/test_lqi.jl` is reorganized into testsets and extended. The gain returned
by `lqi` is pinned against an explicitly augmented plant, which fixes the sign and the
scaling of every channel, and each assembled loop is checked against the poles the
separation principle predicts, `eig(A_aug - B_aug*L)` together with `eig(A - K*C)`,
which is what exposed the discrete-time scaling error. Added coverage for index order,
nonzero `D`, continuous-time `ϵ > 0`, the cross-term argument, load-disturbance
rejection, the `LQGProblem` method with `C1 != I` and `qQ != 0`, and every validation
path. `test/test_augmentation.jl` gains shape, ordering, `neg` and discrete-scaling
tests for `add_output_integrator`. Of the new assertions, 17 in `test_lqi.jl` and 10
in `test_augmentation.jl` fail against the previous implementation.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## master #155 +/- ##
==========================================
- Coverage 91.87% 91.71% -0.17%
==========================================
Files 20 20
Lines 3065 3078 +13
==========================================
+ Hits 2816 2823 +7
- Misses 249 255 +6
Flags with carried forward coverage won't be shown. Click here to find out more. ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
A review of
lqiandlqi_controllerturned up two silent correctness errors, one dimension error, several missing preconditions and a docstring whose stated augmented dynamics did not match the code. The root cause of the first two is the construction used byadd_output_integrator, so that is addressed first and most of the rest follows from it.The augmentation
add_output_integratorbuilt the augmented plant astf(M)*sysfollowed by asimilarity_transformwith a hard-coded state permutation. This was correct only by virtue of properties that the transfer-function conversion does not guarantee: that it yields exactlylength(ind)states, in index order, with unit input gain. It is replaced by an explicit state-space augmentation,with
h = 1, λ = -ϵin continuous time andh = Ts, λ = 1-ϵin discrete time. The integrator state is then a genuine time integral of the selected outputs in both time domains, so the integrator entries ofQ1carry the same meaning for a continuous-time system and for its discretization, and the added states and outputs follow the order in which the indices are given. The input/output behaviour for a scalarindis unchanged.Correctness
The discrete-time integral gain was a factor
Tstoo small. The augmented plant's integrator state was the unscaled running sumΣy, because the conversion placed the factorTsin the added output row rather than in the state dynamics. The integrator assembled bylqi_controller, built fromss(tf(Ts, [1, -(1-ϵ)], Ts)), emittedTs·Σeinstead, soL_imultiplied a signalTstimes too large. The realized loop was exactly the designed loop withL_iscaled byTs; on the example already in the test suite the designed closed-loop pole at 0.647 was realized at 0.972:The controller's integrator is now built from an explicit realization, so it reproduces the state of the plant augmentation rather than only its transfer function. Note that this changes the gain
lqireturns for a discrete-time plant: the integrator columns are now1/Tstimes their previous value, which is what makes the assembled loop the designed loop.A non-ascending
integrator_outputssilently produced the wrong controller. The augmentation ordered the integrator states by output index, since it scannedi ∈ indover1:ny, whilelqi_controllerordered the reference and error channels as given. The two therefore disagreed whenever the indices were not sorted, and the gain was applied to the wrong error signals.integrator_outputs = [2,1]on an asymmetric 2×2 plant produced a closed-loop pole at +1.04. The order given is now honoured throughout — for the integrator states, the reference channels and the integrator entries ofQ1alike.add_output_integratordid not add the documented outputs. For a vectorindit added a single output, equal to the sum of the selected integrals, while still adding one state per index.lqiwas unaffected, since it uses onlyAandB, but the returned system did not match its own documentation and was unusable for anything else.lqi_controller(::LQGProblem, Qi)used the wrong weight.LQGProblem.Q1penalizes the performance outputC1*x, not the state;lqr(::LQGProblem)correspondingly formsC1'Q1*C1 + qQ*C2'C2with the cross termSQ. The method splicedprob.Q1directly into a state penalty, so it threwArgumentError: Q1 must have size 3×3wheneversize(C1,1) != nx, and would have used the wrong weight wheneverC1 != Iwith matching dimensions. It also discardedqQandSQsilently —qQ = 10andqQ = 0produced identical controllers althoughlqr(prob)differs. The augmented weights are now formed the same waylqr(::LQGProblem)forms them, withSQextended by zero rows for the integrator states.Preconditions
lqinow rejects an empty, out-of-range or duplicatedintegrator_outputs, and requiresQ1andQ2to be square and of the right size — only the first dimension was checked, so a non-square weight fell through to aDimensionMismatchfromlqr. Two further conditions govern whether the LQI problem has a solution at all, and neither was checked:length(integrator_outputs) ≤ nu. Integrating more outputs than there are control inputs leaves the augmentation unstabilizable;lqrpreviously returned a non-stabilizing gain, with a closed-loop eigenvalue at exactly 0, without complaint. This is now an error.s = -ϵ, orz = 1-ϵ). This previously surfaced asThe Hamiltonian matrix is not dichotomicfrom the Riccati solver. It is now reported as a warning naming the cause, based on the rank of the Rosenbrock matrix at the integrator pole; a warning rather than an error, since the rank test is numerical.lqi_controlleralso checks thatobshas the shape produced byobserver_predictor(G, K; output_state=true), instead of failing insideconnect.Docstrings
The
lqidocstring stated the augmented dynamics as[A 0; -C 0]withB_aug = [B; 0], and attributed the negated block toadd_output_integrator(...; neg=true). Neither part was correct.negonly ever flipped the sign of the added output row and never the state dynamics, so it was a no-op inlqiand is no longer passed there; the state integrates+Cx, which is what makes the sign bookkeeping inlqi_controllerwork out, since∫(r-y) = -xᵢatr = 0. A nonzeroDfeeds the control signal into the integrator, so the augmented input matrix is[B; D]. The sign convention is now documented where it is relied upon, together with the fact that the controller contains the negative feedback sign and the loop must therefore be closed withpos_feedback = true.No prose under
docs/was changed.Tests
test/test_lqi.jlis reorganized into testsets and extended. Two properties do the load-bearing work:lqiis pinned against an explicitly augmented plant, which fixes the sign and the scaling of every channel including the integrator channels.eig(A_aug - B_aug*L)together witheig(A - K*C). This is what exposed the discrete-time scaling error; the previous discrete assertions,isstable(minreal(Hd))anddcgain ≈ 1, are both invariant to a uniform scaling of the integral gain and so could not see it.Added coverage: index order, nonzero
D, continuous-timeϵ > 0, the cross-term positional argument, load-disturbance rejection (dcgain(G/(1+GC)) = 0, the point of the integrator), a scalarintegrator_outputs, theLQGProblemmethod withC1 != IandqQ != 0, and every validation path.@test 0 ∈ poles(C0)is replaced by a tolerance-based count, matching what the MIMO test already did, and the unusedusing Plotsis removed.test/test_augmentation.jlgains shape, ordering,negand discrete-scaling tests foradd_output_integrator.Of the new assertions, 17 in
test_lqi.jland 10 intest_augmentation.jlfail against the previous implementation.Full suite: 1493 pass, 6 broken, and one pre-existing failure unrelated to this change —
test_uncertainty.jl:42assertsrand(::Diagonal, 100) isa Matrix, butrandon aDiagonalreturns aDiagonalon the LinearAlgebra version in use. It fails identically on master.🤖 Generated with Claude Code