GNU Octave  4.0.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
dSparse.h
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 2004-2015 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 #if !defined (octave_dSparse_h)
25 #define octave_dSparse_h 1
26 
27 #include "dMatrix.h"
28 #include "dNDArray.h"
29 #include "CMatrix.h"
30 #include "dColVector.h"
31 #include "CColVector.h"
32 
33 #include "DET.h"
34 #include "MSparse.h"
35 #include "MSparse-defs.h"
36 
37 #include "Sparse-op-decls.h"
38 
39 #include "MatrixType.h"
40 
41 class PermMatrix;
42 class DiagMatrix;
44 class SparseBoolMatrix;
45 
46 class
47 OCTAVE_API
49 {
50 public:
51 
52  typedef void (*solve_singularity_handler) (double rcond);
53 
54  SparseMatrix (void) : MSparse<double> () { }
55 
57  : MSparse<double> (r, c) { }
58 
59  SparseMatrix (const dim_vector& dv, octave_idx_type nz = 0) :
60  MSparse<double> (dv, nz) { }
61 
62  explicit SparseMatrix (octave_idx_type r, octave_idx_type c, double val)
63  : MSparse<double> (r, c, val) { }
64 
65  SparseMatrix (const SparseMatrix& a) : MSparse<double> (a) { }
66 
67  SparseMatrix (const SparseMatrix& a, const dim_vector& dv)
68  : MSparse<double> (a, dv) { }
69 
71 
72  SparseMatrix (const Sparse<double>& a) : MSparse<double> (a) { }
73 
74  explicit SparseMatrix (const SparseBoolMatrix& a);
75 
76  explicit SparseMatrix (const Matrix& a) : MSparse<double> (a) { }
77 
78  explicit SparseMatrix (const NDArray& a) : MSparse<double> (a) { }
79 
80  SparseMatrix (const Array<double>& a, const idx_vector& r,
81  const idx_vector& c, octave_idx_type nr = -1,
82  octave_idx_type nc = -1, bool sum_terms = true,
83  octave_idx_type nzm = -1)
84  : MSparse<double> (a, r, c, nr, nc, sum_terms, nzm) { }
85 
86  explicit SparseMatrix (const DiagMatrix& a);
87 
88  explicit SparseMatrix (const PermMatrix& a) : MSparse<double>(a) { }
89 
91  octave_idx_type num_nz) : MSparse<double> (r, c, num_nz) { }
92 
94  {
96  return *this;
97  }
98 
99  bool operator == (const SparseMatrix& a) const;
100  bool operator != (const SparseMatrix& a) const;
101 
102  bool is_symmetric (void) const;
103 
104  SparseMatrix max (int dim = -1) const;
105  SparseMatrix max (Array<octave_idx_type>& index, int dim = -1) const;
106  SparseMatrix min (int dim = -1) const;
107  SparseMatrix min (Array<octave_idx_type>& index, int dim = -1) const;
108 
109  // destructive insert/delete/reorder operations
110 
112  octave_idx_type c);
113 
114  SparseMatrix& insert (const SparseMatrix& a,
115  const Array<octave_idx_type>& indx);
116 
117  SparseMatrix concat (const SparseMatrix& rb,
121 
122  friend OCTAVE_API SparseMatrix real (const SparseComplexMatrix& a);
123  friend OCTAVE_API SparseMatrix imag (const SparseComplexMatrix& a);
124 
125  friend OCTAVE_API SparseMatrix atan2 (const double& x, const SparseMatrix& y)
126  GCC_ATTR_DEPRECATED ;
127  friend OCTAVE_API SparseMatrix atan2 (const SparseMatrix& x, const double& y)
128  GCC_ATTR_DEPRECATED ;
129  friend OCTAVE_API SparseMatrix atan2 (const SparseMatrix& x,
130  const SparseMatrix& y)
131  GCC_ATTR_DEPRECATED ;
132 
133  SparseMatrix transpose (void) const
134  {
135  return MSparse<double>::transpose ();
136  }
137  SparseMatrix hermitian (void) const { return transpose (); }
138 
139  // extract row or column i.
140 
141  RowVector row (octave_idx_type i) const;
142 
143  ColumnVector column (octave_idx_type i) const;
144 
145 private:
146  SparseMatrix dinverse (MatrixType &mattyp, octave_idx_type& info,
147  double& rcond, const bool force = false,
148  const bool calccond = true) const;
149 
150  SparseMatrix tinverse (MatrixType &mattyp, octave_idx_type& info,
151  double& rcond, const bool force = false,
152  const bool calccond = true) const;
153 
154 public:
155  SparseMatrix inverse (void) const;
156  SparseMatrix inverse (MatrixType& mattype) const;
157  SparseMatrix inverse (MatrixType& mattype, octave_idx_type& info) const;
158  SparseMatrix inverse (MatrixType& mattype, octave_idx_type& info,
159  double& rcond, int force = 0, int calc_cond = 1) const;
160 
161  DET determinant (void) const;
162  DET determinant (octave_idx_type& info) const;
163  DET determinant (octave_idx_type& info, double& rcond,
164  int calc_cond = 1) const;
165 
166 private:
167  // Diagonal matrix solvers
168  Matrix dsolve (MatrixType &typ, const Matrix& b, octave_idx_type& info,
169  double& rcond, solve_singularity_handler sing_handler,
170  bool calc_cond = false) const;
171 
172  ComplexMatrix dsolve (MatrixType &typ, const ComplexMatrix& b,
173  octave_idx_type& info, double& rcond,
174  solve_singularity_handler sing_handler,
175  bool calc_cond = false) const;
176 
177  SparseMatrix dsolve (MatrixType &typ, const SparseMatrix& b,
178  octave_idx_type& info, double& rcond,
179  solve_singularity_handler sing_handler,
180  bool calc_cond = false) const;
181 
182  SparseComplexMatrix dsolve (MatrixType &typ, const SparseComplexMatrix& b,
183  octave_idx_type& info, double& rcond,
184  solve_singularity_handler sing_handler,
185  bool calc_cond = false) const;
186 
187  // Upper triangular matrix solvers
188  Matrix utsolve (MatrixType &typ, const Matrix& b, octave_idx_type& info,
189  double& rcond, solve_singularity_handler sing_handler,
190  bool calc_cond = false) const;
191 
193  octave_idx_type& info, double& rcond,
194  solve_singularity_handler sing_handler,
195  bool calc_cond = false) const;
196 
198  octave_idx_type& info, double& rcond,
199  solve_singularity_handler sing_handler,
200  bool calc_cond = false) const;
201 
203  octave_idx_type& info, double& rcond,
204  solve_singularity_handler sing_handler,
205  bool calc_cond = false) const;
206 
207  // Lower triangular matrix solvers
208  Matrix ltsolve (MatrixType &typ, const Matrix& b, octave_idx_type& info,
209  double& rcond, solve_singularity_handler sing_handler,
210  bool calc_cond = false) const;
211 
213  octave_idx_type& info, double& rcond,
214  solve_singularity_handler sing_handler,
215  bool calc_cond = false) const;
216 
218  octave_idx_type& info, double& rcond,
219  solve_singularity_handler sing_handler,
220  bool calc_cond = false) const;
221 
223  octave_idx_type& info, double& rcond,
224  solve_singularity_handler sing_handler,
225  bool calc_cond = false) const;
226 
227  // Tridiagonal matrix solvers
228  Matrix trisolve (MatrixType &typ, const Matrix& b, octave_idx_type& info,
229  double& rcond, solve_singularity_handler sing_handler,
230  bool calc_cond = false) const;
231 
232  ComplexMatrix trisolve (MatrixType &typ, const ComplexMatrix& b,
233  octave_idx_type& info, double& rcond,
234  solve_singularity_handler sing_handler,
235  bool calc_cond = false) const;
236 
237  SparseMatrix trisolve (MatrixType &typ, const SparseMatrix& b,
238  octave_idx_type& info, double& rcond,
239  solve_singularity_handler sing_handler,
240  bool calc_cond = false) const;
241 
242  SparseComplexMatrix trisolve (MatrixType &typ, const SparseComplexMatrix& b,
243  octave_idx_type& info, double& rcond,
244  solve_singularity_handler sing_handler,
245  bool calc_cond = false) const;
246 
247  // Banded matrix solvers (umfpack/cholesky)
248  Matrix 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 
252  ComplexMatrix bsolve (MatrixType &typ, const ComplexMatrix& b,
253  octave_idx_type& info, double& rcond,
254  solve_singularity_handler sing_handler,
255  bool calc_cond = false) const;
256 
257  SparseMatrix bsolve (MatrixType &typ, const SparseMatrix& b,
258  octave_idx_type& info, double& rcond,
259  solve_singularity_handler sing_handler,
260  bool calc_cond = false) const;
261 
262  SparseComplexMatrix bsolve (MatrixType &typ, const SparseComplexMatrix& b,
263  octave_idx_type& info, double& rcond,
264  solve_singularity_handler sing_handler,
265  bool calc_cond = false) const;
266 
267  // Full matrix solvers (umfpack/cholesky)
268  void * factorize (octave_idx_type& err, double &rcond, Matrix &Control,
269  Matrix &Info, solve_singularity_handler sing_handler,
270  bool calc_cond = false) const;
271 
272  Matrix fsolve (MatrixType &typ, const Matrix& b, octave_idx_type& info,
273  double& rcond, solve_singularity_handler sing_handler,
274  bool calc_cond = false) const;
275 
276  ComplexMatrix fsolve (MatrixType &typ, const ComplexMatrix& b,
277  octave_idx_type& info, double& rcond,
278  solve_singularity_handler sing_handler,
279  bool calc_cond = false) const;
280 
281  SparseMatrix fsolve (MatrixType &typ, const SparseMatrix& b,
282  octave_idx_type& info, double& rcond,
283  solve_singularity_handler sing_handler,
284  bool calc_cond = false) const;
285 
286  SparseComplexMatrix fsolve (MatrixType &typ, const SparseComplexMatrix& b,
287  octave_idx_type& info, double& rcond,
288  solve_singularity_handler sing_handler,
289  bool calc_cond = false) const;
290 
291 public:
292  // Generic interface to solver with no probing of type
293  Matrix solve (MatrixType &typ, const Matrix& b) const;
294  Matrix solve (MatrixType &typ, const Matrix& b, octave_idx_type& info) const;
295  Matrix solve (MatrixType &typ, const Matrix& b, octave_idx_type& info,
296  double& rcond) const;
297  Matrix solve (MatrixType &typ, const Matrix& b, octave_idx_type& info,
298  double& rcond, solve_singularity_handler sing_handler,
299  bool singular_fallback = true) const;
300 
301  ComplexMatrix solve (MatrixType &typ, const ComplexMatrix& b) const;
302  ComplexMatrix solve (MatrixType &typ, const ComplexMatrix& b,
303  octave_idx_type& info) const;
304  ComplexMatrix solve (MatrixType &typ, const ComplexMatrix& b,
305  octave_idx_type& info, double& rcond) const;
306  ComplexMatrix solve (MatrixType &typ, const ComplexMatrix& b,
307  octave_idx_type& info, double& rcond,
308  solve_singularity_handler sing_handler,
309  bool singular_fallback = true) const;
310 
311  SparseMatrix solve (MatrixType &typ, const SparseMatrix& b) const;
312  SparseMatrix solve (MatrixType &typ, const SparseMatrix& b,
313  octave_idx_type& info) const;
314  SparseMatrix solve (MatrixType &typ, const SparseMatrix& b,
315  octave_idx_type& info, double& rcond) const;
316  SparseMatrix solve (MatrixType &typ, const SparseMatrix& b,
317  octave_idx_type& info, double& rcond,
318  solve_singularity_handler sing_handler,
319  bool singular_fallback = true) const;
320 
321  SparseComplexMatrix solve (MatrixType &typ,
322  const SparseComplexMatrix& b) const;
323  SparseComplexMatrix solve (MatrixType &typ, const SparseComplexMatrix& b,
324  octave_idx_type& info) const;
325  SparseComplexMatrix solve (MatrixType &typ, const SparseComplexMatrix& b,
326  octave_idx_type& info, double& rcond) const;
327  SparseComplexMatrix solve (MatrixType &typ, const SparseComplexMatrix& b,
328  octave_idx_type& info, double& rcond,
329  solve_singularity_handler sing_handler,
330  bool singular_fallabck = true) const;
331 
332  ColumnVector solve (MatrixType &typ, const ColumnVector& b) const;
333  ColumnVector solve (MatrixType &typ, const ColumnVector& b,
334  octave_idx_type& info) const;
335  ColumnVector solve (MatrixType &typ, const ColumnVector& b,
336  octave_idx_type& info, double& rcond) const;
337  ColumnVector solve (MatrixType &typ, const ColumnVector& b,
338  octave_idx_type& info, double& rcond,
339  solve_singularity_handler sing_handler) const;
340 
341  ComplexColumnVector solve (MatrixType &typ,
342  const ComplexColumnVector& b) const;
343  ComplexColumnVector solve (MatrixType &typ, const ComplexColumnVector& b,
344  octave_idx_type& info) const;
345  ComplexColumnVector solve (MatrixType &typ, const ComplexColumnVector& b,
346  octave_idx_type& info, double& rcond) const;
347  ComplexColumnVector solve (MatrixType &typ, const ComplexColumnVector& b,
348  octave_idx_type& info, double& rcond,
349  solve_singularity_handler sing_handler) const;
350 
351  // Generic interface to solver with probing of type
352  Matrix solve (const Matrix& b) const;
353  Matrix solve (const Matrix& b, octave_idx_type& info) const;
354  Matrix solve (const Matrix& b, octave_idx_type& info, double& rcond) const;
355  Matrix solve (const Matrix& b, octave_idx_type& info, double& rcond,
356  solve_singularity_handler sing_handler) const;
357 
358  ComplexMatrix solve (const ComplexMatrix& b) const;
359  ComplexMatrix solve (const ComplexMatrix& b, octave_idx_type& info) const;
360  ComplexMatrix solve (const ComplexMatrix& b, octave_idx_type& info,
361  double& rcond) const;
362  ComplexMatrix solve (const ComplexMatrix& b, octave_idx_type& info,
363  double& rcond,
364  solve_singularity_handler sing_handler) const;
365 
366  SparseMatrix solve (const SparseMatrix& b) const;
367  SparseMatrix solve (const SparseMatrix& b, octave_idx_type& info) const;
368  SparseMatrix solve (const SparseMatrix& b, octave_idx_type& info,
369  double& rcond) const;
370  SparseMatrix solve (const SparseMatrix& b, octave_idx_type& info,
371  double& rcond,
372  solve_singularity_handler sing_handler) const;
373 
374  SparseComplexMatrix solve (const SparseComplexMatrix& b) const;
376  octave_idx_type& info) const;
378  octave_idx_type& info, double& rcond) const;
380  octave_idx_type& info, double& rcond,
381  solve_singularity_handler sing_handler) const;
382 
383  ColumnVector solve (const ColumnVector& b) const;
384  ColumnVector solve (const ColumnVector& b, octave_idx_type& info) const;
385  ColumnVector solve (const ColumnVector& b, octave_idx_type& info,
386  double& rcond) const;
387  ColumnVector solve (const ColumnVector& b, octave_idx_type& info,
388  double& rcond,
389  solve_singularity_handler sing_handler) const;
390 
391  ComplexColumnVector solve (const ComplexColumnVector& b) const;
393  octave_idx_type& info) const;
395  octave_idx_type& info, double& rcond) const;
397  octave_idx_type& info, double& rcond,
398  solve_singularity_handler sing_handler) const;
399 
400  // other operations
401 
402  bool any_element_is_negative (bool = false) const;
403  bool any_element_is_nan (void) const;
404  bool any_element_is_inf_or_nan (void) const;
405  bool any_element_not_one_or_zero (void) const;
406  bool all_elements_are_zero (void) const;
407  bool all_elements_are_int_or_inf_or_nan (void) const;
408  bool all_integers (double& max_val, double& min_val) const;
409  bool too_large_for_float (void) const;
410 
411  SparseBoolMatrix operator ! (void) const;
412 
413  SparseBoolMatrix all (int dim = -1) const;
414  SparseBoolMatrix any (int dim = -1) const;
415 
416  SparseMatrix cumprod (int dim = -1) const;
417  SparseMatrix cumsum (int dim = -1) const;
418  SparseMatrix prod (int dim = -1) const;
419  SparseMatrix sum (int dim = -1) const;
420  SparseMatrix sumsq (int dim = -1) const;
421  SparseMatrix abs (void) const;
422 
423  SparseMatrix diag (octave_idx_type k = 0) const;
424 
425  Matrix matrix_value (void) const;
426 
427  SparseMatrix squeeze (void) const;
428 
429  SparseMatrix reshape (const dim_vector& new_dims) const;
430 
432  bool inv = false) const;
433 
434  SparseMatrix ipermute (const Array<octave_idx_type>& vec) const;
435 
436  // i/o
437 
438  friend OCTAVE_API std::ostream& operator << (std::ostream& os,
439  const SparseMatrix& a);
440  friend OCTAVE_API std::istream& operator >> (std::istream& is,
441  SparseMatrix& a);
442 
443 };
444 
445 // Publish externally used friend functions.
446 
447 extern OCTAVE_API SparseMatrix real (const SparseComplexMatrix& a);
448 extern OCTAVE_API SparseMatrix imag (const SparseComplexMatrix& a);
449 
450 // Other operators.
451 
452 extern OCTAVE_API SparseMatrix operator * (const SparseMatrix& a,
453  const SparseMatrix& b);
454 extern OCTAVE_API Matrix operator * (const Matrix& a,
455  const SparseMatrix& b);
456 extern OCTAVE_API Matrix mul_trans (const Matrix& a,
457  const SparseMatrix& b);
458 extern OCTAVE_API Matrix operator * (const SparseMatrix& a,
459  const Matrix& b);
460 extern OCTAVE_API Matrix trans_mul (const SparseMatrix& a,
461  const Matrix& b);
462 
463 extern OCTAVE_API SparseMatrix operator * (const DiagMatrix&,
464  const SparseMatrix&);
465 extern OCTAVE_API SparseMatrix operator * (const SparseMatrix&,
466  const DiagMatrix&);
467 
468 extern OCTAVE_API SparseMatrix operator + (const DiagMatrix&,
469  const SparseMatrix&);
470 extern OCTAVE_API SparseMatrix operator + (const SparseMatrix&,
471  const DiagMatrix&);
472 extern OCTAVE_API SparseMatrix operator - (const DiagMatrix&,
473  const SparseMatrix&);
474 extern OCTAVE_API SparseMatrix operator - (const SparseMatrix&,
475  const DiagMatrix&);
476 
477 extern OCTAVE_API SparseMatrix operator * (const PermMatrix&,
478  const SparseMatrix&);
479 extern OCTAVE_API SparseMatrix operator * (const SparseMatrix&,
480  const PermMatrix&);
481 
482 extern OCTAVE_API SparseMatrix min (double d, const SparseMatrix& m);
483 extern OCTAVE_API SparseMatrix min (const SparseMatrix& m, double d);
484 extern OCTAVE_API SparseMatrix min (const SparseMatrix& a,
485  const SparseMatrix& b);
486 
487 extern OCTAVE_API SparseMatrix max (double d, const SparseMatrix& m);
488 extern OCTAVE_API SparseMatrix max (const SparseMatrix& m, double d);
489 extern OCTAVE_API SparseMatrix max (const SparseMatrix& a,
490  const SparseMatrix& b);
491 
492 SPARSE_SMS_CMP_OP_DECLS (SparseMatrix, double, OCTAVE_API)
493 SPARSE_SMS_BOOL_OP_DECLS (SparseMatrix, double, OCTAVE_API)
494 
495 SPARSE_SSM_CMP_OP_DECLS (double, SparseMatrix, OCTAVE_API)
496 SPARSE_SSM_BOOL_OP_DECLS (double, SparseMatrix, OCTAVE_API)
497 
498 SPARSE_SMSM_CMP_OP_DECLS (SparseMatrix, SparseMatrix, OCTAVE_API)
499 SPARSE_SMSM_BOOL_OP_DECLS (SparseMatrix, SparseMatrix, OCTAVE_API)
500 
501 SPARSE_FORWARD_DEFS (MSparse, SparseMatrix, Matrix, double)
502 
503 #ifdef USE_64_BIT_IDX_T
504 #define UMFPACK_DNAME(name) umfpack_dl_ ## name
505 #else
506 #define UMFPACK_DNAME(name) umfpack_di_ ## name
507 #endif
508 
509 #endif
OCTAVE_API SparseMatrix operator-(const DiagMatrix &, const SparseMatrix &)
Definition: dSparse.cc:7611
boolMatrix matrix_value(void) const
Definition: boolSparse.cc:248
OCTAVE_API SparseMatrix operator+(const DiagMatrix &, const SparseMatrix &)
Definition: dSparse.cc:7605
#define SPARSE_SMSM_BOOL_OP_DECLS(M1, M2, API)
bool operator!=(const dim_vector &a, const dim_vector &b)
Definition: dim-vector.h:548
SparseMatrix(const Array< double > &a, const idx_vector &r, const 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:80
const octave_base_value const Array< octave_idx_type > & ra_idx
std::istream & operator>>(std::istream &is, SparseBoolMatrix &a)
Definition: boolSparse.cc:279
SparseMatrix transpose(void) const
Definition: dSparse.h:133
static ComplexMatrix ltsolve(const SparseComplexMatrix &, const ColumnVector &, const ComplexMatrix &)
SparseMatrix(const NDArray &a)
Definition: dSparse.h:78
SparseBoolMatrix ipermute(const Array< octave_idx_type > &vec) const
Definition: boolSparse.cc:318
OCTAVE_API Matrix mul_trans(const Matrix &a, const SparseMatrix &b)
Definition: dSparse.cc:7573
OCTAVE_API SparseMatrix imag(const SparseComplexMatrix &a)
Definition: dSparse.cc:661
SparseBoolMatrix & operator=(const SparseBoolMatrix &a)
Definition: boolSparse.h:77
SparseMatrix(const Matrix &a)
Definition: dSparse.h:76
static void transpose(octave_idx_type N, const octave_idx_type *ridx, const octave_idx_type *cidx, octave_idx_type *ridx2, octave_idx_type *cidx2)
Definition: symrcm.cc:382
ComplexNDArray concat(NDArray &ra, ComplexNDArray &rb, const Array< octave_idx_type > &ra_idx)
Definition: CNDArray.cc:664
#define SPARSE_SMS_CMP_OP_DECLS(M, S, API)
SparseBoolMatrix permute(const Array< octave_idx_type > &vec, bool inv=false) const
Definition: boolSparse.cc:312
SparseMatrix(const SparseMatrix &a, const dim_vector &dv)
Definition: dSparse.h:67
SparseMatrix hermitian(void) const
Definition: dSparse.h:137
OCTAVE_API SparseMatrix max(double d, const SparseMatrix &m)
Definition: dSparse.cc:7799
SparseMatrix(octave_idx_type r, octave_idx_type c)
Definition: dSparse.h:56
#define SPARSE_SMS_BOOL_OP_DECLS(M, S, API)
SparseBoolMatrix index(const idx_vector &i, bool resize_ok) const
Definition: boolSparse.cc:293
SparseMatrix(const SparseMatrix &a)
Definition: dSparse.h:65
MSparse< T > transpose(void) const
Definition: MSparse.h:94
SparseMatrix(const MSparse< double > &a)
Definition: dSparse.h:70
F77_RET_T const double const double double * d
#define SPARSE_SSM_CMP_OP_DECLS(S, M, API)
SparseBoolMatrix reshape(const dim_vector &new_dims) const
Definition: boolSparse.cc:306
SparseMatrix(octave_idx_type r, octave_idx_type c, octave_idx_type num_nz)
Definition: dSparse.h:90
SparseMatrix sum(int dim=-1) const
Definition: boolSparse.cc:191
#define SPARSE_SMSM_CMP_OP_DECLS(M1, M2, API)
OCTAVE_API SparseMatrix min(double d, const SparseMatrix &m)
Definition: dSparse.cc:7649
SparseBoolMatrix all(int dim=-1) const
Definition: boolSparse.cc:138
#define SPARSE_SSM_BOOL_OP_DECLS(S, M, API)
SparseMatrix(const dim_vector &dv, octave_idx_type nz=0)
Definition: dSparse.h:59
SparseMatrix(const Sparse< double > &a)
Definition: dSparse.h:72
Definition: DET.h:31
OCTAVE_API SparseMatrix operator*(const SparseMatrix &a, const SparseMatrix &b)
Definition: dSparse.cc:7561
SparseBoolMatrix & insert(const SparseBoolMatrix &a, octave_idx_type r, octave_idx_type c)
Definition: boolSparse.cc:75
SparseMatrix(const PermMatrix &a)
Definition: dSparse.h:88
Definition: dMatrix.h:35
SparseBoolMatrix any(int dim=-1) const
Definition: boolSparse.cc:144
SparseBoolMatrix squeeze(void) const
Definition: boolSparse.cc:287
bool operator==(const dim_vector &a, const dim_vector &b)
Definition: dim-vector.h:519
OCTAVE_API Matrix trans_mul(const SparseMatrix &a, const Matrix &b)
Definition: dSparse.cc:7585
SparseMatrix atan2(const double &x, const SparseMatrix &y)
Definition: dSparse.cc:689
SparseMatrix(void)
Definition: dSparse.h:54
static ComplexMatrix utsolve(const SparseComplexMatrix &, const ColumnVector &, const ComplexMatrix &)
template OCTAVE_API std::ostream & operator<<(std::ostream &, const Array< bool > &)
#define SPARSE_FORWARD_DEFS(B, R, F, T)
Definition: MSparse-defs.h:191
OCTAVE_API SparseMatrix real(const SparseComplexMatrix &a)
Definition: dSparse.cc:640
MSparse< T > & operator=(const MSparse< T > &a)
Definition: MSparse.h:76
SparseBoolMatrix diag(octave_idx_type k=0) const
Definition: boolSparse.cc:242
octave_value operator!(const octave_value &a)
Definition: ov.h:1293
SparseMatrix(octave_idx_type r, octave_idx_type c, double val)
Definition: dSparse.h:62
T abs(T x)
Definition: pr-output.cc:3062
F77_RET_T const double * x