1 | function dfdx = shadowjacobian(f,x) |
---|
2 | % See SDPVAR/jacobian |
---|
3 | |
---|
4 | % Author Johan Löfberg |
---|
5 | % $Id: shadowjacobian.m,v 1.6 2005/05/10 15:04:52 joloef Exp $ |
---|
6 | |
---|
7 | if nargin==1 |
---|
8 | if isa(f,'sdpvar') |
---|
9 | x = recover(depends(f)); |
---|
10 | else |
---|
11 | x = 0; |
---|
12 | end |
---|
13 | else |
---|
14 | if length(getvariables(x))<length(x) |
---|
15 | error('x should be a vector of scalar independant variables'); |
---|
16 | end |
---|
17 | end |
---|
18 | |
---|
19 | if isa(f,'double') |
---|
20 | dfdx = zeros(size(f,1),length(x)); |
---|
21 | end |
---|
22 | |
---|
23 | [n,m]=size(f); |
---|
24 | |
---|
25 | if m>1 |
---|
26 | error('Jacobian only defined for column vectors.') |
---|
27 | end |
---|
28 | |
---|
29 | if n*m==1 |
---|
30 | dfdx = scalar_jacobian(f,x); |
---|
31 | % Argh, fix this (sorts inside scalar_jacobian |
---|
32 | for i = 1:length(x) |
---|
33 | var(i,1)=getvariables(x(i)); |
---|
34 | end |
---|
35 | [i,j]=sort(var); |
---|
36 | dfdx = dfdx(1,j); |
---|
37 | return |
---|
38 | |
---|
39 | else |
---|
40 | dfdx = []; |
---|
41 | AllVars = recover(unique([depends(f) getvariables(x)])); |
---|
42 | |
---|
43 | for i = 1:length(f) |
---|
44 | dfdx = [dfdx;scalar_jacobian(f(i),x,AllVars)]; |
---|
45 | end |
---|
46 | % Argh, fix this (sorts inside scalar_jacobian |
---|
47 | for i = 1:length(x) |
---|
48 | var(i,1)=getvariables(x(i)); |
---|
49 | end |
---|
50 | [i,j]=sort(var); |
---|
51 | dfdx = dfdx(:,j); |
---|
52 | end |
---|
53 | |
---|
54 | function [dfdx,dummy] = scalar_jacobian(f,x,AllVars) |
---|
55 | |
---|
56 | % Aaaaaah, horrible code. FIX! |
---|
57 | |
---|
58 | if isa(f,'double') |
---|
59 | dfdx = zeros(1,length(x)); |
---|
60 | return |
---|
61 | end |
---|
62 | |
---|
63 | if nargin==2 |
---|
64 | AllVars = recover(uniquestripped([depends(f) getvariables(x)])); |
---|
65 | end |
---|
66 | |
---|
67 | exponent_p = exponents(f,AllVars); |
---|
68 | |
---|
69 | coefficients = getbase(f); |
---|
70 | coefficients = coefficients(2:end); |
---|
71 | coefficients = coefficients(:); |
---|
72 | if nnz(exponent_p(1,:))==0 |
---|
73 | exponent_p=exponent_p(2:end,:); |
---|
74 | end |
---|
75 | |
---|
76 | x_variables = getvariables(x); |
---|
77 | AllVars_variables = getvariables(AllVars); |
---|
78 | %AllDeriv = []; |
---|
79 | AllDeriv2 = []; |
---|
80 | for k = 1:length(x) |
---|
81 | wrt = find(ismembc(AllVars_variables,x_variables(k))); |
---|
82 | deriv = exponent_p; |
---|
83 | deriv(:,wrt) = deriv(:,wrt)-1; |
---|
84 | % AllDeriv = [AllDeriv;deriv]; |
---|
85 | keep{k} = find(deriv(:,wrt)~=-1); |
---|
86 | AllDeriv2 = [AllDeriv2;deriv(keep{k},:)]; |
---|
87 | end |
---|
88 | |
---|
89 | if size(AllDeriv2,1)==0 |
---|
90 | dummy = 1; |
---|
91 | else |
---|
92 | dummy = recovermonoms(AllDeriv2,AllVars); |
---|
93 | end |
---|
94 | |
---|
95 | top = 1; |
---|
96 | dfdx=[]; |
---|
97 | for k = 1:length(x) |
---|
98 | wrt = find(ismembc(AllVars_variables,x_variables(k))); |
---|
99 | m = length(keep{k}); |
---|
100 | if m>0 |
---|
101 | poly = sum((coefficients(keep{k}(:)).*exponent_p(keep{k}(:),wrt)).*dummy(top:top+m-1),1); |
---|
102 | top = top + m; |
---|
103 | else |
---|
104 | poly = 0; |
---|
105 | end |
---|
106 | % %else |
---|
107 | % % poly = 0; |
---|
108 | % %end |
---|
109 | % tip = tip + LengthUsed{k};%size(exponent_p,1); |
---|
110 | % end |
---|
111 | %j = find(AllDeriv(top:top+size(exponent_p,1)-1,wrt)~=-1); |
---|
112 | %if ~isempty(j) |
---|
113 | % poly = sum(coefficients(j).*exponent_p(j,wrt).*dummy(top+j-1)); |
---|
114 | %else |
---|
115 | % poly = 0; |
---|
116 | %end |
---|
117 | %top = top + size(exponent_p,1); |
---|
118 | |
---|
119 | dfdx = [dfdx ; poly]; |
---|
120 | end |
---|
121 | dfdx = dfdx'; |
---|