annotate armadillo-3.900.4/include/armadillo_bits/op_inv_meat.hpp @ 84:55a047986812 tip

Update library URI so as not to be document-local
author Chris Cannam
date Wed, 22 Apr 2020 14:21:57 +0100
parents 1ec0e2823891
children
rev   line source
Chris@49 1 // Copyright (C) 2008-2011 NICTA (www.nicta.com.au)
Chris@49 2 // Copyright (C) 2008-2011 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 op_inv
Chris@49 10 //! @{
Chris@49 11
Chris@49 12
Chris@49 13 //! immediate inverse of a matrix, storing the result in a dense matrix
Chris@49 14 template<typename eT>
Chris@49 15 inline
Chris@49 16 void
Chris@49 17 op_inv::apply(Mat<eT>& out, const Mat<eT>& A, const bool slow)
Chris@49 18 {
Chris@49 19 arma_extra_debug_sigprint();
Chris@49 20
Chris@49 21 // no need to check for aliasing, due to:
Chris@49 22 // - auxlib::inv() copies A to out before inversion
Chris@49 23 // - for 2x2 and 3x3 matrices the code is alias safe
Chris@49 24
Chris@49 25 bool status = auxlib::inv(out, A, slow);
Chris@49 26
Chris@49 27 if(status == false)
Chris@49 28 {
Chris@49 29 out.reset();
Chris@49 30 arma_bad("inv(): matrix appears to be singular");
Chris@49 31 }
Chris@49 32 }
Chris@49 33
Chris@49 34
Chris@49 35
Chris@49 36 //! immediate inverse of T1, storing the result in a dense matrix
Chris@49 37 template<typename T1>
Chris@49 38 inline
Chris@49 39 void
Chris@49 40 op_inv::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_inv>& X)
Chris@49 41 {
Chris@49 42 arma_extra_debug_sigprint();
Chris@49 43
Chris@49 44 const strip_diagmat<T1> strip(X.m);
Chris@49 45
Chris@49 46 if(strip.do_diagmat == true)
Chris@49 47 {
Chris@49 48 op_inv::apply_diag(out, strip.M);
Chris@49 49 }
Chris@49 50 else
Chris@49 51 {
Chris@49 52 const uword mode = X.aux_uword_a;
Chris@49 53
Chris@49 54 const bool status = (mode == 0) ? auxlib::inv(out, X.m) : auxlib::inv(out, X.m, true);
Chris@49 55
Chris@49 56 if(status == false)
Chris@49 57 {
Chris@49 58 out.reset();
Chris@49 59 arma_bad("inv(): matrix appears to be singular");
Chris@49 60 }
Chris@49 61 }
Chris@49 62 }
Chris@49 63
Chris@49 64
Chris@49 65
Chris@49 66 template<typename T1>
Chris@49 67 inline
Chris@49 68 void
Chris@49 69 op_inv::apply_diag(Mat<typename T1::elem_type>& out, const Base<typename T1::elem_type, T1>& X)
Chris@49 70 {
Chris@49 71 arma_extra_debug_sigprint();
Chris@49 72
Chris@49 73 typedef typename T1::elem_type eT;
Chris@49 74
Chris@49 75 const diagmat_proxy_check<T1> A(X.get_ref(), out);
Chris@49 76
Chris@49 77 const uword N = A.n_elem;
Chris@49 78
Chris@49 79 out.set_size(N,N);
Chris@49 80
Chris@49 81 for(uword col=0; col<N; ++col)
Chris@49 82 {
Chris@49 83 for(uword row=0; row<col; ++row) { out.at(row,col) = eT(0); }
Chris@49 84
Chris@49 85 out.at(col,col) = eT(1) / A[col];
Chris@49 86
Chris@49 87 for(uword row=col+1; row<N; ++row) { out.at(row,col) = eT(0); }
Chris@49 88 }
Chris@49 89
Chris@49 90 }
Chris@49 91
Chris@49 92
Chris@49 93
Chris@49 94 //! inverse of T1 (triangular matrices)
Chris@49 95 template<typename T1>
Chris@49 96 inline
Chris@49 97 void
Chris@49 98 op_inv_tr::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_inv_tr>& X)
Chris@49 99 {
Chris@49 100 arma_extra_debug_sigprint();
Chris@49 101
Chris@49 102 const bool status = auxlib::inv_tr(out, X.m, X.aux_uword_a);
Chris@49 103
Chris@49 104 if(status == false)
Chris@49 105 {
Chris@49 106 out.reset();
Chris@49 107 arma_bad("inv(): matrix appears to be singular");
Chris@49 108 }
Chris@49 109 }
Chris@49 110
Chris@49 111
Chris@49 112
Chris@49 113 //! inverse of T1 (symmetric positive definite matrices)
Chris@49 114 template<typename T1>
Chris@49 115 inline
Chris@49 116 void
Chris@49 117 op_inv_sympd::apply(Mat<typename T1::elem_type>& out, const Op<T1,op_inv_sympd>& X)
Chris@49 118 {
Chris@49 119 arma_extra_debug_sigprint();
Chris@49 120
Chris@49 121 const bool status = auxlib::inv_sympd(out, X.m, X.aux_uword_a);
Chris@49 122
Chris@49 123 if(status == false)
Chris@49 124 {
Chris@49 125 out.reset();
Chris@49 126 arma_bad("inv(): matrix appears to be singular");
Chris@49 127 }
Chris@49 128 }
Chris@49 129
Chris@49 130
Chris@49 131
Chris@49 132 //! @}