Re: sprand/sprandn fourth parameter implementation

[Felipe G. Nievinski]

Hi Eduardo.  I hope you don't mind me CC'ing the list at large.

On Thu, Mar 6, 2014 at 10:04 AM, Eduardo <address@hidden> wrote:

...

I think I came up with a solution using the clues given in the matlab documentation. This is the algorithm that I coded for the sprand(n,m,density,rc) function (added 4th parameter). The aproach is the following. I assume that I start with a matrix S and I want to get a sparse one as A=U*SV' (svd descomposition scheme). U and V will provide the randomness as they will be Jacobi random rotation matrices.

a) If rc is a vector with the singular values desired for the sparse matrix then:

     1) Create an nxm diagonal matrix with the singular values ordered in the 1st diagonal and then the rest 0s.
         rc = sort(rc, 'descend');
         S = sparse (diag (rc, n, m));

b) If rc is the reciprocal condition number then I do the following:
    
 1)vmin = rc * rand();  % Stablish the minimum single value with a value between 0 and rc
    vmax = vmin / rc;   % Calculate the maximun single value from the minimum and rc.Both vmin, and vmax are < 1
    for (i = 1:min(m, n))
      v(i) = rand() * (vmax - vmin) + vmin; %Generate min(m,n) random values between vmin and vmax
    endfor
    v(1) = vmax;
    v(end) = vmin;
    v = sort (v, 'descend');
    S = sparse (diag (v, n, m)); Then build the sparse diagonal matrix

 c) The last step is  performing random rotations to S using Jacobi rotations (V and U as rotation matrices). That will populate the matrix S with a few new terms in each iteration. Stop  when reached the density of elements desired.

I have done several tests and everything works fine, condition number/singular values are preserved and density is achieved with a maximum error 0.1%. I am concerned about the type and distribution of values I get. Because you end up with lots of negative values and values greater than 1 and lower than -1 as elements. Is that behavior the same as in the Matlab implementation? I havent got a copy to try the function. What do you think about the aproach followed?

Thanks in advance,

Eduardo (Edu159)

It seems that now you have a very good understanding of the problem.

I suggest the following next steps:

- write a few paragraphs describing briefly in plain English the algorithm implemented, and include this in the function documentation; * bonus points if you can cite a book or article!

- update the documentation to describe the new input arguments now accepted;

- try and write an automated test can could be ran without much manual intervention and interpretation; see:

<https://en.wikipedia.org/wiki/Automated_testing>

- send a message to the list with subject starting with "Matlab Run Request", whereby you ask folks to help you compare answers

- have you found a similar implementation in other mathematical FLOSS?

Keep up the good work.

-F.