Arageli: smith.hpp Source File

00001 /*****************************************************************************
00002     
00003     smith.hpp -- Smith Normal Diagonal Form of a matrix.
00004 
00005     This file is a part of the Arageli library.
00006 
00007     Copyright (C) Nikolai Yu. Zolotykh, 1999--2006
00008     Copyright (C) Sergey S. Lyalin, 2005--2006
00009     Copyright (C) Anna Bestaeva, 2006
00010 
00011 *****************************************************************************/
00012 
00019 #ifndef _ARAGELI_smith_hpp_
00020 #define _ARAGELI_smith_hpp_
00021 
00022 #include "config.hpp"
00023 
00024 #include "smith/classic.hpp"
00025 #include "smith/near_optimal.hpp"
00026 #include "smith/storjohann.hpp"
00027 
00028 
00029 namespace Arageli
00030 {
00031 
00032 
00033 namespace ctrl
00034 {
00035 
00036 struct smith_idler
00037 {
00038     void preamble ()
00039     {}
00040 
00041     void conclusion ()
00042     {}
00043 
00044     template <typename Q, typename B, typename P>
00045     void current_matrices (const Q&, const B&, const P&)
00046     {}
00047 
00048     template <typename I, typename J>
00049     void find_smallest_nonzero (const I&, const J&)
00050     {}
00051 
00052     template <typename I, typename J>
00053     void pivot_item (const I&, const J& j)
00054     {}
00055 
00056     void after_pivoting () {}
00057 
00058     template <typename I, typename J>
00059     void nondivisor_entry (const I& k, const I& i, const J& l)
00060     {}
00061 
00062     void pivot_adjustment () {}
00063 
00064     template <typename C, typename I, typename J, typename Q, typename B, typename P>
00065     bool stop (const C& corner, const I& i, const J& j, const Q& q, const B& b, const P& p)
00066     { return false; }
00067 };
00068 
00069 }   // namespace ctrl
00070 
00071 
00073 template <typename T>
00074 struct smith_norm : std::unary_function<T&, void>
00075 {
00076     template <typename P>
00077     void operator() (T& x, P& p, std::size_t i) const
00078     { do_norm(x, p, i, bool_type<type_traits<T>::is_polynom>::value); }
00079 
00080 private:
00081 
00083     template <typename P>
00084     void do_norm (T& x, P& p, std::size_t i, true_type) const
00085     {
00086         ARAGELI_ASSERT_0(!is_null(x));
00087         p.div_row(i, x.leading_coef());
00088         x /= safe_reference(x.leading_coef());
00089     }
00090 
00092     template <typename P>
00093     void do_norm (T& x, P& p, std::size_t i, false_type) const {}
00094 };
00095 
00097 
00104 template
00105 <
00106     typename MA,
00107     typename MB,
00108     typename MP,
00109     typename MQ,
00110     typename Rank,
00111     typename T_det, typename T_factory,
00112     typename Ctrler, typename Norm
00113 >
00114 void smith
00115 (
00116     const MA& A,
00117     MB& B,
00118     MP& P,
00119     MQ& Q,
00120     Rank& rank,
00121     T_det& det,
00122     const T_factory& tfctr,
00123     Ctrler ctrler,
00124     Norm norm
00125 );
00126 
00127 
00129 template
00130 <
00131     typename MA,
00132     typename MB,
00133     typename MP,
00134     typename MQ,
00135     typename Rank,
00136     typename T_det, typename Ctrler
00137 >
00138 inline void smith
00139 (
00140     const MA& A,
00141     MB& B,
00142     MP& P,
00143     MQ& Q,
00144     Rank& rank,
00145     T_det& det,
00146     Ctrler ctrler
00147 )
00148 {
00149     typedef typename MA::element_type T_A;
00150     
00151     return smith
00152     (
00153         A, B, P, Q, rank, det, factory<T_det>(),
00154         ctrler, smith_norm<T_A>()
00155     );
00156 }
00157 
00158 
00160 template
00161 <
00162     typename MA,
00163     typename MB,
00164     typename MP,
00165     typename MQ,
00166     typename Rank,
00167     typename T_det
00168 >
00169 inline void smith
00170 (
00171     const MA& A,
00172     MB& B,
00173     MP& P,
00174     MQ& Q,
00175     Rank& rank,
00176     T_det& det
00177 )
00178 {
00179     typedef typename MA::element_type T_A;
00180     
00181     return smith
00182     (
00183         A, B, P, Q, rank, det, factory<T_det>(),
00184         ctrl::smith_idler(), smith_norm<T_A>()
00185     );
00186 }
00187 
00188 
00190 template <typename M>
00191 M smith (const M& a);
00192 
00193 
00195 
00202 template
00203 <
00204     typename MA,
00205     typename MB,
00206     typename MP,
00207     typename MQ,
00208     typename Rank,
00209     typename T_det
00210 >
00211 void smith_storjohann
00212 (
00213     const MA &A,
00214     MB &B,
00215     MP &P,
00216     MQ &Q,
00217     Rank &rank,
00218     T_det &det
00219 );
00220 
00221 
00223 
00229 template
00230 <
00231     typename MA,
00232     typename MB,
00233     typename MP,
00234     typename MQ,
00235     typename Rank,
00236     typename T_det
00237 >
00238 void smith_near_optimal
00239 (
00240     const MA &A,
00241     MB &B,
00242     MP &P,
00243     MQ &Q,
00244     Rank &rank,
00245     T_det &det
00246 );
00247 
00248 
00249 } // namesapce Arageli
00250 
00251 
00252 #ifdef ARAGELI_INCLUDE_CPP_WITH_EXPORT_TEMPLATE
00253     #define ARAGELI_INCLUDE_CPP_WITH_EXPORT_TEMPLATE_smith
00254     #include "smith.cpp"
00255     #undef  ARAGELI_INCLUDE_CPP_WITH_EXPORT_TEMPLATE_smith
00256 #endif
00257 
00258 #endif  // #ifndef _ARAGELI_smith_hpp_