00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016 #ifndef __NUMERICALINTEGRATION_H__
00017 #define __NUMERICALINTEGRATION_H__
00018
00019
00020
00021 #include "../utilities/vector.h"
00022
00023
00024
00025 #include "../../parsers/tinyxml/tinyxml.h"
00026
00027
00028 namespace OpenNN
00029 {
00030
00033
00034 class NumericalIntegration
00035 {
00036 public:
00037
00038
00039
00040 explicit NumericalIntegration(void);
00041
00042
00043
00044
00045 virtual ~NumericalIntegration(void);
00046
00048
00049 enum NumericalIntegrationMethod{TrapezoidMethod, SimpsonMethod};
00050
00051
00052
00053 const NumericalIntegrationMethod& get_numerical_integration_method(void) const;
00054 std::string write_numerical_integration_method(void) const;
00055
00056 const bool& get_display(void) const;
00057
00058 void set(const NumericalIntegration&);
00059
00060 void set_numerical_integration_method(const NumericalIntegrationMethod&);
00061 void set_numerical_integration_method(const std::string&);
00062
00063 void set_display(const bool&);
00064
00065 void set_default(void);
00066
00067
00068
00069 double calculate_trapezoid_integral(const Vector<double>&, const Vector<double>&) const;
00070 double calculate_Simpson_integral(const Vector<double>&, const Vector<double>&) const;
00071
00072 double calculate_integral(const Vector<double>&, const Vector<double>&) const;
00073
00074
00075
00076 TiXmlElement* to_XML(void) const;
00077 void from_XML(TiXmlElement*);
00078
00079
00080
00081
00082
00083
00084
00092
00093 template<class T>
00094 double calculate_trapezoid_integral(const T& t,
00095 double (T::*f)(const double&) const , const double& a, const double& b, const unsigned int& n) const
00096 {
00097
00098
00099 const double h = (b-a)/(n-1.0);
00100
00101
00102
00103 double sum = (t.*f)(a)/2.0;
00104
00105 for(unsigned int i = 1; i < n-1; i++)
00106 {
00107 sum += (t.*f)(a + i*h);
00108 }
00109
00110 sum += (t.*f)(b)/2.0;
00111
00112
00113
00114 return(h*sum);
00115 }
00116
00117
00118
00119
00120
00128
00129 template<class T>
00130 Vector<double> calculate_trapezoid_integral(const T& t, Vector<double> (T::*f)(const double&) const , const double& a, const double& b, const unsigned int& n) const
00131 {
00132
00133
00134 const double h = (b-a)/(n-1.0);
00135
00136
00137
00138 Vector<double> sum = (t.*f)(a)/2.0;
00139
00140 for(unsigned int i = 1; i < n-1; i++)
00141 {
00142 sum += (t.*f)(a + i*h);
00143 }
00144
00145 sum += (t.*f)(b)/2.0;
00146
00147
00148
00149 return(sum*h);
00150 }
00151
00152
00153
00154
00155
00163
00164 template<class T>
00165 double calculate_Simpson_integral(const T& t, double (T::*f)(const double&) const , const double& a, const double& b, const unsigned int& n) const
00166 {
00167 const double h = (b-a)/(n-1.0);
00168
00169 double sum = (t.*f)(a)/3.0;
00170
00171 for(unsigned int i = 1; i < n-1; i++)
00172 {
00173 if(i%2 != 0)
00174 {
00175 sum += 4.0*(t.*f)(a + i*h)/3.0;
00176 }
00177 else
00178 {
00179 sum += 2.0*(t.*f)(a + i*h)/3.0;
00180 }
00181 }
00182
00183 sum += (t.*f)(b)/3.0;
00184
00185 return(h*sum);
00186 }
00187
00188
00189
00190
00191
00199
00200 template<class T>
00201 Vector<double> calculate_Simpson_integral(const T& t, Vector<double> (T::*f)(const double&) const , const double& a, const double& b, const unsigned int& n) const
00202 {
00203 const double h = (b-a)/(n-1.0);
00204
00205 Vector<double> sum = (t.*f)(a)/3.0;
00206
00207 for(unsigned int i = 1; i < n-1; i++)
00208 {
00209 if(i%2 != 0)
00210 {
00211 sum += (t.*f)(a + i*h)*4.0/3.0;
00212 }
00213 else
00214 {
00215 sum += (t.*f)(a + i*h)*2.0/3.0;
00216 }
00217 }
00218
00219 sum += (t.*f)(b)/3.0;
00220
00221
00222
00223 return(sum*h);
00224 }
00225
00226
00227
00228
00229
00236
00237 template<class T>
00238 double calculate_integral(const T& t, double (T::*f)(const double&) const , const double& a, const double& b, const unsigned int& n) const
00239 {
00240 switch(numerical_integration_method)
00241 {
00242 case TrapezoidMethod:
00243 {
00244 return(calculate_trapezoid_integral(t, &T::f, a, b, n));
00245 }
00246 break;
00247
00248 case SimpsonMethod:
00249 {
00250 return(calculate_Simpson_integral(t, &T::f, a, b, n));
00251 }
00252 break;
00253
00254 default:
00255 {
00256 std::ostringstream buffer;
00257
00258 buffer << "OpenNN Exception: NumericalIntegration class.\n"
00259 << "double calculate_integral(const T&, double (T::*f)(const double&) const, const double&, const double&, const unsigned int&) const method.\n"
00260 << "Unknown numerical integration method.\n";
00261
00262 throw std::logic_error(buffer.str());
00263 }
00264 break;
00265 }
00266 }
00267
00268
00269
00270
00271
00278
00279 template<class T>
00280 Vector<double> calculate_integral(const T& t, Vector<double> (T::*f)(const double&) const , const double& a, const double& b, const unsigned int& n) const
00281 {
00282 switch(numerical_integration_method)
00283 {
00284 case TrapezoidMethod:
00285 {
00286 return(calculate_trapezoid_integral(t, &T::f, a, b, n));
00287 }
00288 break;
00289
00290 case SimpsonMethod:
00291 {
00292 return(calculate_Simpson_integral(t, &T::f, a, b, n));
00293 }
00294 break;
00295
00296 default:
00297 {
00298 std::ostringstream buffer;
00299
00300 buffer << "OpenNN Exception: NumericalIntegration class.\n"
00301 << "double calculate_integral(const T&, Vector<double> (T::*f)(const double&) const, const double&, const double&, const unsigned int&) const method.\n"
00302 << "Unknown numerical integration method.\n";
00303
00304 throw std::logic_error(buffer.str());
00305 }
00306 break;
00307 }
00308 }
00309
00310
00311 private:
00312
00314
00315 NumericalIntegrationMethod numerical_integration_method;
00316
00318
00319 bool display;
00320
00321 };
00322
00323 }
00324
00325 #endif
00326
00327
00328
00329
00330
00331
00332
00333
00334
00335
00336
00337
00338
00339
00340
00341
00342
00343