Line data Source code
1 : module test_simplex
2 :
3 : use num_def
4 : use num_lib
5 : use math_lib
6 :
7 : implicit none
8 :
9 : integer :: num_calls
10 : logical, parameter :: show_details = .false.
11 :
12 : contains
13 :
14 1 : subroutine do_test_simplex
15 : !call test_FR2 ! okay -- escapes from local min
16 1 : call test_FR4 ! okay -- escapes from local min
17 : !call test_FR6 ! okay -- escapes from local min
18 1 : call test_ER ! okay
19 1 : call test_WD ! okay
20 1 : call test_BLE ! okay
21 1 : call test_PS ! okay
22 1 : call test_TR ! okay
23 1 : end subroutine do_test_simplex
24 :
25 6 : subroutine test1_simplex( &
26 6 : n, x_first, x_lower, x_upper, simplex, &
27 : centroid_weight_power, enforce_bounds, &
28 : adaptive_random_search, fcn, str)
29 : integer, intent(in) :: n ! number of dimensions
30 : real(dp), dimension(:) :: x_first, x_lower, x_upper
31 : real(dp) :: simplex(:, :), centroid_weight_power
32 : logical, intent(in) :: enforce_bounds, adaptive_random_search
33 : interface
34 : include 'num_simplex_fcn.dek'
35 : end interface
36 : character(len=*) :: str
37 :
38 60 : real(dp) :: f(n + 1), x_final(n), x_atol, x_rtol, f_final
39 : integer :: iter_max, fcn_calls_max
40 : integer :: lrpar, lipar
41 : logical :: start_from_given_simplex_and_f
42 6 : integer, pointer :: ipar(:) ! (lipar)
43 6 : real(dp), pointer :: rpar(:) ! (lrpar)
44 : integer :: num_iters, num_fcn_calls, &
45 : num_fcn_calls_for_ars, num_accepted_for_ars, ierr
46 : integer :: seed
47 : real(dp) :: alpha, beta, gamma, delta
48 :
49 : include 'formats'
50 :
51 6 : write (*, *) 'testing NM_simplex with '//trim(str)
52 :
53 6 : num_calls = 0
54 6 : lrpar = 0; lipar = 0
55 6 : allocate (rpar(lrpar), ipar(lipar))
56 :
57 6 : x_atol = 1d-10
58 6 : x_rtol = 1d-10
59 :
60 6 : iter_max = 1000
61 6 : fcn_calls_max = iter_max*10
62 6 : seed = 1074698122
63 :
64 6 : start_from_given_simplex_and_f = .false.
65 :
66 6 : alpha = 1d0
67 6 : beta = 2d0
68 6 : gamma = 0.5d0
69 6 : delta = 0.5d0
70 :
71 : call NM_simplex( &
72 : n, x_lower, x_upper, x_first, x_final, f_final, &
73 : simplex, f, start_from_given_simplex_and_f, &
74 : fcn, x_atol, x_rtol, &
75 : iter_max, fcn_calls_max, &
76 : centroid_weight_power, enforce_bounds, &
77 : adaptive_random_search, seed, &
78 : alpha, beta, gamma, delta, &
79 : lrpar, rpar, lipar, ipar, &
80 : num_iters, num_fcn_calls, &
81 6 : num_fcn_calls_for_ars, num_accepted_for_ars, ierr)
82 :
83 6 : if (ierr /= 0) then
84 0 : write (*, *) 'failed in do_simplex'
85 : else
86 6 : if (f_final < 1d-10) then
87 6 : write (*, 1) 'found f_final < 1d-10'
88 : else
89 0 : write (*, 1) 'failed in do_simplex: f_final too large', f_final
90 : end if
91 : if (show_details) then
92 : write (*, 1) 'f_final', f_final
93 : write (*, 1) 'x_final', x_final(1:n)
94 : write (*, 2) 'num_iters', num_iters
95 : write (*, 2) 'num_fcn_calls', num_fcn_calls
96 : if (num_fcn_calls_for_ars > 0) &
97 : write (*, 2) 'num_fcn_calls_for_ars', num_fcn_calls_for_ars
98 : if (num_accepted_for_ars > 0) &
99 : write (*, 2) 'num_accepted_for_ars', num_accepted_for_ars
100 : end if
101 : end if
102 6 : write (*, *)
103 :
104 6 : deallocate (rpar, ipar)
105 :
106 6 : end subroutine test1_simplex
107 :
108 1 : subroutine test_WD
109 : ! testing with 4 dimensional Wood function
110 : integer, parameter :: n = 4
111 15 : real(dp), dimension(n) :: x_first, x_lower, x_upper
112 26 : real(dp) :: simplex(n, n + 1), centroid_weight_power
113 : logical :: enforce_bounds, adaptive_random_search
114 :
115 1 : x_first(1:n) = [-3d0, -1d0, -3d0, -1d0]
116 1 : x_lower(1:n) = [-4d0, -2d0, -4d0, -2d0]
117 1 : x_upper(1:n) = [2d0, 2d0, 2d0, 2d0]
118 :
119 1 : enforce_bounds = .true.
120 1 : adaptive_random_search = .true.
121 1 : centroid_weight_power = 1d0
122 :
123 : call test1_simplex( &
124 : n, x_first, x_lower, x_upper, simplex, &
125 : centroid_weight_power, enforce_bounds, &
126 1 : adaptive_random_search, fcn_WD, 'WD')
127 :
128 1 : end subroutine test_WD
129 :
130 766 : real(dp) function fcn_WD(n, x, lrpar, rpar, lipar, ipar, op_code, ierr)
131 : use const_def, only: dp
132 : integer, intent(in) :: n
133 : real(dp), intent(in) :: x(:) ! (n)
134 : integer, intent(in) :: lrpar, lipar
135 : integer, intent(inout), pointer :: ipar(:) ! (lipar)
136 : real(dp), intent(inout), pointer :: rpar(:) ! (lrpar)
137 : integer, intent(in) :: op_code
138 : integer, intent(out) :: ierr
139 766 : ierr = 0
140 766 : fcn_WD = WD(x)
141 766 : end function fcn_WD
142 :
143 766 : real(dp) function WD(x) ! Wood function
144 : real(dp), intent(in) :: x(:)
145 : integer :: i, n
146 766 : n = size(x, dim=1)
147 766 : WD = 0
148 1532 : do i = 1, n/4
149 : WD = WD + &
150 0 : 100d0*pow2(x(i + 1) - x(i)*x(i)) + pow2(1d0 - x(i)) + &
151 0 : 90d0*pow2(x(i + 3) - x(i + 2)*x(i + 2)) + pow2(1d0 - x(i + 2)) + &
152 1532 : 100d0*pow2(x(i + 1) + x(i + 3) - 2d0) + pow2(x(i + 1) - x(i + 3))/sqrt(10d0)
153 : end do
154 766 : num_calls = num_calls + 1
155 766 : end function WD
156 :
157 1 : subroutine test_ER
158 : ! testing with 4 dimensional extended Rosenbrock
159 : integer, parameter :: n = 4
160 15 : real(dp), dimension(n) :: x_first, x_lower, x_upper
161 26 : real(dp) :: simplex(n, n + 1), centroid_weight_power
162 : logical :: enforce_bounds, adaptive_random_search
163 :
164 5 : x_first(1:n) = 3d0
165 5 : x_lower(1:n) = -2d0
166 5 : x_upper(1:n) = 5d0
167 :
168 1 : enforce_bounds = .true.
169 1 : adaptive_random_search = .true.
170 1 : centroid_weight_power = 1d0
171 :
172 : call test1_simplex( &
173 : n, x_first, x_lower, x_upper, simplex, &
174 : centroid_weight_power, enforce_bounds, &
175 1 : adaptive_random_search, fcn_ER, 'ER')
176 :
177 1 : end subroutine test_ER
178 :
179 686 : real(dp) function fcn_ER(n, x, lrpar, rpar, lipar, ipar, op_code, ierr)
180 : use const_def, only: dp
181 : integer, intent(in) :: n
182 : real(dp), intent(in) :: x(:) ! (n)
183 : integer, intent(in) :: lrpar, lipar
184 : integer, intent(inout), pointer :: ipar(:) ! (lipar)
185 : real(dp), intent(inout), pointer :: rpar(:) ! (lrpar)
186 : integer, intent(in) :: op_code
187 : integer, intent(out) :: ierr
188 686 : ierr = 0
189 686 : fcn_ER = ER(x)
190 686 : end function fcn_ER
191 :
192 686 : real(dp) function ER(x) ! extended Rosenbrock
193 : real(dp), intent(in) :: x(:)
194 : integer :: i, n
195 686 : n = size(x, dim=1)
196 686 : ER = 0
197 2744 : do i = 1, n - 1
198 2744 : ER = ER + 100d0*pow2(x(i + 1) - x(i)*x(i)) + pow2(1D0 - x(i))
199 : end do
200 686 : num_calls = num_calls + 1
201 686 : end function ER
202 :
203 0 : subroutine test_FR2
204 : ! testing with 2 dimensional Freudenstein and Roth function
205 : integer, parameter :: n = 2
206 0 : real(dp), dimension(n) :: x_first, x_lower, x_upper
207 0 : real(dp) :: simplex(n, n + 1), centroid_weight_power
208 : logical :: enforce_bounds, adaptive_random_search
209 :
210 : ! NOTE --- this function has a local minimum
211 : ! and the starting values lead to the false-minimum first.
212 :
213 : ! true min = 0 at (5,4)
214 : ! false min = 48.98... at (11.41..., -0.8968...)
215 : ! starting at (0.5, -2) leads to the local min.
216 :
217 0 : x_first(1:n) = [0.5d0, -2d0]
218 0 : x_lower(1:n) = -2d0
219 0 : x_upper(1:n) = 6d0
220 :
221 0 : enforce_bounds = .false.
222 0 : centroid_weight_power = 1d0
223 0 : adaptive_random_search = .true.
224 :
225 : call test1_simplex( &
226 : n, x_first, x_lower, x_upper, simplex, &
227 : centroid_weight_power, enforce_bounds, &
228 0 : adaptive_random_search, fcn_FR, 'FR2')
229 :
230 0 : end subroutine test_FR2
231 :
232 1 : subroutine test_FR4
233 : ! testing with 4 dimensional Freudenstein and Roth function
234 : integer, parameter :: n = 4
235 15 : real(dp), dimension(n) :: x_first, x_lower, x_upper
236 26 : real(dp) :: simplex(n, n + 1), centroid_weight_power
237 : logical :: enforce_bounds, adaptive_random_search
238 :
239 : ! = 0 at (5,4)
240 : ! = 48.98... at (11.41..., -0.8968...)
241 : ! starting at (0.5, -2) leads to the bad local min.
242 :
243 1 : x_first(1:n) = [0.5d0, -2d0, 0.5d0, -2d0]
244 5 : x_lower(1:n) = -2d0
245 5 : x_upper(1:n) = 6d0
246 :
247 1 : enforce_bounds = .false.
248 1 : centroid_weight_power = 1d0
249 1 : adaptive_random_search = .true.
250 :
251 : call test1_simplex( &
252 : n, x_first, x_lower, x_upper, simplex, &
253 : centroid_weight_power, enforce_bounds, &
254 1 : adaptive_random_search, fcn_FR, 'FR4')
255 :
256 1 : end subroutine test_FR4
257 :
258 0 : subroutine test_FR6
259 : ! testing with 6 dimensional Freudenstein and Roth function
260 : integer, parameter :: n = 6
261 0 : real(dp), dimension(n) :: x_first, x_lower, x_upper
262 0 : real(dp) :: simplex(n, n + 1), centroid_weight_power
263 : logical :: enforce_bounds, adaptive_random_search
264 :
265 : ! = 0 at (5,4)
266 : ! = 48.98... at (11.41..., -0.8968...)
267 : ! starting at (0.5, -2) leads to the bad local min.
268 :
269 0 : x_first(1:n) = [0.5d0, -2d0, 0.5d0, -2d0, 0.5d0, -2d0]
270 0 : x_lower(1:n) = -2d0
271 0 : x_upper(1:n) = 6d0
272 :
273 0 : enforce_bounds = .false.
274 0 : centroid_weight_power = 1d0
275 0 : adaptive_random_search = .true.
276 :
277 : call test1_simplex( &
278 : n, x_first, x_lower, x_upper, simplex, &
279 : centroid_weight_power, enforce_bounds, &
280 0 : adaptive_random_search, fcn_FR, 'FR6')
281 :
282 0 : end subroutine test_FR6
283 :
284 7422 : real(dp) function fcn_FR(n, x, lrpar, rpar, lipar, ipar, op_code, ierr)
285 : use const_def, only: dp
286 : integer, intent(in) :: n
287 : real(dp), intent(in) :: x(:) ! (n)
288 : integer, intent(in) :: lrpar, lipar
289 : integer, intent(inout), pointer :: ipar(:) ! (lipar)
290 : real(dp), intent(inout), pointer :: rpar(:) ! (lrpar)
291 : integer, intent(in) :: op_code
292 : integer, intent(out) :: ierr
293 7422 : ierr = 0
294 7422 : fcn_FR = FR(x)
295 7422 : end function fcn_FR
296 :
297 7422 : real(dp) function FR(x) ! Freudenstein and Roth function
298 : real(dp), intent(in) :: x(:)
299 : integer :: i, n
300 7422 : n = size(x, dim=1)
301 7422 : FR = 0
302 22266 : do i = 1, n/2
303 : FR = FR + &
304 0 : pow2(-13d0 + x(2*i - 1) + ((5d0 - x(2*i))*x(2*i) - 2d0)*x(2*i)) + &
305 22266 : pow2(-29d0 + x(2*i - 1) + ((1d0 + x(2*i))*x(2*i) - 14d0)*x(2*i))
306 : end do
307 7422 : num_calls = num_calls + 1
308 7422 : end function FR
309 :
310 1 : subroutine test_BLE
311 : ! testing with 4 dimensional Beale function
312 : integer, parameter :: n = 4
313 15 : real(dp), dimension(n) :: x_first, x_lower, x_upper
314 26 : real(dp) :: simplex(n, n + 1), centroid_weight_power
315 : logical :: enforce_bounds, adaptive_random_search
316 :
317 5 : x_first(1:n) = 3d0
318 5 : x_lower(1:n) = -2d0
319 5 : x_upper(1:n) = 5d0
320 :
321 1 : enforce_bounds = .true.
322 1 : adaptive_random_search = .true.
323 1 : centroid_weight_power = 1d0
324 :
325 : call test1_simplex( &
326 : n, x_first, x_lower, x_upper, simplex, &
327 : centroid_weight_power, enforce_bounds, &
328 1 : adaptive_random_search, fcn_BLE, 'BLE')
329 :
330 1 : end subroutine test_BLE
331 :
332 612 : real(dp) function fcn_BLE(n, x, lrpar, rpar, lipar, ipar, op_code, ierr)
333 : use const_def, only: dp
334 : integer, intent(in) :: n
335 : real(dp), intent(in) :: x(:) ! (n)
336 : integer, intent(in) :: lrpar, lipar
337 : integer, intent(inout), pointer :: ipar(:) ! (lipar)
338 : real(dp), intent(inout), pointer :: rpar(:) ! (lrpar)
339 : integer, intent(in) :: op_code
340 : integer, intent(out) :: ierr
341 612 : ierr = 0
342 612 : fcn_BLE = BLE(x)
343 612 : end function fcn_BLE
344 :
345 612 : real(dp) function BLE(x) ! Beale function
346 : real(dp), intent(in) :: x(:)
347 : integer :: i, n
348 612 : n = size(x, dim=1)
349 612 : BLE = 0
350 1836 : do i = 1, n/2
351 : BLE = BLE + &
352 0 : pow2(1.5d0 - x(2*i - 1)*(1d0 - x(2*i))) + &
353 : pow2(2.25d0 - x(2*i - 1)*(1d0 - x(2*i)*x(2*i))) + &
354 1836 : pow2(2.625d0 - x(2*i - 1)*(1d0 - x(2*i)*x(2*i)*x(2*i)))
355 : end do
356 612 : num_calls = num_calls + 1
357 612 : end function BLE
358 :
359 1 : subroutine test_PS
360 : ! testing with 4 dimensional Powell singular function
361 : integer, parameter :: n = 4
362 15 : real(dp), dimension(n) :: x_first, x_lower, x_upper
363 26 : real(dp) :: simplex(n, n + 1), centroid_weight_power
364 : logical :: enforce_bounds, adaptive_random_search
365 :
366 1 : x_first(1:n) = [3d0, -1d0, 0d0, 1d0]
367 5 : x_lower(1:n) = -1d0
368 5 : x_upper(1:n) = 3d0
369 :
370 1 : enforce_bounds = .true.
371 1 : adaptive_random_search = .true.
372 1 : centroid_weight_power = 1d0
373 :
374 : call test1_simplex( &
375 : n, x_first, x_lower, x_upper, simplex, &
376 : centroid_weight_power, enforce_bounds, &
377 1 : adaptive_random_search, fcn_PS, 'PS')
378 :
379 1 : end subroutine test_PS
380 :
381 1350 : real(dp) function fcn_PS(n, x, lrpar, rpar, lipar, ipar, op_code, ierr)
382 : use const_def, only: dp
383 : integer, intent(in) :: n
384 : real(dp), intent(in) :: x(:) ! (n)
385 : integer, intent(in) :: lrpar, lipar
386 : integer, intent(inout), pointer :: ipar(:) ! (lipar)
387 : real(dp), intent(inout), pointer :: rpar(:) ! (lrpar)
388 : integer, intent(in) :: op_code
389 : integer, intent(out) :: ierr
390 1350 : ierr = 0
391 1350 : fcn_PS = PS(x)
392 1350 : end function fcn_PS
393 :
394 1350 : real(dp) function PS(x) ! Powell singular function
395 : real(dp), intent(in) :: x(:)
396 : integer :: i, n
397 1350 : n = size(x, dim=1)
398 1350 : PS = 0
399 2700 : do i = 1, n/4
400 : PS = PS + &
401 0 : pow2(x(i) + 10d0*x(i + 1)) + &
402 0 : 5d0*pow2(x(i + 2) - x(i + 3)) + &
403 : pow4(x(i + 1) - 2d0*x(i + 2)) + &
404 2700 : 10d0*pow4(x(i) - x(i + 3))
405 : end do
406 1350 : num_calls = num_calls + 1
407 1350 : end function PS
408 :
409 1 : subroutine test_TR
410 : ! testing with 4 dimensional Trigonometric function
411 : integer, parameter :: n = 4
412 15 : real(dp), dimension(n) :: x_first, x_lower, x_upper
413 26 : real(dp) :: simplex(n, n + 1), centroid_weight_power
414 : logical :: enforce_bounds, adaptive_random_search
415 :
416 5 : x_first(1:n) = 0.01d0
417 5 : x_lower(1:n) = -0.1d0
418 5 : x_upper(1:n) = 0.1d0
419 :
420 1 : enforce_bounds = .true.
421 1 : adaptive_random_search = .true.
422 1 : centroid_weight_power = 1d0
423 :
424 : call test1_simplex( &
425 : n, x_first, x_lower, x_upper, simplex, &
426 : centroid_weight_power, enforce_bounds, &
427 1 : adaptive_random_search, fcn_TR, 'TR')
428 :
429 1 : end subroutine test_TR
430 :
431 266 : real(dp) function fcn_TR(n, x, lrpar, rpar, lipar, ipar, op_code, ierr)
432 : use const_def, only: dp
433 : integer, intent(in) :: n
434 : real(dp), intent(in) :: x(:) ! (n)
435 : integer, intent(in) :: lrpar, lipar
436 : integer, intent(inout), pointer :: ipar(:) ! (lipar)
437 : real(dp), intent(inout), pointer :: rpar(:) ! (lrpar)
438 : integer, intent(in) :: op_code
439 : integer, intent(out) :: ierr
440 266 : ierr = 0
441 266 : fcn_TR = TR(x)
442 266 : end function fcn_TR
443 :
444 266 : real(dp) function TR(x) ! Trigonometric function
445 : real(dp), intent(in) :: x(:)
446 : integer :: j, i, n
447 266 : real(dp) :: sum_cos
448 266 : n = size(x, dim=1)
449 266 : sum_cos = 0
450 1330 : do j = 1, n
451 1330 : sum_cos = sum_cos + cos(x(j))
452 : end do
453 : TR = 0
454 1330 : do i = 1, n
455 1330 : TR = TR + pow2(n - sum_cos + i*(1d0 - cos(x(i))) - sin(x(i)))
456 : end do
457 266 : num_calls = num_calls + 1
458 266 : end function TR
459 :
460 : end module test_simplex
|
