00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 template<typename T1>
00029 inline
00030 void
00031 eig_sym(Col<typename T1::pod_type>& eigval, const Base<typename T1::elem_type,T1>& X, const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0)
00032 {
00033 arma_extra_debug_sigprint();
00034
00035 typedef typename T1::elem_type eT;
00036
00037
00038
00039
00040 const unwrap<T1> tmp(X.get_ref());
00041 const Mat<eT>& A = tmp.M;
00042
00043 auxlib::eig_sym(eigval, A);
00044 }
00045
00046
00047
00048
00049 template<typename T1>
00050 inline
00051 Col<typename T1::pod_type>
00052 eig_sym(const Base<typename T1::elem_type,T1>& X, const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0)
00053 {
00054 arma_extra_debug_sigprint();
00055
00056 Col<typename T1::pod_type> out;
00057 eig_sym(out, X);
00058
00059 return out;
00060 }
00061
00062
00063
00064 template<typename T1>
00065 inline
00066 void
00067 eig_sym
00068 (
00069 Col<typename T1::pod_type>& eigval,
00070 Mat<typename T1::elem_type>& eigvec,
00071 const Base<typename T1::elem_type,T1>& X,
00072 const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0
00073 )
00074 {
00075 arma_extra_debug_sigprint();
00076
00077 typedef typename T1::elem_type eT;
00078
00079 const unwrap<T1> tmp(X.get_ref());
00080 const Mat<eT>& A = tmp.M;
00081
00082 auxlib::eig_sym(eigval, eigvec, A);
00083 }
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094 template<typename T1>
00095 inline
00096 void
00097 eig_gen
00098 (
00099 Col< std::complex<typename T1::pod_type> >& eigval,
00100 Mat<typename T1::elem_type>& l_eigvec,
00101 Mat<typename T1::elem_type>& r_eigvec,
00102 const Base<typename T1::elem_type,T1>& X,
00103 const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0
00104 )
00105 {
00106 arma_extra_debug_sigprint();
00107
00108 typedef typename T1::elem_type eT;
00109
00110 const unwrap<T1> tmp(X.get_ref());
00111 const Mat<eT>& A = tmp.M;
00112
00113 auxlib::eig_gen(eigval, l_eigvec, r_eigvec, A, 'b');
00114 }
00115
00116
00117
00118
00119
00120
00121 template<typename eT, typename T1>
00122 inline
00123 void
00124 eig_gen
00125 (
00126 Col< std::complex<eT> >& eigval,
00127 Mat< std::complex<eT> >& eigvec,
00128 const Base<eT, T1>& X,
00129 const char side = 'r',
00130 const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0
00131 )
00132 {
00133 arma_extra_debug_sigprint();
00134
00135
00136
00137 const unwrap<T1> tmp(X.get_ref());
00138 const Mat<eT>& A = tmp.M;
00139
00140 Mat<eT> dummy_eigvec;
00141 Mat<eT> tmp_eigvec;
00142
00143 switch(side)
00144 {
00145 case 'r':
00146 auxlib::eig_gen(eigval, dummy_eigvec, tmp_eigvec, A, side);
00147 break;
00148
00149 case 'l':
00150 auxlib::eig_gen(eigval, tmp_eigvec, dummy_eigvec, A, side);
00151 break;
00152
00153 default:
00154 arma_stop("eig_gen(): parameter 'side' is invalid");
00155 return;
00156 }
00157
00158
00159 const u32 n = A.n_rows;
00160
00161 if(n > 0)
00162 {
00163 eigvec.set_size(n,n);
00164
00165 for(u32 j=0; j<n; ++j)
00166 {
00167 if( (j < n-1) && (eigval[j] == std::conj(eigval[j+1])) )
00168 {
00169
00170
00171
00172 for(u32 i=0; i<n; ++i)
00173 {
00174 eigvec.at(i,j) = std::complex<eT>( tmp_eigvec.at(i,j), tmp_eigvec.at(i,j+1) );
00175 eigvec.at(i,j+1) = std::complex<eT>( tmp_eigvec.at(i,j), -tmp_eigvec.at(i,j+1) );
00176 }
00177
00178 ++j;
00179 }
00180 else
00181 {
00182
00183
00184 for(u32 i=0; i<n; ++i)
00185 {
00186 eigvec.at(i,j) = std::complex<eT>(tmp_eigvec.at(i,j), eT(0));
00187 }
00188
00189 }
00190 }
00191 }
00192
00193 }
00194
00195
00196
00197
00198
00199
00200 template<typename T, typename T1>
00201 inline
00202 void
00203 eig_gen
00204 (
00205 Col< std::complex<T> >& eigval,
00206 Mat< std::complex<T> >& eigvec,
00207 const Base<std::complex<T>, T1>& X,
00208 const char side = 'r',
00209 const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0
00210 )
00211 {
00212 arma_extra_debug_sigprint();
00213
00214
00215 typedef typename std::complex<T> eT;
00216
00217 const unwrap<T1> tmp(X.get_ref());
00218 const Mat<eT>& A = tmp.M;
00219
00220 Mat<eT> dummy_eigvec;
00221
00222 switch(side)
00223 {
00224 case 'r':
00225 auxlib::eig_gen(eigval, dummy_eigvec, eigvec, A, side);
00226 break;
00227
00228 case 'l':
00229 auxlib::eig_gen(eigval, eigvec, dummy_eigvec, A, side);
00230 break;
00231
00232 default:
00233 arma_stop("eig_gen(): parameter 'side' is invalid");
00234 }
00235 }
00236
00237
00238
00239
00240