Numerical derivative of matrix eqn

[spasmous]

A = [1 2; 3 4; 5 6];
b = [1;2;3];
db = [1e-2;1e-2;1e-2];

x = pinv(A)*b;
dx = pinv(A)*db;

From matrix calculus, I know that dx/db = pinv(A) however the
numerical calculation of dx*pinv(db) is:

-0.33333 -0.33333 -0.33333
0.33333 0.33333 0.33333

which is not pinv(A). What am doing wrong here :(

[Greg Heath]

clear all, clc

A = [1 2; 3 4; 5 6] ;
b = [1;2;3] ;
db = [1e-2;1e-2;1e-2] ;

pinvA = pinv(A)
x = pinvA*b
dx = pinvA*db
err = A*x-b

% pinvA =
% -1.3333 -0.3333 0.6667
% 1.0833 0.3333 -0.4167
% x =
% 0.0000
% 0.5000
% dx =
% -0.0100
% 0.0100
% err =
% 1.0e-015 *
% 0.2220
% 0.8882
% 0.8882
%

% >From matrix calculus, I know that dx/db = pinv(A)
% -----SNIP
%
% Incorrect. You know that

pinvA1 = x/b
pinvA2 = dx/db

% pinvA1 =
% 0 0 0.0000
% 0 0 0.1667
% pinvA2 =
% -1.0000 0 0
% 1.0000 0 0
%
% are minimizers of

norm(x-pinvA1*b) % = 0
norm(dx-pinvA2*db) % = 0

% However, neither is pinvA.

Hope this helps.

Greg

[Peter Spellucci]

In article <34d21c55-0c55-44cb...@p35g2000yqh.googlegroups.com>,

you assume that b*pinv(b) = id

hth
peter

[Greg Heath]

On Nov 27, 11:02 am, spellu...@fb04633.mathematik.tu-darmstadt.de
(Peter Spellucci) wrote:
> In article <34d21c55-0c55-44cb-b5fc-3c96d846b...@p35g2000yqh.googlegroups.com>,

>
>  spasmous <spasm...@gmail.com> writes:
> >A = [1 2; 3 4; 5 6];
> >b = [1;2;3];
> >db = [1e-2;1e-2;1e-2];
>
> >x = pinv(A)*b;
> >dx = pinv(A)*db;
>
> >From matrix calculus, I know that dx/db = pinv(A) however the
> >numerical calculation of dx*pinv(db) is:
>
> >  -0.33333  -0.33333  -0.33333
> >   0.33333   0.33333   0.33333
>
> >which is not pinv(A). What am doing wrong here :(
>
>  you assume that b*pinv(b) = id

Thanks for getting to the heart of the matter.

On Nov 27, 11:02 am, spellu...@fb04633.mathematik.tu-darmstadt.de
(Peter Spellucci) wrote:
> In article <34d21c55-0c55-44cb-b5fc-3c96d846b...@p35g2000yqh.googlegroups.com>,

>
> spasmous <spasm...@gmail.com> writes:
> >A = [1 2; 3 4; 5 6];
> >b = [1;2;3];
> >db = [1e-2;1e-2;1e-2];
>
> >x = pinv(A)*b;
> >dx = pinv(A)*db;
>
> >From matrix calculus, I know that dx/db = pinv(A) however the
> >numerical calculation of dx*pinv(db) is:
>
> > -0.33333 -0.33333 -0.33333
> > 0.33333 0.33333 0.33333
>
> >which is not pinv(A). What am doing wrong here :(
>
> you assume that b*pinv(b) = id

Thanks for getting to the heart of the matter.

Examples:

clear all, clc

A = [1 2; 3 4; 5 6] ;
b = [1;2;3] ;
db = [1e-2;1e-2;1e-2] ;

pinvA = pinv(A)
pinvb = pinv(b)
pinvdb = pinv(db)

ApinvA = A*pinvA
bpinvb = b*pinvb
dbpinvdb = db*pinvdb

pinvAA = pinvA*A
pinvbb = pinvb*b
pinvdbdb = pinvdb*db

% pinvA =
% -1.3333 -0.3333 0.6667
% 1.0833 0.3333 -0.4167

% pinvb =
% 0.0714 0.1429 0.2143
% pinvdb =
% 33.3333 33.3333 33.3333
%
% ApinvA =
% 0.8333 0.3333 -0.1667
% 0.3333 0.3333 0.3333
% -0.1667 0.3333 0.8333
% bpinvb =
% 0.0714 0.1429 0.2143
% 0.1429 0.2857 0.4286
% 0.2143 0.4286 0.6429
% dbpinvdb =
% 0.3333 0.3333 0.3333
% 0.3333 0.3333 0.3333
% 0.3333 0.3333 0.3333
%
% pinvAA =
% 1.0000 0.0000
% 0 1.0000
% pinvbb =
% 1.0000
% pinvdbdb =
% 1

Hope this helps.

Greg

[spasmous]

Yes, its clear I made an incorrect assumption about b*pinv(b). Thanks,
I should have noticed that myself! Thanks for your replies BTW, I'm
really happy this group is still alive to post questions. It looks
almost overrun with spam lately.