fn_misc.hpp

Go to the documentation of this file.
00001 // Copyright (C) 2009 NICTA
00002 // 
00003 // Authors:
00004 // - Conrad Sanderson (conradsand at ieee dot org)
00005 // 
00006 // This file is part of the Armadillo C++ library.
00007 // It is provided without any warranty of fitness
00008 // for any purpose. You can redistribute this file
00009 // and/or modify it under the terms of the GNU
00010 // Lesser General Public License (LGPL) as published
00011 // by the Free Software Foundation, either version 3
00012 // of the License or (at your option) any later version.
00013 // (see http://www.opensource.org/licenses for more info)
00014 
00015 
00016 //! \addtogroup fn_misc
00017 //! @{
00018 
00019 //! \brief
00020 //! Generate a vector with 'num' elements.
00021 //! The values of the elements linearly increase from 'start' upto (and including) 'end'.
00022 
00023 template<typename eT>
00024 inline
00025 Mat<eT>
00026 linspace(const eT start, const eT end, const u32 num, const u32 dim = 0)
00027   {
00028   arma_extra_debug_sigprint();
00029   
00030   arma_debug_check( (num < 2), "linspace(): num must be >= 2");
00031   
00032   Mat<eT> x;
00033   
00034   if(dim == 0)
00035     {
00036     x.set_size(num,1);  // column vector
00037     }
00038   else
00039     {
00040     x.set_size(1,num);  // row vector
00041     }
00042   
00043   
00044   const eT delta = (end-start)/(num-1);
00045   
00046   x[0] = start;
00047   
00048   for(u32 i=1; i<num; ++i)
00049     {
00050     x[i] = x[i-1] + delta;
00051     }
00052   
00053   return x; 
00054   }
00055 
00056 
00057 
00058 inline
00059 mat
00060 linspace(const double start, const double end, const u32 num, const u32 dim = 0)
00061   {
00062   arma_extra_debug_sigprint();
00063   return linspace<double>(start, end, num, dim);
00064   }
00065 
00066 
00067 
00068 //
00069 // real
00070 
00071 template<typename T, typename T1>
00072 inline
00073 Mat<T>
00074 real(const Base<std::complex<T>, T1>& X)
00075   {
00076   arma_extra_debug_sigprint();
00077   
00078   typedef std::complex<T> eT;
00079 
00080   const unwrap<T1> A_tmp(X.get_ref());
00081   const Mat<eT>& A = A_tmp.M;
00082   
00083   Mat<T> out(A.n_rows, A.n_cols);
00084   
00085   const eT* A_mem = A.mem;
00086   T* out_mem = out.memptr();
00087   
00088   for(u32 i=0; i<out.n_elem; ++i)
00089     {
00090     out_mem[i] = std::real(A_mem[i]);
00091     }
00092   
00093   return out;
00094   }
00095 
00096 
00097 
00098 template<typename T, typename T1>
00099 inline
00100 Cube<T>
00101 real(const BaseCube<std::complex<T>, T1>& X)
00102   {
00103   arma_extra_debug_sigprint();
00104   
00105   typedef std::complex<T> eT;
00106 
00107   const unwrap_cube<T1> A_tmp(X.get_ref());
00108   const Cube<eT>& A   = A_tmp.M;
00109   
00110   Cube<T> out(A.n_rows, A.n_cols, A.n_slices);
00111   
00112   const eT* A_mem = A.mem;
00113   T* out_mem = out.memptr();
00114   
00115   for(u32 i=0; i<out.n_elem; ++i)
00116     {
00117     out_mem[i] = std::real(A_mem[i]);
00118     }
00119   
00120   return out;
00121   }
00122 
00123 
00124 
00125 //
00126 // imag
00127 
00128 template<typename T, typename T1>
00129 inline
00130 Mat<T>
00131 imag(const Base<std::complex<T>,T1>& X)
00132   {
00133   arma_extra_debug_sigprint();
00134   
00135   typedef std::complex<T> eT;
00136 
00137   const unwrap<T1> A_tmp(X.get_ref());
00138   const Mat<eT>& A = A_tmp.M;
00139   
00140   Mat<T> out(A.n_rows, A.n_cols);
00141   
00142   const eT* A_mem = A.mem;
00143   T* out_mem = out.memptr();
00144   
00145   for(u32 i=0; i<out.n_elem; ++i)
00146     {
00147     out_mem[i] = std::imag(A_mem[i]);
00148     }
00149   
00150   return out;
00151   }
00152 
00153 
00154 
00155 template<typename T, typename T1>
00156 inline
00157 Cube<T>
00158 imag(const BaseCube<std::complex<T>,T1>& X)
00159   {
00160   arma_extra_debug_sigprint();
00161   
00162   typedef std::complex<T> eT;
00163 
00164   const unwrap_cube<T1> A_tmp(X.get_ref());
00165   const Cube<eT>& A   = A_tmp.M;
00166   
00167   Cube<T> out(A.n_rows, A.n_cols, A.n_slices);
00168   
00169   const eT* A_mem = A.mem;
00170   T* out_mem = out.memptr();
00171   
00172   for(u32 i=0; i<out.n_elem; ++i)
00173     {
00174     out_mem[i] = std::imag(A_mem[i]);
00175     }
00176   
00177   return out;
00178   }
00179 
00180 
00181 
00182 //
00183 // log_add
00184 
00185 template<typename eT>
00186 inline
00187 eT
00188 log_add(eT log_a, eT log_b)
00189   {
00190   if(log_a < log_b)
00191     {
00192     std::swap(log_a, log_b);
00193     }
00194   
00195   const eT negdelta = log_b - log_a;
00196   
00197   if( (negdelta < Math<eT>::log_min()) || (arma_isfinite(negdelta) == false) )
00198     {
00199     return log_a;
00200     }
00201   else
00202     {
00203     #if defined(ARMA_HAVE_LOG1P)
00204       return (log_a + log1p(std::exp(negdelta)));
00205     #else
00206       return (log_a + std::log(1.0 + std::exp(negdelta)));
00207     #endif
00208     }
00209   }
00210 
00211 
00212 
00213 template<typename eT>
00214 inline
00215 eT
00216 trunc_log(const eT x)
00217   {
00218   if(std::numeric_limits<eT>::is_iec559)
00219     {
00220     if(x == std::numeric_limits<eT>::infinity())
00221       {
00222       return Math<eT>::log_max();
00223       }
00224     if(x <= 0)
00225       {
00226       return Math<eT>::log_min();
00227       }
00228     else
00229       {
00230       return std::log(x);
00231       }
00232     }
00233   else
00234     {
00235     return std::log(x);
00236     }
00237   }
00238 
00239 
00240 
00241 template<typename eT>
00242 inline
00243 eT
00244 trunc_exp(const eT x)
00245   {
00246   if(std::numeric_limits<eT>::is_iec559 && (x >= Math<eT>::log_max() ))
00247     {
00248     return std::numeric_limits<eT>::max();
00249     }
00250   else
00251     {
00252     return std::exp(x);
00253     }
00254   }
00255 
00256 
00257 
00258 //
00259 // log
00260 
00261 template<typename T1>
00262 inline
00263 const Op<T1, op_log>
00264 log(const Base<typename T1::elem_type,T1>& A)
00265   {
00266   arma_extra_debug_sigprint();
00267   
00268   return Op<T1, op_log>(A.get_ref());
00269   }
00270 
00271 
00272 
00273 template<typename T1>
00274 inline
00275 const OpCube<T1, op_log>
00276 log(const BaseCube<typename T1::elem_type,T1>& A)
00277   {
00278   arma_extra_debug_sigprint();
00279   
00280   return OpCube<T1, op_log>(A.get_ref());
00281   }
00282 
00283 
00284 
00285 //
00286 // trunc_log
00287 
00288 template<typename T1>
00289 inline
00290 const Op<T1, op_trunc_log>
00291 trunc_log(const Base<typename T1::elem_type,T1>& A)
00292   {
00293   arma_extra_debug_sigprint();
00294   
00295   return Op<T1, op_trunc_log>(A.get_ref());
00296   }
00297 
00298 
00299 
00300 template<typename T1>
00301 inline
00302 const OpCube<T1, op_trunc_log>
00303 trunc_log(const BaseCube<typename T1::elem_type,T1>& A)
00304   {
00305   arma_extra_debug_sigprint();
00306   
00307   return OpCube<T1, op_trunc_log>(A.get_ref());
00308   }
00309 
00310 
00311 
00312 //
00313 // log10
00314 
00315 template<typename T1>
00316 inline
00317 const Op<T1, op_log10>
00318 log10(const Base<typename T1::elem_type,T1>& A)
00319   {
00320   arma_extra_debug_sigprint();
00321   
00322   return Op<T1, op_log10>(A.get_ref());
00323   }
00324 
00325 
00326 
00327 template<typename T1>
00328 inline
00329 const OpCube<T1, op_log10>
00330 log10(const BaseCube<typename T1::elem_type,T1>& A)
00331   {
00332   arma_extra_debug_sigprint();
00333   
00334   return OpCube<T1, op_log10>(A.get_ref());
00335   }
00336 
00337 
00338 
00339 //
00340 // exp
00341 
00342 template<typename T1>
00343 inline
00344 const Op<T1, op_exp>
00345 exp(const Base<typename T1::elem_type,T1>& A)
00346   {
00347   arma_extra_debug_sigprint();
00348   
00349   return Op<T1, op_exp>(A.get_ref());
00350   }
00351 
00352 
00353 
00354 template<typename T1>
00355 inline
00356 const OpCube<T1, op_exp>
00357 exp(const BaseCube<typename T1::elem_type,T1>& A)
00358   {
00359   arma_extra_debug_sigprint();
00360   
00361   return OpCube<T1, op_exp>(A.get_ref());
00362   }
00363 
00364 
00365 
00366 //
00367 // trunc_exp
00368 
00369 template<typename T1>
00370 inline
00371 const Op<T1, op_trunc_exp>
00372 trunc_exp(const Base<typename T1::elem_type,T1>& A)
00373   {
00374   arma_extra_debug_sigprint();
00375   
00376   return Op<T1, op_trunc_exp>(A.get_ref());
00377   }
00378 
00379 
00380 
00381 template<typename T1>
00382 inline
00383 const OpCube<T1, op_trunc_exp>
00384 trunc_exp(const BaseCube<typename T1::elem_type,T1>& A)
00385   {
00386   arma_extra_debug_sigprint();
00387   
00388   return OpCube<T1, op_trunc_exp>(A.get_ref());
00389   }
00390 
00391 
00392 
00393 //
00394 // abs
00395 
00396 template<typename T1>
00397 inline
00398 Mat<typename T1::pod_type>
00399 abs(const Base<typename T1::elem_type,T1>& X)
00400   {
00401   arma_extra_debug_sigprint();
00402   
00403   const unwrap<T1> A_tmp(X.get_ref());
00404 
00405   // if T1 is a complex matrix,
00406   // pod_type is the underlying type used by std::complex;
00407   // otherwise pod_type is the same as elem_type
00408   
00409   typedef typename T1::elem_type  in_eT;
00410   typedef typename T1::pod_type  out_eT;
00411 
00412   const Mat<in_eT>& A = A_tmp.M;
00413   
00414   Mat<out_eT> out(A.n_rows, A.n_cols);
00415   
00416   const in_eT* A_mem   = A.mem;
00417   out_eT*      out_mem = out.memptr();
00418   
00419   for(u32 i=0; i<out.n_elem; ++i)
00420     {
00421     out_mem[i] = std::abs(A_mem[i]);
00422     }
00423   
00424   return out;
00425   }
00426 
00427 
00428 
00429 template<typename T1>
00430 inline
00431 Cube<typename T1::pod_type>
00432 abs(const BaseCube<typename T1::elem_type,T1>& X)
00433   {
00434   arma_extra_debug_sigprint();
00435   
00436   const unwrap_cube<T1> A_tmp(X.get_ref());
00437 
00438   // if T1 is a complex matrix,
00439   // pod_type is the underlying type used by std::complex;
00440   // otherwise pod_type is the same as elem_type
00441   
00442   typedef typename T1::elem_type  in_eT;
00443   typedef typename T1::pod_type  out_eT;
00444 
00445   const Cube<in_eT>& A = A_tmp.M;
00446   
00447   Cube<out_eT> out(A.n_rows, A.n_cols, A.n_slices);
00448   
00449   const in_eT* A_mem   = A.mem;
00450   out_eT*      out_mem = out.memptr();
00451   
00452   for(u32 i=0; i<out.n_elem; ++i)
00453     {
00454     out_mem[i] = std::abs(A_mem[i]);
00455     }
00456   
00457   return out;
00458   }
00459 
00460 
00461 
00462 //
00463 // fabs
00464 
00465 template<typename T1>
00466 inline
00467 Mat<typename T1::pod_type>
00468 fabs(const Base<typename T1::elem_type,T1>& A)
00469   {
00470   arma_extra_debug_sigprint();
00471   
00472   return abs(A);
00473   }
00474 
00475 
00476 
00477 template<typename T1>
00478 inline
00479 Cube<typename T1::pod_type>
00480 fabs(const BaseCube<typename T1::elem_type,T1>& A)
00481   {
00482   arma_extra_debug_sigprint();
00483   
00484   return abs(A);
00485   }
00486 
00487 
00488 
00489 //
00490 // square
00491 
00492 template<typename T1>
00493 inline
00494 const Op<T1, op_square>
00495 square(const Base<typename T1::elem_type,T1>& A)
00496   {
00497   arma_extra_debug_sigprint();
00498   
00499   return Op<T1, op_square>(A.get_ref());
00500   }
00501 
00502 
00503 
00504 template<typename T1>
00505 inline
00506 const OpCube<T1, op_square>
00507 square(const BaseCube<typename T1::elem_type,T1>& A)
00508   {
00509   arma_extra_debug_sigprint();
00510   
00511   return OpCube<T1, op_square>(A.get_ref());
00512   }
00513 
00514 
00515 
00516 //
00517 // sqrt
00518 
00519 template<typename T1>
00520 inline
00521 const Op<T1, op_sqrt>
00522 sqrt(const Base<typename T1::elem_type,T1>& A)
00523   {
00524   arma_extra_debug_sigprint();
00525   
00526   return Op<T1, op_sqrt>(A.get_ref());
00527   }
00528 
00529 
00530 
00531 template<typename T1>
00532 inline
00533 const OpCube<T1, op_sqrt>
00534 sqrt(const BaseCube<typename T1::elem_type,T1>& A)
00535   {
00536   arma_extra_debug_sigprint();
00537   
00538   return OpCube<T1, op_sqrt>(A.get_ref());
00539   }
00540 
00541 
00542 
00543 // pow
00544 
00545 template<typename T1>
00546 inline
00547 const Op<T1, op_pow>
00548 pow(const Base<typename T1::elem_type,T1>& A, const typename T1::elem_type exponent)
00549   {
00550   arma_extra_debug_sigprint();
00551   
00552   return Op<T1, op_pow>(A.get_ref(), exponent);
00553   }
00554 
00555 
00556 
00557 template<typename T1>
00558 inline
00559 const OpCube<T1, op_pow>
00560 pow(const BaseCube<typename T1::elem_type,T1>& A, const typename T1::elem_type exponent)
00561   {
00562   arma_extra_debug_sigprint();
00563   
00564   return OpCube<T1, op_pow>(A.get_ref(), exponent);
00565   }
00566 
00567 
00568 
00569 // pow, specialised handling (non-complex exponent for complex matrices)
00570 
00571 template<typename T1>
00572 inline
00573 const Op<T1, op_pow>
00574 pow(const Base<typename T1::elem_type,T1>& A, const typename T1::elem_type::value_type exponent)
00575   {
00576   arma_extra_debug_sigprint();
00577   
00578   return Op<T1, op_pow>(A.get_ref(), eT(exponent));
00579   }
00580 
00581 
00582 
00583 template<typename T1>
00584 inline
00585 const OpCube<T1, op_pow>
00586 pow(const BaseCube<typename T1::elem_type,T1>& A, const typename T1::elem_type::value_type exponent)
00587   {
00588   arma_extra_debug_sigprint();
00589   
00590   return OpCube<T1, op_pow>(A.get_ref(), eT(exponent));
00591   }
00592 
00593 
00594 
00595 #if defined(ARMA_GOOD_COMPILER)
00596 
00597 
00598 // pow_s32  (integer exponent)
00599 
00600 template<typename T1>
00601 inline
00602 const Op<T1, op_pow_s32>
00603 pow(const Base<typename T1::elem_type,T1>& A, const s32 exponent)
00604   {
00605   arma_extra_debug_sigprint();
00606   
00607   if(exponent >= 0)
00608     {
00609     return Op<T1, op_pow_s32>(A.get_ref(), exponent, 0);
00610     }
00611   else
00612     {
00613     return Op<T1, op_pow_s32>(A.get_ref(), -exponent, 1);
00614     }
00615   }
00616 
00617 
00618 
00619 template<typename T1>
00620 inline
00621 const OpCube<T1, op_pow_s32>
00622 pow(const BaseCube<typename T1::elem_type,T1>& A, const s32 exponent)
00623   {
00624   arma_extra_debug_sigprint();
00625   
00626   if(exponent >= 0)
00627     {
00628     return OpCube<T1, op_pow_s32>(A.get_ref(), exponent, 0);
00629     }
00630   else
00631     {
00632     return OpCube<T1, op_pow_s32>(A.get_ref(), -exponent, 1);
00633     }
00634   }
00635 
00636 
00637 
00638 #endif
00639 
00640 
00641 
00642 // conj
00643 
00644 template<typename T, typename T1>
00645 inline
00646 const Op<T1, op_conj>
00647 conj(const Base<std::complex<T>,T1>& A)
00648   {
00649   arma_extra_debug_sigprint();
00650 
00651   return Op<T1, op_conj>(A.get_ref());
00652   }
00653 
00654 
00655 
00656 template<typename T, typename T1>
00657 inline
00658 const OpCube<T1, op_conj>
00659 conj(const BaseCube<std::complex<T>,T1>& A)
00660   {
00661   arma_extra_debug_sigprint();
00662 
00663   return OpCube<T1, op_conj>(A.get_ref());
00664   }
00665 
00666 
00667 
00668 template<typename T1>
00669 inline
00670 const T1&
00671 conj(const Op<T1, op_conj>& A)
00672   {
00673   arma_extra_debug_sigprint();
00674   
00675   return A.m;
00676   }
00677 
00678 
00679 
00680 template<typename T1>
00681 inline
00682 const T1&
00683 conj(const OpCube<T1, op_conj>& A)
00684   {
00685   arma_extra_debug_sigprint();
00686   
00687   return A.m;
00688   }
00689 
00690 
00691 //! the conjugate of the transpose of a complex matrix is the same as the hermitian transpose
00692 template<typename T1>
00693 inline
00694 const Op<T1, op_htrans>
00695 conj(const Op<T1, op_trans>& A)
00696   {
00697   arma_extra_debug_sigprint();
00698   
00699   arma_type_check< is_complex<typename T1::elem_type>::value == false >::apply();
00700 
00701   return Op<T1, op_htrans>(A.m);
00702   }
00703 
00704 
00705 //! @}