1 | function sys = sdpvar(varargin) |
---|
2 | %SDPVAR Create symbolic decision variable |
---|
3 | % |
---|
4 | % You can create a sdpvar variable by: |
---|
5 | % X = SDPVAR(n) Symmetric nxn matrix |
---|
6 | % X = SDPVAR(n,n) Symmetric nxn matrix |
---|
7 | % X = SDPVAR(n,m) Full nxm matrix (n~=m) |
---|
8 | % |
---|
9 | % Definition of multiple scalars can be simplified |
---|
10 | % SDPVAR x y z w |
---|
11 | % |
---|
12 | % The parametrizations supported are |
---|
13 | % X = SDPVAR(n,n,'full') Full nxn matrix |
---|
14 | % X = SDPVAR(n,n,'symmetric') Symmetric nxn matrix |
---|
15 | % X = SDPVAR(n,n,'toeplitz') Symmetric Toeplitz |
---|
16 | % X = SDPVAR(n,n,'hankel') Symmetric Hankel |
---|
17 | % X = SDPVAR(n,n,'skew') Skew-symmetric |
---|
18 | % |
---|
19 | % The letters 'sy','f','ha', 't' and 'sk' are searched for in the third argument |
---|
20 | % hence sdpvar(n,n,'toeplitz') gives the same result as sdpvar(n,n,'t') |
---|
21 | % |
---|
22 | % Only square Toeplitz and Hankel matries are supported |
---|
23 | % |
---|
24 | % A scalar is defined as a 1x1 matrix |
---|
25 | % |
---|
26 | % Higher-dimensional matrices are also supported, although this currently |
---|
27 | % is an experimental feature with limited use. The type flag applies to |
---|
28 | % the lowest level slice. |
---|
29 | % |
---|
30 | % X = SDPVAR(n,n,n,'full') Full nxnxn matrix |
---|
31 | % |
---|
32 | % In addition to the matrix type, a fourth argument |
---|
33 | % can be used to obtain a complex matrix. All the |
---|
34 | % matrix types above apply to a complex matrix, and |
---|
35 | % in addition a Hermitian type is added |
---|
36 | % |
---|
37 | % X = SDPVAR(n,n,'hermitian','complex') Complex Hermitian nxn matrix (X=X'=conj(X.')) |
---|
38 | % |
---|
39 | % The other types are obtained as above |
---|
40 | % X = SDPVAR(n,n,'symmetric','complex') Complex symmetric nxn matrix (X=X.') |
---|
41 | % X = SDPVAR(n,n,'full','complex') Complex full nxn matrix |
---|
42 | % ... and the same for Toeplitz, Hankel and skew-symmetric |
---|
43 | % |
---|
44 | % See also @SDPVAR/SET, INTVAR, BINVAR, methods('sdpvar'), SEE |
---|
45 | |
---|
46 | % Author Johan Löfberg |
---|
47 | % $Id: sdpvar.m,v 1.46 2006/08/11 15:38:20 joloef Exp $ |
---|
48 | |
---|
49 | superiorto('double'); |
---|
50 | if nargin==0 |
---|
51 | return |
---|
52 | end |
---|
53 | |
---|
54 | if isstruct(varargin{1}) |
---|
55 | sys = class(varargin{1},'sdpvar'); |
---|
56 | return |
---|
57 | end |
---|
58 | |
---|
59 | %if nargin == 1 |
---|
60 | |
---|
61 | % To speed up dualization, we keep track of primal SDP cones |
---|
62 | % [0 0] : Nothing known (cleared in some operator, or none-cone to start with) |
---|
63 | % [1 0] : Primal cone |
---|
64 | % [1 1] : Primal cone + DOUBLE |
---|
65 | % [1 2 x] : Primal cone + SDPVAR |
---|
66 | % [-1 1] : -Primal cone + DOUBLE |
---|
67 | % [-1 2 x] : -Primal cone + SDPVAR |
---|
68 | |
---|
69 | conicinfo = [0 0]; |
---|
70 | |
---|
71 | if ischar(varargin{1}) |
---|
72 | switch varargin{1} |
---|
73 | case 'clear' |
---|
74 | disp('Obsolete comand'); |
---|
75 | return |
---|
76 | case 'nvars' |
---|
77 | sys = yalmip('nvars');%THIS IS OBSAOLETE AND SHOULD NOT BE USED |
---|
78 | return |
---|
79 | otherwise |
---|
80 | n = length(varargin); |
---|
81 | varnames = varargin; |
---|
82 | for k = 1:n |
---|
83 | varcmd{k}='(1,1)'; |
---|
84 | lp=findstr(varargin{k},'('); |
---|
85 | rp=findstr(varargin{k},')'); |
---|
86 | if isempty(lp) & isempty(rp) |
---|
87 | if ~isvarname(varargin{k}) |
---|
88 | error('Not a valid variable name.') |
---|
89 | end |
---|
90 | else |
---|
91 | if (~isempty(lp))&(~isempty(rp)) |
---|
92 | if min(lp)<max(rp) |
---|
93 | varnames{k} = varargin{k}(1:lp-1); |
---|
94 | varcmd{k}=varargin{k}(lp:rp); |
---|
95 | else |
---|
96 | error('Not a valid variable name.') |
---|
97 | end |
---|
98 | else |
---|
99 | error('Not a valid variable name.') |
---|
100 | end |
---|
101 | end |
---|
102 | end |
---|
103 | for k = 1:n |
---|
104 | if isequal(varnames{k},'i') | isequal(varnames{k},'j') |
---|
105 | if length(dbstack) == 1 |
---|
106 | assignin('caller',varnames{k},eval(['sdpvar' varcmd{k}])); |
---|
107 | else |
---|
108 | error(['Due to a bug in MATLAB, use ' varnames{k} ' = sdpvar' varcmd{k} ' instead.']); |
---|
109 | end |
---|
110 | else |
---|
111 | assignin('caller',varnames{k},eval(['sdpvar' varcmd{k}])); |
---|
112 | end |
---|
113 | end |
---|
114 | return |
---|
115 | end |
---|
116 | end |
---|
117 | |
---|
118 | % ************************************************************************* |
---|
119 | % Maybe new NDSDPVAR syntax |
---|
120 | % ************************************************************************* |
---|
121 | if nargin > 2 & isa(varargin{3},'double') & ~isempty(varargin{3}) |
---|
122 | sys = ndsdpvar(varargin{:}); |
---|
123 | return |
---|
124 | end |
---|
125 | |
---|
126 | |
---|
127 | % Supported matrix types |
---|
128 | % - symm |
---|
129 | % - full |
---|
130 | % - skew |
---|
131 | % - hank |
---|
132 | % - toep |
---|
133 | switch nargin |
---|
134 | case 1 %Bug in MATLAB 5.3!! sdpvar called from horzcat!!!!???? |
---|
135 | if isempty(varargin{1}) |
---|
136 | sys = varargin{1}; |
---|
137 | return |
---|
138 | end |
---|
139 | if isa(varargin{1},'sdpvar') |
---|
140 | sys = varargin{1}; |
---|
141 | sys.typeflag = 0; |
---|
142 | return |
---|
143 | end |
---|
144 | n = varargin{1}; |
---|
145 | m = varargin{1}; |
---|
146 | if sum(n.*m)==0 |
---|
147 | sys = zeros(n,m); |
---|
148 | return |
---|
149 | end |
---|
150 | if (n==m) |
---|
151 | matrix_type = 'symm'; |
---|
152 | nvar = sum(n.*(n+1)/2); |
---|
153 | conicinfo = [1 0]; |
---|
154 | else |
---|
155 | matrix_type = 'full'; |
---|
156 | nvar = sum(n.*m); |
---|
157 | end |
---|
158 | case 2 |
---|
159 | n = varargin{1}; |
---|
160 | m = varargin{2}; |
---|
161 | if length(n)~=length(m) |
---|
162 | error('The dimensions must have the same lengths') |
---|
163 | end |
---|
164 | if sum(n.*m)==0 |
---|
165 | sys = zeros(n,m); |
---|
166 | return |
---|
167 | end |
---|
168 | if (n==m) |
---|
169 | matrix_type = 'symm'; |
---|
170 | nvar = sum(n.*(n+1)/2); |
---|
171 | conicinfo = [1 0]; |
---|
172 | else |
---|
173 | matrix_type = 'full'; |
---|
174 | nvar = sum(n.*m); |
---|
175 | end |
---|
176 | case {3,4} |
---|
177 | n = varargin{1}; |
---|
178 | m = varargin{2}; |
---|
179 | if sum(n.*m)==0 |
---|
180 | sys = zeros(n,m); |
---|
181 | return |
---|
182 | end |
---|
183 | |
---|
184 | % Check for complex or real |
---|
185 | if (nargin == 4) |
---|
186 | if isempty(varargin{4}) |
---|
187 | varargin{4} = 'real'; |
---|
188 | else |
---|
189 | if ~ischar(varargin{4}) |
---|
190 | help sdpvar |
---|
191 | error('Fourth argument should be ''complex'' or ''real''') |
---|
192 | end |
---|
193 | end |
---|
194 | index_cmrl = strmatch(varargin{4},{'real','complex'}); |
---|
195 | if isempty(index_cmrl) |
---|
196 | error('Fourth argument should be ''complex'' or ''real''. See help above') |
---|
197 | end |
---|
198 | else |
---|
199 | if ~ischar(varargin{3}) |
---|
200 | help sdpvar |
---|
201 | error('Third argument should be ''symmetric'', ''full'', ''hermitian'',...See help above') |
---|
202 | end |
---|
203 | index_cmrl = 1; |
---|
204 | end; |
---|
205 | |
---|
206 | if isempty(varargin{3}) |
---|
207 | if n==m |
---|
208 | index_type = 7; %Default symmetric |
---|
209 | else |
---|
210 | index_type = 4; |
---|
211 | end |
---|
212 | else |
---|
213 | if ~isempty(strmatch(varargin{3},{'complex','real'})) |
---|
214 | % User had third argument as complex or real |
---|
215 | error(['Third argument should be ''symmetric'', ''full'', ''toeplitz''... Maybe you meant sdpvar(n,n,''full'',''' varargin{3} ''')']) |
---|
216 | end |
---|
217 | index_type = strmatch(varargin{3},{'toeplitz','hankel','symmetric','full','rhankel','skew','hermitian'}); |
---|
218 | end |
---|
219 | |
---|
220 | if isempty(index_type) |
---|
221 | error(['Matrix type "' varargin{3} '" not supported']) |
---|
222 | else |
---|
223 | switch index_type+100*(index_cmrl-1) |
---|
224 | case 1 |
---|
225 | if n~=m |
---|
226 | error('Toeplitz matrix must be square') |
---|
227 | else |
---|
228 | matrix_type = 'toep'; |
---|
229 | nvar = n; |
---|
230 | end |
---|
231 | |
---|
232 | case 2 |
---|
233 | if n~=m |
---|
234 | error('Hankel matrix must be square') |
---|
235 | else |
---|
236 | matrix_type = 'hank'; |
---|
237 | nvar = n; |
---|
238 | end |
---|
239 | |
---|
240 | case 3 |
---|
241 | if n~=m |
---|
242 | error('Symmetric matrix must be square') |
---|
243 | else |
---|
244 | matrix_type = 'symm'; |
---|
245 | nvar = sum(n.*(n+1)/2); |
---|
246 | end |
---|
247 | |
---|
248 | case 4 |
---|
249 | matrix_type = 'full'; |
---|
250 | nvar = sum(n.*m); |
---|
251 | if nvar==1 |
---|
252 | matrix_type = 'symm'; |
---|
253 | end |
---|
254 | |
---|
255 | case 5 |
---|
256 | if n~=m |
---|
257 | error('Hankel matrix must be square') |
---|
258 | else |
---|
259 | matrix_type = 'rhankel'; |
---|
260 | nvar = 2*n-1; |
---|
261 | end |
---|
262 | |
---|
263 | case 6 |
---|
264 | if n~=m |
---|
265 | error('Skew symmetric matrix must be square') |
---|
266 | else |
---|
267 | matrix_type = 'skew'; |
---|
268 | nvar = (n*(n+1)/2)-n; |
---|
269 | end |
---|
270 | |
---|
271 | case 7 |
---|
272 | if n~=m |
---|
273 | error('Symmetric matrix must be square') |
---|
274 | else |
---|
275 | matrix_type = 'symm'; |
---|
276 | nvar = n*(n+1)/2; |
---|
277 | end |
---|
278 | |
---|
279 | case 101 |
---|
280 | if n~=m |
---|
281 | error('Toeplitz matrix must be square') |
---|
282 | else |
---|
283 | matrix_type = 'toep complex'; |
---|
284 | nvar = 2*n; |
---|
285 | end |
---|
286 | |
---|
287 | case 102 |
---|
288 | if n~=m |
---|
289 | error('Hankel matrix must be square') |
---|
290 | else |
---|
291 | matrix_type = 'hank complex'; |
---|
292 | nvar = (2*n); |
---|
293 | end |
---|
294 | |
---|
295 | case 103 |
---|
296 | if n~=m |
---|
297 | error('Symmetric matrix must be square') |
---|
298 | else |
---|
299 | matrix_type = 'symm complex'; |
---|
300 | nvar = 2*n*(n+1)/2; |
---|
301 | end |
---|
302 | |
---|
303 | case 104 |
---|
304 | matrix_type = 'full complex'; |
---|
305 | nvar = 2*n*m; |
---|
306 | if nvar==1 |
---|
307 | matrix_type = 'symm complex'; |
---|
308 | end |
---|
309 | |
---|
310 | case 105 |
---|
311 | if n~=m |
---|
312 | error('Hankel matrix must be square') |
---|
313 | else |
---|
314 | matrix_type = 'rhankel complex'; |
---|
315 | nvar = 2*(2*n-1); |
---|
316 | end |
---|
317 | |
---|
318 | case 106 |
---|
319 | if n~=m |
---|
320 | error('Skew symmetric matrix must be square') |
---|
321 | else |
---|
322 | matrix_type = 'skew complex'; |
---|
323 | nvar = 2*((n*(n+1)/2)-n); |
---|
324 | end |
---|
325 | |
---|
326 | case 107 |
---|
327 | if n~=m |
---|
328 | error('Hermitian matrix must be square') |
---|
329 | else |
---|
330 | matrix_type = 'herm complex'; |
---|
331 | nvar = n*(n+1)/2+(n*(n+1)/2-n); |
---|
332 | end |
---|
333 | |
---|
334 | |
---|
335 | otherwise |
---|
336 | error('Bug! Report!'); |
---|
337 | end |
---|
338 | |
---|
339 | end |
---|
340 | |
---|
341 | case 5 % Fast version for internal use |
---|
342 | sys.basis = varargin{5}; |
---|
343 | sys.lmi_variables=varargin{4}; |
---|
344 | sys.dim(1) = varargin{1}; |
---|
345 | sys.dim(2) = varargin{2}; |
---|
346 | sys.typeflag = 0; |
---|
347 | sys.savedata = []; |
---|
348 | sys.extra = []; |
---|
349 | sys.extra.expanded = []; |
---|
350 | sys.conicinfo = 0; |
---|
351 | % Find zero-variables |
---|
352 | constants = find(sys.lmi_variables==0); |
---|
353 | if ~isempty(constants); |
---|
354 | sys.lmi_variables(constants)=[]; |
---|
355 | sys.basis(:,1) = sys.basis(:,1) + sum(sys.basis(:,1+constants),2); |
---|
356 | sys.basis(:,1+constants)=[]; |
---|
357 | end |
---|
358 | if isempty(sys.lmi_variables) |
---|
359 | sys = full(reshape(sys.basis(:,1),sys.dim(1),sys.dim(2))); |
---|
360 | else |
---|
361 | sys = class(sys,'sdpvar'); |
---|
362 | end |
---|
363 | return |
---|
364 | case 6 % Fast version for internal use |
---|
365 | sys.basis = varargin{5}; |
---|
366 | sys.lmi_variables=varargin{4}; |
---|
367 | sys.dim(1) = varargin{1}; |
---|
368 | sys.dim(2) = varargin{2}; |
---|
369 | sys.typeflag = varargin{6}; |
---|
370 | sys.savedata = []; |
---|
371 | sys.extra = []; |
---|
372 | sys.extra.expanded = []; |
---|
373 | sys.conicinfo = 0; |
---|
374 | % Find zero-variables |
---|
375 | constants = find(sys.lmi_variables==0); |
---|
376 | if ~isempty(constants); |
---|
377 | sys.lmi_variables(constants)=[]; |
---|
378 | sys.basis(:,1) = sys.basis(:,1) + sum(sys.basis(:,1+constants),2); |
---|
379 | sys.basis(:,1+constants)=[]; |
---|
380 | end |
---|
381 | if isempty(sys.lmi_variables) |
---|
382 | sys = full(reshape(sys.basis(:,1),sys.dim(1),sys.dim(2))); |
---|
383 | else |
---|
384 | sys = class(sys,'sdpvar'); |
---|
385 | end |
---|
386 | return |
---|
387 | case 7 % Fast version for internal use |
---|
388 | sys.basis = varargin{5}; |
---|
389 | sys.lmi_variables=varargin{4}; |
---|
390 | sys.dim(1) = varargin{1}; |
---|
391 | sys.dim(2) = varargin{2}; |
---|
392 | sys.typeflag = varargin{6}; |
---|
393 | sys.savedata = []; |
---|
394 | sys.extra = varargin{7}; |
---|
395 | sys.extra.expanded = []; |
---|
396 | sys.conicinfo = varargin{7}; |
---|
397 | % Find zero-variables |
---|
398 | constants = find(sys.lmi_variables==0); |
---|
399 | if ~isempty(constants); |
---|
400 | sys.lmi_variables(constants)=[]; |
---|
401 | sys.basis(:,1) = sys.basis(:,1) + sum(sys.basis(:,1+constants),2); |
---|
402 | sys.basis(:,1+constants)=[]; |
---|
403 | end |
---|
404 | if isempty(sys.lmi_variables) |
---|
405 | sys = full(reshape(sys.basis(:,1),sys.dim(1),sys.dim(2))); |
---|
406 | else |
---|
407 | sys = class(sys,'sdpvar'); |
---|
408 | end |
---|
409 | return |
---|
410 | |
---|
411 | otherwise |
---|
412 | error('Wrong number of arguments in sdpvar creation'); |
---|
413 | end |
---|
414 | |
---|
415 | if isempty(n) | isempty(m) |
---|
416 | error('Size must be integer valued') |
---|
417 | end; |
---|
418 | if ~((abs((n-ceil(n)))+ abs((m-ceil(m))))==0) |
---|
419 | error('Size must be integer valued') |
---|
420 | end |
---|
421 | |
---|
422 | [mt,variabletype] = yalmip('monomtable'); |
---|
423 | lmi_variables = (1:nvar)+size(mt,1); |
---|
424 | |
---|
425 | for blk = 1:length(n) |
---|
426 | switch matrix_type |
---|
427 | |
---|
428 | case 'full' |
---|
429 | basis{blk} = [spalloc(n(blk)*m(blk),1,0) speye(n(blk)*m(blk))];%speye(nvar)]; |
---|
430 | |
---|
431 | case 'full complex' |
---|
432 | basis = [spalloc(n*m,1,0) speye(nvar/2) speye(nvar/2)*sqrt(-1)]; |
---|
433 | |
---|
434 | case 'symm' |
---|
435 | if 0 |
---|
436 | basis = spalloc(n^2,1+nvar,n^2); |
---|
437 | l = 2; |
---|
438 | an_empty = spalloc(n,n,2); |
---|
439 | for i=1:n |
---|
440 | temp = an_empty; |
---|
441 | temp(i,i)=1; |
---|
442 | basis(:,l)=temp(:); |
---|
443 | l = l+1; |
---|
444 | for j=i+1:n, |
---|
445 | temp = an_empty; |
---|
446 | temp(i,j)=1; |
---|
447 | temp(j,i)=1; |
---|
448 | basis(:,l)=temp(:); |
---|
449 | l = l+1; |
---|
450 | end |
---|
451 | end |
---|
452 | else |
---|
453 | % Hrm...fast but completely f*d up |
---|
454 | |
---|
455 | Y = reshape(1:n(blk)^2,n(blk),n(blk)); |
---|
456 | Y = tril(Y); |
---|
457 | Y = (Y+Y')-diag(sparse(diag(Y))); |
---|
458 | [uu,oo,pp] = unique(Y(:)); |
---|
459 | if 1 |
---|
460 | basis{blk} = sparse(1:n(blk)^2,pp+1,1); |
---|
461 | else |
---|
462 | basis{blk} = lazybasis(n^2,1+(n*(n+1)/2),1:n(blk)^2,pp+1,ones(n(blk)^2,1)); |
---|
463 | end |
---|
464 | |
---|
465 | end |
---|
466 | |
---|
467 | case 'symm complex' |
---|
468 | basis = spalloc(n^2,1+nvar,2); |
---|
469 | l = 2; |
---|
470 | an_empty = spalloc(n,n,2); |
---|
471 | for i=1:n |
---|
472 | temp = an_empty; |
---|
473 | temp(i,i)=1; |
---|
474 | basis(:,l)=temp(:); |
---|
475 | l = l+1; |
---|
476 | for j=i+1:n, |
---|
477 | temp = an_empty; |
---|
478 | temp(i,j)=1; |
---|
479 | temp(j,i)=1; |
---|
480 | basis(:,l)=temp(:); |
---|
481 | l = l+1; |
---|
482 | end |
---|
483 | end |
---|
484 | for i=1:n |
---|
485 | temp = an_empty; |
---|
486 | temp(i,i)=sqrt(-1); |
---|
487 | basis(:,l)=temp(:); |
---|
488 | l = l+1; |
---|
489 | for j=i+1:n, |
---|
490 | temp = an_empty; |
---|
491 | temp(i,j)=sqrt(-1); |
---|
492 | temp(j,i)=sqrt(-1); |
---|
493 | basis(:,l)=temp(:); |
---|
494 | l = l+1; |
---|
495 | end |
---|
496 | end |
---|
497 | |
---|
498 | case 'herm complex' |
---|
499 | basis = spalloc(n^2,1+nvar,2); |
---|
500 | l = 2; |
---|
501 | an_empty = spalloc(n,n,2); |
---|
502 | for i=1:n |
---|
503 | temp = an_empty; |
---|
504 | temp(i,i)=1; |
---|
505 | basis(:,l)=temp(:); |
---|
506 | l = l+1; |
---|
507 | for j=i+1:n, |
---|
508 | temp = an_empty; |
---|
509 | temp(i,j)=1; |
---|
510 | temp(j,i)=1; |
---|
511 | basis(:,l)=temp(:); |
---|
512 | l = l+1; |
---|
513 | end |
---|
514 | end |
---|
515 | for i=1:n |
---|
516 | for j=i+1:n, |
---|
517 | temp = an_empty; |
---|
518 | temp(i,j)=sqrt(-1); |
---|
519 | temp(j,i)=-sqrt(-1); |
---|
520 | basis(:,l)=temp(:); |
---|
521 | l = l+1; |
---|
522 | end |
---|
523 | end |
---|
524 | |
---|
525 | case 'skew' |
---|
526 | basis = spalloc(n^2,1+nvar,2); |
---|
527 | l = 2; |
---|
528 | an_empty = spalloc(n,n,2); |
---|
529 | for i=1:n |
---|
530 | for j=i+1:n, |
---|
531 | temp = an_empty; |
---|
532 | temp(i,j)=1; |
---|
533 | temp(j,i)=-1; |
---|
534 | basis(:,l)=temp(:); |
---|
535 | l = l+1; |
---|
536 | end |
---|
537 | end |
---|
538 | |
---|
539 | case 'skew complex' |
---|
540 | basis = spalloc(n^2,1+nvar,2); |
---|
541 | l = 2; |
---|
542 | an_empty = spalloc(n,n,2); |
---|
543 | for i=1:n |
---|
544 | for j=i+1:n, |
---|
545 | temp = an_empty; |
---|
546 | temp(i,j)=1; |
---|
547 | temp(j,i)=-1; |
---|
548 | basis(:,l)=temp(:); |
---|
549 | l = l+1; |
---|
550 | end |
---|
551 | end |
---|
552 | for i=1:n |
---|
553 | for j=i+1:n, |
---|
554 | temp = an_empty; |
---|
555 | temp(i,j)=sqrt(-1); |
---|
556 | temp(j,i)=-sqrt(-1); |
---|
557 | basis(:,l)=temp(:); |
---|
558 | l = l+1; |
---|
559 | end |
---|
560 | end |
---|
561 | |
---|
562 | case 'toep' |
---|
563 | basis = spalloc(n^2,1+nvar,2); |
---|
564 | an_empty = spalloc(n,1,1); |
---|
565 | for i=1:n, |
---|
566 | v = an_empty; |
---|
567 | v(i)=1; |
---|
568 | temp = sparse(toeplitz(v)); |
---|
569 | basis(:,i+1) = temp(:); |
---|
570 | end |
---|
571 | |
---|
572 | % Notice, complex Toeplitz not Hermitian |
---|
573 | case 'toep complex' |
---|
574 | basis = spalloc(n^2,1+nvar,2); |
---|
575 | an_empty = spalloc(n,1,1); |
---|
576 | for i=1:n, |
---|
577 | v = an_empty; |
---|
578 | v(i)=1; |
---|
579 | temp = sparse(toeplitz(v)); |
---|
580 | basis(:,i+1) = temp(:); |
---|
581 | end |
---|
582 | for i=1:n, |
---|
583 | v = an_empty; |
---|
584 | v(i)=sqrt(-1); |
---|
585 | temp = sparse(toeplitz(v)); |
---|
586 | basis(:,n+i+1) = temp(:); |
---|
587 | end |
---|
588 | |
---|
589 | case 'hank' |
---|
590 | basis = spalloc(n^2,1+nvar,2); |
---|
591 | an_empty = spalloc(n,1,1); |
---|
592 | for i=1:n, |
---|
593 | v = an_empty; |
---|
594 | v(i)=1; |
---|
595 | temp = sparse(hankel(v)); |
---|
596 | basis(:,i+1) = temp(:); |
---|
597 | end |
---|
598 | |
---|
599 | case 'hank complex' |
---|
600 | basis = spalloc(n^2,1+nvar,2); |
---|
601 | an_empty = spalloc(n,1,1); |
---|
602 | for i=1:n, |
---|
603 | v = an_empty; |
---|
604 | v(i)=1; |
---|
605 | temp = sparse(hankel(v)); |
---|
606 | basis(:,i+1) = temp(:); |
---|
607 | end |
---|
608 | for i=1:n, |
---|
609 | v = an_empty; |
---|
610 | v(i)=sqrt(-1); |
---|
611 | temp = sparse(hankel(v)); |
---|
612 | basis(:,n+i+1) = temp(:); |
---|
613 | end |
---|
614 | |
---|
615 | case 'rhankel' |
---|
616 | basis = spalloc(n^2,1+nvar,2); |
---|
617 | an_empty = spalloc(2*n-1,1,1); |
---|
618 | for i=1:nvar, |
---|
619 | v = an_empty; |
---|
620 | v(i)=1; |
---|
621 | temp = sparse(hankel(v(1:n),[v(n);v(n+1:2*n-1)])); |
---|
622 | basis(:,i+1) = temp(:); |
---|
623 | end |
---|
624 | |
---|
625 | case 'rhankel complex' |
---|
626 | basis = spalloc(n^2,1+nvar,2); |
---|
627 | an_empty = spalloc(2*n-1,1,1); |
---|
628 | for i=1:nvar/2, |
---|
629 | v = an_empty; |
---|
630 | v(i)=1; |
---|
631 | temp = sparse(hankel(v(1:n),[v(n);v(n+1:2*n-1)])); |
---|
632 | basis(:,i+1) = temp(:); |
---|
633 | end |
---|
634 | for i=1:nvar/2, |
---|
635 | v = an_empty; |
---|
636 | v(i)=sqrt(-1); |
---|
637 | temp = sparse(hankel(v(1:n),[v(n);v(n+1:2*n-1)])); |
---|
638 | basis(:,nvar/2+i+1) = temp(:); |
---|
639 | end |
---|
640 | |
---|
641 | otherwise |
---|
642 | error('Bug! Report') |
---|
643 | end |
---|
644 | |
---|
645 | end |
---|
646 | |
---|
647 | % Update monomtable and pre-calculated variable type |
---|
648 | n_mt = size(mt,1); |
---|
649 | m_mt = size(mt,2); |
---|
650 | if min(lmi_variables)>m_mt % New variables |
---|
651 | if size(mt,1)~=size(mt,2) |
---|
652 | mt(size(mt,1),size(mt,1))=0; |
---|
653 | end |
---|
654 | fill=spalloc(size(mt,1),length(lmi_variables),0); |
---|
655 | mt=[mt fill;fill' speye(length(lmi_variables))]; |
---|
656 | else |
---|
657 | mt(lmi_variables,lmi_variables) = speye(length(lmi_variables)); |
---|
658 | end |
---|
659 | variabletype(1,size(mt,1)) = 0; |
---|
660 | yalmip('setmonomtable',mt,variabletype); |
---|
661 | |
---|
662 | % Create an object |
---|
663 | if isa(basis,'cell') |
---|
664 | top = 1; |
---|
665 | for blk = 1:length(n) |
---|
666 | sys{blk}.basis=basis{blk}; |
---|
667 | nn = size(sys{blk}.basis,2)-1; |
---|
668 | sys{blk}.lmi_variables = lmi_variables(top:top+nn-1); |
---|
669 | top = top + nn; |
---|
670 | sys{blk}.dim(1) = n(blk); |
---|
671 | sys{blk}.dim(2) = m(blk); |
---|
672 | sys{blk}.typeflag = 0; |
---|
673 | sys{blk}.savedata = []; |
---|
674 | sys{blk}.extra = []; |
---|
675 | sys{blk}.extra.expanded = []; |
---|
676 | sys{blk}.conicinfo = conicinfo; |
---|
677 | sys{blk} = class(sys{blk},'sdpvar'); |
---|
678 | end |
---|
679 | if length(n)==1 |
---|
680 | sys = sys{1}; |
---|
681 | end |
---|
682 | else |
---|
683 | sys.basis=basis; |
---|
684 | sys.lmi_variables = lmi_variables; |
---|
685 | sys.dim(1) = n; |
---|
686 | sys.dim(2) = m; |
---|
687 | sys.typeflag = 0; |
---|
688 | sys.savedata = []; |
---|
689 | sys.extra = []; |
---|
690 | sys.extra.expanded = []; |
---|
691 | sys.conicinfo = conicinfo; |
---|
692 | sys = class(sys,'sdpvar'); |
---|
693 | if ~isreal(basis) |
---|
694 | % Add internal information about complex pairs |
---|
695 | complex_elements = find(any(imag(basis),2)); |
---|
696 | complex_pairs = []; |
---|
697 | for i = 1:length(complex_elements) |
---|
698 | complex_pairs = [complex_pairs;lmi_variables(find(basis(complex_elements(i),:))-1)]; |
---|
699 | end |
---|
700 | complex_pairs = uniquesafe(complex_pairs,'rows'); |
---|
701 | yalmip('addcomplexpair',complex_pairs); |
---|
702 | end |
---|
703 | end |
---|
704 | |
---|
705 | % Typeflags |
---|
706 | % 0 Standard variable |
---|
707 | % 1 Inequality (LMI) |
---|
708 | % 2 Inequality (element) |
---|
709 | % 3 Equality |
---|
710 | % 4 Cone |
---|
711 | % 5 norm object (osbolete) |
---|
712 | % 6 logdet object |
---|
713 | % 8 KYP object |
---|