Safety in RL #101 added condensed matrices

Project.toml

@@ -15,6 +15,7 @@ Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 FileIO = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549"
 Flux = "587475ba-b771-5e3f-ad9e-33799f191a9c"
+GLPK = "60bf3e95-4087-53dc-ae20-288a0d20c6a6"
 Glob = "c27321d9-0574-5035-807b-f59d2c89b15c"
 GraphPlot = "a2cc645c-3eea-5389-862e-a155d0052231"
 Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
@@ -32,6 +33,7 @@ NLSolvers = "337daf1e-9722-11e9-073e-8b9effe078ba"
 NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec"
 PlotlyBase = "a03496cd-edff-5a9b-9e67-9cda94a718b5"
 PlotlyJS = "f0f68f2c-4968-5e81-91da-67840de0976a"
+Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 Polyhedra = "67491407-f73d-577b-9b50-8179a7c68029"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 ReinforcementLearning = "158674fc-8238-5cab-b5ba-03dfc80d1318"

examples/scripts/Polyhedra_minimal_robust_set.jl

@@ -0,0 +1,63 @@
+# https://juliapolyhedra.github.io/Polyhedra.jl/stable/generated/Minimal%20Robust%20Positively%20Invariant%20Set/
+
+
+A = [1 1; 0 1] - [1; 1] * [1.17 1.03]
+
+using Polyhedra
+Wv = vrep([[x, y] for x in [-1.0, 1.0] for y in [-1.0, 1.0]])
+
+using GLPK
+using JuMP
+lib = DefaultLibrary{Float64}(GLPK.Optimizer)
+W = polyhedron(Wv, lib)
+
+
+function Fs(s::Integer, verbose=1)
+    @assert s ≥ 1
+    F = W
+    A_W = W
+    for i in 1:(s-1)
+        A_W = A * A_W
+        F += A_W
+        if verbose ≥ 1
+            println("Number of points after adding A^$i * W: ", npoints(F))
+        end
+        removevredundancy!(F)
+        if verbose ≥ 1
+            println("Number of points after removing redundant ones: ", npoints(F))
+        end
+    end
+    return F
+end
+
+F_1 = Fs(1)
+F_2 = Fs(2)
+F_3 = Fs(3)
+
+using Plots
+plot()
+for i in 10:-1:1
+    plot!(Fs(i, 0))
+end
+plot!()
+display(p)
+
+using PlotlyJS
+layout = Layout(
+            plot_bgcolor="#f1f3f7",
+            title = "DoubleIntegrator",
+            xaxis_title = "x1",
+            yaxis_title = "x2",
+            width = 1200,
+            height = 850,
+            margin=attr(l=100, r=200, b=80, t=100, pad=10)
+        )
+
+traces = []
+
+for i in 10:-1:1
+    push!(traces, Fs(i, 0))
+end
+
+p = PlotlyBase.Plot(traces, layout)#, config = PlotConfig(scrollZoom=true))
+display(p)

examples/scripts/Safeguard_test.jl

@@ -1,4 +1,5 @@
 using ElectricGrid
+using LinearAlgebra
 
 CM = [ 0. 0. 1.
         0. 0. 2.
@@ -81,7 +82,7 @@ omega_x = [
 
 A_sqare = A_P10^2
 
-N = 4  # maximum number of iteration
+N = 3  # maximum number of iteration
 
 m = 1 # number of inputs to the system,
 n = size(B_P10)[1]
@@ -91,33 +92,35 @@ n = size(B_P10)[1]
 global A_N = 1* Matrix(I, size(A_P10)[1], size(A_P10)[1])
 global B_N = zeros(n, m*N)
 
+using SparseArrays
+
+W_x_cal = sparse(W_x)
+W_x_sp = sparse(W_x)
+W_u_cal = sparse(W_u)
+
+W_x_cal = blockdiag(W_x_cal, W_x_cal)
+W_u_cal = cat(W_u,  reverse(W_u); dims=(1,2))
+
 for ii  = 1:N
 
     global A_N = [A_N; A_P10^ii]
 
+    global W_x_cal = blockdiag(W_x_cal, W_x_sp)
+    #Wglobal W_u_cal = blockdiag(W_u_cal, W_u_cal)
+
 
 
     B_N_newline = A_P10^(ii-1)*B_P10
-    println("First:")
-    println(B_N_newline)
-    println("")
+
 
     for jj in 2:ii
-        println("First:")
-        println(A_P10^(ii-jj)*B_P10)
-        println("")
+
         B_N_newline = hcat(B_N_newline, A_P10^(ii-jj)*B_P10)
 
-        println("")
-        println(jj)
-        println("")
-        println(B_N_newline)
-        println("")
+
     end
     B_N_newline = hcat(B_N_newline, zeros(n, m*(N-ii)))
-    println("Last:")
-    println(B_N_newline)
-    println("")
+
 
 
 
@@ -132,3 +135,8 @@ println(A_N)
 println("")
 println(B_N)
 println("")
+println("")
+println(W_x_cal)
+println("")
+println(W_u_cal)
+println("")