Does anyone have experience in moving from Matlab’s LMI toolbox to the various JuliaOpt packages? I’m taking a non-linear dynamics course (AAE666 at Purdue) and we’ve just been introduced to linear matrix inequalities in the context of quadratic stability. As such, I’m still trying to grapple with the toolbox and LMIs in general and I’m not sure how to translate such a problem to a JuliaOPT package. Any resources or thoughts would be much appreciated.
Did you have a look at JuMP, especially the Section Semidefinite Constraints or the Example Example minimal Ellipse?
I have no experience so far in LMI’s but I think this would be a good starting point.
For example solving the Lyapunov Inequality for a continuous system could look like this:
using JuMP
using MosekTools
using LinearAlgebra
A = [-0.5 1; 0 -0.3]; n = size(A,1); Q = [1 0; 0 1]
model = Model(with_optimizer(Mosek.Optimizer))
@variable(model, P[i=1:n, j=1:n], PSD)
@SDconstraint(model, P ⪰ 0)
@SDconstraint(model, A'*P + P*A ⪯ -Q)
@SDconstraint(model, A'*P + P*A ⪰ -Q)
optimize!(model)
P = value.(P)
# Test result:
A'*P + P*A ≈ -Q
I wanted to try the example in the above post with a different optimizer. But, I got a StactOverflowError
. Any ideas what I am doing wrong?
using JuMP
using SCS
using LinearAlgebra
A = [-0.5 1; 0 -0.3]
n = size(A,1)
Q = [1 0; 0 1.]
model = Model(with_optimizer(SCS.Optimizer))
@variable(model, P[1:n, 1:n], PSD)
@SDconstraint(model, P ⪰ 0)
@SDconstraint(model, A'*P + P*A ≥ -Q)
@SDconstraint(model, A'*P + P*A ≤ -Q)
optimize!(model)
P = value.(P)
# Test result:
A'*P + P*A ≈ -Q
ERROR: LoadError: StackOverflowError:
Stacktrace:
[1] mutable_operate!(::typeof(MutableArithmetics.sub_mul), ::Matrix{AffExpr}, ::Matrix{VariableRef}, ::Matrix{Float64}) (repeats 79984 times)
@ MutableArithmetics ~/.julia/packages/MutableArithmetics/bPWR4/src/linear_algebra.jl:132
How about using Convex.jl?
julia> using Convex, SCS, LinearAlgebra
julia> A = [-0.5 1; 0 -0.3]
2×2 Matrix{Float64}:
-0.5 1.0
0.0 -0.3
julia> n = size(A,1)
2
julia> eigvals(A)
2-element Vector{Float64}:
-0.5
-0.3
julia> X = Semidefinite(n)
Variable
size: (2, 2)
sign: real
vexity: affine
id: 404…993
julia> constraint = A'*X + X*A ⪯ -Matrix{Float64}(I, n, n)
sdp constraint (affine)
└─ + (affine; real)
├─ 2×2 Matrix{Float64}
└─ - (affine; real)
└─ + (affine; real)
├─ …
└─ …
julia> problem = satisfy(constraint)
minimize
└─ 0
subject to
└─ sdp constraint (affine)
└─ + (affine; real)
├─ 2×2 Matrix{Float64}
└─ - (affine; real)
└─ …
status: `solve!` not called yet
julia> solve!(problem,() -> SCS.Optimizer(verbose=false))
julia> problem.status
OPTIMAL::TerminationStatusCode = 1
julia> X.value
2×2 Matrix{Float64}:
1.83577 2.88281
2.88281 11.9403
julia> eigvals(X.value)
2-element Vector{Float64}:
1.071158392495276
12.704883497869783
Someone asked this yesterday on StackOverflow: optimization - Julia JuMP StackOverflow Error when defining a constraint on linear matrix inequalities - Stack Overflow
It’s a bug: https://github.com/jump-dev/MutableArithmetics.jl/issues/86
As a work-around, split it into an expression first:
@expression(model, ex, A' * P + P * A)
@SDconstraint(model, ex <= -Q)
@zdenek_hurak Thank you for pointing out Convex.jl and providing me with the example.
@odow Thank you for the workaround. It works
If you update your packages, is this fixed now? If not, can you provide the output of ] st -m
.
The problem is fixed after the update. After the update, MutableArithmetics
package is upgraded from v.0.2.14 to v.0.2.15. The problem is fixed already, but in case you are interested, here is the output of ]st -m
(@dev-env) pkg> st -m
Status `~/.julia/environments/dev-env/Manifest.toml`
[c3fe647b] AbstractAlgebra v0.15.1
[621f4979] AbstractFFTs v1.0.1
[537997a7] AbstractPlotting v0.15.27
[1520ce14] AbstractTrees v0.3.4
[79e6a3ab] Adapt v3.3.0
[27a7e980] Animations v0.4.1
[c7e460c6] ArgParse v1.1.3
[ec485272] ArnoldiMethod v0.1.0
[4fba245c] ArrayInterface v3.1.7
[4c555306] ArrayLayouts v0.6.5
[bf4720bc] AssetRegistry v0.1.0
[13072b0f] AxisAlgorithms v1.0.0
[aae01518] BandedMatrices v0.16.8
[6e4b80f9] BenchmarkTools v0.5.0
[9e28174c] BinDeps v1.0.2
[b99e7846] BinaryProvider v0.5.10
[ad839575] Blink v0.12.5
[8e7c35d0] BlockArrays v0.15.2
[a74b3585] Blosc v0.7.0
[764a87c0] BoundaryValueDiffEq v2.7.1
[e1450e63] BufferedStreams v1.0.0
[fa961155] CEnum v0.4.1
[15fcbb24] CGAL v0.4.0
[5d742f6a] CSVFiles v1.0.0
[49dc2e85] Calculus v0.5.1
[324d7699] CategoricalArrays v0.9.5
[c45f814d] Causal v0.3.0 `~/.julia/dev/Causal`
[d360d2e6] ChainRulesCore v0.9.37
[fe2786b3] ChaoticCommunications v0.1.0 `~/.julia/dev/ChaoticCommunications`
[aaaa29a8] Clustering v0.14.2
[da1fd8a2] CodeTracking v1.0.5
[523fee87] CodecBzip2 v0.7.2
[944b1d66] CodecZlib v0.6.0
[a2cac450] ColorBrewer v0.4.0
[35d6a980] ColorSchemes v3.11.0
[3da002f7] ColorTypes v0.10.12
[c3611d14] ColorVectorSpace v0.8.7
[5ae59095] Colors v0.12.7
[861a8166] Combinatorics v1.0.2
[38540f10] CommonSolve v0.2.0
[bbf7d656] CommonSubexpressions v0.3.0
[34da2185] Compat v3.27.0
[a81c6b42] Compose v0.9.2
[8f4d0f93] Conda v1.5.1
[187b0558] ConstructionBase v1.1.0
[d38c429a] Contour v0.5.7
[f65535da] Convex v0.14.6
[a8cc5b0e] Crayons v4.0.4
[667455a9] Cubature v1.5.1
[1f15a43c] CxxWrap v0.11.2
[717857b8] DSP v0.6.10
[9a962f9c] DataAPI v1.6.0
[a93c6f00] DataFrames v0.22.7
[864edb3b] DataStructures v0.18.9
[e2d170a0] DataValueInterfaces v1.0.0
[e7dc6d0d] DataValues v0.4.13
[bcd4f6db] DelayDiffEq v5.29.1
[2b5f629d] DiffEqBase v6.58.0
[459566f4] DiffEqCallbacks v2.16.1
[5a0ffddc] DiffEqFinancial v2.4.0
[c894b116] DiffEqJump v6.14.0
[77a26b50] DiffEqNoiseProcess v5.7.1
[055956cb] DiffEqPhysics v3.9.0
[163ba53b] DiffResults v1.0.3
[b552c78f] DiffRules v1.0.2
[0c46a032] DifferentialEquations v6.16.0
[667095cb] DigitalCommunications v0.1.0 `~/.julia/dev/DigitalCommunications`
[c619ae07] DimensionalPlotRecipes v1.2.0
[b4f34e82] Distances v0.10.2
[31c24e10] Distributions v0.24.15
[ffbed154] DocStringExtensions v0.8.4
[e30172f5] Documenter v0.26.3
[35a29f4d] DocumenterTools v0.1.9
[497a8b3b] DoubleFloats v1.1.19
[78f8389a] DynamicalNetworks v0.1.0 `~/.julia/dev/DynamicalNetworks`
[da5c29d0] EllipsisNotation v1.1.0
[7876af07] Example v0.5.3
[89b67f3b] ExcelFiles v1.0.0
[c04bee98] ExcelReaders v0.11.0
[d4d017d3] ExponentialUtilities v1.8.3
[e2ba6199] ExprTools v0.1.3
[8f5d6c58] EzXML v1.1.0
[c87230d0] FFMPEG v0.4.0
[7a1cc6ca] FFTW v1.3.2
[9aa1b823] FastClosures v0.3.2
[5789e2e9] FileIO v1.6.5
[1a297f60] FillArrays v0.11.7
[6a86dc24] FiniteDiff v2.8.0
[53c48c17] FixedPointNumbers v0.8.4
[59287772] Formatting v0.4.2
[f6369f11] ForwardDiff v0.10.18
[b405c40b] FractalTools v0.1.0 `~/.julia/dev/FractalTools`
[b38be410] FreeType v3.0.1
[663a7486] FreeTypeAbstraction v0.8.4
[069b7b12] FunctionWrappers v1.1.2
[de31a74c] FunctionalCollections v0.5.0
[fb4132e2] FuzzyCompletions v0.4.1
[f7f18e0c] GLFW v3.4.1
[e9467ef8] GLMakie v0.1.30
[60bf3e95] GLPK v0.14.8
[28b8d3ca] GR v0.57.3
[0337cf30] GRUtils v0.6.1
[14197337] GenericLinearAlgebra v0.2.5
[fd0ad045] GeometricalPredicates v0.4.0
[5c1252a2] GeometryBasics v0.3.12
[a2cc645c] GraphPlot v0.4.4
[a2bd30eb] Graphics v1.1.0
[3955a311] GridLayoutBase v0.5.4
[42e2da0e] Grisu v1.0.0
[f67ccb44] HDF5 v0.15.4
[cd3eb016] HTTP v0.8.19
[9fb69e20] Hiccup v0.2.2
[7073ff75] IJulia v1.23.2
[b5f81e59] IOCapture v0.1.1
[615f187c] IfElse v0.1.0
[a09fc81d] ImageCore v0.8.21
[82e4d734] ImageIO v0.5.3
[9b13fd28] IndirectArrays v0.5.1
[d25df0c9] Inflate v0.1.2
[83e8ac13] IniFile v0.5.0
[a98d9a8b] Interpolations v0.12.10
[8197267c] IntervalSets v0.5.3
[d8418881] Intervals v1.5.0
[41ab1584] InvertedIndices v1.0.0
[b6b21f68] Ipopt v0.6.5
[f1662d9f] Isoband v0.1.1
[c8e1da08] IterTools v1.3.0
[1c8ee90f] IterableTables v1.0.0
[42fd0dbc] IterativeSolvers v0.9.0
[82899510] IteratorInterfaceExtensions v1.0.0
[033835bb] JLD2 v0.2.4
[692b3bcd] JLLWrappers v1.2.0
[97c1335a] JSExpr v0.5.2
[682c06a0] JSON v0.21.1
[7d188eb4] JSONSchema v0.3.3
[4076af6c] JuMP v0.21.6
[aa1ae85d] JuliaInterpreter v0.8.11
[5ab0869b] KernelDensity v0.6.2
[b964fa9f] LaTeXStrings v1.2.1
[2ee39098] LabelledArrays v1.6.0
[23fbe1c1] Latexify v0.14.12
[50d2b5c4] Lazy v0.15.1
[093fc24a] LightGraphs v1.3.5
[d3d80556] LineSearches v7.1.1
[e6f89c97] LoggingExtras v0.4.6
[6f1432cf] LoweredCodeUtils v1.2.9
[23992714] MAT v0.10.1
[1914dd2f] MacroTools v0.5.6
[ee78f7c6] Makie v0.12.0
[dbb5928d] MappedArrays v0.3.0
[7eb4fadd] Match v1.1.0
[b8f27783] MathOptInterface v0.9.20
[fdba3010] MathProgBase v0.7.8
[739be429] MbedTLS v1.0.3
[442fdcdd] Measures v0.3.1
[7269a6da] MeshIO v0.4.6
[e1d29d7a] Missings v0.4.5
[78c3b35d] Mocking v0.7.1
[961ee093] ModelingToolkit v5.14.2
[66fc600b] ModernGL v1.1.2
[e94cdb99] MosaicViews v0.2.4
[6405355b] Mosek v1.1.3
[1ec41992] MosekTools v0.9.4
[99f44e22] MsgPack v1.1.0
[46d2c3a1] MuladdMacro v0.2.2
[f9640e96] MultiScaleArrays v1.8.1
[ffc61752] Mustache v1.0.10
[d8a4904e] MutableArithmetics v0.2.15
[a975b10e] Mux v0.7.5
[eb2d6272] MyPackage v0.1.3 `~/.julia/dev/MyPackage`
[d41bc354] NLSolversBase v7.8.0
[2774e3e8] NLsolve v4.5.1
[77ba4419] NaNMath v0.3.5
[b8a86587] NearestNeighbors v0.4.8
[f09324ee] Netpbm v1.0.0
[8913a72c] NonlinearSolve v0.3.8
[4d1e1d77] Nullables v1.0.0
[47be7bcc] ORCA v0.5.0
[510215fc] Observables v0.3.3
[6fe1bfb0] OffsetArrays v1.6.2
[429524aa] Optim v1.3.0
[bac558e1] OrderedCollections v1.4.0
[1dea7af3] OrdinaryDiffEq v5.52.3
[90014a1f] PDMats v0.11.0
[f57f5aa1] PNGFiles v0.3.6
[9b87118b] PackageCompiler v1.2.5
[19eb6ba3] Packing v0.4.1
[5432bcbf] PaddedViews v0.5.8
[65888b18] ParameterizedFunctions v5.10.0
[d96e819e] Parameters v0.12.2
[69de0a69] Parsers v1.1.0
[fa939f87] Pidfile v1.2.0
[14b8a8f1] PkgTemplates v0.7.16
[eebad327] PkgVersion v0.1.1
[a83a3aaf] PkgVersionHelper v0.1.2
[ccf2f8ad] PlotThemes v2.0.1
[995b91a9] PlotUtils v1.0.10
[a03496cd] PlotlyBase v0.4.3
[f0f68f2c] PlotlyJS v0.14.1
[91a5bcdd] Plots v1.11.2
[c3e4b0f8] Pluto v0.14.1
[e409e4f3] PoissonRandom v0.4.0
[647866c9] PolygonOps v0.1.1
[67491407] Polyhedra v0.6.13
[f27b6e38] Polynomials v1.2.1
[2dfb63ee] PooledArrays v1.2.1
[85a6dd25] PositiveFactorizations v0.2.4
[08abe8d2] PrettyTables v0.11.1
[92933f4c] ProgressMeter v1.5.0
[438e738f] PyCall v1.92.3
[d330b81b] PyPlot v2.9.0
[1fd47b50] QuadGK v2.4.1
[be4d8f0f] Quadmath v0.5.5
[1a8c2f83] Query v1.0.0
[2aef5ad7] QueryOperators v0.9.3
[74087812] Random123 v1.3.1
[fb686558] RandomExtensions v0.4.3
[e6cf234a] RandomNumbers v1.4.0
[c84ed2f1] Ratios v0.4.0
[3cdcf5f2] RecipesBase v1.1.1
[01d81517] RecipesPipeline v0.3.2
[731186ca] RecursiveArrayTools v2.11.1
[f2c3362d] RecursiveFactorization v0.1.0
[189a3867] Reexport v1.0.0
[ae029012] Requires v1.1.3
[ae5879a3] ResettableStacks v1.1.0
[295af30f] Revise v3.1.14
[79098fc4] Rmath v0.7.0
[7e49a35a] RuntimeGeneratedFunctions v0.5.2
[c946c3f1] SCS v0.7.1
[1bc83da4] SafeTestsets v0.0.1
[322a6be2] Sass v0.1.0
[0bca4576] SciMLBase v1.11.0
[6c6a2e73] Scratch v1.0.3
[efcf1570] Setfield v0.7.0
[65257c39] ShaderAbstractions v0.2.5
[992d4aef] Showoff v0.3.2
[73760f76] SignedDistanceFields v0.4.0
[699a6c99] SimpleTraits v0.9.3
[b85f4697] SoftGlobalScope v1.1.0
[a2af1166] SortingAlgorithms v0.3.1
[47a9eef4] SparseDiffTools v1.13.0
[276daf66] SpecialFunctions v1.3.0
[aedffcd0] Static v0.2.4
[90137ffa] StaticArrays v0.12.5
[2913bbd2] StatsBase v0.33.5
[4c63d2b9] StatsFuns v0.9.7
[9672c7b4] SteadyStateDiffEq v1.6.2
[789caeaf] StochasticDiffEq v6.33.1
[69024149] StringEncodings v0.3.3
[09ab397b] StructArrays v0.5.1
[856f2bd8] StructTypes v1.5.2
[c3572dad] Sundials v4.4.1
[d1185830] SymbolicUtils v0.10.3
[0c5d862f] Symbolics v0.1.18
[d1efa939] TableIOInterface v0.1.6
[5e66a065] TableShowUtils v0.2.5
[3783bdb8] TableTraits v1.0.1
[382cd787] TableTraitsUtils v1.0.1
[bd369af6] Tables v1.4.1
[c5d3f3f7] TetGen v1.1.0
[e0df1984] TextParse v1.0.1
[b718987f] TextWrap v1.0.1
[731e570b] TiffImages v0.3.0
[f269a46b] TimeZones v1.5.3
[a759f4b9] TimerOutputs v0.5.8
[3bb67fe8] TranscodingStreams v0.9.5
[a2a6695c] TreeViews v0.3.0
[9f332540] Triangulation v0.1.0 `~/.julia/dev/Triangulation`
[30578b45] URIParser v0.4.1
[3a884ed6] UnPack v1.0.2
[1cfade01] UnicodeFun v0.4.1
[1986cc42] Unitful v1.7.0
[81def892] VersionParsing v1.2.0
[19fa3120] VertexSafeGraphs v0.1.2
[72f80fcb] VoronoiDelaunay v0.4.0
[ea10d353] WeakRefStrings v0.6.2
[0f1e0344] WebIO v0.8.15
[104b5d7c] WebSockets v1.5.9
[cc8bc4a8] Widgets v0.6.2
[1b915085] WinReg v0.3.1
[efce3f68] WoodburyMatrices v0.5.3
[fdbf4ff8] XLSX v0.7.6
[c2297ded] ZMQ v1.2.1
[a5390f91] ZipFile v0.9.3
[700de1a5] ZygoteRules v0.2.1
[ae81ac8f] ASL_jll v0.1.2+0
[0b7ba130] Blosc_jll v1.14.3+1
[6e34b625] Bzip2_jll v1.0.6+5
[8fcd9439] CGAL_jll v5.2.1+0
[83423d85] Cairo_jll v1.16.0+6
[7bc98958] Cubature_jll v1.0.4+0
[5ae413db] EarCut_jll v2.1.5+1
[2e619515] Expat_jll v2.2.7+6
[b22a6f82] FFMPEG_jll v4.3.1+4
[f5851436] FFTW_jll v3.3.9+7
[a3f928ae] Fontconfig_jll v2.13.1+14
[d7e528f0] FreeType2_jll v2.10.1+5
[559328eb] FriBidi_jll v1.0.5+6
[0656b61e] GLFW_jll v3.3.3+0
[e8aa6df9] GLPK_jll v5.0.0+0
[d2c73de3] GR_jll v0.57.2+0
[78b55507] Gettext_jll v0.20.1+7
[7746bdde] Glib_jll v2.59.0+4
[0234f1f7] HDF5_jll v1.12.0+1
[1d5cc7b8] IntelOpenMP_jll v2018.0.3+2
[9cc047cb] Ipopt_jll v3.13.4+1
[aacddb02] JpegTurbo_jll v2.0.1+3
[c1c5ebd0] LAME_jll v3.100.0+3
[dd4b983a] LZO_jll v2.10.0+3
[dd192d2f] LibVPX_jll v1.9.0+1
[e9f186c6] Libffi_jll v3.2.1+4
[d4300ac3] Libgcrypt_jll v1.8.5+4
[7e76a0d4] Libglvnd_jll v1.3.0+3
[7add5ba3] Libgpg_error_jll v1.36.0+3
[94ce4f54] Libiconv_jll v1.16.0+7
[4b2f31a3] Libmount_jll v2.34.0+3
[89763e89] Libtiff_jll v4.1.0+2
[38a345b3] Libuuid_jll v2.34.0+7
[5ced341a] Lz4_jll v1.9.2+2
[d00139f3] METIS_jll v5.1.0+5
[856f044c] MKL_jll v2021.1.1+1
[d7ed1dd3] MUMPS_seq_jll v5.2.1+4
[e7412a2a] Ogg_jll v1.3.4+2
[656ef2d0] OpenBLAS32_jll v0.3.12+1
[458c3c95] OpenSSL_jll v1.1.1+6
[efe28fd5] OpenSpecFun_jll v0.5.3+4
[91d4177d] Opus_jll v1.3.1+3
[2f80f16e] PCRE_jll v8.42.0+4
[30392449] Pixman_jll v0.40.0+0
[ea2cea3b] Qt5Base_jll v5.15.2+0
[f50d1b31] Rmath_jll v0.3.0+0
[af6e375f] SCS_GPU_jll v2.1.2+0
[f4f2fc5b] SCS_jll v2.1.2+1
[fb77eaff] Sundials_jll v5.2.0+1
[b47fdcd6] TetGen_jll v1.5.1+1
[a2964d1f] Wayland_jll v1.17.0+4
[2381bf8a] Wayland_protocols_jll v1.18.0+4
[02c8fc9c] XML2_jll v2.9.11+0
[aed1982a] XSLT_jll v1.1.33+4
[4f6342f7] Xorg_libX11_jll v1.6.9+4
[0c0b7dd1] Xorg_libXau_jll v1.0.9+4
[935fb764] Xorg_libXcursor_jll v1.2.0+4
[a3789734] Xorg_libXdmcp_jll v1.1.3+4
[1082639a] Xorg_libXext_jll v1.3.4+4
[d091e8ba] Xorg_libXfixes_jll v5.0.3+4
[a51aa0fd] Xorg_libXi_jll v1.7.10+4
[d1454406] Xorg_libXinerama_jll v1.1.4+4
[ec84b674] Xorg_libXrandr_jll v1.5.2+4
[ea2f1a96] Xorg_libXrender_jll v0.9.10+4
[14d82f49] Xorg_libpthread_stubs_jll v0.1.0+3
[c7cfdc94] Xorg_libxcb_jll v1.13.0+3
[cc61e674] Xorg_libxkbfile_jll v1.1.0+4
[12413925] Xorg_xcb_util_image_jll v0.4.0+1
[2def613f] Xorg_xcb_util_jll v0.4.0+1
[975044d2] Xorg_xcb_util_keysyms_jll v0.4.0+1
[0d47668e] Xorg_xcb_util_renderutil_jll v0.3.9+1
[c22f9ab0] Xorg_xcb_util_wm_jll v0.4.1+1
[35661453] Xorg_xkbcomp_jll v1.4.2+4
[33bec58e] Xorg_xkeyboard_config_jll v2.27.0+4
[c5fb5394] Xorg_xtrans_jll v1.4.0+3
[8f1865be] ZeroMQ_jll v4.3.2+6
[3161d3a3] Zstd_jll v1.4.8+0
[28df3c45] boost_jll v1.71.0+3
[9a68df92] isoband_jll v0.2.2+0
[0ac62f75] libass_jll v0.14.0+4
[e9ad47b2] libcgal_julia_jll v0.16.106+0
[3eaa8342] libcxxwrap_julia_jll v0.8.6+0
[f638f0a6] libfdk_aac_jll v0.1.6+4
[b53b4c65] libpng_jll v1.6.37+6
[a9144af2] libsodium_jll v1.0.20+0
[f27f6e37] libvorbis_jll v1.3.6+6
[1270edf5] x264_jll v2020.7.14+2
[dfaa095f] x265_jll v3.0.0+3
[d8fb68d0] xkbcommon_jll v0.9.1+5
[0dad84c5] ArgTools
[56f22d72] Artifacts
[2a0f44e3] Base64
[ade2ca70] Dates
[8bb1440f] DelimitedFiles
[8ba89e20] Distributed
[f43a241f] Downloads
[7b1f6079] FileWatching
[9fa8497b] Future
[b77e0a4c] InteractiveUtils
[4af54fe1] LazyArtifacts
[b27032c2] LibCURL
[76f85450] LibGit2
[8f399da3] Libdl
[37e2e46d] LinearAlgebra
[56ddb016] Logging
[d6f4376e] Markdown
[a63ad114] Mmap
[ca575930] NetworkOptions
[44cfe95a] Pkg
[de0858da] Printf
[3fa0cd96] REPL
[9a3f8284] Random
[ea8e919c] SHA
[9e88b42a] Serialization
[1a1011a3] SharedArrays
[6462fe0b] Sockets
[2f01184e] SparseArrays
[10745b16] Statistics
[4607b0f0] SuiteSparse
[fa267f1f] TOML
[a4e569a6] Tar
[8dfed614] Test
[cf7118a7] UUIDs
[4ec0a83e] Unicode
[e66e0078] CompilerSupportLibraries_jll
[781609d7] GMP_jll
[deac9b47] LibCURL_jll
[29816b5a] LibSSH2_jll
[3a97d323] MPFR_jll
[c8ffd9c3] MbedTLS_jll
[14a3606d] MozillaCACerts_jll
[4536629a] OpenBLAS_jll
[bea87d4a] SuiteSparse_jll
[83775a58] Zlib_jll
[8e850ede] nghttp2_jll
[3f19e933] p7zip_jll
(@dev-env) pkg>
Great! I had fixed this a while ago, it just needed a new release.
Dear Friends
I need to solve an LMI as follows
As I am a beginner to Julia I do not know how to solve it
(build it in Julia context) and which commands are appropriate.
here is the problem
I need to find P
I appreciate your help with this.
Well, I think that the reason why you were directed by @rafael.guerra from your original post to this post was not exactly to replicate your question here but rather to have a look at what had been discussed here. The discussion here turns out to be relevant for your need. Have you at least tried and struggled? Where exactly have you struggled?
But let me restate here for convenience: you have two options: either use JuMP or Convex to formulate your LMI feasibility problem. For this particular problem both do the job. The documentation for both contains a sufficient description how to define and solve semidefinite problems, in this case just LMI feasibility problems. Examples of usage for both packages are also above in this thread. It should not be difficult to massage those code snippets so that they address your problem.
I am eager to help you once you encounter problems.
At minimum, please, paste here a code for generating those matrices A
, B
, Q
and R
(inside triple quotes).
Two notes. First, as this LMI corresponds to the standard discrete-time LQ-optimal regulation (aka LQR) problem - it is just the discrete-time algebraic Riccati equation (ARE) rewritten into an LMI format -, it is quite vital to specify that the matrices Q
and R
are at least positive semidefinite (R
even positive definite), right? Without it the problem does not have much use, right? This is an important part of the problem statement.
Second, since the modeling (and solving) packages that I have listed above can combine individual LMIs, it may be more convenient to split that single LMI into two: the first will contain the left upper 2x2 block, the second will be just the constraint that P
is positive (semi?)definite. The softwares allow imposing constraints directly on the variables (see the examples above in this thread), hence you will only have to encode a single 2x2 matrix LMI constraint.
Well, although it is a different problem from what you originally posted, the approach to its solution remains the same - just assemble a large matrix from individual blocks. Below I show the code for the original problem:
using Convex
using SCS
n = 3
m = 2
A = rand(n,n)
B = rand(n,m)
Q = diagm([1.0, 10.0, 100.0])
R = diagm([1.0, 10.0])
P = Semidefinite(n)
constraint = ([R+B'*P*B B'*P*A; A'*P*B A'*P*A+Q-P] ⪰ 0)
problem = satisfy(constraint)
solve!(problem,SCS.Optimizer)
P.value
Note: the inequalities in the code are nonstrict while your problem asked for strict inequalities. This can be accomplished by including some nonzero term on the right hand side. Typically an identity matrix is used.
What worries me, however, is that in your new problem the matrix is not symmetric. And symmetry is typically assumed in LMIs. Is your matrix correct?
Thank you so much for your reply
Yes I have checked it but it seems to be strange obtaining
If this be the only result (I mean the obtained none-symmetric matrix)
Then this means that it can not be solved by the current
LMI methods? As they are for symmetric matrices?
the negative definiteness of Q will not affect the solution or the solver?
I am afraid I do not quite get what you mean. So, you started with a matrix defining a positive (or negative, it does not matter, just multiply by -1) definite form and you factorized it and now you want to enforce the (semi)definiteness constraint on the factor? Why? What is actually your original problem?
Aha, it was a mistake I made while I was posting here the first problem was not my problem I am really sorry for this mistake.
In my problem, all the matrices except Q are semi-positive definite (positive definite) and I need to find P. (being positive definite or semi positive definite does not hurt the result which I am going to make and it depends on the rest part that which one is better I can take them positive definite)
If we multiply Q by a -1 how will problem change?
So this asymmetric Matrix inequality can be solved by the LMI solvers offered by Julia?
My original problem is proving (negative definiteness of \Lambda)
So if P is positive definite Must I write :
P = Definite(n)
In the constraints
Must I put
constraint = (\Lambda<0)?
. I also thought of developing a numerical method if it is possible and if this matrix is not an LMI
but it is not my main problem during my study.
Thank you so much
I am afraid I am not able to follow you. Now you start talking about some \Lambda
variable, which appeared in none of your previously stated problems…
A few general remarks:
-
I can hardly think of a need to analyze positive (semi) definiteness of a matrix while not assuming a symmetry. At least in control theory (which is where your problems are obviously coming from), these matrices define quadratic matrix forms and even if you somehow start with a nonsymmetrix matrix
P
, you can symmetrize it through(P+P')/2
with no impact on the value of the quadratic form. Thats why you will barely encounter a discussion of nonsymmetric matrices in the LMI/SDP literature. -
Whether you can use the function
Positive()
, you can easily check by yourself. You do not need discussion forum for that, do you? Just read the documentation: for Convex.jl package at Basic Types · Convex.jl. For JuMP. jl package at Variables · JuMP. -
In order to increase your chances of being helped by somebody here, please read and follow the instructions at Please read: make it easier to help you. In particular (but not only) the item 2 on inserting the code into your posts.
My dear friend thank you so much for your notes.
As Julia is ratherly young programming language (but at the same time really effective) all sources
does not have examples that can cover almost every corner and this made me ask.
sorry if they are too primitive maybe.
Again I really appreciate your help.