Mercurial > octave-dspies
view liboctave/util/nz-iterators.h @ 19009:8d47ce2053f2 draft
Added safety checks to Array::xelem
There's no reason to have a method which never checks invariants, ever. Added
debugging checks to Array::xelem to help catch and debug out-of-bounds errors
and reference overlap
* configure.ac: Added configuration option for uniqueness-checking with xelem
* jit-typeinfo.cc (octave_jit_paren_scalar): Call const Array::xelem rather
than Array::xelem
* Array-util.h, Array-util.cc (check_out_of_range): Extract common pattern to
method
(check_index): Methods to check index is in-bounds
(compute_index): Added bool parameter check. does not check bounds when check
is false and BOUNDS_CHECKING is off
* Array.h, Array.cc (xelem): Use methods from Array-util.h to compute indices
(is_unique): Check if this is the only reference to data
* CmplxQR.cc, dbleQR.cc, fCmplxQR.cc, floatQR.cc
(form): Move second assignment to after the call to xelem
* lo-array-gripes.h, lo-array-gripes.cc (gripe_modifying_nonunique): Added
error message for when non-const xelem is called on non-unique array
author | David Spies <dnspies@gmail.com> |
---|---|
date | Mon, 14 Jul 2014 13:07:59 -0600 |
parents | 2e0613dadfee |
children | 3fb030666878 |
line wrap: on
line source
/* 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_nz_iterators_h) #define octave_nz_iterators_h 1 #include "interp-idx.h" #include "oct-inttypes.h" #include "Array.h" #include "DiagArray2.h" #include "PermMatrix.h" #include "Sparse.h" #include "direction.h" // This file contains generic column-major iterators over // the nonzero elements of any array or matrix. If you have a matrix mat // of type M, you can construct the proper iterator type using // M::iter_type iter(mat) and iter will iterate efficiently (forwards // or backwards) over the nonzero elements of mat. // // The parameter T is the element-type except for PermMatrix where // the element type is always bool // // Use a dir_handler to indicate which direction you intend to iterate. // (see step-dir.h). begin() resets the iterator to the beginning or // end of the matrix (for dir_handler<FORWARD> and <BACKWARD> respectively). // finished(dirc) indicates whether the iterators has finished traversing // the nonzero elements. step(dirc) steps from one element to the next. // // You can, for instance, use a for-loop as follows: // // typedef M::iter_type iter_t; // typedef M::element_type T; // // iter_t iter(mat); // dir_handler<1> dirc; // // for(iter.begin (dirc); !iter.finished (dirc); iter.step (dirc)) // { // octave_idx_type row = iter.row(); // octave_idx_type col = iter.col(); // double doub_index = iter.interp_index (); // T elem = iter.data(); // // ... Do something with these // } // // Note that array_iter for indexing over full matrices also includes // a iter.flat_index () method which returns an octave_idx_type. // // The other iterators to not have a flat_index() method because they // risk overflowing octave_idx_type. It is recommended you take care // to implement your function in a way that accounts for this problem. // // FIXME: I'd like to add in these // default no-parameter versions of // begin() and step() to each of the // classes. But the C++ compiler complains // because apparently I'm not allowed to overload // templated methods with non-templated ones. Any // ideas for work-arounds? // //#define INCLUDE_DEFAULT_STEPS \ // void begin (void) \ // { \ // dir_handler<FORWARD> dirc; \ // begin (dirc); \ // } \ // void step (void) \ // { \ // dir_handler<FORWARD> dirc; \ // step (dirc); \ // } \ // bool finished (void) const \ // { \ // dir_handler<FORWARD> dirc; \ // return finished (dirc); \ // } // A generic method for checking if some element of a matrix with // element type T is zero. template<typename T> bool is_zero (T t) { return t == static_cast<T> (0); } // An iterator over full arrays. When the number of dimensions exceeds // 2, calls to iter.col() may exceed mat.cols() up to mat.dims().numel(1) // // This mimics the behavior of the "find" method (both in Octave and Matlab) // on many-dimensional matrices. template<typename T> class array_iterator { private: const Array<T>& mat; //Actual total number of columns = mat.dims().numel(1) //can be different from length of row dimension const octave_idx_type totcols; const octave_idx_type numels; octave_idx_type coli; octave_idx_type rowj; octave_idx_type my_idx; template<direction dir> void step_once (dir_handler<dir> dirc) { my_idx += dir; rowj += dir; if (dirc.is_ended (rowj, mat.rows ())) { rowj = dirc.begin (mat.rows ()); coli += dir; } } template<direction dir> void move_to_nz (dir_handler<dir> dirc) { while (!finished (dirc) && is_zero (data ())) { step_once (dirc); } } public: array_iterator (const Array<T>& arg_mat) : mat (arg_mat), totcols (arg_mat.dims ().numel (1)), numels ( totcols * arg_mat.rows ()) { dir_handler<FORWARD> dirc; begin (dirc); } template<direction dir> void begin (dir_handler<dir> dirc) { coli = dirc.begin (totcols); rowj = dirc.begin (mat.rows ()); my_idx = dirc.begin (mat.numel ()); move_to_nz (dirc); } octave_idx_type col (void) const { return coli; } octave_idx_type row (void) const { return rowj; } double interp_idx (void) const { return to_interp_idx (my_idx); } octave_idx_type flat_idx (void) const { return my_idx; } T data (void) const { return mat.elem (my_idx); } template<direction dir> void step (dir_handler<dir> dirc) { step_once (dirc); move_to_nz (dirc); } template<direction dir> bool finished (dir_handler<dir> dirc) const { return dirc.is_ended (my_idx, numels); } }; template<typename T> class sparse_iterator { private: const Sparse<T>& mat; octave_idx_type coli; octave_idx_type my_idx; template<direction dir> void adjust_col (dir_handler<dir> dirc) { while (!finished (dirc) && dirc.is_ended (my_idx, mat.cidx (coli), mat.cidx (coli + 1))) coli += dir; } public: sparse_iterator (const Sparse<T>& arg_mat) : mat (arg_mat) { dir_handler<FORWARD> dirc; begin (dirc); } template<direction dir> void begin (dir_handler<dir> dirc) { coli = dirc.begin (mat.cols ()); my_idx = dirc.begin (mat.nnz ()); adjust_col (dirc); } double interp_idx (void) const { return to_interp_idx (row (), col (), mat.dims ()); } octave_idx_type col (void) const { return coli; } octave_idx_type row (void) const { return mat.ridx (my_idx); } T data (void) const { return mat.data (my_idx); } template<direction dir> void step (dir_handler<dir> dirc) { my_idx += dir; adjust_col (dirc); } template<direction dir> bool finished (dir_handler<dir> dirc) const { return dirc.is_ended (coli, mat.cols ()); } }; template<typename T> class diag_iterator { private: const DiagArray2<T>& mat; octave_idx_type my_idx; template <direction dir> void move_to_nz (dir_handler<dir> dirc) { while (!finished (dirc) && is_zero (data ())) { my_idx += dir; } } public: diag_iterator (const DiagArray2<T>& arg_mat) : mat (arg_mat) { dir_handler<FORWARD> dirc; begin (dirc); } template<direction dir> void begin (dir_handler<dir> dirc) { my_idx = dirc.begin (mat.diag_length ()); move_to_nz (dirc); } double interp_idx (void) const { return to_interp_idx (row (), col (), mat.dims ()); } octave_idx_type col (void) const { return my_idx; } octave_idx_type row (void) const { return my_idx; } T data (void) const { return mat.dgelem (my_idx); } template<direction dir> void step (dir_handler<dir> dirc) { my_idx += dir; move_to_nz (dirc); } template<direction dir> bool finished (dir_handler<dir> dirc) const { return dirc.is_ended (my_idx, mat.diag_length ()); } }; class perm_iterator { private: const PermMatrix& mat; octave_idx_type my_idx; public: perm_iterator (const PermMatrix& arg_mat) : mat (arg_mat) { dir_handler<FORWARD> dirc; begin (dirc); } template<direction dir> void begin (dir_handler<dir> dirc) { my_idx = dirc.begin (mat.cols ()); } octave_idx_type interp_idx (void) const { return to_interp_idx (row (), col (), mat.dims ()); } octave_idx_type col (void) const { return my_idx; } octave_idx_type row (void) const { return mat.perm_elem (my_idx); } bool data (void) const { return true; } template<direction dir> void step (dir_handler<dir>) { my_idx += dir; } template<direction dir> bool finished (dir_handler<dir> dirc) const { return dirc.is_ended (my_idx, mat.rows ()); } }; #endif