Once more pinv / dgelss

[thombos]

How did you allocate your matrices A and B ? I see you are using 2D arrays. Remember that the matrix needs to be stored in a contiguous array and it needs to be column major in your case.

Thomas

[Julien Langou]

Wait. I did not run your code. But if A is the input matrix of DGELSS
(5-by-2), B is the input right-hand side (5-by-5), then in output of DGELSS,
you will find in the array A the right singular vector of A (stored rowwise); in
the array B, you will find in the first 2 rows of B the 2-by-5 pseudo inverse
and in the last 3 rows the norm of residuals for each right-hand side.
(You need to take the sqrt of the sum of the squares of the 3 elements in
each column to get the norm of the residual for each right-hand side.)

Where are you reading the pseudo-inverse?
Try to print the whole 5-by-5 matrix B and understand where the pseudo-
inverse is.

Julien.

[Julien Langou]

Oh yes. Good point from Thomas there. You should provide your data
in column major format so stored by column (as you have just mentioned).
Still I am not sure you are reading the pseudo-inverse correctly in the array
B. (Once more it will be given back to you by column and you just care
about the first two rows of B.)
-j

[thombos]

But I can see it from here (possibly tousands of miles away). :D

[thombos]

Try setting your NRHS parameter to be _r as well.

Thomas

What I mean to say is that some numbers in your result are correct, in particular the first column and one further one, indicating an error with your matrix format specifications to me, so the NRHS parameter is the number of columns of B, which should equal the number of rows in A.

[thombos]

Just on a hunch, the first time you call gelss, just hand it NULL instead of A and B.

Are you actually using NRHS in your call, because I'm seeing

dgelss_(&_r, &_c, &_c ....

which should be
dgelss_(&_r, &_c, &NRHS ....
or
dgelss_(&_r, &_c, &_r .....

[thombos]

Are you actually using NRHS in your call, because I'm seeing

dgelss_(&_r, &_c, &_c ....

which should be
dgelss_(&_r, &_c, &NRHS ....
or
dgelss_(&_r, &_c, &_r .....

[thombos]

Your result is BTW correct already. Its like Julie said, its all a matter of interpreting the output matrix correctly, check this:

B={0.510716 -0.229171 -0.802211 0.122618 0.0520118 0.0715534 0.397695 0.144067 -0.637568 -0.446299 0.770416 -0.575517 0.571603 0.245704 0.15571 -0.369687 0.775273 0.0850571 0.688318 -0.21754 -0.228787 0.525211 0.0416721 -0.210524 0.852375}

Thomas (Cambridge, MA)

If you read the documentation for dgelss you will see how to calculate the indices of the solution.