Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

debug: fallback if MovingHorizonEstimator arrival covariance update fails #151

Merged
merged 17 commits into from
Jan 23, 2025

Conversation

franckgaga
Copy link
Member

@franckgaga franckgaga commented Jan 22, 2025

context: See #98 (comment) point 1.

@franckgaga franckgaga changed the title added: fallback if MovingHorizonEstimator arrival covariance update fails debug: fallback if MovingHorizonEstimator arrival covariance update fails Jan 22, 2025
@codecov-commenter
Copy link

codecov-commenter commented Jan 22, 2025

Codecov Report

Attention: Patch coverage is 93.61702% with 3 lines in your changes missing coverage. Please review.

Project coverage is 98.89%. Comparing base (3ed42e4) to head (c2bb875).
Report is 7 commits behind head on main.

Files with missing lines Patch % Lines
src/controller/execute.jl 75.00% 1 Missing ⚠️
src/estimator/mhe/execute.jl 96.15% 1 Missing ⚠️
src/general.jl 85.71% 1 Missing ⚠️
Additional details and impacted files
@@           Coverage Diff           @@
##             main     #151   +/-   ##
=======================================
  Coverage   98.88%   98.89%           
=======================================
  Files          24       24           
  Lines        3766     3789   +23     
=======================================
+ Hits         3724     3747   +23     
  Misses         42       42           

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

@franckgaga
Copy link
Member Author

franckgaga commented Jan 22, 2025

As a fallback for the MHE, I first tried to implement the method of

Nicholas J. Higham, “Computing a Nearest Symmetric Positive Semidefinite Matrix,” Linear Algebra and its Applications, vol. 103, pp. 103–118, 1988.

described here on page 5-3. The eigenvalues of the resulting matrix are indeed all >= 0, but the cholesky call was still failing with IsPosDefException, probably because of the eigenvalues at 0.

I finally opt to keep the arrival covariance of the last iteration $k-1$ (or, in other words, to simply not update/correct it) and emit a warning. This is motivated by the fact that, for the MHE, the estimation horizon reduces the impact of the arrival covariance estimate in the total cost, so an approximate value is okay (compared to UKF/KF). As I mentioned in the issue, I used the MHE with a constant arrival covariance in the past and it was working well with a sufficiently long horizon. Moreover, the arrival covariance is always only a rough approximation when the MHE constraint are actives!

I'll add the Joseph form for the KalmanFilter and UnscentedKalmanFilter and default to that form in the MHE when I will have more time available.

Copy link
Member

@baggepinnen baggepinnen left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This approach is fine by me, the warning indicates what's being done. There are factorizations that can handle semi-definiteness, such as https://github.com/timholy/PositiveFactorizations.jl, but in practice I don't think it'll matter much as you say in #98 (comment)

if err isa PosDefException
@warn("Arrival covariance is not positive definite: keeping the old one")
else
rethrow(err)
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

yes true thanks! That's what I normally do but the copilot wrote this line for me. The more I use the copilot, the more I find mistake like that.

@baggepinnen
Copy link
Member

baggepinnen commented Jan 23, 2025

Keep in mind that you can get PosDepException for other reasons as well, this one has happened to me a few times and is a sign of other problems

julia> cholesky([NaN;;])
ERROR: PosDefException: matrix is not Hermitian; Factorization failed.
Stacktrace:
 [1] checkpositivedefinite(info::Int64)

we don't want to "hide" this issue by the try catch. any(isnan, P), or even better any(!isfinite, P), catches that

if err isa PosDefException
@warn("Arrival covariance is not positive definite: keeping the old one")
else
rethrow(err)
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

yes true thanks! That's what I normally do but the copilot wrote this line for me. The more I use the copilot, the more I find mistake like that.

@JuliaControl JuliaControl deleted a comment from baggepinnen Jan 23, 2025
@JuliaControl JuliaControl deleted a comment from baggepinnen Jan 23, 2025
@JuliaControl JuliaControl deleted a comment from baggepinnen Jan 23, 2025
the `Union` type was overly complex and `getinfo` should not be used in "normal" operation, only for troubleshoothing. The computational performance are not that important (there is already many allocations in this function)
The full output without ellipsis can still be shown in the REPL with the following logger:
`debuglogger = ConsoleLogger(stderr, Logging.Debug, show_limited=false)`
@franckgaga
Copy link
Member Author

franckgaga commented Jan 23, 2025

Keep in mind that you can get PosDepException for other reasons as well, this one has happened to me a few times and is a sign of other problems

julia> cholesky([NaN;;])
ERROR: PosDefException: matrix is not Hermitian; Factorization failed.
Stacktrace:
 [1] checkpositivedefinite(info::Int64)

we don't want to "hide" this issue by the try catch. any(isnan, P), or even better any(!isfinite, P), catches that

FYI, It does not generate this exception when the matrix is a Hermitian:

julia> cholesky(Hermitian([1 NaN;NaN NaN], :L))
Cholesky{Float64, Matrix{Float64}}
L factor:
2×2 LowerTriangular{Float64, Matrix{Float64}}:
   1.0      
 NaN    NaN

No exception at all is thrown as you can see. But the NaNs are propagated correctly so It will presumably crash later in the code.

@franckgaga franckgaga merged commit 35ee483 into main Jan 23, 2025
3 of 4 checks passed
@franckgaga franckgaga deleted the try_mhe branch January 23, 2025 23:28
@baggepinnen
Copy link
Member

I just had a related idea, when I write optimization over positive definite matrices, I usually optimize the triangular part of a cholesky factor instead of over the matrix itself, example here, see the functions triangular, invtriangular and

    n2 = length(params) ÷ 2
    Qchol = triangular(params[1:n2])
    Rchol = triangular(params[n2+1:2n2])
    Q = Qchol'Qchol
    R = Rchol'Rchol

This would be equivalent to using the "square-root form" of a Kalman filter, but for optimiztion-based approaches it's much simpler to set up. Not sure if it's relevant in the current case of non-pos-def arrival covariance in MHE, but if it is, it could perhaps give a boost to performance.

@franckgaga
Copy link
Member Author

Oh interesting!
Yep I'll give it a try! :)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants