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