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

Compute the 1D integral of the input scalar integrand getFunc on the finite or infinite interval [a, b] and estimate its absolute error via the requested non-adaptive Gauss-Kronrod (GK) extension rule. More...

Detailed Description

Compute the 1D integral of the input scalar integrand getFunc on the finite or infinite interval [a, b] and estimate its absolute error via the requested non-adaptive Gauss-Kronrod (GK) extension rule.

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 abserr.
  2. On exit, it must generate an input scalar of the same type and kind as abserr, representing the corresponding function value.
    The following illustrates the general interface of getFunc:
    function getFunc(x) result(func)
    use pm_kind, only: RK => RKG
    real(RK) , intent(in) :: x
    real(RK) :: func
    end function
    This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
    Definition: pm_kind.F90:268
    integer, parameter RK
    The default real kind in the ParaMonte library: real64 in Fortran, c_double in C-Fortran Interoperati...
    Definition: pm_kind.F90:543
    where RKG refers to any desired real kind supported by the processor that is desired for the output abserr.
[in]lb: The input scalar argument that can be either,
  1. a value of type real of the same kind as the kind for the output abserr, 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 the kind for the output abserr, representing the upper limit of integration, or
  2. the constant pinf, representing positive infinity ( \(+\infty\)) as the upper limit of integration.
[out]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 can 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 abserr, 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 can 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 abserr, 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 can 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 abserr, 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 can be present if and only if nodeK and weightK optional input arguments are present and qrule is missing.)
[out]abserr: The output scalar argument of type real of kind any supported by the processor (e.g., RK, RK32, RK64, or RK128), representing the estimated absolute error in the resulting integral.
[out]intAbsFunc: The output scalar argument of the same type and kind as abserr, representing the integral of the absolute value of the input integrand function getFunc() over the specified range [lb, ub].
This output is not of primary use to the end users, but is required for computation of global error and the stopping rule in Global Adaptive Quadrature algorithms.
[out]smoothness: The output scalar argument of the same type and kind as abserr, representing a measure of the smoothness of the input integrand function getFunc() over the specified range [lb, ub].
This output is not of primary use to the end users, but is required for computation of global error and the stopping rule in Global Adaptive Quadrature algorithms.
Returns
quadGK : The output scalar of the same type and kind as the output abserr, containing the integral of the specified integrand getFunc() over the specified range [lb, ub].


Possible calling interfaces

quadGK = getQuadGK(getFunc, lb, ub, qrule, abserr, intAbsFunc, smoothness)
quadGK = getQuadGK(getFunc, lb, ub, nodeK, weightK, weightG, abserr, intAbsFunc, smoothness)
Compute the 1D integral of the input scalar integrand getFunc on the finite or infinite interval [a,...
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.
Also, the conditions size(nodeK) == size(weightK) and ssize(weightG) == size(nodeK) / 2 must hold for the corresponding procedure arguments, if present.
These conditions are verified only if the library is built with the preprocessor macro CHECK_ENABLED=1.
If the target function contains points of difficulties, singularities, or discontinuities, user must ensure the abscissas of the specified Gauss-Kronrod rule do not match such points.
Particularly, computing an integrand at its singularities can lead to undefined values that can lead to unexpected segmentation fault or propagation of NaN values within the computational flow or other strange errors that can be extremely difficult to debug.
A simple check can be added within the target integrand implementations to ensure no such difficulty point matches an input value at which the function must be evaluated.
Alternatively, one should consider using the adaptive integration routines isFailedQuad or getQuadErr while setting their input help arguments to the points of difficulties of the integrand.
Remarks
The procedures under discussion are impure.
See also
Robert Piessens, Maria Branders, 1974, A Note on the Optimal Addition of Abscissas to Quadrature Formulas of Gauss and Lobatto, Mathematics of Computation, 28, 125, 135-139.
Laurie, 1997, Calculation of gauss-kronrod quadrature rules.
Kronrod, 1965, Nodes and weights of quadrature formulas.
FORTRAN90 version by John Burkardt.


Example usage

1program example
2
3 use pm_kind, only: SK, IK, RK => RKH ! All other real kinds are also supported.
4 use pm_quadpack, only: getQuadGK, ninf, pinf
5 use pm_quadpack, only: setNodeWeightGK
6 use pm_quadpack, only: GK15, GK21
7 use pm_quadpack, only: GK31, GK41
8 use pm_quadpack, only: GK51, GK61
9 use pm_io, only: display_type
10
11 implicit none
12
13 real(RK) :: integral, abserr, intAbsFunc, smoothness, truth
14
15 type(display_type) :: disp
16 disp = display_type(file = "main.out.F90")
17
18 call disp%skip()
19 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
20 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
21 call disp%show("! Gauss-Kronrod quadrature over a finite interval.")
22 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
23 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
24 call disp%skip()
25
26 call disp%skip()
27 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
28 call disp%show("! Compute the Gauss-Kronrod quadrature of the probability density function of the Beta distribution.")
29 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
30 call disp%skip()
31
32 truth = 1._RK
33
34 call disp%skip()
35 call disp%show("integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK15, abserr, intAbsFunc, smoothness)")
36 integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK15, abserr, intAbsFunc, smoothness)
37 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
38 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
39 call disp%skip()
40
41 call disp%skip()
42 call disp%show("integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK21, abserr, intAbsFunc, smoothness)")
43 integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK21, abserr, intAbsFunc, smoothness)
44 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
45 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
46 call disp%skip()
47
48 call disp%skip()
49 call disp%show("integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK41, abserr, intAbsFunc, smoothness)")
50 integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK41, abserr, intAbsFunc, smoothness)
51 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
52 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
53 call disp%skip()
54
55 call disp%skip()
56 call disp%show("integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK51, abserr, intAbsFunc, smoothness)")
57 integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK51, abserr, intAbsFunc, smoothness)
58 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
59 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
60 call disp%skip()
61
62 call disp%skip()
63 call disp%show("integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK61, abserr, intAbsFunc, smoothness)")
64 integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK61, abserr, intAbsFunc, smoothness)
65 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
66 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
67 call disp%skip()
68
69 call disp%skip()
70 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
71 call disp%show("! Compute the Gauss-Kronrod quadrature of the probability density function of the Beta distribution for non-default node counts.")
72 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
73 call disp%skip()
74
75 block
76 integer(IK) , parameter :: MAX_NPG = 30_IK
77 integer(IK) :: npg
78 real(RK) :: nodeK(MAX_NPG+1), weightK(MAX_NPG+1), weightG((MAX_NPG+1)/2)
79 do npg = 1_IK, MAX_NPG
80 call disp%skip()
81 call disp%show("npg ! The number of points in the Gauss-Legendre quadrature rule.")
82 call disp%show( npg )
83 call disp%show("call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights")
84 call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2))
85 call disp%show("integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.")
86 integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness)
87 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
88 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
89 call disp%skip()
90 end do
91 end block
92
93 call disp%skip()
94 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
95 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
96 call disp%show("! Gauss-Kronrod quadrature over a positive semi-infinite interval.")
97 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
98 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
99 call disp%skip()
100
101 call disp%skip()
102 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
103 call disp%show("! Compute the Gauss-Kronrod quadrature of the probability density function of the LogNormal distribution on the semi-infinite range `(0,+infinity)`.")
104 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
105 call disp%skip()
106
107 truth = 1._RK
108
109 call disp%skip()
110 call disp%show("integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK15, abserr, intAbsFunc, smoothness)")
111 integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK15, abserr, intAbsFunc, smoothness)
112 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
113 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
114 call disp%skip()
115
116 call disp%skip()
117 call disp%show("integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK21, abserr, intAbsFunc, smoothness)")
118 integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK21, abserr, intAbsFunc, smoothness)
119 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
120 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
121 call disp%skip()
122
123 call disp%skip()
124 call disp%show("integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK41, abserr, intAbsFunc, smoothness)")
125 integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK41, abserr, intAbsFunc, smoothness)
126 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
127 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
128 call disp%skip()
129
130 call disp%skip()
131 call disp%show("integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK51, abserr, intAbsFunc, smoothness)")
132 integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK51, abserr, intAbsFunc, smoothness)
133 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
134 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
135 call disp%skip()
136
137 call disp%skip()
138 call disp%show("integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK61, abserr, intAbsFunc, smoothness)")
139 integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK61, abserr, intAbsFunc, smoothness)
140 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
141 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
142 call disp%skip()
143
144 call disp%skip()
145 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
146 call disp%show("! Compute the Gauss-Kronrod quadrature of the probability density function of the LogNormal distribution for non-default node counts.")
147 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
148 call disp%skip()
149
150 block
151 integer(IK) , parameter :: MAX_NPG = 30_IK
152 integer(IK) :: npg
153 real(RK) :: nodeK(MAX_NPG+1), weightK(MAX_NPG+1), weightG((MAX_NPG+1)/2)
154 do npg = 1_IK, MAX_NPG
155 call disp%skip()
156 call disp%show("npg ! The number of points in the Gauss-Legendre quadrature rule.")
157 call disp%show( npg )
158 call disp%show("call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights")
159 call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2))
160 call disp%show("integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.")
161 integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness)
162 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
163 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
164 call disp%skip()
165 end do
166 end block
167
168 call disp%skip()
169 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
170 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
171 call disp%show("! Gauss-Kronrod quadrature over a negative semi-infinite interval.")
172 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
173 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
174 call disp%skip()
175
176 call disp%skip()
177 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
178 call disp%show("! Compute the Gauss-Kronrod quadrature of the probability density function of the Exponential distribution on the semi-infinite range `(-infinity,0)`.")
179 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
180 call disp%skip()
181
182 truth = 1._RK
183
184 call disp%skip()
185 call disp%show("integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK15, abserr, intAbsFunc, smoothness)")
186 integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK15, abserr, intAbsFunc, smoothness)
187 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
188 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
189 call disp%skip()
190
191 call disp%skip()
192 call disp%show("integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK21, abserr, intAbsFunc, smoothness)")
193 integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK21, abserr, intAbsFunc, smoothness)
194 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
195 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
196 call disp%skip()
197
198 call disp%skip()
199 call disp%show("integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK41, abserr, intAbsFunc, smoothness)")
200 integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK41, abserr, intAbsFunc, smoothness)
201 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
202 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
203 call disp%skip()
204
205 call disp%skip()
206 call disp%show("integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK51, abserr, intAbsFunc, smoothness)")
207 integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK51, abserr, intAbsFunc, smoothness)
208 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
209 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
210 call disp%skip()
211
212 call disp%skip()
213 call disp%show("integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK61, abserr, intAbsFunc, smoothness)")
214 integral = getQuadGK(getNegExpPDF, ninf, 0._RK, GK61, abserr, intAbsFunc, smoothness)
215 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
216 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
217 call disp%skip()
218
219 call disp%skip()
220 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
221 call disp%show("! Compute the Gauss-Kronrod quadrature of the probability density function of the Exponential distribution for non-default node counts.")
222 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
223 call disp%skip()
224
225 block
226 integer(IK) , parameter :: MAX_NPG = 30_IK
227 integer(IK) :: npg
228 real(RK) :: nodeK(MAX_NPG+1), weightK(MAX_NPG+1), weightG((MAX_NPG+1)/2)
229 do npg = 1_IK, MAX_NPG
230 call disp%skip()
231 call disp%show("npg ! The number of points in the Gauss-Legendre quadrature rule.")
232 call disp%show( npg )
233 call disp%show("call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights")
234 call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2))
235 call disp%show("integral = getQuadGK(getNegExpPDF, ninf, 0._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.")
236 integral = getQuadGK(getNegExpPDF, ninf, 0._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness)
237 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
238 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
239 call disp%skip()
240 end do
241 end block
242
243 call disp%skip()
244 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
245 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
246 call disp%show("! Gauss-Kronrod quadrature over a full infinite interval (-infinity, +infinity).")
247 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
248 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
249 call disp%skip()
250
251 call disp%skip()
252 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
253 call disp%show("! Compute the Gauss-Kronrod quadrature of the probability density function of the Normal distribution on the infinite range `(-infinity,+infinity)`.")
254 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
255 call disp%skip()
256
257 truth = 1._RK
258
259 call disp%skip()
260 call disp%show("integral = getQuadGK(getNormPDF, ninf, pinf, GK15, abserr, intAbsFunc, smoothness)")
261 integral = getQuadGK(getNormPDF, ninf, pinf, GK15, abserr, intAbsFunc, smoothness)
262 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
263 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
264 call disp%skip()
265
266 call disp%skip()
267 call disp%show("integral = getQuadGK(getNormPDF, ninf, pinf, GK21, abserr, intAbsFunc, smoothness)")
268 integral = getQuadGK(getNormPDF, ninf, pinf, GK21, abserr, intAbsFunc, smoothness)
269 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
270 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
271 call disp%skip()
272
273 call disp%skip()
274 call disp%show("integral = getQuadGK(getNormPDF, ninf, pinf, GK41, abserr, intAbsFunc, smoothness)")
275 integral = getQuadGK(getNormPDF, ninf, pinf, GK41, abserr, intAbsFunc, smoothness)
276 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
277 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
278 call disp%skip()
279
280 call disp%skip()
281 call disp%show("integral = getQuadGK(getNormPDF, ninf, pinf, GK51, abserr, intAbsFunc, smoothness)")
282 integral = getQuadGK(getNormPDF, ninf, pinf, GK51, abserr, intAbsFunc, smoothness)
283 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
284 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
285 call disp%skip()
286
287 call disp%skip()
288 call disp%show("integral = getQuadGK(getNormPDF, ninf, pinf, GK61, abserr, intAbsFunc, smoothness)")
289 integral = getQuadGK(getNormPDF, ninf, pinf, GK61, abserr, intAbsFunc, smoothness)
290 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
291 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
292 call disp%skip()
293
294 call disp%skip()
295 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
296 call disp%show("! Compute the Gauss-Kronrod quadrature of the probability density function of the Normal distribution for non-default node counts.")
297 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
298 call disp%skip()
299
300 block
301 integer(IK) , parameter :: MAX_NPG = 30_IK
302 integer(IK) :: npg
303 real(RK) :: nodeK(MAX_NPG+1), weightK(MAX_NPG+1), weightG((MAX_NPG+1)/2)
304 do npg = 1_IK, MAX_NPG
305 call disp%skip()
306 call disp%show("npg ! The number of points in the Gauss-Legendre quadrature rule.")
307 call disp%show( npg )
308 call disp%show("call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights")
309 call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2))
310 call disp%show("integral = getQuadGK(getNormPDF, ninf, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.")
311 integral = getQuadGK(getNormPDF, ninf, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness)
312 call disp%show("[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]")
313 call disp%show( [truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)] )
314 call disp%skip()
315 end do
316 end block
317
318contains
319
320 function getBetaPDF(x) result(pdf)
321 use pm_distBeta, only: getBetaLogPDF
322 real(RK) , intent(in) :: x
323 real(RK) :: pdf
324 pdf = exp(getBetaLogPDF(x, alpha = 2._RK, beta = 2._RK))
325 end function
326
327 function getLogNormPDF(x) result(pdf)
329 real(RK) , intent(in) :: x
330 real(RK) :: pdf
331 pdf = exp(getLogNormLogPDF(x))
332 end function
333
334 function getNegExpPDF(x) result(pdf)
335 use pm_distExp, only: getExpLogPDF
336 real(RK) , intent(in) :: x
337 real(RK) :: pdf
338 pdf = exp(getExpLogPDF(-x))
339 end function
340
341 function getNormPDF(x) result(pdf)
342 use pm_distNorm, only: getNormLogPDF
343 real(RK) , intent(in) :: x
344 real(RK) :: pdf
345 pdf = exp(getNormLogPDF(x))
346 end function
347
348end program example
Generate and return the natural logarithm of the Probability Density Function (PDF) of the Beta distr...
Generate and return the natural logarithm of the Probability Density Function (PDF) of the Exponentia...
Definition: pm_distExp.F90:221
Generate the natural logarithm of probability density function (PDF) of the univariate Lognormal dist...
Generate the natural logarithm of probability density function (PDF) of the univariate Normal distrib...
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
Return the Kronrod nodes and weights of the extension to the -point Gauss-Legendre quadrature,...
This module contains classes and procedures for computing various statistical quantities related to t...
Definition: pm_distBeta.F90:99
This module contains classes and procedures for computing various statistical quantities related to t...
Definition: pm_distExp.F90:112
This module contains classes and procedures for computing various statistical quantities related to t...
This module contains classes and procedures for computing various statistical quantities related to t...
This module contains classes and procedures for input/output (IO) or generic display operations on st...
Definition: pm_io.F90:252
type(display_type) disp
This is a scalar module variable an object of type display_type for general display.
Definition: pm_io.F90:11393
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RKH
The scalar integer constant of intrinsic default kind, representing the highest-precision real kind t...
Definition: pm_kind.F90:858
type(GK31_type), parameter GK31
The scalar constant object of type GK31_type that indicates the use of 15-point Gauss-Legendre quadra...
type(GK51_type), parameter GK51
The scalar constant object of type GK51_type that indicates the use of 25-point Gauss-Legendre quadra...
type(GK61_type), parameter GK61
The scalar constant object of type GK61_type that indicates the use of 30-point Gauss-Legendre quadra...
type(GK21_type), parameter GK21
The scalar constant object of type GK21_type that indicates the use of 10-point Gauss-Legendre quadra...
type(GK15_type), parameter GK15
The scalar constant object of type GK15_type that indicates the use of 7-point Gauss-Legendre quadrat...
type(GK41_type), parameter GK41
The scalar constant object of type GK41_type that indicates the use of 20-point Gauss-Legendre quadra...
Generate and return an object of type display_type.
Definition: pm_io.F90:10282

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! Gauss-Kronrod quadrature over a finite interval.
5!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
6!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
7
8
9!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
10! Compute the Gauss-Kronrod quadrature of the probability density function of the Beta distribution.
11!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
12
13
14integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK15, abserr, intAbsFunc, smoothness)
15[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
16+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000000, +0.962964972193617926527988971292463659E-32, +1.00000000000000000000000000000000000, +0.332350865347540181865392442941994453, +0.00000000000000000000000000000000000
17
18
19integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK21, abserr, intAbsFunc, smoothness)
20[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
21+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999904, +0.962964972193617926527988971292463541E-32, +0.999999999999999999999999999999999904, +0.348820029290895558486501301340103054, +0.962964972193617926527988971292463659E-34
22
23
24integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK41, abserr, intAbsFunc, smoothness)
25[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
26+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999904, +0.962964972193617926527988971292463541E-32, +0.999999999999999999999999999999999904, +0.365617398082792305848065507834403699, +0.962964972193617926527988971292463659E-34
27
28
29integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK51, abserr, intAbsFunc, smoothness)
30[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
31+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999807, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999807, +0.369298545656546792050550199588929904, +0.192592994438723585305597794258492732E-33
32
33
34integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, GK61, abserr, intAbsFunc, smoothness)
35[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
36+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999904, +0.962964972193617926527988971292463541E-32, +0.999999999999999999999999999999999904, +0.372002569341007109923787616661336924, +0.962964972193617926527988971292463659E-34
37
38
39!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
40! Compute the Gauss-Kronrod quadrature of the probability density function of the Beta distribution for non-default node counts.
41!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
42
43
44npg ! The number of points in the Gauss-Legendre quadrature rule.
45+1
46call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
47integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
48[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
49+1.00000000000000000000000000000000000, +0.480421111969383347900839915541048870, +0.270837007627134129622102202568011808, +0.480421111969383347900839915541048870, +0.270837007627134129622102202568011808, +0.519578888030616652099160084458951130
50
51
52npg ! The number of points in the Gauss-Legendre quadrature rule.
53+2
54call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
55integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
56[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
57+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000019, +0.962964972193617926527988971292463897E-32, +1.00000000000000000000000000000000019, +0.205050505050505050505050505050505040, +0.192592994438723585305597794258492732E-33
58
59
60npg ! The number of points in the Gauss-Legendre quadrature rule.
61+3
62call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
63integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
64[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
65+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999807, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999807, +0.313218923737658024129334619196987163, +0.192592994438723585305597794258492732E-33
66
67
68npg ! The number of points in the Gauss-Legendre quadrature rule.
69+4
70call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
71integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
72[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
73+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999807, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999807, +0.315930942664627962287706737178712867, +0.192592994438723585305597794258492732E-33
74
75
76npg ! The number of points in the Gauss-Legendre quadrature rule.
77+5
78call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
79integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
80[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
81+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999904, +0.962964972193617926527988971292463541E-32, +0.999999999999999999999999999999999904, +0.321609174541580455838518696237276792, +0.962964972193617926527988971292463659E-34
82
83
84npg ! The number of points in the Gauss-Legendre quadrature rule.
85+6
86call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
87integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
88[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
89+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999711, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999711, +0.337728935330107303671417159619312911, +0.288889491658085377958396691387739098E-33
90
91
92npg ! The number of points in the Gauss-Legendre quadrature rule.
93+7
94call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
95integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
96[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
97+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999422, +0.962964972193617926527988971292463066E-32, +0.999999999999999999999999999999999422, +0.332350865347540181865392442941993586, +0.577778983316170755916793382775478196E-33
98
99
100npg ! The number of points in the Gauss-Legendre quadrature rule.
101+8
102call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
103integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
104[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
105+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000000, +0.962964972193617926527988971292463659E-32, +1.00000000000000000000000000000000000, +0.345422184218279383037014278751830416, +0.00000000000000000000000000000000000
106
107
108npg ! The number of points in the Gauss-Legendre quadrature rule.
109+9
110call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
111integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
112[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
113+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000019, +0.962964972193617926527988971292463897E-32, +1.00000000000000000000000000000000019, +0.348302550445394272753570649136973960, +0.192592994438723585305597794258492732E-33
114
115
116npg ! The number of points in the Gauss-Legendre quadrature rule.
117+10
118call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
119integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
120[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
121+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999807, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999807, +0.348820029290895558486501301340102573, +0.192592994438723585305597794258492732E-33
122
123
124npg ! The number of points in the Gauss-Legendre quadrature rule.
125+11
126call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
127integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
128[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
129+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000019, +0.962964972193617926527988971292463897E-32, +1.00000000000000000000000000000000019, +0.354999841253936686479701825785108681, +0.192592994438723585305597794258492732E-33
130
131
132npg ! The number of points in the Gauss-Legendre quadrature rule.
133+12
134call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
135integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
136[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
137+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000039, +0.962964972193617926527988971292464015E-32, +1.00000000000000000000000000000000039, +0.354157860946384654881202770232920338, +0.385185988877447170611195588516985464E-33
138
139
140npg ! The number of points in the Gauss-Legendre quadrature rule.
141+13
142call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
143integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
144[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
145+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000000, +0.962964972193617926527988971292463659E-32, +1.00000000000000000000000000000000000, +0.358244419047881521537644154374837196, +0.00000000000000000000000000000000000
146
147
148npg ! The number of points in the Gauss-Legendre quadrature rule.
149+14
150call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
151integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
152[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
153+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999807, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999807, +0.360043536617096848677395286636972682, +0.192592994438723585305597794258492732E-33
154
155
156npg ! The number of points in the Gauss-Legendre quadrature rule.
157+15
158call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
159integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
160[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
161+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999807, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999807, +0.359770623875105141649427242169437024, +0.192592994438723585305597794258492732E-33
162
163
164npg ! The number of points in the Gauss-Legendre quadrature rule.
165+16
166call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
167integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
168[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
169+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999904, +0.962964972193617926527988971292463541E-32, +0.999999999999999999999999999999999904, +0.363157479470576703962170302649563992, +0.962964972193617926527988971292463659E-34
170
171
172npg ! The number of points in the Gauss-Legendre quadrature rule.
173+17
174call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
175integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
176[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
177+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000000, +0.962964972193617926527988971292463659E-32, +1.00000000000000000000000000000000000, +0.363109511055304158303972009711228538, +0.00000000000000000000000000000000000
178
179
180npg ! The number of points in the Gauss-Legendre quadrature rule.
181+18
182call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
183integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
184[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
185+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999807, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999807, +0.364796933312968640518518780790230271, +0.192592994438723585305597794258492732E-33
186
187
188npg ! The number of points in the Gauss-Legendre quadrature rule.
189+19
190call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
191integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
192[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
193+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999807, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999807, +0.366088266354431031746529998266692271, +0.192592994438723585305597794258492732E-33
194
195
196npg ! The number of points in the Gauss-Legendre quadrature rule.
197+20
198call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
199integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
200[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
201+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000000, +0.962964972193617926527988971292463659E-32, +1.00000000000000000000000000000000000, +0.365617398082792305848065507834402592, +0.00000000000000000000000000000000000
202
203
204npg ! The number of points in the Gauss-Legendre quadrature rule.
205+21
206call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
207integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
208[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
209+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999711, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999711, +0.367821045564605980960034310589386864, +0.288889491658085377958396691387739098E-33
210
211
212npg ! The number of points in the Gauss-Legendre quadrature rule.
213+22
214call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
215integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
216[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
217+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000058, +0.962964972193617926527988971292464253E-32, +1.00000000000000000000000000000000058, +0.367991165630253018476809914215071252, +0.577778983316170755916793382775478196E-33
218
219
220npg ! The number of points in the Gauss-Legendre quadrature rule.
221+23
222call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
223integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
224[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
225+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999615, +0.962964972193617926527988971292463303E-32, +0.999999999999999999999999999999999615, +0.368768410353483515935716732597477937, +0.385185988877447170611195588516985464E-33
226
227
228npg ! The number of points in the Gauss-Legendre quadrature rule.
229+24
230call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
231integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
232[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
233+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999711, +0.962964972193617926527988971292463422E-32, +0.999999999999999999999999999999999711, +0.369761279627099167963894850246886530, +0.288889491658085377958396691387739098E-33
234
235
236npg ! The number of points in the Gauss-Legendre quadrature rule.
237+25
238call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
239integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
240[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
241+1.00000000000000000000000000000000000, +0.999999999999999999999999999999999615, +0.962964972193617926527988971292463303E-32, +0.999999999999999999999999999999999615, +0.369298545656546792050550199588930818, +0.385185988877447170611195588516985464E-33
242
243
244npg ! The number of points in the Gauss-Legendre quadrature rule.
245+26
246call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
247integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
248[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
249+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000000, +0.962964972193617926527988971292463659E-32, +1.00000000000000000000000000000000000, +0.370827437016519009582483758754977133, +0.00000000000000000000000000000000000
250
251
252npg ! The number of points in the Gauss-Legendre quadrature rule.
253+27
254call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
255integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
256[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
257+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000039, +0.962964972193617926527988971292464015E-32, +1.00000000000000000000000000000000039, +0.371074860800886503248166477335971387, +0.385185988877447170611195588516985464E-33
258
259
260npg ! The number of points in the Gauss-Legendre quadrature rule.
261+28
262call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
263integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
264[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
265+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000077, +0.962964972193617926527988971292464371E-32, +1.00000000000000000000000000000000077, +0.371427050294648513107276864275321340, +0.770371977754894341222391177033970927E-33
266
267
268npg ! The number of points in the Gauss-Legendre quadrature rule.
269+29
270call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
271integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
272[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
273+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000019, +0.962964972193617926527988971292463897E-32, +1.00000000000000000000000000000000019, +0.372223979361284634235934277787081296, +0.192592994438723585305597794258492732E-33
274
275
276npg ! The number of points in the Gauss-Legendre quadrature rule.
277+30
278call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
279integral = getQuadGK(getBetaPDF, 0._RK, +1._RK, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
280[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
281+1.00000000000000000000000000000000000, +1.00000000000000000000000000000000000, +0.962964972193617926527988971292463659E-32, +1.00000000000000000000000000000000000, +0.372002569341007109923787616661335335, +0.00000000000000000000000000000000000
282
283
284!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
285!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
286! Gauss-Kronrod quadrature over a positive semi-infinite interval.
287!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
288!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
289
290
291!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
292! Compute the Gauss-Kronrod quadrature of the probability density function of the LogNormal distribution on the semi-infinite range `(0,+infinity)`.
293!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
294
295
296integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK15, abserr, intAbsFunc, smoothness)
297[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
298+1.00000000000000000000000000000000000, +0.999997299690519869241720281259818052, +0.256665312152048430438215466789267068E-1, +0.999997299690519869241720281259818052, +0.405482311281071585135996295867852563, +0.270030948013075827971874018194840906E-5
299
300
301integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK21, abserr, intAbsFunc, smoothness)
302[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
303+1.00000000000000000000000000000000000, +1.00000009798761157060346918719858444, +0.295193558106936309656912906523764044E-2, +1.00000009798761157060346918719858444, +0.426989872574134944871374539334134092, +0.979876115706034691871985844378680066E-7
304
305
306integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK41, abserr, intAbsFunc, smoothness)
307[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
308+1.00000000000000000000000000000000000, +0.999999999795081386330293090324938485, +0.125441756876698348170916154866295533E-5, +0.999999999795081386330293090324938485, +0.446857633223355518536040613616829822, +0.204918613669706909675061515079233957E-9
309
310
311integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK51, abserr, intAbsFunc, smoothness)
312[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
313+1.00000000000000000000000000000000000, +0.999999999955606206754649785466436600, +0.391224786873945559603990633870625383E-7, +0.999999999955606206754649785466436600, +0.451151051357497766162941844671249309, +0.443937932453502145335634001947330452E-10
314
315
316integral = getQuadGK(getLogNormPDF, 0._RK, pinf, GK61, abserr, intAbsFunc, smoothness)
317[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
318+1.00000000000000000000000000000000000, +0.999999999995853818011903866479753249, +0.225157476189573037457067065221885697E-10, +0.999999999995853818011903866479753249, +0.454073139160562458663076636221205096, +0.414618198809613352024675133634889360E-11
319
320
321!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
322! Compute the Gauss-Kronrod quadrature of the probability density function of the LogNormal distribution for non-default node counts.
323!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
324
325
326npg ! The number of points in the Gauss-Legendre quadrature rule.
327+1
328call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
329integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
330[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
331+1.00000000000000000000000000000000000, +0.440654390519982570923254671966087876, +0.263157810538753065658119728890170155, +0.440654390519982570923254671966087876, +0.263157810538753065658119728890170155, +0.559345609480017429076745328033912172
332
333
334npg ! The number of points in the Gauss-Legendre quadrature rule.
335+2
336call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
337integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
338[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
339+1.00000000000000000000000000000000000, +1.00115159535561551677250223806146437, +0.248267677849831406984651517189212306, +1.00115159535561551677250223806146437, +0.248267677849831406984651517189212306, +0.115159535561551677250223806146437201E-2
340
341
342npg ! The number of points in the Gauss-Legendre quadrature rule.
343+3
344call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
345integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
346[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
347+1.00000000000000000000000000000000000, +1.00040478970647333319876225416181032, +0.386544757028673834265268612489689942, +1.00040478970647333319876225416181032, +0.386544757028673834265268612489689942, +0.404789706473333198762254161810319913E-3
348
349
350npg ! The number of points in the Gauss-Legendre quadrature rule.
351+4
352call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
353integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
354[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
355+1.00000000000000000000000000000000000, +1.00002610414724307748322615021718026, +0.386231062970553098636444648476170314, +1.00002610414724307748322615021718026, +0.386231062970553098636444648476170314, +0.261041472430774832261502171802586866E-4
356
357
358npg ! The number of points in the Gauss-Legendre quadrature rule.
359+5
360call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
361integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
362[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
363+1.00000000000000000000000000000000000, +0.999985810059247339669805714912328236, +0.394991742073170665588114892341666176, +0.999985810059247339669805714912328236, +0.394991742073170665588114892341666176, +0.141899407526603301942850876717643897E-4
364
365
366npg ! The number of points in the Gauss-Legendre quadrature rule.
367+6
368call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
369integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
370[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
371+1.00000000000000000000000000000000000, +0.999991153146646905609618866850569932, +0.201055388449807341789299677668162966E-2, +0.999991153146646905609618866850569932, +0.413788278103359200140993657807761858, +0.884685335309439038113314943006817385E-5
372
373
374npg ! The number of points in the Gauss-Legendre quadrature rule.
375+7
376call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
377integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
378[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
379+1.00000000000000000000000000000000000, +0.999997299690519869241720281259817474, +0.256665312152048430438215466789497668E-1, +0.999997299690519869241720281259817474, +0.405482311281071585135996295867851504, +0.270030948013075827971874018252618804E-5
380
381
382npg ! The number of points in the Gauss-Legendre quadrature rule.
383+8
384call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
385integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
386[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
387+1.00000000000000000000000000000000000, +0.999999381352227797311100131928905262, +0.183053043078258033477296544446724386E-1, +0.999999381352227797311100131928905262, +0.423097534055426140755284335767787928, +0.618647772202688899868071094737841243E-6
388
389
390npg ! The number of points in the Gauss-Legendre quadrature rule.
391+9
392call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
393integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
394[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
395+1.00000000000000000000000000000000000, +0.999999990612983549879343252846014940, +0.818358877620782547271452401299804686E-2, +0.999999990612983549879343252846014940, +0.425704502695006040033705707820210107, +0.938701645012065674715398505995817499E-8
396
397
398npg ! The number of points in the Gauss-Legendre quadrature rule.
399+10
400call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
401integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
402[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
403+1.00000000000000000000000000000000000, +1.00000009798761157060346918719858444, +0.295193558106936309656912906523764233E-2, +1.00000009798761157060346918719858444, +0.426989872574134944871374539334133514, +0.979876115706034691871985844378680066E-7
404
405
406npg ! The number of points in the Gauss-Legendre quadrature rule.
407+11
408call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
409integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
410[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
411+1.00000000000000000000000000000000000, +1.00000007745615047036466107235347550, +0.899333271967751714425000608569671298E-3, +1.00000007745615047036466107235347550, +0.434133277977773222279492925965094215, +0.774561504703646610723534755037693165E-7
412
413
414npg ! The number of points in the Gauss-Legendre quadrature rule.
415+12
416call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
417integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
418[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
419+1.00000000000000000000000000000000000, +1.00000004799397022376443202052945016, +0.227785342107320067761047547616597876E-3, +1.00000004799397022376443202052945016, +0.432341032766134749110172349751338347, +0.479939702237644320205294501602025633E-7
420
421
422npg ! The number of points in the Gauss-Legendre quadrature rule.
423+13
424call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
425integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
426[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
427+1.00000000000000000000000000000000000, +1.00000002468784697448240228843729433, +0.391085523137819425490135854698818053E-4, +1.00000002468784697448240228843729433, +0.438129154649774798484788272162522590, +0.246878469744824022884372943262141470E-7
428
429
430npg ! The number of points in the Gauss-Legendre quadrature rule.
431+14
432call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
433integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
434[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
435+1.00000000000000000000000000000000000, +1.00000001214523492936835710483860060, +0.869361016703036103633369759211822842E-6, +1.00000001214523492936835710483860060, +0.439848436180187079104872956531366532, +0.121452349293683571048386005983282358E-7
436
437
438npg ! The number of points in the Gauss-Legendre quadrature rule.
439+15
440call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
441integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
442[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
443+1.00000000000000000000000000000000000, +1.00000000524709950801412578397433880, +0.453855310215365294604835174030162360E-5, +1.00000000524709950801412578397433880, +0.439938137444660635481863189048735278, +0.524709950801412578397433879702286054E-8
444
445
446npg ! The number of points in the Gauss-Legendre quadrature rule.
447+16
448call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
449integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
450[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
451+1.00000000000000000000000000000000000, +1.00000000207728443433581482672936140, +0.634663524822756477103426987347064793E-5, +1.00000000207728443433581482672936140, +0.443801817673865015557406858740782638, +0.207728443433581482672936140386060941E-8
452
453
454npg ! The number of points in the Gauss-Legendre quadrature rule.
455+17
456call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
457integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
458[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
459+1.00000000000000000000000000000000000, +1.00000000060138193848116631223235291, +0.527551248353347731943956964631672420E-5, +1.00000000060138193848116631223235291, +0.443283372221765898366160921339535976, +0.601381938481166312232352911385048352E-9
460
461
462npg ! The number of points in the Gauss-Legendre quadrature rule.
463+18
464call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
465integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
466[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
467+1.00000000000000000000000000000000000, +1.00000000002337684456012439615534067, +0.358328444953570690596202537321935925E-5, +1.00000000002337684456012439615534067, +0.445857503888062707775638352030904815, +0.233768445601243961553406690044477554E-10
468
469
470npg ! The number of points in the Gauss-Legendre quadrature rule.
471+19
472call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
473integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
474[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
475+1.00000000000000000000000000000000000, +0.999999999827357937381159052944207481, +0.219175345657139401480060276808435034E-5, +0.999999999827357937381159052944207481, +0.447116124776096284469552631143427524, +0.172642062618840947055792519112665077E-9
476
477
478npg ! The number of points in the Gauss-Legendre quadrature rule.
479+20
480call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
481integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
482[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
483+1.00000000000000000000000000000000000, +0.999999999795081386330293090324938581, +0.125441756876698348170916154825542999E-5, +0.999999999795081386330293090324938581, +0.446857633223355518536040613616828666, +0.204918613669706909675061418782736738E-9
484
485
486npg ! The number of points in the Gauss-Legendre quadrature rule.
487+21
488call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
489integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
490[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
491+1.00000000000000000000000000000000000, +0.999999999823694334447531529435443256, +0.680808135708412824895859552233285701E-6, +0.999999999823694334447531529435443256, +0.449341868923503714056919684482635832, +0.176305665552468470564556744368424523E-9
492
493
494npg ! The number of points in the Gauss-Legendre quadrature rule.
495+22
496call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
497integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
498[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
499+1.00000000000000000000000000000000000, +0.999999999863908587433285676154034471, +0.355423836548274635452003188878673935E-6, +0.999999999863908587433285676154034471, +0.449228245653034840442224144352105582, +0.136091412566714323845965528887489281E-9
500
501
502npg ! The number of points in the Gauss-Legendre quadrature rule.
503+23
504call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
505integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
506[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
507+1.00000000000000000000000000000000000, +0.999999999903139510554949594591368601, +0.178094484871936054857683960232484674E-6, +0.999999999903139510554949594591368601, +0.450553065813618873446119310989872083, +0.968604894450504054086313994102149628E-10
508
509
510npg ! The number of points in the Gauss-Legendre quadrature rule.
511+24
512call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
513integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
514[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
515+1.00000000000000000000000000000000000, +0.999999999932930631263633182006277053, +0.855734277235847496242840188933614528E-7, +0.999999999932930631263633182006277053, +0.451528553045033627003930017346024233, +0.670693687363668179937229473940755511E-10
516
517
518npg ! The number of points in the Gauss-Legendre quadrature rule.
519+25
520call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
521integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
522[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
523+1.00000000000000000000000000000000000, +0.999999999955606206754649785466436407, +0.391224786873945559603990633870624924E-7, +0.999999999955606206754649785466436407, +0.451151051357497766162941844671250223, +0.443937932453502145335635927877274839E-10
524
525
526npg ! The number of points in the Gauss-Legendre quadrature rule.
527+26
528call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
529integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
530[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
531+1.00000000000000000000000000000000000, +0.999999999970905117727407763552248053, +0.166381656789620823670390268592916473E-7, +0.999999999970905117727407763552248053, +0.452917709799272825368192692193337627, +0.290948822725922364477519471835892967E-10
532
533
534npg ! The number of points in the Gauss-Legendre quadrature rule.
535+27
536call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
537integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
538[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
539+1.00000000000000000000000000000000000, +0.999999999981605731220288947315054768, +0.635900931335095252218769876840873286E-8, +0.999999999981605731220288947315054768, +0.452975362039077355596520008820995626, +0.183942687797110526849452318238009307E-10
540
541
542npg ! The number of points in the Gauss-Legendre quadrature rule.
543+28
544call setNodeWeightGK(nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2)) ! Compute the Gauss-Kronrod nodes and weights
545integral = getQuadGK(getLogNormPDF, 0._RK, pinf, nodeK(1:npg+1), weightK(1:npg+1), weightG(1:(npg+1)/2), abserr, intAbsFunc, smoothness) ! Compute the (2*npg + 1)-points Gauss-Kronrod quadrature.
546[truth, integral, abserr, intAbsFunc, smoothness, abs((integral - truth) / truth)]
547+1.0000000000000000000000000