Mercurial > hg > octave-lyh

comparison libinterp/octave-value/ov-re-sparse.cc @ 15195:2fc554ffbc28

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
split libinterp from src * libinterp: New directory. Move all files from src directory here except Makefile.am, main.cc, main-cli.cc, mkoctfile.in.cc, mkoctfilr.in.sh, octave-config.in.cc, octave-config.in.sh. * libinterp/Makefile.am: New file, extracted from src/Makefile.am. * src/Makefile.am: Delete everything except targets and definitions needed to build and link main and utility programs. * Makefile.am (SUBDIRS): Include libinterp in the list. * autogen.sh: Run config-module.sh in libinterp/dldfcn directory, not src/dldfcn directory. * configure.ac (AC_CONFIG_SRCDIR): Use libinterp/octave.cc, not src/octave.cc. (DL_LDFLAGS, LIBOCTINTERP): Use libinterp, not src. (AC_CONFIG_FILES): Include libinterp/Makefile in the list. * find-docstring-files.sh: Look in libinterp, not src. * gui/src/Makefile.am (liboctgui_la_CPPFLAGS): Find header files in libinterp, not src.
author John W. Eaton <jwe@octave.org>
date Sat, 18 Aug 2012 16:23:39 -0400
parents src/octave-value/ov-re-sparse.cc@62a35ae7d6a2
children 9020dddc925a
comparison
equal deleted inserted replaced
15194:0f0b795044c3 15195:2fc554ffbc28
1 /*
2
3 Copyright (C) 2004-2012 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 <climits>
29
30 #include <iostream>
31 #include <vector>
32
33 #include "lo-specfun.h"
34 #include "lo-mappers.h"
35 #include "oct-locbuf.h"
36
37 #include "mxarray.h"
38 #include "ov-base.h"
39 #include "ov-scalar.h"
40 #include "gripes.h"
41
42 #include "ls-hdf5.h"
43
44 #include "ov-re-sparse.h"
45
46 #include "ov-base-sparse.h"
47 #include "ov-base-sparse.cc"
48
49 #include "ov-bool-sparse.h"
50
51 template class OCTINTERP_API octave_base_sparse<SparseMatrix>;
52
53 DEFINE_OCTAVE_ALLOCATOR (octave_sparse_matrix);
54
55 DEFINE_OV_TYPEID_FUNCTIONS_AND_DATA (octave_sparse_matrix, "sparse matrix", "double");
56
57 idx_vector
58 octave_sparse_matrix::index_vector (void) const
59 {
60 if (matrix.numel () == matrix.nnz ())
61 return idx_vector (array_value ());
62 else
63 {
64 std::string nm = type_name ();
65 error ("%s type invalid as index value", nm.c_str ());
66 return idx_vector ();
67 }
68 }
69
70 octave_base_value *
71 octave_sparse_matrix::try_narrowing_conversion (void)
72 {
73 octave_base_value *retval = 0;
74
75 if (Vsparse_auto_mutate)
76 {
77 // Don't use numel, since it can overflow for very large matrices
78 // Note that for the second test, this means it becomes approximative
79 // since it involves a cast to double to avoid issues of overflow
80 if (matrix.rows () == 1 && matrix.cols () == 1)
81 {
82 // Const copy of the matrix, so the right version of () operator used
83 const SparseMatrix tmp (matrix);
84
85 retval = new octave_scalar (tmp (0));
86 }
87 else if (matrix.cols () > 0 && matrix.rows () > 0
88 && (double (matrix.byte_size ()) > double (matrix.rows ())
89 * double (matrix.cols ()) * sizeof (double)))
90 retval = new octave_matrix (matrix.matrix_value ());
91 }
92
93 return retval;
94 }
95
96 double
97 octave_sparse_matrix::double_value (bool) const
98 {
99 double retval = lo_ieee_nan_value ();
100
101 if (numel () > 0)
102 {
103 if (numel () > 1)
104 gripe_implicit_conversion ("Octave:array-to-scalar",
105 "real sparse matrix", "real scalar");
106
107 retval = matrix (0, 0);
108 }
109 else
110 gripe_invalid_conversion ("real sparse matrix", "real scalar");
111
112 return retval;
113 }
114
115 Complex
116 octave_sparse_matrix::complex_value (bool) const
117 {
118 double tmp = lo_ieee_nan_value ();
119
120 Complex retval (tmp, tmp);
121
122 // FIXME -- maybe this should be a function, valid_as_scalar()
123 if (rows () > 0 && columns () > 0)
124 {
125 if (numel () > 1)
126 gripe_implicit_conversion ("Octave:array-to-scalar",
127 "real sparse matrix", "complex scalar");
128
129 retval = matrix (0, 0);
130 }
131 else
132 gripe_invalid_conversion ("real sparse matrix", "complex scalar");
133
134 return retval;
135 }
136
137 Matrix
138 octave_sparse_matrix::matrix_value (bool) const
139 {
140 return matrix.matrix_value ();
141 }
142
143 boolNDArray
144 octave_sparse_matrix::bool_array_value (bool warn) const
145 {
146 NDArray m = matrix.matrix_value ();
147
148 if (m.any_element_is_nan ())
149 gripe_nan_to_logical_conversion ();
150 else if (warn && m.any_element_not_one_or_zero ())
151 gripe_logical_conversion ();
152
153 return boolNDArray (m);
154 }
155
156 charNDArray
157 octave_sparse_matrix::char_array_value (bool) const
158 {
159 charNDArray retval (dims (), 0);
160 octave_idx_type nc = matrix.cols ();
161 octave_idx_type nr = matrix.rows ();
162
163 for (octave_idx_type j = 0; j < nc; j++)
164 for (octave_idx_type i = matrix.cidx (j); i < matrix.cidx (j+1); i++)
165 retval(matrix.ridx (i) + nr * j) = static_cast<char>(matrix.data (i));
166
167 return retval;
168 }
169
170 ComplexMatrix
171 octave_sparse_matrix::complex_matrix_value (bool) const
172 {
173 return ComplexMatrix (matrix.matrix_value ());
174 }
175
176 ComplexNDArray
177 octave_sparse_matrix::complex_array_value (bool) const
178 {
179 return ComplexNDArray (ComplexMatrix (matrix.matrix_value ()));
180 }
181
182 NDArray
183 octave_sparse_matrix::array_value (bool) const
184 {
185 return NDArray (matrix.matrix_value ());
186 }
187
188 SparseBoolMatrix
189 octave_sparse_matrix::sparse_bool_matrix_value (bool warn) const
190 {
191 if (matrix.any_element_is_nan ())
192 gripe_nan_to_logical_conversion ();
193 else if (warn && matrix.any_element_not_one_or_zero ())
194 gripe_logical_conversion ();
195
196 return mx_el_ne (matrix, 0.0);
197 }
198
199 octave_value
200 octave_sparse_matrix::convert_to_str_internal (bool, bool, char type) const
201 {
202 octave_value retval;
203 dim_vector dv = dims ();
204 octave_idx_type nel = dv.numel ();
205
206 if (nel == 0)
207 {
208 char s = '\0';
209 retval = octave_value (&s, type);
210 }
211 else
212 {
213 octave_idx_type nr = matrix.rows ();
214 octave_idx_type nc = matrix.cols ();
215 charNDArray chm (dv, static_cast<char> (0));
216
217 bool warned = false;
218
219 for (octave_idx_type j = 0; j < nc; j++)
220 for (octave_idx_type i = matrix.cidx (j);
221 i < matrix.cidx (j+1); i++)
222 {
223 octave_quit ();
224
225 double d = matrix.data (i);
226
227 if (xisnan (d))
228 {
229 gripe_nan_to_character_conversion ();
230 return retval;
231 }
232 else
233 {
234 int ival = NINT (d);
235
236 if (ival < 0 || ival > UCHAR_MAX)
237 {
238 // FIXME -- is there something
239 // better we could do?
240
241 ival = 0;
242
243 if (! warned)
244 {
245 ::warning ("range error for conversion to character value");
246 warned = true;
247 }
248 }
249
250 chm (matrix.ridx (i) + j * nr) =
251 static_cast<char> (ival);
252 }
253 }
254
255 retval = octave_value (chm, type);
256 }
257
258 return retval;
259 }
260
261 bool
262 octave_sparse_matrix::save_binary (std::ostream& os, bool&save_as_floats)
263 {
264 dim_vector d = this->dims ();
265 if (d.length () < 1)
266 return false;
267
268 // Ensure that additional memory is deallocated
269 matrix.maybe_compress ();
270
271 int nr = d(0);
272 int nc = d(1);
273 int nz = nnz ();
274
275 int32_t itmp;
276 // Use negative value for ndims to be consistent with other formats
277 itmp= -2;
278 os.write (reinterpret_cast<char *> (&itmp), 4);
279
280 itmp= nr;
281 os.write (reinterpret_cast<char *> (&itmp), 4);
282
283 itmp= nc;
284 os.write (reinterpret_cast<char *> (&itmp), 4);
285
286 itmp= nz;
287 os.write (reinterpret_cast<char *> (&itmp), 4);
288
289 save_type st = LS_DOUBLE;
290 if (save_as_floats)
291 {
292 if (matrix.too_large_for_float ())
293 {
294 warning ("save: some values too large to save as floats --");
295 warning ("save: saving as doubles instead");
296 }
297 else
298 st = LS_FLOAT;
299 }
300 else if (matrix.nnz () > 8192) // FIXME -- make this configurable.
301 {
302 double max_val, min_val;
303 if (matrix.all_integers (max_val, min_val))
304 st = get_save_type (max_val, min_val);
305 }
306
307 // add one to the printed indices to go from
308 // zero-based to one-based arrays
309 for (int i = 0; i < nc+1; i++)
310 {
311 octave_quit ();
312 itmp = matrix.cidx (i);
313 os.write (reinterpret_cast<char *> (&itmp), 4);
314 }
315
316 for (int i = 0; i < nz; i++)
317 {
318 octave_quit ();
319 itmp = matrix.ridx (i);
320 os.write (reinterpret_cast<char *> (&itmp), 4);
321 }
322
323 write_doubles (os, matrix.data (), st, nz);
324
325 return true;
326 }
327
328 bool
329 octave_sparse_matrix::load_binary (std::istream& is, bool swap,
330 oct_mach_info::float_format fmt)
331 {
332 int32_t nz, nc, nr, tmp;
333 char ctmp;
334
335 if (! is.read (reinterpret_cast<char *> (&tmp), 4))
336 return false;
337
338 if (swap)
339 swap_bytes<4> (&tmp);
340
341 if (tmp != -2) {
342 error ("load: only 2D sparse matrices are supported");
343 return false;
344 }
345
346 if (! is.read (reinterpret_cast<char *> (&nr), 4))
347 return false;
348 if (! is.read (reinterpret_cast<char *> (&nc), 4))
349 return false;
350 if (! is.read (reinterpret_cast<char *> (&nz), 4))
351 return false;
352
353 if (swap)
354 {
355 swap_bytes<4> (&nr);
356 swap_bytes<4> (&nc);
357 swap_bytes<4> (&nz);
358 }
359
360 SparseMatrix m (static_cast<octave_idx_type> (nr),
361 static_cast<octave_idx_type> (nc),
362 static_cast<octave_idx_type> (nz));
363
364 for (int i = 0; i < nc+1; i++)
365 {
366 octave_quit ();
367 if (! is.read (reinterpret_cast<char *> (&tmp), 4))
368 return false;
369 if (swap)
370 swap_bytes<4> (&tmp);
371 m.xcidx (i) = tmp;
372 }
373
374 for (int i = 0; i < nz; i++)
375 {
376 octave_quit ();
377 if (! is.read (reinterpret_cast<char *> (&tmp), 4))
378 return false;
379 if (swap)
380 swap_bytes<4> (&tmp);
381 m.xridx (i) = tmp;
382 }
383
384 if (! is.read (reinterpret_cast<char *> (&ctmp), 1))
385 return false;
386
387 read_doubles (is, m.xdata (), static_cast<save_type> (ctmp), nz, swap, fmt);
388
389 if (error_state || ! is)
390 return false;
391
392 if (! m.indices_ok ())
393 return false;
394
395 matrix = m;
396
397 return true;
398 }
399
400 #if defined (HAVE_HDF5)
401
402 bool
403 octave_sparse_matrix::save_hdf5 (hid_t loc_id, const char *name,
404 bool save_as_floats)
405 {
406 dim_vector dv = dims ();
407 int empty = save_hdf5_empty (loc_id, name, dv);
408 if (empty)
409 return (empty > 0);
410
411 // Ensure that additional memory is deallocated
412 matrix.maybe_compress ();
413
414 #if HAVE_HDF5_18
415 hid_t group_hid = H5Gcreate (loc_id, name, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
416 #else
417 hid_t group_hid = H5Gcreate (loc_id, name, 0);
418 #endif
419 if (group_hid < 0)
420 return false;
421
422 hid_t space_hid = -1, data_hid = -1;
423 bool retval = true;
424 SparseMatrix m = sparse_matrix_value ();
425 octave_idx_type tmp;
426 hsize_t hdims[2];
427
428 space_hid = H5Screate_simple (0, hdims, 0);
429 if (space_hid < 0)
430 {
431 H5Gclose (group_hid);
432 return false;
433 }
434 #if HAVE_HDF5_18
435 data_hid = H5Dcreate (group_hid, "nr", H5T_NATIVE_IDX, space_hid,
436 H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
437 #else
438 data_hid = H5Dcreate (group_hid, "nr", H5T_NATIVE_IDX, space_hid,
439 H5P_DEFAULT);
440 #endif
441 if (data_hid < 0)
442 {
443 H5Sclose (space_hid);
444 H5Gclose (group_hid);
445 return false;
446 }
447
448 tmp = m.rows ();
449 retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT,
450 &tmp) >= 0;
451 H5Dclose (data_hid);
452 if (!retval)
453 {
454 H5Sclose (space_hid);
455 H5Gclose (group_hid);
456 return false;
457 }
458 #if HAVE_HDF5_18
459 data_hid = H5Dcreate (group_hid, "nc", H5T_NATIVE_IDX, space_hid,
460 H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
461 #else
462 data_hid = H5Dcreate (group_hid, "nc", H5T_NATIVE_IDX, space_hid,
463 H5P_DEFAULT);
464 #endif
465 if (data_hid < 0)
466 {
467 H5Sclose (space_hid);
468 H5Gclose (group_hid);
469 return false;
470 }
471
472 tmp = m.cols ();
473 retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT,
474 &tmp) >= 0;
475 H5Dclose (data_hid);
476 if (!retval)
477 {
478 H5Sclose (space_hid);
479 H5Gclose (group_hid);
480 return false;
481 }
482
483 #if HAVE_HDF5_18
484 data_hid = H5Dcreate (group_hid, "nz", H5T_NATIVE_IDX, space_hid,
485 H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
486 #else
487 data_hid = H5Dcreate (group_hid, "nz", H5T_NATIVE_IDX, space_hid,
488 H5P_DEFAULT);
489 #endif
490 if (data_hid < 0)
491 {
492 H5Sclose (space_hid);
493 H5Gclose (group_hid);
494 return false;
495 }
496
497 tmp = m.nnz ();
498 retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT,
499 &tmp) >= 0;
500 H5Dclose (data_hid);
501 if (!retval)
502 {
503 H5Sclose (space_hid);
504 H5Gclose (group_hid);
505 return false;
506 }
507
508 H5Sclose (space_hid);
509
510 hdims[0] = m.cols () + 1;
511 hdims[1] = 1;
512
513 space_hid = H5Screate_simple (2, hdims, 0);
514
515 if (space_hid < 0)
516 {
517 H5Gclose (group_hid);
518 return false;
519 }
520
521 #if HAVE_HDF5_18
522 data_hid = H5Dcreate (group_hid, "cidx", H5T_NATIVE_IDX, space_hid,
523 H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
524 #else
525 data_hid = H5Dcreate (group_hid, "cidx", H5T_NATIVE_IDX, space_hid,
526 H5P_DEFAULT);
527 #endif
528 if (data_hid < 0)
529 {
530 H5Sclose (space_hid);
531 H5Gclose (group_hid);
532 return false;
533 }
534
535 octave_idx_type * itmp = m.xcidx ();
536 retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT,
537 itmp) >= 0;
538 H5Dclose (data_hid);
539 if (!retval)
540 {
541 H5Sclose (space_hid);
542 H5Gclose (group_hid);
543 return false;
544 }
545
546 H5Sclose (space_hid);
547
548 hdims[0] = m.nnz ();
549 hdims[1] = 1;
550
551 space_hid = H5Screate_simple (2, hdims, 0);
552
553 if (space_hid < 0)
554 {
555 H5Gclose (group_hid);
556 return false;
557 }
558 #if HAVE_HDF5_18
559 data_hid = H5Dcreate (group_hid, "ridx", H5T_NATIVE_IDX, space_hid,
560 H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
561 #else
562 data_hid = H5Dcreate (group_hid, "ridx", H5T_NATIVE_IDX, space_hid,
563 H5P_DEFAULT);
564 #endif
565 if (data_hid < 0)
566 {
567 H5Sclose (space_hid);
568 H5Gclose (group_hid);
569 return false;
570 }
571
572 itmp = m.xridx ();
573 retval = H5Dwrite (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL, H5P_DEFAULT,
574 itmp) >= 0;
575 H5Dclose (data_hid);
576 if (!retval)
577 {
578 H5Sclose (space_hid);
579 H5Gclose (group_hid);
580 return false;
581 }
582
583 hid_t save_type_hid = H5T_NATIVE_DOUBLE;
584
585 if (save_as_floats)
586 {
587 if (m.too_large_for_float ())
588 {
589 warning ("save: some values too large to save as floats --");
590 warning ("save: saving as doubles instead");
591 }
592 else
593 save_type_hid = H5T_NATIVE_FLOAT;
594 }
595 #if HAVE_HDF5_INT2FLOAT_CONVERSIONS
596 // hdf5 currently doesn't support float/integer conversions
597 else
598 {
599 double max_val, min_val;
600
601 if (m.all_integers (max_val, min_val))
602 save_type_hid
603 = save_type_to_hdf5 (get_save_type (max_val, min_val));
604 }
605 #endif /* HAVE_HDF5_INT2FLOAT_CONVERSIONS */
606
607 #if HAVE_HDF5_18
608 data_hid = H5Dcreate (group_hid, "data", save_type_hid, space_hid,
609 H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
610 #else
611 data_hid = H5Dcreate (group_hid, "data", save_type_hid, space_hid,
612 H5P_DEFAULT);
613 #endif
614 if (data_hid < 0)
615 {
616 H5Sclose (space_hid);
617 H5Gclose (group_hid);
618 return false;
619 }
620
621 double * dtmp = m.xdata ();
622 retval = H5Dwrite (data_hid, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL,
623 H5P_DEFAULT, dtmp) >= 0;
624 H5Dclose (data_hid);
625 H5Sclose (space_hid);
626 H5Gclose (group_hid);
627
628 return retval;
629 }
630
631 bool
632 octave_sparse_matrix::load_hdf5 (hid_t loc_id, const char *name)
633 {
634 octave_idx_type nr, nc, nz;
635 hid_t group_hid, data_hid, space_hid;
636 hsize_t rank;
637
638 dim_vector dv;
639 int empty = load_hdf5_empty (loc_id, name, dv);
640 if (empty > 0)
641 matrix.resize (dv);
642 if (empty)
643 return (empty > 0);
644
645 #if HAVE_HDF5_18
646 group_hid = H5Gopen (loc_id, name, H5P_DEFAULT);
647 #else
648 group_hid = H5Gopen (loc_id, name);
649 #endif
650 if (group_hid < 0) return false;
651
652 #if HAVE_HDF5_18
653 data_hid = H5Dopen (group_hid, "nr", H5P_DEFAULT);
654 #else
655 data_hid = H5Dopen (group_hid, "nr");
656 #endif
657 space_hid = H5Dget_space (data_hid);
658 rank = H5Sget_simple_extent_ndims (space_hid);
659
660 if (rank != 0)
661 {
662 H5Dclose (data_hid);
663 H5Gclose (group_hid);
664 return false;
665 }
666
667 if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
668 H5P_DEFAULT, &nr) < 0)
669 {
670 H5Dclose (data_hid);
671 H5Gclose (group_hid);
672 return false;
673 }
674
675 H5Dclose (data_hid);
676
677 #if HAVE_HDF5_18
678 data_hid = H5Dopen (group_hid, "nc", H5P_DEFAULT);
679 #else
680 data_hid = H5Dopen (group_hid, "nc");
681 #endif
682 space_hid = H5Dget_space (data_hid);
683 rank = H5Sget_simple_extent_ndims (space_hid);
684
685 if (rank != 0)
686 {
687 H5Dclose (data_hid);
688 H5Gclose (group_hid);
689 return false;
690 }
691
692 if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
693 H5P_DEFAULT, &nc) < 0)
694 {
695 H5Dclose (data_hid);
696 H5Gclose (group_hid);
697 return false;
698 }
699
700 H5Dclose (data_hid);
701
702 #if HAVE_HDF5_18
703 data_hid = H5Dopen (group_hid, "nz", H5P_DEFAULT);
704 #else
705 data_hid = H5Dopen (group_hid, "nz");
706 #endif
707 space_hid = H5Dget_space (data_hid);
708 rank = H5Sget_simple_extent_ndims (space_hid);
709
710 if (rank != 0)
711 {
712 H5Dclose (data_hid);
713 H5Gclose (group_hid);
714 return false;
715 }
716
717 if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
718 H5P_DEFAULT, &nz) < 0)
719 {
720 H5Dclose (data_hid);
721 H5Gclose (group_hid);
722 return false;
723 }
724
725 H5Dclose (data_hid);
726
727 SparseMatrix m (static_cast<octave_idx_type> (nr),
728 static_cast<octave_idx_type> (nc),
729 static_cast<octave_idx_type> (nz));
730
731 #if HAVE_HDF5_18
732 data_hid = H5Dopen (group_hid, "cidx", H5P_DEFAULT);
733 #else
734 data_hid = H5Dopen (group_hid, "cidx");
735 #endif
736 space_hid = H5Dget_space (data_hid);
737 rank = H5Sget_simple_extent_ndims (space_hid);
738
739 if (rank != 2)
740 {
741 H5Sclose (space_hid);
742 H5Dclose (data_hid);
743 H5Gclose (group_hid);
744 return false;
745 }
746
747 OCTAVE_LOCAL_BUFFER (hsize_t, hdims, rank);
748 OCTAVE_LOCAL_BUFFER (hsize_t, maxdims, rank);
749
750 H5Sget_simple_extent_dims (space_hid, hdims, maxdims);
751
752 if (static_cast<int> (hdims[0]) != nc + 1
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_idx_type *itmp = m.xcidx ();
762 if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
763 H5P_DEFAULT, itmp) < 0)
764 {
765 H5Sclose (space_hid);
766 H5Dclose (data_hid);
767 H5Gclose (group_hid);
768 return false;
769 }
770
771 H5Sclose (space_hid);
772 H5Dclose (data_hid);
773
774 #if HAVE_HDF5_18
775 data_hid = H5Dopen (group_hid, "ridx", H5P_DEFAULT);
776 #else
777 data_hid = H5Dopen (group_hid, "ridx");
778 #endif
779 space_hid = H5Dget_space (data_hid);
780 rank = H5Sget_simple_extent_ndims (space_hid);
781
782 if (rank != 2)
783 {
784 H5Sclose (space_hid);
785 H5Dclose (data_hid);
786 H5Gclose (group_hid);
787 return false;
788 }
789
790 H5Sget_simple_extent_dims (space_hid, hdims, maxdims);
791
792 if (static_cast<int> (hdims[0]) != nz || static_cast<int> (hdims[1]) != 1)
793 {
794 H5Sclose (space_hid);
795 H5Dclose (data_hid);
796 H5Gclose (group_hid);
797 return false;
798 }
799
800 itmp = m.xridx ();
801 if (H5Dread (data_hid, H5T_NATIVE_IDX, H5S_ALL, H5S_ALL,
802 H5P_DEFAULT, itmp) < 0)
803 {
804 H5Sclose (space_hid);
805 H5Dclose (data_hid);
806 H5Gclose (group_hid);
807 return false;
808 }
809
810 H5Sclose (space_hid);
811 H5Dclose (data_hid);
812
813 #if HAVE_HDF5_18
814 data_hid = H5Dopen (group_hid, "data", H5P_DEFAULT);
815 #else
816 data_hid = H5Dopen (group_hid, "data");
817 #endif
818 space_hid = H5Dget_space (data_hid);
819 rank = H5Sget_simple_extent_ndims (space_hid);
820
821 if (rank != 2)
822 {
823 H5Sclose (space_hid);
824 H5Dclose (data_hid);
825 H5Gclose (group_hid);
826 return false;
827 }
828
829 H5Sget_simple_extent_dims (space_hid, hdims, maxdims);
830
831 if (static_cast<int> (hdims[0]) != nz || static_cast<int> (hdims[1]) != 1)
832 {
833 H5Sclose (space_hid);
834 H5Dclose (data_hid);
835 H5Gclose (group_hid);
836 return false;
837 }
838
839 double *dtmp = m.xdata ();
840 bool retval = false;
841 if (H5Dread (data_hid, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL,
842 H5P_DEFAULT, dtmp) >= 0
843 && m.indices_ok ())
844 {
845 retval = true;
846 matrix = m;
847 }
848
849 H5Sclose (space_hid);
850 H5Dclose (data_hid);
851 H5Gclose (group_hid);
852
853 return retval;
854 }
855
856 #endif
857
858 mxArray *
859 octave_sparse_matrix::as_mxArray (void) const
860 {
861 mwSize nz = nzmax ();
862 mwSize nr = rows ();
863 mwSize nc = columns ();
864 mxArray *retval = new mxArray (mxDOUBLE_CLASS, nr, nc, nz, mxREAL);
865 double *pr = static_cast<double *> (retval->get_data ());
866 mwIndex *ir = retval->get_ir ();
867 mwIndex *jc = retval->get_jc ();
868
869 for (mwIndex i = 0; i < nz; i++)
870 {
871 pr[i] = matrix.data (i);
872 ir[i] = matrix.ridx (i);
873 }
874
875 for (mwIndex i = 0; i < nc + 1; i++)
876 jc[i] = matrix.cidx (i);
877
878 return retval;
879 }
880
881 octave_value
882 octave_sparse_matrix::map (unary_mapper_t umap) const
883 {
884 switch (umap)
885 {
886 case umap_imag:
887 return SparseMatrix (matrix.rows (), matrix.cols (), 0.0);
888
889 case umap_real:
890 case umap_conj:
891 return matrix;
892
893 // Mappers handled specially.
894 #define ARRAY_METHOD_MAPPER(UMAP, FCN) \
895 case umap_ ## UMAP: \
896 return octave_value (matrix.FCN ())
897
898 ARRAY_METHOD_MAPPER (abs, abs);
899
900 #define ARRAY_MAPPER(UMAP, TYPE, FCN) \
901 case umap_ ## UMAP: \
902 return octave_value (matrix.map<TYPE> (FCN))
903
904 ARRAY_MAPPER (acos, Complex, rc_acos);
905 ARRAY_MAPPER (acosh, Complex, rc_acosh);
906 ARRAY_MAPPER (angle, double, ::arg);
907 ARRAY_MAPPER (arg, double, ::arg);
908 ARRAY_MAPPER (asin, Complex, rc_asin);
909 ARRAY_MAPPER (asinh, double, ::asinh);
910 ARRAY_MAPPER (atan, double, ::atan);
911 ARRAY_MAPPER (atanh, Complex, rc_atanh);
912 ARRAY_MAPPER (erf, double, ::erf);
913 ARRAY_MAPPER (erfinv, double, ::erfinv);
914 ARRAY_MAPPER (erfcinv, double, ::erfcinv);
915 ARRAY_MAPPER (erfc, double, ::erfc);
916 ARRAY_MAPPER (gamma, double, xgamma);
917 ARRAY_MAPPER (lgamma, Complex, rc_lgamma);
918 ARRAY_MAPPER (cbrt, double, ::cbrt);
919 ARRAY_MAPPER (ceil, double, ::ceil);
920 ARRAY_MAPPER (cos, double, ::cos);
921 ARRAY_MAPPER (cosh, double, ::cosh);
922 ARRAY_MAPPER (exp, double, ::exp);
923 ARRAY_MAPPER (expm1, double, ::expm1);
924 ARRAY_MAPPER (fix, double, ::fix);
925 ARRAY_MAPPER (floor, double, ::floor);
926 ARRAY_MAPPER (log, Complex, rc_log);
927 ARRAY_MAPPER (log2, Complex, rc_log2);
928 ARRAY_MAPPER (log10, Complex, rc_log10);
929 ARRAY_MAPPER (log1p, Complex, rc_log1p);
930 ARRAY_MAPPER (round, double, xround);
931 ARRAY_MAPPER (roundb, double, xroundb);
932 ARRAY_MAPPER (signum, double, ::signum);
933 ARRAY_MAPPER (sin, double, ::sin);
934 ARRAY_MAPPER (sinh, double, ::sinh);
935 ARRAY_MAPPER (sqrt, Complex, rc_sqrt);
936 ARRAY_MAPPER (tan, double, ::tan);
937 ARRAY_MAPPER (tanh, double, ::tanh);
938 ARRAY_MAPPER (isnan, bool, xisnan);
939 ARRAY_MAPPER (isna, bool, octave_is_NA);
940 ARRAY_MAPPER (isinf, bool, xisinf);
941 ARRAY_MAPPER (finite, bool, xfinite);
942
943 default: // Attempt to go via dense matrix.
944 return octave_base_sparse<SparseMatrix>::map (umap);
945 }
946 }