[37] | 1 | function varargout = geomean(varargin) |
---|
| 2 | %GEOMEAN (overloaded) |
---|
| 3 | % |
---|
| 4 | % t = GEOMEAN(X) |
---|
| 5 | % |
---|
| 6 | % For Hermitian matrix X, returns det(X)^(1/length(X)) |
---|
| 7 | % |
---|
| 8 | % For real vector X, returns prod(X)^(1/length(X)) |
---|
| 9 | % |
---|
| 10 | % This concave function is monotonically growing in det(P) |
---|
| 11 | % for P>0, so it can be used for maximizing det(P), |
---|
| 12 | % or to add lower bound constraints on the determinant. |
---|
| 13 | % |
---|
| 14 | % When GEOMEAN is used in a problem, the domain constraint |
---|
| 15 | % set(X>0) is automatically added to the problem. |
---|
| 16 | % |
---|
| 17 | % If you only use GEOMEAN as the objective function, you |
---|
| 18 | % may want to consider GEOMEAN2 which may results in a |
---|
| 19 | % slightly smaller SDP model. |
---|
| 20 | % |
---|
| 21 | % See also SDPVAR, GEOMEAN2, SUMK, SUMABSK |
---|
| 22 | |
---|
| 23 | % Author Johan Löfberg |
---|
| 24 | % $Id: geomean.m,v 1.1 2006/08/10 18:00:20 joloef Exp $ |
---|
| 25 | |
---|
| 26 | switch class(varargin{1}) |
---|
| 27 | |
---|
| 28 | case 'sdpvar' % Overloaded operator for SDPVAR objects. Pass on args and save them. |
---|
| 29 | X = varargin{1}; |
---|
| 30 | [n,m] = size(X); |
---|
| 31 | if is(varargin{1},'hermitian') | min(n,m)==1 |
---|
| 32 | varargout{1} = yalmip('addextendedvariable',mfilename,varargin{:}); |
---|
| 33 | else |
---|
| 34 | % Create one variable for each column |
---|
| 35 | y = []; |
---|
| 36 | for i = 1:m |
---|
| 37 | index = (1+n*(i-1)):i*n; |
---|
| 38 | x = extsubsref(X,index); |
---|
| 39 | y = [y yalmip('addextendedvariable',mfilename,x)]; |
---|
| 40 | end |
---|
| 41 | varargout{1} = y; |
---|
| 42 | end |
---|
| 43 | |
---|
| 44 | case 'char' % YALMIP send 'model' when it wants the epigraph or hypograph |
---|
| 45 | if isequal(varargin{1},'graph') |
---|
| 46 | t = varargin{2}; % Second arg is the extended operator variable |
---|
| 47 | X = varargin{3}; % Third arg and above are the args user used when defining t. |
---|
| 48 | |
---|
| 49 | % Extend if not power of 2 |
---|
| 50 | n = length(X); |
---|
| 51 | m=2^ceil(log2(n)); |
---|
| 52 | X0=X; |
---|
| 53 | F = set([]); |
---|
| 54 | if n ~= m |
---|
| 55 | d=m-n; |
---|
| 56 | if size(X,1)==size(X,2) |
---|
| 57 | % model determinant |
---|
| 58 | % Convert to a real problem |
---|
| 59 | if is(X,'complex') |
---|
| 60 | X = [real(X) imag(X);-imag(X) real(X)]; |
---|
| 61 | n = length(X); |
---|
| 62 | % We will now get (detX)^2, so we need to pad |
---|
| 63 | % differently to get (detX)^(1/n) |
---|
| 64 | d=2*m-n; |
---|
| 65 | end |
---|
| 66 | |
---|
| 67 | D = tril(sdpvar(n,n)); |
---|
| 68 | delta = diag(D); |
---|
| 69 | F = set([X D;D' diag(delta)] > 0); |
---|
| 70 | p = 2^ceil(log2(n)); |
---|
| 71 | if 1 |
---|
| 72 | x = [delta;ones(d,1)*t]; |
---|
| 73 | else |
---|
| 74 | % More efficient in detset, but not tested |
---|
| 75 | % sufficiently yet |
---|
| 76 | x = delta; |
---|
| 77 | end |
---|
| 78 | |
---|
| 79 | elseif size(X,1)>size(X,2) |
---|
| 80 | x = [X;ones(d,1)*t]; |
---|
| 81 | |
---|
| 82 | else |
---|
| 83 | x = [X ones(1,d)*t]; |
---|
| 84 | end |
---|
| 85 | else |
---|
| 86 | if size(X,1)==size(X,2) |
---|
| 87 | % model determinant |
---|
| 88 | % Convert to a real problem |
---|
| 89 | if is(X,'complex') |
---|
| 90 | X = [real(X) imag(X);-imag(X) real(X)]; |
---|
| 91 | n = length(X); |
---|
| 92 | % We will now get (detX)^2, so we need to pad |
---|
| 93 | % differently to get (detX)^(1/n) |
---|
| 94 | d=2*m-n; |
---|
| 95 | end |
---|
| 96 | |
---|
| 97 | D = tril(sdpvar(n,n)); |
---|
| 98 | delta = diag(D); |
---|
| 99 | F = set([X D;D' diag(delta)] > 0); |
---|
| 100 | p = 2^ceil(log2(n)); |
---|
| 101 | x = delta; |
---|
| 102 | |
---|
| 103 | elseif size(X,1)>size(X,2) |
---|
| 104 | x = X; |
---|
| 105 | |
---|
| 106 | elseif size(X,1)<size(X,2) |
---|
| 107 | x = X; |
---|
| 108 | end |
---|
| 109 | end |
---|
| 110 | |
---|
| 111 | varargout{1} = F + detset(t,x); |
---|
| 112 | |
---|
| 113 | if (min(size(X))>1) & ishermitian(X0) |
---|
| 114 | varargout{2} = struct('convexity','concave','monotonicity','none','definiteness','positive'); |
---|
| 115 | else |
---|
| 116 | varargout{2} = struct('convexity','concave','monotonicity','increasing','definiteness','positive'); |
---|
| 117 | end |
---|
| 118 | varargout{3} = X0; %We have altered the original input, lucky we saved it! |
---|
| 119 | elseif isequal(varargin{1},'milp') |
---|
| 120 | |
---|
| 121 | varargout{1} = []; |
---|
| 122 | varargout{2} = []; |
---|
| 123 | varargout{3} = []; |
---|
| 124 | |
---|
| 125 | else |
---|
| 126 | end |
---|
| 127 | otherwise |
---|
| 128 | end |
---|