changeset 10484:6a9571b57745

rewrite sparse null assignment (part 1)
author Jaroslav Hajek <highegg@gmail.com>
date Fri, 02 Apr 2010 16:04:52 +0200
parents 7b5f706f3a83
children b4e14e628fc9
files liboctave/ChangeLog liboctave/Sparse.cc
diffstat 2 files changed, 78 insertions(+), 239 deletions(-) [+]
line wrap: on
line diff
--- a/liboctave/ChangeLog	Fri Apr 02 07:31:59 2010 +0200
+++ b/liboctave/ChangeLog	Fri Apr 02 16:04:52 2010 +0200
@@ -1,3 +1,8 @@
+2010-04-02  Jaroslav Hajek  <highegg@gmail.com>
+
+	* Sparse.cc (Sparse<T>::maybe_delete_elements): Rewrite. Optimize for
+	sparse column vectors.
+
 2010-04-01  Jaroslav Hajek  <highegg@gmail.com>
 
 	* mx-inlines.cc: Declare all loops as throw (). Ditto for
--- a/liboctave/Sparse.cc	Fri Apr 02 07:31:59 2010 +0200
+++ b/liboctave/Sparse.cc	Fri Apr 02 16:04:52 2010 +0200
@@ -1219,260 +1219,94 @@
 
 template <class T>
 void
-Sparse<T>::maybe_delete_elements (idx_vector& idx_arg)
+Sparse<T>::maybe_delete_elements (idx_vector& idx)
 {
+  Sparse<T> retval;
+
+  assert (ndims () == 2);
+
+  // FIXME: please don't fix the shadowed member warning yet because
+  // Sparse<T>::idx will eventually go away.
+
   octave_idx_type nr = dim1 ();
   octave_idx_type nc = dim2 ();
-  octave_idx_type orig_nr = nr;
-  octave_idx_type orig_nc = nc;
-  octave_idx_type orig_nnz = nnz ();
-
-  if (nr == 0 && nc == 0)
-    return;
-
-  octave_idx_type n;
-  if (nr == 1)
-    n = nc;
+  octave_idx_type nz = nnz ();
+
+  octave_idx_type nel = numel (); // Can throw.
+
+  const dim_vector idx_dims = idx.orig_dimensions ();
+
+  if (idx.extent (nel) > nel)
+    gripe_del_index_out_of_range (true, idx.extent (nel), nel);
   else if (nc == 1)
-    n = nr;
-  else
-    {
-      // Reshape to row vector for Matlab compatibility.
-
-      n = nr * nc;
-      nr = 1;
-      nc = n;
-    }
-
-  if (idx_arg.is_colon_equiv (n, 1))
     {
-      // Either A(:) = [] or A(idx) = [] with idx enumerating all
-      // elements, so we delete all elements and return [](0x0).  To
-      // preserve the orientation of the vector, you have to use
-      // A(idx,:) = [] (delete rows) or A(:,idx) (delete columns).
-
-      resize (0, 0);
-      return;
-    }
-
-  idx_arg.sort (true);
-
-  octave_idx_type num_to_delete = idx_arg.length (n);
-
-  if (num_to_delete != 0)
-    {
-      octave_idx_type new_n = n;
-      octave_idx_type new_nnz = orig_nnz;
-
-      octave_idx_type iidx = 0;
-      octave_idx_type kk = idx_arg.elem (iidx);
-
-      const Sparse<T> tmp (*this);
-
-      if (orig_nr == 1)
+      // Sparse column vector.
+      const Sparse<T> tmp = *this; // constant copy to prevent COW.
+
+      octave_idx_type lb, ub;
+
+      if (idx.is_cont_range (nel, lb, ub))
         {
-          for (octave_idx_type i = 0; i < n; i++)
-            {
-              octave_quit ();
-
-              if (i == kk)
-                {
-                  iidx++;
-                  kk = idx_arg.elem (iidx);
-                  new_n--;
-
-                  if (tmp.cidx(i) != tmp.cidx(i + 1))
-                    new_nnz--;
-
-                  if (iidx == num_to_delete)
-                    break;
-                }
-            }
-        }
-      else if (orig_nc == 1)
-        {
-          octave_idx_type next_ridx = -1;
-          octave_idx_type next_ridx_val = -1;
-          if (orig_nnz > 0)
-            {
-              next_ridx = 0;
-              next_ridx_val = tmp.ridx (0);
-            }
-
-          for (octave_idx_type i = 0; i < n; i++)
-            {
-              octave_quit ();
-
-              if (i == kk)
-                {
-                  iidx++;
-                  kk = idx_arg.elem (iidx);
-                  new_n--;
-
-                  while (i > next_ridx_val)
-                    {
-                      next_ridx++;
-                      if (next_ridx >= orig_nnz)
-                        {
-                          next_ridx = -1;
-                          next_ridx_val = n;
-                          break;
-                        }
-                      else
-                        next_ridx_val = tmp.ridx (next_ridx);
-                    }
-
-                  if (i == next_ridx_val)
-                    new_nnz--;
-
-                  if (iidx == num_to_delete)
-                    break;
-                }
-            }
+          // Special-case a contiguous range.
+          // Look-up indices first.
+          octave_idx_type li = lblookup (tmp.ridx (), nz, lb);
+          octave_idx_type ui = lblookup (tmp.ridx (), nz, ub);
+          // Copy data and adjust indices.
+          octave_idx_type nz_new = nz - (ui - li);
+          *this = Sparse<T> (nr - (ub - lb), 1, nz_new);
+          copy_or_memcpy (li, tmp.data (), data ());
+          copy_or_memcpy (li, tmp.ridx (), xridx ());
+          copy_or_memcpy (nz - ui, tmp.data () + ui, xdata () + li);
+          mx_inline_sub (nz - ui, xridx () + li, tmp.ridx () + ui, ub - lb);
+          xcidx(1) = nz_new;
         }
       else
         {
-          for (octave_idx_type i = 0; i < n; i++)
+          OCTAVE_LOCAL_BUFFER (octave_idx_type, ridx_new, nz);
+          OCTAVE_LOCAL_BUFFER (T, data_new, nz);
+          idx_vector sidx = idx.sorted (true);
+          const octave_idx_type *sj = sidx.raw ();
+          octave_idx_type sl = sidx.length (nel), nz_new = 0, j = 0;
+          for (octave_idx_type i = 0; i < nz; i++)
             {
-              octave_quit ();
-
-              if (i == kk)
+              octave_idx_type r = tmp.ridx(i);
+              for (;j < sl && sj[j] < r; j++) ;
+              if (j == sl || sj[j] > r)
                 {
-                  iidx++;
-                  kk = idx_arg.elem (iidx);
-                  new_n--;
-
-                  if (tmp.elem (i) != T ())
-                    new_nnz--;
-
-                  if (iidx == num_to_delete)
-                    break;
-                }
-            }
-        }
-
-      if (new_n > 0)
-        {
-          rep->count--;
-
-          if (nr == 1)
-            rep = new typename Sparse<T>::SparseRep (1, new_n, new_nnz);
-          else
-            rep = new typename Sparse<T>::SparseRep (new_n, 1, new_nnz);
-
-          octave_idx_type ii = 0;
-          octave_idx_type jj = 0;
-          iidx = 0;
-          kk = idx_arg.elem (iidx);
-
-          if (orig_nr == 1)
-            {
-              for (octave_idx_type i = 0; i < n; i++)
-                {
-                  octave_quit ();
-
-                  if (iidx < num_to_delete && i == kk)
-                    kk = idx_arg.elem (++iidx);
-                  else
-                    {
-                      if (tmp.cidx(i) != tmp.cidx(i + 1))
-                        {
-                          data(ii) = tmp.data (tmp.cidx(i));
-                          ridx(ii++) = jj;
-                        }
-                      jj++;
-                    }
+                  data_new[nz_new] = tmp.data(i);
+                  ridx_new[nz_new++] = r - j;
                 }
             }
-          else if (orig_nc == 1)
-            {
-              octave_idx_type next_ridx = -1;
-              octave_idx_type next_ridx_val = n;
-              if (orig_nnz > 0)
-                {
-                  next_ridx = 0;
-                  next_ridx_val = tmp.ridx (0);
-                }
-
-              for (octave_idx_type i = 0; i < n; i++)
-                {
-                  octave_quit ();
-
-                  if (iidx < num_to_delete && i == kk)
-                    kk = idx_arg.elem (++iidx);
-                  else
-                    {
-                      while (i > next_ridx_val)
-                        {
-                          next_ridx++;
-                          if (next_ridx >= orig_nnz)
-                            {
-                              next_ridx = -1;
-                              next_ridx_val = n;
-                              break;
-                            }
-                          else
-                            next_ridx_val = tmp.ridx (next_ridx);
-                        }
-
-                      if (i == next_ridx_val)
-                        {
-                          data(ii) = tmp.data(next_ridx);
-                          ridx(ii++) = jj;
-                        }
-                      jj++;
-                    }
-                }
-            }
-          else
-            {
-              for (octave_idx_type i = 0; i < n; i++)
-                {
-                  octave_quit ();
-
-                  if (iidx < num_to_delete && i == kk)
-                    kk = idx_arg.elem (++iidx);
-                  else
-                    {
-                      T el = tmp.elem (i);
-                      if (el != T ())
-                        {
-                          data(ii) = el;
-                          ridx(ii++) = jj;
-                        }
-                      jj++;
-                    }
-                }
-            }
-
-          dimensions.resize (2);
-
-          if (nr == 1)
-            {
-              ii = 0;
-              cidx(0) = 0;
-              for (octave_idx_type i = 0; i < new_n; i++)
-                {
-                  octave_quit ();
-                  if (ridx(ii) == i)
-                    ridx(ii++) = 0;
-                  cidx(i+1) = ii;
-                }
-
-              dimensions(0) = 1;
-              dimensions(1) = new_n;
-            }
-          else
-            {
-              cidx(0) = 0;
-              cidx(1) = new_nnz;
-              dimensions(0) = new_n;
-              dimensions(1) = 1;
-            }
+
+          *this = Sparse<T> (nr - sl, 1, nz_new);
+          copy_or_memcpy (nz_new, ridx_new, ridx ());
+          copy_or_memcpy (nz_new, data_new, xdata ());
+          xcidx(1) = nz_new;
+        }
+    }
+  else if (nr == 1)
+    {
+      // Sparse row vector.
+      octave_idx_type lb, ub;
+      if (idx.is_cont_range (nc, lb, ub))
+        {
+          const Sparse<T> tmp = *this;
+          octave_idx_type lbi = tmp.cidx(lb), ubi = tmp.cidx(ub), new_nz = nz - (ubi - lbi);
+          *this = Sparse<T> (1, nc - (ub - lb), new_nz);
+          copy_or_memcpy (lbi, tmp.data (), data ());
+          copy_or_memcpy (nz - ubi, tmp.data () + ubi, xdata () + lbi);
+          fill_or_memset (new_nz, static_cast<octave_idx_type> (0), ridx ());
+          copy_or_memcpy (lb, tmp.cidx () + 1, cidx () + 1);
+          mx_inline_sub (nc - ub, xcidx () + 1, tmp.cidx () + ub + 1, ubi - lbi);
         }
       else
-        (*current_liboctave_error_handler)
-          ("A(idx) = []: index out of range");
+        *this = index (idx.complement (nc));
+    }
+  else
+    {
+      *this = index (idx_vector::colon);
+      maybe_delete_elements (idx);
+      *this = transpose (); // We want a row vector.
     }
 }