multiparametric

A = randn(15,3);
b = rand(15,1);
E = randn(15,2);

z = sdpvar(3,1);
x = [0.1;0.2];

F = [A*z <= b+E*x];
obj = (z-1)'*(z-1);

sol = optimize(F,obj);
value(z)

x = sdpvar(2,1);
F = [A*z <= b+E*x, -1 <= x <= 1];
sol = solvemp(F,obj,[],x);

xx = [0.1;0.2];
[i,j] = isinside(sol{1}.Pn,xx)
sol{1}.Fi{j}*xx + sol{1}.Gi{j}

[sol,diagnostics,aux,Valuefunction,OptimalSolution] = solvemp(F,obj,[],x);

assign(x,[0.1;0.2]);
value(OptimalSolution)


plot(Valuefunction);
figure
plot(OptimalSolution(1));



A = [2 -1;1 0];
B = [1;0];
C = [0.5 0.5];

N = 5;
x = sdpvar(2,N+1,'full');
u = sdpvar(1,N);

Model = [];
Objective = 0;
for i = 1:N
    Objective = Objective + norm(C*x(:,i+1),1) + norm(u(i),1);
    Model = [Model, x(:,i+1) == A*x(:,i) + B*u(i)];
    Model = [Model, -1 <= u(i) <= 1, -5 <= x(:,i+1) <= 5];  
end
Model = [Model, -5 <= x(:,1) <= 5];
[sol,diagnostics,aux,Valuefunction,OptimalSolution] = solvemp(Model,Objective,[],x(:,1),u(1));
plot(Valuefunction);
plot(OptimalSolution);

N = 3;
x = sdpvar(2,N+1,'full');
u = sdpvar(1,N);

Model = [];
Objective = 0;
for i = 1:N
    Objective = Objective + norm(C*x(:,i+1),1) + norm(u(i),1);
    Model = [Model, x(:,i+1) == A*x(:,i) + B*u(i)];
    Model = [Model, -1 <= u(i) <= 1, -5 <= x(:,i+1) <= 5];
    Model = [Model, ismember(u(i),[-1:1/3:1])];
end
Model = [Model, -5 <= x(:,1) <= 5];
[sol,diagnostics,aux,Valuefunction,OptimalSolution] = solvemp(Model,Objective,[],x(:,1),u(1));
plot(Valuefunction);
plot(OptimalSolution);