[dongxuelian2]

Related to ControlSystems.jl#1081, specifically the place_knvd rank and iteration tolerances.

place_knvd can treat a rank-deficient Float32 input matrix as full rank and either return the wrong closed-loop poles or throw from its full-rank shortcut. It also stops on any increase in cond(X), although KNV method 0 optimizes a different quantity. This change corrects the rank precision and stops on the determinant objective, with an additional consistency check before projecting complex feedback to real feedback.

Reproduction

using ControlSystemsBase, LinearAlgebra
A = Float32[0 1; 0 0]
B = Float32[1 1 1; 2 2 2]
p = Float32[-1, -2]
K = place(A, B, p)
eigvals(A - B*K)

On master 2737b535, this produces approximately -0.8 ± 0.4im. The computed second singular value is rounding noise, but 100eps(Float64) counts it as an independent input direction. The corrected result is approximately [-2, -1], with K ≈ fill(1/3, 3, 2). Replacing B with ones(Float32,2,2) instead exposes a SingularException on master; it now follows the rank-one reduction. Both problems reproduce with input scales 1e-6, 1, and 1e6.

The existing complex-pole example with state dimension five is a second regression: master stops after six sweeps when the condition number increases, leaving a maximum absolute pole error of about 0.0221. The continued iteration assigns these poles to near working precision. Its unit-test tolerance is tightened from 0.025 to 1e-10.

Numerical choices

• Rank: use max(100,size(B)...)*eps(eltype(svdvals(B)))*σmax. Precision follows the SVD itself: promoting A does not improve a Float32 factorization of B. Keep the existing safety factor rather than unconditionally replacing it with the matrix dimension. For example, with A=[0.0 1;0 0], B=1e6*[1.0 1;2 2+32eps()], and poles [-1.0,-2], the dimension-only cutoff retains a marginal direction, generates enormous gains, and gives an error around 0.03; the existing safety factor gives a well-scaled rank-one solution. Account for dimensions above that floor, and reject zero/empty numerical rank explicitly.
• Stopping: KNV method 0 maximizes abs(det(X)) one column at a time; it does not monotonically minimize cond(X). Monitor logabsdet(X)[1], require a finite value, and compare the absolute change with sqrt(eps(real(T))). Differences of logarithms measure relative determinant changes without determinant underflow. Keep condition-number logging when verbose.
• Real feedback: for real A,B, determinant stagnation alone can precede accurate real feedback. This includes real target poles with init=:rand, since random initialization uses complex eigenvectors. Before stopping, also require norm(imag(M)) ≤ n*eps(real(T))*norm(M) for M=X*Λ/X, with finite norm. Otherwise the final real(...) projection can change the assigned poles. This is a numerical consistency check, not an eigenvalue-error guarantee.

The column updates, admissible subspaces, gain recovery, public keywords, and 200-sweep cap/warning are unchanged. acker and the outer place gain-realness assertion are unchanged. No dependency is added.

References: KNV (1985), Tits–Yang's determinant interpretation, and SciPy's determinant-based stopping documentation. No convergence theorem for real KNV is asserted for this implementation's complex extension.

Validation

Follow-up correction on Linux x86_64 (demeter2), Julia 1.13.1: all 267 synthesis tests pass. The new regression tests cover real target poles with random initialization for Float32 and Float64 using three seeds. All six pole-assignment checks failed before extending the realness condition to real targets. An additional 240-case comparison covering both precisions, real/conjugate targets, and default/random initialization confirmed that all eight regressions identified during review were resolved, with no previously passing case failing the fixed reporting thresholds (1e-3 for Float32 and 1e-9 for Float64).

The following full-suite, stress-test, and benchmark results were obtained for the original implementation at bd03f252, before the follow-up correction.

Windows x86_64, AMD Ryzen 7 7840H, Julia 1.10.12 / 1.13.2, OpenBLAS ILP64 with one BLAS thread.

• Julia 1.10.12: the complete synthesis test file passes 255 tests.
• Julia 1.10.12: GROUP=ControlSystems julia --project -e 'using Pkg; Pkg.test()' passes 1156 tests.
• Julia 1.13.2: julia --project=lib/ControlSystemsBase -e 'using Pkg; Pkg.test()' passes 6454 tests, with 6 existing broken tests, zero failures, and 5 Aqua checks passing (10m10s test execution).
• git diff --check passes.

The regressions cover Float32/Float64, mixed precision, rank-one and rank-two redundant inputs, square singular inputs, uniform input scaling, weak independent directions, the near-dependent Float64 gain regression, zero/no-input rejection, an uncontrollable example, and conjugate/repeated poles. Existing tolerances were not relaxed.

Separately, MersenneTwister(1081) generated 100 systems of orders 4/6/8 with two inputs, a dependent third column, and real/conjugate targets, for both types and three scales: 1800 cases. Pole error uses one-to-one Hungarian matching and is normalized by the largest target magnitude. Reporting thresholds, fixed across before/after runs, are 1e-3 for Float32 and 1e-9 for Float64.

Group (300 cases each)	Master successes	New successes	Median sweeps, master → new
Float32, real targets	300	300	6.5 → 10
Float64, real targets	300	300	6.5 → 22.5
Float32, dependent columns	0	300	22 → 10
Float64, dependent columns	300	300	6 → 22
Float32, conjugate targets	234	300	5 → 61.5
Float64, conjugate targets	66	288	5 → 183.5

No previously passing case failed the same reporting threshold after the change. The reference was the independent Varga Schur assignment in the existing MatrixPencils dependency, evaluated in Float64 on the supplied input values, with the appropriate feedback sign. Neither exact gain agreement nor exactness of the floating-point reference was assumed. Eight weak-controllability and two repeated-pole adversarial cases were investigated separately.

BenchmarkTools 1.8.0, Julia 1.10.12 (100 samples, one evaluation/sample; warmed implementations, one BLAS thread, no other Julia jobs running). Fixed-seed real-target examples, median milliseconds:

State dimension	Type	Master	New	Master/new sweeps
4	Float32	0.100	0.135	7 / 14
4	Float64	0.098	0.370	7 / 37
8	Float32	0.170	0.847	3 / 18
8	Float64	0.185	2.007	3 / 43
16	Float32	0.777	28.403	3 / 135
16	Float64	1.258	10.788	4 / 39

The extra sweeps have a real cost; this is not a speedup claim. The Float32 example with state dimension 16 is already beyond reliable working precision (cond(X) about 8e6, relative pole error about 0.058 after the change versus 0.111 before). It is retained in the performance comparison rather than excluded. Neither stopping rule fixes that precision limitation.

Limits and compatibility

These remain heuristic rank/stagnation decisions. Weak physical input directions near the rank cutoff may be discarded. Uniform scaling is tested away from overflow/underflow; invariance under arbitrary state/input coordinate scaling is not claimed.

Continuing past condition-number increases takes more sweeps. In the complex stress cases, 55/300 Float32 and 137/300 Float64 runs reached the existing cap, including some accurate assignments that did not satisfy the realness check. Twelve Float64 complex cases still missed the reporting threshold (worst relative error 4.55e-4); the Schur reference assigned them successfully. The existing sum-of-basis initialization can also fail for repeated poles. This PR does not redesign that initialization or claim to repair all complex/ill-conditioned assignments.

The initial top-level Julia 1.13 invocation failed during standard-library Markdown precompilation with StackOverflowError, before the package tests. The successful top-level run above uses Julia 1.10.12; the Base suite uses Julia 1.13.2. The follow-up correction was additionally tested on Linux as described above; macOS execution has not been performed locally. Upstream CI is separate.

AI assistance: OpenAI Codex assisted with the research, implementation, tests, numerical audit, and this description. All reported local test, stress, and benchmark results were actually executed.

---

Follow-up correction: host demeter2, Codex session 01a15eb9-1234-7a03-a4be-7f8814d35ca2, working directory /tmp/controlsystems-review-1091.

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 92.29%. Comparing base (2737b53) to head (88ab540).
⚠️ Report is 1 commits behind head on master.

Additional details and impacted files

@@            Coverage Diff             @@
##           master    #1091      +/-   ##
==========================================
+ Coverage   92.23%   92.29%   +0.05%     
==========================================
  Files          42       42              
  Lines        6157     6177      +20     
==========================================
+ Hits         5679     5701      +22     
+ Misses        478      476       -2     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.