make a contraction metric for the system defined by xdot in matlab

main.m
% Define system dynamics: xdot = f(x, u)
% Here, f is a function that takes x and u as input and returns xdot
f = @(x, u) [x(2); -x(1) + u];

% Define the state space
% Here, we assume x1 and x2 varies between -1 and 1, and u varies between -0.5 and 0.5
x1range = [-1, 1];
x2range = [-1, 1];
urange = [-0.5, 0.5];
N = 51;   % number of grid points in each dimension (increase N for higher accuracy)
[X1, X2, U] = ndgrid(linspace(x1range(1), x1range(2), N), ...
                     linspace(x2range(1), x2range(2), N), ...
                     linspace(urange(1), urange(2), N));

% Compute the contraction metric
% Here, we use the HJI toolbox for MATLAB (https://github.com/HJReachability/hj_matlab)
% to solve the HJI equation and obtain the contraction metric
schemeData.grid = {X1, X2, U};
schemeData.dynSys = f;
schemeData.accuracy = 'high';  % increase accuracy for higher precision
[K, ~] = HJIPDE_solve(@minWithDeriv, schemeData, 'zero');

% Plot the contraction metric over the state space
S = sqrt(K);  % the contraction metric is the square root of the solution to HJI equation
slice(X1, X2, U, S, [0], [0], [0]);  % plot the contraction metric at u=0
xlabel('x_1');
ylabel('x_2');
zlabel('S(x)');
1210 chars
29 lines