GNU Octave  3.8.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ov-bool-sparse.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 2004-2013 David Bateman
4 Copyright (C) 1998-2004 Andy Adler
5 
6 This file is part of Octave.
7 
8 Octave is free software; you can redistribute it and/or modify it
9 under the terms of the GNU General Public License as published by the
10 Free Software Foundation; either version 3 of the License, or (at your
11 option) any later version.
12 
13 Octave is distributed in the hope that it will be useful, but WITHOUT
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
16 for more details.
17 
18 You should have received a copy of the GNU General Public License
19 along with Octave; see the file COPYING. If not, see
20 <http://www.gnu.org/licenses/>.
21 
22 */
23 
24 #ifdef HAVE_CONFIG_H
25 #include <config.h>
26 #endif
27 
28 #include <iostream>
29 #include <limits>
30 #include <vector>
31 
32 #include "dim-vector.h"
33 
34 #include "mxarray.h"
35 #include "ov-base.h"
36 #include "ov-scalar.h"
37 #include "ov-bool.h"
38 #include "ov-bool-mat.h"
39 #include "gripes.h"
40 #include "ops.h"
41 #include "oct-locbuf.h"
42 
43 #include "ov-re-sparse.h"
44 #include "ov-cx-sparse.h"
45 #include "ov-bool-sparse.h"
46 
47 #include "ov-base-sparse.h"
48 #include "ov-base-sparse.cc"
49 
51 
53 
55  "sparse bool matrix", "logical");
56 
57 static octave_base_value *
59 {
61 
62  return
63  new octave_sparse_matrix (SparseMatrix (v.sparse_bool_matrix_value ()));
64 }
65 
68 {
71 }
72 
75 {
76  octave_base_value *retval = 0;
77 
79  {
80  // Don't use numel, since it can overflow for very large matrices
81  // Note that for the second test, this means it becomes approximative
82  // since it involves a cast to double to avoid issues of overflow
83  if (matrix.rows () == 1 && matrix.cols () == 1)
84  {
85  // Const copy of the matrix, so the right version of () operator used
86  const SparseBoolMatrix tmp (matrix);
87 
88  retval = new octave_bool (tmp (0));
89  }
90  else if (matrix.cols () > 0 && matrix.rows () > 0
91  && (double (matrix.byte_size ()) > double (matrix.rows ())
92  * double (matrix.cols ()) * sizeof (bool)))
93  retval = new octave_bool_matrix (matrix.matrix_value ());
94  }
95 
96  return retval;
97 }
98 
99 double
101 {
102  double retval = lo_ieee_nan_value ();
103 
104  if (numel () > 0)
105  {
106  if (numel () > 1)
107  gripe_implicit_conversion ("Octave:array-to-scalar",
108  "bool sparse matrix", "real scalar");
109 
110  retval = matrix (0, 0);
111  }
112  else
113  gripe_invalid_conversion ("bool sparse matrix", "real scalar");
114 
115  return retval;
116 }
117 
118 Complex
120 {
121  double tmp = lo_ieee_nan_value ();
122 
123  Complex retval (tmp, tmp);
124 
125  if (rows () > 0 && columns () > 0)
126  {
127  if (numel () > 1)
128  gripe_implicit_conversion ("Octave:array-to-scalar",
129  "bool sparse matrix", "complex scalar");
130 
131  retval = matrix (0, 0);
132  }
133  else
134  gripe_invalid_conversion ("bool sparse matrix", "complex scalar");
135 
136  return retval;
137 }
138 
141  char type) const
142 {
144  return tmp.convert_to_str (pad, force, type);
145 }
146 
147 // FIXME: These are inefficient ways of creating full matrices
148 
149 Matrix
151 {
152  return Matrix (matrix.matrix_value ());
153 }
154 
157 {
158  return ComplexMatrix (matrix.matrix_value ());
159 }
160 
163 {
165 }
166 
167 NDArray
169 {
170  return NDArray (Matrix (matrix.matrix_value ()));
171 }
172 
175 {
176  charNDArray retval (dims (), 0);
177  octave_idx_type nc = matrix.cols ();
178  octave_idx_type nr = matrix.rows ();
179 
180  for (octave_idx_type j = 0; j < nc; j++)
181  for (octave_idx_type i = matrix.cidx (j); i < matrix.cidx (j+1); i++)
182  retval(matrix.ridx (i) + nr * j) = static_cast<char>(matrix.data (i));
183 
184  return retval;
185 }
186 
189 {
190  return matrix.matrix_value ();
191 }
192 
195 {
196  return boolNDArray (matrix.matrix_value ());
197 }
198 
199 
202 {
203  return SparseMatrix (this->matrix);
204 }
205 
208 {
209  return SparseComplexMatrix (this->matrix);
210 }
211 
212 bool
213 octave_sparse_bool_matrix::save_binary (std::ostream& os, bool&)
214 {
215  dim_vector d = this->dims ();
216  if (d.length () < 1)
217  return false;
218 
219  // Ensure that additional memory is deallocated
221 
222  int nr = d(0);
223  int nc = d(1);
224  int nz = nnz ();
225 
226  int32_t itmp;
227  // Use negative value for ndims to be consistent with other formats
228  itmp= -2;
229  os.write (reinterpret_cast<char *> (&itmp), 4);
230 
231  itmp= nr;
232  os.write (reinterpret_cast<char *> (&itmp), 4);
233 
234  itmp= nc;
235  os.write (reinterpret_cast<char *> (&itmp), 4);
236 
237  itmp= nz;
238  os.write (reinterpret_cast<char *> (&itmp), 4);
239 
240  // add one to the printed indices to go from
241  // zero-based to one-based arrays
242  for (int i = 0; i < nc+1; i++)
243  {
244  octave_quit ();
245  itmp = matrix.cidx (i);
246  os.write (reinterpret_cast<char *> (&itmp), 4);
247  }
248 
249  for (int i = 0; i < nz; i++)
250  {
251  octave_quit ();
252  itmp = matrix.ridx (i);
253  os.write (reinterpret_cast<char *> (&itmp), 4);
254  }
255 
256  OCTAVE_LOCAL_BUFFER (char, htmp, nz);
257 
258  for (int i = 0; i < nz; i++)
259  htmp[i] = (matrix.data (i) ? 1 : 0);
260 
261  os.write (htmp, nz);
262 
263  return true;
264 }
265 
266 bool
267 octave_sparse_bool_matrix::load_binary (std::istream& is, bool swap,
268  oct_mach_info::float_format /* fmt */)
269 {
270  int32_t nz, nc, nr, tmp;
271  if (! is.read (reinterpret_cast<char *> (&tmp), 4))
272  return false;
273 
274  if (swap)
275  swap_bytes<4> (&tmp);
276 
277  if (tmp != -2)
278  {
279  error ("load: only 2-D sparse matrices are supported");
280  return false;
281  }
282 
283  if (! is.read (reinterpret_cast<char *> (&nr), 4))
284  return false;
285  if (! is.read (reinterpret_cast<char *> (&nc), 4))
286  return false;
287  if (! is.read (reinterpret_cast<char *> (&nz), 4))
288  return false;
289 
290  if (swap)
291  {
292  swap_bytes<4> (&nr);
293  swap_bytes<4> (&nc);
294  swap_bytes<4> (&nz);
295  }
296 
297  SparseBoolMatrix m (static_cast<octave_idx_type> (nr),
298  static_cast<octave_idx_type> (nc),
299  static_cast<octave_idx_type> (nz));
300 
301  for (int i = 0; i < nc+1; i++)
302  {
303  octave_quit ();
304  if (! is.read (reinterpret_cast<char *> (&tmp), 4))
305  return false;
306  if (swap)
307  swap_bytes<4> (&tmp);
308  m.cidx (i) = tmp;
309  }
310 
311  for (int i = 0; i < nz; i++)
312  {
313  octave_quit ();
314  if (! is.read (reinterpret_cast<char *> (&tmp), 4))
315  return false;
316  if (swap)
317  swap_bytes<4> (&tmp);
318  m.ridx (i) = tmp;
319  }
320 
321  if (error_state || ! is)
322  return false;
323 
324  OCTAVE_LOCAL_BUFFER (char, htmp, nz);
325 
326  if (! is.read (htmp, nz))
327  return false;
328 
329  for (int i = 0; i < nz; i++)
330  m.data(i) = (htmp[i] ? 1 : 0);
331 
332  if (! m.indices_ok ())
333  return false;
334 
335  matrix = m;
336 
337  return true;
338 }
339 
340 #if defined (HAVE_HDF5)
341 
342 bool
343 octave_sparse_bool_matrix::save_hdf5 (hid_t loc_id, const char *name, bool)
344 {
345  dim_vector dv = dims ();
346  int empty = save_hdf5_empty (loc_id, name, dv);
347  if (empty)
348  return (empty > 0);
349 
350  // Ensure that additional memory is deallocated
352 #if HAVE_HDF5_18
353  hid_t group_hid = H5Gcreate (loc_id, name, H5P_DEFAULT, H5P_DEFAULT,
354  H5P_DEFAULT);
355 #else
356  hid_t group_hid = H5Gcreate (loc_id, name, 0);
357 #endif
358  if (group_hid < 0)
359  return false;
360 
361  hid_t space_hid = -1, data_hid = -1;
362  bool retval = true;
364  octave_idx_type tmp;
365  hsize_t hdims[2];
366 
367  space_hid = H5Screate_simple (0, hdims, 0);
368  if (space_hid < 0)
369  {
370  H5Gclose (group_hid);
371  return false;
372  }
373 #if HAVE_HDF5_18
374  data_hid = H5Dcreate (group_hid, "nr", H5T_NATIVE_IDX, space_hid,
375  H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
376 #else
377  data_hid = H5Dcreate (group_hid, "nr", H5T_NATIVE_IDX, space_hid,
378  H5P_DEFAULT);
379 #endif
380  if (data_hid < 0)
381  {
382  H5Sclose (space_hid);
383  H5Gclose (group_hid);
384  return false;
385  }
386 
387  tmp = m.rows ();
388  retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL,
389  H5S_ALL, H5P_DEFAULT, &tmp) >= 0;
390  H5Dclose (data_hid);
391  if (!retval)
392  {
393  H5Sclose (space_hid);
394  H5Gclose (group_hid);
395  return false;
396  }
397 
398 #if HAVE_HDF5_18
399  data_hid = H5Dcreate (group_hid, "nc", H5T_NATIVE_IDX, space_hid,
400  H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
401 #else
402  data_hid = H5Dcreate (group_hid, "nc", H5T_NATIVE_IDX, space_hid,
403  H5P_DEFAULT);
404 #endif
405  if (data_hid < 0)
406  {
407  H5Sclose (space_hid);
408  H5Gclose (group_hid);
409  return false;
410  }
411 
412  tmp = m.cols ();
413  retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
414  H5P_DEFAULT, &tmp) >= 0;
415  H5Dclose (data_hid);
416  if (!retval)
417  {
418  H5Sclose (space_hid);
419  H5Gclose (group_hid);
420  return false;
421  }
422 
423 #if HAVE_HDF5_18
424  data_hid = H5Dcreate (group_hid, "nz", H5T_NATIVE_IDX, space_hid,
425  H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
426 #else
427  data_hid = H5Dcreate (group_hid, "nz", H5T_NATIVE_IDX, space_hid,
428  H5P_DEFAULT);
429 #endif
430  if (data_hid < 0)
431  {
432  H5Sclose (space_hid);
433  H5Gclose (group_hid);
434  return false;
435  }
436 
437  tmp = m.nnz ();
438  retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
439  H5P_DEFAULT, &tmp) >= 0;
440  H5Dclose (data_hid);
441  if (!retval)
442  {
443  H5Sclose (space_hid);
444  H5Gclose (group_hid);
445  return false;
446  }
447 
448  H5Sclose (space_hid);
449 
450  hdims[0] = m.cols () + 1;
451  hdims[1] = 1;
452 
453  space_hid = H5Screate_simple (2, hdims, 0);
454 
455  if (space_hid < 0)
456  {
457  H5Gclose (group_hid);
458  return false;
459  }
460 
461 #if HAVE_HDF5_18
462  data_hid = H5Dcreate (group_hid, "cidx", H5T_NATIVE_IDX, space_hid,
463  H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
464 #else
465  data_hid = H5Dcreate (group_hid, "cidx", H5T_NATIVE_IDX, space_hid,
466  H5P_DEFAULT);
467 #endif
468  if (data_hid < 0)
469  {
470  H5Sclose (space_hid);
471  H5Gclose (group_hid);
472  return false;
473  }
474 
475  octave_idx_type * itmp = m.xcidx ();
476  retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
477  H5P_DEFAULT, itmp) >= 0;
478  H5Dclose (data_hid);
479  if (!retval)
480  {
481  H5Sclose (space_hid);
482  H5Gclose (group_hid);
483  return false;
484  }
485 
486  H5Sclose (space_hid);
487 
488  hdims[0] = m.nnz ();
489  hdims[1] = 1;
490 
491  space_hid = H5Screate_simple (2, hdims, 0);
492 
493  if (space_hid < 0)
494  {
495  H5Gclose (group_hid);
496  return false;
497  }
498 
499 #if HAVE_HDF5_18
500  data_hid = H5Dcreate (group_hid, "ridx", H5T_NATIVE_IDX, space_hid,
501  H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
502 #else
503  data_hid = H5Dcreate (group_hid, "ridx", H5T_NATIVE_IDX, space_hid,
504  H5P_DEFAULT);
505 #endif
506  if (data_hid < 0)
507  {
508  H5Sclose (space_hid);
509  H5Gclose (group_hid);
510  return false;
511  }
512 
513  itmp = m.xridx ();
514  retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
515  H5P_DEFAULT, itmp) >= 0;
516  H5Dclose (data_hid);
517  if (!retval)
518  {
519  H5Sclose (space_hid);
520  H5Gclose (group_hid);
521  return false;
522  }
523 
524 #if HAVE_HDF5_18
525  data_hid = H5Dcreate (group_hid, "data", H5T_NATIVE_HBOOL, space_hid,
526  H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
527 #else
528  data_hid = H5Dcreate (group_hid, "data", H5T_NATIVE_HBOOL, space_hid,
529  H5P_DEFAULT);
530 #endif
531  if (data_hid < 0)
532  {
533  H5Sclose (space_hid);
534  H5Gclose (group_hid);
535  return false;
536  }
537 
538  OCTAVE_LOCAL_BUFFER (hbool_t, htmp, m.nnz ());
539  for (int i = 0; i < m.nnz (); i++)
540  htmp[i] = m.xdata(i);
541 
542  retval = H5Dwrite (data_hid, H5T_NATIVE_HBOOL, H5S_ALL, H5S_ALL,
543  H5P_DEFAULT, htmp) >= 0;
544  H5Dclose (data_hid);
545  H5Sclose (space_hid);
546  H5Gclose (group_hid);
547 
548  return retval;
549 }
550 
551 bool
552 octave_sparse_bool_matrix::load_hdf5 (hid_t loc_id, const char *name)
553 {
554  octave_idx_type nr, nc, nz;
555  hid_t group_hid, data_hid, space_hid;
556  hsize_t rank;
557 
558  dim_vector dv;
559  int empty = load_hdf5_empty (loc_id, name, dv);
560  if (empty > 0)
561  matrix.resize (dv);
562  if (empty)
563  return (empty > 0);
564 
565 #if HAVE_HDF5_18
566  group_hid = H5Gopen (loc_id, name, H5P_DEFAULT);
567 #else
568  group_hid = H5Gopen (loc_id, name);
569 #endif
570  if (group_hid < 0 ) return false;
571 
572 #if HAVE_HDF5_18
573  data_hid = H5Dopen (group_hid, "nr", H5P_DEFAULT);
574 #else
575  data_hid = H5Dopen (group_hid, "nr");
576 #endif
577  space_hid = H5Dget_space (data_hid);
578  rank = H5Sget_simple_extent_ndims (space_hid);
579 
580  if (rank != 0)
581  {
582  H5Dclose (data_hid);
583  H5Gclose (group_hid);
584  return false;
585  }
586 
587  if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT, &nr)
588  < 0)
589  {
590  H5Dclose (data_hid);
591  H5Gclose (group_hid);
592  return false;
593  }
594 
595  H5Dclose (data_hid);
596 
597 #if HAVE_HDF5_18
598  data_hid = H5Dopen (group_hid, "nc", H5P_DEFAULT);
599 #else
600  data_hid = H5Dopen (group_hid, "nc");
601 #endif
602  space_hid = H5Dget_space (data_hid);
603  rank = H5Sget_simple_extent_ndims (space_hid);
604 
605  if (rank != 0)
606  {
607  H5Dclose (data_hid);
608  H5Gclose (group_hid);
609  return false;
610  }
611 
612  if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT, &nc)
613  < 0)
614  {
615  H5Dclose (data_hid);
616  H5Gclose (group_hid);
617  return false;
618  }
619 
620  H5Dclose (data_hid);
621 
622 #if HAVE_HDF5_18
623  data_hid = H5Dopen (group_hid, "nz", H5P_DEFAULT);
624 #else
625  data_hid = H5Dopen (group_hid, "nz");
626 #endif
627  space_hid = H5Dget_space (data_hid);
628  rank = H5Sget_simple_extent_ndims (space_hid);
629 
630  if (rank != 0)
631  {
632  H5Dclose (data_hid);
633  H5Gclose (group_hid);
634  return false;
635  }
636 
637  if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT, &nz)
638  < 0)
639  {
640  H5Dclose (data_hid);
641  H5Gclose (group_hid);
642  return false;
643  }
644 
645  H5Dclose (data_hid);
646 
647  SparseBoolMatrix m (static_cast<octave_idx_type> (nr),
648  static_cast<octave_idx_type> (nc),
649  static_cast<octave_idx_type> (nz));
650 
651 #if HAVE_HDF5_18
652  data_hid = H5Dopen (group_hid, "cidx", H5P_DEFAULT);
653 #else
654  data_hid = H5Dopen (group_hid, "cidx");
655 #endif
656  space_hid = H5Dget_space (data_hid);
657  rank = H5Sget_simple_extent_ndims (space_hid);
658 
659  if (rank != 2)
660  {
661  H5Sclose (space_hid);
662  H5Dclose (data_hid);
663  H5Gclose (group_hid);
664  return false;
665  }
666 
667  OCTAVE_LOCAL_BUFFER (hsize_t, hdims, rank);
668  OCTAVE_LOCAL_BUFFER (hsize_t, maxdims, rank);
669 
670  H5Sget_simple_extent_dims (space_hid, hdims, maxdims);
671 
672  if (static_cast<int> (hdims[0]) != nc + 1
673  || static_cast<int> (hdims[1]) != 1)
674  {
675  H5Sclose (space_hid);
676  H5Dclose (data_hid);
677  H5Gclose (group_hid);
678  return false;
679  }
680 
681  octave_idx_type *itmp = m.xcidx ();
682  if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT, itmp)
683  < 0)
684  {
685  H5Sclose (space_hid);
686  H5Dclose (data_hid);
687  H5Gclose (group_hid);
688  return false;
689  }
690 
691  H5Sclose (space_hid);
692  H5Dclose (data_hid);
693 
694 #if HAVE_HDF5_18
695  data_hid = H5Dopen (group_hid, "ridx", H5P_DEFAULT);
696 #else
697  data_hid = H5Dopen (group_hid, "ridx");
698 #endif
699  space_hid = H5Dget_space (data_hid);
700  rank = H5Sget_simple_extent_ndims (space_hid);
701 
702  if (rank != 2)
703  {
704  H5Sclose (space_hid);
705  H5Dclose (data_hid);
706  H5Gclose (group_hid);
707  return false;
708  }
709 
710  H5Sget_simple_extent_dims (space_hid, hdims, maxdims);
711 
712  if (static_cast<int> (hdims[0]) != nz
713  || static_cast<int> (hdims[1]) != 1)
714  {
715  H5Sclose (space_hid);
716  H5Dclose (data_hid);
717  H5Gclose (group_hid);
718  return false;
719  }
720 
721  itmp = m.xridx ();
722  if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
723  H5P_DEFAULT, itmp) < 0)
724  {
725  H5Sclose (space_hid);
726  H5Dclose (data_hid);
727  H5Gclose (group_hid);
728  return false;
729  }
730 
731  H5Sclose (space_hid);
732  H5Dclose (data_hid);
733 
734 #if HAVE_HDF5_18
735  data_hid = H5Dopen (group_hid, "data", H5P_DEFAULT);
736 #else
737  data_hid = H5Dopen (group_hid, "data");
738 #endif
739  space_hid = H5Dget_space (data_hid);
740  rank = H5Sget_simple_extent_ndims (space_hid);
741 
742  if (rank != 2)
743  {
744  H5Sclose (space_hid);
745  H5Dclose (data_hid);
746  H5Gclose (group_hid);
747  return false;
748  }
749 
750  H5Sget_simple_extent_dims (space_hid, hdims, maxdims);
751 
752  if (static_cast<int> (hdims[0]) != nz
753  || static_cast<int> (hdims[1]) != 1)
754  {
755  H5Sclose (space_hid);
756  H5Dclose (data_hid);
757  H5Gclose (group_hid);
758  return false;
759  }
760 
761  OCTAVE_LOCAL_BUFFER (hbool_t, htmp, nz);
762  bool retval = false;
763  if (H5Dread (data_hid, H5T_NATIVE_HBOOL, H5S_ALL, H5S_ALL,
764  H5P_DEFAULT, htmp) >= 0
765  && m.indices_ok ())
766  {
767  retval = true;
768 
769  for (int i = 0; i < nz; i++)
770  m.xdata(i) = htmp[i];
771 
772  matrix = m;
773  }
774 
775  H5Sclose (space_hid);
776  H5Dclose (data_hid);
777  H5Gclose (group_hid);
778 
779  return retval;
780 }
781 
782 #endif
783 
784 mxArray *
786 {
787  mwSize nz = nzmax ();
788  mxArray *retval = new mxArray (mxLOGICAL_CLASS, rows (), columns (),
789  nz, mxREAL);
790  bool *pr = static_cast<bool *> (retval->get_data ());
791  mwIndex *ir = retval->get_ir ();
792  mwIndex *jc = retval->get_jc ();
793 
794  for (mwIndex i = 0; i < nz; i++)
795  {
796  pr[i] = matrix.data (i);
797  ir[i] = matrix.ridx (i);
798  }
799 
800  for (mwIndex i = 0; i < columns () + 1; i++)
801  jc[i] = matrix.cidx (i);
802 
803  return retval;
804 }