source: proiecte/pmake3d/make3d_original/Make3dSingleImageStanford_version0.1/third_party/opt/yalmip/@sdpvar/double.m @ 37

Last change on this file since 37 was 37, checked in by (none), 14 years ago

Added original make3d

File size: 9.2 KB
Line 
1function sys=double(X)
2%DOUBLE Returns current numerical value
3%bajs i bajs
4% Author Johan Löfberg
5% $Id: double.m,v 1.43 2006/09/21 14:34:41 joloef Exp $
6
7solution = yalmip('getsolution');
8lmi_variables = X.lmi_variables;
9opt_variables = solution.variables;
10
11% Definition of nonlinear variables
12[mt,variabletype] = yalmip('monomtable');
13nonlinears = lmi_variables(find(variabletype(X.lmi_variables)));
14
15if ~isempty(solution.values)
16    if max(lmi_variables) <= length(solution.values) & isempty(nonlinears)
17        if ~any(isnan(solution.values(lmi_variables(:))))
18            % Yihoo, we can do this really fast by
19            % re-using the old values
20            sys = X.basis*[1;solution.values(lmi_variables(:))];
21            if X.typeflag==10
22                sys = eig(full(reshape(sys,X.dim(1),X.dim(2))));
23            else
24                sys = full(reshape(sys,X.dim(1),X.dim(2)));
25            end
26            return
27        end
28    end
29end
30
31% Okey, we could not do it really fast...
32
33% Definition of nonlinear variables
34allextended   = yalmip('extvariables');
35allevaluators = yalmip('evalVariables');
36allextended = union(allextended,allevaluators);
37
38if isempty(nonlinears) & isempty(allextended) & all(ismembc(lmi_variables,solution.variables))
39
40    % speed up code for simple linear case
41    values = solution.values;
42    if isempty(solution.values)
43        values = sparse(solution.variables,ones(length(solution.variables),1),solution.optvar,size(mt,1),1);
44        %values = zeros(size(mt,1),1);
45        %values(solution.variables) = solution.optvar;
46        yalmip('setvalues',values);
47    else
48        if any(isnan(solution.values(lmi_variables(:))))
49            values = sparse(solution.variables,ones(length(solution.variables),1),solution.optvar,size(mt,1),1);
50            %            values = zeros(size(mt,1),1);
51            %            values(solution.variables) = solution.optvar;
52            yalmip('setvalues',values);
53        end
54    end
55
56    sys = X.basis*[1;values(lmi_variables(:))];
57    if X.typeflag==10
58        sys = eig(full(reshape(sys,X.dim(1),X.dim(2))));
59    else
60        sys = full(reshape(sys,X.dim(1),X.dim(2)));
61    end
62
63    return
64end
65
66% All double values
67values(size(mt,1),1)=nan;values(:)=nan;
68values(solution.variables) = solution.optvar;
69
70% Evaluate the extended operators
71if ~isempty(allextended)
72    extended_variables = find(ismembc(X.lmi_variables,allextended));
73    if ~isempty(extended_variables)
74        extstructE = lmi_variables(extended_variables);
75        extstructE = yalmip('extstruct',extstructE);
76        if ~isa(extstructE,'cell')
77            extstructE = {extstructE};
78        end
79        for i = 1:length(extended_variables)
80            extvar = lmi_variables(extended_variables(i));
81            %extstruct = yalmip('extstruct',extvar);
82            extstruct = extstructE{i};%yalmip('extstruct',extvar);
83         %   if ~isequal(extstructE{i},extstruct)
84         %       1
85         %   end
86            for k = 1:length(extstruct.arg)
87                if isa(extstruct.arg{k},'sdpvar') | isa(extstruct.arg{k},'constraint')
88                    extstruct.arg{k} = double(extstruct.arg{k});
89                end
90            end
91
92            switch extstruct.fcn
93               
94                case 'sort'
95                    [w,loc] = sort(extstruct.arg{1});
96                    if extstruct.arg{2}.isthisloc
97                        val = loc(extstruct.arg{2}.i);
98                    else
99                        val = w(extstruct.arg{2}.i);
100                    end
101                   
102                case 'semivar'
103                    % A bit messy, since it not really is a nonlinear
104                    % operator. semivar(1) does not return 1, but a new
105                    % semivar variable...
106                    val = values(getvariables(extstruct.var));
107
108                case 'geomean' % Not 100% MATLAB consistent (Hermitian case differ)
109                    val = extstruct.arg{1};
110                    if ~any(any(isnan(val)))
111                        [n,m] = size(val);
112                        if n == m
113                            if issymmetric(val)
114                                val = max(0,real(det(val)))^(1/n);
115                            else
116                                val = geomean(val);
117                            end
118                        else
119                            val = geomean(val);
120                        end
121                    else
122                        val = nan;
123                    end
124
125                case 'pwf'
126                    % Has to be placed here due to the case when
127                    % all functions are double, since in this case,
128                    % we cannot determine in pwf if we want the double or
129                    % create a pw constant function...
130                    n = length(extstruct.arg)/2;
131                    i = 1;
132                    val = nan;
133                    while i<=n
134                        if min(checkset(extstruct.arg{2*i}))>=0
135                            val = extstruct.arg{2*i-1};
136                            break
137                        end
138                        i = i + 1;
139                    end
140
141                case 'mpower'     
142                    val = extstruct.arg{1};
143                case {'max','min'} % Cannot take several inputs, put everything in a vector
144                    val = feval(extstruct.fcn,[extstruct.arg{:}]);
145                case {'or','and'}
146                    try
147                        val = feval(extstruct.fcn,extstruct.arg{:});
148                    catch
149                        val = nan;
150                    end
151
152                otherwise
153                    try
154                        if ismembc(extvar,allevaluators)
155                            % Evaluator variables (used for exp, log, ...)
156                            % have an appended argument that we need to
157                            % remove. The appended variables is used
158                            % internally to convert from general format
159                            % f(a'x+c) to f(z), z==a'z+b
160                            val = feval(extstruct.fcn,extstruct.arg{1:end-1});
161                        else
162                            val = feval(extstruct.fcn,extstruct.arg{:});
163                        end
164                    catch
165                        val = nan;
166                    end
167            end
168            values(extvar) = full(val);
169        end
170    end
171end
172
173if ~isempty(nonlinears)
174    mt_t = mt'; %Working columnwise is faster
175    use_these = find(ismember(lmi_variables,nonlinears));
176    all_extended_variables  = yalmip('extvariables');
177    all_evaluator_variables = yalmip('evalVariables');
178    all_extended_variables = union(all_extended_variables,all_evaluator_variables);
179
180    for i = use_these
181        monom_i = mt_t(:,lmi_variables(i));
182        used_in_monom = find(monom_i);
183
184        if ~isempty(all_extended_variables)
185            extended_variables = find(ismember(used_in_monom,all_extended_variables));
186            if ~isempty(extended_variables)
187                for ii = 1:length(extended_variables)
188                    extvar = used_in_monom(extended_variables(ii));
189                    extstruct = yalmip('extstruct',extvar);
190                    for k = 1:length(extstruct.arg)
191                        if isa(extstruct.arg{k},'sdpvar')
192                            extstruct.arg{k} = double(extstruct.arg{k});
193                        end
194                    end
195
196                    switch extstruct.fcn
197                       
198                        case 'sort'
199                            w = sort(extstruct.arg{1});
200                            val = w(extstruct.arg{2});
201
202                        case 'semivar'
203                            val = values(getvariables(extstruct.var));
204
205                        case 'mpower'      %
206                            val = extstruct.arg{1};
207                        case {'max','min'} % Cannot take several inputs, put everything in a vector
208                            val = feval(extstruct.fcn,[extstruct.arg{:}]);
209                        otherwise          % VAL = OPERATOR(X1,X2,...)
210                            if ismembc(extvar,all_evaluator_variables)
211                                % Evaluator variables (used for exp, log, ...)
212                                % have an appended argument that we need to
213                                % remove. The appended variables is used
214                                % internally to convert from genreal format
215                                % f(a'x+c) to f(z), z==a'z+b
216                                val = feval(extstruct.fcn,extstruct.arg{1:end-1});
217                            else
218                                val = feval(extstruct.fcn,extstruct.arg{:});
219                            end
220                    end
221                    values(extvar) = val;
222                end
223            end
224        end
225        % This code is a bit shaky due to the 0^0 bug in linux 6.5
226        %the_product = prod(values(used_in_monom).^monom_i(used_in_monom));
227        the_product = 1;
228        for j = 1:length(used_in_monom)
229            the_product = the_product*values(used_in_monom(j))^monom_i(used_in_monom(j));
230        end
231        values(lmi_variables(i)) = the_product;
232    end
233end
234sys = X.basis*[1;values(lmi_variables(:))];
235
236if X.typeflag==10
237    sys = eig(full(reshape(sys,X.dim(1),X.dim(2))));
238else
239    sys = full(reshape(sys,X.dim(1),X.dim(2)));
240end
Note: See TracBrowser for help on using the repository browser.