ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_quadPack::getQuadErr 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 Gauss-Kronrod (GK) extension rules. 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 Gauss-Kronrod (GK) extension rules.

The resulting integral hopefully satisfies the accuracy condition abs(truth - result) <= max(abstol, reltol * abs(i)),
where truth is the true value of the integral and abstol, reltol are the user-specified absolute and relative tolerances, respectively.
This interface combines and significantly extends the functionalities of the QAG, QAGI, QAGP, QAGS, QAWC, QAWO, QAWS subroutines of the venerable QUADPACK library.

Specifically,

  1. When the input argument help is missing, this interface performs adaptive global Gauss-Kronrod (GK) with the user-specified predefined (7-15 (GK15), 10-21 (GK21), 15-31 (GK31), 20-41 (GK41), 25-51 (GK51), 30-61 (GK61)) Gauss-Kronrod extension rules,
    or with user-specified arbitrarily-defined Gauss-Kronrod extension rules via the triple (nodeK, weightK, weightG) input arguments.
    This interface replicates and extends the QAG/QAGE functionalities of the QuadPack library to all Gauss-Kronrod extension rules mentioned above and all real kinds supported by the processor.
  2. When the input argument help is present and set to weps, this interface accelerates the convergence of adaptive global Gauss-Kronrod (GK) with predefined or arbitrary Gauss-Kronrod rules are mentioned above, via the Epsilon extrapolation algorithm of Wynn (1961)).
    This extrapolation is particularly helpful with,
    1. integration of integrands with singularities at or within the integration limits, or
    2. integration of integrands with semi or doubly infinite integration limits,
    possibly leading to orders of magnitude speedup in computing the integral within the requested tolerance.
    This interface replicates and extends the QAGS/QAGSE functionalities of the QuadPack library to all Gauss-Kronrod extension rules and all real kinds supported by the processor.
  3. When the input argument help is present and set to a vector of points of difficulties within the integration domain, this interface accelerates the convergence of adaptive global Gauss-Kronrod (GK) with predefined or arbitrary Gauss-Kronrod rules are mentioned above via the Epsilon algorithm of Wynn (1961)).
    Simultaneously, the interface carefully takes into account the specified points of integration difficulty (such as singularities or discontinuities) specified in help.
    This interface replicates and extends the QAGP/QAGPE functionalities of the QuadPack library to all Gauss-Kronrod extension rules and all real kinds supported by the processor.
  4. When the integration limits are set to,
    1. (lb = ninf, ub), implying the semi-infinite range \((-\infty, \ms{ub})\), or
    2. (lb, ub = pinf), implying the semi-infinite range \((\ms{lb}, +\infty)\), or
    3. (lb = ninf, ub = pinf), implying the fully-infinite range \((-\infty, +\infty)\),
    this interface first applies the following corresponding change of variables,
    1. \(x = \ms{ub} - \frac{1 - t}{t}\),
    2. \(x = \ms{lb} + \frac{1 - t}{t}\),
    3. \(x = \frac{1 - t}{t}\),
    to the integrand to transform the semi/fully-infinite integration range to the finite range \((0, 1)\).
    Then it calls the user-specified interface to perform the integration.
    The functionality of this interface significantly extends the functionalities of the QAGI/QAGIE routines of the QuadPack library.
    Specifically, this interface
    1. computes integrals with arbitrary machine precision by accepting all real kinds supported by the processor.
    2. allows the use of all predefined or arbitrary Gauss-Kronrod rules for semi or fully infinite integrations.
    3. allows explicit specification of singularities or discontinuities via the input help argument for semi or fully infinite integrations.
Parameters
getFunc: The input function to be integrated (i.e., the integrand).
  1. On entry, it must take an input scalar of the same type and kind as integral.
  2. On exit, it must generate an input scalar of the same type and kind as integral, representing the corresponding function value.
The following illustrates the general interface of getFunc:
function getFunc(x) result(func)
real(RKG) , intent(in) :: x
real(RKG) :: func
end function
where RKG must match the user-specified kind type parameter for the integral output argument below.
[in]lb: The input scalar argument that can be either,
  1. a value of type real of the same kind as integral, representing the lower limit of integration, or
  2. the constant ninf, representing negative infinity ( \(-\infty\)) as the lower limit of integration.
[in]ub: The input scalar argument that can be either,
  1. a value of type real of the same kind as integral, representing the upper limit of integration, or
  2. the constant pinf, representing positive infinity ( \(+\infty\)) as the upper limit of integration.
[in]abstol: The input scalar argument of the same type and kind as integral, representing the absolute tolerance of integration.
If the estimated integration error reaches a value below this threshold, the integration is assumed to have converged.
This argument can be set to any non-negative value, including abstol = 0..
[in]reltol: The input scalar argument of the same type and kind as integral, representing the relative tolerance of integration.
If the relative accuracy of integration reaches a value below this threshold, the integration is assumed to have converged.
This argument can be set to any positive value significantly larger (e.g., \(\times10000\)) than epsilon(real(0., kind(reltol)).
A good rule of thumb is to set reltol = epsilon(real(0, kind(integral)))**(2./3.).
The integration result is generally orders of magnitude more precise than the specified reltol.
[in]qrule: The input scalar constant argument that can be either,
  1. GK15 of type GK15_type, or
  2. GK21 of type GK21_type, or
  3. GK31 of type GK31_type, or
  4. GK41 of type GK41_type, or
  5. GK51 of type GK51_type, or
  6. GK61 of type GK61_type.
The specified objects are empty and merely serve to differentiate the multitude of orders of Gauss-Kronrod quadrature rules.
For example, specifying GK15 dictates the use of 15-points Gauss-Kronrod quadrature rules for computing the integral and estimating its error.
(optional. It must be present if and only if nodeK, weightK and weightG optional input arguments are missing.)
[in]nodeK: The input contiguous vector argument of the same type and kind as integral, of size \(n + 1\), where \(n\) is the number of points in the Gauss rule to be used for the integration.
It contains the nodes of the \(n\)-points Gauss-Legendre quadrature rule and its \(n+1\)-points Kronrod extension rule.
The procedures under the generic interface setNodeWeightGK return this vector.
(optional. It must be present if and only if weightK and weightG optional input arguments are present and qrule is missing.)
[in]weightK: The input contiguous vector argument of the same type and kind as integral, of size \(n + 1\), where \(n\) is the number of points in the Gauss rule to be used for the integration.
It contains the Kronrod optimal extension weights for the \(2n+1\)-points Gauss-Legendre-Kronrod quadrature method.
The procedures under the generic interface setNodeWeightGK return this vector.
(optional. It must be present if and only if nodeK and weightG optional input arguments are present and qrule is missing.)
[in]weightG: The input contiguous vector argument of the same type and kind as integral, of size \((n + 1)/2\), where \(n\) is the number of points in the Gauss rule to be used for the integration.
It contains the weights for the \(n\)-points Gauss-Legendre quadrature method.
The procedures under the generic interface setNodeWeightGK return this vector.
(optional. It must be present if and only if nodeK and weightK optional input arguments are present and qrule is missing.)
[in]help: The input scalar constant argument that can be any of the following:
  1. The scalar constant weps of an object of type weps_type, implying the use of the Epsilon extrapolation algorithm of Wynn (1961)) for computing the integral and estimating its error.
    The specified objects are empty and merely serve to differentiate the multitude of different extrapolation methods to accelerate the integration convergence.
    Specify this argument as weps if you suspect the integrand contains integrable singularities or the integration range is semi/fully-infinite.
    In such cases, the extrapolation method can lead to orders of magnitude faster convergence of integration compared with no acceleration.
    When the input argument help is set to wesp, the integration routine corresponds to the QAGI and QAGS routines of QuadPack software.
    However, unlike QAGI and QAGS which use Gauss-Kronrod 7-15 and 10-21 quadrature rules respectively, this generic interface allows the extrapolation and acceleration of any Gauss-Kronrod quadrature rule with the Epsilon method of Wynn (1961).
    See examples below for an illustration of the speedup and appropriate use cases.
  2. A vector of non-zero length of the same type and kind as the input integral, containing the points of integration difficulty within the specified integration range (lb, ub).
    These points can be (integrable) singularities or discontinuities in the integrand that cannot not be modeled well by polynomials.
    Specifying the difficulty points, when they are known, can lead to significant computational speedups and convergence of the integral.
    On input, the specified help points must be in ascending order.
    If manual sorting is impossible, help can be sorted automatically via setSorted before calling this procedure.
    Additionally, all elements of help must be unique.
    If uniqueness is unknown at runtime, the unique elements of help can be automatically obtained via getUnique before calling this procedure.
    **When both the lower and upper limits of integration are infinities ( \((\ms{lb}, \ms{ub}) = (-\infty, +\infty)\)), the specified difficulty points in help cannot contain -1**.
    In such cases, use a simple change of variable to shift the singularity to values other than -1.
  3. A scalar constant of type wcauchy_type implying that the Cauchy Principal Value of the input function must be computed.
    In such a case, the Cauchy singularity must be stored in the cs component of the input object of wcauchy_type before calling this generic interface.
    If the integrand has more than one Cauchy pole, split the integration range such that each separate range contains only one Cauchy singularity.
    Then, pass each integrand separately to this generic interface and sum the results of the two integrations.
(optional. If missing, the procedure attempts to compute the integral without convergence-accelerating extrapolations or explicit assumptions about the singularities or discontinuities.)
[out]integral: The output scalar of type real of kind any supported by the processor (e.g., RK, RK32, RK64, or RK128), containing the integral of the specified integrand getFunc() over the specified range (lb, ub).
[out]abserr: The output scalar argument of the same type and kind as integral, representing the estimated absolute error in the resulting integral,
which should equal or exceed abs(truth - integral) where truth is the true value of the integral.
[out]sinfo: The output contiguous matrix argument of the same type and kind as integral, of size (1:4, 1:nintmax) where nintmax is the maximum number of subintervals the procedure is allowed to create to refine the computed integral.
  1. The first row, sinfo(1,1:nint), contains the lower bound of each of nint <= nintmax final surviving subintervals formed in the integration process.
  2. The second row, sinfo(2,1:nint), contains the upper bound of each of nint <= nintmax final surviving subintervals formed in the integration process.
  3. The third row, sinfo(3,1:nint), contains the estimated integral of each of nint <= nintmax final surviving subintervals formed in the integration process.
  4. The fourth row, sinfo(4,1:nint), contains the estimated integral error of each of nint <= nintmax final surviving subintervals formed in the integration process.
[out]sindex: The output contiguous vector argument of type integer of default kind IK of size (1:nintmax), the first k elements of which contain pointers to the error estimates over the subintervals,
such that sinfo(3, sindex(1)), ..., sinfo(sindex(k)) form a decreasing sequence with k = nint if nint <= (nintmax / 2 + 2), and k = nintmax - nint + 1 otherwise.
[out]neval: The output scalar argument of type integer of default kind IK, containing the number of function calls made during the integration.
[out]nint: The output scalar argument of type integer of default kind IK, containing the number of final surviving subintervals formed during the integration.
Returns
err : The output scalar of type integer of default kind IK, that is set to zero if the integration converges without any errors.
Otherwise, a non-zero value of err indicates the occurrence of an error of varying severities,
  1. If err == 1, then the maximum number of subdivisions allowed has reached.
    One can allow more subdivisions by increasing the value of limit (and taking the according dimension adjustments into account).
    However, if this yields no improvement, it is advised to analyze the integrand in order to determine the integration difficulties.
    If the position of a local difficulty can be determined (i.e., a singularity or discontinuity within the interval), one will probably gain from splitting up the interval at this point and calling the integrator on the subranges.
    If possible, an appropriate special-purpose integrator should be used which is designed for handling the type of difficulty involved.
  2. If err == 2, then the occurrence of roundoff error is detected, which prevents the requested tolerance from being achieved.
  3. If err == 3, then an extremely bad integrand behavior occurs at some points of the integration interval.
  4. If err == 4, then the algorithm has not converged on return due to the detection of roundoff errors in the extrapolation method.
    The returned result is the best that could be obtained and is not necessarily wrong.
    In such cases, check whether the output abserr is satisfactory for the specific task, even though it may be larger than the requested tolerances.
    Increasing the requested absolute and relative tolerances may help convergence.
  5. If err == 5, then the algorithm has not converged.
    The integral is likely divergent, or slowly convergent.
    Note that divergence can be the underlying cause for any of the error numbers listed here.


Possible calling interfaces

err = getQuadErr(getFunc, lb, ub, abstol, reltol, qrule, integral, abserr, sinfo, sindex, neval, nint) ! extends QAG/QAGI routines of QuadPack.
err = getQuadErr(getFunc, lb, ub, abstol, reltol, qrule, help, integral, abserr, sinfo, sindex, neval, nint) ! extends QAGS/QAGI/QAWC routines of QuadPack.
err = getQuadErr(getFunc, lb, ub, abstol, reltol, qrule, help(:), integral, abserr, sinfo, sindex, neval, nint) ! extends QAGS/QAGI/QAGP routines of QuadPack.
err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK(:), weightK(:), weightG(:), integral, abserr, sinfo, sindex, neval, nint) ! extends QAG/QAGI/QAWC routines of QuadPack.
err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK(:), weightK(:), weightG(:), help, integral, abserr, sinfo, sindex, neval, nint) ! extends QAGS/QAGI routines of QuadPack.
err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK(:), weightK(:), weightG(:), help(:), integral, abserr, sinfo, sindex, neval, nint) ! extends QAGS/QAGI/QAGP routines of QuadPack.
!
Compute the 1D integral of the input scalar (potentially singular) integrand getFunc on a finite or s...
This module contains classes and procedures for non-adaptive and adaptive global numerical quadrature...
Warning
The condition lb < ub must hold for the corresponding procedure argument.
The condition abstol > 0. must hold for the corresponding procedure arguments.
The condition reltol > 0. must hold for the corresponding procedure arguments.
The condition size(sindex) > 0_IK must hold for the corresponding procedure arguments.
The condition size(sinfo, 1) == 4_IK must hold for the corresponding procedure arguments.
The condition size(sinfo, 2) == size(sindex) must hold for the corresponding procedure arguments.
The condition size(help) > 0 must hold for the corresponding procedure arguments when help contains a set of points of difficulties within the domain of integration.
The condition all(help > lb) must hold for the corresponding procedure arguments when help contains a set of points of difficulties within the domain of integration.
The condition all(help < ub) must hold for the corresponding procedure arguments when help contains a set of points of difficulties within the domain of integration.
The condition all(help /= -1.) .and. lb == ninf .and. ub == pinf must hold for the corresponding procedure arguments when help contains a set of points of difficulties within the domain of integration.
The element values of the input help must be in ascending order when help contains a set of points of difficulties within the domain of integration.
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.
Note
semi/fully-infinite intervals
If the integration range is semi or fully infinite, the recommended quadrature rule is GK21.
Higher quadrature rules are also possible and available for semi/fully-infinite intervals although benefits could be minimal due to the introduction of the singularity in the transformed finite-range 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.
5 use pm_quadTest, only: int1_type
6 use pm_quadTest, only: int2_type
7 use pm_quadTest, only: int3_type
8 use pm_quadTest, only: int4_type
9 use pm_quadTest, only: int5_type
10 use pm_quadTest, only: int6_type
11 use pm_quadTest, only: int7_type
12 use pm_quadTest, only: int8_type
13 use pm_quadTest, only: int9_type
23 use pm_io, only: display_type
24
25 implicit none
26
27 type(display_type) :: disp
28 disp = display_type(file = "main.out.F90")
29
30 call disp%skip()
31 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
32 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
33 call disp%show("! Compare integration via the Adaptive Global Gauss-Kronrod Quadrature method with (QAG) and without (QAGS) extrapolation acceleration.")
34 call disp%show("! Note the significant efficiency improvement in the integration via QAGS when the integrand has integrable singularity at some points.")
35 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
36 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
37 call disp%skip()
38
39 call test_getQuadErr(disp, int1_type(-4._RKH, +4._RKH))
40 call test_getQuadErr(disp, int2_type(+2._RKH, +3._RKH))
41 call test_getQuadErr(disp, int3_type(ub = +3._RKH))
43 call test_getQuadErr(disp, intSinCos_type(-4_IK, +4_IK, a = 10._RKH, b = -10._RKH))
44 call test_getQuadErr(disp, intNormPDF_type(-3._RKH, +3._RKH))
45 call test_getQuadErr(disp, intLogNormPDF_type(lb = exp(-6._RKG), ub = exp(+6._RKG)))
46
47 call disp%skip()
48 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
49 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
50 call disp%show("! Compare integration via:")
51 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with automatic break point handling and extrapolation acceleration (QAGS),")
52 call disp%show("! with,")
53 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with specified break point handling and extrapolation acceleration (QAGP).")
54 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
55 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
56 call disp%skip()
57
58 call test_getQuadErr(disp, int5_type(0._RKG, 3._RKG))
59 call test_getQuadErr(disp, int5_type(-2._RKG, +5._RKG))
60
61 call disp%skip()
62 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
63 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
64 call disp%show("! Compare integration via:")
65 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a semi-infinite range `(lb,+Inf)`,")
66 call disp%show("! with,")
67 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a semi-infinite range `(lb,+Inf)`.")
68 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
69 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
70 call disp%skip()
71
76
77 call disp%skip()
78 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
79 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
80 call disp%show("! Compare integration via:")
81 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a semi-infinite range `(-Inf, ub)`,")
82 call disp%show("! with,")
83 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a semi-infinite range `(-Inf, ub)`.")
84 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
85 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
86 call disp%skip()
87
88 call test_getQuadErr(disp, int8_type())
90 call test_getQuadErr(disp, intDoncker2_type(ub = -1._RKG))
91
92 call disp%skip()
93 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
94 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
95 call disp%show("! Compare integration via:")
96 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method without extrapolation on a full-infinite range `(-Inf, +Inf)`,")
97 call disp%show("! with,")
98 call disp%show("! the Adaptive Global Gauss-Kronrod Quadrature method with extrapolation on a full-infinite range `(-Inf, +Inf)`.")
99 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
100 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
101 call disp%skip()
102
106
107 call disp%skip()
108 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
109 call disp%show("! Compute the Cauchy Principal Value over a finite interval.")
110 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
111 call disp%skip()
112
115
116end program example
Generate and return an IEEE-compliant negative infinity.
Definition: pm_except.F90:1719
Generate and return an IEEE-compliant positive infinity.
Definition: pm_except.F90:1189
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 procedures and generic interfaces for testing for exceptional cases at runtime.
Definition: pm_except.F90:46
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
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
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 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
9
10!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
11! int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)
12!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
13
14!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
15! Using arbitrary Gauss-Kronrod nodes and weights (here: 35-71).
16!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
17
18integrand%desc
19"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
20[abstol, reltol]
21+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
22[lb, ub]
23-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
24
25call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
26err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK, weightK, weightG, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
27if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
28getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
290.592319847946765693983261139952398973, 0.592319847946765693983261139952399069, 0.420679694161706175066609104098715408E-29, 0.962964972193617926527988971292463659E-34 >= 1 (unbiased)? TRUE
30[numFuncEval, numInterval]
31+213, +2
32
33integrand%desc
34"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
35[abstol, reltol]
36+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
37[lb, ub]
38-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
39
40call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
41err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK, weightK, weightG, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
42if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
43getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
440.592319847946765693983261139952398973, 0.592319847946765693983261139952399069, 0.420679694161706175066609104098715408E-29, 0.962964972193617926527988971292463659E-34 >= 1 (unbiased)? TRUE
45[numFuncEval, numInterval]
46+213, +2
47
48!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
49! Using predefined Gauss-Kronrod 30-61 nodes and weights.
50!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
51
52integrand%desc
53"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
54[abstol, reltol]
55+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
56[lb, ub]
57-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
58
59call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
60err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK61, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
61if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
62getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
630.592319847946765693983261139952398973, 0.592319847946765693983261139952395699, 0.749057713854147205376412430387292458E-25, 0.327408090545830095019516250239437644E-32 >= 1 (unbiased)? TRUE
64[numFuncEval, numInterval]
65+183, +2
66
67integrand%desc
68"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
69[abstol, reltol]
70+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
71[lb, ub]
72-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
73
74call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
75err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK61, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
76if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
77getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
780.592319847946765693983261139952398973, 0.592319847946765693983261139952395699, 0.749057713854147205376412430387292458E-25, 0.327408090545830095019516250239437644E-32 >= 1 (unbiased)? TRUE
79[numFuncEval, numInterval]
80+183, +2
81
82!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
83! Using predefined Gauss-Kronrod 25-51 nodes and weights.
84!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
85
86integrand%desc
87"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
88[abstol, reltol]
89+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
90[lb, ub]
91-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
92
93call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
94err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK51, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
95if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
96getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
970.592319847946765693983261139952398973, 0.592319847946765693983261139952399069, 0.282354400096319256248938688151661529E-30, 0.962964972193617926527988971292463659E-34 >= 1 (unbiased)? TRUE
98[numFuncEval, numInterval]
99+357, +4
100
101integrand%desc
102"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
103[abstol, reltol]
104+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
105[lb, ub]
106-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
107
108call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
109err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK51, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
110if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
111getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1120.592319847946765693983261139952398973, 0.592319847946765693983261139952399069, 0.282354400096319256248938688151661529E-30, 0.962964972193617926527988971292463659E-34 >= 1 (unbiased)? TRUE
113[numFuncEval, numInterval]
114+357, +4
115
116!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
117! Using predefined Gauss-Kronrod 20-41 nodes and weights.
118!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
119
120integrand%desc
121"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
122[abstol, reltol]
123+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
124[lb, ub]
125-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
126
127call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
128err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK41, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
129if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
130getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1310.592319847946765693983261139952398973, 0.592319847946765693983261139952369891, 0.130763857408511882147171505007517289E-23, 0.290815421602472613811452669330324025E-31 >= 1 (unbiased)? TRUE
132[numFuncEval, numInterval]
133+287, +4
134
135integrand%desc
136"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
137[abstol, reltol]
138+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
139[lb, ub]
140-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
141
142call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
143err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK41, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
144if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
145getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1460.592319847946765693983261139952398973, 0.592319847946765693983261139952369891, 0.130763857408511882147171505007517289E-23, 0.290815421602472613811452669330324025E-31 >= 1 (unbiased)? TRUE
147[numFuncEval, numInterval]
148+287, +4
149
150!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
151! Using predefined Gauss-Kronrod 15-31 nodes and weights.
152!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
153
154integrand%desc
155"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
156[abstol, reltol]
157+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
158[lb, ub]
159-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
160
161call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
162err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK31, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
163if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
164getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1650.592319847946765693983261139952398973, 0.592319847946765693983261139952398973, 0.908014202174603824005803987572878882E-26, 0.00000000000000000000000000000000000 >= 1 (unbiased)? TRUE
166[numFuncEval, numInterval]
167+341, +6
168
169integrand%desc
170"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
171[abstol, reltol]
172+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
173[lb, ub]
174-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
175
176call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
177err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK31, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
178if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
179getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1800.592319847946765693983261139952398973, 0.592319847946765693983261139952398973, 0.908014202174603824005803987572878882E-26, 0.00000000000000000000000000000000000 >= 1 (unbiased)? TRUE
181[numFuncEval, numInterval]
182+341, +6
183
184!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
185! Using predefined Gauss-Kronrod 10-21 nodes and weights.
186!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
187
188integrand%desc
189"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
190[abstol, reltol]
191+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
192[lb, ub]
193-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
194
195call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
196err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK21, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
197if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
198getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
1990.592319847946765693983261139952398973, 0.592319847946765693983261139950914562, 0.242973532478881507761272478019151270E-22, 0.148441050463646203374289499924733273E-29 >= 1 (unbiased)? TRUE
200[numFuncEval, numInterval]
201+399, +10
202
203integrand%desc
204"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
205[abstol, reltol]
206+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
207[lb, ub]
208-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
209
210call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
211err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK21, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
212if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
213getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
2140.592319847946765693983261139952398973, 0.592319847946765693983261139950914466, 0.242973532478881507761272477921271703E-22, 0.148450680113368139553554779814446198E-29 >= 1 (unbiased)? TRUE
215[numFuncEval, numInterval]
216+399, +10
217
218!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
219! Using predefined Gauss-Kronrod 7-15 nodes and weights.
220!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
221
222integrand%desc
223"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
224[abstol, reltol]
225+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
226[lb, ub]
227-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
228
229call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
230err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK15, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
231if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
232getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
2330.592319847946765693983261139952398973, 0.592319847946765693983261139952415054, 0.311644046752404524231338211060753717E-22, 0.160815150356334193730174158205841431E-31 >= 1 (unbiased)? TRUE
234[numFuncEval, numInterval]
235+765, +26
236
237integrand%desc
238"int1_type: an algebraic integrand of the form f(x) = x**2 / (x**2 + 1) / (x**2 + 4) for x in (lb, ub)"
239[abstol, reltol]
240+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
241[lb, ub]
242-4.00000000000000000000000000000000000, +4.00000000000000000000000000000000000
243
244call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
245err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK15, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
246if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
247getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
2480.592319847946765693983261139952398973, 0.592319847946765693983261139952415054, 0.311644046752404524231537458497527643E-22, 0.160815150356334193730174158205841431E-31 >= 1 (unbiased)? TRUE
249[numFuncEval, numInterval]
250+765, +26
251
252
253!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
254! 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
255!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
256
257!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
258! Using arbitrary Gauss-Kronrod nodes and weights (here: 35-71).
259!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
260
261integrand%desc
262"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"
263[abstol, reltol]
264+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
265[lb, ub]
266+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
267
268call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
269err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK, weightK, weightG, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
270if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
271
272 **** INTEGRATION FAILED: ERR = 3 ****
273
274getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
2750.942809041582063365867792482806465323, 0.942809041582063353517270922122486063, 0.275996032154812596404096004412307632E-15, 0.123505215606839792599902659423179129E-16 >= 1 (unbiased)? TRUE
276[numFuncEval, numInterval]
277+14271, +101
278
279integrand%desc
280"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"
281[abstol, reltol]
282+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
283[lb, ub]
284+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
285
286call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
287err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK, weightK, weightG, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
288if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
289getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
2900.942809041582063365867792482806465323, 0.942809041582063365867792482806360745, 0.780290516968488605865629463438283303E-30, 0.104577995980226906820939602282361553E-30 >= 1 (unbiased)? TRUE
291[numFuncEval, numInterval]
292+781, +6
293
294!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
295! Using predefined Gauss-Kronrod 30-61 nodes and weights.
296!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
297
298integrand%desc
299"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"
300[abstol, reltol]
301+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
302[lb, ub]
303+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
304
305call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
306err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK61, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
307if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
308
309 **** INTEGRATION FAILED: ERR = 3 ****
310
311getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3120.942809041582063365867792482806465323, 0.942809041582063354181389418572630975, 0.395526291759776583849399257993266119E-15, 0.116864030642338343481331615363661636E-16 >= 1 (unbiased)? TRUE
313[numFuncEval, numInterval]
314+12261, +101
315
316integrand%desc
317"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"
318[abstol, reltol]
319+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
320[lb, ub]
321+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
322
323call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
324err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK61, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
325if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
326getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3270.942809041582063365867792482806465323, 0.942809041582063365867792482806368256, 0.683127351274152557078955376234873720E-30, 0.970668691971166869940212883062803369E-31 >= 1 (unbiased)? TRUE
328[numFuncEval, numInterval]
329+671, +6
330
331!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
332! Using predefined Gauss-Kronrod 25-51 nodes and weights.
333!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
334
335integrand%desc
336"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"
337[abstol, reltol]
338+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
339[lb, ub]
340+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
341
342call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
343err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK51, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
344if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
345
346 **** INTEGRATION FAILED: ERR = 3 ****
347
348getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3490.942809041582063365867792482806465323, 0.942809041582063353612391313951514454, 0.393612625583696433467216358595637327E-15, 0.122554011688549508692539077120199527E-16 >= 1 (unbiased)? TRUE
350[numFuncEval, numInterval]
351+10251, +101
352
353integrand%desc
354"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"
355[abstol, reltol]
356+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
357[lb, ub]
358+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
359
360call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
361err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK51, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
362if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
363getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3640.942809041582063365867792482806465323, 0.942809041582063365867792482806399167, 0.654527291600002104661074103787487549E-30, 0.661556935897015515524728423277922534E-31 >= 1 (unbiased)? TRUE
365[numFuncEval, numInterval]
366+561, +6
367
368!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
369! Using predefined Gauss-Kronrod 20-41 nodes and weights.
370!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
371
372integrand%desc
373"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"
374[abstol, reltol]
375+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
376[lb, ub]
377+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
378
379call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
380err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK41, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
381if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
382
383 **** INTEGRATION FAILED: ERR = 3 ****
384
385getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
3860.942809041582063365867792482806465323, 0.942809041582063354504151452306557739, 0.393427741808206315477620382019741615E-15, 0.113636410304999075837096309795641304E-16 >= 1 (unbiased)? TRUE
387[numFuncEval, numInterval]
388+8241, +101
389
390integrand%desc
391"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"
392[abstol, reltol]
393+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
394[lb, ub]
395+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
396
397call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
398err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK41, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
399if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
400getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
4010.942809041582063365867792482806465323, 0.942809041582063365867792482806474856, 0.311326575510196675646498834418853501E-30, 0.953335322471681747262709081579539023E-32 >= 1 (unbiased)? TRUE
402[numFuncEval, numInterval]
403+451, +6
404
405!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
406! Using predefined Gauss-Kronrod 15-31 nodes and weights.
407!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
408
409integrand%desc
410"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"
411[abstol, reltol]
412+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
413[lb, ub]
414+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
415
416call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
417err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK31, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
418if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
419
420 **** INTEGRATION FAILED: ERR = 3 ****
421
422getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
4230.942809041582063365867792482806465323, 0.942809041582063360076317716910846285, 0.279164711945302751623913154970201256E-15, 0.579147476589561903767056045130549251E-17 >= 1 (unbiased)? TRUE
424[numFuncEval, numInterval]
425+6293, +102
426
427integrand%desc
428"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"
429[abstol, reltol]
430+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
431[lb, ub]
432+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
433
434call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
435err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK31, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
436if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
437getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
4380.942809041582063365867792482806465323, 0.942809041582063365867792482806439901, 0.205400428568898703728420047576682499E-30, 0.254222752659115132603389088421210406E-31 >= 1 (unbiased)? TRUE
439[numFuncEval, numInterval]
440+341, +6
441
442!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
443! Using predefined Gauss-Kronrod 10-21 nodes and weights.
444!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
445
446integrand%desc
447"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"
448[abstol, reltol]
449+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
450[lb, ub]
451+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
452
453call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
454err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK21, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
455if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
456
457 **** INTEGRATION FAILED: ERR = 3 ****
458
459getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
4600.942809041582063365867792482806465323, 0.942809041582063358489159110539769991, 0.274189625236721766830515122145815584E-15, 0.737863337226669533225571680774701636E-17 >= 1 (unbiased)? TRUE
461[numFuncEval, numInterval]
462+4263, +102
463
464integrand%desc
465"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"
466[abstol, reltol]
467+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
468[lb, ub]
469+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
470
471call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
472err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK21, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
473if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
474getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
4750.942809041582063365867792482806465323, 0.942809041582063365867792482806451456, 0.747866619686254516023298199584791888E-27, 0.138666955995880981420030411866114767E-31 >= 1 (unbiased)? TRUE
476[numFuncEval, numInterval]
477+399, +10
478
479!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
480! Using predefined Gauss-Kronrod 7-15 nodes and weights.
481!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
482
483integrand%desc
484"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"
485[abstol, reltol]
486+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
487[lb, ub]
488+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
489
490call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
491err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK15, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
492if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
493
494 **** INTEGRATION FAILED: ERR = 3 ****
495
496getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
4970.942809041582063365867792482806465323, 0.942809041582063352985651916754259394, 0.107881013975510488055028051825168280E-14, 0.128821405660522059292977722684510948E-16 >= 1 (unbiased)? TRUE
498[numFuncEval, numInterval]
499+3465, +116
500
501integrand%desc
502"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"
503[abstol, reltol]
504+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
505[lb, ub]
506+0.00000000000000000000000000000000000, +0.666666666666666666666666666666666635
507
508call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
509err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK15, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
510if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
511getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
5120.942809041582063365867792482806465323, 0.942809041582063365867792482806467056, 0.894595922874628788043750076873935104E-26, 0.173333694994851226775038014832643459E-32 >= 1 (unbiased)? TRUE
513[numFuncEval, numInterval]
514+585, +20
515
516
517!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
518! 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
519!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
520
521!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
522! Using arbitrary Gauss-Kronrod nodes and weights (here: 35-71).
523!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
524
525integrand%desc
526"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"
527[abstol, reltol]
528+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
529[lb, ub]
530+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
531
532call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
533err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK, weightK, weightG, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
534if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
535getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
536-3.12249862669012531099145873075559453, -3.12249862669012531099145703288302519, 0.159554531774406733373916406327365902E-21, 0.169787256934451856030425956397998728E-23 >= 1 (unbiased)? TRUE
537[numFuncEval, numInterval]
538+22791, +161
539
540integrand%desc
541"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"
542[abstol, reltol]
543+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
544[lb, ub]
545+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
546
547call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
548err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK, weightK, weightG, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
549if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
550getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
551-3.12249862669012531099145873075559453, -3.12249862669012531099145873075586917, 0.104192809991349459650328406693844568E-29, 0.274637610069619832645782454612610636E-30 >= 1 (unbiased)? TRUE
552[numFuncEval, numInterval]
553+1065, +8
554
555!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
556! Using predefined Gauss-Kronrod 30-61 nodes and weights.
557!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
558
559integrand%desc
560"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"
561[abstol, reltol]
562+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
563[lb, ub]
564+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
565
566call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
567err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK61, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
568if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
569getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
570-3.12249862669012531099145873075559453, -3.12249862669012531099145675992893943, 0.158821823031096796049113138417448082E-21, 0.197082665510145750102784224886897408E-23 >= 1 (unbiased)? TRUE
571[numFuncEval, numInterval]
572+19581, +161
573
574integrand%desc
575"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"
576[abstol, reltol]
577+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
578[lb, ub]
579+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
580
581call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
582err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK61, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
583if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
584getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
585-3.12249862669012531099145873075559453, -3.12249862669012531099145873075560878, 0.295052467480124532688175820804010865E-30, 0.142518815884655453126142367751284622E-31 >= 1 (unbiased)? TRUE
586[numFuncEval, numInterval]
587+915, +8
588
589!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
590! Using predefined Gauss-Kronrod 25-51 nodes and weights.
591!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
592
593integrand%desc
594"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"
595[abstol, reltol]
596+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
597[lb, ub]
598+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
599
600call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
601err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK51, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
602if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
603getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
604-3.12249862669012531099145873075559453, -3.12249862669012531099145637935769671, 0.157863935952385241531199129294364333E-21, 0.235139789781688314990001714829882296E-23 >= 1 (unbiased)? TRUE
605[numFuncEval, numInterval]
606+16371, +161
607
608integrand%desc
609"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"
610[abstol, reltol]
611+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
612[lb, ub]
613+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
614
615call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
616err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK51, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
617if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
618getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
619-3.12249862669012531099145873075559453, -3.12249862669012531099145873075562303, 0.278104283969516857181283214909263505E-30, 0.285037631769310906252284735502569243E-31 >= 1 (unbiased)? TRUE
620[numFuncEval, numInterval]
621+765, +8
622
623!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
624! Using predefined Gauss-Kronrod 20-41 nodes and weights.
625!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
626
627integrand%desc
628"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"
629[abstol, reltol]
630+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
631[lb, ub]
632+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
633
634call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
635err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK41, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
636if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
637getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
638-3.12249862669012531099145873075559453, -3.12249862669012531099145581781496689, 0.156352167288494027740417861647194691E-21, 0.291294062764122468285662650331711561E-23 >= 1 (unbiased)? TRUE
639[numFuncEval, numInterval]
640+13161, +161
641
642integrand%desc
643"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"
644[abstol, reltol]
645+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
646[lb, ub]
647+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
648
649call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
650err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK41, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
651if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
652getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
653-3.12249862669012531099145873075559453, -3.12249862669012531099145873075552443, 0.445660189131206376397153295914152182E-30, 0.701038499756953850512375971100913544E-31 >= 1 (unbiased)? TRUE
654[numFuncEval, numInterval]
655+615, +8
656
657!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
658! Using predefined Gauss-Kronrod 15-31 nodes and weights.
659!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
660
661integrand%desc
662"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"
663[abstol, reltol]
664+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
665[lb, ub]
666+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
667
668call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
669err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK31, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
670if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
671getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
672-3.12249862669012531099145873075559453, -3.12249862669012531099145489154805592, 0.153862855021951093878837768152750949E-21, 0.383920753861232441434641228203070917E-23 >= 1 (unbiased)? TRUE
673[numFuncEval, numInterval]
674+9951, +161
675
676integrand%desc
677"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"
678[abstol, reltol]
679+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
680[lb, ub]
681+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
682
683call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
684err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK31, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
685if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
686getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
687-3.12249862669012531099145873075559453, -3.12249862669012531099145873075556603, 0.182963344716787406040317904545568095E-30, 0.285037631769310906252284735502569243E-31 >= 1 (unbiased)? TRUE
688[numFuncEval, numInterval]
689+465, +8
690
691!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
692! Using predefined Gauss-Kronrod 10-21 nodes and weights.
693!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
694
695integrand%desc
696"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"
697[abstol, reltol]
698+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
699[lb, ub]
700+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
701
702call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
703err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK21, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
704if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
705getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
706-3.12249862669012531099145873075559453, -3.12249862669012531099145729306953740, 0.166762278786010609240106982579981836E-21, 0.143768605713088174426049560344490814E-23 >= 1 (unbiased)? TRUE
707[numFuncEval, numInterval]
708+7959, +190
709
710integrand%desc
711"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"
712[abstol, reltol]
713+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
714[lb, ub]
715+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
716
717call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
718err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK21, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
719if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
720getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
721-3.12249862669012531099145873075559453, -3.12249862669012531099145873075553329, 0.638742935999207969456768107027576605E-25, 0.612445722315141001271800985742006887E-31 >= 1 (unbiased)? TRUE
722[numFuncEval, numInterval]
723+567, +14
724
725!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
726! Using predefined Gauss-Kronrod 7-15 nodes and weights.
727!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
728
729integrand%desc
730"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"
731[abstol, reltol]
732+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
733[lb, ub]
734+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
735
736call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
737err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK15, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
738if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
739getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
740-3.12249862669012531099145873075559453, -3.12249862669012531099145821665001954, 0.168984966934533206681021840309021027E-21, 0.514105574994657236239111087685608621E-24 >= 1 (unbiased)? TRUE
741[numFuncEval, numInterval]
742+9345, +312
743
744integrand%desc
745"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"
746[abstol, reltol]
747+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
748[lb, ub]
749+0.00000000000000000000000000000000000, +3.00000000000000000000000000000000000
750
751call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
752err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK15, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
753if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
754getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
755-3.12249862669012531099145873075559453, -3.12249862669012531099145873075569198, 0.736077254850655834241784402267592237E-24, 0.974520551859941341646324838947973223E-31 >= 1 (unbiased)? TRUE
756[numFuncEval, numInterval]
757+825, +28
758
759
760!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
761! 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
762!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
763
764!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
765! Using arbitrary Gauss-Kronrod nodes and weights (here: 35-71).
766!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
767
768integrand%desc
769"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"
770[abstol, reltol]
771+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
772[lb, ub]
773+0.00000000000000000000000000000000000, +1.00000000000000000000000000000000000
774
775call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
776err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK, weightK, weightG, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
777if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
778getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
779-0.189275187882093321180367135892330336, -0.189275187882093321180367157318588076, 0.709016787220663753280987987987449705E-23, 0.214262577393144583924325489311691073E-25 >= 1 (unbiased)? TRUE
780[numFuncEval, numInterval]
781+7313, +52
782
783integrand%desc
784"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"
785[abstol, reltol]
786+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
787[lb, ub]
788+0.00000000000000000000000000000000000, +1.00000000000000000000000000000000000
789
790call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
791err = getQuadErr(getFunc, lb, ub, abstol, reltol, nodeK, weightK, weightG, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
792if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
793getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
794-0.189275187882093321180367135892330336, -0.189275187882093321180367124314328969, 0.875190015574302935869778694464445775E-23, 0.115780013668100956229775935139168619E-25 >= 1 (unbiased)? TRUE
795[numFuncEval, numInterval]
796+3053, +22
797
798!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
799! Using predefined Gauss-Kronrod 30-61 nodes and weights.
800!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
801
802integrand%desc
803"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"
804[abstol, reltol]
805+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
806[lb, ub]
807+0.00000000000000000000000000000000000, +1.00000000000000000000000000000000000
808
809call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
810err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK61, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
811if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
812getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
813-0.189275187882093321180367135892330336, -0.189275187882093321180367149919369956, 0.551012606285227548119424458552174544E-23, 0.140270396193317866081626615262131133E-25 >= 1 (unbiased)? TRUE
814[numFuncEval, numInterval]
815+6405, +53
816
817integrand%desc
818"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"
819[abstol, reltol]
820+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
821[lb, ub]
822+0.00000000000000000000000000000000000, +1.00000000000000000000000000000000000
823
824call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
825err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK61, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
826if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
827getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
828-0.189275187882093321180367135892330336, -0.189275187882093321180367089569735215, 0.151044926215718500490728968093616592E-23, 0.463225951208172227937270578607855780E-25 >= 1 (unbiased)? TRUE
829[numFuncEval, numInterval]
830+2623, +22
831
832!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
833! Using predefined Gauss-Kronrod 25-51 nodes and weights.
834!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
835
836integrand%desc
837"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"
838[abstol, reltol]
839+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
840[lb, ub]
841+0.00000000000000000000000000000000000, +1.00000000000000000000000000000000000
842
843call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
844err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK51, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
845if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
846getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
847-0.189275187882093321180367135892330336, -0.189275187882093321180367156609609617, 0.990374035889893040148759005915834788E-23, 0.207172792803981521991901235694034717E-25 >= 1 (unbiased)? TRUE
848[numFuncEval, numInterval]
849+5355, +53
850
851integrand%desc
852"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"
853[abstol, reltol]
854+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
855[lb, ub]
856+0.00000000000000000000000000000000000, +1.00000000000000000000000000000000000
857
858call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
859err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK51, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
860if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
861getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
862-0.189275187882093321180367135892330336, -0.189275187882093321180367136737390768, 0.206538524316503455750529471463185887E-23, 0.845060431740830626700484040037212854E-27 >= 1 (unbiased)? TRUE
863[numFuncEval, numInterval]
864+2295, +23
865
866!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
867! Using predefined Gauss-Kronrod 20-41 nodes and weights.
868!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
869
870integrand%desc
871"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"
872[abstol, reltol]
873+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
874[lb, ub]
875+0.00000000000000000000000000000000000, +1.00000000000000000000000000000000000
876
877call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
878err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK41, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
879if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
880getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
881-0.189275187882093321180367135892330336, -0.189275187882093321180367151548235028, 0.957392833849128474709628758212286240E-23, 0.156559046920400077057876054566418227E-25 >= 1 (unbiased)? TRUE
882[numFuncEval, numInterval]
883+4387, +54
884
885integrand%desc
886"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"
887[abstol, reltol]
888+0.00000000000000000000000000000000000, +0.559579552072376591733071579965869345E-22
889[lb, ub]
890+0.00000000000000000000000000000000000, +1.00000000000000000000000000000000000
891
892call setNodeWeightGK(order = 35_IK, nodeK = nodeK, weightK = weightK, weightG = weightG) ! Compute the 35-71 Gauss-Kronrod nodes and weights.
893err = getQuadErr(getFunc, lb, ub, abstol, reltol, GK41, weps, integral, abserr, sinfo, sindex, numFuncEval, numInterval)
894if (err /= 0_IK) call disp%show(getStrUpper(SK_' **** integration failed: err = '//getStr(err)//SK_' ****'), tmsize = 1_IK, bmsize = 1_IK)
895getStr([truth, integral, abserr, abs(integral - truth)])//SK_' >= 1 (unbiased)? '//getStr(abs(integral - truth) <= abserr)
896-0.1892751878820933211803671358