Hi c., thank you for your very interesting suggestions: it was exactly what I hoped to see. I tried the iterative solution, and although the result is not what I expected (the one I erroneously called c, was the x0 matrix, indeed), it seems it works.
Il giorno venerdì 31 maggio 2013, c. <
address@hidden> ha scritto:
>
> On 31 May 2013, at 00:25, andrea console <
address@hidden> wrote:
>
>> Thank you all.
>> In this example the matrices size is:
>> A: 1072452 X 323136 (sparse)
>> b: 1072452 X 1
>> c: 323136 X 1
> what is c?
>> This is the file (67MB):
https://docs.google.com/file/d/0B4nOH7xkcKiabkdSNktLT0xhYjg/edit?usp=sharing
>>
>> I didn't save it as blocks because, before computation, I have to remove some rows (and corresponding elements in b array), so the block structure is not preserved
>>
>
>
> The matrix is very large but very sparse and "very rectangular".
> As the original block structure is lost, the best I can suggest at this point is
>
> AA = A' * A;
> bb = A' * b;
>
> %% then either
> x = AA \ bb; %% this will take veeeeery long
> %% or
> x = pcg (AA, bb, [], []);
> %% this will be inaccurate but you can refine the solution by an additional step
> x = pcg (AA, bb, [], [], x); %% repeat ad libitum
>
> Alternatively you may want to use a singular value decomposition, unfortunately
> I just found a small bug in svds.m that creates a full size vector from the sparse matrix input
> at an intermediate step of the algorithm, which I fixed with this changeset:
>
>
http://hg.savannah.gnu.org/hgweb/octave/rev/bd5e1b3a1ef9
>
> so you'll also need to apply the same change to svds.m before you can use it:
>
> --- a/scripts/sparse/svds.m
> +++ b/scripts/sparse/svds.m
> @@ -136,8 +136,9 @@
>
> error ("svds: SIGMA must be a positive real value or the string 'L'");
> endif
>
> - [m, n] = size (A);
> - max_a = max (abs (A(:)));
> + [m, n] = size (A);
> + max_a = max (abs (nonzeros (A)));
> +
> if (max_a == 0)
> s = zeros (k, 1); # special case of zero matrix
> else
>
> c.
>
>
>
>
>
>
>
>
>