Mercurial > octave-dspies
diff liboctave/util/find.h @ 19006:2e0613dadfee draft
All calls to "find" use the same generic implementation (bug #42408, 42421)
* find.cc: Rewrite.
Move generic "find" logic to find.h
(Ffind) : Changed calls to find_nonzero_elem_idx to find_templated
Added unit test for bug #42421
* Array.cc (and .h) (Array::find): Deleted function. Replaced with find::find(Array)
from find.h
* Array.h: Added typedef for array_iterator (in nz-iterators.h) as
Array::iter_type
* DiagArray2.h: Added typedef for diag_iterator (in nz-iterators.h) as
DiagArray2::iter_type
* PermMatrix.h: Added typedef for perm_iterator (in nz-iterators.h) as
PermMatrix::iter_type
Also added typedef for bool as PermMatrix::element_type
(not octave_idx_type)
Added an nnz() function (which is an alias for perm_length) and a
perm_elem(i) function for retrieving the ith element of the permutation
* Sparse.h: Added typedef for sparse_iterator (in nz-iterators.h) as
Sparse::iter_type
Added a short comment documenting the the argument to the numel
function
* idx-vector.cc (idx_vector::idx_mask_rep::as_array): Changed Array.find to
find::find(Array) (in find.h)
* (new file) find.h
* (new file) interp-idx.h: Simple methods for converting between interpreter
index type and internal octave_idx_type/row-col pair
* (new file) min-with-nnz.h: Fast methods for taking an arbitrary matrix M and
an octave_idx_type n and finding min(M.nnz(), n)
* (new file) nz-iterators.h: Iterators for traversing (in column-major order)
the nonzero elements of any array or matrix backwards or forwards
* (new file) direction.h: Generic methods for simplifying code has to deal with
a "backwards or forwards" template argument
* build-sparse-tests.sh: Removed 5-return-value calls to "find" in unit-tests;
Admittedly this commit breaks this "feature" which was undocumented and only
partially supported to begin with (ie never worked for full matrices,
permutation matrices, or diagonal matrices)
author | David Spies <dnspies@gmail.com> |
---|---|
date | Tue, 17 Jun 2014 16:41:11 -0600 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/liboctave/util/find.h Tue Jun 17 16:41:11 2014 -0600 @@ -0,0 +1,144 @@ +/* + +Copyright (C) 2014 David Spies + +This file is part of Octave. + +Octave is free software; you can redistribute it and/or modify it +under the terms of the GNU General Public License as published by the +Free Software Foundation; either version 3 of the License, or (at your +option) any later version. + +Octave is distributed in the hope that it will be useful, but WITHOUT +ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +for more details. + +You should have received a copy of the GNU General Public License +along with Octave; see the file COPYING. If not, see +<http://www.gnu.org/licenses/>. + +*/ +#if !defined (octave_find_h) +#define octave_find_h 1 + +#include "Array.h" + +#include "interp-idx.h" +#include "min-with-nnz.h" +#include "nz-iterators.h" +#include "direction.h" + +namespace find +{ + // This struct mostly matches the signature of ffind_result (see find.cc) + // except without a get_list function. Internal calls to "find" from + // liboctave should call the "find" method which in turn constructs one + // of these as a value to be filled. + // It should not be used anywhere else + + struct find_result + { + Array<octave_idx_type> res; + template<typename ELT_T> + void + add (octave_idx_type place, const array_iterator<ELT_T>& iter) + { + res.xelem (place) = iter.flat_idx (); + } + + find_result (void) { } + find_result (const dim_vector& nnz) : res (nnz) { } + }; + + // This is the generic "find" method for all matrix types. All calls to + // "find" get routed here one way or another as this method contains logic + // for determining the proper output dimensions and handling forward/reversed, + // n_to_find etc. + // It would be unwise to try and duplicate this logic elsewhere as it's fairly + // nuanced and unintuitive (but Matlab compatible). + + template<typename R, direction dir, typename MT> + R + find_to_R (dir_handler<dir> dirc, const MT& v, octave_idx_type n_to_find) + { + typedef typename MT::iter_type iter_t; + + octave_idx_type numres; + if (n_to_find == -1) + numres = v.nnz (); + else + numres = min_with_nnz (v, n_to_find); + + const dim_vector& dv = v.dims (); + octave_idx_type col_len = dv.elem (0); + octave_idx_type cols = dv.numel (1); + + // Fixup return dimensions, for Matlab compatibility. + // find (zeros (0,0)) -> zeros (0,0) + // find (zeros (1,0)) -> zeros (1,0) + // find (zeros (0,1)) -> zeros (0,1) + // find (zeros (0,X)) -> zeros (0,1) + // find (zeros (1,1)) -> zeros (0,0) !!!! WHY? + // find (zeros (0,1,0)) -> zeros (0,0) + // find (zeros (0,1,0,1)) -> zeros (0,0) etc + + dim_vector res_dims; + + if ((cols == 1 && col_len == 1 && numres == 0) + || (col_len == 0 && cols == 0)) + res_dims = dim_vector (0, 0); + else if (col_len == 1) + res_dims = dim_vector (1, numres); + else + res_dims = dim_vector (numres, 1); + + iter_t iterator (v); + R resvec (res_dims); + octave_idx_type count = dirc.begin (numres); + + for (iterator.begin (dirc); + !iterator.finished (dirc) && !dirc.is_ended (count, numres); + iterator.step (dirc), count += dir) + { + resvec.add (count, iterator); + } + + return resvec; + } + + // The method call used for internal liboctave calls to find. It constructs + // an find_result and fills that rather than a res_type_of (see "find.cc"). + // Note that this method only supports full arrays because octave_idx_type + // will generally overflow if used to represent the nonzero indices of sparse + // arrays. Instead consider using the sparse-iterator class from + // nz-iterators.h to avoid this problem + + template<typename ELT_T> + Array<octave_idx_type> + find (const Array<ELT_T>& v, octave_idx_type n_to_find = -1, direction dir = + FORWARD) + { + find_result resvec; + + switch (dir) + { + case FORWARD: + dir_handler<FORWARD> dirc_forward; + resvec = find_to_R<find_result> (dirc_forward, v, n_to_find); + break; + case BACKWARD: + dir_handler<BACKWARD> dirc_back; + resvec = find_to_R<find_result> (dirc_back, v, n_to_find); + break; + default: + liboctave_fatal( + "find: unknown direction: %d; must be FORWARD or BACKWARD", + dir); + } + + return resvec.res; + } +} + +#endif