GNU Octave  4.4.1
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
det.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 1996-2018 John W. Eaton
4 
5 This file is part of Octave.
6 
7 Octave is free software: you can redistribute it and/or modify it
8 under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11 
12 Octave is distributed in the hope that it will be useful, but
13 WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16 
17 You should have received a copy of the GNU General Public License
18 along with Octave; see the file COPYING. If not, see
19 <https://www.gnu.org/licenses/>.
20 
21 */
22 
23 #if defined (HAVE_CONFIG_H)
24 # include "config.h"
25 #endif
26 
27 #include "DET.h"
28 
29 #include "defun.h"
30 #include "error.h"
31 #include "errwarn.h"
32 #include "ovl.h"
33 #include "ops.h"
34 
35 #include "ov-re-mat.h"
36 #include "ov-cx-mat.h"
37 #include "ov-flt-re-mat.h"
38 #include "ov-flt-cx-mat.h"
39 #include "ov-re-diag.h"
40 #include "ov-cx-diag.h"
41 #include "ov-flt-re-diag.h"
42 #include "ov-flt-cx-diag.h"
43 #include "ov-perm.h"
44 
45 #define MAYBE_CAST(VAR, CLASS) \
46  const CLASS *VAR = (arg.type_id () == CLASS::static_type_id () \
47  ? dynamic_cast<const CLASS *> (&arg.get_rep ()) \
48  : nullptr)
49 
50 DEFUN (det, args, nargout,
51  doc: /* -*- texinfo -*-
52 @deftypefn {} {} det (@var{A})
53 @deftypefnx {} {[@var{d}, @var{rcond}] =} det (@var{A})
54 Compute the determinant of @var{A}.
55 
56 Return an estimate of the reciprocal condition number if requested.
57 
58 Programming Notes: Routines from @sc{lapack} are used for full matrices and
59 code from @sc{umfpack} is used for sparse matrices.
60 
61 The determinant should not be used to check a matrix for singularity.
62 For that, use any of the condition number functions: @code{cond},
63 @code{condest}, @code{rcond}.
64 @seealso{cond, condest, rcond}
65 @end deftypefn */)
66 {
67  if (args.length () != 1)
68  print_usage ();
69 
70  octave_value arg = args(0);
71 
72  if (arg.isempty ())
73  return ovl (1.0);
74 
75  if (arg.rows () != arg.columns ())
76  err_square_matrix_required ("det", "A");
77 
79 
80  bool isfloat = arg.is_single_type ();
81 
82  if (arg.is_diag_matrix ())
83  {
84  if (nargout <= 1)
85  retval.resize (1);
86 
87  if (arg.iscomplex ())
88  {
89  if (isfloat)
90  {
92  .determinant ().value ();
93  if (nargout > 1)
95  }
96  else
97  {
99  .determinant ().value ();
100  if (nargout > 1)
102  }
103  }
104  else
105  {
106  if (isfloat)
107  {
109  .determinant ().value ();
110  if (nargout > 1)
112  }
113  else
114  {
116  if (nargout > 1)
117  retval(1) = arg.diag_matrix_value ().rcond ();
118  }
119  }
120  }
121  else if (arg.is_perm_matrix ())
122  {
123  if (nargout <= 1)
124  retval.resize (1);
125 
126  retval(0) = static_cast<double> (arg.perm_matrix_value ().determinant ());
127  if (nargout > 1)
128  retval(1) = 1.0;
129  }
130  else if (arg.is_single_type ())
131  {
132  if (arg.isreal ())
133  {
134  octave_idx_type info;
135  float rcond = 0.0;
136  // Always compute rcond, so we can detect singular matrices.
138 
140  MatrixType mtype = (rep ? rep -> matrix_type () : MatrixType ());
141  FloatDET det = m.determinant (mtype, info, rcond);
142  retval(0) = (info == -1 ? 0.0f : det.value ());
143  retval(1) = rcond;
144  if (rep)
145  rep->matrix_type (mtype);
146  }
147  else if (arg.iscomplex ())
148  {
149  octave_idx_type info;
150  float rcond = 0.0;
151  // Always compute rcond, so we can detect singular matrices.
153 
155  MatrixType mtype = (rep ? rep -> matrix_type () : MatrixType ());
156  FloatComplexDET det = m.determinant (mtype, info, rcond);
157  retval(0) = (info == -1 ? FloatComplex (0.0) : det.value ());
158  retval(1) = rcond;
159  if (rep)
160  rep->matrix_type (mtype);
161  }
162  }
163  else
164  {
165  if (arg.isreal ())
166  {
167  octave_idx_type info;
168  double rcond = 0.0;
169  // Always compute rcond, so we can detect singular matrices.
170  if (arg.issparse ())
171  {
173 
174  DET det = m.determinant (info, rcond);
175  retval(0) = (info == -1 ? 0.0 : det.value ());
176  retval(1) = rcond;
177  }
178  else
179  {
180  Matrix m = arg.matrix_value ();
181 
182  MAYBE_CAST (rep, octave_matrix);
183  MatrixType mtype = (rep ? rep -> matrix_type ()
184  : MatrixType ());
185  DET det = m.determinant (mtype, info, rcond);
186  retval(0) = (info == -1 ? 0.0 : det.value ());
187  retval(1) = rcond;
188  if (rep)
189  rep->matrix_type (mtype);
190  }
191  }
192  else if (arg.iscomplex ())
193  {
194  octave_idx_type info;
195  double rcond = 0.0;
196  // Always compute rcond, so we can detect singular matrices.
197  if (arg.issparse ())
198  {
200 
201  ComplexDET det = m.determinant (info, rcond);
202  retval(0) = (info == -1 ? Complex (0.0) : det.value ());
203  retval(1) = rcond;
204  }
205  else
206  {
208 
210  MatrixType mtype = (rep ? rep -> matrix_type ()
211  : MatrixType ());
212  ComplexDET det = m.determinant (mtype, info, rcond);
213  retval(0) = (info == -1 ? Complex (0.0) : det.value ());
214  retval(1) = rcond;
215  if (rep)
216  rep->matrix_type (mtype);
217  }
218  }
219  else
220  err_wrong_type_arg ("det", arg);
221  }
222 
223  return retval;
224 }
225 
226 /*
227 %!assert (det ([1, 2; 3, 4]), -2, 10*eps)
228 %!assert (det (single ([1, 2; 3, 4])), single (-2), 10*eps ("single"))
229 %!assert (det (eye (2000)), 1)
230 %!error det ()
231 %!error det (1, 2)
232 %!error <must be a square matrix> det ([1, 2; 3, 4; 5, 6])
233 */
T value() const
Definition: DET.h:69
FloatDET determinant(void) const
Definition: fDiagMatrix.cc:307
ComplexDET determinant(void) const
Definition: CMatrix.cc:1340
bool isempty(void) const
Definition: ov.h:529
DET determinant(void) const
Definition: dMatrix.cc:1035
OCTINTERP_API void print_usage(void)
Definition: defun.cc:54
#define DEFUN(name, args_name, nargout_name, doc)
Macro to define a builtin function.
Definition: defun.h:53
DET determinant(void) const
Definition: dSparse.cc:992
DiagMatrix diag_matrix_value(bool force=false) const
Definition: ov.h:891
FloatDiagMatrix float_diag_matrix_value(bool force=false) const
Definition: ov.h:894
FloatMatrix float_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:837
void err_square_matrix_required(const char *fcn, const char *name)
Definition: errwarn.cc:118
bool is_perm_matrix(void) const
Definition: ov.h:574
octave_value arg
Definition: pr-output.cc:3244
DET determinant(void) const
Definition: dDiagMatrix.cc:307
octave_value resize(const dim_vector &dv, bool fill=false) const
Definition: ov.h:511
FloatComplexMatrix float_complex_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:856
octave_idx_type determinant(void) const
Definition: PermMatrix.cc:117
OCTAVE_EXPORT octave_value_list return the number of command line arguments passed to Octave If called with the optional argument the function xample nargout(@histc)
Definition: ov-usr-fcn.cc:997
FloatComplexDET determinant(void) const
Definition: fCMatrix.cc:1341
OCTAVE_EXPORT octave_value_list isfloat
Definition: data.cc:3212
double rcond(void) const
Definition: CDiagMatrix.cc:485
octave_idx_type columns(void) const
Definition: ov.h:474
ComplexDiagMatrix complex_diag_matrix_value(bool force=false) const
Definition: ov.h:897
Definition: DET.h:34
bool is_single_type(void) const
Definition: ov.h:651
double rcond(void) const
Definition: dDiagMatrix.cc:321
octave_idx_type rows(void) const
Definition: ov.h:472
bool issparse(void) const
Definition: ov.h:730
octave_value retval
Definition: data.cc:6246
#define MAYBE_CAST(VAR, CLASS)
Definition: det.cc:45
float rcond(void) const
Definition: dMatrix.h:36
void err_wrong_type_arg(const char *name, const char *s)
Definition: errwarn.cc:162
SparseComplexMatrix sparse_complex_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:885
float rcond(void) const
Definition: fDiagMatrix.cc:321
bool isreal(void) const
Definition: ov.h:703
FloatComplexDiagMatrix float_complex_diag_matrix_value(bool force=false) const
Definition: ov.h:901
FloatComplexDET determinant(void) const
OCTAVE_EXPORT octave_value_list isa nd deftypefn *return ovl(args(0).isinteger())
ComplexDET determinant(void) const
Definition: CSparse.cc:1052
SparseMatrix sparse_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:881
bool iscomplex(void) const
Definition: ov.h:710
bool is_diag_matrix(void) const
Definition: ov.h:571
FloatDET determinant(void) const
Definition: fMatrix.cc:1039
PermMatrix perm_matrix_value(void) const
Definition: ov.h:904
std::complex< float > FloatComplex
Definition: oct-cmplx.h:32
ComplexMatrix complex_matrix_value(bool frc_str_conv=false) const
Definition: ov.h:852
std::complex< double > Complex
Definition: oct-cmplx.h:31
ComplexDET determinant(void) const
Definition: CDiagMatrix.cc:471
Matrix matrix_value(bool frc_str_conv=false) const
Definition: ov.h:834