ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_quadPack::isFailedQuad Interface Reference

Compute the 1D integral of the input scalar (potentially singular) integrand getFunc on a finite or semi/fully-infinite interval (a, b) and estimate its absolute error via the requested adaptive global quadrature rule. More...

Detailed Description

Compute the 1D integral of the input scalar (potentially singular) integrand getFunc on a finite or semi/fully-infinite interval (a, b) and estimate its absolute error via the requested adaptive global quadrature rule.

This interface provides a simple nimble powerful wrapper for the higher-performance but lower-level getQuadErr interface, much like the \(\ms{integral}\) routine of MATLAB.
The algorithm uses an Adaptive Global Quadrature with Gauss-Kronrod 10-21 quadrature rules combined with the Epsilon extrapolation method of Wynn (1961) to evaluate the integral.
See the documentation of getQuadErr for further details on the integration methodologies and implementations.

Parameters
getFunc: See the description of the corresponding argument in the documentation of getQuadErr.
[in]lb: The input scalar of type real of the same kind as integral, representing the lower limit of integration.
Set lb = -huge(lb) or to the IEEE-compliant negative infinity (lb =getInfNeg(lb)) to imply \(-\infty\) as the lower bound of integration.
[in]ub: The input scalar of type real of the same kind as integral, representing the upper limit of integration.
Set ub = huge(ub) or to the IEEE-compliant positive infinity (ub =getInfPos(lb)) to imply \(+\infty\) as the upper bound of integration.
[out]integral: See the description of the corresponding argument in the documentation of getQuadErr.
[in]abstol: See the description of the corresponding argument in the documentation of getQuadErr.
(optional, default = 0)
[in]reltol: See the description of the corresponding argument in the documentation of getQuadErr.
(optional, default = epsilon(real(0, kind(integral)))**(2./3.))
[in]help: See the description of the corresponding argument in the documentation of getQuadErr.
(optional. If missing, the procedure attempts to compute the integral without explicit assumptions about the singularities or discontinuities.)
[out]abserr: See the description of the corresponding argument in the documentation of getQuadErr.
(optional. If missing, the estimated error of the integral will not be returned.)
[out]neval: See the description of the corresponding argument in the documentation of getQuadErr.
(optional. If missing, the number of function evaluations will not be returned.)
[out]nint: See the description of the corresponding argument in the documentation of getQuadErr.
(optional. If missing, the number of subintervals formed will not be returned.)
[out]msg: The output scalar argument of type character of default kind SK of arbitrary length type parameter that is set to a diagnostic message if the integration fails to converge.
A length type parameter of 127 is sufficient to capture all error messages.
If msg has shorter length parameter, the output message will be trimmed from the end, otherwise padded with blanks as necessary.
(optional. If missing, no diagnostic message will be returned.)
Returns
failed : The output scalar of type logical of default kind LK, that is set to .true. if and only if the integration fails to converge within the requested tolerances.
Otherwise, it is set .false. if the integration succeeds with no errors.
See the description of the output argument err of getQuadErr for information on the kinds of integration failures that can happen.


Possible calling interfaces

use pm_kind, only: LK
logical(LK) :: failed
failed = isFailedQuad(getFunc, lb, ub , integral, abserr = abserr, abstol = abstol, reltol = reltol, neval = neval, nint = nint, msg = msg)
failed = isFailedQuad(getFunc, lb, ub, help , integral, abserr = abserr, abstol = abstol, reltol = reltol, neval = neval, nint = nint, msg = msg)
Compute the 1D integral of the input scalar (potentially singular) integrand getFunc on a finite or s...
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541
This module contains classes and procedures for non-adaptive and adaptive global numerical quadrature...
Warning
All conditions in the warning section of the documentation of getQuadErr also apply to this interface.
These conditions are verified only if the library is built with the preprocessor macro CHECK_ENABLED=1.
If the target function contains points of difficulties, singularities, or discontinuities, user must ensure the abscissas of the specified Gauss-Kronrod rule do not match such points.
Particularly, computing an integrand at its singularities can lead to undefined values that can lead to unexpected segmentation fault or propagation of NaN values within the computational flow or other strange errors that can be extremely difficult to debug.
A simple check can be added within the target integrand implementations to ensure no such difficulty point matches an input value at which the function must be evaluated.
Alternatively, one should consider using the adaptive integration routines isFailedQuad or getQuadErr while setting their input help arguments to the points of difficulties of the integrand.
Remarks
The procedures under discussion are impure.
See also
getQuadGK
getQuadErr
getQuadRomb
isFailedQuad


Example usage

1program example
2
3 use pm_kind, only: SK, IK, RKH, RKG => RKH ! All other real kinds are also supported.
4 use pm_quadTest, only: int1_type
5 use pm_quadTest, only: int2_type
6 use pm_quadTest, only: int3_type
7 use pm_quadTest, only: int4_type
8 use pm_quadTest, only: int5_type
9 use pm_quadTest, only: int6_type
10 use pm_quadTest, only: int7_type
11 use pm_quadTest, only: int8_type
12 use pm_quadTest, only: int9_type
24 use pm_io, only: display_type
25
26 implicit none
27
28 type(display_type) :: disp
29 disp = display_type(file = "main.out.F90")
30
31 call disp%skip()
32 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
33 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
34 call disp%show("! Compare integration via the Adaptive Global Gauss-Kronrod Quadrature method with (QAG) and without (QAGS) extrapolation acceleration.")
35 call disp%show("! Note the significant efficiency improvement in the integration via QAGS when the integrand has integrable singularity at some points.")
36 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
37 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
38 call disp%skip()
39
40 call test_isFailedQuad(disp, int1_type(-4._RKH, +4._RKH))
41 call test_isFailedQuad(disp, int2_type(+2._RKH, +3._RKH))
42 call test_isFailedQuad(disp, int3_type(ub = +3._RKH))
44 call test_isFailedQuad(disp, intSinCos_type(-4_IK, +4_IK, a = 10._RKH, b = -10._RKH))
45 call test_isFailedQuad(disp, intNormPDF_type(-3._RKH, +3._RKH))
46 call test_isFailedQuad(disp, intLogNormPDF_type(lb = exp(-6._RKG), ub = exp(+6._RKG)))
47
48 call disp%skip()
49 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
50 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
51 call disp%show("! Compare integration via:")
52 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with automatic break point handling and extrapolation acceleration (QAGS),")
53 call disp%show("! with,")
54 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with specified break point handling and extrapolation acceleration (QAGP).")
55 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
56 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
57 call disp%skip()
58
59 call test_isFailedQuad(disp, int5_type(0._RKG, 3._RKG))
60 call test_isFailedQuad(disp, int5_type(-2._RKG, +5._RKG))
61
62 call disp%skip()
63 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
64 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
65 call disp%show("! Compare integration via:")
66 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a semi-infinite range `(lb, +Inf)`,")
67 call disp%show("! with,")
68 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a semi-infinite range `(lb, +Inf)`.")
69 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
70 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
71 call disp%skip()
72
74 !call disp%show("call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = +1.0_RKG, beta = 3._RKG))")
75 ! call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = +1.0_RKG, beta = 3._RKG))
76 !call disp%show("call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = +0.5_RKG, beta = 3._RKG))")
77 ! call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = +0.5_RKG, beta = 3._RKG))
78 !call disp%show("call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = +0.0_RKG, beta = 3._RKG))")
79 ! call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = +0.0_RKG, beta = 3._RKG))
80 !call disp%show("call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = -0.5_RKG, beta = 3._RKG))")
81 ! call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = -0.5_RKG, beta = 3._RKG))
82 !call disp%show("call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = -1.0_RKG, beta = 3._RKG))")
83 ! call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = -1.0_RKG, beta = 3._RKG))
84 !call disp%show("call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = -1.5_RKG, beta = 3._RKG))")
85 ! call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = -1.5_RKG, beta = 3._RKG))
86 !call disp%show("call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = -2.0_RKG, beta = 3._RKG))")
87 ! call test_isFailedQuad(disp, intGamUpp_type(lb = 1._RKG, alpha = -2.0_RKG, beta = 3._RKG))
92
93 call disp%skip()
94 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
95 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
96 call disp%show("! Compare integration via:")
97 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a semi-infinite range `(-Inf, ub)`,")
98 call disp%show("! with,")
99 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a semi-infinite range `(-Inf, ub)`.")
100 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
101 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
102 call disp%skip()
103
104 call test_isFailedQuad(disp, int8_type())
106 call test_isFailedQuad(disp, intDoncker2_type(ub = -1._RKG))
107
108 call disp%skip()
109 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
110 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
111 call disp%show("! Compare integration via:")
112 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a full-infinite range `(-Inf, +Inf)`,")
113 call disp%show("! with,")
114 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a full-infinite range `(-Inf, +Inf)`.")
115 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
116 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
117 call disp%skip()
118
122 call test_isFailedQuad(disp, intGenExpGammaPDF_type(kappa = 5._RKH, invOmega = 0.1_RKH, logSigma = 2._RKH))
124
125 call disp%skip()
126 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
127 call disp%show("! Compute the Cauchy Principal Value over a finite interval.")
128 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
129 call disp%skip()
130
133
134end program example
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11726
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11508
Run the adaptive global quadrature methods for the specified input integrand object.
This module contains classes and procedures for input/output (IO) or generic display operations on st...
Definition: pm_io.F90:252
type(display_type) disp
This is a scalar module variable an object of type display_type for general display.
Definition: pm_io.F90:11393
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RKH
The scalar integer constant of intrinsic default kind, representing the highest-precision real kind t...
Definition: pm_kind.F90:858
This module contains a collection of interesting or challenging integrands for testing or examining t...
Definition: pm_quadTest.F90:54
Generate and return an object of type display_type.
Definition: pm_io.F90:10282
This is the derived type for generating test integrand objects of the algebraic form as described bel...
This is the derived type for generating test integrand objects of algebraic form as described below.
This is the derived type for generating test integrand objects of algebraic form as described below.
This is the derived type for generating test integrand objects of the following algebraic form.
This is the derived type for generating test integrand objects of the following algebraic form.
This is the derived type for generating test integrand objects of the following algebraic form.
This is the derived type for generating test integrand objects of the following algebraic form.
This is the derived type for generating test integrand objects of the following algebraic form.
This is the derived type for generating test integrand objects of the algebraic form as described bel...
This is the derived type for generating test integrand objects of the algebraic form as described bel...
This is the derived type for generating test integrand objects of algebraic form as described below.
This is the derived type for generating test integrand objects of algebraic form as described below.
This is the derived type for generating test integrand objects of the following algebraic form.
This is the derived type for generating test integrand objects of the Probability Density Function of...
This is the derived type for generating test integrand objects of the Probability Density Function of...
This is the derived type for generating test integrand objects of the Probability Density Function of...
This is the derived type for generating test integrand objects of the sum of five Probability Density...
This is the derived type for generating test integrand objects of the trigonometric form as described...

Example Unix compile command via Intel ifort compiler
1#!/usr/bin/env sh
2rm main.exe
3ifort -fpp -standard-semantics -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example Windows Batch compile command via Intel ifort compiler
1del main.exe
2set PATH=..\..\..\lib;%PATH%
3ifort /fpp /standard-semantics /O3 /I:..\..\..\include main.F90 ..\..\..\lib\libparamonte*.lib /exe:main.exe
4main.exe

Example Unix / MinGW compile command via GNU gfortran compiler
1#!/usr/bin/env sh
2rm main.exe
3gfortran -cpp -ffree-line-length-none -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example output
1
2!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4! Compare integration via the Adaptive Global Gauss-Kronrod Quadrature method with (QAG) and without (QAGS) extrapolation acceleration.
5! Note the significant efficiency improvement in the integration via QAGS when the integrand has integrable singularity at some points.
6!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
7!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
8
9integrand%desc
10"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
11
12!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
13! int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)
14!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
15
16
17! Regular adaptive global quadrature.
18
19if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
20getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
210.592319847946765693983261139952398973, 0.592319847946765693983261139951656719, 0.122995853428115837177070158359807283E-22, 0.742253400566840697767773899072230989E-30 (unbiased)? TRUE
22numFuncEval
23+441
24
25
26! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
27
28if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
29getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
300.592319847946765693983261139952398973, 0.592319847946765693983261139951656623, 0.122995853428115837177070158261927716E-22, 0.742349697064060059560426697969360235E-30 (unbiased)? TRUE
31numFuncEval
32+441
33
34integrand%desc
35"int2_type: an algebraic integrand of the form f(x) = 1 / sqrt(a - b * x) for x in (0, a / b) with a > 0 and b > 0 with a singularity at the upper bound of integration"
36
37!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
38! int2_type: an algebraic integrand of the form f(x) = 1 / sqrt(a - b * x) for x in (0, a / b) with a > 0 and b > 0 with a singularity at the upper bound of integration
39!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
40
41
42! Regular adaptive global quadrature.
43
44if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
45
46 ******** Integration did NOT converge. ********
47
48getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
490.942809041582063365867792482806465323, 0.942809041582063358489159110539769991, 0.274189625236721766830515122145815584E-15, 0.737863337226669533225571680774701636E-17 (unbiased)? TRUE
50numFuncEval
51+4263
52
53
54! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
55
56if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
57getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
580.942809041582063365867792482806465323, 0.942809041582063365867792482806451456, 0.747866619686254516023298199584791888E-27, 0.138666955995880981420030411866114767E-31 (unbiased)? TRUE
59numFuncEval
60+399
61
62integrand%desc
63"int3_type: an algebraic integrand of the form f(x) = log(x) / sqrt(x) for x in (0, ub) with a singularity at the lower limit of integration"
64
65!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
66! int3_type: an algebraic integrand of the form f(x) = log(x) / sqrt(x) for x in (0, ub) with a singularity at the lower limit of integration
67!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
68
69
70! Regular adaptive global quadrature.
71
72if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
73getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
74-3.12249862669012531099145873075559453, -3.12249862669012531099145800385027096, 0.871742302516757365494019872303469907E-22, 0.726905323565232669836577649328480277E-24 (unbiased)? TRUE
75numFuncEval
76+8127
77
78
79! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
80
81if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
82getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
83-3.12249862669012531099145873075559453, -3.12249862669012531099145873075553329, 0.638742935999207969456768107027576605E-25, 0.612445722315141001271800985742006887E-31 (unbiased)? TRUE
84numFuncEval
85+567
86
87integrand%desc
88"int4_type: an algebraic integrand of the form f(x) = log(x) / (1. + log(x)**2)**2 for x in (0, 1) with a singularity at the lower limit"
89
90!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
91! int4_type: an algebraic integrand of the form f(x) = log(x) / (1. + log(x)**2)**2 for x in (0, 1) with a singularity at the lower limit
92!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
93
94
95! Regular adaptive global quadrature.
96
97if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
98getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
99-0.189275187882093321180367135892330336, -0.189275187882093321180367139193492826, 0.439350295463286972817507956603694862E-23, 0.330116248974192359134514995223265581E-26 (unbiased)? TRUE
100numFuncEval
101+2457
102
103
104! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
105
106if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
107getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
108-0.189275187882093321180367135892330336, -0.189275187882093321180367392427969773, 0.530509501355781083785828386429659476E-23, 0.256535639436815431080742250597631077E-24 (unbiased)? TRUE
109numFuncEval
110+1029
111
112integrand%desc
113"intSinCos_typer(): a highly oscillatory integrand of the form f(x) = cos(a * sin(b * x)) for x in (lb, ub)"
114
115!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
116! intSinCos_typer(): a highly oscillatory integrand of the form f(x) = cos(a * sin(b * x)) for x in (lb, ub)
117!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
118
119
120! Regular adaptive global quadrature.
121
122if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
123getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
124-6.18103992684276605416490727881801341, -6.18103992684276605416490727881796180, 0.204459144829722446056783263048240039E-21, 0.516149225095779208619002088612760521E-31 (unbiased)? TRUE
125numFuncEval
126+40509
127
128
129! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
130
131if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
132getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
133-6.18103992684276605416490727881801341, -6.18103992684276605416490727881796180, 0.204459144829722446056783263048240039E-21, 0.516149225095779208619002088612760521E-31 (unbiased)? TRUE
134numFuncEval
135+40509
136
137integrand%desc
138"intNormPDF_type: f(x) = exp(-0.5 * (log(x) - mu)**2 / sigma**2) / (sigma * sqrt(2 * acos(-1.))) for x in (lb, ub)"
139
140!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
141! intNormPDF_type: f(x) = exp(-0.5 * (log(x) - mu)**2 / sigma**2) / (sigma * sqrt(2 * acos(-1.))) for x in (lb, ub)
142!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
143
144
145! Regular adaptive global quadrature.
146
147if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
148getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1490.997300203936739810946696370464810073, 0.997300203936739810946696370464810073, 0.863216418769238025171827058295059452E-23, 0.00000000000000000000000000000000000 (unbiased)? TRUE
150numFuncEval
151+147
152
153
154! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
155
156if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
157getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1580.997300203936739810946696370464810073, 0.997300203936739810946696370464810073, 0.863216418769238025171827058295059452E-23, 0.00000000000000000000000000000000000 (unbiased)? TRUE
159numFuncEval
160+147
161
162integrand%desc
163"intLogNormPDF_type: f(x) = exp(-0.5 * (log(x) - mu)**2 / sigma**2) / (x * sigma * sqrt(2 * acos(-1.))) for x in (0 <= lb, ub)"
164
165!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
166! intLogNormPDF_type: f(x) = exp(-0.5 * (log(x) - mu)**2 / sigma**2) / (x * sigma * sqrt(2 * acos(-1.))) for x in (0 <= lb, ub)
167!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
168
169
170! Regular adaptive global quadrature.
171
172if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
173getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1740.999999998026824709924603718598271714, 0.999999998026824709924603718597611697, 0.287685113478269648540938767224899488E-23, 0.660016191941505726842283640923854592E-30 (unbiased)? TRUE
175numFuncEval
176+1113
177
178
179! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
180
181if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
182getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1830.999999998026824709924603718598271714, 0.999999998026824709924603718597611601, 0.287685113534716101717891864391506212E-23, 0.660112488438725088634936439820983838E-30 (unbiased)? TRUE
184numFuncEval
185+1113
186
187
188!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
189!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
190! Compare integration via:
191! the Adaptive Global Gauss-Kronrod Quadrature method with automatic break point handling and extrapolation acceleration (QAGS),
192! with,
193! the Adaptive Global Gauss-Kronrod Quadrature method with specified break point handling and extrapolation acceleration (QAGP).
194!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
195!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
196
197integrand%desc
198"int5_type: an algebraic integrand of the form f(x) = x**3 log(abs((x**2 - 1) * (x**2 - 2))) for x in (0., 3.) with 4 possible singularities: [-sqrt(2.), -1., 1., sqrt(2.)]"
199
200!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
201! int5_type: an algebraic integrand of the form f(x) = x**3 log(abs((x**2 - 1) * (x**2 - 2))) for x in (0., 3.) with 4 possible singularities: [-sqrt(2.), -1., 1., sqrt(2.)]
202!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
203
204
205! Regular adaptive global quadrature.
206
207if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
208getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
20952.7407483834714449977291997202299809, 52.7407483834714449977291765299254724, 0.172884396198263085498243239907807437E-20, 0.231903045085095141212385933523689371E-22 (unbiased)? TRUE
210numFuncEval
211+10101
212
213
214! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
215
216if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
217getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
21852.7407483834714449977291997202299809, 52.7407483834714449977292388198386134, 0.157545601084532647841636140405395188E-20, 0.390996086324092930047679517794757520E-22 (unbiased)? TRUE
219numFuncEval
220+9975
221
222
223! Assisted adaptive global quadrature by specifying points of difficulties of the integrand.
224
225break = integrand%break
226break
227+1.00000000000000000000000000000000000, +1.41421356237309504880168872420969798
228if (isFailedQuad(getFunc, lb, ub, break, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
229getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
23052.7407483834714449977291997202299809, 52.7407483834714449977291997201071898, 0.145935898758650954111037868386060001E-21, 0.122791130278308118836119374489798694E-27 (unbiased)? TRUE
231numFuncEval
232+2163
233
234integrand%desc
235"int5_type: an algebraic integrand of the form f(x) = x**3 log(abs((x**2 - 1) * (x**2 - 2))) for x in (-2., 5.) with 4 possible singularities: [-sqrt(2.), -1., 1., sqrt(2.)]"
236
237!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
238! int5_type: an algebraic integrand of the form f(x) = x**3 log(abs((x**2 - 1) * (x**2 - 2))) for x in (-2., 5.) with 4 possible singularities: [-sqrt(2.), -1., 1., sqrt(2.)]
239!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
240
241
242! Regular adaptive global quadrature.
243
244if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
245getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
246808.362636933094758891687400938039628, 808.362636933094758891687268473760062, 0.255379463719878359511890040859317480E-19, 0.132464279565495221439626229092202295E-21 (unbiased)? TRUE
247numFuncEval
248+19047
249
250
251! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
252
253if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
254getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
255808.362636933094758891687400938039628, 808.362636933094758891681251389789229, 0.334323235794388144484139093037256979E-20, 0.614954825039840811287942669490845986E-20 (unbiased)? FALSE
256numFuncEval
257+17955
258
259
260! Assisted adaptive global quadrature by specifying points of difficulties of the integrand.
261
262break = integrand%break
263break
264-1.41421356237309504880168872420969798, -1.00000000000000000000000000000000000, +1.00000000000000000000000000000000000, +1.41421356237309504880168872420969798
265if (isFailedQuad(getFunc, lb, ub, break, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
266getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
267808.362636933094758891687400938039628, 808.362636933094758891687400937280645, 0.990216586020630180919255548708256071E-21, 0.758982798435765983281759345872700701E-27 (unbiased)? TRUE
268numFuncEval
269+3045
270
271
272!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
273!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
274! Compare integration via:
275! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a semi-infinite range `(lb, +Inf)`,
276! with,
277! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a semi-infinite range `(lb, +Inf)`.
278!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
279!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
280
281integrand%desc
282"intGamUpp_type: an algebraic integrand of the form f(x; lb, alpha, beta) = (x / lb)**alpha * exp(-beta * (x - lb)) for x in (lb, +Inf), lb > 0."
283
284!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
285! intGamUpp_type: an algebraic integrand of the form f(x; lb, alpha, beta) = (x / lb)**alpha * exp(-beta * (x - lb)) for x in (lb, +Inf), lb > 0.
286!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
287
288
289! Regular adaptive global quadrature.
290
291if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
292getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
2932.00000000000000000000000000000000039, 2.00000000000000000000000000000000039, 0.503988260578339162149833493686791185E-24, 0.00000000000000000000000000000000000 (unbiased)? TRUE
294numFuncEval
295+483
296
297
298! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
299
300if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
301getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3022.00000000000000000000000000000000039, 2.00000000000000000000000000000000039, 0.503988260578338140171031197806738174E-24, 0.00000000000000000000000000000000000 (unbiased)? TRUE
303numFuncEval
304+483
305
306integrand%desc
307"intLogNormPDF_type: f(x) = exp(-0.5 * (log(x) - mu)**2 / sigma**2) / (x * sigma * sqrt(2 * acos(-1.))) for x in (0 <= lb, ub)"
308
309!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
310! intLogNormPDF_type: f(x) = exp(-0.5 * (log(x) - mu)**2 / sigma**2) / (x * sigma * sqrt(2 * acos(-1.))) for x in (0 <= lb, ub)
311!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
312
313
314! Regular adaptive global quadrature.
315
316if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
317getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3181.00000000000000000000000000000000000, 1.00000000000000000000000005509718469, 0.185875326636098607607964145008442373E-22, 0.550971846938518212456720189377984095E-25 (unbiased)? TRUE
319numFuncEval
320+735
321
322
323! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
324
325if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
326getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3271.00000000000000000000000000000000000, 1.00000000000000000000000005509718469, 0.185875326636098606689609183428530258E-22, 0.550971846938518212456720189377984095E-25 (unbiased)? TRUE
328numFuncEval
329+735
330
331integrand%desc
332"intDoncker1_type: f(x) = 1 / (1 + x) / sqrt(x) for x in (0 <= lb, ub) with a square-root singularity at 0"
333
334!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
335! intDoncker1_type: f(x) = 1 / (1 + x) / sqrt(x) for x in (0 <= lb, ub) with a square-root singularity at 0
336!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
337
338
339! Regular adaptive global quadrature.
340
341if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
342
343 ******** Integration did NOT converge. ********
344
345getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3463.14159265358979323846264338327950280, 3.14159265358979319764120376296608615, 0.114874526811124567747394666325693956E-14, 0.408214396203134166514689154707924291E-16 (unbiased)? TRUE
347numFuncEval
348+8463
349
350
351! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
352
353if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
354getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3553.14159265358979323846264338327950280, 3.14159265358979323846264338346019817, 0.105499447763580851205123272015818952E-22, 0.180695369614276997099759184920380085E-27 (unbiased)? TRUE
356numFuncEval
357+2121
358
359integrand%desc
360"int6_type: an algebraic integrand of the form f(x) = log(x) / (1 + 100 * x**2) for x in (0, +Inf) with a singularity at the lower limit"
361
362!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
363! int6_type: an algebraic integrand of the form f(x) = log(x) / (1 + 100 * x**2) for x in (0, +Inf) with a singularity at the lower limit
364!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
365
366
367! Regular adaptive global quadrature.
368
369if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
370getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
371-0.361689220620773240624502327513089541, -0.361689220620773240624502327908133835, 0.118682246318398965349923011739062857E-22, 0.395044294632030724453129549211143496E-27 (unbiased)? TRUE
372numFuncEval
373+6657
374
375
376! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
377
378if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
379getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
380-0.361689220620773240624502327513089541, -0.361689220620773240624502327510642406, 0.120941504830908329737177699048668687E-23, 0.244713473558703155578925197329697327E-29 (unbiased)? TRUE
381numFuncEval
382+1659
383
384integrand%desc
385"int7_type: an algebraic integrand of the form f(x) = -log(abs((1 - x**2) * (1 - 2 * x**2)) / x**4) / x**5 for x in (1./3., +Inf) with two singularities at [1 / sqrt(2), 1]"
386
387!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
388! int7_type: an algebraic integrand of the form f(x) = -log(abs((1 - x**2) * (1 - 2 * x**2)) / x**4) / x**5 for x in (1./3., +Inf) with two singularities at [1 / sqrt(2), 1]
389!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
390
391
392! Regular adaptive global quadrature.
393
394if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
395getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
39652.7407483834714449977291997202299809, 52.7407483834714449977292221649683846, 0.168283512274674281393657004006737318E-20, 0.224447384036677384031379327649902600E-22 (unbiased)? TRUE
397numFuncEval
398+10143
399
400
401! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
402
403if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
404getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
40552.7407483834714449977291997202299809, 52.7407483834714449977306706202548049, 0.169942945536630143394077258382982232E-20, 0.147090002482391233646949135657205724E-20 (unbiased)? TRUE
406numFuncEval
407+9681
408
409
410! Assisted adaptive global quadrature by specifying points of difficulties of the integrand.
411
412break = integrand%break
413break
414+0.707106781186547524400844362104849088, +1.00000000000000000000000000000000000
415if (isFailedQuad(getFunc, lb, ub, break, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
416getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
41752.7407483834714449977291997202299809, 52.7407483834714449977291997199544774, 0.550213208024888852525213619612864128E-21, 0.275503508172616333885316422295596819E-27 (unbiased)? TRUE
418numFuncEval
419+2499
420
421
422!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
423!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
424! Compare integration via:
425! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a semi-infinite range `(-Inf, ub)`,
426! with,
427! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a semi-infinite range `(-Inf, ub)`.
428!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
429!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
430
431integrand%desc
432"int8_type: an algebraic integrand of the form f(x) = log(abs((1 - x**2) * (1 - 2 * x**2)) / x**4) / x**5 for x in (-Inf, -1./3.) with two singularities at [-1, -1 / sqrt(2)]"
433
434!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
435! int8_type: an algebraic integrand of the form f(x) = log(abs((1 - x**2) * (1 - 2 * x**2)) / x**4) / x**5 for x in (-Inf, -1./3.) with two singularities at [-1, -1 / sqrt(2)]
436!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
437
438
439! Regular adaptive global quadrature.
440
441if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
442getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
44352.7407483834714449977291997202299809, 52.7407483834714449977292221649683846, 0.168283512274674281393657004006737318E-20, 0.224447384036677384031379327649902600E-22 (unbiased)? TRUE
444numFuncEval
445+10143
446
447
448! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
449
450if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
451getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
45252.7407483834714449977291997202299809, 52.7407483834714449977306706202548049, 0.169942945536630143394077258382982232E-20, 0.147090002482391233646949135657205724E-20 (unbiased)? TRUE
453numFuncEval
454+9681
455
456
457! Assisted adaptive global quadrature by specifying points of difficulties of the integrand.
458
459break = integrand%break
460break
461-1.00000000000000000000000000000000000, -0.707106781186547524400844362104849088
462if (isFailedQuad(getFunc, lb, ub, break, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
463getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
46452.7407483834714449977291997202299809, 52.7407483834714449977291997199544774, 0.550213208024888852525213619612864128E-21, 0.275503508172616333885316422295596819E-27 (unbiased)? TRUE
465numFuncEval
466+2499
467
468integrand%desc
469"intDoncker2_type: f(x) = exp(x) / sqrt(-x) for x in (lb, ub <= 0) with a square-root singularity at 0"
470
471!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
472! intDoncker2_type: f(x) = exp(x) / sqrt(-x) for x in (lb, ub <= 0) with a square-root singularity at 0
473!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
474
475
476! Regular adaptive global quadrature.
477
478if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
479
480 ******** Integration did NOT converge. ********
481
482getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
4831.77245385090551602729816748334114514, 1.77245385090551600687177993068626502, 0.574360549312137963665679582228894917E-15, 0.204263875526548801162109406012147667E-16 (unbiased)? TRUE
484numFuncEval
485+4431
486
487
488! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
489
490if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
491getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
4921.77245385090551602729816748334114514, 1.77245385090551602729816748322058404, 0.316674820796741424790461948518443858E-23, 0.120561095995702138441865857630079607E-27 (unbiased)? TRUE
493numFuncEval
494+1575
495
496integrand%desc
497"intDoncker2_type: f(x) = exp(x) / sqrt(-x) for x in (lb, ub <= 0) with a square-root singularity at 0"
498
499!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
500! intDoncker2_type: f(x) = exp(x) / sqrt(-x) for x in (lb, ub <= 0) with a square-root singularity at 0
501!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
502
503
504! Regular adaptive global quadrature.
505
506if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
507getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
5080.278805585280661976499232611077439155, 0.278805585280661976499232611077447244, 0.136847162139771908855491210817521382E-24, 0.808890576642639058283510735885669474E-32 (unbiased)? TRUE
509numFuncEval
510+441
511
512
513! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
514
515if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
516getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
5170.278805585280661976499232611077439155, 0.278805585280661976499232611077447292, 0.136847162139771908864182157536569162E-24, 0.813705401503607147916150680742131792E-32 (unbiased)? TRUE
518numFuncEval
519+441
520
521
522!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
523!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
524! Compare integration via:
525! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a full-infinite range `(-Inf, +Inf)`,
526! with,
527! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a full-infinite range `(-Inf, +Inf)`.
528!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
529!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
530
531integrand%desc
532"int9_type: an algebraic piecewise integrand of the form f(x) = log(abs((1 - x**2) * (1 - 2 * x**2)) / x**4) / x**5 for x in (1./3., +Inf) and 1 / (acos(-1) * sqrt(-(x+10) * (x+9))) for x in (-10, 9), otherwise 0, with two singularities at [-10, -9, 1/3., 1 / sqrt(2), 1]"
533
534!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
535! int9_type: an algebraic piecewise integrand of the form f(x) = log(abs((1 - x**2) * (1 - 2 * x**2)) / x**4) / x**5 for x in (1./3., +Inf) and 1 / (acos(-1) * sqrt(-(x+10) * (x+9))) for x in (-10, 9), otherwise 0, with two singularities at [-10, -9, 1/3., 1 / sqrt(2), 1]
536!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
537
538
539! Regular adaptive global quadrature.
540
541if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
542
543 ******** Integration did NOT converge. ********
544
545getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
54653.7407483834714449977291997202299809, 53.7407483834714446968169568664147327, 0.297103187254200474475915711578675253E-14, 0.300912242853815248241860288447560029E-15 (unbiased)? TRUE
547numFuncEval
548+19593
549
550
551! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
552
553if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
554
555 ******** Integration did NOT converge. ********
556
557getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
55853.7407483834714449977291997202299809, 53.7407483834714786408368937861800291, 0.597341489232741586490157071452958427E-12, 0.336431076940659500481100214330126955E-13 (unbiased)? TRUE
559numFuncEval
560+16653
561
562
563! Assisted adaptive global quadrature by specifying points of difficulties of the integrand.
564
565break = integrand%break
566break
567-10.0000000000000000000000000000000000, -9.00000000000000000000000000000000000, +0.333333333333333333333333333333333317, +0.707106781186547524400844362104849088, +1.00000000000000000000000000000000000
568if (isFailedQuad(getFunc, lb, ub, break, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
569getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
57053.7407483834714449977291997202299809, 53.7407483834714449977291998296210022, 0.511241450110203608013944520330859404E-22, 0.109391021252254961020398451633034320E-24 (unbiased)? TRUE
571numFuncEval
572+5796
573
574integrand%desc
575"intNormPDF_type: f(x) = exp(-0.5 * (log(x) - mu)**2 / sigma**2) / (sigma * sqrt(2 * acos(-1.))) for x in (lb, ub)"
576
577!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
578! intNormPDF_type: f(x) = exp(-0.5 * (log(x) - mu)**2 / sigma**2) / (sigma * sqrt(2 * acos(-1.))) for x in (lb, ub)
579!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
580
581
582! Regular adaptive global quadrature.
583
584if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
585getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
5861.00000000000000000000000000000000000, 0.999999999999999999999999999999999133, 0.386431949672064481044988439592369946E-24, 0.866668474974256133875190074163217293E-33 (unbiased)? TRUE
587numFuncEval
588+483
589
590
591! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
592
593if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
594getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
5951.00000000000000000000000000000000000, 0.999999999999999999999999999999999133, 0.386431949672064466119465327815932357E-24, 0.866668474974256133875190074163217293E-33 (unbiased)? TRUE
596numFuncEval
597+483
598
599integrand%desc
600"intGenExpGammaPDF_type: f(x) = GenExpGamma(x; kappa = 1., omega = 1., logSigma = 0.) for x in (-Inf, Inf)"
601
602!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
603! intGenExpGammaPDF_type: f(x) = GenExpGamma(x; kappa = 1., omega = 1., logSigma = 0.) for x in (-Inf, Inf)
604!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
605
606
607! Regular adaptive global quadrature.
608
609if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
610getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
6111.00000000000000000000000000000000000, 1.00000000000000000000000000000013597, 0.148409338409973194577277535014189193E-22, 0.135970654073738851225752042746495869E-30 (unbiased)? TRUE
612numFuncEval
613+651
614
615
616! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
617
618if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
619getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
6201.00000000000000000000000000000000000, 1.00000000000000000000000000000013597, 0.148409338409973194495421223865688864E-22, 0.135970654073738851225752042746495869E-30 (unbiased)? TRUE
621numFuncEval
622+651
623
624integrand%desc
625"intGenExpGammaPDF_type: f(x) = GenExpGamma(x; kappa = 5., omega = 10., logSigma = 2.) for x in (-Inf, Inf)"
626
627!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
628! intGenExpGammaPDF_type: f(x) = GenExpGamma(x; kappa = 5., omega = 10., logSigma = 2.) for x in (-Inf, Inf)
629!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
630
631
632! Regular adaptive global quadrature.
633
634if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
635getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
6361.00000000000000000000000000000000000, 1.00000000000000000000000000000000578, 0.132417436079502091744367222811402828E-22, 0.577778983316170755916793382775478196E-32 (unbiased)? TRUE
637numFuncEval
638+903
639
640
641! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
642
643if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
644getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
6451.00000000000000000000000000000000000, 1.00000000000000000000000000000000578, 0.132417436079839015749459180202287325E-22, 0.577778983316170755916793382775478196E-32 (unbiased)? TRUE
646numFuncEval
647+903
648
649integrand%desc
650"intPentaGammaInf_type: f(x) = sum of five Gamma PDFs with five break points in the integration range x in (-Inf, +Inf)"
651
652!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
653! intPentaGammaInf_type: f(x) = sum of five Gamma PDFs with five break points in the integration range x in (-Inf, +Inf)
654!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
655
656
657! Regular adaptive global quadrature.
658
659if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
660
661 ******** Integration did NOT converge. ********
662
663getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
6645.00000000000000000000000000000000000, 4.99999999999999999999972410297847334, 0.200738257888944093088270474559436884E-19, 0.275897021526661457740995046287992818E-21 (unbiased)? TRUE
665numFuncEval
666+21945
667
668
669! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
670
671if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
672getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
6735.00000000000000000000000000000000000, 4.99999999999999999999999459738972784, 0.120968281919312175939861440033309516E-21, 0.540261027215530879807578922336193197E-23 (unbiased)? TRUE
674numFuncEval
675+11025
676
677
678! Assisted adaptive global quadrature by specifying points of difficulties of the integrand.
679
680break = integrand%break
681break
682-9.00000000000000000000000000000000000, -5.00000000000000000000000000000000000, +2.00000000000000000000000000000000000, +5.00000000000000000000000000000000000, +7.00000000000000000000000000000000000
683if (isFailedQuad(getFunc, lb, ub, break, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
684getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
6855.00000000000000000000000000000000000, 4.99999999999999999999999999996475779, 0.454086489790053140200838733525734823E-22, 0.352422068663531514279007291757730680E-28 (unbiased)? TRUE
686numFuncEval
687+5502
688
689
690!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
691! Compute the Cauchy Principal Value over a finite interval.
692!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
693
694integrand%desc
695"intCauchy1_type: an integrand of the form w(x) * f(x) with Cauchy weight w(x) 1 / (x - cs) and f(x) = 1 ~,~ x \in (lb < cs, cs < ub)"
696
697!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
698! intCauchy1_type: an integrand of the form w(x) * f(x) with Cauchy weight w(x) 1 / (x - cs) and f(x) = 1 ~,~ x \in (lb < cs, cs < ub)
699!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
700
701
702! Regular adaptive global quadrature.
703
704if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
705
706 ******** Integration did NOT converge. ********
707
708getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
709-0.405465108108164381978013115464349228, -2.04356971055286136189443548589605091, 7.24146052036062818462419495305868050, 1.63810460244469697991642237043170168 (unbiased)? TRUE
710numFuncEval
711+4389
712
713
714! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
715
716if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
717
718 ******** Integration did NOT converge. ********
719
720getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
721-0.405465108108164381978013115464349228, -0.405465108108164381978013115464718525, 0.948013123403926992764805620933964846E-25, 0.369297066836252474823483770490659813E-30 (unbiased)? TRUE
722numFuncEval
723+1155
724
725
726! Assisted adaptive global quadrature by specifying Cauchy weight of the integrand.
727
728integrand%wcauchy%cs
729+1.00000000000000000000000000000000000
730if (isFailedQuad(getFuncUnweighted, lb, ub, integrand%wcauchy, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
731getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
732-0.405465108108164381978013115464349228, -0.405465108108164381978013115464349036, 0.192592994438723585305597794258492732E-33, 0.192592994438723585305597794258492732E-33 (unbiased)? TRUE
733numFuncEval
734+25
735
736integrand%desc
737"intCauchy2_type: an integrand of the form w(x) * f(x) = 1 / (x + 2.) / (x - 3.) for x in (-3., 2.)"
738
739!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
740! intCauchy2_type: an integrand of the form w(x) * f(x) = 1 / (x + 2.) / (x - 3.) for x in (-3., 2.)
741!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
742
743
744! Regular adaptive global quadrature.
745
746if (isFailedQuad(getFunc, lb, ub, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
747
748 ******** Integration did NOT converge. ********
749
750getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
751-0.635610766069589123929388320259411106, -0.313635651749756707119425316680125425, 1.43562881981449143882504794728653286, 0.321975114319832416809963003579285681 (unbiased)? TRUE
752numFuncEval
753+4347
754
755
756! Assisted adaptive global quadrature by the Wynn Epsilon extrapolation.
757
758if (isFailedQuad(getFunc, lb, ub, weps, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
759
760 ******** Integration did NOT converge. ********
761
762getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
763-0.635610766069589123929388320259411106, -0.635610766069589123929388320259460313, 0.597574167278168521319843806854314588E-24, 0.492075100790938760455802364330448930E-31 (unbiased)? TRUE
764numFuncEval
765+1281
766
767
768! Assisted adaptive global quadrature by specifying Cauchy weight of the integrand.
769
770integrand%wcauchy%cs
771-2.00000000000000000000000000000000000
772if (isFailedQuad(getFuncUnweighted, lb, ub, integrand%wcauchy, integral, abserr, neval = numFuncEval)) call disp%show(' ******** Integration did NOT converge. ********', tmsize = 1_IK, bmsize = 1_IK)
773getStr([truth, integral, abserr, abs(integral - truth)])//SK_' (unbiased)? '//getStr(abs(integral - truth) <= abserr)
774-0.635610766069589123929388320259411106, -0.635610766069589123929388320259411106, 0.485640976828224618892203332811773371E-24, 0.00000000000000000000000000000000000 (unbiased)? TRUE
775numFuncEval
776+771
777
778
type(weps_type), parameter weps
The scalar constant object of type weps_type that indicates the use of Epsilon extrapolation method o...
Test:
test_pm_quadPack


Final Remarks


If you believe this algorithm or its documentation can be improved, we appreciate your contribution and help to edit this page's documentation and source file on GitHub.
For details on the naming abbreviations, see this page.
For details on the naming conventions, see this page.
This software is distributed under the MIT license with additional terms outlined below.

  1. If you use any parts or concepts from this library to any extent, please acknowledge the usage by citing the relevant publications of the ParaMonte library.
  2. If you regenerate any parts/ideas from this library in a programming environment other than those currently supported by this ParaMonte library (i.e., other than C, C++, Fortran, MATLAB, Python, R), please also ask the end users to cite this original ParaMonte library.

This software is available to the public under a highly permissive license.
Help us justify its continued development and maintenance by acknowledging its benefit to society, distributing it, and contributing to it.

Author:
Amir Shahmoradi, Oct 16, 2009, 11:14 AM, Michigan

Definition at line 4989 of file pm_quadPack.F90.


The documentation for this interface was generated from the following file: