fn_norm.hpp

Go to the documentation of this file.
00001 // Copyright (C) 2010 NICTA and the authors listed below
00002 // http://nicta.com.au
00003 // 
00004 // Authors:
00005 // - Conrad Sanderson (conradsand at ieee dot org)
00006 // 
00007 // This file is part of the Armadillo C++ library.
00008 // It is provided without any warranty of fitness
00009 // for any purpose. You can redistribute this file
00010 // and/or modify it under the terms of the GNU
00011 // Lesser General Public License (LGPL) as published
00012 // by the Free Software Foundation, either version 3
00013 // of the License or (at your option) any later version.
00014 // (see http://www.opensource.org/licenses for more info)
00015 
00016 
00017 //! \addtogroup fn_norm
00018 //! @{
00019 
00020 
00021 
00022 template<typename T1>
00023 arma_hot
00024 inline
00025 typename T1::pod_type
00026 norm_unwrap(const Base<typename T1::elem_type,T1>& X, const u32 k)
00027   {
00028   arma_extra_debug_sigprint();
00029   
00030   typedef typename T1::elem_type eT;
00031   typedef typename T1::pod_type   T;
00032   
00033   const unwrap<T1>   tmp(X.get_ref());
00034   const Mat<eT>& A = tmp.M;
00035 
00036   arma_debug_check(    (A.n_elem == 0),                      "norm(): given object has no elements"  );
00037   arma_debug_check( !( (A.n_rows == 1) || (A.n_cols == 1) ), "norm(): given object must be a vector" );
00038   arma_debug_check(    (k == 0),                             "norm(): k must be greater than zero"   );
00039 
00040   const eT* A_mem = A.memptr();
00041   const u32 N     = A.n_elem;
00042 
00043   if(k==1)
00044     {
00045     T acc = T(0);
00046     
00047     for(u32 i=0; i<N; ++i)
00048       {
00049       acc += std::abs(A_mem[i]);
00050       }
00051     
00052     return acc;
00053     }
00054   else
00055   if(k==2)
00056     {
00057     if(is_complex<eT>::value == false)
00058       {
00059       eT acc = eT(0);
00060       
00061       for(u32 i=0; i<N; ++i)
00062         {
00063         const eT tmp = A_mem[i];
00064         acc += tmp*tmp;
00065         }
00066       
00067       return std::sqrt(access::tmp_real(acc));
00068       }
00069     else
00070       {
00071       T acc = T(0);
00072       
00073       for(u32 i=0; i<N; ++i)
00074         {
00075         const T tmp = std::abs(A_mem[i]);
00076         acc += tmp*tmp;
00077         }
00078       
00079       return std::sqrt(acc);
00080       }
00081     }
00082   else
00083     {
00084     T acc = T(0);
00085     
00086     for(u32 i=0; i<N; ++i)
00087       {
00088       acc += std::pow(std::abs(A_mem[i]), int(k));
00089       }
00090     
00091     return std::pow(acc, T(1)/T(k));
00092     }
00093   
00094   }
00095 
00096 
00097 
00098 template<typename T1>
00099 arma_hot
00100 inline
00101 typename T1::pod_type
00102 norm_unwrap(const Base<typename T1::elem_type,T1>& X, const char* method)
00103   {
00104   arma_extra_debug_sigprint();
00105   
00106   typedef typename T1::elem_type eT;
00107   typedef typename T1::pod_type   T;
00108   
00109   const unwrap<T1>   tmp(X.get_ref());
00110   const Mat<eT>& A = tmp.M;
00111 
00112   arma_debug_check(    (A.n_elem == 0),                      "norm(): given object has no elements"  );
00113   arma_debug_check( !( (A.n_rows == 1) || (A.n_cols == 1) ), "norm(): given object must be a vector" );
00114   
00115   const eT* A_mem = A.memptr();
00116   const u32 N     = A.n_elem;
00117   
00118   const char sig = method[0];
00119   
00120   if( (sig == 'i') || (sig == 'I') || (sig == '+') )   // max norm
00121     {
00122     T max_val = std::abs(A_mem[0]);
00123     
00124     for(u32 i=1; i<N; ++i)
00125       {
00126       const T tmp_val = std::abs(A_mem[i]);
00127       
00128       if(tmp_val > max_val)
00129         {
00130         max_val = tmp_val; 
00131         }
00132       }
00133     
00134     return max_val;
00135     }
00136   else
00137   if(sig == '-')   // min norm
00138     {
00139     T min_val = std::abs(A_mem[0]);
00140     
00141     for(u32 i=1; i<N; ++i)
00142       {
00143       const T tmp_val = std::abs(A_mem[i]);
00144       
00145       if(tmp_val < min_val)
00146         {
00147         min_val = tmp_val; 
00148         }
00149       }
00150     
00151     return min_val;
00152     }
00153   else
00154     {
00155     arma_stop("norm(): unknown norm type");
00156     
00157     return T(0);
00158     }
00159   
00160   }
00161 
00162 
00163 
00164 template<typename T1>
00165 arma_hot
00166 inline
00167 typename T1::pod_type
00168 norm_proxy(const Base<typename T1::elem_type,T1>& X, const u32 k)
00169   {
00170   arma_extra_debug_sigprint();
00171   
00172   typedef typename T1::elem_type eT;
00173   typedef typename T1::pod_type   T;
00174   
00175   const Proxy<T1> A(X.get_ref());
00176   
00177   arma_debug_check(    (A.n_elem == 0),                      "norm(): given object has no elements"  );
00178   arma_debug_check( !( (A.n_rows == 1) || (A.n_cols == 1) ), "norm(): given object must be a vector" );
00179   arma_debug_check(    (k == 0),                             "norm(): k must be greater than zero"   );
00180   
00181   const u32 N = A.n_elem;
00182   
00183   if(k==1)
00184     {
00185     T acc = T(0);
00186     
00187     for(u32 i=0; i<N; ++i)
00188       {
00189       acc += std::abs(A[i]);
00190       }
00191     
00192     return acc;
00193     }
00194   else
00195   if(k==2)
00196     {
00197     if(is_complex<eT>::value == false)
00198       {
00199       eT acc = eT(0);
00200       
00201       for(u32 i=0; i<N; ++i)
00202         {
00203         const eT tmp = A[i];
00204         acc += tmp*tmp;
00205         }
00206       
00207       return std::sqrt(access::tmp_real(acc));
00208       }
00209     else
00210       {
00211       T acc = T(0);
00212       
00213       for(u32 i=0; i<N; ++i)
00214         {
00215         const T tmp = std::abs(A[i]);
00216         acc += tmp*tmp;
00217         }
00218       
00219       return std::sqrt(acc);
00220       }
00221     }
00222   else
00223     {
00224     T acc = T(0);
00225     
00226     for(u32 i=0; i<N; ++i)
00227       {
00228       acc += std::pow(std::abs(A[i]), int(k));
00229       }
00230     
00231     return std::pow(acc, T(1)/T(k));
00232     }
00233   
00234   }
00235 
00236 
00237 
00238 template<typename T1>
00239 arma_hot
00240 inline
00241 typename T1::pod_type
00242 norm_proxy(const Base<typename T1::elem_type,T1>& X, const char* method)
00243   {
00244   arma_extra_debug_sigprint();
00245   
00246   typedef typename T1::elem_type eT;
00247   typedef typename T1::pod_type   T;
00248   
00249   const Proxy<T1> A(X.get_ref());
00250   
00251   arma_debug_check(    (A.n_elem == 0),                      "norm(): given object has no elements"  );
00252   arma_debug_check( !( (A.n_rows == 1) || (A.n_cols == 1) ), "norm(): given object must be a vector" );
00253   
00254   const u32 N = A.n_elem;
00255   
00256   const char sig = method[0];
00257   
00258   if( (sig == 'i') || (sig == 'I') || (sig == '+') )   // max norm
00259     {
00260     T max_val = std::abs(A[0]);
00261     
00262     for(u32 i=1; i<N; ++i)
00263       {
00264       const T tmp_val = std::abs(A[i]);
00265       
00266       if(tmp_val > max_val)
00267         {
00268         max_val = tmp_val; 
00269         }
00270       }
00271     
00272     return max_val;
00273     }
00274   else
00275   if(sig == '-')   // min norm
00276     {
00277     T min_val = std::abs(A[0]);
00278     
00279     for(u32 i=1; i<N; ++i)
00280       {
00281       const T tmp_val = std::abs(A[i]);
00282       
00283       if(tmp_val < min_val)
00284         {
00285         min_val = tmp_val; 
00286         }
00287       }
00288     
00289     return min_val;
00290     }
00291   else
00292     {
00293     arma_stop("norm(): unknown norm type");
00294     
00295     return T(0);
00296     }
00297   
00298   }
00299 
00300 
00301 
00302 template<typename T1>
00303 arma_inline
00304 arma_warn_unused
00305 typename T1::pod_type
00306 norm(const Base<typename T1::elem_type,T1>& X, const u32 k)
00307   {
00308   arma_extra_debug_sigprint();
00309   
00310   return (is_Mat<T1>::value == true) ? norm_unwrap(X, k) : norm_proxy(X, k);
00311   }
00312 
00313 
00314 
00315 template<typename T1>
00316 arma_inline
00317 arma_warn_unused
00318 typename T1::pod_type
00319 norm(const Base<typename T1::elem_type,T1>& X, const char* method)
00320   {
00321   arma_extra_debug_sigprint();
00322   
00323   return (is_Mat<T1>::value == true) ? norm_unwrap(X, method) : norm_proxy(X, method);
00324   }
00325 
00326 
00327 
00328 //! @}