Model Predictive control of LTI system?

[adnan jaffar]

Dear Dr. Johan,
I am looking for Model predictive control Yalmip code if it is already tested and working. Could you please refer me any Model predictive control which could solve the control problem of attached plant model.
In attached plant, all state space matrices and disturbance input are known.
Your help would be highly appreciated

[Johan Löfberg]

https://yalmip.github.io/example/standardmpc/

[adnan jaffar]

Dear Dr. Johan,
Thank you so much for your quick reply.
I have followed the code and i got following warning once i run the following code. Can you please comment on following where i am doing something wrong. Thanks


                                [Code]
A=eye(7)+0.05*[-0.350312629 0.001652708 -22.38023628 0 -0.007 0 0;119.9592691 -1.236693495 16746.7579 0 0 0 -0.001872252;1 -0.010309278 0 0 0 0 0;0 0 0 0 1 0 0;0.000212974    0 0 -25.98376792 -0.055356012 -7.48E-06    0;0    0 0 0 0 -2 0;0 0 0 0 0 0 -10];
B=0.05*[0 0;0 0;0 0;0 0;0 0;2 0;0 10];
C=[1 0 0 0 0 0 0;0 1 0 0 0 0 0;0 0 1 0 0 0 0;0 0 0 1 0 0 0;0 0 0 0 1 0 0;0 0 0 0 0 1 0;0 0 0 0 0 0 1;0 3776 0 0 0 66.08    0];
E=0.05*[0.007;0;0;0;0.000212974;0;0];
N = 10;
U = sdpvar(N,2);
W = sdpvar(N,1);
x = sdpvar(7,1);

Y = [];
xk = x;
for k = 1:N
 xk = A*xk + B*[U(k);U(k+2)]+E*W(k);
 Y = [Y;C*xk];
end
F = [70<=Y(1)<= 123,-0.72<=Y(2)<=-1.26,-0.5<=Y(3)<=0.5,0<=Y(4)<=47403,-15000<=Y(5)<=15000,-0.017<=Y(6)<= 1.57,-0.139<=Y(7)<=0.139,-0.0865<=U(1)<=1.6395,15000<=U(2)<=48903];
objective = norm(Y-1,1) + norm(U,1)*0.01;
G = [3.98 <= W <= 4.02]
%G=[W==1]
[Frobust,h] = robustify(F + G,objective,[],W);
xk = [0;0;0;0;0;0;0];
ops = sdpsettings;
for i = 1:25
    optimize([Frobust, x == xk(:,end)],h,ops);
    xk = [xk A*xk(:,end) + B*[value(U(1));value(U(2))] + E*(-1+2*rand(1))];
end
plot(C*xk)


                                                                                                           [Warning]





+++++++++++++++++++++++++++++++++++++++
|   ID|                     Constraint|
+++++++++++++++++++++++++++++++++++++++
|   #1|   Element-wise inequality 20x1|
+++++++++++++++++++++++++++++++++++++++
***** Starting YALMIP robustification module. *********************
 - Detected 10 uncertain variables
 - Detected 10 independent group(s) of uncertain variables
 - Using possibly conservative approach to deal with uncertainty dependent auxilliary variables.
 - (change robust.auxreduce to 'projection', 'enumeration' for exact solution.)
 - Eliminating uncertainty using explicit maximization of inf-norm
***** Derivation of robust counterpart done ***********************
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
Exiting: One or more of the residuals, duality gap, or total relative error
 has stalled:
         the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6.
>> 

[Johan Löfberg]

You said w was known, so why are you setting up (an intractable) robust optimization problem?

The problem is simply infeasible already in the initial state

>> optimize([Frobust, x == xk(:,end)],h,ops);

presolving:

presolving (1 rounds: 1 fast, 0 medium, 0 exhaustive):

 7 deleted vars, 13 deleted constraints, 0 added constraints, 17 tightened bounds, 0 added holes, 1 changed sides, 0 changed coefficients

 0 implications, 0 cliques

presolving detected infeasibility

Presolving Time: 0.01

SCIP Status        : problem is solved [infeasible]

Solving Time (sec) : 0.06

Solving Nodes      : 0

Primal Bound       : +1.00000000000000e+20 (0 solutions)

Dual Bound         : -1.00000000000000e+20

Gap                : infinite

[adnan jaffar]

Dr. Johan,
Sorry i got your point now. I have revised the code but again i am getting

 "the primal appears to be infeasible and the dual unbounded since
         the dual objective > 1e+10
         and the primal objective > -1e+6 "

Here is the Code

A=eye(7)+0.05*[-0.350312629 0.001652708 -22.38023628 0 -0.007 0 0;119.9592691 -1.236693495 16746.7579 0 0 0 -0.001872252;1 -0.010309278 0 0 0 0 0;0 0 0 0 1 0 0;0.000212974    0 0 -25.98376792 -0.055356012 -7.48E-06    0;0    0 0 0 0 -2 0;0 0 0 0 0 0 -10];

Bb=0.05*[0 0;0 0;0 0;0 0;0 0;2 0;0 10];

C=[1 0 0 0 0 0 0;0 1 0 0 0 0 0;0 0 1 0 0 0 0;0 0 0 1 0 0 0;0 0 0 0 1 0 0;0 0 0 0 0 1 0;0 0 0 0 0 0 1;0 3776 0 0 0 66.08    0];
E=0.05*[0.007;0;0;0;0.000212974;0;0];

B=[Bb E];

nx = 7; % Number of states
nu = 2; % Number of inputs

% MPC data
Q = eye(7);
%R = eye(2);
R=diag([2,5]);
N = 30;

% Initial state
x0 = [0;0;0;0;0;0;0];
 xmin = [669.3;-41.25;-0.5;0;-15000;-1;-8];
    xmax = [1173.7;-72.19;0.5;47403;15000;90;8];

    umin = [-5; -1500];
    umax = [94; 48903];
u = sdpvar(repmat(nu,1,N),repmat(1,1,N));

constraints = [];
objective = 0;
x = x0;

for k = 1:N

 x = A*x + B*[u{k};4];
 objective = objective + norm(Q*x,1) + norm(R*u{k},1);
 constraints = [constraints, umin <= u{k}<= umax, xmin<=x<=xmax];
end
optimize(constraints,objective);
value(u{1})



Values obtained are
153.4074
 -915.6646

[Johan Löfberg]

Do you understand what the word "infeasible" means in the context of optimization? Your problem is simply infeasible. Trivially infeasible even, as you cannot even take one step and get into/stay in the feasible region

optimize([xmin <= A*x0 + B*[u{k};4] <= xmax, umin <= u{1}<= umax])

To me, xmax(2) looks odd (x has to be very negative, while you initial state is 0)