GNU Octave  9.1.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
dSparse.h
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////
2 //
3 // Copyright (C) 1998-2024 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 #if ! defined (octave_dSparse_h)
27 #define octave_dSparse_h 1
28 
29 #include "octave-config.h"
30 
31 #include "mx-fwd.h"
32 
33 #include "CColVector.h"
34 #include "CMatrix.h"
35 #include "DET.h"
36 #include "MSparse.h"
37 #include "MatrixType.h"
38 #include "Sparse-op-decls.h"
39 #include "dColVector.h"
40 #include "dMatrix.h"
41 #include "dNDArray.h"
42 
43 class
45 {
46 public:
47 
48  // Corresponding dense matrix type for this sparse matrix type.
50 
51  typedef void (*solve_singularity_handler) (double rcond);
52 
53  SparseMatrix () : MSparse<double> () { }
54 
56  : MSparse<double> (r, c) { }
57 
58  SparseMatrix (const dim_vector& dv, octave_idx_type nz = 0) :
59  MSparse<double> (dv, nz) { }
60 
61  explicit SparseMatrix (octave_idx_type r, octave_idx_type c, double val)
62  : MSparse<double> (r, c, val) { }
63 
64  SparseMatrix (const SparseMatrix& a) : MSparse<double> (a) { }
65 
66  SparseMatrix (const SparseMatrix& a, const dim_vector& dv)
67  : MSparse<double> (a, dv) { }
68 
69  SparseMatrix (const MSparse<double>& a) : MSparse<double> (a) { }
70 
71  SparseMatrix (const Sparse<double>& a) : MSparse<double> (a) { }
72 
73  explicit OCTAVE_API SparseMatrix (const SparseBoolMatrix& a);
74 
75  explicit SparseMatrix (const Matrix& a) : MSparse<double> (a) { }
76 
77  explicit SparseMatrix (const NDArray& a) : MSparse<double> (a) { }
78 
80  const octave::idx_vector& c, octave_idx_type nr = -1,
81  octave_idx_type nc = -1, bool sum_terms = true,
82  octave_idx_type nzm = -1)
83  : MSparse<double> (a, r, c, nr, nc, sum_terms, nzm) { }
84 
85  explicit OCTAVE_API SparseMatrix (const DiagMatrix& a);
86 
87  explicit SparseMatrix (const PermMatrix& a) : MSparse<double> (a) { }
88 
90  octave_idx_type num_nz) : MSparse<double> (r, c, num_nz) { }
91 
92  SparseMatrix& operator = (const SparseMatrix& a)
93  {
95  return *this;
96  }
97 
98  ~SparseMatrix () = default;
99 
100  OCTAVE_API bool operator == (const SparseMatrix& a) const;
101  OCTAVE_API bool operator != (const SparseMatrix& a) const;
102 
103  OCTAVE_API bool issymmetric () const;
104 
105  OCTAVE_API SparseMatrix max (int dim = -1) const;
106  OCTAVE_API SparseMatrix max (Array<octave_idx_type>& index, int dim = -1) const;
107  OCTAVE_API SparseMatrix min (int dim = -1) const;
108  OCTAVE_API SparseMatrix min (Array<octave_idx_type>& index, int dim = -1) const;
109 
110  // destructive insert/delete/reorder operations
111 
113  insert (const SparseMatrix& a, octave_idx_type r, octave_idx_type c);
114 
116  insert (const SparseMatrix& a, const Array<octave_idx_type>& indx);
117 
119  concat (const SparseMatrix& rb, const Array<octave_idx_type>& ra_idx);
122 
125 
127  {
128  return MSparse<double>::transpose ();
129  }
130  SparseMatrix hermitian () const { return transpose (); }
131 
132  // extract row or column i.
133 
134  OCTAVE_API RowVector row (octave_idx_type i) const;
135 
136  OCTAVE_API ColumnVector column (octave_idx_type i) const;
137 
138 private:
140  dinverse (MatrixType& mattype, octave_idx_type& info, double& rcond,
141  const bool force = false, const bool calccond = true) const;
142 
144  tinverse (MatrixType& mattype, octave_idx_type& info, double& rcond,
145  const bool force = false, const bool calccond = true) const;
146 
147 public:
148  OCTAVE_API SparseMatrix inverse () const;
149  OCTAVE_API SparseMatrix inverse (MatrixType& mattype) const;
151  inverse (MatrixType& mattype, octave_idx_type& info) const;
153  inverse (MatrixType& mattype, octave_idx_type& info, double& rcond,
154  bool force = false, bool calc_cond = true) const;
155 
156  OCTAVE_API DET determinant () const;
157  OCTAVE_API DET determinant (octave_idx_type& info) const;
158  OCTAVE_API DET determinant (octave_idx_type& info, double& rcond,
159  bool calc_cond = true) const;
160 
161 private:
162  // Diagonal matrix solvers
164  dsolve (MatrixType& typ, const Matrix& b, octave_idx_type& info,
165  double& rcond, solve_singularity_handler sing_handler,
166  bool calc_cond = false) const;
167 
169  dsolve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info,
170  double& rcond, solve_singularity_handler sing_handler,
171  bool calc_cond = false) const;
172 
174  dsolve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info,
175  double& rcond, solve_singularity_handler sing_handler,
176  bool calc_cond = false) const;
177 
179  dsolve (MatrixType& typ, const SparseComplexMatrix& b, octave_idx_type& info,
180  double& rcond, solve_singularity_handler sing_handler,
181  bool calc_cond = false) const;
182 
183  // Upper triangular matrix solvers
185  utsolve (MatrixType& typ, const Matrix& b, octave_idx_type& info,
186  double& rcond, solve_singularity_handler sing_handler,
187  bool calc_cond = false) const;
188 
190  utsolve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info,
191  double& rcond, solve_singularity_handler sing_handler,
192  bool calc_cond = false) const;
193 
195  utsolve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info,
196  double& rcond, solve_singularity_handler sing_handler,
197  bool calc_cond = false) const;
198 
200  utsolve (MatrixType& typ, const SparseComplexMatrix& b, octave_idx_type& info,
201  double& rcond, solve_singularity_handler sing_handler,
202  bool calc_cond = false) const;
203 
204  // Lower triangular matrix solvers
206  ltsolve (MatrixType& typ, const Matrix& b, octave_idx_type& info,
207  double& rcond, solve_singularity_handler sing_handler,
208  bool calc_cond = false) const;
209 
211  ltsolve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info,
212  double& rcond, solve_singularity_handler sing_handler,
213  bool calc_cond = false) const;
214 
216  ltsolve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info,
217  double& rcond, solve_singularity_handler sing_handler,
218  bool calc_cond = false) const;
219 
221  ltsolve (MatrixType& typ, const SparseComplexMatrix& b, octave_idx_type& info,
222  double& rcond, solve_singularity_handler sing_handler,
223  bool calc_cond = false) const;
224 
225  // Tridiagonal matrix solvers
227  trisolve (MatrixType& typ, const Matrix& b, octave_idx_type& info,
228  double& rcond, solve_singularity_handler sing_handler,
229  bool calc_cond = false) const;
230 
232  trisolve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info,
233  double& rcond, solve_singularity_handler sing_handler,
234  bool calc_cond = false) const;
235 
237  trisolve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info,
238  double& rcond, solve_singularity_handler sing_handler,
239  bool calc_cond = false) const;
240 
242  trisolve (MatrixType& typ, const SparseComplexMatrix& b, octave_idx_type& info,
243  double& rcond, solve_singularity_handler sing_handler,
244  bool calc_cond = false) const;
245 
246  // Banded matrix solvers (umfpack/cholesky)
248  bsolve (MatrixType& typ, const Matrix& b, octave_idx_type& info,
249  double& rcond, solve_singularity_handler sing_handler,
250  bool calc_cond = false) const;
251 
253  bsolve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info,
254  double& rcond, solve_singularity_handler sing_handler,
255  bool calc_cond = false) const;
256 
258  bsolve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info,
259  double& rcond, solve_singularity_handler sing_handler,
260  bool calc_cond = false) const;
261 
263  bsolve (MatrixType& typ, const SparseComplexMatrix& b, octave_idx_type& info,
264  double& rcond, solve_singularity_handler sing_handler,
265  bool calc_cond = false) const;
266 
267  // Full matrix solvers (umfpack/cholesky)
268  OCTAVE_API void *
269  factorize (octave_idx_type& err, double& rcond, Matrix& Control,
270  Matrix& Info, solve_singularity_handler sing_handler,
271  bool calc_cond = false) const;
272 
274  fsolve (MatrixType& typ, const Matrix& b, octave_idx_type& info,
275  double& rcond, solve_singularity_handler sing_handler,
276  bool calc_cond = false) const;
277 
279  fsolve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info,
280  double& rcond, solve_singularity_handler sing_handler,
281  bool calc_cond = false) const;
282 
284  fsolve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info,
285  double& rcond, solve_singularity_handler sing_handler,
286  bool calc_cond = false) const;
287 
289  fsolve (MatrixType& typ, const SparseComplexMatrix& b, octave_idx_type& info,
290  double& rcond, solve_singularity_handler sing_handler,
291  bool calc_cond = false) const;
292 
293 public:
294  // Generic interface to solver with no probing of type
295  OCTAVE_API Matrix solve (MatrixType& typ, const Matrix& b) const;
297  solve (MatrixType& typ, const Matrix& b, octave_idx_type& info) const;
299  solve (MatrixType& typ, const Matrix& b, octave_idx_type& info,
300  double& rcond) const;
302  solve (MatrixType& typ, const Matrix& b, octave_idx_type& info,
303  double& rcond, solve_singularity_handler sing_handler,
304  bool singular_fallback = true) const;
305 
307  solve (MatrixType& typ, const ComplexMatrix& b) const;
309  solve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info) const;
311  solve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info,
312  double& rcond) const;
314  solve (MatrixType& typ, const ComplexMatrix& b, octave_idx_type& info,
315  double& rcond, solve_singularity_handler sing_handler,
316  bool singular_fallback = true) const;
317 
318  OCTAVE_API SparseMatrix solve (MatrixType& typ, const SparseMatrix& b) const;
320  solve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info) const;
322  solve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info,
323  double& rcond) const;
325  solve (MatrixType& typ, const SparseMatrix& b, octave_idx_type& info,
326  double& rcond, solve_singularity_handler sing_handler,
327  bool singular_fallback = true) const;
328 
330  solve (MatrixType& typ, const SparseComplexMatrix& b) const;
332  solve (MatrixType& typ, const SparseComplexMatrix& b,
333  octave_idx_type& info) const;
335  solve (MatrixType& typ, const SparseComplexMatrix& b,
336  octave_idx_type& info, double& rcond) const;
338  solve (MatrixType& typ, const SparseComplexMatrix& b,
339  octave_idx_type& info, double& rcond,
340  solve_singularity_handler sing_handler,
341  bool singular_fallabck = true) const;
342 
343  OCTAVE_API ColumnVector solve (MatrixType& typ, const ColumnVector& b) const;
345  solve (MatrixType& typ, const ColumnVector& b, octave_idx_type& info) const;
347  solve (MatrixType& typ, const ColumnVector& b, octave_idx_type& info,
348  double& rcond) const;
350  solve (MatrixType& typ, const ColumnVector& b, octave_idx_type& info,
351  double& rcond, solve_singularity_handler sing_handler) const;
352 
354  solve (MatrixType& typ, const ComplexColumnVector& b) const;
356  solve (MatrixType& typ, const ComplexColumnVector& b,
357  octave_idx_type& info) const;
359  solve (MatrixType& typ, const ComplexColumnVector& b,
360  octave_idx_type& info, double& rcond) const;
362  solve (MatrixType& typ, const ComplexColumnVector& b,
363  octave_idx_type& info, double& rcond,
364  solve_singularity_handler sing_handler) const;
365 
366  // Generic interface to solver with probing of type
367  OCTAVE_API Matrix solve (const Matrix& b) const;
368  OCTAVE_API Matrix solve (const Matrix& b, octave_idx_type& info) const;
370  solve (const Matrix& b, octave_idx_type& info, double& rcond) const;
372  solve (const Matrix& b, octave_idx_type& info, double& rcond,
373  solve_singularity_handler sing_handler) const;
374 
377  solve (const ComplexMatrix& b, octave_idx_type& info) const;
379  solve (const ComplexMatrix& b, octave_idx_type& info, double& rcond) const;
381  solve (const ComplexMatrix& b, octave_idx_type& info, double& rcond,
382  solve_singularity_handler sing_handler) const;
383 
384  OCTAVE_API SparseMatrix solve (const SparseMatrix& b) const;
386  solve (const SparseMatrix& b, octave_idx_type& info) const;
388  solve (const SparseMatrix& b, octave_idx_type& info, double& rcond) const;
390  solve (const SparseMatrix& b, octave_idx_type& info, double& rcond,
391  solve_singularity_handler sing_handler) const;
392 
393  OCTAVE_API SparseComplexMatrix solve (const SparseComplexMatrix& b) const;
395  solve (const SparseComplexMatrix& b, octave_idx_type& info) const;
397  solve (const SparseComplexMatrix& b, octave_idx_type& info,
398  double& rcond) const;
400  solve (const SparseComplexMatrix& b, octave_idx_type& info,
401  double& rcond, solve_singularity_handler sing_handler) const;
402 
403  OCTAVE_API ColumnVector solve (const ColumnVector& b) const;
405  solve (const ColumnVector& b, octave_idx_type& info) const;
407  solve (const ColumnVector& b, octave_idx_type& info, double& rcond) const;
409  solve (const ColumnVector& b, octave_idx_type& info, double& rcond,
410  solve_singularity_handler sing_handler) const;
411 
412  OCTAVE_API ComplexColumnVector solve (const ComplexColumnVector& b) const;
414  solve (const ComplexColumnVector& b, octave_idx_type& info) const;
416  solve (const ComplexColumnVector& b, octave_idx_type& info,
417  double& rcond) const;
419  solve (const ComplexColumnVector& b, octave_idx_type& info,
420  double& rcond, solve_singularity_handler sing_handler) const;
421 
422  // other operations
423 
424  OCTAVE_API bool any_element_is_negative (bool = false) const;
425  OCTAVE_API bool any_element_is_nan () const;
426  OCTAVE_API bool any_element_is_inf_or_nan () const;
427  OCTAVE_API bool any_element_not_one_or_zero () const;
428  OCTAVE_API bool all_elements_are_zero () const;
429  OCTAVE_API bool all_elements_are_int_or_inf_or_nan () const;
430  OCTAVE_API bool all_integers (double& max_val, double& min_val) const;
431  OCTAVE_API bool too_large_for_float () const;
432 
434 
435  OCTAVE_API SparseBoolMatrix all (int dim = -1) const;
436  OCTAVE_API SparseBoolMatrix any (int dim = -1) const;
437 
438  OCTAVE_API SparseMatrix cumprod (int dim = -1) const;
439  OCTAVE_API SparseMatrix cumsum (int dim = -1) const;
440  OCTAVE_API SparseMatrix prod (int dim = -1) const;
441  OCTAVE_API SparseMatrix sum (int dim = -1) const;
442  OCTAVE_API SparseMatrix sumsq (int dim = -1) const;
443  OCTAVE_API SparseMatrix abs () const;
444 
445  OCTAVE_API SparseMatrix diag (octave_idx_type k = 0) const;
446 
447  OCTAVE_API Matrix matrix_value () const;
448 
449  OCTAVE_API SparseMatrix squeeze () const;
450 
451  OCTAVE_API SparseMatrix reshape (const dim_vector& new_dims) const;
452 
454  permute (const Array<octave_idx_type>& vec, bool inv = false) const;
455 
456  OCTAVE_API SparseMatrix ipermute (const Array<octave_idx_type>& vec) const;
457 
458  // i/o
459 
460  friend OCTAVE_API std::ostream& operator << (std::ostream& os,
461  const SparseMatrix& a);
462  friend OCTAVE_API std::istream& operator >> (std::istream& is,
463  SparseMatrix& a);
464 
465 };
466 
467 // Publish externally used friend functions.
468 
471 
472 // Other operators.
473 
475  const SparseMatrix& b);
476 extern OCTAVE_API Matrix operator * (const Matrix& a,
477  const SparseMatrix& b);
478 extern OCTAVE_API Matrix mul_trans (const Matrix& a,
479  const SparseMatrix& b);
480 extern OCTAVE_API Matrix operator * (const SparseMatrix& a,
481  const Matrix& b);
482 extern OCTAVE_API Matrix trans_mul (const SparseMatrix& a,
483  const Matrix& b);
484 
486  const SparseMatrix&);
488  const DiagMatrix&);
489 
491  const SparseMatrix&);
493  const DiagMatrix&);
495  const SparseMatrix&);
497  const DiagMatrix&);
498 
500  const SparseMatrix&);
502  const PermMatrix&);
503 
504 extern OCTAVE_API SparseMatrix min (double d, const SparseMatrix& m);
505 extern OCTAVE_API SparseMatrix min (const SparseMatrix& m, double d);
506 extern OCTAVE_API SparseMatrix min (const SparseMatrix& a,
507  const SparseMatrix& b);
508 
509 extern OCTAVE_API SparseMatrix max (double d, const SparseMatrix& m);
510 extern OCTAVE_API SparseMatrix max (const SparseMatrix& m, double d);
511 extern OCTAVE_API SparseMatrix max (const SparseMatrix& a,
512  const SparseMatrix& b);
513 
516 
519 
522 
524 
525 #endif
template std::ostream & operator<<(std::ostream &, const Array< bool > &)
ComplexNDArray concat(NDArray &ra, ComplexNDArray &rb, const Array< octave_idx_type > &ra_idx)
Definition: CNDArray.cc:418
#define SPARSE_FORWARD_DEFS(B, R, F, T)
Definition: MSparse.h:194
#define SPARSE_SMSM_CMP_OP_DECLS(M1, M2, API)
#define SPARSE_SMSM_BOOL_OP_DECLS(M1, M2, API)
#define SPARSE_SSM_CMP_OP_DECLS(S, M, API)
#define SPARSE_SMS_BOOL_OP_DECLS(M, S, API)
#define SPARSE_SSM_BOOL_OP_DECLS(S, M, API)
#define SPARSE_SMS_CMP_OP_DECLS(M, S, API)
std::istream & operator>>(std::istream &is, SparseBoolMatrix &a)
Definition: boolSparse.cc:279
MSparse< T > transpose() const
Definition: MSparse.h:98
MSparse< T > & operator=(const MSparse< T > &a)
Definition: MSparse.h:80
Definition: dMatrix.h:42
SparseMatrix transpose() const
Definition: dSparse.h:126
SparseMatrix(const Sparse< double > &a)
Definition: dSparse.h:71
ComplexMatrix solve(const ComplexMatrix &b) const
SparseMatrix(const SparseMatrix &a, const dim_vector &dv)
Definition: dSparse.h:66
~SparseMatrix()=default
SparseMatrix(const SparseMatrix &a)
Definition: dSparse.h:64
SparseMatrix(octave_idx_type r, octave_idx_type c)
Definition: dSparse.h:55
SparseMatrix(const MSparse< double > &a)
Definition: dSparse.h:69
Matrix dense_matrix_type
Definition: dSparse.h:49
SparseMatrix(const Array< double > &a, const octave::idx_vector &r, const octave::idx_vector &c, octave_idx_type nr=-1, octave_idx_type nc=-1, bool sum_terms=true, octave_idx_type nzm=-1)
Definition: dSparse.h:79
SparseMatrix hermitian() const
Definition: dSparse.h:130
SparseMatrix(octave_idx_type r, octave_idx_type c, octave_idx_type num_nz)
Definition: dSparse.h:89
SparseMatrix()
Definition: dSparse.h:53
SparseMatrix(const NDArray &a)
Definition: dSparse.h:77
SparseMatrix(const PermMatrix &a)
Definition: dSparse.h:87
SparseMatrix(const Matrix &a)
Definition: dSparse.h:75
SparseMatrix(const dim_vector &dv, octave_idx_type nz=0)
Definition: dSparse.h:58
SparseMatrix(octave_idx_type r, octave_idx_type c, double val)
Definition: dSparse.h:61
Definition: Sparse.h:49
Definition: DET.h:39
Vector representing the dimensions (size) of an Array.
Definition: dim-vector.h:94
SparseMatrix operator+(const DiagMatrix &, const SparseMatrix &)
Definition: dSparse.cc:7555
Matrix trans_mul(const SparseMatrix &a, const Matrix &b)
Definition: dSparse.cc:7535
SparseMatrix operator*(const SparseMatrix &a, const SparseMatrix &b)
Definition: dSparse.cc:7511
SparseMatrix real(const SparseComplexMatrix &a)
Definition: dSparse.cc:554
SparseMatrix min(double d, const SparseMatrix &m)
Definition: dSparse.cc:7599
SparseMatrix operator-(const DiagMatrix &, const SparseMatrix &)
Definition: dSparse.cc:7561
SparseMatrix max(double d, const SparseMatrix &m)
Definition: dSparse.cc:7749
Matrix mul_trans(const Matrix &a, const SparseMatrix &b)
Definition: dSparse.cc:7523
SparseMatrix imag(const SparseComplexMatrix &a)
Definition: dSparse.cc:575
bool operator!=(const dim_vector &a, const dim_vector &b)
Definition: dim-vector.h:536
bool operator==(const dim_vector &a, const dim_vector &b)
Definition: dim-vector.h:520
octave::idx_vector idx_vector
Definition: idx-vector.h:1022
F77_RET_T const F77_DBLE const F77_DBLE F77_DBLE * d
bool too_large_for_float(double x)
Definition: lo-utils.cc:56
#define OCTAVE_API
Definition: main.cc:55
T octave_idx_type m
Definition: mx-inlines.cc:781
T * r
Definition: mx-inlines.cc:781
const octave_base_value const Array< octave_idx_type > & ra_idx
octave_value operator!(const octave_value &a)
Definition: ov.h:1633
template int8_t abs(int8_t)