math_util.h

説明を見る。
00001 // -*-c++-*-
00002 
00008 /*
00009  *Copyright:
00010 
00011  Copyright (C) Hidehisa AKIYAMA
00012 
00013  This code is free software; you can redistribute it and/or
00014  modify it under the terms of the GNU Lesser General Public
00015  License as published by the Free Software Foundation; either
00016  version 2.1 of the License, or (at your option) any later version.
00017 
00018  This library is distributed in the hope that it will be useful,
00019  but WITHOUT ANY WARRANTY; without even the implied warranty of
00020  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00021  Lesser General Public License for more details.
00022 
00023  You should have received a copy of the GNU Lesser General Public
00024  License along with this library; if not, write to the Free Software
00025  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00026 
00027  *EndCopyright:
00028  */
00029 
00031 
00032 #ifndef RCSC_MATH_UTIL_H
00033 #define RCSC_MATH_UTIL_H
00034 
00035 #include <algorithm>
00036 #include <cmath>
00037 
00038 namespace rcsc {
00039 
00041 static const double EPS = 1.0e-10;
00042 
00043 /*-------------------------------------------------------------------*/
00051 template < typename  T >
00052 const T &
00053 bound( const T & low, const T & x, const T & high )
00054 {
00055     return std::min( std::max( low, x ), high );
00056 }
00057 
00058 /*-------------------------------------------------------------------*/
00066 template < typename  T >
00067 const T &
00068 min_max( const T & low, const T & x, const T & high )
00069 {
00070     return std::min( std::max( low, x ), high );
00071 }
00072 
00073 /*-------------------------------------------------------------------*/
00078 template < typename T >
00079 T
00080 square( const T & x )
00081 {
00082     return x * x;
00083 }
00084 
00085 /*-------------------------------------------------------------------*/
00092 inline
00093 double
00094 sign( const double & x )
00095 {
00096     return x > 0.0 ? 1.0 : -1.0;
00097 }
00098 
00099 /*-------------------------------------------------------------------*/
00109 inline
00110 int
00111 quadratic_formula( const double & a,
00112                    const double & b,
00113                    const double & c,
00114                    double * sol1,
00115                    double * sol2 )
00116 {
00117     double d = b * b - 4.0 * a * c;
00118 
00119     // ignore small noise
00120     if ( std::fabs( d ) < 0.001 )
00121     {
00122         if ( sol1 ) *sol1 = -b / ( 2.0 * a );
00123         return 1;
00124     }
00125 
00126     if ( d < 0.0 )
00127     {
00128         return 0;
00129     }
00130 
00131     d = std::sqrt( d );
00132     if ( sol1 ) *sol1 = (-b + d) / (2.0 * a);
00133     if ( sol2 ) *sol2 = (-b - d) / (2.0 * a);
00134     return 2;
00135 }
00136 
00137 /*-------------------------------------------------------------------*/
00148 inline
00149 double
00150 calc_sum_geom_series( const double & first_term,
00151                       const double & r,
00152                       const int len )
00153 {
00154     // sum     = f + fr + fr^2 + ... + fr^(n-1)
00155     // sum * r =     fr + fr^2 + ... + fr^(n-1) + fr^n
00156     // sum * ( r - 1 ) = fr^n - f
00157     // sum = f * ( r^n - 1.0 ) / ( r - 1 )
00158     return first_term * ( ( std::pow( r, len ) - 1.0 ) / ( r - 1.0 ) );
00159 }
00160 
00161 
00162 /*-------------------------------------------------------------------*/
00172 inline
00173 double
00174 calc_sum_inf_geom_series( const double & first_term,
00175                           const double & r )
00176 {
00177     if ( r < 0.0 || 1.0 <= r )
00178     {
00179         return 0.0;
00180     }
00181 
00182     // limit(n->inf, 0<r<1)  sum = f * ( 1 - r^n ) / ( 1 - r )
00183     return first_term / ( 1.0 - r );
00184 }
00185 
00186 /*-------------------------------------------------------------------*/
00196 inline
00197 double
00198 calc_first_term_geom_series( const double & sum,
00199                              const double & r,
00200                              const int len )
00201 {
00202     // sum = f * ( 1 - r^n ) / ( 1 - r )
00203     // f   = sum * ( 1 - r ) / ( 1 - r^n )
00204     return sum * ( 1.0 - r ) / ( 1.0 - std::pow( r, len ) );
00205 }
00206 
00207 /*-------------------------------------------------------------------*/
00218 inline
00219 double
00220 calc_first_term_inf_geom_series( const double & sum,
00221                                  const double & r )
00222 {
00223     // limit(n->inf, 0<r<1) f = sum * ( 1 - r ) / ( 1 - r^n )
00224     return sum * ( 1.0 - r );
00225 }
00226 
00227 /*-------------------------------------------------------------------*/
00239 inline
00240 double
00241 calc_first_term_geom_series_last( const double & last_term,
00242                                   const double & sum,
00243                                   const double & r )
00244 {
00245     if ( std::fabs( last_term ) < 0.001 )
00246     {
00247         return sum * ( 1.0 - r );
00248     }
00249 
00250     // l + (l * 1/r) + ... + (l * 1/r^(n-1))               = sum
00251     //     (l * 1/r) + ... + (l * 1/r^(n-1)) + (l * 1/r^n) = sum * (1/r)
00252     // l*(1/r^n) - l = sum * (1/r - 1)
00253     // (1/r^n) = sum * (1/r - 1) / l + 1
00254     double inverse = 1.0 / r;
00255     double tmp = 1.0 + sum * (inverse - 1.0) / last_term;
00256     if ( tmp < 0.001 )
00257     {
00258         return last_term;
00259     }
00260 
00261     //double len = std::log( tmp ) / std::log( inverse );
00262     //return last_term * std::pow( inverse, len );
00263     return last_term * std::pow( inverse, std::log( tmp ) / std::log( inverse ) );
00264 }
00265 
00266 /*-------------------------------------------------------------------*/
00276 inline
00277 double
00278 calc_length_geom_series( const double & first_term,
00279                          const double & sum,
00280                          const double & r )
00281 {
00282     if ( first_term <= EPS
00283          || sum < 0.0
00284          || r <= EPS )
00285     {
00286         // cannot take the zero first term
00287         // cannot take the negative sum
00288         // cannot take the negative ratio
00289         return -1.0;
00290     }
00291 
00292     if ( sum <= EPS )
00293     {
00294         // already there
00295         return 0.0;
00296     }
00297 
00298     // f + fr + fr^2 + ... + fr^(n-1)        = sum
00299     //     fr + fr^2 + ... + fr^(n-1) + fr^n = sum * r
00300     // fr^n - f = sum * ( r - 1 )
00301     // r^n = 1 + sum * ( r - 1 ) / f
00302 
00303     double tmp = 1.0 + sum * ( r - 1.0 ) / first_term;
00304     if ( tmp <= EPS )
00305     {
00306         return -1.0;
00307     }
00308     return std::log( tmp ) / std::log( r );
00309 }
00310 
00311 } // end namespace
00312 
00313 #endif

librcscに対してThu May 1 15:41:20 2008に生成されました。  doxygen 1.5.0