GNU Octave 7.1.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
ov-base-sparse.cc
Go to the documentation of this file.
1////////////////////////////////////////////////////////////////////////
2//
3// Copyright (C) 1998-2022 The Octave Project Developers
4//
5// See the file COPYRIGHT.md in the top-level directory of this
6// distribution or <https://octave.org/copyright/>.
7//
8// This file is part of Octave.
9//
10// Octave is free software: you can redistribute it and/or modify it
11// under the terms of the GNU General Public License as published by
12// the Free Software Foundation, either version 3 of the License, or
13// (at your option) any later version.
14//
15// Octave is distributed in the hope that it will be useful, but
16// WITHOUT ANY WARRANTY; without even the implied warranty of
17// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18// GNU General Public License for more details.
19//
20// You should have received a copy of the GNU General Public License
21// along with Octave; see the file COPYING. If not, see
22// <https://www.gnu.org/licenses/>.
23//
24////////////////////////////////////////////////////////////////////////
25
26// This file should not include config.h. It is only included in other
27// C++ source files that should have included config.h before including
28// this file.
29
30#include <iomanip>
31#include <istream>
32#include <ostream>
33#include <sstream>
34
35#include "ovl.h"
36#include "ov-base.h"
37#include "quit.h"
38#include "pr-output.h"
39
40#include "byte-swap.h"
41#include "ls-oct-text.h"
42#include "ls-utils.h"
43#include "ls-hdf5.h"
44
45#include "boolSparse.h"
46#include "ov-base-sparse.h"
48#include "pager.h"
49#include "utils.h"
50
51#include "lo-array-errwarn.h"
52
53template <typename T>
56 bool resize_ok)
57{
58 octave_value retval;
59
60 octave_idx_type n_idx = idx.length ();
61
62 // If we catch an indexing error in index_vector, we flag an error in
63 // index k. Ensure it is the right value before each idx_vector call.
64 // Same variable as used in the for loop in the default case.
65
66 octave_idx_type k = 0;
67
68 try
69 {
70 switch (n_idx)
71 {
72 case 0:
73 retval = matrix;
74 break;
75
76 case 1:
77 {
78 octave::idx_vector i = idx (0).index_vector ();
79
80 retval = octave_value (matrix.index (i, resize_ok));
81 }
82 break;
83
84 case 2:
85 {
86 octave::idx_vector i = idx (0).index_vector ();
87
88 k = 1;
89 octave::idx_vector j = idx (1).index_vector ();
90
91 retval = octave_value (matrix.index (i, j, resize_ok));
92 }
93 break;
94
95 default:
96 error ("sparse indexing needs 1 or 2 indices");
97 }
98 }
99 catch (octave::index_exception& ie)
100 {
101 // Rethrow to allow more info to be reported later.
102 ie.set_pos_if_unset (n_idx, k+1);
103 throw;
104 }
105
106 return retval;
107}
108
109template <typename T>
111octave_base_sparse<T>::subsref (const std::string& type,
112 const std::list<octave_value_list>& idx)
113{
114 octave_value retval;
115
116 switch (type[0])
117 {
118 case '(':
119 retval = do_index_op (idx.front ());
120 break;
121
122 case '{':
123 case '.':
124 {
125 std::string nm = type_name ();
126 error ("%s cannot be indexed with %c", nm.c_str (), type[0]);
127 }
128 break;
129
130 default:
132 }
133
134 return retval.next_subsref (type, idx);
135}
136
137template <typename T>
139octave_base_sparse<T>::subsasgn (const std::string& type,
140 const std::list<octave_value_list>& idx,
141 const octave_value& rhs)
142{
143 octave_value retval;
144
145 switch (type[0])
146 {
147 case '(':
148 {
149 if (type.length () != 1)
150 {
151 std::string nm = type_name ();
152 error ("in indexed assignment of %s, last lhs index must be ()",
153 nm.c_str ());
154 }
155
156 retval = numeric_assign (type, idx, rhs);
157 }
158 break;
159
160 case '{':
161 case '.':
162 {
163 if (! isempty ())
165 std::string nm = type_name ();
166 error ("%s cannot be indexed with %c", nm.c_str (), type[0]);
167 }
168
170
171 retval = tmp.subsasgn (type, idx, rhs);
172 }
173 break;
174
175 default:
178
179 return retval;
180}
181
182template <typename MT>
183void
185{
186 octave_idx_type len = idx.length ();
187
188 // If we catch an indexing error in index_vector, we flag an error in
189 // index k. Ensure it is the right value before each idx_vector call.
190 // Same variable as used in the for loop in the default case.
191
192 octave_idx_type k = 0;
193
194 try
195 {
196 switch (len)
197 {
198 case 1:
199 {
200 octave::idx_vector i = idx (0).index_vector ();
201
202 matrix.delete_elements (i);
203
204 break;
205 }
206
207 case 2:
208 {
209 octave::idx_vector i = idx (0).index_vector ();
210
211 k = 1;
212 octave::idx_vector j = idx (1).index_vector ();
214 matrix.delete_elements (i, j);
216 break;
218
219 default:
220 error ("sparse indexing needs 1 or 2 indices");
221 }
222 }
224 {
225 // Rethrow to allow more info to be reported later.
226 ie.set_pos_if_unset (len, k+1);
227 throw;
228 }
230 // Invalidate the matrix type
231 typ.invalidate_type ();
233
234template <typename T>
237{
238 T retval (matrix);
239 retval.resize (dv);
240 return retval;
241}
242
243template <typename T>
244bool
246{
247 bool retval = false;
248 dim_vector dv = matrix.dims ();
249 octave_idx_type nel = dv.numel ();
250 octave_idx_type nz = nnz ();
251
252 if (nel > 0)
253 {
254 T t1 (matrix.reshape (dim_vector (nel, 1)));
255
256 if (t1.any_element_is_nan ())
258
259 if (nel > 1)
261
262 if (nz == nel)
263 {
264 SparseBoolMatrix t2 = t1.all ();
265
266 retval = t2(0);
267 }
268 }
269
270 return retval;
271}
272
273template <typename T>
274bool
276{
277 dim_vector dv = dims ();
278
279 return (dv.all_ones () || dv.any_zero ());
280}
281
282template <typename T>
283void
284octave_base_sparse<T>::print (std::ostream& os, bool pr_as_read_syntax)
285{
286 print_raw (os, pr_as_read_syntax);
287 newline (os);
288}
289
290template <typename T>
291void
293 const std::string& prefix) const
294{
295 matrix.print_info (os, prefix);
296}
297
298template <typename T>
299void
301 bool pr_as_read_syntax) const
302{
303 octave::preserve_stream_state stream_state (os);
304
305 octave_idx_type nr = matrix.rows ();
306 octave_idx_type nc = matrix.cols ();
307 octave_idx_type nz = nnz ();
308
309 // FIXME: this should probably all be handled by a
310 // separate octave_print_internal function that can handle format
311 // compact, loose, etc.
312
313 os << "Compressed Column Sparse (rows = " << nr
314 << ", cols = " << nc
315 << ", nnz = " << nz;
316
317 // Avoid calling numel here since it can easily overflow
318 // octave_idx_type even when there is no real problem storing the
319 // sparse array.
320
321 double dnr = nr;
322 double dnc = nc;
323 double dnel = dnr * dnc;
324
325 if (dnel > 0)
326 {
327 double pct = (nz / dnel * 100);
328
329 int prec = 2;
330
331 // Display at least 2 significant figures and up to 4 as we
332 // approach 100%. Avoid having limited precision of the display
333 // result in reporting 100% for matrices that are not actually
334 // 100% full.
335
336 if (pct == 100)
337 prec = 3;
338 else
339 {
340 if (pct > 99.9)
341 prec = 4;
342 else if (pct > 99)
343 prec = 3;
344
345 if (pct > 99.99)
346 pct = 99.99;
347 }
348
349 os << " [" << std::setprecision (prec) << pct << "%]";
350 }
351
352 os << ")\n";
353
354 // add one to the printed indices to go from
355 // zero-based to one-based arrays
356
357 if (nz != 0)
358 {
359 for (octave_idx_type j = 0; j < nc; j++)
360 {
361 octave_quit ();
362
363 // FIXME: is there an easy way to get the max row
364 // and column indices so we can set the width appropriately
365 // and line up the columns here? Similarly, we should look
366 // at all the nonzero values and display them with the same
367 // formatting rules that apply to columns of a matrix.
368
369 for (octave_idx_type i = matrix.cidx (j); i < matrix.cidx (j+1); i++)
370 {
371 os << "\n";
372 os << " (" << matrix.ridx (i)+1 << ", " << j+1 << ") -> ";
373
374 octave_print_internal (os, matrix.data (i), pr_as_read_syntax);
375 }
376 }
377 }
378}
379
380template <typename MT>
383{
384 return float_display_format ();
385 // return make_format (this->matrix);
386}
387
388template <typename MT>
389std::string
392 octave_idx_type j) const
393{
394 std::ostringstream buf;
395 octave_print_internal (buf, fmt, this->matrix(i, j));
396 return buf.str ();
397}
398
399template <typename T>
400bool
402{
403 dim_vector dv = this->dims ();
404
405 // Ensure that additional memory is deallocated
406 matrix.maybe_compress ();
407
408 os << "# nnz: " << nnz () << "\n";
409 os << "# rows: " << dv(0) << "\n";
410 os << "# columns: " << dv(1) << "\n";
411
412 os << this->matrix;
413
414 return true;
415}
416
417template <typename T>
418bool
420{
421 octave_idx_type nz = 0;
422 octave_idx_type nr = 0;
423 octave_idx_type nc = 0;
424
425 if (! extract_keyword (is, "nnz", nz, true)
426 || ! extract_keyword (is, "rows", nr, true)
427 || ! extract_keyword (is, "columns", nc, true))
428 error ("load: failed to extract number of rows and columns");
429
430 T tmp (nr, nc, nz);
431
432 is >> tmp;
433
434 if (! is)
435 error ("load: failed to load matrix constant");
436
437 matrix = tmp;
438
439 return true;
440}
441
442template <typename T>
445{
446 octave_idx_type nr = matrix.rows ();
447 octave_idx_type nc = matrix.cols ();
448
449 octave_idx_type i = n % nr;
450 octave_idx_type j = n / nr;
451
452 return (i < nr && j < nc) ? octave_value (matrix(i, j)) : octave_value ();
453}
454
455template <typename T>
458{
459 if (umap == umap_xtolower || umap == umap_xtoupper)
460 return matrix;
461
462 // Try the map on the dense value.
463 // FIXME: We should probably be smarter about this, especially for the
464 // cases that are expected to return sparse matrices.
465 octave_value retval = this->full_value ().map (umap);
466
467 // Sparsify the result if possible.
468
469 switch (umap)
470 {
471 case umap_xisalnum:
472 case umap_xisalpha:
473 case umap_xisascii:
474 case umap_xiscntrl:
475 case umap_xisdigit:
476 case umap_xisgraph:
477 case umap_xislower:
478 case umap_xisprint:
479 case umap_xispunct:
480 case umap_xisspace:
481 case umap_xisupper:
482 case umap_xisxdigit:
483 // FIXME: intentionally skip this step for string mappers.
484 // Is this wanted?
485 break;
486
487 default:
488 {
489 switch (retval.builtin_type ())
490 {
491 case btyp_double:
492 retval = retval.sparse_matrix_value ();
493 break;
494
495 case btyp_complex:
496 retval = retval.sparse_complex_matrix_value ();
497 break;
498
499 case btyp_bool:
500 retval = retval.sparse_bool_matrix_value ();
501 break;
502
503 default:
504 break;
505 }
506 }
507 }
508
509 return retval;
510}
OCTAVE_API SparseBoolMatrix all(int dim=-1) const
Definition: boolSparse.cc:140
Vector representing the dimensions (size) of an Array.
Definition: dim-vector.h:94
octave_idx_type numel(int n=0) const
Number of elements that a matrix with this dimensions would have.
Definition: dim-vector.h:335
bool all_ones(void) const
Definition: dim-vector.h:324
bool any_zero(void) const
Definition: dim-vector.h:316
void set_pos_if_unset(octave_idx_type nd_arg, octave_idx_type dim_arg)
OCTINTERP_API bool print_as_scalar(void) const
OCTINTERP_API float_display_format get_edit_display_format(void) const
OCTINTERP_API octave_value subsref(const std::string &type, const std::list< octave_value_list > &idx)
OCTINTERP_API void print(std::ostream &os, bool pr_as_read_syntax=false)
OCTINTERP_API bool load_ascii(std::istream &is)
OCTINTERP_API void print_raw(std::ostream &os, bool pr_as_read_syntax=false) const
OCTINTERP_API bool save_ascii(std::ostream &os)
OCTINTERP_API octave_value fast_elem_extract(octave_idx_type n) const
OCTINTERP_API octave_value subsasgn(const std::string &type, const std::list< octave_value_list > &idx, const octave_value &rhs)
OCTINTERP_API void delete_elements(const octave_value_list &idx)
OCTINTERP_API octave_value map(octave_base_value::unary_mapper_t umap) const
OCTINTERP_API bool is_true(void) const
OCTINTERP_API octave_value resize(const dim_vector &dv, bool=false) const
OCTINTERP_API void print_info(std::ostream &os, const std::string &prefix) const
OCTINTERP_API octave_value do_index_op(const octave_value_list &idx, bool resize_ok=false)
OCTINTERP_API std::string edit_display(const float_display_format &fmt, octave_idx_type i, octave_idx_type j) const
octave_idx_type length(void) const
Definition: ovl.h:113
SparseMatrix sparse_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:945
builtin_type_t builtin_type(void) const
Definition: ov.h:735
static OCTINTERP_API octave_value empty_conv(const std::string &type, const octave_value &rhs=octave_value())
OCTINTERP_API octave_value next_subsref(const std::string &type, const std::list< octave_value_list > &idx, std::size_t skip=1)
SparseBoolMatrix sparse_bool_matrix_value(bool warn=false) const
Definition: ov.h:952
OCTINTERP_API octave_value subsasgn(const std::string &type, const std::list< octave_value_list > &idx, const octave_value &rhs)
SparseComplexMatrix sparse_complex_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:949
void error(const char *fmt,...)
Definition: error.cc:980
#define panic_impossible()
Definition: error.h:411
void warn_array_as_logical(const dim_vector &dv)
Definition: errwarn.cc:286
std::string extract_keyword(std::istream &is, const char *keyword, const bool next_only)
Definition: ls-oct-text.cc:84
void err_nan_to_logical_conversion(void)
return octave_value(v1.char_array_value() . concat(v2.char_array_value(), ra_idx),((a1.is_sq_string()||a2.is_sq_string()) ? '\'' :'"'))
@ btyp_double
Definition: ov-base.h:76
@ btyp_bool
Definition: ov-base.h:88
@ btyp_complex
Definition: ov-base.h:78
void octave_print_internal(std::ostream &os, const float_display_format &fmt, bool d, bool pr_as_read_syntax)
Definition: pr-output.cc:1762
F77_RET_T len
Definition: xerbla.cc:61