<br><br><div class="gmail_quote">On Mon, Apr 29, 2013 at 6:10 AM, Jordi Gutiérrez Hermoso <span dir="ltr"><<a href="mailto:jordigh@octave.org" target="_blank">jordigh@octave.org</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

On 29 April 2013 06:25, Miroslaw Kwasniak <<a href="mailto:Miroslaw.Kwasniak@pwr.wroc.pl">Miroslaw.Kwasniak@pwr.wroc.pl</a>> wrote:<br>
> it's something wrong whith sparse matrices A(n,n) when n is a multiple<br>
> of 65536=2^16.<br>
><br>
> Demonstration code ======================================<br>
><br>
> for i=1:3;<br>
>   for n=i*2^16+(-1:1);<br>
>     A=spdiags(ones(n,1),0,n,n);<br>
>     t=trace(A);<br>
>     printf("n=%8d trace=%8d %s\n",n,t,["ERR";"ok"]((t==n)+1,:));<br>
>   endfor;<br>
> endfor<br>
><br>
> Results ======================================<br>
><br>
> n=   65535 trace=   65535 ok<br>
> n=   65536 trace=       0 ERR<br>
> n=   65537 trace=   65537 ok<br>
> n=  131071 trace=  131071 ok<br>
> n=  131072 trace=       0 ERR<br>
> n=  131073 trace=  131073 ok<br>
> n=  196607 trace=  196607 ok<br>
> n=  196608 trace=       0 ERR<br>
> n=  196609 trace=  196609 ok<br>
<br>
Confirmed. The problem is that the numel function is limited to<br>
returning octave_idx_type, which ordinarily of size 2^32, and<br>
certainly is so for Debian. This makes sense, since you can only index<br>
that many elements in a matrix. You're hitting the indexing limit. To<br>
get 64-bit indexing, you would need to recompile all of Octave's<br>
Fortran dependencies with -fdefault-integer-8.<br>
<br>
I'm not sure exactly what the bug is here. For instance, you can't<br>
index your matrix A either, and this is checked for correctly:<br>
<br>
    A(end)<br>
<br>
Perhaps the best thing  to do would be to forbid creation of sparse<br>
matrices where numel(A) > std::numeric_limits<int>::max().<br>
<br>
Your matrix is simply too large to be indexed, and this breaks<br>
assumptions elsewhere in our code.<br>
<br>
- Jordi G. H.<br>
</blockquote></div><br>I'm confused - this is a diagonal sparse matrix so you should be able to trace() (or any other op)<div>up to n = 2^32, not n^2 = 2^32. The limit on sparse matrices should be number of non-zeros < 2^32</div>

<div><div><br></div>-- <br>Ed Meyer<br>
</div>