Chris@49
|
1 // Copyright (C) 2008-2013 NICTA (www.nicta.com.au)
|
Chris@49
|
2 // Copyright (C) 2008-2013 Conrad Sanderson
|
Chris@49
|
3 //
|
Chris@49
|
4 // This Source Code Form is subject to the terms of the Mozilla Public
|
Chris@49
|
5 // License, v. 2.0. If a copy of the MPL was not distributed with this
|
Chris@49
|
6 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
|
Chris@49
|
7
|
Chris@49
|
8
|
Chris@49
|
9 //! \addtogroup fn_det
|
Chris@49
|
10 //! @{
|
Chris@49
|
11
|
Chris@49
|
12
|
Chris@49
|
13
|
Chris@49
|
14 //! determinant of mat
|
Chris@49
|
15 template<typename T1>
|
Chris@49
|
16 inline
|
Chris@49
|
17 arma_warn_unused
|
Chris@49
|
18 typename T1::elem_type
|
Chris@49
|
19 det
|
Chris@49
|
20 (
|
Chris@49
|
21 const Base<typename T1::elem_type,T1>& X,
|
Chris@49
|
22 const bool slow = false,
|
Chris@49
|
23 const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0
|
Chris@49
|
24 )
|
Chris@49
|
25 {
|
Chris@49
|
26 arma_extra_debug_sigprint();
|
Chris@49
|
27 arma_ignore(junk);
|
Chris@49
|
28
|
Chris@49
|
29 return auxlib::det(X, slow);
|
Chris@49
|
30 }
|
Chris@49
|
31
|
Chris@49
|
32
|
Chris@49
|
33
|
Chris@49
|
34 //! determinant of diagmat
|
Chris@49
|
35 template<typename T1>
|
Chris@49
|
36 inline
|
Chris@49
|
37 arma_warn_unused
|
Chris@49
|
38 typename T1::elem_type
|
Chris@49
|
39 det
|
Chris@49
|
40 (
|
Chris@49
|
41 const Op<T1, op_diagmat>& X,
|
Chris@49
|
42 const bool slow = false
|
Chris@49
|
43 )
|
Chris@49
|
44 {
|
Chris@49
|
45 arma_extra_debug_sigprint();
|
Chris@49
|
46 arma_ignore(slow);
|
Chris@49
|
47
|
Chris@49
|
48 typedef typename T1::elem_type eT;
|
Chris@49
|
49
|
Chris@49
|
50 const diagmat_proxy<T1> A(X.m);
|
Chris@49
|
51
|
Chris@49
|
52 const uword N = A.n_elem;
|
Chris@49
|
53
|
Chris@49
|
54 eT val1 = eT(1);
|
Chris@49
|
55 eT val2 = eT(1);
|
Chris@49
|
56
|
Chris@49
|
57 uword i,j;
|
Chris@49
|
58 for(i=0, j=1; j<N; i+=2, j+=2)
|
Chris@49
|
59 {
|
Chris@49
|
60 val1 *= A[i];
|
Chris@49
|
61 val2 *= A[j];
|
Chris@49
|
62 }
|
Chris@49
|
63
|
Chris@49
|
64
|
Chris@49
|
65 if(i < N)
|
Chris@49
|
66 {
|
Chris@49
|
67 val1 *= A[i];
|
Chris@49
|
68 }
|
Chris@49
|
69
|
Chris@49
|
70 return val1 * val2;
|
Chris@49
|
71 }
|
Chris@49
|
72
|
Chris@49
|
73
|
Chris@49
|
74
|
Chris@49
|
75 //! determinant of a triangular matrix
|
Chris@49
|
76 template<typename T1>
|
Chris@49
|
77 inline
|
Chris@49
|
78 arma_warn_unused
|
Chris@49
|
79 typename T1::elem_type
|
Chris@49
|
80 det
|
Chris@49
|
81 (
|
Chris@49
|
82 const Op<T1, op_trimat>& X,
|
Chris@49
|
83 const bool slow = false
|
Chris@49
|
84 )
|
Chris@49
|
85 {
|
Chris@49
|
86 arma_extra_debug_sigprint();
|
Chris@49
|
87 arma_ignore(slow);
|
Chris@49
|
88
|
Chris@49
|
89 typedef typename T1::elem_type eT;
|
Chris@49
|
90
|
Chris@49
|
91 const Proxy<T1> P(X.m);
|
Chris@49
|
92
|
Chris@49
|
93 const uword N = P.get_n_rows();
|
Chris@49
|
94
|
Chris@49
|
95 arma_debug_check( (N != P.get_n_cols()), "det(): matrix is not square" );
|
Chris@49
|
96
|
Chris@49
|
97 eT val1 = eT(1);
|
Chris@49
|
98 eT val2 = eT(1);
|
Chris@49
|
99
|
Chris@49
|
100 uword i,j;
|
Chris@49
|
101 for(i=0, j=1; j<N; i+=2, j+=2)
|
Chris@49
|
102 {
|
Chris@49
|
103 val1 *= P.at(i,i);
|
Chris@49
|
104 val2 *= P.at(j,j);
|
Chris@49
|
105 }
|
Chris@49
|
106
|
Chris@49
|
107 if(i < N)
|
Chris@49
|
108 {
|
Chris@49
|
109 val1 *= P.at(i,i);
|
Chris@49
|
110 }
|
Chris@49
|
111
|
Chris@49
|
112 return val1 * val2;
|
Chris@49
|
113 }
|
Chris@49
|
114
|
Chris@49
|
115
|
Chris@49
|
116
|
Chris@49
|
117 //! determinant of inv(A), without doing the inverse operation
|
Chris@49
|
118 template<typename T1>
|
Chris@49
|
119 inline
|
Chris@49
|
120 arma_warn_unused
|
Chris@49
|
121 typename T1::elem_type
|
Chris@49
|
122 det
|
Chris@49
|
123 (
|
Chris@49
|
124 const Op<T1,op_inv>& in,
|
Chris@49
|
125 const bool slow = false,
|
Chris@49
|
126 const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0
|
Chris@49
|
127 )
|
Chris@49
|
128 {
|
Chris@49
|
129 arma_extra_debug_sigprint();
|
Chris@49
|
130 arma_ignore(junk);
|
Chris@49
|
131
|
Chris@49
|
132 typedef typename T1::elem_type eT;
|
Chris@49
|
133
|
Chris@49
|
134 eT tmp = det(in.m, slow);
|
Chris@49
|
135 arma_warn( (tmp == eT(0)), "det(): warning: denominator is zero" );
|
Chris@49
|
136
|
Chris@49
|
137 return eT(1) / tmp;
|
Chris@49
|
138 }
|
Chris@49
|
139
|
Chris@49
|
140
|
Chris@49
|
141
|
Chris@49
|
142 //! determinant of trans(A)
|
Chris@49
|
143 template<typename T1>
|
Chris@49
|
144 inline
|
Chris@49
|
145 arma_warn_unused
|
Chris@49
|
146 typename T1::elem_type
|
Chris@49
|
147 det
|
Chris@49
|
148 (
|
Chris@49
|
149 const Op<T1,op_htrans>& in,
|
Chris@49
|
150 const bool slow = false,
|
Chris@49
|
151 const typename arma_blas_type_only<typename T1::elem_type>::result* junk1 = 0,
|
Chris@49
|
152 const typename arma_not_cx<typename T1::elem_type>::result* junk2 = 0
|
Chris@49
|
153 )
|
Chris@49
|
154 {
|
Chris@49
|
155 arma_extra_debug_sigprint();
|
Chris@49
|
156 arma_ignore(junk1);
|
Chris@49
|
157 arma_ignore(junk2);
|
Chris@49
|
158
|
Chris@49
|
159 return auxlib::det(in.m, slow); // bypass op_htrans
|
Chris@49
|
160 }
|
Chris@49
|
161
|
Chris@49
|
162
|
Chris@49
|
163
|
Chris@49
|
164 template<typename T>
|
Chris@49
|
165 arma_inline
|
Chris@49
|
166 arma_warn_unused
|
Chris@49
|
167 const typename arma_scalar_only<T>::result &
|
Chris@49
|
168 det(const T& x)
|
Chris@49
|
169 {
|
Chris@49
|
170 return x;
|
Chris@49
|
171 }
|
Chris@49
|
172
|
Chris@49
|
173
|
Chris@49
|
174
|
Chris@49
|
175 //! @}
|