IT++ Logo Newcom Logo

matfunc.cpp

Go to the documentation of this file.
00001 
00033 #include <itpp/base/matfunc.h>
00034 #include <itpp/base/schur.h>
00035 #include <itpp/base/converters.h>
00036 #include <limits>
00037 
00038 
00039 namespace itpp {
00040 
00041   // Square root of a square matrix (based on Octave sqrtm implementation)
00042   cmat sqrtm(const mat& A)
00043   {
00044     return sqrtm(to_cmat(A));
00045   }
00046 
00047   // Square root of the complex square matrix A
00048   cmat sqrtm(const cmat& A)
00049   {
00050     cmat U, T;
00051     schur(A, U, T);
00052 
00053     int n = U.rows();
00054     cmat R(n, n);
00055 
00056     R.zeros();
00057     for (int j = 0; j < n; j++)
00058       R(j, j) = std::sqrt(T(j, j));
00059 
00060     const double fudge = std::sqrt(std::numeric_limits<double>::min());
00061 
00062     for (int p = 0; p < n - 1; p++) {
00063       for (int i = 0; i < n - (p + 1); i++) {
00064         const int j = i + p + 1;
00065         std::complex<double> s = T(i, j);
00066 
00067         for (int k = i + 1; k < j; k++)
00068           s -= R(i, k) * R(k, j);
00069 
00070         const std::complex<double> d = R(i, i) + R(j, j) + fudge;
00071         const std::complex<double> conj_d = conj(d);
00072 
00073         R(i, j) = (s * conj_d) / (d * conj_d);
00074       }
00075     }
00076 
00077     return U * R * U.H();
00078   }
00079 
00080 
00081   // ---------------------- Instantiations ------------------------------------
00082 
00083   template int length(const vec &v);
00084   template int length(const cvec &v);
00085   template int length(const svec &v);
00086   template int length(const ivec &v);
00087   template int length(const bvec &v);
00088 
00089 
00090   template double sum(const vec &v);
00091   template std::complex<double> sum(const cvec &v);
00092   template short sum(const svec &v);
00093   template int sum(const ivec &v);
00094   template bin sum(const bvec &v);
00095 
00096   template double sum_sqr(const vec &v);
00097   template std::complex<double> sum_sqr(const cvec &v);
00098   template short sum_sqr(const svec &v);
00099   template int sum_sqr(const ivec &v);
00100   template bin sum_sqr(const bvec &v);
00101 
00102   template vec cumsum(const vec &v);
00103   template cvec cumsum(const cvec &v);
00104   template svec cumsum(const svec &v);
00105   template ivec cumsum(const ivec &v);
00106   template bvec cumsum(const bvec &v);
00107 
00108   template double prod(const vec &v);
00109   template std::complex<double> prod(const cvec &v);
00110   template short prod(const svec &v);
00111   template int prod(const ivec &v);
00112   template bin prod(const bvec &v);
00113 
00114   template vec cross(const vec &v1, const vec &v2);
00115   template ivec cross(const ivec &v1, const ivec &v2);
00116   template svec cross(const svec &v1, const svec &v2);
00117 
00118   template vec reverse(const vec &in);
00119   template cvec reverse(const cvec &in);
00120   template svec reverse(const svec &in);
00121   template ivec reverse(const ivec &in);
00122   template bvec reverse(const bvec &in);
00123 
00124   template vec repeat(const vec &v, int norepeats);
00125   template cvec repeat(const cvec &v, int norepeats);
00126   template svec repeat(const svec &v, int norepeats);
00127   template ivec repeat(const ivec &v, int norepeats);
00128   template bvec repeat(const bvec &v, int norepeats);
00129 
00130   template vec apply_function(float (*f)(float), const vec &data);
00131   template vec apply_function(double (*f)(double), const vec &data);
00132   template cvec apply_function(std::complex<double> (*f)(std::complex<double>), const cvec &data);
00133   template svec apply_function(short (*f)(short), const svec &data);
00134   template ivec apply_function(int (*f)(int), const ivec &data);
00135   template bvec apply_function(bin (*f)(bin), const bvec &data);
00136 
00137 
00138   template ivec zero_pad(const ivec &v, int n);
00139   template vec zero_pad(const vec &v, int n);
00140   template cvec zero_pad(const cvec &v, int n);
00141   template bvec zero_pad(const bvec &v, int n);
00142 
00143   template ivec zero_pad(const ivec &v);
00144   template vec zero_pad(const vec &v);
00145   template cvec zero_pad(const cvec &v);
00146   template bvec zero_pad(const bvec &v);
00147 
00148   template mat  zero_pad(const mat &, int, int);
00149   template cmat zero_pad(const cmat &, int, int);
00150   template imat zero_pad(const imat &, int, int);
00151   template bmat zero_pad(const bmat &, int, int);
00152 
00153   template vec sum(const mat &m, int dim);
00154   template cvec sum(const cmat &m, int dim);
00155   template svec sum(const smat &m, int dim);
00156   template ivec sum(const imat &m, int dim);
00157   template bvec sum(const bmat &m, int dim);
00158 
00159   template vec sum_sqr(const mat & m, int dim);
00160   template cvec sum_sqr(const cmat &m, int dim);
00161   template svec sum_sqr(const smat &m, int dim);
00162   template ivec sum_sqr(const imat &m, int dim);
00163   template bvec sum_sqr(const bmat &m, int dim);
00164 
00165   template mat cumsum(const mat &m, int dim);
00166   template cmat cumsum(const cmat &m, int dim);
00167   template smat cumsum(const smat &m, int dim);
00168   template imat cumsum(const imat &m, int dim);
00169   template bmat cumsum(const bmat &m, int dim);
00170 
00171   template vec prod(const mat &m, int dim);
00172   // Template instantiation of product
00173   template cvec prod(const cmat &v, int dim);
00174   template svec prod(const smat &m, int dim);
00175   template ivec prod(const imat &m, int dim);
00176 
00177   template vec diag(const mat &in);
00178   template cvec diag(const cmat &in);
00179 
00180   template void diag(const vec &in, mat &m);
00181   template void diag(const cvec &in, cmat &m);
00182 
00183   template mat diag(const vec &v, const int K);
00184   template cmat diag(const cvec &v, const int K);
00185 
00186   template mat bidiag(const vec &, const vec &);
00187   template cmat bidiag(const cvec &, const cvec &);
00188 
00189   template void bidiag(const vec &, const vec &, mat &);
00190   template void bidiag(const cvec &, const cvec &, cmat &);
00191 
00192   template void bidiag(const mat &, vec &, vec &);
00193   template void bidiag(const cmat &, cvec &, cvec &);
00194 
00195   template mat tridiag(const vec &main, const vec &, const vec &);
00196   template cmat tridiag(const cvec &main, const cvec &, const cvec &);
00197 
00198   template void tridiag(const vec &main, const vec &, const vec &, mat &);
00199   template void tridiag(const cvec &main, const cvec &, const cvec &, cmat &);
00200 
00201   template void tridiag(const mat &m, vec &, vec &, vec &);
00202   template void tridiag(const cmat &m, cvec &, cvec &, cvec &);
00203 
00204   template double trace(const mat &in);
00205   template std::complex<double> trace(const cmat &in);
00206   template short trace(const smat &in);
00207   template int trace(const imat &in);
00208   template bin trace(const bmat &in);
00209 
00210   template void transpose(const mat &m, mat &out);
00211   template void transpose(const cmat &m, cmat &out);
00212   template void transpose(const smat &m, smat &out);
00213   template void transpose(const imat &m, imat &out);
00214   template void transpose(const bmat &m, bmat &out);
00215 
00216   template mat transpose(const mat &m);
00217   template cmat transpose(const cmat &m);
00218   template smat transpose(const smat &m);
00219   template imat transpose(const imat &m);
00220   template bmat transpose(const bmat &m);
00221 
00222 
00223   template void hermitian_transpose(const mat &m, mat &out);
00224   template void hermitian_transpose(const cmat &m, cmat &out);
00225   template void hermitian_transpose(const smat &m, smat &out);
00226   template void hermitian_transpose(const imat &m, imat &out);
00227   template void hermitian_transpose(const bmat &m, bmat &out);
00228 
00229   template mat hermitian_transpose(const mat &m);
00230   template cmat hermitian_transpose(const cmat &m);
00231   template smat hermitian_transpose(const smat &m);
00232   template imat hermitian_transpose(const imat &m);
00233   template bmat hermitian_transpose(const bmat &m);
00234 
00235   template mat repeat(const mat &m, int norepeats);
00236   template cmat repeat(const cmat &m, int norepeats);
00237   template smat repeat(const smat &m, int norepeats);
00238   template imat repeat(const imat &m, int norepeats);
00239   template bmat repeat(const bmat &m, int norepeats);
00240 
00241   template mat apply_function(float (*f)(float), const mat &data);
00242   template mat apply_function(double (*f)(double), const mat &data);
00243   template cmat apply_function(std::complex<double> (*f)(std::complex<double>), const cmat &data);
00244   template smat apply_function(short (*f)(short), const smat &data);
00245   template imat apply_function(int (*f)(int), const imat &data);
00246   template bmat apply_function(bin (*f)(bin), const bmat &data);
00247 
00248   template  vec rvectorize(const  mat &m);
00249   template cvec rvectorize(const cmat &m);
00250   template  ivec rvectorize(const  imat &m);
00251   template  bvec rvectorize(const  bmat &m);
00252 
00253   template  vec cvectorize(const  mat &m);
00254   template cvec cvectorize(const cmat &m);
00255   template  ivec cvectorize(const  imat &m);
00256   template  bvec cvectorize(const  bmat &m);
00257 
00258   template  mat reshape(const  mat &m, int rows, int cols);
00259   template cmat reshape(const cmat &m, int rows, int cols);
00260   template  imat reshape(const  imat &m, int rows, int cols);
00261   template  bmat reshape(const  bmat &m, int rows, int cols);
00262 
00263   template  mat reshape(const  vec &m, int rows, int cols);
00264   template cmat reshape(const cvec &m, int rows, int cols);
00265   template  imat reshape(const  ivec &m, int rows, int cols);
00266   template  bmat reshape(const  bvec &m, int rows, int cols);
00267 
00268   template vec upsample(const vec &v, int usf);
00269   template cvec upsample(const cvec &v, int usf);
00270   template svec upsample(const svec &v, int usf);
00271   template ivec upsample(const ivec &v, int usf);
00272   template bvec upsample(const bvec &v, int usf);
00273 
00274   template mat upsample(const mat &v, int usf);
00275   template cmat upsample(const cmat &v, int usf);
00276   template smat upsample(const smat &v, int usf);
00277   template imat upsample(const imat &v, int usf);
00278   template bmat upsample(const bmat &v, int usf);
00279 
00280   template void upsample(const vec &v, int usf,  vec & u);
00281   template void upsample(const cvec &v, int usf,  cvec & u);
00282   template void upsample(const svec &v, int usf,  svec & u);
00283   template void upsample(const ivec &v, int usf,  ivec & u);
00284   template void upsample(const bvec &v, int usf,  bvec & u);
00285 
00286   template void upsample(const mat &v, int usf,  mat & u);
00287   template void upsample(const cmat &v, int usf,  cmat & u);
00288   template void upsample(const smat &v, int usf,  smat & u);
00289   template void upsample(const imat &v, int usf,  imat & u);
00290   template void upsample(const bmat &v, int usf,  bmat & u);
00291  
00292   template vec lininterp(const vec &v, int usf);
00293   template cvec lininterp(const cvec &v, int usf);
00294 
00295   template mat lininterp(const mat &v, int usf);
00296   template cmat lininterp(const cmat &v, int usf);
00297 
00298   template void lininterp(const vec &v, int usf,  vec & u);
00299   template void lininterp(const cvec &v, int usf,  cvec & u);
00300 
00301   template void lininterp(const mat &v, int usf,  mat & u);
00302   template void lininterp(const cmat &v, int usf,  cmat & u);
00303 
00304   template mat lininterp(const mat &m, const double f_base, const double f_ups, const int nrof_samples, const double t_start);
00305   template cmat lininterp(const cmat &m, const double f_base, const double f_ups, const int nrof_samples, const double t_start);
00306 
00307   template vec lininterp(const vec &v, const double f_base, const double f_ups, const int nrof_samples, const double t_start);
00308   template cvec lininterp(const cvec &v, const double f_base, const double f_ups, const int nrof_samples, const double t_start);
00309 
00310 } // namespace itpp
SourceForge Logo

Generated on Thu Apr 19 14:18:29 2007 for IT++ by Doxygen 1.5.1