[37] | 1 | function [base,v] = coefficients(p,x) |
---|
| 2 | %COEFFICIENTS Extract coefficients and monomials from polynomials |
---|
| 3 | % |
---|
| 4 | % [c,v] = COEFFICIENTS(p,x) extracts the coefficents |
---|
| 5 | % of a scalar polynomial p(x) = c'*v(x) |
---|
| 6 | % |
---|
| 7 | % c = COEFFICIENTS(p,x) extracts the all coefficents |
---|
| 8 | % of a matrix polynomial. |
---|
| 9 | |
---|
| 10 | % |
---|
| 11 | % INPUT |
---|
| 12 | % p : SDPVAR object |
---|
| 13 | % x : SDPVAR object |
---|
| 14 | % |
---|
| 15 | % OUTPUT |
---|
| 16 | % c : SDPVAR object |
---|
| 17 | % v : SDPVAR object |
---|
| 18 | % |
---|
| 19 | % EXAMPLE |
---|
| 20 | % sdpvar x y s t |
---|
| 21 | % p = x^2+x*y*(s+t)+s^2+t^2; % define p(x,y), parameterized with s and t |
---|
| 22 | % [c,v] = coefficients(p,[x y]); |
---|
| 23 | % sdisplay([c v]) |
---|
| 24 | % |
---|
| 25 | % See also SDPVAR |
---|
| 26 | |
---|
| 27 | % Author Johan Löfberg |
---|
| 28 | % $Id: coefficients.m,v 1.14 2006/08/30 14:11:14 joloef Exp $ |
---|
| 29 | |
---|
| 30 | if isa(p,'double') |
---|
| 31 | base = p(:); |
---|
| 32 | v = 1; |
---|
| 33 | return |
---|
| 34 | end |
---|
| 35 | |
---|
| 36 | if nargout>1 & (max(size(p))>1) |
---|
| 37 | error('For matrix inputs, only the coefficients can be returned. Request feature if you need this...'); |
---|
| 38 | end |
---|
| 39 | |
---|
| 40 | if nargin==1 |
---|
| 41 | allvar = depends(p); |
---|
| 42 | xvar = allvar; |
---|
| 43 | x = recover(xvar); |
---|
| 44 | else |
---|
| 45 | xvar = intersect(depends(x),depends(p)); |
---|
| 46 | end |
---|
| 47 | |
---|
| 48 | % Try to debug this! |
---|
| 49 | p = p(:); |
---|
| 50 | base = []; |
---|
| 51 | for i = 1:length(p) |
---|
| 52 | allvar = depends(p(i)); |
---|
| 53 | t = setdiff(allvar,xvar); |
---|
| 54 | if 1%isa(p(i),'sdpvar') |
---|
| 55 | [exponent_p,p_base] = getexponentbase(p(i),recover(depends(p(i)))); |
---|
| 56 | ParametricIndicies = find(ismember(allvar,t)); |
---|
| 57 | % FIX : don't define it here, wait until sparser below. Speed!! |
---|
| 58 | tempbase = parameterizedbase(p(i),[],recover(t),ParametricIndicies,exponent_p,p_base); |
---|
| 59 | [i,j,k] = unique(full(exponent_p(:,find(~ismember(allvar,t)))),'rows'); |
---|
| 60 | V = sparse(1:length(k),k,1,length(tempbase),max(k))'; |
---|
| 61 | base = [base;V*tempbase]; |
---|
| 62 | if nargout == 2 |
---|
| 63 | keepthese = j(1:max(k)); |
---|
| 64 | v = recovermonoms(exponent_p(keepthese,find(~ismember(allvar,t))),recover(xvar)); |
---|
| 65 | end |
---|
| 66 | end |
---|
| 67 | if isa(p,'ncvar') |
---|
| 68 | [exponent_p,ordered_list] = exponents(p,recover(depends(p(i)))); |
---|
| 69 | ParametricIndicies = find(ismember(allvar,t)); |
---|
| 70 | |
---|
| 71 | for j = 1:length(ParametricIndicies) |
---|
| 72 | a = find(ordered_list(:,1) == ParametricIndicies(j)) |
---|
| 73 | b = []; |
---|
| 74 | for k = 1:length(a) |
---|
| 75 | b = [b ordered_list(a(k),2:end)] |
---|
| 76 | end |
---|
| 77 | b = b(find(b)); |
---|
| 78 | basetemp = []; |
---|
| 79 | for k = 1:length(b) |
---|
| 80 | basetemp = [basetemp ncvar(struct(recover(t((k)))))]; |
---|
| 81 | end |
---|
| 82 | base = [base;sum(basetemp)]; |
---|
| 83 | end |
---|
| 84 | end |
---|
| 85 | end |
---|
| 86 | |
---|
| 87 | |
---|
| 88 | function p_base_parametric = parameterizedbase(p,z, params,ParametricIndicies,exponent_p,p_base) |
---|
| 89 | |
---|
| 90 | % Check for linear parameterization |
---|
| 91 | parametric_basis = exponent_p(:,ParametricIndicies); |
---|
| 92 | if all(sum(parametric_basis,2)==0) |
---|
| 93 | p_base_parametric = full(p_base(:)); |
---|
| 94 | return |
---|
| 95 | end |
---|
| 96 | if all(sum(parametric_basis,2)<=1) |
---|
| 97 | p_base_parametric = full(p_base(:)); |
---|
| 98 | n = length(p_base_parametric); |
---|
| 99 | ii = []; |
---|
| 100 | vars = []; |
---|
| 101 | js = sum(parametric_basis,1); |
---|
| 102 | for i = 1:size(parametric_basis,2) |
---|
| 103 | if js(i) |
---|
| 104 | j = find(parametric_basis(:,i)); |
---|
| 105 | ii = [ii j(:)']; |
---|
| 106 | vars = [vars repmat(i,1,js(i))]; |
---|
| 107 | end |
---|
| 108 | end |
---|
| 109 | k = setdiff1D(1:n,ii); |
---|
| 110 | if isempty(k) |
---|
| 111 | p_base_parametric = p_base_parametric.*sparse(ii,repmat(1,1,n),params(vars)); |
---|
| 112 | else |
---|
| 113 | pp = params(vars); % Must do this, bug in ML 6.1 (x=sparse(1);x([1 1]) gives different result in 6.1 and 7.0!) |
---|
| 114 | p_base_parametric = p_base_parametric.*sparse([ii k(:)'],repmat(1,1,n),[pp(:)' ones(1,1,length(k))]); |
---|
| 115 | end |
---|
| 116 | else |
---|
| 117 | % Bummer, nonlinear parameterization sucks... |
---|
| 118 | for i = 1:length(p_base) |
---|
| 119 | j = find(exponent_p(i,ParametricIndicies)); |
---|
| 120 | if ~isempty(j) |
---|
| 121 | temp = p_base(i); |
---|
| 122 | for k = 1:length(j) |
---|
| 123 | if exponent_p(i,ParametricIndicies(j(k)))==1 |
---|
| 124 | temp = temp*params(j(k)); |
---|
| 125 | else |
---|
| 126 | temp = temp*params(j(k))^exponent_p(i,ParametricIndicies(j(k))); |
---|
| 127 | end |
---|
| 128 | end |
---|
| 129 | xx{i} = temp; |
---|
| 130 | else |
---|
| 131 | xx{i} = p_base(i); |
---|
| 132 | end |
---|
| 133 | end |
---|
| 134 | p_base_parametric = stackcell(sdpvar(1,1),xx)'; |
---|
| 135 | end |
---|
| 136 | |
---|
| 137 | |
---|
| 138 | |
---|
| 139 | % |
---|
| 140 | % |
---|
| 141 | % |
---|
| 142 | % |
---|
| 143 | % |
---|
| 144 | % function [base,v] = coefficients(p,x) |
---|
| 145 | % %COEFFICIENTS Extract coefficients and monomials from polynomials |
---|
| 146 | % % |
---|
| 147 | % % [c,v] = COEFFICIENTS(p,x) extracts the coefficents |
---|
| 148 | % % of a polynomial p(x) = c'*v(x) |
---|
| 149 | % % |
---|
| 150 | % % INPUT |
---|
| 151 | % % p : SDPVAR object |
---|
| 152 | % % x : SDPVAR object |
---|
| 153 | % % |
---|
| 154 | % % OUTPUT |
---|
| 155 | % % c : SDPVAR object |
---|
| 156 | % % v : SDPVAR object |
---|
| 157 | % % |
---|
| 158 | % % EXAMPLE |
---|
| 159 | % % sdpvar x y s t |
---|
| 160 | % % p = x^2+x*y*(s+t)+s^2+t^2; % define p(x,y), parameterized with s and t |
---|
| 161 | % % [c,v] = coefficients(p,[x y]); |
---|
| 162 | % % sdisplay([c v]) |
---|
| 163 | % % |
---|
| 164 | % % See also SDPVAR |
---|
| 165 | % |
---|
| 166 | % % Author Johan Löfberg |
---|
| 167 | % % $Id: coefficients.m,v 1.14 2006/08/30 14:11:14 joloef Exp $ |
---|
| 168 | % |
---|
| 169 | % |
---|
| 170 | % if length(p) > 1%size(p,2) > 1 |
---|
| 171 | % error('Coefficents can only be applied to column vectors'); |
---|
| 172 | % end |
---|
| 173 | % |
---|
| 174 | % allvar = depends(p); |
---|
| 175 | % if nargin==1 |
---|
| 176 | % xvar = allvar; |
---|
| 177 | % x = recover(xvar); |
---|
| 178 | % else |
---|
| 179 | % xvar = depends(x); |
---|
| 180 | % end |
---|
| 181 | % |
---|
| 182 | % pvar = recover(depends(p)); |
---|
| 183 | % |
---|
| 184 | % base = []; |
---|
| 185 | % for i = 1:length(p) |
---|
| 186 | % [bi{i},vi{i}] = coefficientsi(p(i),xvar,pvar,allvar); |
---|
| 187 | % end |
---|
| 188 | % |
---|
| 189 | % % Fix the lengths of the basis to use same basis for all elements |
---|
| 190 | % if length(bi)>1 |
---|
| 191 | % allvars = []; |
---|
| 192 | % for i = 1:length(bi) |
---|
| 193 | % bivar{i} = getvariables(vi{i}); |
---|
| 194 | % if isequal(vi{i}(1),1) |
---|
| 195 | % bivar{i} = [0 bivar{i}]; |
---|
| 196 | % end |
---|
| 197 | % allvars = unique([allvars bivar{i}]); |
---|
| 198 | % end |
---|
| 199 | % v = recover(allvars); |
---|
| 200 | % c = zeros(length(p),length(allvars))'; |
---|
| 201 | % ci = []; |
---|
| 202 | % cj = []; |
---|
| 203 | % cv = []; |
---|
| 204 | % for i = 1:length(bi) |
---|
| 205 | % index = find(ismember(allvars,bivar{i})); |
---|
| 206 | % ci = [ci index]; |
---|
| 207 | % cj = [cj repmat(i,1,length(index))]; |
---|
| 208 | % cv = [cv bi{i}']; |
---|
| 209 | % end |
---|
| 210 | % base = sparse(ci,cj,cv); |
---|
| 211 | % else |
---|
| 212 | % base = bi{1}; |
---|
| 213 | % v = vi{1}; |
---|
| 214 | % end |
---|
| 215 | % |
---|
| 216 | % |
---|
| 217 | % function [base,v] = coefficientsi(p,xvar,pvar,allvar) |
---|
| 218 | % |
---|
| 219 | % % Try to debug this! |
---|
| 220 | % t = setdiff(allvar,xvar); |
---|
| 221 | % [exponent_p,p_base] = getexponentbase(p,pvar); |
---|
| 222 | % ParametricIndicies = find(ismember(allvar,t)); |
---|
| 223 | % % FIX : don't define it here, wait until sparser below. Speed!! |
---|
| 224 | % tempbase = parameterizedbase(p,[],recover(t),ParametricIndicies,exponent_p,p_base); |
---|
| 225 | % [i,j,k] = unique(full(exponent_p(:,find(~ismember(allvar,t)))),'rows'); |
---|
| 226 | % %V = sparse(max(k),length(tempbase)); |
---|
| 227 | % %for i = 1:max(k) |
---|
| 228 | % % V(i,find(k==i)) = 1; |
---|
| 229 | % %end |
---|
| 230 | % V = sparse(1:length(k),k,1,length(tempbase),max(k))'; |
---|
| 231 | % base = V*tempbase; |
---|
| 232 | % if nargout == 2 |
---|
| 233 | % keepthese = j(1:max(k)); |
---|
| 234 | % v = recovermonoms(exponent_p(keepthese,find(~ismember(allvar,t))),recover(xvar)); |
---|
| 235 | % end |
---|
| 236 | % |
---|
| 237 | % |
---|
| 238 | % function p_base_parametric = parameterizedbase(p,z, params,ParametricIndicies,exponent_p,p_base) |
---|
| 239 | % |
---|
| 240 | % % Check for linear parameterization |
---|
| 241 | % parametric_basis = exponent_p(:,ParametricIndicies); |
---|
| 242 | % if all(sum(parametric_basis,2)==0) |
---|
| 243 | % p_base_parametric = full(p_base(:)); |
---|
| 244 | % return |
---|
| 245 | % end |
---|
| 246 | % if all(sum(parametric_basis,2)<=1) |
---|
| 247 | % p_base_parametric = full(p_base(:)); |
---|
| 248 | % n = length(p_base_parametric); |
---|
| 249 | % ii = []; |
---|
| 250 | % vars = []; |
---|
| 251 | % js = sum(parametric_basis,1); |
---|
| 252 | % for i = 1:size(parametric_basis,2) |
---|
| 253 | % if js(i) |
---|
| 254 | % j = find(parametric_basis(:,i)); |
---|
| 255 | % ii = [ii j(:)']; |
---|
| 256 | % vars = [vars repmat(i,1,js(i))]; |
---|
| 257 | % end |
---|
| 258 | % end |
---|
| 259 | % k = setdiff1D(1:n,ii); |
---|
| 260 | % if isempty(k) |
---|
| 261 | % p_base_parametric = p_base_parametric.*sparse(ii,repmat(1,1,n),params(vars)); |
---|
| 262 | % else |
---|
| 263 | % pp = params(vars); % Must do this, bug in ML 6.1 (x=sparse(1);x([1 1]) gives different result in 6.1 and 7.0!) |
---|
| 264 | % p_base_parametric = p_base_parametric.*sparse([ii k(:)'],repmat(1,1,n),[pp(:)' ones(1,1,length(k))]); |
---|
| 265 | % end |
---|
| 266 | % else |
---|
| 267 | % % Bummer, nonlinear parameterization sucks... |
---|
| 268 | % for i = 1:length(p_base) |
---|
| 269 | % j = find(exponent_p(i,ParametricIndicies)); |
---|
| 270 | % if ~isempty(j) |
---|
| 271 | % temp = p_base(i); |
---|
| 272 | % for k = 1:length(j) |
---|
| 273 | % if exponent_p(i,ParametricIndicies(j(k)))==1 |
---|
| 274 | % temp = temp*params(j(k)); |
---|
| 275 | % else |
---|
| 276 | % temp = temp*params(j(k))^exponent_p(i,ParametricIndicies(j(k))); |
---|
| 277 | % end |
---|
| 278 | % end |
---|
| 279 | % xx{i} = temp; |
---|
| 280 | % else |
---|
| 281 | % xx{i} = p_base(i); |
---|
| 282 | % end |
---|
| 283 | % end |
---|
| 284 | % p_base_parametric = stackcell(sdpvar(1,1),xx)'; |
---|
| 285 | % end |
---|