$search
00001 // Copyright (C) 2011 NICTA (www.nicta.com.au) 00002 // Copyright (C) 2011 Conrad Sanderson 00003 // 00004 // This file is part of the Armadillo C++ library. 00005 // It is provided without any warranty of fitness 00006 // for any purpose. You can redistribute this file 00007 // and/or modify it under the terms of the GNU 00008 // Lesser General Public License (LGPL) as published 00009 // by the Free Software Foundation, either version 3 00010 // of the License or (at your option) any later version. 00011 // (see http://www.opensource.org/licenses for more info) 00012 00013 00016 00017 00018 00019 template<typename T1> 00020 inline 00021 void 00022 op_symmat::apply 00023 ( 00024 Mat<typename T1::elem_type>& out, 00025 const Op<T1,op_symmat>& in, 00026 const typename arma_not_cx<typename T1::elem_type>::result* junk 00027 ) 00028 { 00029 arma_extra_debug_sigprint(); 00030 arma_ignore(junk); 00031 00032 typedef typename T1::elem_type eT; 00033 00034 const unwrap<T1> tmp(in.m); 00035 const Mat<eT>& A = tmp.M; 00036 00037 arma_debug_check( (A.is_square() == false), "symmatu()/symmatl(): given matrix must be square" ); 00038 00039 const uword N = A.n_rows; 00040 const bool upper = (in.aux_uword_a == 0); 00041 00042 if(&out != &A) 00043 { 00044 out.copy_size(A); 00045 00046 if(upper) 00047 { 00048 // upper triangular: copy the diagonal and the elements above the diagonal 00049 00050 for(uword i=0; i<N; ++i) 00051 { 00052 const eT* A_data = A.colptr(i); 00053 eT* out_data = out.colptr(i); 00054 00055 arrayops::copy( out_data, A_data, i+1 ); 00056 } 00057 } 00058 else 00059 { 00060 // lower triangular: copy the diagonal and the elements below the diagonal 00061 00062 for(uword i=0; i<N; ++i) 00063 { 00064 const eT* A_data = A.colptr(i); 00065 eT* out_data = out.colptr(i); 00066 00067 arrayops::copy( &out_data[i], &A_data[i], N-i ); 00068 } 00069 } 00070 } 00071 00072 00073 if(upper) 00074 { 00075 // reflect elements across the diagonal from upper triangle to lower triangle 00076 00077 for(uword col=1; col < N; ++col) 00078 { 00079 const eT* coldata = out.colptr(col); 00080 00081 for(uword row=0; row < col; ++row) 00082 { 00083 out.at(col,row) = coldata[row]; 00084 } 00085 } 00086 } 00087 else 00088 { 00089 // reflect elements across the diagonal from lower triangle to upper triangle 00090 00091 for(uword col=0; col < N; ++col) 00092 { 00093 const eT* coldata = out.colptr(col); 00094 00095 for(uword row=(col+1); row < N; ++row) 00096 { 00097 out.at(col,row) = coldata[row]; 00098 } 00099 } 00100 } 00101 } 00102 00103 00104 00105 template<typename T1> 00106 inline 00107 void 00108 op_symmat::apply 00109 ( 00110 Mat<typename T1::elem_type>& out, 00111 const Op<T1,op_symmat>& in, 00112 const typename arma_cx_only<typename T1::elem_type>::result* junk 00113 ) 00114 { 00115 arma_extra_debug_sigprint(); 00116 arma_ignore(junk); 00117 00118 typedef typename T1::elem_type eT; 00119 00120 const unwrap<T1> tmp(in.m); 00121 const Mat<eT>& A = tmp.M; 00122 00123 arma_debug_check( (A.is_square() == false), "symmatu()/symmatl(): given matrix must be square" ); 00124 00125 const uword N = A.n_rows; 00126 const bool upper = (in.aux_uword_a == 0); 00127 00128 if(&out != &A) 00129 { 00130 out.copy_size(A); 00131 00132 if(upper) 00133 { 00134 // upper triangular: copy the diagonal and the elements above the diagonal 00135 00136 for(uword i=0; i<N; ++i) 00137 { 00138 const eT* A_data = A.colptr(i); 00139 eT* out_data = out.colptr(i); 00140 00141 arrayops::copy( out_data, A_data, i+1 ); 00142 } 00143 } 00144 else 00145 { 00146 // lower triangular: copy the diagonal and the elements below the diagonal 00147 00148 for(uword i=0; i<N; ++i) 00149 { 00150 const eT* A_data = A.colptr(i); 00151 eT* out_data = out.colptr(i); 00152 00153 arrayops::copy( &out_data[i], &A_data[i], N-i ); 00154 } 00155 } 00156 } 00157 00158 00159 if(upper) 00160 { 00161 // reflect elements across the diagonal from upper triangle to lower triangle 00162 00163 for(uword col=1; col < N; ++col) 00164 { 00165 const eT* coldata = out.colptr(col); 00166 00167 for(uword row=0; row < col; ++row) 00168 { 00169 out.at(col,row) = std::conj(coldata[row]); 00170 } 00171 } 00172 } 00173 else 00174 { 00175 // reflect elements across the diagonal from lower triangle to upper triangle 00176 00177 for(uword col=0; col < N; ++col) 00178 { 00179 const eT* coldata = out.colptr(col); 00180 00181 for(uword row=(col+1); row < N; ++row) 00182 { 00183 out.at(col,row) = std::conj(coldata[row]); 00184 } 00185 } 00186 } 00187 } 00188 00189 00190