Object-oriented Scientific Computing Library: Version 0.910
gsl_inte.h
00001 /*
00002   -------------------------------------------------------------------
00003   
00004   Copyright (C) 2006-2012, Andrew W. Steiner
00005   
00006   This file is part of O2scl.
00007   
00008   O2scl is free software; you can redistribute it and/or modify
00009   it under the terms of the GNU General Public License as published by
00010   the Free Software Foundation; either version 3 of the License, or
00011   (at your option) any later version.
00012   
00013   O2scl is distributed in the hope that it will be useful,
00014   but WITHOUT ANY WARRANTY; without even the implied warranty of
00015   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00016   GNU General Public License for more details.
00017   
00018   You should have received a copy of the GNU General Public License
00019   along with O2scl. If not, see <http://www.gnu.org/licenses/>.
00020 
00021   -------------------------------------------------------------------
00022 */
00023 /* 
00024  * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Brian Gough
00025  * 
00026  * This program is free software; you can redistribute it and/or modify
00027  * it under the terms of the GNU General Public License as published by
00028  * the Free Software Foundation; either version 3 of the License, or (at
00029  * your option) any later version.
00030  * 
00031  * This program is distributed in the hope that it will be useful, but
00032  * WITHOUT ANY WARRANTY; without even the implied warranty of
00033  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00034  * General Public License for more details.
00035  * 
00036  * You should have received a copy of the GNU General Public License
00037  * along with this program; if not, write to the Free Software
00038  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 
00039  * 02110-1301, USA.
00040  */
00041 #ifndef O2SCL_GSL_INTE_H
00042 #define O2SCL_GSL_INTE_H
00043 
00044 #include <gsl/gsl_machine.h>
00045 
00046 #ifndef DOXYGENP
00047 namespace o2scl {
00048 #endif
00049   
00050   /** \brief GSL integration base
00051       
00052       This base class does not perform any actual integration, but
00053       just provides functions to be used in the integration
00054       classes based on GSL.
00055   */
00056   class gsl_inte {
00057     
00058   public:
00059 
00060     gsl_inte() {
00061     }
00062     
00063   protected:
00064 
00065     /** \brief QUADPACK's nonlinear rescaling of the absolute-error 
00066         estimate.
00067           
00068         The values \f$ \rho_{\mathrm{abs}} \f$ (stored in
00069         <tt>result_abs</tt>) and \f$ \rho_{\mathrm{abs}} \f$ (stored
00070         in <tt>result_asc</tt>) are assumed to be
00071         \f{eqnarray*}{
00072         \rho_{\mathrm{abs}} &=& \int_a^b |f|\,dx, \\
00073         \rho_{\mathrm{asc}} &=& \int_a^b |f - \mu(f)|\, dx, \qquad
00074         \mu(f) = \frac{1}{b-a}\int_a^b f\, dx,
00075         \f}
00076         all of which are computed from the best (i.e., finest-grid)
00077         approximation of the integrals.  The rescaled error, \f$
00078         \sigma_\mathrm{err}, \f$ is computed from the raw error, \c
00079         err, by
00080         \f[
00081         \sigma_\mathrm{err} =
00082         \rho_\mathrm{asc} \cdot 
00083         \min \left\{1, \; 
00084         \left(\frac{200 |\mathrm{err}|}{\rho_\mathrm{asc}} \right)^{3/2}
00085         \right\},
00086         \f]
00087         or
00088         \f[
00089         \sigma_\mathrm{err} = 
00090         50\cdot \epsilon_\mathrm{mach} \cdot \rho_\mathrm{abs},
00091         \f]
00092         whichever of the two is greater. The value \f$
00093         \epsilon_\mathrm{mach} \f$ denotes "machine epsilon." (In the
00094         case that the second value underflows, the first value is
00095         automatically accepted.)
00096 
00097         This function is used in \ref gsl_inte_qng and \ref
00098         gsl_inte_kronrod::gauss_kronrod_base().
00099     */
00100     double rescale_error(double err, const double result_abs, 
00101                          const double result_asc) {
00102 
00103       err = fabs(err);
00104       
00105       if (result_asc != 0 && err != 0) {
00106         
00107         double scale = pow((200 * err / result_asc), 1.5);
00108         
00109         if (scale < 1) {
00110           err = result_asc * scale;
00111         } else {
00112           err = result_asc;
00113         }
00114       }
00115 
00116       if (result_abs > GSL_DBL_MIN / (50 * GSL_DBL_EPSILON)) {
00117 
00118         double min_err = 50 * GSL_DBL_EPSILON * result_abs;
00119         
00120         if (min_err > err) {
00121           err = min_err;
00122         }
00123       }
00124       
00125       return err;
00126     }
00127     
00128   };
00129   
00130   
00131 #ifndef DOXYGENP
00132 }
00133 #endif
00134 
00135 #endif
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).

Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads.