Line data Source code
1 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3 : !!!! !!!!
4 : !!!! ParaMonte: Parallel Monte Carlo and Machine Learning Library. !!!!
5 : !!!! !!!!
6 : !!!! Copyright (C) 2012-present, The Computational Data Science Lab !!!!
7 : !!!! !!!!
8 : !!!! This file is part of the ParaMonte library. !!!!
9 : !!!! !!!!
10 : !!!! LICENSE !!!!
11 : !!!! !!!!
12 : !!!! https://github.com/cdslaborg/paramonte/blob/main/LICENSE.md !!!!
13 : !!!! !!!!
14 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
15 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16 :
17 : !> \brief
18 : !> This include file contains procedure implementation of [pm_mathCumSum](@ref pm_mathCumSum).
19 : !>
20 : !> \finmain
21 : !>
22 : !> \author
23 : !> \AmirShahmoradi, Sunday 3:33 AM, September 19, 2021, Dallas, TX
24 :
25 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
26 :
27 : !%%%%%%%%%%%%%%%%%%%%
28 : #if getCumPropExp_ENABLED
29 : !%%%%%%%%%%%%%%%%%%%%
30 :
31 : #if Def_ENABLED
32 : type(sequence_type), parameter :: control = sequence_type()
33 : #elif !(Sel_ENABLED || Seq_ENABLED)
34 : #error "Unrecognized interface."
35 : #endif
36 5992 : if (0_IK < size(array, 1, IK)) then
37 5992 : if (present(direction) .and. present(action)) then
38 2401 : if (same_type_as(direction, forward)) then
39 1193 : if (same_type_as(action, nothing)) then
40 588 : call setCumPropExp(cumPropExp, array, maxArray, control, forward, nothing)
41 605 : elseif (same_type_as(action, reverse)) then
42 605 : call setCumPropExp(cumPropExp, array, maxArray, control, forward, reverse)
43 : else
44 0 : error stop "@getCumPropExp(): Unrecognized action."
45 : end if
46 1208 : elseif (same_type_as(direction, backward)) then
47 1208 : if (same_type_as(action, nothing)) then
48 600 : call setCumPropExp(cumPropExp, array, maxArray, control, backward, nothing)
49 608 : elseif (same_type_as(action, reverse)) then
50 608 : call setCumPropExp(cumPropExp, array, maxArray, control, backward, reverse)
51 : else
52 0 : error stop "@getCumPropExp(): Unrecognized action."
53 : end if
54 : else
55 0 : error stop "@getCumPropExp(): Unrecognized direction."
56 : end if
57 3591 : elseif (present(direction)) then
58 1201 : if (same_type_as(direction, forward)) then
59 0 : call setCumPropExp(cumPropExp, array, maxArray, control, forward, nothing)
60 1201 : elseif (same_type_as(direction, backward)) then
61 1201 : call setCumPropExp(cumPropExp, array, maxArray, control, backward, nothing)
62 : else
63 0 : error stop "@getCumPropExp(): Unrecognized direction."
64 : end if
65 2390 : elseif (present(action)) then
66 1211 : if (same_type_as(action, nothing)) then
67 0 : call setCumPropExp(cumPropExp, array, maxArray, control, forward, nothing)
68 1211 : elseif (same_type_as(action, reverse)) then
69 1211 : call setCumPropExp(cumPropExp, array, maxArray, control, forward, reverse)
70 : else
71 0 : error stop "@getCumPropExp(): Unrecognized action."
72 : end if
73 : else
74 1179 : call setCumPropExp(cumPropExp, array, maxArray, control)
75 : end if
76 : end if
77 :
78 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
79 : #elif setCumPropExp_ENABLED && Old_ENABLED
80 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
81 :
82 : real(RKC), parameter :: LOGTINY = log(tiny(0._RKC))
83 : integer(IK) :: lenArray
84 1453 : lenArray = size(array, kind = IK)
85 1453 : CHECK_ASSERTION(__LINE__, 0_IK < lenArray, SK_"@setCumPropExp(): The condition `0 < size(array)` must hold. size(array) = "//getStr(lenArray)) ! fpp
86 22653 : CHECK_ASSERTION(__LINE__, maxval(array, 1) == maxArray, SK_"@setCumPropExp(): The condition `maxval(array, 1) == maxArray` must hold. maxval(array, 1), maxArray = "//getStr([maxval(array, 1), maxArray])) ! fpp
87 : #if For_ENABLED && Non_ENABLED
88 : block
89 : integer(IK) :: i
90 : real(RKC) :: cumPropExpInv
91 : #if Seq_ENABLED
92 428 : array(1) = exp(array(1) - maxArray)
93 : #elif Sel_ENABLED
94 : real(RKC) :: exponent
95 449 : exponent = array(1) - maxArray
96 449 : array(1) = 0._RKC
97 449 : if (LOGTINY < exponent) array(1) = exp(exponent)
98 : #endif
99 4676 : do i = 2, lenArray
100 : #if Seq_ENABLED
101 2279 : array(i) = array(i - 1) + exp(array(i) - maxArray)
102 : #elif Sel_ENABLED
103 1948 : exponent = array(i) - maxArray
104 1948 : array(i) = array(i - 1)
105 2397 : if (LOGTINY < exponent) array(i) = array(i) + exp(exponent)
106 : #endif
107 : end do
108 877 : cumPropExpInv = 1._RKC / array(lenArray)
109 4676 : do i = 1, lenArray - 1
110 4676 : array(i) = array(i) * cumPropExpInv
111 : end do
112 877 : array(lenArray) = 1._RKC
113 : end block
114 : #elif For_ENABLED && Rev_ENABLED
115 303 : call setCumPropExp(array, maxArray, control)
116 303 : call setReversed(array)
117 : #elif Bac_ENABLED && Non_ENABLED
118 272 : call setReversed(array)
119 272 : call setCumPropExp(array, maxArray, control)
120 : #elif Bac_ENABLED && Rev_ENABLED
121 : block
122 : integer(IK) :: i
123 : real(RKC) :: cumPropExpInv
124 : #if Seq_ENABLED
125 1 : array(lenArray) = exp(array(lenArray) - maxArray)
126 : #elif Sel_ENABLED
127 : real(RKC) :: exponent
128 0 : exponent = array(lenArray) - maxArray
129 0 : array(lenArray) = 0._RKC
130 0 : if (LOGTINY < exponent) array(lenArray) = exp(exponent)
131 : #endif
132 4 : do i = lenArray - 1, 1, -1
133 : #if Seq_ENABLED
134 4 : array(i) = array(i + 1) + exp(array(i) - maxArray)
135 : #elif Sel_ENABLED
136 0 : exponent = array(i) - maxArray
137 0 : array(i) = array(i + 1)
138 0 : if (LOGTINY < exponent) array(i) = array(i) + exp(exponent)
139 : #endif
140 : end do
141 1 : cumPropExpInv = 1._RKC / array(1)
142 1 : array(1) = 1._RKC
143 4 : do i = 2, lenArray
144 4 : array(i) = array(i) * cumPropExpInv
145 : end do
146 : end block
147 : #else
148 : #error "Unrecognized interface."
149 : #endif
150 :
151 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
152 : #elif setCumPropExp_ENABLED && New_ENABLED
153 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
154 :
155 : #if Sel_ENABLED
156 : real(RKC) :: exponent
157 : #elif !Seq_ENABLED
158 : #error "Unrecognized interface."
159 : #endif
160 : integer(IK) :: i, j
161 : integer(IK) :: lenArray
162 : real(RKC), parameter :: LOGTINY = log(tiny(0._RKC))
163 : real(RKC) :: cumPropExpInv
164 7522 : lenArray = size(array, kind = IK)
165 22566 : CHECK_ASSERTION(__LINE__, size(array, 1, IK) == size(cumPropExp, 1, IK), SK_"@setCumPropExp(): The condition `size(array, 1) == size(cumPropExp, 1)` must hold. size(array), size(cumPropExp) = "//getStr([size(array, 1, IK), size(cumPropExp,1, IK)]))
166 120282 : CHECK_ASSERTION(__LINE__, maxval(array, 1) == maxArray, SK_"@setCumPropExp(): The condition `maxval(array, 1) == maxArray` must hold. maxval(array, 1), maxArray = "//getStr([maxval(array, 1), maxArray]))
167 7522 : CHECK_ASSERTION(__LINE__, 0_IK < lenArray, SK_"@getCumPropExp()/setCumPropExp(): The condition `0 < size(array)` must hold. size(array) = "//getStr(lenArray))
168 :
169 : #if For_ENABLED && Non_ENABLED
170 :
171 : #if Seq_ENABLED
172 1308 : cumPropExp(1) = exp(array(1) - maxArray)
173 : #elif Sel_ENABLED
174 755 : exponent = array(1) - maxArray
175 755 : cumPropExp(1) = 0._RKC
176 755 : if (LOGTINY < exponent) cumPropExp(1) = exp(exponent)
177 : #endif
178 11439 : do i = 2, lenArray
179 : #if Seq_ENABLED
180 7290 : cumPropExp(i) = cumPropExp(i - 1) + exp(array(i) - maxArray)
181 : #elif Sel_ENABLED
182 3394 : exponent = array(i) - maxArray
183 4149 : if (LOGTINY < exponent) then
184 1713 : cumPropExp(i) = cumPropExp(i - 1) + exp(exponent)
185 : else
186 1681 : cumPropExp(i) = cumPropExp(i - 1)
187 : end if
188 : #endif
189 : end do
190 2063 : cumPropExpInv = 1._RKC / cumPropExp(lenArray)
191 2063 : cumPropExp(lenArray) = 1._RKC
192 11439 : do j = lenArray - 1, 1, -1
193 11439 : cumPropExp(j) = cumPropExp(j) * cumPropExpInv
194 : end do
195 :
196 : #elif For_ENABLED && Rev_ENABLED
197 :
198 : #if Seq_ENABLED
199 1349 : cumPropExp(lenArray) = exp(array(1) - maxArray)
200 : #elif Sel_ENABLED
201 732 : exponent = array(1) - maxArray
202 732 : cumPropExp(lenArray) = 0._RKC
203 732 : if (LOGTINY < exponent) cumPropExp(lenArray) = exp(exponent)
204 : #endif
205 11550 : do i = 2, lenArray
206 9469 : j = lenArray - i + 1_IK
207 : #if Seq_ENABLED
208 7497 : cumPropExp(j) = cumPropExp(j + 1) + exp(array(i) - maxArray)
209 : #elif Sel_ENABLED
210 3321 : cumPropExp(j) = cumPropExp(j + 1)
211 3321 : exponent = array(i) - maxArray
212 4053 : if (LOGTINY < exponent) cumPropExp(j) = cumPropExp(j) + exp(exponent)
213 : #endif
214 : end do
215 2081 : cumPropExpInv = 1._RKC / cumPropExp(1)
216 2081 : cumPropExp(1) = 1._RKC
217 11550 : do j = 2, lenArray
218 11550 : cumPropExp(j) = cumPropExp(j) * cumPropExpInv
219 : end do
220 :
221 : #elif Bac_ENABLED && Non_ENABLED
222 :
223 : #if Seq_ENABLED
224 1379 : cumPropExp(1) = exp(array(lenArray) - maxArray)
225 : #elif Sel_ENABLED
226 732 : cumPropExp(1) = 0._RKC
227 732 : exponent = array(lenArray) - maxArray
228 732 : if (LOGTINY < exponent) cumPropExp(1) = exp(exponent)
229 : #endif
230 11481 : do i = 2, lenArray
231 : #if Seq_ENABLED
232 7513 : cumPropExp(i) = cumPropExp(i - 1) + exp(array(lenArray - i + 1) - maxArray)
233 : #elif Sel_ENABLED
234 3236 : cumPropExp(i) = cumPropExp(i - 1)
235 3236 : exponent = array(lenArray - i + 1) - maxArray
236 3968 : if (LOGTINY < exponent) cumPropExp(i) = cumPropExp(i) + exp(exponent)
237 : #endif
238 : end do
239 2111 : cumPropExpInv = 1._RKC / cumPropExp(lenArray)
240 2111 : cumPropExp(lenArray) = 1._RKC
241 11481 : do j = lenArray - 1, 1, -1
242 11481 : cumPropExp(j) = cumPropExp(j) * cumPropExpInv
243 : end do
244 :
245 : #elif Bac_ENABLED && Rev_ENABLED
246 :
247 : #if Seq_ENABLED
248 646 : cumPropExp(lenArray) = exp(array(lenArray) - maxArray)
249 : #elif Sel_ENABLED
250 621 : cumPropExp(lenArray) = 0._RKC
251 621 : exponent = array(lenArray) - maxArray
252 621 : if (LOGTINY < exponent) cumPropExp(lenArray) = exp(exponent)
253 : #endif
254 6866 : do i = lenArray - 1, 1, -1
255 : #if Seq_ENABLED
256 3468 : cumPropExp(i) = cumPropExp(i + 1) + exp(array(i) - maxArray)
257 : #elif Sel_ENABLED
258 2777 : cumPropExp(i) = cumPropExp(i + 1)
259 2777 : exponent = array(i) - maxArray
260 3398 : if (LOGTINY < exponent) cumPropExp(i) = cumPropExp(i) + exp(exponent)
261 : #endif
262 : end do
263 1267 : cumPropExpInv = 1._RKC / cumPropExp(1)
264 1267 : cumPropExp(1) = 1._RKC
265 6866 : do j = 2, lenArray
266 6866 : cumPropExp(j) = cumPropExp(j) * cumPropExpInv
267 : end do
268 : #else
269 : #error "Unrecognized interface."
270 : #endif
271 :
272 : #else
273 : !%%%%%%%%%%%%%%%%%%%%%%%%
274 : #error "Unrecognized interface."
275 : !%%%%%%%%%%%%%%%%%%%%%%%%
276 : #endif
|
