// @(#)root/mathmore:$Id: RootFinder.h 25486 2008-09-22 12:43:03Z moneta $
// Authors: L. Moneta, A. Zsenei 08/2005
/**********************************************************************
* *
* Copyright (c) 2004 ROOT Foundation, CERN/PH-SFT *
* *
* This library is free software; you can redistribute it and/or *
* modify it under the terms of the GNU General Public License *
* as published by the Free Software Foundation; either version 2 *
* of the License, or (at your option) any later version. *
* *
* This library is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
* General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this library (see file COPYING); if not, write *
* to the Free Software Foundation, Inc., 59 Temple Place, Suite *
* 330, Boston, MA 02111-1307 USA, or contact the author. *
* *
**********************************************************************/
// Header file for class RootFinder
//
// Created by: moneta at Sun Nov 14 16:59:55 2004
//
// Last update: Sun Nov 14 16:59:55 2004
//
#ifndef ROOT_Math_RootFinder
#define ROOT_Math_RootFinder
#ifndef ROOT_Math_IFunctionfwd
#include "Math/IFunctionfwd.h"
#endif
#ifndef ROOT_Math_IRootFinderMethod
#include "Math/IRootFinderMethod.h"
#endif
/**
@defgroup RootFinders One-dimensional Root-Finding algorithms
@ingroup NumAlgo
*/
namespace ROOT {
namespace Math {
//_____________________________________________________________________________________
/**
Class to find the Root of one dimensional functions.
The class is templated on the type of Root solver algorithms.
The possible types of Root-finding algorithms are:
- Root Bracketing Algorithms which they do not require function derivatives
- Roots::Bisection
- Roots::FalsePos
- Roots::Brent
- Root Finding Algorithms using Derivatives
- Roots::Newton
- Roots::Secant
- Roots::Steffenson
This class does not cupport copying
@ingroup RootFinders
*/
class RootFinder {
public:
enum EType { kBRENT, // Methods from MathCore
kGSL_BISECTION, kGSL_FALSE_POS, kGSL_BRENT, // GSL Normal
kGSL_NEWTON, kGSL_SECANT, kGSL_STEFFENSON // GSL Derivatives
};
/**
Construct a Root-Finder algorithm
*/
RootFinder(RootFinder::EType type = RootFinder::kBRENT);
virtual ~RootFinder();
private:
// usually copying is non trivial, so we make this unaccessible
RootFinder(const RootFinder & ) {}
RootFinder & operator = (const RootFinder & rhs)
{
if (this == &rhs) return *this; // time saving self-test
return *this;
}
public:
int SetMethod(RootFinder::EType type = RootFinder::kBRENT);
/**
Provide to the solver the function and the initial search interval [xlow, xup]
for algorithms not using derivatives (bracketing algorithms)
The templated function f must be of a type implementing the \a operator() method,
double operator() ( double x )
Returns non zero if interval is not valid (i.e. does not contains a root)
*/
int SetFunction( const IGenFunction & f, double xlow, double xup) {
return fSolver->SetFunction( f, xlow, xup);
}
/**
Provide to the solver the function and an initial estimate of the root,
for algorithms using derivatives.
The templated function f must be of a type implementing the \a operator()
and the \a Gradient() methods.
double operator() ( double x )
Returns non zero if starting point is not valid
*/
int SetFunction( const IGradFunction & f, double xstart) {
return fSolver->SetFunction( f, xstart);
}
template
int Solve(Function f, Derivative d, double start,
int maxIter = 100, double absTol = 1E-3, double relTol = 1E-6);
template
int Solve(Function f, double min, double max,
int maxIter = 100, double absTol = 1E-3, double relTol = 1E-6);
/**
Compute the roots iterating until the estimate of the Root is within the required tolerance returning
the iteration Status
*/
int Solve( int maxIter = 100, double absTol = 1E-3, double relTol = 1E-6) {
return fSolver->Solve( maxIter, absTol, relTol );
}
/**
Return the number of iteration performed to find the Root.
*/
int Iterations() const {
return fSolver->Iterations();
}
/**
Perform a single iteration and return the Status
*/
int Iterate() {
return fSolver->Iterate();
}
/**
Return the current and latest estimate of the Root
*/
double Root() const {
return fSolver->Root();
}
/**
Return the current and latest estimate of the lower value of the Root-finding interval (for bracketing algorithms)
*/
/* double XLower() const { */
/* return fSolver->XLower(); */
/* } */
/**
Return the current and latest estimate of the upper value of the Root-finding interval (for bracketing algorithms)
*/
/* double XUpper() const { */
/* return fSolver->XUpper(); */
/* } */
/**
Get Name of the Root-finding solver algorithm
*/
const char * Name() const {
return fSolver->Name();
}
#ifdef LATER
/**
Test convertgence Status of current iteration using interval values (for bracketing algorithms)
*/
static int TestInterval( double xlow, double xup, double epsAbs, double epsRel) {
return GSLRootHelper::TestInterval(xlow, xup, epsAbs, epsRel);
}
/**
Test convergence Status of current iteration using last Root estimates (for algorithms using function derivatives)
*/
static int TestDelta( double r1, double r0, double epsAbs, double epsRel) {
return GSLRootHelper::TestDelta(r1, r0, epsAbs, epsRel);
}
/**
Test function residual
*/
static int TestResidual(double f, double epsAbs) {
return GSLRootHelper::TestResidual(f, epsAbs);
}
#endif
protected:
private:
IRootFinderMethod* fSolver; // type of algorithm to be used
};
} // namespace Math
} // namespace ROOT
#include "Math/WrappedFunction.h"
#include "Math/Functor.h"
template
int ROOT::Math::RootFinder::Solve(Function f, Derivative d, double start,
int maxIter, double absTol, double relTol)
{
ROOT::Math::GradFunctor1D wf(f, d);
if (fSolver) fSolver->SetFunction(wf, start);
return Solve(maxIter, absTol, relTol);
}
template
int ROOT::Math::RootFinder::Solve(Function f, double min, double max,
int maxIter, double absTol, double relTol)
{
ROOT::Math::WrappedFunction wf(f);
if (fSolver) fSolver->SetFunction(wf, min, max);
return Solve(maxIter, absTol, relTol);
}
#endif /* ROOT_Math_RootFinder */