13 #ifndef HEX_GAUSSKRONROD
14 #define HEX_GAUSSKRONROD
18 #include <gsl/gsl_integration.h>
61 gsl_integration_workspace * Workspace;
70 : Integrand(f), Result(special::constant::
Nan), AbsErr(special::constant::
Nan),
71 Ok(false), EpsAbs(0.), EpsRel(1
e-5), Limit(limit),
72 Workspace(gsl_integration_workspace_alloc(Limit)) {}
77 gsl_integration_workspace_free(Workspace);
92 static double eval (
double x,
void * ptr_radf)
137 if (std::isnan(a) or std::isnan(b))
141 Status =
"Some of the bounds is not finite.";
159 int err = GSL_SUCCESS, err1 = GSL_SUCCESS, err2 = GSL_SUCCESS;
160 double Result1, Result2, AbsErr1, AbsErr2;
163 if (std::isfinite(a) and std::isfinite(b))
165 err = gsl_integration_qag
168 EpsAbs, EpsRel, Limit,
173 Status = (err != GSL_SUCCESS ? gsl_strerror(err) :
"");
175 else if (std::isfinite(a) and not std::isfinite(b))
177 err = gsl_integration_qagiu
180 EpsAbs, EpsRel, Limit,
184 Status = (err != GSL_SUCCESS ? gsl_strerror(err) :
"");
186 else if (not std::isfinite(a) and std::isfinite(b))
188 err = gsl_integration_qagil
191 EpsAbs, EpsRel, Limit,
195 Status = (err != GSL_SUCCESS ? gsl_strerror(err) :
"");
199 err1 = gsl_integration_qagiu
202 EpsAbs, EpsRel, Limit,
207 err2 = gsl_integration_qagil
210 EpsAbs, EpsRel, Limit,
215 Result = Result1 + Result2;
216 AbsErr = AbsErr1 + AbsErr2;
218 Status = (err1 != GSL_SUCCESS or err2 != GSL_SUCCESS ? std::string(gsl_strerror(err1)) + std::string(
" & ") + std::string(gsl_strerror(err2)) : std::string());
221 Ok = (err == GSL_SUCCESS and err1 == GSL_SUCCESS and err2 == GSL_SUCCESS);
229 bool ok()
const {
return Ok; }
231 double error()
const {
return AbsErr; }
232 std::string
const &
status()
const {
return Status; }
237 double epsabs ()
const {
return EpsAbs; }
238 double epsrel ()
const {
return EpsRel; }
~GaussKronrod()
Definition: gausskronrod.h:75
double epsrel() const
Definition: gausskronrod.h:238
static double eval(double x, void *ptr_radf)
Evaluates radial function supplied in a pointer.
Definition: gausskronrod.h:92
double epsabs() const
Definition: gausskronrod.h:237
double F(double k, int l, double r, double sigma)
Definition: hydrogen.cpp:163
double operator()(double x) const
Evaluate the function.
Definition: gausskronrod.h:81
const double e
Definition: special.h:44
double result() const
Definition: gausskronrod.h:230
double error() const
Definition: gausskronrod.h:231
bool ok() const
Definition: gausskronrod.h:229
bool integrate(double a, double b)
Compute the integral.
Definition: gausskronrod.h:121
Numerical integrator.
Definition: gausskronrod.h:35
GaussKronrod(Functor f, size_t limit=1000)
Definition: gausskronrod.h:69
std::string const & status() const
Definition: gausskronrod.h:232
void setEpsAbs(double epsabs)
Definition: gausskronrod.h:234
Base class for radial functions.
Definition: special.h:536
const double Nan
Definition: special.h:68
void setEpsRel(double epsrel)
Definition: gausskronrod.h:235