16 integer(IK) :: info, infos(
3)
17 real(RKG) :: alpha, beta, betaInc, betaInv
18 real(RKG),
allocatable :: alphas(:), betas(:), betaIncs(:), betaInvs(:)
20 type(display_type) :: disp
30 call disp%show(
"alpha = getUnifRand(0.1_RKG, 10._RKG); beta = getUnifRand(0.1_RKG, 10._RKG)")
34 call disp%show(
"call setBetaInv(betaInv, betaInc, alpha, beta, getLogBeta(alpha, beta), signed, info)")
36 call disp%show(
"if (info /= 0) error stop 'Beta inversion failed.'")
37 if (info
/= 0)
error stop 'Beta inversion failed.'
40 call disp%show(
"getBetaInc(betaInv, alpha, beta, signed)")
49 call disp%show(
"alpha = getUnifRand(0.1_RKG, 10._RKG); beta = getUnifRand(0.1_RKG, 10._RKG)")
53 call disp%show(
"call setBetaInv(betaInv, betaInc, alpha, beta, getLogBeta(alpha, beta), signed, info)")
55 call disp%show(
"if (info /= 0) error stop 'Beta inversion failed.'")
56 if (info
/= 0)
error stop 'Beta inversion failed.'
59 call disp%show(
"getBetaInc(betaInv, alpha, beta, signed)")
66 call disp%show(
"betaInc = 1._RKG - epsilon(1._RKG)")
67 betaInc
= 1._RKG - epsilon(
1._RKG)
68 call disp%show(
"alpha = getUnifRand(0.1_RKG, 10._RKG); beta = getUnifRand(0.1_RKG, 10._RKG)")
72 call disp%show(
"call setBetaInv(betaInv, betaInc, alpha, beta, getLogBeta(alpha, beta), signed, info)")
74 call disp%show(
"if (info /= 0) error stop 'Beta inversion failed.'")
75 if (info
/= 0)
error stop 'Beta inversion failed.'
78 call disp%show(
"getBetaInc(betaInv, alpha, beta, signed)")
85 call disp%show(
"betaIncs = [0._RKG, .5_RKG, 1._RKG]")
86 betaIncs
= [
0._RKG, .
5_RKG,
1._RKG]
87 call disp%show(
"alpha = getUnifRand(0.1_RKG, 10._RKG); beta = getUnifRand(0.1_RKG, 10._RKG)")
91 call disp%show(
"call setResized(betaInvs, size(betaIncs, 1, IK))")
93 call disp%show(
"call setBetaInv(betaInvs, betaIncs, alpha, beta, getLogBeta(alpha, beta), signed, infos)")
95 call disp%show(
"if (any(infos /= 0)) error stop 'Beta inversion failed.'")
96 if (
any(infos
/= 0))
error stop 'Beta inversion failed.'
99 call disp%show(
"getBetaInc(betaInvs, alpha, beta, signed)")
106 call disp%show(
"betaIncs = [0._RKG, .5_RKG, 1._RKG]")
107 betaIncs
= [
0._RKG, .
5_RKG,
1._RKG]
108 call disp%show(
"alphas = [.1_RKG, 1._RKG, 10._RKG]")
109 alphas
= [.
1_RKG,
1._RKG,
10._RKG]
112 call disp%show(
"call setResized(betaInvs, size(betaIncs, 1, IK))")
114 call disp%show(
"call setBetaInv(betaInvs, betaIncs, alphas, beta, getLogBeta(alpha, beta), signed, infos)")
116 call disp%show(
"if (any(infos /= 0)) error stop 'Beta inversion failed.'")
117 if (
any(infos
/= 0))
error stop 'Beta inversion failed.'
120 call disp%show(
"getBetaInc(betaInvs, alpha = [.1_RKG, 1._RKG, 10._RKG], beta = 3._RKG, signed = signed)")
121 call disp%show(
getBetaInc(betaInvs, alpha
= [.
1_RKG,
1._RKG,
10._RKG], beta
= 3._RKG, signed
= signed) )
132 integer(IK) ,
parameter :: NP
= 1000
133 real(RKG) ,
parameter :: alpha(
*)
= [
0.1_RKG,
10._RKG,
1._RKG,
0.1_RKG,
10._RKG], beta(
*)
= [
0.1_RKG,
0.1_RKG,
1._RKG,
10._RKG,
10._RKG]
134 real(RKG) :: betaInv(
max(
size(alpha),
size(beta)))
135 real(RKG) :: betaInvRef(
size(betaInv))
136 real(RKG) :: betaInc(NP)
137 integer(IK) :: infos(
size(betaInv))
138 integer :: fileUnit, i, j
140 call setLinSpace(betaInc,
0._RKG,
1._RKG, fopen
= .true._LK, lopen
= .true._LK)
141 open(newunit
= fileUnit, file
= "setBetaInv.RK.txt")
144 if (
any(infos
/= 0))
error stop "Beta inversion failed."
145 write(fileUnit,
"(*(g0,:,','))") betaInc(i), betaInv
151 character(
*),
parameter :: RKSTR
= "RKS"
152 real(RKG) :: betaInv(
max(
size(alpha),
size(beta)))
153 open(newunit
= fileUnit, file
= "setBetaInv."//RKSTR
//".abserr.txt")
155 call setBetaInv(betaInv,
real(betaInc(i), RKG),
real(alpha, RKG),
real(beta, RKG),
getLogBeta(
real(alpha, RKG),
real(beta, RKG)), signed, infos)
156 if (
any(infos
/= 0))
error stop "Beta inversion failed."
157 betaInvRef
= getBetaInv(betaInc(i), alpha, beta, signed
= .true._LK)
158 write(fileUnit,
"(*(g0,:,','))") betaInc(i), abs(betaInv
- merge(
1 + betaInvRef, betaInvRef, betaInvRef
< 0))
165 character(
*),
parameter :: RKSTR
= "RKD"
166 real(RKG) :: betaInv(
max(
size(alpha),
size(beta)))
167 open(newunit
= fileUnit, file
= "setBetaInv."//RKSTR
//".abserr.txt")
169 call setBetaInv(betaInv,
real(betaInc(i), RKG),
real(alpha, RKG),
real(beta, RKG),
getLogBeta(
real(alpha, RKG),
real(beta, RKG)), signed, infos)
170 if (
any(infos
/= 0))
error stop "Beta inversion failed."
171 betaInvRef
= getBetaInv(betaInc(i), alpha, beta, signed
= .true._LK)
172 write(fileUnit,
"(*(g0,:,','))") betaInc(i), abs(betaInv
- merge(
1 + betaInvRef, betaInvRef, betaInvRef
< 0))
179 character(
*),
parameter :: RKSTR
= "RKH"
180 real(RKG) :: betaInv(
max(
size(alpha),
size(beta)))
181 open(newunit
= fileUnit, file
= "setBetaInv."//RKSTR
//".abserr.txt")
183 call setBetaInv(betaInv,
real(betaInc(i), RKG),
real(alpha, RKG),
real(beta, RKG),
getLogBeta(
real(alpha, RKG),
real(beta, RKG)), signed, infos)
184 if (
any(infos
/= 0))
error stop "Beta inversion failed."
185 betaInvRef
= getBetaInv(betaInc(i), alpha, beta, signed
= .true._LK)
186 write(fileUnit,
"(*(g0,:,','))") betaInc(i), abs(betaInv
- merge(
1 + betaInvRef, betaInvRef, betaInvRef
< 0))
Allocate or resize (shrink or expand) an input allocatable scalar string or array of rank 1....
Return the linSpace output argument with size(linSpace) elements of evenly-spaced values over the int...
Generate and return a scalar or a contiguous array of rank 1 of length s1 of randomly uniformly distr...
This is a generic method of the derived type display_type with pass attribute.
This is a generic method of the derived type display_type with pass attribute.
Generate and return the regularized Incomplete Beta Function as defined in the details section of pm...
Generate and return the regularized Inverse Incomplete Beta Function as defined in the details secti...
Generate and return the natural logarithm of the Beta Function as defined in the details section of ...
This module contains procedures and generic interfaces for resizing allocatable arrays of various typ...
This module contains procedures and generic interfaces for generating arrays with linear or logarithm...
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...
type(display_type) disp
This is a scalar module variable an object of type display_type for general display.
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
integer, parameter RK
The default real kind in the ParaMonte library: real64 in Fortran, c_double in C-Fortran Interoperati...
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
integer, parameter RKD
The double precision real kind in Fortran mode. On most platforms, this is an 64-bit real kind.
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
integer, parameter RKH
The scalar integer constant of intrinsic default kind, representing the highest-precision real kind t...
integer, parameter RKS
The single-precision real kind in Fortran mode. On most platforms, this is an 32-bit real kind.
Generate and return an object of type display_type.
7+3.3586226449696928,
+1.6777439673037353
9if (info
/= 0)
error stop 'Beta inversion failed.'
21+9.5141715968095379,
+5.0634513784404760
23if (info
/= 0)
error stop 'Beta inversion failed.'
32betaInc
= 1._RKG - epsilon(
1._RKG)
35+0.67793192402639557,
+0.60257550338669730
37if (info
/= 0)
error stop 'Beta inversion failed.'
46betaIncs
= [
0._RKG, .
5_RKG,
1._RKG]
49+6.6260890908919174,
+8.2940262014142299
52if (
any(infos
/= 0))
error stop 'Beta inversion failed.'
54+0.0000000000000000,
+0.44154198194786182,
+1.0000000000000000
56+0.0000000000000000,
+0.49999999999999978,
+1.0000000000000000
61betaIncs
= [
0._RKG, .
5_RKG,
1._RKG]
62alphas
= [.
1_RKG,
1._RKG,
10._RKG]
66if (
any(infos
/= 0))
error stop 'Beta inversion failed.'
68+0.0000000000000000,
+0.22994604328204467E-2,
+1.0000000000000000
69getBetaInc(betaInvs, alpha
= [.
1_RKG,
1._RKG,
10._RKG], beta
= 3._RKG, signed
= signed)
70+0.0000000000000000,
+0.68825309020540924E-2,
+1.0000000000000000
3import matplotlib.pyplot
as plt
13label = [
r"$\alpha, \beta = .1, .1$"
14 ,
r"$\alpha, \beta = 10, .1$"
15 ,
r"$\alpha, \beta = 1., 1.$"
16 ,
r"$\alpha, \beta = .1, 10$"
17 ,
r"$\alpha, \beta = 10, 10$"
20parent = os.path.basename(os.path.dirname(__file__))
21pattern = parent +
"*.txt"
22fileList = glob.glob(pattern)
25 df = pd.read_csv(file, delimiter =
",")
27 fig = plt.figure(figsize = 1.25 * np.array([6.4, 4.8]), dpi = 200)
30 for i
in range(1,len(df.values[0,:]+1)):
32 plt.plot( df.values[:, 0]
37 plt.xticks(fontsize = fontsize - 2)
38 plt.yticks(fontsize = fontsize - 2)
40 nbit = file.split(
".")[1][2:]
41 ax.set_xlabel(
"X : Regularized Incomplete Beta Function", fontsize = fontsize)
42 ax.set_ylabel(
"{}-bits Absolute Error: X - BetaInc(BetaInv(X))".format(nbit), fontsize = fontsize)
44 ax.set_xlabel(
"x", fontsize = fontsize)
45 ax.set_ylabel(
"Regularized Inverse Beta Function", fontsize = fontsize)
53 plt.grid(visible =
True, which =
"both", axis =
"both", color =
"0.85", linestyle =
"-")
54 ax.tick_params(axis =
"y", which =
"minor")
55 ax.tick_params(axis =
"x", which =
"minor")
58 plt.savefig(file.replace(
".txt",
".png"))