# HG changeset patch
# User Jaroslav Hajek <highegg@gmail.com>
# Date 1238069100 -3600
# Node ID 484756d558d606e3a19d7c18bae42f3a623da36a
# Parent  15774d99f4cce686eb825e82fe2d81c5a5de264a
add Array<T>::find

diff --git a/liboctave/Array.cc b/liboctave/Array.cc
--- a/liboctave/Array.cc
+++ b/liboctave/Array.cc
@@ -2487,6 +2487,69 @@
   return idx;
 }
 
+template <class T>
+Array<octave_idx_type> 
+Array<T>::find (octave_idx_type n, bool backward) const
+{
+  Array<octave_idx_type> retval;
+  const T *src = data ();
+  octave_idx_type nel = nelem ();
+  const T zero = T ();
+  if (n < 0)
+    {
+      // We want all elements, which means we'll almost surely need
+      // to resize. So count first, then allocate array of exact size.
+      octave_idx_type cnt = 0;
+      for (octave_idx_type i = 0; i < nel; i++)
+        cnt += src[i] != zero;
+
+      retval = Array<octave_idx_type> (cnt);
+      octave_idx_type *dest = retval.fortran_vec ();
+      for (octave_idx_type i = 0; i < nel; i++)
+        if (src[i] != zero) *dest++ = i;
+    }
+  else
+    {
+      // We want a fixed max number of elements, usually small. So be
+      // optimistic, alloc the array in advance, and then resize if
+      // needed.
+      retval = Array<octave_idx_type> (n);
+      if (backward)
+        {
+          // Do the search as a series of successive single-element searches.
+          octave_idx_type k = 0, l = nel - 1;
+          for (; k < n; k++)
+            {
+              for (;l >= 0 && src[l] == zero; l--) ;
+              if (l >= 0)
+                retval(k) = l--;
+              else
+                break;
+            }
+          if (k < n)
+            retval.resize (k);
+        }
+      else
+        {
+          // Do the search as a series of successive single-element searches.
+          octave_idx_type k = 0, l = 0;
+          for (; k < n; k++)
+            {
+              for (;l != nel && src[l] == zero; l++) ;
+              if (l != n)
+                retval(k) = l++;
+              else
+                break;
+            }
+          if (k < n)
+            retval.resize (k);
+        }
+    }
+
+  return retval;
+}
+
+
 #define INSTANTIATE_ARRAY_SORT(T) template class octave_sort<T>;
 
 #define NO_INSTANTIATE_ARRAY_SORT(T) \
@@ -2520,7 +2583,9 @@
 template <> Array<octave_idx_type>  \
 Array<T>::lookup (const Array<T>&, sortmode, bool, bool) const \
 { return Array<octave_idx_type> (); } \
-
+template <> Array<octave_idx_type> \
+Array<T>::find (octave_idx_type, bool) const\
+{ return Array<octave_idx_type> (); } \
 
 
 template <class T>
diff --git a/liboctave/Array.h b/liboctave/Array.h
--- a/liboctave/Array.h
+++ b/liboctave/Array.h
@@ -579,6 +579,10 @@
   Array<octave_idx_type> lookup (const Array<T>& values, sortmode mode = UNSORTED, 
                                  bool linf = false, bool rinf = false) const;
 
+  // Find indices of (at most n) nonzero elements. If n is specified, backward
+  // specifies search from backward.
+  Array<octave_idx_type> find (octave_idx_type n = -1, bool backward = false) const;
+
   Array<T> diag (octave_idx_type k = 0) const;
 
   template <class U, class F>
diff --git a/liboctave/ChangeLog b/liboctave/ChangeLog
--- a/liboctave/ChangeLog
+++ b/liboctave/ChangeLog
@@ -1,3 +1,8 @@
+2009-03-26  Jaroslav Hajek  <highegg@gmail.com>
+
+	* Array.cc (Array<T>::find): New method.
+	* Array.h: Declare it.
+
 2009-03-25  John W. Eaton  <jwe@octave.org>
 
 	* EIG.cc (EIG::init (const Matrix&, bool),