All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
inte_qaws_gsl.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2014, Jerry Gagelman
5  and Andrew W. Steiner
6 
7  This file is part of O2scl.
8 
9  O2scl is free software; you can redistribute it and/or modify
10  it under the terms of the GNU General Public License as published by
11  the Free Software Foundation; either version 3 of the License, or
12  (at your option) any later version.
13 
14  O2scl is distributed in the hope that it will be useful,
15  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  GNU General Public License for more details.
18 
19  You should have received a copy of the GNU General Public License
20  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
21 
22  -------------------------------------------------------------------
23 */
24 /*
25  * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Brian Gough
26  *
27  * This program is free software; you can redistribute it and/or modify
28  * it under the terms of the GNU General Public License as published by
29  * the Free Software Foundation; either version 3 of the License, or (at
30  * your option) any later version.
31  *
32  * This program is distributed in the hope that it will be useful, but
33  * WITHOUT ANY WARRANTY; without even the implied warranty of
34  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
35  * General Public License for more details.
36  *
37  * You should have received a copy of the GNU General Public License
38  * along with this program; if not, write to the Free Software
39  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
40  * 02110-1301, USA.
41  */
42 #ifndef GSL_INTE_QAWS_H
43 #define GSL_INTE_QAWS_H
44 
45 /** \file inte_qaws_gsl.h
46  \brief File defining \ref o2scl::inte_qaws_gsl
47 */
48 
49 #include <o2scl/inte_qawc_gsl.h>
50 
51 #ifndef DOXYGEN_NO_O2NS
52 namespace o2scl {
53 #endif
54 
55  /** \brief Adaptive integration with with algebraic-logarithmic
56  singularities at the end-points (GSL)
57 
58  This class computes the weighted integral
59  \f[
60  \int_a^b f(x)(x - a)^\alpha (b - x)^\beta \log^\mu(x - a)
61  \log^\nu(b - x)~dx
62  \f]
63  where the parameters of the weight function must satisfy
64  \f[
65  \alpha > -1, \quad \beta > -1, \quad
66  \mu \in \{0, 1\}, \quad \nu \in \{0, 1\},
67  \f]
68  and which are set by \ref set_weight(). Note that setting \f$
69  \mu=0 \f$ or \f$ \nu=0 \f$ removes the respective factor \f$
70  \log^mu(\ldots) \f$ or \f$ \log^\nu(\ldots) \f$ from the weight.
71 
72  The adaptive refinement algorithm described for \ref
73  inte_qag_gsl is used. When a subinterval contains one of the
74  endpoints, a special 25-point modified Clenshaw-Curtis rule is
75  used to control the singularities. For subintervals which do not
76  include the endpoints, a Gauss-Kronrod integration rule is used.
77 
78  See \ref gslinte_subsect in the User's guide for general
79  information about the GSL integration classes.
80  */
81  template<class func_t=funct11> class inte_qaws_gsl :
82  public inte_cheb_gsl<func_t> {
83 
84 #ifndef DOXYGEN_INTERNAL
85 
86  protected:
87 
88  /** \name Data from \c gsl_integration_qaws_table
89  */
90  //@{
91  double alpha;
92  double beta;
93  int mu;
94  int nu;
95 
96  double ri[25];
97  double rj[25];
98  double rg[25];
99  double rh[25];
100  //@}
101 
102  /** \brief Set the array values \c ri, \c rj, \c rg, \c rh from the
103  current values \c alpha and \c beta.
104 
105  This is the function from the GSL source code \c integration/qmomo.c
106  that initializes \c gsl_integration_qaws_table.
107  */
109 
110  const double alpha_p1 = this->alpha + 1.0;
111  const double beta_p1 = this->beta + 1.0;
112 
113  const double alpha_p2 = this->alpha + 2.0;
114  const double beta_p2 = this->beta + 2.0;
115 
116  const double r_alpha = pow (2.0, alpha_p1);
117  const double r_beta = pow (2.0,beta_p1);
118 
119  size_t i;
120 
121  double an, anm1;
122 
123  this->ri[0] = r_alpha / alpha_p1;
124  this->ri[1] = this->ri[0] * this->alpha / alpha_p2;
125 
126  an = 2.0;
127  anm1 = 1.0;
128 
129  for (i = 2; i < 25; i++) {
130  this->ri[i] = -(r_alpha + an * (an - alpha_p2) * this->ri[i - 1])
131  / (anm1 * (an + alpha_p1));
132  anm1 = an;
133  an = an + 1.0;
134  }
135 
136  rj[0] = r_beta / beta_p1;
137  rj[1] = rj[0] * this->beta / beta_p2;
138 
139  an = 2.0;
140  anm1 = 1.0;
141 
142  for (i = 2; i < 25; i++) {
143  rj[i] = (-(r_beta + an * (an - beta_p2) * rj[i - 1])
144  / (anm1 * (an + beta_p1)));
145  anm1 = an;
146  an = an + 1.0;
147  }
148 
149  this->rg[0] = -this->ri[0] / alpha_p1;
150  this->rg[1] = -this->rg[0] - 2.0 * r_alpha / (alpha_p2 * alpha_p2);
151 
152  an = 2.0;
153  anm1 = 1.0;
154 
155  for (i = 2; i < 25; i++) {
156  this->rg[i] = (-(an * (an - alpha_p2) *
157  this->rg[i - 1] - an * this->ri[i - 1]
158  + anm1 * this->ri[i])
159  / (anm1 * (an + alpha_p1)));
160  anm1 = an;
161  an = an + 1.0;
162  }
163 
164  this->rh[0] = -this->rj[0] / beta_p1;
165  this->rh[1] = -this->rh[0] - 2.0 * r_beta / (beta_p2 * beta_p2);
166 
167  an = 2.0;
168  anm1 = 1.0;
169 
170  for (i = 2; i < 25; i++) {
171  this->rh[i] = (-(an * (an - beta_p2) *
172  this->rh[i - 1] - an * this->rj[i - 1]
173  + anm1 * this->rj[i])
174  / (anm1 * (an + beta_p1)));
175  anm1 = an;
176  an = an + 1.0;
177  }
178 
179  for (i = 1; i < 25; i += 2) {
180  this->rj[i] *= -1;
181  this->rh[i] *= -1;
182  }
183 
184  return;
185  }
186 
187  /** \brief True if algebraic-logarithmic singularity is present at the
188  right endpoint in the definition \c f_trans.
189  */
190  bool fn_qaws_R;
191 
192  /** \brief True if algebraic-logarithmic singularity is present at the
193  left endpoint in the definition \c f_trans.
194  */
195  bool fn_qaws_L;
196 
197  /// Left endpoint in definition of \c f_trans
199 
200  /// Right endpoint in definition of \c f_trans.
202 
203  /** \brief Weighted integrand.
204  */
205  virtual double transform(double t, func_t &func) {
206 
207  double factor = 1.0,y;
208 
209  if (fn_qaws_L) {
210  if (alpha != 0.0) {
211  factor *= pow(t - left_endpoint,alpha);
212  }
213  if (mu == 1) {
214  factor *= log(t - left_endpoint);
215  }
216  }
217 
218  if (fn_qaws_R) {
219  if (beta != 0.0) {
220  factor *= pow(right_endpoint - t,beta);
221  }
222  if (nu == 1) {
223  factor *= log(right_endpoint - t);
224  }
225  }
226 
227  return func(t)*factor;
228  }
229 
230  /** \brief Clenshaw-Curtis 25-point integration and error estimator
231  for functions with an algebraic-logarithmic singularity at the
232  endpoint(s).
233  */
234  void qc25s(func_t &func, double a, double b, double a1, double b1,
235  double &result, double &abserr, int &err_reliable) {
236 
237  // Transformed function object for inte_cheb_series()
238  funct11 fmp=
239  std::bind(std::mem_fn<double(double,func_t &)>
241  this,std::placeholders::_1,func);
242 
243  this->left_endpoint = a;
244  this->right_endpoint = b;
245 
246  if (a1 == a && (this->alpha != 0.0 || this->mu != 0)) {
247 
248  double cheb12[13], cheb24[25];
249 
250  double factor = pow(0.5 * (b1 - a1),this->alpha + 1.0);
251 
252  // weighted_function.function = &fn_qaws_R;
253  this->fn_qaws_R = true;
254  this->fn_qaws_L = false;
255 
256  this->inte_cheb_series(fmp,a1,b1,cheb12,cheb24);
257 
258  if (this->mu == 0) {
259  double res12 = 0,res24 = 0;
260  double u = factor;
261 
262  this->compute_result(this->ri,cheb12,cheb24,res12,res24);
263 
264  result = u * res24;
265  abserr = fabs(u * (res24 - res12));
266 
267  } else {
268 
269  double res12a = 0,res24a = 0;
270  double res12b = 0,res24b = 0;
271 
272  double u = factor * log(b1 - a1);
273  double v = factor;
274 
275  this->compute_result(this->ri,cheb12,cheb24,res12a,res24a);
276  this->compute_result(this->rg,cheb12,cheb24,res12b,res24b);
277 
278  result = u * res24a + v * res24b;
279  abserr = fabs(u*(res24a - res12a)) + fabs(v*(res24b - res12b));
280  }
281 
282  err_reliable = 0;
283  return;
284 
285  } else if (b1 == b && (this->beta != 0.0 || this->nu != 0)) {
286 
287  double cheb12[13], cheb24[25];
288  double factor = pow(0.5 * (b1 - a1), this->beta + 1.0);
289 
290  // weighted_function.function = &fn_qaws_L;
291  this->fn_qaws_L = true;
292  this->fn_qaws_R = false;
293 
294  this->inte_cheb_series(fmp,a1,b1,cheb12,cheb24);
295 
296  if (this->nu == 0) {
297 
298  double res12 = 0, res24 = 0;
299  double u = factor;
300 
301  this->compute_result(this->rj, cheb12,cheb24,res12,res24);
302 
303  result = u * res24;
304  abserr = fabs(u * (res24 - res12));
305 
306  } else {
307 
308  double res12a = 0, res24a = 0;
309  double res12b = 0, res24b = 0;
310 
311  double u = factor * log(b1 - a1);
312  double v = factor;
313 
314  this->compute_result(this->rj,cheb12,cheb24,res12a,res24a);
315  this->compute_result(this->rh,cheb12,cheb24,res12b,res24b);
316 
317  result = u * res24a + v * res24b;
318  abserr = fabs(u*(res24a - res12a)) + fabs(v*(res24b - res12b));
319  }
320 
321  err_reliable = 0;
322  return;
323 
324  } else {
325 
326  double resabs, resasc;
327 
328  // weighted_function.function = &fn_qaws;
329  this->fn_qaws_R = true;
330  this->fn_qaws_L = true;
331 
332  this->gauss_kronrod(func,a1,b1,&result,&abserr,&resabs,&resasc);
333 
334  if (abserr == resasc) {
335  err_reliable = 0;
336  } else {
337  err_reliable = 1;
338  }
339 
340  return;
341  }
342  }
343 
344  /** \brief Compute the 13-point and 25-point approximations from
345  the Chebyshev moments and coefficients.
346  */
347  void compute_result(double *r, double *cheb12, double *cheb24,
348  double &result12, double &result24) {
349 
350  result12=0.0;
351  result24=0.0;
352 
353  size_t i;
354  for (i = 0; i < 13; i++) {
355  result12 += r[i] * cheb12[i];
356  }
357 
358  for (i = 0; i < 25; i++) {
359  result24 += r[i] * cheb24[i];
360  }
361  }
362 
363 #endif
364 
365  public:
366 
367  /** \brief Initialize the adptive workspace as with the constructor
368  \ref inte_qag_gsl::inte_qag_gsl.
369 
370  The default paramters \f$ \alpha, \beta, \mu, \nu \f$ of the weight
371  function are all zero.
372  */
373  inte_qaws_gsl() : inte_cheb_gsl<func_t>() {
374  set_weight(0.0,0.0,0,0);
375  }
376 
377  ~inte_qaws_gsl() {}
378 
379  /** \brief Sets the exponents of singularites of the weight function.
380 
381  The parameters determine the exponents of the weight function
382  \f[
383  W(x) = (x-a)^\alpha (b-x)^\beta \log^\mu(x-a) \log^\nu(b-x),
384  \f]
385  and must satsify
386  \f[
387  \alpha > -1, \quad \beta > -1, \quad
388  \mu \in \{0, 1\}, \quad \nu \in \{0, 1\}.
389  \f]
390  In order for the adaptive algorithm to run quickly, a table of
391  Chebyshev weights for the particular parameters are computed in
392  advance.
393  */
394  int set_weight(double u_alpha, double u_beta, int u_mu, int u_nu) {
395 
396  if (u_alpha < -1.0) {
397  std::string estr=((std::string)"Variable alpha must be ")+
398  "greater than -1.0 in inte_qaws_gsl().";
399  O2SCL_ERR(estr.c_str(),exc_einval);
400  }
401  if (u_beta < -1.0) {
402  std::string estr=((std::string)"Variable beta must be ")+
403  "greater than -1.0 in inte_qaws_gsl().";
404  O2SCL_ERR(estr.c_str(),exc_einval);
405  }
406  if (u_mu != 0 && u_mu != 1) {
407  std::string estr=((std::string)"Variable mu must be 0 or 1 ")+
408  "in inte_qaws_gsl().";
409  O2SCL_ERR(estr.c_str(),exc_einval);
410  }
411  if (u_nu != 0 && u_nu != 1) {
412  std::string estr=((std::string)"Variable nu must be 0 or 1 ")+
413  "in inte_qaws_gsl().";
414  O2SCL_ERR(estr.c_str(),exc_einval);
415  }
416 
417  this->alpha = u_alpha;
418  this->beta = u_beta;
419  this->mu = u_mu;
420  this->nu = u_nu;
421 
423  return success;
424  }
425 
426  /** \brief Returns the current values (via reference) of the
427  weight-function's parameters.
428  */
429  void get_weight(double &u_alpha, double &u_beta, int &u_mu, int &u_nu) {
430  u_alpha = this->alpha;
431  u_beta = this->beta;
432  u_mu = this->mu;
433  u_nu = this->nu;
434  }
435 
436  /** \brief Integrate the function \c func on the interval (\c a, \c b)
437  returning the \c result and error estimate \c abserr.
438  */
439  virtual int integ_err(func_t &func, double a, double b,
440  double &result, double &abserr) {
441 
442  double area, errsum;
443  double result0, abserr0;
444  double tolerance;
445  this->last_iter = 0;
446  int roundoff_type1 = 0, roundoff_type2 = 0, error_type = 0;
447 
448  /* Initialize results */
449 
450  this->w->initialise(a,b);
451 
452  result = 0;
453  abserr = 0;
454 
455  size_t limit=this->w->limit;
456 
457  if (b <= a) {
458  std::string estr="Integration limits, a="+dtos(a);
459  estr+=" and b="+dtos(b)+", must satisfy a < b";
460  estr+=" in inte_qaws_gsl::gsl_qaws().";
461  O2SCL_ERR(estr.c_str(),exc_einval);
462  }
463 
464  double dbl_eps=std::numeric_limits<double>::epsilon();
465 
466  if (this->tol_abs <= 0 && (this->tol_rel < 50 * dbl_eps ||
467  this->tol_rel < 0.5e-28)) {
468  this->last_iter=0;
469  std::string estr="Tolerance cannot be achieved with given ";
470  estr+="value of tol_abs, "+dtos(this->tol_abs)+", and tol_rel, "+
471  dtos(this->tol_rel)+", in inte_qaws_gsl::integ_err().";
472  O2SCL_ERR(estr.c_str(),exc_ebadtol);
473  }
474 
475  /* perform the first integration */
476 
477  {
478  double area1, area2;
479  double error1, error2;
480  int err_reliable1, err_reliable2;
481  double a1 = a;
482  double b1 = 0.5 * (a + b);
483  double a2 = b1;
484  double b2 = b;
485 
486  this->qc25s(func, a, b, a1, b1, area1, error1, err_reliable1);
487  this->qc25s(func, a, b, a2, b2, area2, error2, err_reliable2);
488 
489  this->last_iter = 2;
490 
491  if (error1 > error2) {
492  this->w->append_interval(a1, b1, area1, error1);
493  this->w->append_interval(a2, b2, area2, error2);
494  } else {
495  this->w->append_interval(a2, b2, area2, error2);
496  this->w->append_interval(a1, b1, area1, error1);
497  }
498 
499  result0 = area1 + area2;
500  abserr0 = error1 + error2;
501  }
502 
503  /* Test on accuracy */
504 
505  tolerance = GSL_MAX_DBL (this->tol_abs, this->tol_rel * fabs (result0));
506 
507  /* Test on accuracy, use 0.01 relative error as an extra safety
508  margin on the first iteration (ignored for subsequent iterations) */
509 
510  if (abserr0 < tolerance && abserr0 < 0.01 * fabs(result0)) {
511  result = result0;
512  abserr = abserr0;
513  return success;
514  } else if (limit == 1) {
515  result = result0;
516  abserr = abserr0;
517 
518  std::string estr = "A maximum of 1 iteration was insufficient ";
519  estr += "in inte_qaws_gsl::gsl_qaws().";
520  O2SCL_CONV_RET(estr.c_str(), exc_emaxiter, this->err_nonconv);
521  }
522 
523  area = result0;
524  errsum = abserr0;
525 
526  do {
527  double a1, b1, a2, b2;
528  double a_i, b_i, r_i, e_i;
529  double area1 = 0, area2 = 0, area12 = 0;
530  double error1 = 0, error2 = 0, error12 = 0;
531  int err_reliable1, err_reliable2;
532 
533  /* Bisect the subinterval with the largest error estimate */
534  this->w->retrieve(&a_i,&b_i,&r_i,&e_i);
535 
536  a1 = a_i;
537  b1 = 0.5 * (a_i + b_i);
538  a2 = b1;
539  b2 = b_i;
540 
541  qc25s(func, a, b, a1, b1, area1, error1, err_reliable1);
542  qc25s(func, a, b, a2, b2, area2, error2, err_reliable2);
543 
544  area12 = area1 + area2;
545  error12 = error1 + error2;
546 
547  errsum += (error12 - e_i);
548  area += area12 - r_i;
549 
550  if (err_reliable1 && err_reliable2) {
551 
552  double delta = r_i - area12;
553 
554  if (fabs(delta) <= 1.0e-5 * fabs (area12)
555  && error12 >= 0.99 * e_i) {
556  roundoff_type1++;
557  }
558  if (this->last_iter >= 10 && error12 > e_i) {
559  roundoff_type2++;
560  }
561  }
562 
563  tolerance = GSL_MAX_DBL (this->tol_abs, this->tol_rel * fabs (area));
564 
565  if (errsum > tolerance) {
566  if (roundoff_type1 >= 6 || roundoff_type2 >= 20) {
567  // round off error
568  error_type = 2;
569  }
570 
571  /* set error flag in the case of bad integrand behaviour at
572  a point of the integration range */
573  if (this->w->subinterval_too_small (a1, a2, b2)) {
574  error_type = 3;
575  }
576  }
577 
578  this->w->update(a1, b1, area1, error1, a2, b2, area2, error2);
579  this->w->retrieve(&a_i,&b_i,&r_i,&e_i);
580  this->last_iter++;
581 
582  } while (this->last_iter < this->w->limit
583  && !error_type && errsum > tolerance);
584 
585  result = this->w->sum_results();
586  abserr = errsum;
587  this->interror = abserr;
588 
589  if (errsum <= tolerance) {
590  return success;
591  } else if (error_type == 2) {
592  std::string estr="Round-off error prevents tolerance ";
593  estr+="from being achieved in inte_qaws_gsl::gsl_qaws().";
594  O2SCL_CONV_RET(estr.c_str(),exc_eround,this->err_nonconv);
595  } else if (error_type == 3) {
596  std::string estr="Bad integrand behavior ";
597  estr+=" in inte_qaws_gsl::gsl_qaws().";
598  O2SCL_CONV_RET(estr.c_str(),exc_esing,this->err_nonconv);
599  } else if (this->last_iter == limit) {
600  std::string estr="Maximum number of subdivisions ("+itos(limit);
601  estr+=") reached in inte_qaws_gsl::gsl_qaws().";
602  O2SCL_CONV_RET(estr.c_str(),exc_emaxiter,this->err_nonconv);
603  } else {
604  std::string estr="Could not integrate function in ";
605  estr+="inte_qaws_gsl::gsl_qaws().";
606  O2SCL_ERR(estr.c_str(),exc_efailed);
607  }
608 
609  // No return statement needed since the above if statement
610  // always forces a return
611  }
612 
613  /// Return string denoting type ("inte_qaws_gsl")
614  const char *type() { return "inte_qaws_gsl"; }
615 
616  };
617 
618 #ifndef DOXYGEN_NO_O2NS
619 }
620 #endif
621 
622 #endif
virtual int integ_err(func_t &func, double a, double b, double &result, double &abserr)
Integrate the function func on the interval (a, b) returning the result and error estimate abserr...
bool fn_qaws_L
True if algebraic-logarithmic singularity is present at the left endpoint in the definition f_trans...
bool fn_qaws_R
True if algebraic-logarithmic singularity is present at the right endpoint in the definition f_trans...
void get_weight(double &u_alpha, double &u_beta, int &u_mu, int &u_nu)
Returns the current values (via reference) of the weight-function's parameters.
#define O2SCL_CONV_RET(d, n, b)
Set a "convergence" error and return the error value.
Definition: err_hnd.h:292
std::function< double(double)> funct11
One-dimensional function typedef.
Definition: funct.h:44
double left_endpoint
Left endpoint in definition of f_trans.
void inte_cheb_series(func2_t &f, double a, double b, double *cheb12, double *cheb24)
Compute Chebyshev series expansion using a FFT method.
invalid argument supplied by user
Definition: err_hnd.h:59
apparent singularity detected
Definition: err_hnd.h:93
exceeded max number of iterations
Definition: err_hnd.h:73
Adaptive integration with with algebraic-logarithmic singularities at the end-points (GSL) ...
Definition: inte_qaws_gsl.h:81
size_t last_iter
The most recent number of iterations taken.
Definition: inte.h:63
bool err_nonconv
If true, call the error handler if the routine does not converge or reach the desired tolerance (defa...
Definition: inte.h:81
double right_endpoint
Right endpoint in definition of f_trans.
int update(double a1, double b1, double area1, double error1, double a2, double b2, double area2, double error2)
Determine which new subinterval to add to the workspace stack and perform update. ...
generic failure
Definition: err_hnd.h:61
Chebyshev integration base class (GSL)
Definition: inte_qawc_gsl.h:46
int subinterval_too_small(double a1, double a2, double b2)
Test whether the proposed subdivision falls before floating-point precision.
double tol_abs
The maximum absolute uncertainty in the value of the integral (default )
Definition: inte.h:73
virtual void gauss_kronrod(func_t &func, double a, double b, double *result, double *abserr, double *resabs, double *resasc)
Integration wrapper for internal transformed function type.
int set_weight(double u_alpha, double u_beta, int u_mu, int u_nu)
Sets the exponents of singularites of the weight function.
int initialise(double a, double b)
Initialize the workspace for an integration with limits a and b.
std::string dtos(double x, int prec=6, bool auto_prec=false)
Convert a double to a string.
inte_qaws_gsl()
Initialize the adptive workspace as with the constructor inte_qag_gsl::inte_qag_gsl.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
Definition: err_hnd.h:273
virtual double transform(double t, func_t &func)
Weighted integrand.
void compute_result(double *r, double *cheb12, double *cheb24, double &result12, double &result24)
Compute the 13-point and 25-point approximations from the Chebyshev moments and coefficients.
double interror
The uncertainty for the last integration computation.
Definition: inte.h:113
int retrieve(double *a, double *b, double *r, double *e) const
Retrieve the ith result from the workspace stack.
Integrate a function with a singularity (GSL) [abstract base].
void append_interval(double a1, double b1, double area1, double error1)
Push a new interval to the workspace stack.
inte_workspace_gsl * w
The integration workspace.
double sum_results()
Add up all of the contributions to construct the final result.
void qc25s(func_t &func, double a, double b, double a1, double b1, double &result, double &abserr, int &err_reliable)
Clenshaw-Curtis 25-point integration and error estimator for functions with an algebraic-logarithmic ...
double tol_rel
The maximum relative uncertainty in the value of the integral (default )
Definition: inte.h:68
size_t limit
Maximum number of subintervals allocated.
const char * type()
Return string denoting type ("inte_qaws_gsl")
std::string itos(int x)
Convert an integer to a string.
user specified an invalid tolerance
Definition: err_hnd.h:77
Success.
Definition: err_hnd.h:47
failed because of roundoff error
Definition: err_hnd.h:87
void initialise_qaws_table()
Set the array values ri, rj, rg, rh from the current values alpha and beta.

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).
Hosted at Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads..