Merge pull request #213 from mldiego/master

CORA updates, ARCH-COMP 2024

code/nnv/engine/nncs/NonLinearODE.m

@@ -347,14 +347,12 @@ function set_output_mat(obj, output_mat)
                 for ik=1:Nz
                     try
                         Z1 = zonotope(Z{ik}); % ensure it's a zonotope
-                        Z1 = Z1.Z;
+                        % Z1 = Z1.Z;
                     catch
                         Z1 = zonotope(Z);
-                        Z1 = Z1.Z; % get c and V 
+                        % Z1 = Z1.Z; % get c and V 
                     end
-                    c = Z1(:,1); % center vector
-                    V = Z1(:, 2:size(Z1, 2)); % generators
-                    Zz = Zono(c, V);
+                    Zz = Zono(Z1.c, Z1.G);
                     S = [S Zz.toStar];
                 end
             end
@@ -402,14 +400,12 @@ function get_interval_sets(obj)
                         for iz=1:Nz
                             try
                                 Z2 = zonotope(Z1{iz}); % ensure it's a zonotope
-                                Z2 = Z2.Z;
+                                % Z2 = Z2.Z;
                             catch
                                 Z2 = zonotope(Z1); % ensure it's a zonotope
-                                Z2 = Z2.Z; % get c and V 
+                                % Z2 = Z2.Z; % get c and V 
                             end
-                            c = Z2(:,1); % center vector
-                            V = Z2(:, 2:size(Z2, 2)); % generators
-                            Zz = Zono(c, V);
+                            Zz = Zono(Z2.c, Z2.G);
                             Ss = [Ss Zz.toStar];
                         end
                     end
@@ -435,14 +431,12 @@ function get_point_sets(obj)
                         for iz=1:Nz
                             try
                                 Z2 = zonotope(Z1{iz}); % ensure it's a zonotope
-                                Z2 = Z2.Z;
+                                % Z2 = Z2.Z;
                             catch
                                 Z2 = zonotope(Z1); % ensure it's a zonotope
-                                Z2 = Z2.Z; % get c and V 
+                                % Z2 = Z2.Z; % get c and V 
                             end
-                            c = Z2(:,1); % center vector
-                            V = Z2(:, 2:size(Z2, 2)); % generators
-                            Zz = Zono(c, V);
+                            Zz = Zono(Z2.c, Z2.G);
                             Ss = [Ss Zz.toStar];
                         end
                     end

code/nnv/examples/Submission/ARCH-COMP2024/README.md

@@ -0,0 +1,47 @@
+# ARCH-COMP 2024
+ARCH-COMP AINNCS Category 2024 Model Files
+
+Event info: https://cps-vo.org/group/ARCH/FriendlyCompetition
+
+
+## Benchmarks 
+
+We plan to reuse benchmarks from 2023 but  feel free to propose new benchmarks or modifications to existing ones: https://github.com/verivital/ARCH-COMP2024/issues/1
+
+
+
+### [2023 Benchmarks](https://github.com/verivital/ARCH-COMP2023)
+
+- Adaptive Cruise Controller (ACC)
+
+- Airplane
+
+- Attitude Control 
+
+- Double Pendulum
+
+- Single Pendulum
+
+- QUAD
+
+- TORA with heterogeneous and sigmoid controller
+
+- TORA with ReLU controller (benchmark 9)
+
+- Unicycle (benchmark 10)
+
+- VCAS
+
+- 2D Spacecraft Docking
+
+
+### Competition History
+
+Prior year reports:
+- 2023: https://easychair.org/publications/paper/Vfq4b
+- 2022: https://easychair.org/publications/paper/C1J8
+- 2021: https://easychair.org/publications/paper/Jq4h
+- 2020: https://easychair.org/publications/paper/Jvwg
+- 2019: https://easychair.org/publications/paper/BFKs
+
+Repeatability archives: https://gitlab.com/goranf/ARCH-COMP/

code/nnv/examples/Submission/ARCH-COMP2024/benchmarks/ACC/Specifications.txt

@@ -0,0 +1,21 @@
+Parameters:
+
+a_lead = -2
+T_gap = 1.4
+D_default = 10
+
+Initial states:
+x_lead = x(1) = [90, 110]
+v_lead = x(2) = [32, 32.2]
+γ_lead = x(3) = 0
+x_ego  = x(4) = [10, 11]
+v_ego  = x(5) = [30, 30.2]
+γ_ego  = x(6) = 0
+
+t = 5 seconds
+control period = 0.1 s
+
+Safety Property:
+
+x_lead - x_ego >= D_safe, where D_safe = D_default + T_gap * v_ego;
+

code/nnv/examples/Submission/ARCH-COMP2024/benchmarks/ACC/dynamicsACC.m

@@ -0,0 +1,21 @@
+function [dx]=dynamicsACC(x,a_ego)
+
+mu=0.0001; % friction parameter
+
+% x1 = lead_car position
+% x2 = lead_car velocity
+% x3 = lead_car internal state
+
+% x4 = ego_car position
+% x5 = ego_car velocity
+% x6 = ego_car internal state
+
+% lead car dynamics
+a_lead  = -2; 
+dx(1,1) = x(2);
+dx(2,1) = x(3);
+dx(3,1) = -2 * x(3) + 2 * a_lead - mu*x(2)^2;
+% ego car dyanmics
+dx(4,1) = x(5); 
+dx(5,1) = x(6);
+dx(6,1) = -2 * x(6) + 2 * a_ego - mu*x(5)^2;

code/nnv/examples/Submission/ARCH-COMP2024/benchmarks/ACC/reach.m

@@ -0,0 +1,77 @@
+function rT = reach()
+
+%% Reachability analysis of ACC
+
+% Load components
+net = load_NN_from_mat('controller_5_20.mat');
+reachStep = 0.02;
+controlPeriod = 0.1;
+output_mat = [0 0 0 0 1 0;1 0 0 -1 0 0; 0 1 0 0 -1 0]; % feedback: relative distance, relative velocity and ego-car velocity
+plant = NonLinearODE(6,1,@dynamicsACC, reachStep, controlPeriod, output_mat);
+nncs = NonlinearNNCS(net,plant);
+
+%% Reachability analysis
+
+tF = 5; % seconds
+time = 0:controlPeriod:5;
+steps = length(time);
+input_ref = [30;1.4];
+
+% Initial set
+lb = [90; 32; 0; 10; 30; 0];
+ub = [110; 32.2; 0; 11; 30.2; 0];
+init_set = Star(lb,ub);
+
+% Input set
+lb = 0;
+ub = 0;
+input_set = Star(lb,ub);
+
+% Reachabilty analysis
+reachPRM.ref_input = input_ref;
+reachPRM.numSteps = 50;
+reachPRM.init_set = init_set;
+reachPRM.numCores = 1;
+reachPRM.reachMethod = 'approx-star';
+[R,rT] = nncs.reach(reachPRM);
+% disp("Time to compute ACC reach sets: " +string(rT));
+
+% Save results
+if is_codeocean
+    save('/results/logs/acc.mat', 'R','rT','-v7.3');
+else
+    save('acc.mat', 'R','rT','-v7.3');
+end
+
+%% Visualize results
+
+% Transform reach set into actual distance vs safe distance
+t_gap = 1.4;
+D_default = 10;
+outAll = [];
+safe_dis = [];
+% Transfrom intermediate reachsets from cora to NNV
+nncs.plant.get_interval_sets();
+% Get intermediate reach sets
+allsets = nncs.plant.intermediate_reachSet;
+for i=1:length(allsets)
+    outAll = [outAll allsets(i).affineMap(output_mat,[])];
+    safe_dis = [safe_dis allsets(i).affineMap([0 0 0 0 t_gap 0], D_default)];
+end
+times = reachStep:reachStep:tF;
+
+% Plot results
+f = figure;
+Star.plotRanges_2D(outAll,2,times,'b');
+hold on;
+Star.plotRanges_2D(safe_dis,1,times,'r');
+xlabel('Time (s)');
+ylabel('Distance (m)');
+% Save figure
+if is_codeocean
+    exportgraphics(f,'/results/logs/acc.pdf', 'ContentType', 'vector');
+else
+    exportgraphics(f,'acc.pdf','ContentType', 'vector');
+end
+
+end