GNU Octave  4.2.1
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-base-sparse.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 2004-2017 David Bateman
4 Copyright (C) 1998-2004 Andy Adler
5 Copyright (C) 2010 VZLU Prague
6 
7 This file is part of Octave.
8 
9 Octave is free software; you can redistribute it and/or modify it
10 under the terms of the GNU General Public License as published by the
11 Free Software Foundation; either version 3 of the License, or (at your
12 option) any later version.
13 
14 Octave is distributed in the hope that it will be useful, but WITHOUT
15 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
16 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
17 for more details.
18 
19 You should have received a copy of the GNU General Public License
20 along with Octave; see the file COPYING. If not, see
21 <http://www.gnu.org/licenses/>.
22 
23 */
24 
25 // This file should not include config.h. It is only included in other
26 // C++ source files that should have included config.h before including
27 // this file.
28 
29 #include <iomanip>
30 #include <iostream>
31 
32 #include "ovl.h"
33 #include "ov-base.h"
34 #include "quit.h"
35 #include "pr-output.h"
36 
37 #include "byte-swap.h"
38 #include "ls-oct-text.h"
39 #include "ls-utils.h"
40 #include "ls-hdf5.h"
41 
42 #include "boolSparse.h"
43 #include "ov-base-sparse.h"
45 #include "pager.h"
46 #include "utils.h"
47 
48 #include "lo-array-errwarn.h"
49 
50 template <typename T>
53  bool resize_ok)
54 {
56 
57  octave_idx_type n_idx = idx.length ();
58 
59  // If we catch an indexing error in index_vector, we flag an error in
60  // index k. Ensure it is the right value befor each idx_vector call.
61  // Same variable as used in the for loop in the default case.
62 
63  octave_idx_type k = 0;
64 
65  try
66  {
67  switch (n_idx)
68  {
69  case 0:
70  retval = matrix;
71  break;
72 
73  case 1:
74  {
75  idx_vector i = idx (0).index_vector ();
76 
77  retval = octave_value (matrix.index (i, resize_ok));
78  }
79  break;
80 
81  case 2:
82  {
83  idx_vector i = idx (0).index_vector ();
84 
85  k = 1;
86  idx_vector j = idx (1).index_vector ();
87 
88  retval = octave_value (matrix.index (i, j, resize_ok));
89  }
90  break;
91 
92  default:
93  error ("sparse indexing needs 1 or 2 indices");
94  }
95  }
96  catch (octave::index_exception& e)
97  {
98  // Rethrow to allow more info to be reported later.
99  e.set_pos_if_unset (n_idx, k+1);
100  throw;
101  }
102 
103  return retval;
104 }
105 
106 template <typename T>
109  const std::list<octave_value_list>& idx)
110 {
112 
113  switch (type[0])
114  {
115  case '(':
116  retval = do_index_op (idx.front ());
117  break;
118 
119  case '{':
120  case '.':
121  {
122  std::string nm = type_name ();
123  error ("%s cannot be indexed with %c", nm.c_str (), type[0]);
124  }
125  break;
126 
127  default:
128  panic_impossible ();
129  }
130 
131  return retval.next_subsref (type, idx);
132 }
133 
134 template <typename T>
137  const std::list<octave_value_list>& idx,
138  const octave_value& rhs)
139 {
141 
142  switch (type[0])
143  {
144  case '(':
145  {
146  if (type.length () != 1)
147  {
148  std::string nm = type_name ();
149  error ("in indexed assignment of %s, last lhs index must be ()",
150  nm.c_str ());
151  }
152 
153  retval = numeric_assign (type, idx, rhs);
154  }
155  break;
156 
157  case '{':
158  case '.':
159  {
160  if (! is_empty ())
161  {
162  std::string nm = type_name ();
163  error ("%s cannot be indexed with %c", nm.c_str (), type[0]);
164  }
165 
167 
168  retval = tmp.subsasgn (type, idx, rhs);
169  }
170  break;
171 
172  default:
173  panic_impossible ();
174  }
175 
176  return retval;
177 }
178 
179 template <typename T>
180 void
182 {
183 
184  octave_idx_type len = idx.length ();
185 
186  // If we catch an indexing error in index_vector, we flag an error in
187  // index k. Ensure it is the right value befor each idx_vector call.
188  // Same variable as used in the for loop in the default case.
189 
190  octave_idx_type k = 0;
191 
192  try
193  {
194  switch (len)
195  {
196  case 1:
197  {
198  idx_vector i = idx (0).index_vector ();
199 
200  matrix.assign (i, rhs);
201 
202  break;
203  }
204 
205  case 2:
206  {
207  idx_vector i = idx (0).index_vector ();
208 
209  k = 1;
210  idx_vector j = idx (1).index_vector ();
211 
212  matrix.assign (i, j, rhs);
213 
214  break;
215  }
216 
217  default:
218  error ("sparse indexing needs 1 or 2 indices");
219  }
220  }
221  catch (octave::index_exception& e)
222  {
223  // Rethrow to allow more info to be reported later.
224  e.set_pos_if_unset (len, k+1);
225  throw;
226  }
227 
228  // Invalidate matrix type.
229  typ.invalidate_type ();
230 }
231 
232 template <typename MT>
233 void
235 {
236  octave_idx_type len = idx.length ();
237 
238  // If we catch an indexing error in index_vector, we flag an error in
239  // index k. Ensure it is the right value befor each idx_vector call.
240  // Same variable as used in the for loop in the default case.
241 
242  octave_idx_type k = 0;
243 
244  try
245  {
246  switch (len)
247  {
248  case 1:
249  {
250  idx_vector i = idx (0).index_vector ();
251 
252  matrix.delete_elements (i);
253 
254  break;
255  }
256 
257  case 2:
258  {
259  idx_vector i = idx (0).index_vector ();
260 
261  k = 1;
262  idx_vector j = idx (1).index_vector ();
263 
264  matrix.delete_elements (i, j);
265 
266  break;
267  }
268 
269  default:
270  error ("sparse indexing needs 1 or 2 indices");
271  }
272  }
273  catch (octave::index_exception& e)
274  {
275  // Rethrow to allow more info to be reported later.
276  e.set_pos_if_unset (len, k+1);
277  throw;
278  }
279 
280  // Invalidate the matrix type
281  typ.invalidate_type ();
282 }
283 
284 template <typename T>
287 {
288  T retval (matrix);
289  retval.resize (dv);
290  return retval;
291 }
292 
293 template <typename T>
294 bool
296 {
297  bool retval = false;
298  dim_vector dv = matrix.dims ();
299  octave_idx_type nel = dv.numel ();
300  octave_idx_type nz = nnz ();
301 
302  if (nel > 0)
303  {
304  T t1 (matrix.reshape (dim_vector (nel, 1)));
305 
306  if (t1.any_element_is_nan ())
308 
309  if (nel > 1)
311 
312  if (nz == nel)
313  {
314  SparseBoolMatrix t2 = t1.all ();
315 
316  retval = t2(0);
317  }
318  }
319 
320  return retval;
321 }
322 
323 template <typename T>
324 bool
326 {
327  dim_vector dv = dims ();
328 
329  return (dv.all_ones () || dv.any_zero ());
330 }
331 
332 template <typename T>
333 void
334 octave_base_sparse<T>::print (std::ostream& os, bool pr_as_read_syntax)
335 {
336  print_raw (os, pr_as_read_syntax);
337  newline (os);
338 }
339 
340 template <typename T>
341 void
343  const std::string& prefix) const
344 {
345  matrix.print_info (os, prefix);
346 }
347 
348 template <typename T>
349 void
351  bool pr_as_read_syntax) const
352 {
353  octave_preserve_stream_state stream_state (os);
354 
355  octave_idx_type nr = matrix.rows ();
356  octave_idx_type nc = matrix.cols ();
357  octave_idx_type nz = nnz ();
358 
359  // FIXME: this should probably all be handled by a
360  // separate octave_print_internal function that can handle format
361  // compact, loose, etc.
362 
363  os << "Compressed Column Sparse (rows = " << nr
364  << ", cols = " << nc
365  << ", nnz = " << nz;
366 
367  // Avoid calling numel here since it can easily overflow
368  // octave_idx_type even when there is no real problem storing the
369  // sparse array.
370 
371  double dnr = nr;
372  double dnc = nc;
373  double dnel = dnr * dnc;
374 
375  if (dnel > 0)
376  {
377  double pct = (nz / dnel * 100);
378 
379  int prec = 2;
380 
381  // Display at least 2 significant figures and up to 4 as we
382  // approach 100%. Avoid having limited precision of the display
383  // result in reporting 100% for matrices that are not actually
384  // 100% full.
385 
386  if (pct == 100)
387  prec = 3;
388  else
389  {
390  if (pct > 99.9)
391  prec = 4;
392  else if (pct > 99)
393  prec = 3;
394 
395  if (pct > 99.99)
396  pct = 99.99;
397  }
398 
399  os << " [" << std::setprecision (prec) << pct << "%]";
400  }
401 
402  os << ")\n";
403 
404  // add one to the printed indices to go from
405  // zero-based to one-based arrays
406 
407  if (nz != 0)
408  {
409  for (octave_idx_type j = 0; j < nc; j++)
410  {
411  octave_quit ();
412 
413  // FIXME: is there an easy way to get the max row
414  // and column indices so we can set the width appropriately
415  // and line up the columns here? Similarly, we should look
416  // at all the nonzero values and display them with the same
417  // formatting rules that apply to columns of a matrix.
418 
419  for (octave_idx_type i = matrix.cidx (j); i < matrix.cidx (j+1); i++)
420  {
421  os << "\n";
422  os << " (" << matrix.ridx (i)+1 << ", " << j+1 << ") -> ";
423 
424  octave_print_internal (os, matrix.data (i), pr_as_read_syntax);
425  }
426  }
427  }
428 }
429 
430 template <typename T>
431 bool
433 {
434  dim_vector dv = this->dims ();
435 
436  // Ensure that additional memory is deallocated
437  matrix.maybe_compress ();
438 
439  os << "# nnz: " << nnz () << "\n";
440  os << "# rows: " << dv(0) << "\n";
441  os << "# columns: " << dv(1) << "\n";
442 
443  os << this->matrix;
444 
445  return true;
446 }
447 
448 template <typename T>
449 bool
451 {
452  octave_idx_type nz = 0;
453  octave_idx_type nr = 0;
454  octave_idx_type nc = 0;
455 
456  if (! extract_keyword (is, "nnz", nz, true)
457  || ! extract_keyword (is, "rows", nr, true)
458  || ! extract_keyword (is, "columns", nc, true))
459  error ("load: failed to extract number of rows and columns");
460 
461  T tmp (nr, nc, nz);
462 
463  is >> tmp;
464 
465  if (! is)
466  error ("load: failed to load matrix constant");
467 
468  matrix = tmp;
469 
470  return true;
471 }
472 
473 template <typename T>
476 {
477  octave_idx_type nr = matrix.rows ();
478  octave_idx_type nc = matrix.cols ();
479 
480  octave_idx_type i = n % nr;
481  octave_idx_type j = n / nr;
482 
483  return (i < nr && j < nc) ? octave_value (matrix(i,j)) : octave_value ();
484 }
485 
486 template <typename T>
489 {
490  if (umap == umap_xtolower || umap == umap_xtoupper)
491  return matrix;
492 
493  // Try the map on the dense value.
494  // FIXME: We should probably be smarter about this, especially for the
495  // cases that are expected to return sparse matrices.
496  octave_value retval = this->full_value ().map (umap);
497 
498  // Sparsify the result if possible.
499 
500  switch (umap)
501  {
502  case umap_xisalnum:
503  case umap_xisalpha:
504  case umap_xisascii:
505  case umap_xiscntrl:
506  case umap_xisdigit:
507  case umap_xisgraph:
508  case umap_xislower:
509  case umap_xisprint:
510  case umap_xispunct:
511  case umap_xisspace:
512  case umap_xisupper:
513  case umap_xisxdigit:
514  case umap_xtoascii:
515  // FIXME: intentionally skip this step for string mappers.
516  // Is this wanted?
517  break;
518 
519  default:
520  {
521  switch (retval.builtin_type ())
522  {
523  case btyp_double:
524  retval = retval.sparse_matrix_value ();
525  break;
526 
527  case btyp_complex:
528  retval = retval.sparse_complex_matrix_value ();
529  break;
530 
531  case btyp_bool:
532  retval = retval.sparse_bool_matrix_value ();
533  break;
534 
535  default:
536  break;
537  }
538  }
539  }
540 
541  return retval;
542 }
octave_value subsref(const std::string &type, const std::list< octave_value_list > &idx)
octave_value subsasgn(const std::string &type, const std::list< octave_value_list > &idx, const octave_value &rhs)
void set_pos_if_unset(octave_idx_type nd_arg, octave_idx_type dim_arg)
octave_idx_type length(void) const
Definition: ovl.h:96
bool load_ascii(std::istream &is)
for large enough k
Definition: lu.cc:606
the sparsity preserving column transformation such that that defines the pivoting threshold can be given in which case it defines the then the first element defines the pivoting tolerance for the unsymmetric the values defined such that for full matrix
Definition: lu.cc:138
void error(const char *fmt,...)
Definition: error.cc:570
octave_value map(octave_base_value::unary_mapper_t umap) const
Definition: ov.h:1413
bool print_as_scalar(void) const
bool all_ones(void) const
Definition: dim-vector.h:377
void err_nan_to_logical_conversion(void)
i e
Definition: data.cc:2724
octave_idx_type numel(int n=0) const
Number of elements that a matrix with this dimensions would have.
Definition: dim-vector.h:389
octave_value subsasgn(const std::string &type, const std::list< octave_value_list > &idx, const octave_value &rhs)
Definition: ov.cc:1540
void assign(const octave_value_list &idx, const T &rhs)
OCTAVE_EXPORT octave_value_list isdir nd deftypefn *std::string nm
Definition: utils.cc:941
SparseBoolMatrix all(int dim=-1) const
Definition: boolSparse.cc:138
std::string extract_keyword(std::istream &is, const char *keyword, const bool next_only)
Definition: ls-oct-text.cc:80
octave_value map(octave_base_value::unary_mapper_t umap) const
SparseBoolMatrix sparse_bool_matrix_value(bool warn=false) const
Definition: ov.h:841
double tmp
Definition: data.cc:6300
octave_value retval
Definition: data.cc:6294
#define panic_impossible()
Definition: error.h:40
void print_raw(std::ostream &os, bool pr_as_read_syntax=false) const
idx type
Definition: ov.cc:3129
void print(std::ostream &os, bool pr_as_read_syntax=false)
void print_info(std::ostream &os, const std::string &prefix) const
bool any_zero(void) const
Definition: dim-vector.h:359
SparseComplexMatrix sparse_complex_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:838
the exceeded dimensions are set to if fewer subscripts than dimensions are the exceeding dimensions are merged into the final requested dimension For consider the following dims
Definition: sub2ind.cc:255
bool save_ascii(std::ostream &os)
void delete_elements(const octave_value_list &idx)
octave_value fast_elem_extract(octave_idx_type n) const
octave_value resize(const dim_vector &dv, bool=false) const
=val(i)}if ode{val(i)}occurs in table i
Definition: lookup.cc:239
void octave_print_internal(std::ostream &, char, bool)
Definition: pr-output.cc:1722
octave_value do_index_op(const octave_value_list &idx, bool resize_ok=false)
void warn_array_as_logical(const dim_vector &dv)
Definition: errwarn.cc:276
SparseMatrix sparse_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:834
write the output to stdout if nargout is
Definition: load-save.cc:1576
Vector representing the dimensions (size) of an Array.
Definition: dim-vector.h:87
If this string is the system will ring the terminal sometimes it is useful to be able to print the original representation of the string
Definition: utils.cc:854
builtin_type_t builtin_type(void) const
Definition: ov.h:619
dim_vector dv
Definition: sub2ind.cc:263
bool is_true(void) const
octave_value next_subsref(const std::string &type, const std::list< octave_value_list > &idx, size_t skip=1)
Definition: ov.cc:1462
return octave_value(v1.char_array_value().concat(v2.char_array_value(), ra_idx),((a1.is_sq_string()||a2.is_sq_string())? '\'': '"'))
static octave_value empty_conv(const std::string &type, const octave_value &rhs=octave_value())
Definition: ov.cc:2882