-
Notifications
You must be signed in to change notification settings - Fork 9
/
Copy pathscript_generate_feasible_region_modified.m
177 lines (137 loc) · 4.89 KB
/
script_generate_feasible_region_modified.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
%% Script - An attempt (with failure) to use yalmip plot to get the feasible region
clear all;close all;clc;
%% Call yalmpi and solvers - This is an example - I am running Matlab on a mac.
addpath(genpath('~/Documents/MATLAB/yalmip'))
addpath(genpath('~/Documents/MATLAB/cvx/sedumi'))
addpath(genpath('~/Documents/MATLAB/cvx/sdpt3'))
addpath(genpath('~/mosek/8/toolbox/r2014a'));
%% Define Constants
flagc=1;
disp('*******************************************************************');
disp('* *');
disp('* Dr. Lu''s procedure *');
disp('* *');
disp('*******************************************************************');
disp('');
if flagc == 1
disp('A Constrained problem will be solved');
disp('');
disp('LMIs 11, 15, 16, 20, 21, 23 and 25 will be used');
elseif flagc == 2
disp('A Constrained problem will be solved');
disp('');
disp('LMIs 11, 15, 16, 20, 21, 23, 24, 25 and 26 will be used');
else
disp('An unconstrained problem will be solved');
disp('');
disp('LMIs 11, 15 and 16 will be used');
end
disp('');
disp('Feasible Region using plot - lmilab is recommended but mosek works');
disp('');
disp('*******************************************************************');
disp('');
opts = sdpsettings('verbose',1, 'warning',1,'solver','mosek');
model_parameters_jianbo
nq=nx;
nr=nu;
umax=1;
xmax=0.7;
rk=1;
N=1;
define_lmi_variables_jianbo
xk=sdpvar(2,1);
B{2,1}=input('Enter value for B{2,1} (default =[0;1]) => ');
if isempty(B{2,1})
B{2,1}=[0;1];
end
%% LMI 11
pre_build_lmi11_jianbo
actual_build_lmi11_jianbo
%% LMI 15
pre_build_lmi15_jianbo
actual_build_lmi15_jianbo
%% LMI 16
pre_build_lmi16_jianbo
actual_build_lmi16_jianbo
%% Big LMI
if flagc == 1 % Constrained on the input
pre_build_lmi20_jianbo
actual_build_lmi20_jianbo
pre_build_lmi21_jianbo
actual_build_lmi21_jianbo
pre_build_lmi23_jianbo
actual_build_lmi23_jianbo
pre_build_lmi25_jianbo
actual_build_lmi25_jianbo
%LMIs_jianbo=[biglmi11;biglmi15;biglmi16;biglmi20;biglmi21;biglmi23;biglmi25,biglmi25_1];
LMIs_jianbo=[biglmi20;biglmi21;biglmi23;biglmi25;biglmi25_1];
elseif flagc == 2 % Constrained on the input and states
pre_build_lmi20_jianbo
actual_build_lmi20_jianbo
pre_build_lmi21_jianbo
actual_build_lmi21_jianbo
pre_build_lmi23_jianbo
actual_build_lmi23_jianbo
pre_build_lmi24_jianbo
actual_build_lmi24_jianbo
pre_build_lmi25_jianbo
actual_build_lmi25_jianbo
pre_build_lmi26_jianbo
actual_build_lmi26_jianbo
LMIs_jianbo=[biglmi11;biglmi15;biglmi16;biglmi20;biglmi21;biglmi23;biglmi24;biglmi25,biglmi25_1;biglmi26,biglmi26_1];
else % unconstrained
LMIs_jianbo=[biglmi11;biglmi15;biglmi16];
end
%% Plot
npoints=100;
tfig=0;
tfig=tfig+1;
figure(tfig);
p=plot(LMIs_jianbo,xk,'b',npoints,opts);
p=cell2mat(p)';
fill(p(:,1),p(:,2),'blue','LineWidth',2,'FaceAlpha',.3);
if flagc == 1
title(sprintf('Feasible Region - Input constraint with u_{max} = %g',umax));
else
title(sprintf('Feasible Region - Input and State constraints with u_{max} = %g and x_{max} = %g',umax,xmax));
end
xlabel('x_1(0)');ylabel('x_2(0)');
% Doing it again just to save the points
if flagc == 1
s=sprintf('points_feasreg_jianbo_u_constraint_N_%d_umax_%s=p;',N,strrep(num2str(umax),'.','_'));
elseif flagc == 2
s=sprintf('points_feasreg_jianbo_u_and_x_constraints_N_%d_umax_%s_xmax_%s=p;',N,strrep(num2str(umax),'.','_'),strrep(num2str(xmax),'.','_'));
else
s=sprintf('points_feasreg_jianbo_no_constraints_N_%d=p;',N);
end
eval(s);
if exist('data','dir') ~= 7 % Please see help for exist
mkdir('data')
end
if flagc == 1
s=sprintf('save data/pfeasreg_jianbo_u_constraint_N_%d_umax_%s points_feasreg_jianbo_u_constraint_N_%d_umax_%s',N,strrep(num2str(umax),'.','_'),N,strrep(num2str(umax),'.','_'));
elseif flagc == 2
s=sprintf('save data/pfeasreg_jianbo_u_and_x_constraints_N_%d_umax_%s_xmax_%s.mat points_feasreg_jianbo_u_and_x_constraints_N_%d_umax_%s_xmax_%s',...
N,strrep(num2str(umax),'.','_'),strrep(num2str(xmax),'.','_'),N,strrep(num2str(umax),'.','_'),strrep(num2str(xmax),'.','_'));
else
s=sprintf('save data/pfeasreg_jianbo_no_constraints_N_%d.mat points_feasreg_jianbo_no_constraints_N_%d',N.N);
end
eval(s)
disp('Feasible region saved to dir data');
%% Images
% if exist('images','dir') ~= 7 % Please see help for exist
% mkdir('images')
% end
%
% for i=1:tfig
% if flagc == 1
% s=sprintf('images/feasregp_jianbo_%d_u_constraint_N_%d_umax_%g.png',i,N,umax);
% elseif flagc == 2
% s=sprintf('images/feasregp_jianbo_%d_u_and_x_constraints_N_%d_umax_%g_xmax_%g.png',i,N,umax,xmax);
% else
% s=sprintf('images/feasregp_jianbo_%d_no_constraints_N_%d.png',i,N);
% end
% figure(i);
% print(s,'-dpng');
% end