1INTEG(3) User Contributed Perl Documentation INTEG(3)
2
6 PDL::GSL::INTEG - PDL interface to numerical integration routines in
7 GSL
8
10 This is an interface to the numerical integration package present in
11 the GNU Scientific Library, which is an implementation of QUADPACK.
12 Functions are named gslinteg_{algorithm} where {algorithm} is the
14 QUADPACK naming convention. The available functions are:
15 gslinteg_qng: Non-adaptive Gauss-Kronrod integration
17 gslinteg_qag: Adaptive integration
18 gslinteg_qags: Adaptive integration with singularities
19 gslinteg_qagp: Adaptive integration with known singular points
20 gslinteg_qagi: Adaptive integration on infinite interval of the form
21 (-\infty,\infty)
22 gslinteg_qagiu: Adaptive integration on infinite interval of the form
23 (la,\infty)
24 gslinteg_qagil: Adaptive integration on infinite interval of the form
25 (-\infty,lb)
26 gslinteg_qawc: Adaptive integration for Cauchy principal values
27 gslinteg_qaws: Adaptive integration for singular functions
28 gslinteg_qawo: Adaptive integration for oscillatory functions
29 gslinteg_qawf: Adaptive integration for Fourier integrals
30 Each algorithm computes an approximation to the integral, I, of the
32 function f(x)w(x), where w(x) is a weight function (for general
33 integrands w(x)=1). The user provides absolute and relative error
34 bounds (epsabs,epsrel) which specify the following accuracy
35 requirement:
36 |RESULT - I| <= max(epsabs, epsrel |I|)
38 The routines will fail to converge if the error bounds are too
40 stringent, but always return the best approximation obtained up to that
41 stage
42 All functions return the result, and estimate of the absolute error and
44 an error flag (which is zero if there were no problems). You are
45 responsible for checking for any errors, no warnings are issued unless
46 the option {Warn => 'y'} is specified in which case the reason of
47 failure will be printed.
48 You can nest integrals up to 20 levels. If you find yourself in the
50 unlikely situation that you need more, you can change the value of
51 'max_nested_integrals' in the first line of the file 'FUNC.c' and
52 recompile.
53
55 Throughout this documentation we strive to use the same variables that
56 are present in the original GSL documentation (see See Also).
57 Oftentimes those variables are called "a" and "b". Since good Perl
58 coding practices discourage the use of Perl variables $a and $b, here
59 we refer to Parameters "a" and "b" as $pa and $pb, respectively, and
60 Limits (of domain or integration) as $la and $lb.
61 Please check the GSL documentation for more information.
63
65 use PDL;
66 use PDL::GSL::INTEG;
67 my $la = 1.2;
69 my $lb = 3.7;
70 my $epsrel = 0;
71 my $epsabs = 1e-6;
72 # Non adaptive integration
74 my ($res,$abserr,$ierr,$neval) = gslinteg_qng(\&myf,$la,$lb,$epsrel,$epsabs);
75 # Warnings on
76 my ($res,$abserr,$ierr,$neval) = gslinteg_qng(\&myf,$la,$lb,$epsrel,$epsabs,{Warn=>'y'});
77 # Adaptive integration with warnings on
79 my $limit = 1000;
80 my $key = 5;
81 my ($res,$abserr,$ierr) = gslinteg_qag(\&myf,$la,$lb,$epsrel,
82 $epsabs,$limit,$key,{Warn=>'y'});
83 sub myf{
85 my ($x) = @_;
86 return exp(-$x**2);
87 }
88
90 qng_meat
91 Signature: (double a(); double b(); double epsabs();
92 double epsrel(); double [o] result(); double [o] abserr();
93 int [o] neval(); int [o] ierr(); int gslwarn(); SV* function)
94 info not available
96 qng_meat does not process bad values. It will set the bad-value flag
98 of all output ndarrays if the flag is set for any of the input
99 ndarrays.
100 qag_meat
102 Signature: (double a(); double b(); double epsabs();double epsrel(); int limit();
103 int key(); double [o] result(); double [o] abserr();int n();int [o] ierr();int gslwarn();; SV* function)
104 info not available
106 qag_meat does not process bad values. It will set the bad-value flag
108 of all output ndarrays if the flag is set for any of the input
109 ndarrays.
110 qags_meat
112 Signature: (double a(); double b(); double epsabs();double epsrel(); int limit();
113 double [o] result(); double [o] abserr();int n();int [o] ierr();int gslwarn();; SV* function)
114 info not available
116 qags_meat does not process bad values. It will set the bad-value flag
118 of all output ndarrays if the flag is set for any of the input
119 ndarrays.
120 qagp_meat
122 Signature: (double pts(l); double epsabs();double epsrel();int limit();
123 double [o] result(); double [o] abserr();int n();int [o] ierr();int gslwarn();; SV* function)
124 info not available
126 qagp_meat does not process bad values. It will set the bad-value flag
128 of all output ndarrays if the flag is set for any of the input
129 ndarrays.
130 qagi_meat
132 Signature: (double epsabs();double epsrel(); int limit();
133 double [o] result(); double [o] abserr(); int n(); int [o] ierr();int gslwarn();; SV* function)
134 info not available
136 qagi_meat does not process bad values. It will set the bad-value flag
138 of all output ndarrays if the flag is set for any of the input
139 ndarrays.
140 qagiu_meat
142 Signature: (double a(); double epsabs();double epsrel();int limit();
143 double [o] result(); double [o] abserr();int n();int [o] ierr();int gslwarn();; SV* function)
144 info not available
146 qagiu_meat does not process bad values. It will set the bad-value flag
148 of all output ndarrays if the flag is set for any of the input
149 ndarrays.
150 qagil_meat
152 Signature: (double b(); double epsabs();double epsrel();int limit();
153 double [o] result(); double [o] abserr();int n();int [o] ierr();int gslwarn();; SV* function)
154 info not available
156 qagil_meat does not process bad values. It will set the bad-value flag
158 of all output ndarrays if the flag is set for any of the input
159 ndarrays.
160 qawc_meat
162 Signature: (double a(); double b(); double c(); double epsabs();double epsrel();int limit();
163 double [o] result(); double [o] abserr();int n();int [o] ierr();int gslwarn();; SV* function)
164 info not available
166 qawc_meat does not process bad values. It will set the bad-value flag
168 of all output ndarrays if the flag is set for any of the input
169 ndarrays.
170 qaws_meat
172 Signature: (double a(); double b();double epsabs();double epsrel();int limit();
173 double [o] result(); double [o] abserr();int n();
174 double alpha(); double beta(); int mu(); int nu();int [o] ierr();int gslwarn();; SV* function)
175 info not available
177 qaws_meat does not process bad values. It will set the bad-value flag
179 of all output ndarrays if the flag is set for any of the input
180 ndarrays.
181 qawo_meat
183 Signature: (double a(); double b();double epsabs();double epsrel();int limit();
184 double [o] result(); double [o] abserr();int n();
185 int sincosopt(); double omega(); double L(); int nlevels();int [o] ierr();int gslwarn();; SV* function)
186 info not available
188 qawo_meat does not process bad values. It will set the bad-value flag
190 of all output ndarrays if the flag is set for any of the input
191 ndarrays.
192 qawf_meat
194 Signature: (double a(); double epsabs();int limit();
195 double [o] result(); double [o] abserr();int n();
196 int sincosopt(); double omega(); int nlevels();int [o] ierr();int gslwarn();; SV* function)
197 info not available
199 qawf_meat does not process bad values. It will set the bad-value flag
201 of all output ndarrays if the flag is set for any of the input
202 ndarrays.
203 gslinteg_qng
205 Non-adaptive Gauss-Kronrod integration
206 This function applies the Gauss-Kronrod 10-point, 21-point, 43-point
208 and 87-point integration rules in succession until an estimate of the
209 integral of f over ($la,$lb) is achieved within the desired absolute
210 and relative error limits, $epsabs and $epsrel. It is meant for fast
211 integration of smooth functions. It returns an array with the result,
212 an estimate of the absolute error, an error flag and the number of
213 function evaluations performed.
214 Usage:
216 ($res,$abserr,$ierr,$neval) = gslinteg_qng($function_ref,$la,$lb,
218 $epsrel,$epsabs,[{Warn => $warn}]);
219 Example:
221 my ($res,$abserr,$ierr,$neval) = gslinteg_qng(\&f,0,1,0,1e-9);
223 # with warnings on
224 my ($res,$abserr,$ierr,$neval) = gslinteg_qng(\&f,0,1,0,1e-9,{Warn => 'y'});
225 sub f{
227 my ($x) = @_;
228 return ($x**2.6)*log(1.0/$x);
229 }
230 gslinteg_qag
232 Adaptive integration
233 This function applies an integration rule adaptively until an estimate
235 of the integral of f over ($la,$lb) is achieved within the desired
236 absolute and relative error limits, $epsabs and $epsrel. On each
237 iteration the adaptive integration strategy bisects the interval with
238 the largest error estimate; the maximum number of allowed subdivisions
239 is given by the parameter $limit. The integration rule is determined
240 by the value of $key, which has to be one of (1,2,3,4,5,6) and
241 correspond to the 15, 21, 31, 41, 51 and 61 point Gauss-Kronrod rules
242 respectively. It returns an array with the result, an estimate of the
243 absolute error and an error flag.
244 Please check the GSL documentation for more information.
246 Usage:
248 ($res,$abserr,$ierr) = gslinteg_qag($function_ref,$la,$lb,$epsrel,
250 $epsabs,$limit,$key,[{Warn => $warn}]);
251 Example:
253 my ($res,$abserr,$ierr) = gslinteg_qag(\&f,0,1,0,1e-10,1000,1);
255 # with warnings on
256 my ($res,$abserr,$ierr) = gslinteg_qag(\&f,0,1,0,1e-10,1000,1,{Warn => 'y'});
257 sub f{
259 my ($x) = @_;
260 return ($x**2.6)*log(1.0/$x);
261 }
262 gslinteg_qags
264 Adaptive integration with singularities
265 This function applies the Gauss-Kronrod 21-point integration rule
267 adaptively until an estimate of the integral of f over ($la,$lb) is
268 achieved within the desired absolute and relative error limits, $epsabs
269 and $epsrel. The algorithm is such that it accelerates the convergence
270 of the integral in the presence of discontinuities and integrable
271 singularities. The maximum number of allowed subdivisions done by the
272 adaptive algorithm must be supplied in the parameter $limit.
273 Please check the GSL documentation for more information.
275 Usage:
277 ($res,$abserr,$ierr) = gslinteg_qags($function_ref,$la,$lb,$epsrel,
279 $epsabs,$limit,[{Warn => $warn}]);
280 Example:
282 my ($res,$abserr,$ierr) = gslinteg_qags(\&f,0,1,0,1e-10,1000);
284 # with warnings on
285 ($res,$abserr,$ierr) = gslinteg_qags(\&f,0,1,0,1e-10,1000,{Warn => 'y'});
286 sub f{
288 my ($x) = @_;
289 return ($x)*log(1.0/$x);
290 }
291 gslinteg_qagp
293 Adaptive integration with known singular points
294 This function applies the adaptive integration algorithm used by
296 gslinteg_qags taking into account the location of singular points until
297 an estimate of the integral of f over ($la,$lb) is achieved within the
298 desired absolute and relative error limits, $epsabs and $epsrel.
299 Singular points are supplied in the ndarray $points, whose endpoints
300 determine the integration range. So, for example, if the function has
301 singular points at x_1 and x_2 and the integral is desired from a to b
302 (a < x_1 < x_2 < b), $points = pdl(a,x_1,x_2,b). The maximum number of
303 allowed subdivisions done by the adaptive algorithm must be supplied in
304 the parameter $limit.
305 Please check the GSL documentation for more information.
307 Usage:
309 ($res,$abserr,$ierr) = gslinteg_qagp($function_ref,$points,$epsabs,
311 $epsrel,$limit,[{Warn => $warn}])
312 Example:
314 my $points = pdl(0,1,sqrt(2),3);
316 my ($res,$abserr,$ierr) = gslinteg_qagp(\&f,$points,0,1e-3,1000);
317 # with warnings on
318 ($res,$abserr,$ierr) = gslinteg_qagp(\&f,$points,0,1e-3,1000,{Warn => 'y'});
319 sub f{
321 my ($x) = @_;
322 my $x2 = $x**2;
323 my $x3 = $x**3;
324 return $x3 * log(abs(($x2-1.0)*($x2-2.0)));
325 }
326 gslinteg_qagi
328 Adaptive integration on infinite interval
329 This function estimates the integral of the function f over the
331 infinite interval (-\infty,+\infty) within the desired absolute and
332 relative error limits, $epsabs and $epsrel. After a transformation,
333 the algorithm of gslinteg_qags with a 15-point Gauss-Kronrod rule is
334 used. The maximum number of allowed subdivisions done by the adaptive
335 algorithm must be supplied in the parameter $limit.
336 Please check the GSL documentation for more information.
338 Usage:
340 ($res,$abserr,$ierr) = gslinteg_qagi($function_ref,$epsabs,
342 $epsrel,$limit,[{Warn => $warn}]);
343 Example:
345 my ($res,$abserr,$ierr) = gslinteg_qagi(\&myfn,1e-7,0,1000);
347 # with warnings on
348 ($res,$abserr,$ierr) = gslinteg_qagi(\&myfn,1e-7,0,1000,{Warn => 'y'});
349 sub myfn{
351 my ($x) = @_;
352 return exp(-$x - $x*$x) ;
353 }
354 gslinteg_qagiu
356 Adaptive integration on infinite interval
357 This function estimates the integral of the function f over the
359 infinite interval (la,+\infty) within the desired absolute and relative
360 error limits, $epsabs and $epsrel. After a transformation, the
361 algorithm of gslinteg_qags with a 15-point Gauss-Kronrod rule is used.
362 The maximum number of allowed subdivisions done by the adaptive
363 algorithm must be supplied in the parameter $limit.
364 Please check the GSL documentation for more information.
366 Usage:
368 ($res,$abserr,$ierr) = gslinteg_qagiu($function_ref,$la,$epsabs,
370 $epsrel,$limit,[{Warn => $warn}]);
371 Example:
373 my $alfa = 1;
375 my ($res,$abserr,$ierr) = gslinteg_qagiu(\&f,99.9,1e-7,0,1000);
376 # with warnings on
377 ($res,$abserr,$ierr) = gslinteg_qagiu(\&f,99.9,1e-7,0,1000,{Warn => 'y'});
378 sub f{
380 my ($x) = @_;
381 if (($x==0) && ($alfa == 1)) {return 1;}
382 if (($x==0) && ($alfa > 1)) {return 0;}
383 return ($x**($alfa-1))/((1+10*$x)**2);
384 }
385 gslinteg_qagil
387 Adaptive integration on infinite interval
388 This function estimates the integral of the function f over the
390 infinite interval (-\infty,lb) within the desired absolute and relative
391 error limits, $epsabs and $epsrel. After a transformation, the
392 algorithm of gslinteg_qags with a 15-point Gauss-Kronrod rule is used.
393 The maximum number of allowed subdivisions done by the adaptive
394 algorithm must be supplied in the parameter $limit.
395 Please check the GSL documentation for more information.
397 Usage:
399 ($res,$abserr,$ierr) = gslinteg_qagl($function_ref,$lb,$epsabs,
401 $epsrel,$limit,[{Warn => $warn}]);
402 Example:
404 my ($res,$abserr,$ierr) = gslinteg_qagil(\&myfn,1.0,1e-7,0,1000);
406 # with warnings on
407 ($res,$abserr,$ierr) = gslinteg_qagil(\&myfn,1.0,1e-7,0,1000,{Warn => 'y'});
408 sub myfn{
410 my ($x) = @_;
411 return exp($x);
412 }
413 gslinteg_qawc
415 Adaptive integration for Cauchy principal values
416 This function computes the Cauchy principal value of the integral of f
418 over (la,lb), with a singularity at c, I = \int_{la}^{lb} dx f(x)/(x -
419 c). The integral is estimated within the desired absolute and relative
420 error limits, $epsabs and $epsrel. The maximum number of allowed
421 subdivisions done by the adaptive algorithm must be supplied in the
422 parameter $limit.
423 Please check the GSL documentation for more information.
425 Usage:
427 ($res,$abserr,$ierr) = gslinteg_qawc($function_ref,$la,$lb,$c,$epsabs,$epsrel,$limit)
429 Example:
431 my ($res,$abserr,$ierr) = gslinteg_qawc(\&f,-1,5,0,0,1e-3,1000);
433 # with warnings on
434 ($res,$abserr,$ierr) = gslinteg_qawc(\&f,-1,5,0,0,1e-3,1000,{Warn => 'y'});
435 sub f{
437 my ($x) = @_;
438 return 1.0 / (5.0 * $x * $x * $x + 6.0) ;
439 }
440 gslinteg_qaws
442 Adaptive integration for singular functions
443 The algorithm in gslinteg_qaws is designed for integrands with
445 algebraic-logarithmic singularities at the end-points of an integration
446 region. Specifically, this function computes the integral given by I =
447 \int_{la}^{lb} dx f(x) (x-la)^alpha (lb-x)^beta log^mu (x-la) log^nu
448 (lb-x). The integral is estimated within the desired absolute and
449 relative error limits, $epsabs and $epsrel. The maximum number of
450 allowed subdivisions done by the adaptive algorithm must be supplied in
451 the parameter $limit.
452 Please check the GSL documentation for more information.
454 Usage:
456 ($res,$abserr,$ierr) =
458 gslinteg_qawc($function_ref,$alpha,$beta,$mu,$nu,$la,$lb,
459 $epsabs,$epsrel,$limit,[{Warn => $warn}]);
460 Example:
462 my ($res,$abserr,$ierr) = gslinteg_qaws(\&f,0,0,1,0,0,1,0,1e-7,1000);
464 # with warnings on
465 ($res,$abserr,$ierr) = gslinteg_qaws(\&f,0,0,1,0,0,1,0,1e-7,1000,{Warn => 'y'});
466 sub f{
468 my ($x) = @_;
469 if($x==0){return 0;}
470 else{
471 my $u = log($x);
472 my $v = 1 + $u*$u;
473 return 1.0/($v*$v);
474 }
475 }
476 gslinteg_qawo
478 Adaptive integration for oscillatory functions
479 This function uses an adaptive algorithm to compute the integral of f
481 over (la,lb) with the weight function sin(omega*x) or cos(omega*x) --
482 which of sine or cosine is used is determined by the parameter $opt
483 ('cos' or 'sin'). The integral is estimated within the desired
484 absolute and relative error limits, $epsabs and $epsrel. The maximum
485 number of allowed subdivisions done by the adaptive algorithm must be
486 supplied in the parameter $limit.
487 Please check the GSL documentation for more information.
489 Usage:
491 ($res,$abserr,$ierr) = gslinteg_qawo($function_ref,$omega,$sin_or_cos,
493 $la,$lb,$epsabs,$epsrel,$limit,[opt])
494 Example:
496 my $PI = 3.14159265358979323846264338328;
498 my ($res,$abserr,$ierr) = PDL::GSL::INTEG::gslinteg_qawo(\&f,10*$PI,'sin',0,1,0,1e-7,1000);
499 # with warnings on
500 ($res,$abserr,$ierr) = PDL::GSL::INTEG::gslinteg_qawo(\&f,10*$PI,'sin',0,1,0,1e-7,1000,{Warn => 'y'});
501 sub f{
503 my ($x) = @_;
504 if($x==0){return 0;}
505 else{ return log($x);}
506 }
507 gslinteg_qawf
509 Adaptive integration for Fourier integrals
510 This function attempts to compute a Fourier integral of the function f
512 over the semi-infinite interval [la,+\infty). Specifically, it attempts
513 tp compute I = \int_{la}^{+\infty} dx f(x)w(x), where w(x) is
514 sin(omega*x) or cos(omega*x) -- which of sine or cosine is used is
515 determined by the parameter $opt ('cos' or 'sin'). The integral is
516 estimated within the desired absolute error limit $epsabs. The maximum
517 number of allowed subdivisions done by the adaptive algorithm must be
518 supplied in the parameter $limit.
519 Please check the GSL documentation for more information.
521 Usage:
523 gslinteg_qawf($function_ref,$omega,$sin_or_cos,$la,$epsabs,$limit,[opt])
525 Example:
527 my ($res,$abserr,$ierr) = gslinteg_qawf(\&f,$PI/2.0,'cos',0,1e-7,1000);
529 # with warnings on
530 ($res,$abserr,$ierr) = gslinteg_qawf(\&f,$PI/2.0,'cos',0,1e-7,1000,{Warn => 'y'});
531 sub f{
533 my ($x) = @_;
534 if ($x == 0){return 0;}
535 return 1.0/sqrt($x)
536 }
537
539 Feedback is welcome. Log bugs in the PDL bug database (the database is
540 always linked from <http://pdl.perl.org>).
541
543 PDL
544 The GSL documentation for numerical integration is online at
546 <https://www.gnu.org/software/gsl/doc/html/integration.html>
547
549 This file copyright (C) 2003,2005 Andres Jordan <ajordan@eso.org> All
550 rights reserved. There is no warranty. You are allowed to redistribute
551 this software documentation under certain conditions. For details, see
552 the file COPYING in the PDL distribution. If this file is separated
553 from the PDL distribution, the copyright notice should be included in
554 the file.
555 The GSL integration routines were written by Brian Gough. QUADPACK was
557 written by Piessens, Doncker-Kapenga, Uberhuber and Kahaner.
558perl v5.34.0 2021-08-16 INTEG(3)