Yet Another eXchange Tool 0.11.3
Loading...
Searching...
No Matches
ftest_common.f90
1!>
2!! @file ftest_common.f90
3!!
4!! @copyright Copyright (C) 2016 Jörg Behrens <behrens@dkrz.de>
5!! Moritz Hanke <hanke@dkrz.de>
6!! Thomas Jahns <jahns@dkrz.de>
7!!
8!! @author Jörg Behrens <behrens@dkrz.de>
9!! Moritz Hanke <hanke@dkrz.de>
10!! Thomas Jahns <jahns@dkrz.de>
11!!
12
13!
14! Maintainer: Jörg Behrens <behrens@dkrz.de>
15! Moritz Hanke <hanke@dkrz.de>
16! Thomas Jahns <jahns@dkrz.de>
17! URL: https://dkrz-sw.gitlab-pages.dkrz.de/yaxt/
18!
19! Redistribution and use in source and binary forms, with or without
20! modification, are permitted provided that the following conditions are
21! met:
22!
23! Redistributions of source code must retain the above copyright notice,
24! this list of conditions and the following disclaimer.
25!
26! Redistributions in binary form must reproduce the above copyright
27! notice, this list of conditions and the following disclaimer in the
28! documentation and/or other materials provided with the distribution.
29!
30! Neither the name of the DKRZ GmbH nor the names of its contributors
31! may be used to endorse or promote products derived from this software
32! without specific prior written permission.
33!
34! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
35! IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
36! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
37! PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
38! OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
39! EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
40! PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
41! PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
42! LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
43! NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
44! SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
45!
46#include "fc_feature_defs.inc"
47MODULE ftest_common
48#include "xt_slice_c_loc.inc"
49 USE mpi
50 USE xt_core, ONLY: i2, i4, i8
51 USE iso_c_binding, ONLY: c_int, c_double, &
52 c_int16_t, c_int32_t, c_int64_t, c_size_t, c_loc, c_ptr
53 IMPLICIT NONE
54 PRIVATE
55 INTEGER, PUBLIC, PARAMETER :: dp = selected_real_kind(12, 307)
56
57 TYPE timer
58 CHARACTER(len=20) :: label = 'undef'
59 INTEGER :: istate = -1
60 REAL(dp) :: t0 = 0.0_dp
61 REAL(dp) :: dt_work = 0.0_dp
62 END TYPE timer
63
64 INTERFACE test_abort
65 MODULE PROCEDURE test_abort_cmsl_f
66 MODULE PROCEDURE test_abort_msl_f
67 END INTERFACE test_abort
68
69 INTERFACE icmp
70 MODULE PROCEDURE icmp_2d
71 MODULE PROCEDURE icmp_2d_i2
72 MODULE PROCEDURE icmp_2d_i8
73 MODULE PROCEDURE icmp_3d
74 MODULE PROCEDURE icmp_3d_i2
75 MODULE PROCEDURE icmp_3d_i8
76 END INTERFACE icmp
77
78 INTERFACE
79 SUBROUTINE posix_exit(code) BIND(c, name="exit")
80 IMPORT :: c_int
81 INTEGER(c_int), VALUE, INTENT(in) :: code
82 END SUBROUTINE posix_exit
83 END INTERFACE
84 PUBLIC :: posix_exit
85
86 INTERFACE
87 ! this function is meant as an escape hatch to prevent warnings
88 ! about elision of non-pure calls
89 PURE FUNCTION pure_print(s) RESULT(is_true)
90 INTEGER :: is_true ! returns non-zero value to prevent
91 ! optimizing out of the call
92 CHARACTER(len=*), INTENT(in) :: s
93 END FUNCTION pure_print
94 END INTERFACE
95
96 INTERFACE cmp_arrays
97
98 PURE FUNCTION cmp_dbl_arrays(asize, a, b)
99 IMPORT :: c_double
100 INTEGER :: cmp_dbl_arrays
101 INTEGER, INTENT(in) :: asize
102 REAL(c_double), INTENT(in) :: a(asize), b(asize)
103 END FUNCTION cmp_dbl_arrays
104#define cmp_dbl_arrays(asize,a,b) (Cmp_dbl_arrays(asize,a,b) /= 0)
105
106 PURE FUNCTION cmp_int16_arrays(asize, a, b)
107 IMPORT :: c_int16_t
108 INTEGER :: cmp_int16_arrays
109 INTEGER, INTENT(in) :: asize
110 INTEGER(c_int16_t), INTENT(in) :: a(asize), b(asize)
111 END FUNCTION cmp_int16_arrays
112#define cmp_int16_arrays(asize,a,b) (Cmp_int16_arrays(asize,a,b) /= 0)
113
114 PURE FUNCTION cmp_int32_arrays(asize, a, b)
115 IMPORT :: c_int32_t
116 INTEGER :: cmp_int32_arrays
117 INTEGER, INTENT(in) :: asize
118 INTEGER(c_int32_t), INTENT(in) :: a(asize), b(asize)
119 END FUNCTION cmp_int32_arrays
120#define cmp_int32_arrays(asize,a,b) (Cmp_int32_arrays(asize,a,b) /= 0)
121
122 PURE FUNCTION cmp_int64_arrays(asize, a, b)
123 IMPORT :: c_int64_t
124 INTEGER :: cmp_int64_arrays
125 INTEGER, INTENT(in) :: asize
126 INTEGER(c_int64_t), INTENT(in) :: a(asize), b(asize)
127 END FUNCTION cmp_int64_arrays
128#define cmp_int64_arrays(asize,a,b) (Cmp_int64_arrays(asize,a,b) /= 0)
129
130 MODULE PROCEDURE cmp_dbl_arrays_a1d_a1d
131 MODULE PROCEDURE cmp_dbl_arrays_a2d_a2d
132 MODULE PROCEDURE cmp_dbl_arrays_a3d_a3d
133 MODULE PROCEDURE cmp_i2_arrays_a1d_a1d
134 MODULE PROCEDURE cmp_i4_arrays_a1d_a1d
135 MODULE PROCEDURE cmp_i8_arrays_a1d_a1d
136 MODULE PROCEDURE cmp_i2_arrays_a2d_a2d
137 MODULE PROCEDURE cmp_i4_i2_arrays_a2d_a2d
138 MODULE PROCEDURE cmp_i4_arrays_a2d_a2d
139 MODULE PROCEDURE cmp_i4_i8_arrays_a2d_a2d
140 MODULE PROCEDURE cmp_i8_arrays_a2d_a2d
141 MODULE PROCEDURE cmp_i2_arrays_a3d_a3d
142 MODULE PROCEDURE cmp_i4_i2_arrays_a3d_a3d
143 MODULE PROCEDURE cmp_i4_arrays_a3d_a3d
144 MODULE PROCEDURE cmp_i4_i8_arrays_a3d_a3d
145 MODULE PROCEDURE cmp_i8_arrays_a3d_a3d
146 END INTERFACE cmp_arrays
147
148 INTERFACE id_map
149 MODULE PROCEDURE id_map_i2, id_map_i4, id_map_i8
150 END INTERFACE id_map
151
152 INTERFACE icbrt
153 MODULE PROCEDURE icbrt_i2, icbrt_i4, icbrt_i8
154 END INTERFACE icbrt
155
156 INTERFACE random_fill
157 MODULE PROCEDURE random_fill_i2, random_fill_i4, random_fill_i8
158 END INTERFACE random_fill
159
160 INTERFACE crc32
161 FUNCTION memcrc(a, n) BIND(c, name="Xt_memcrc") RESULT(crcval)
162 IMPORT :: c_ptr, c_size_t, c_int32_t
163 TYPE(c_ptr), VALUE, INTENT(in) :: a
164 INTEGER(c_size_t), VALUE, INTENT(in) :: n
165 INTEGER(c_int32_t) :: crcval
166 END function memcrc
167 MODULE PROCEDURE memcrc_i2, memcrc_i4, memcrc_i8
168 END INTERFACE crc32
169
170 INTERFACE permute
171 MODULE PROCEDURE permute_i2
172 MODULE PROCEDURE permute_i4
173 MODULE PROCEDURE permute_i8
174 END INTERFACE permute
175
176 REAL(dp) :: sync_dt_sum = 0.0_dp
177 LOGICAL, PARAMETER :: debug = .false.
178 LOGICAL :: verbose = .false.
179
180 PUBLIC :: init_mpi, finish_mpi
181 PUBLIC :: timer, treset, tstart, tstop, treport, mysync
182 PUBLIC :: id_map, icbrt, icmp, factorize, regular_deco
183 PUBLIC :: test_abort, set_verbose, get_verbose
184 PUBLIC :: cmp_arrays, crc32
185 PUBLIC :: random_fill, permute
186 PUBLIC :: run_randomized_tests, init_fortran_random
187 CHARACTER(len=*), PARAMETER :: filename = 'ftest_common.f90'
188CONTAINS
189
190 SUBROUTINE init_mpi
191 CHARACTER(len=*), PARAMETER :: context = 'init_mpi: '
192 INTEGER :: ierror
193#ifndef _OPENMP
194 CALL mpi_init(ierror)
195#else
196 INTEGER :: th_provided
197 th_provided = mpi_thread_single
198 CALL mpi_init_thread(mpi_thread_multiple, th_provided, ierror)
199#endif
200 IF (ierror /= mpi_success) CALL test_abort(context//'MPI_INIT failed', &
201 filename, __line__)
202 END SUBROUTINE init_mpi
203
204 SUBROUTINE finish_mpi
205 CHARACTER(len=*), PARAMETER :: context = 'finish_mpi: '
206 INTEGER :: ierror
207 CALL mpi_finalize(ierror)
208 IF (ierror /= mpi_success) CALL test_abort(context//'MPI_FINALIZE failed', &
209 filename, &
210 __line__)
211 END SUBROUTINE finish_mpi
212
213 SUBROUTINE set_verbose(verb)
214 LOGICAL, INTENT(in) :: verb
215 verbose = verb
216 END SUBROUTINE set_verbose
217
218 SUBROUTINE get_verbose(verb)
219 LOGICAL, INTENT(out) :: verb
220 verb = verbose
221 END SUBROUTINE get_verbose
222
223 PURE SUBROUTINE treset(t, label)
224 TYPE(timer), INTENT(inout) :: t
225 CHARACTER(len=*), INTENT(in) :: label
226 t%label = label
227 t%istate = 0
228 t%t0 = 0.0_dp
229 t%dt_work = 0.0_dp
230 END SUBROUTINE treset
231
232 SUBROUTINE tstart(t)
233 TYPE(timer), INTENT(inout) :: t
234 IF (debug) WRITE(0,*) 'tstart: ',t%label
235 CALL mysync
236 t%istate = 1
237 t%t0 = work_time()
238 END SUBROUTINE tstart
239
240 SUBROUTINE tstop(t)
241 TYPE(timer), INTENT(inout) :: t
242 REAL(dp) :: t1
243 IF (debug) WRITE(0,*) 'tstop: ',t%label
244 t1 = work_time()
245 t%dt_work = t%dt_work + (t1 - t%t0)
246 t%istate = 0
247 CALL mysync
248
249 END SUBROUTINE tstop
250
251 SUBROUTINE treport(t,extra_label,comm)
252 TYPE(timer), INTENT(in) :: t
253 CHARACTER(len=*), INTENT(in) :: extra_label
254 INTEGER, INTENT(in) :: comm
255
256 CHARACTER(len=*), PARAMETER :: context = 'treport: '
257 REAL(dp) :: work_sum, work_max, work_avg, e
258 REAL(dp) :: sbuf(1)
259 REAL(dp), ALLOCATABLE :: rbuf(:)
260 INTEGER :: nprocs, rank, ierror
261
262 sbuf(1) = t%dt_work
263 CALL mpi_comm_rank(comm, rank, ierror)
264 IF (ierror /= mpi_success) &
265 CALL test_abort(context//'mpi_comm_rank failed', filename, __line__)
266 CALL mpi_comm_size(comm, nprocs, ierror)
267 IF (ierror /= mpi_success) &
268 CALL test_abort(context//'mpi_comm_size failed', filename, __line__)
269 ALLOCATE(rbuf(0:nprocs-1))
270 rbuf = -1.0_dp
271 CALL mpi_gather(sbuf, 1, mpi_double_precision, &
272 & rbuf, 1, mpi_double_precision, &
273 & 0, comm, ierror)
274 IF (ierror /= mpi_success) CALL test_abort(context//'MPI_GATHER failed', &
275 filename, __line__)
276
277 IF (rank == 0) THEN
278 IF (cmp_dbl_arrays(1, rbuf, sbuf)) &
279 CALL test_abort(context//'internal error (1)', &
280 filename, __line__)
281 IF (any(rbuf < 0.0_dp)) CALL test_abort(context//'internal error (2)', &
282 filename, __line__)
283 work_sum = sum(rbuf)
284 work_max = maxval(rbuf)
285 work_avg = work_sum / real(nprocs, dp)
286 e = work_avg / (work_max + 1.e-20_dp)
287
288 IF (verbose) WRITE(0,'(A,I4,2X,A16,3F18.8)') &
289 'nprocs, label, wmax, wavg, e =', &
290 nprocs, extra_label//':'//t%label, &
291 work_max, work_avg, e
292 ENDIF
293
294 END SUBROUTINE treport
295
296 SUBROUTINE mysync
297 CHARACTER(len=*), PARAMETER :: context = 'mysync: '
298 INTEGER :: ierror
299 REAL(dp) :: t0, dt
300
301 t0 = mpi_wtime()
302
303 CALL mpi_barrier(mpi_comm_world, ierror)
304 IF (ierror /= mpi_success) CALL test_abort(context//'MPI_BARRIER failed', &
305 filename, __line__)
306
307 dt = (mpi_wtime() - t0)
308 sync_dt_sum = sync_dt_sum + dt
309
310 END SUBROUTINE mysync
311
312 REAL(dp) FUNCTION work_time()
313 work_time = mpi_wtime() - sync_dt_sum
314 RETURN
315 END FUNCTION work_time
316
317 PURE SUBROUTINE id_map_i2(map, ofs)
318 INTEGER(i2), INTENT(out) :: map(:,:)
319 INTEGER(i2), OPTIONAL, INTENT(in) :: ofs
320
321 INTEGER :: i, j, m, n
322 INTEGER :: ofs_
323
324 m = SIZE(map, 1)
325 n = SIZE(map, 2)
326 IF (PRESENT(ofs)) THEN
327 ofs_ = int(ofs) - 1
328 ELSE
329 ofs_ = 0
330 END IF
331 DO j = 1, n
332 DO i = 1, m
333 map(i,j) = int((j - 1) * m + i + ofs_, i2)
334 ENDDO
335 ENDDO
336
337 END SUBROUTINE id_map_i2
338
339 PURE SUBROUTINE id_map_i4(map, ofs)
340 INTEGER(i4), INTENT(out) :: map(:,:)
341 INTEGER(i4), OPTIONAL, INTENT(in) :: ofs
342
343 INTEGER :: i, j, m, n
344 INTEGER(i4) :: ofs_
345
346 m = SIZE(map, 1)
347 n = SIZE(map, 2)
348 IF (PRESENT(ofs)) THEN
349 ofs_ = ofs - 1_i4
350 ELSE
351 ofs_ = 0_i4
352 END IF
353 DO j = 1, n
354 DO i = 1, m
355 map(i,j) = int((j - 1) * m + i, i4) + ofs_
356 ENDDO
357 ENDDO
358
359 END SUBROUTINE id_map_i4
360
361 PURE SUBROUTINE id_map_i8(map, ofs)
362 INTEGER(i8), INTENT(out) :: map(:,:)
363 INTEGER(i8), OPTIONAL, INTENT(in) :: ofs
364
365 INTEGER :: i, j, m, n
366 INTEGER(i8) :: ofs_
367
368 m = SIZE(map, 1)
369 n = SIZE(map, 2)
370 IF (PRESENT(ofs)) THEN
371 ofs_ = ofs - 1_i8
372 ELSE
373 ofs_ = 0_i8
374 END IF
375 DO j = 1, n
376 DO i = 1, m
377 map(i,j) = int((j - 1) * m + i, i8) + ofs_
378 ENDDO
379 ENDDO
380
381 END SUBROUTINE id_map_i8
382
383 SUBROUTINE test_abort_msl_f(msg, source, line)
384 CHARACTER(*), INTENT(in) :: msg
385 CHARACTER(*), INTENT(in) :: source
386 INTEGER, INTENT(in) :: line
387 CALL test_abort_cmsl_f(mpi_comm_world, msg, source, line)
388 END SUBROUTINE test_abort_msl_f
389
390 SUBROUTINE test_abort_cmsl_f(comm, msg, source, line)
391 INTEGER, INTENT(in):: comm
392 CHARACTER(*), INTENT(in) :: msg
393 CHARACTER(*), INTENT(in) :: source
394 INTEGER, INTENT(in) :: line
395
396 INTERFACE
397 SUBROUTINE c_abort() bind(c, name='abort')
398 END SUBROUTINE c_abort
399 END INTERFACE
400
401 INTEGER :: ierror
402 LOGICAL :: flag
403
404 WRITE (0, '(3a,i0,2a)') 'Fatal error in ', source, ', line ', line, &
405 ': ', trim(msg)
406 FLUSH(0)
407 CALL mpi_initialized(flag, ierror)
408 IF (ierror == mpi_success .AND. flag) &
409 CALL mpi_abort(comm, 1, ierror)
410 CALL c_abort
411 END SUBROUTINE test_abort_cmsl_f
412
413 SUBROUTINE icmp_2d_i2(label, f,g, rank)
414 CHARACTER(len=*), PARAMETER :: context = 'ftest_common::icmp_2d_i8: '
415 CHARACTER(len=*), INTENT(in) :: label
416 INTEGER, INTENT(in) :: f(:,:)
417 INTEGER(i2), INTENT(in) :: g(:,:)
418 INTEGER, INTENT(in) :: rank
419
420 INTEGER :: i, j, n1, n2
421 LOGICAL :: mismatch_found
422
423 n1 = SIZE(f,1)
424 n2 = SIZE(f,2)
425 IF (SIZE(g,1) /= n1 .OR. SIZE(g,2) /= n2) &
426 CALL test_abort(context//'shape mismatch error', filename, __line__)
427
428 mismatch_found = .false.
429 DO j = 1, n2
430 DO i = 1, n1
431 mismatch_found = mismatch_found .OR. f(i,j) /= int(g(i,j))
432 END DO
433 END DO
434 IF (mismatch_found) THEN
435 DO j = 1, n2
436 DO i = 1, n1
437 IF (f(i,j) /= int(g(i,j))) THEN
438 WRITE(0,'(2a,4(a,i0))') context, label, ' test failed: i=', &
439 i, ', j=', j, ', f(i,j)=', f(i,j), ', g(i,j)=', g(i,j)
440 CALL test_abort(context//label//' test failed', filename, __line__)
441 ENDIF
442 ENDDO
443 ENDDO
444 END IF
445 IF (verbose) WRITE(0,*) rank,':',context//label//' passed'
446 END SUBROUTINE icmp_2d_i2
447
448 SUBROUTINE icmp_2d(label, f,g, rank)
449 CHARACTER(len=*), PARAMETER :: context = 'ftest_common::icmp_2d: '
450 CHARACTER(len=*), INTENT(in) :: label
451 INTEGER, INTENT(in) :: f(:,:), g(:,:)
452 INTEGER, INTENT(in) :: rank
453
454 INTEGER :: i, j, n1, n2
455 LOGICAL :: mismatch_found
456
457 n1 = SIZE(f,1)
458 n2 = SIZE(f,2)
459 IF (SIZE(g,1) /= n1 .OR. SIZE(g,2) /= n2) &
460 CALL test_abort(context//'shape mismatch error', filename, __line__)
461
462 mismatch_found = .false.
463 DO j = 1, n2
464 DO i = 1, n1
465 mismatch_found = mismatch_found .OR. f(i,j) /= g(i,j)
466 END DO
467 END DO
468 IF (mismatch_found) THEN
469 DO j = 1, n2
470 DO i = 1, n1
471 IF (f(i,j) /= g(i,j)) THEN
472 WRITE(0,'(2a,4(a,i0))') context, label, ' test failed: i=', &
473 i, ', j=', j, ', f(i,j)=', f(i,j), ', g(i,j)=', g(i,j)
474 CALL test_abort(context//label//' test failed', filename, __line__)
475 ENDIF
476 ENDDO
477 ENDDO
478 END IF
479 IF (verbose) WRITE(0,*) rank,':',context//label//' passed'
480 END SUBROUTINE icmp_2d
481
482 SUBROUTINE icmp_2d_i8(label, f,g, rank)
483 CHARACTER(len=*), PARAMETER :: context = 'ftest_common::icmp_2d_i8: '
484 CHARACTER(len=*), INTENT(in) :: label
485 INTEGER, INTENT(in) :: f(:,:)
486 INTEGER(i8), INTENT(in) :: g(:,:)
487 INTEGER, INTENT(in) :: rank
488
489 INTEGER :: i, j, n1, n2
490 LOGICAL :: mismatch_found
491
492 n1 = SIZE(f,1)
493 n2 = SIZE(f,2)
494 IF (SIZE(g,1) /= n1 .OR. SIZE(g,2) /= n2) &
495 CALL test_abort(context//'shape mismatch error', filename, __line__)
496
497 mismatch_found = .false.
498 DO j = 1, n2
499 DO i = 1, n1
500 mismatch_found = mismatch_found .OR. int(f(i,j),i8) /= g(i,j)
501 END DO
502 END DO
503 IF (mismatch_found) THEN
504 DO j = 1, n2
505 DO i = 1, n1
506 IF (int(f(i,j), i8) /= g(i,j)) THEN
507 WRITE(0,'(2a,4(a,i0))') context, label, ' test failed: i=', &
508 i, ', j=', j, ', f(i,j)=', f(i,j), ', g(i,j)=', g(i,j)
509 CALL test_abort(context//label//' test failed', filename, __line__)
510 ENDIF
511 ENDDO
512 ENDDO
513 END IF
514 IF (verbose) WRITE(0,*) rank,':',context//label//' passed'
515 END SUBROUTINE icmp_2d_i8
516
517 SUBROUTINE icmp_3d_i2(label, f,g, rank)
518 CHARACTER(len=*), PARAMETER :: context = 'ftest_common::icmp_3d_i8: '
519 CHARACTER(len=*), INTENT(in) :: label
520 INTEGER, INTENT(in) :: f(:,:,:)
521 INTEGER(i2), INTENT(in) :: g(:,:,:)
522 INTEGER, INTENT(in) :: rank
523
524 INTEGER :: i, j, k, n1, n2, n3
525 LOGICAL :: mismatch_found
526
527 n1 = SIZE(f,1)
528 n2 = SIZE(f,2)
529 n3 = SIZE(f,3)
530 IF (SIZE(g,1) /= n1 .OR. SIZE(g,2) /= n2 .OR. SIZE(g,3) /= n3) &
531 CALL test_abort(context//label//' shape mismatch', filename, __line__)
532
533 mismatch_found = .false.
534 DO k = 1, n3
535 DO j = 1, n2
536 DO i = 1, n1
537 mismatch_found = mismatch_found &
538 .OR. f(i,j,k) /= int(g(i,j,k))
539 END DO
540 END DO
541 END DO
542 IF (mismatch_found) THEN
543 DO k = 1, n3
544 DO j = 1, n2
545 DO i = 1, n1
546 IF (f(i,j,k) /= int(g(i,j,k))) THEN
547 WRITE(0,*) context,label,&
548 ' test failed: i, j, k, f(i,j,k), g(i,j,k) =', &
549 i, j, k, f(i,j,k), g(i,j,k)
550 CALL test_abort(context//label//' test failed', &
551 filename, __line__)
552 ENDIF
553 ENDDO
554 ENDDO
555 ENDDO
556 END IF
557 IF (verbose) WRITE(0,*) rank,':',context//label//' passed'
558 END SUBROUTINE icmp_3d_i2
559
560
561 SUBROUTINE icmp_3d(label, f,g, rank)
562 CHARACTER(len=*), PARAMETER :: context = 'ftest_common::icmp_3d: '
563 CHARACTER(len=*), INTENT(in) :: label
564 INTEGER, INTENT(in) :: f(:,:,:), g(:,:,:)
565 INTEGER, INTENT(in) :: rank
566
567 INTEGER :: i1, i2, i3, n1, n2, n3
568 LOGICAL :: mismatch_found
569
570 n1 = SIZE(f,1)
571 n2 = SIZE(f,2)
572 n3 = SIZE(f,3)
573 IF (SIZE(g,1) /= n1 .OR. SIZE(g,2) /= n2 .OR. SIZE(g,3) /= n3) &
574 CALL test_abort(context//label//' shape mismatch', filename, __line__)
575
576 mismatch_found = .false.
577 DO i3 = 1, n3
578 DO i2 = 1, n2
579 DO i1 = 1, n1
580 mismatch_found = mismatch_found .OR. f(i1,i2,i3) /= g(i1,i2,i3)
581 END DO
582 END DO
583 END DO
584 IF (mismatch_found) THEN
585 DO i3 = 1, n3
586 DO i2 = 1, n2
587 DO i1 = 1, n1
588 IF (f(i1,i2,i3) /= g(i1,i2,i3)) THEN
589 WRITE(0,*) context,label,&
590 ' test failed: i1, i2, i3, f(i1,i2,i3), g(i1,i2,i3) =', &
591 i1, i2, i3, f(i1,i2,i3), g(i1,i2,i3)
592 CALL test_abort(context//label//' test failed', &
593 filename, __line__)
594 ENDIF
595 ENDDO
596 ENDDO
597 ENDDO
598 END IF
599 IF (verbose) WRITE(0,*) rank,':',context//label//' passed'
600 END SUBROUTINE icmp_3d
601
602 SUBROUTINE icmp_3d_i8(label, f,g, rank)
603 CHARACTER(len=*), PARAMETER :: context = 'ftest_common::icmp_3d_i8: '
604 CHARACTER(len=*), INTENT(in) :: label
605 INTEGER, INTENT(in) :: f(:,:,:)
606 INTEGER(i8), INTENT(in) :: g(:,:,:)
607 INTEGER, INTENT(in) :: rank
608
609 INTEGER :: i1, i2, i3, n1, n2, n3
610 LOGICAL :: mismatch_found
611
612 n1 = SIZE(f,1)
613 n2 = SIZE(f,2)
614 n3 = SIZE(f,3)
615 IF (SIZE(g,1) /= n1 .OR. SIZE(g,2) /= n2 .OR. SIZE(g,3) /= n3) &
616 CALL test_abort(context//label//' shape mismatch', filename, __line__)
617
618 mismatch_found = .false.
619 DO i3 = 1, n3
620 DO i2 = 1, n2
621 DO i1 = 1, n1
622 mismatch_found = mismatch_found &
623 .OR. int(f(i1,i2,i3), i8) /= g(i1,i2,i3)
624 END DO
625 END DO
626 END DO
627 IF (mismatch_found) THEN
628 DO i3 = 1, n3
629 DO i2 = 1, n2
630 DO i1 = 1, n1
631 IF (int(f(i1,i2,i3), i8) /= g(i1,i2,i3)) THEN
632 WRITE(0,*) context,label,&
633 ' test failed: i1, i2, i3, f(i1,i2,i3), g(i1,i2,i3) =', &
634 i1, i2, i3, f(i1,i2,i3), g(i1,i2,i3)
635 CALL test_abort(context//label//' test failed', &
636 filename, __line__)
637 ENDIF
638 ENDDO
639 ENDDO
640 ENDDO
641 END IF
642 IF (verbose) WRITE(0,*) rank,':',context//label//' passed'
643 END SUBROUTINE icmp_3d_i8
644
645
646 PURE FUNCTION cmp_dbl_arrays_a1d_a1d(a, b) RESULT(differ)
647 DOUBLE PRECISION, INTENT(in) :: a(:), b(:)
648 LOGICAL :: differ
649 INTEGER :: asize, bsize
650 INTEGER(c_int) :: asize_c
651
652 asize = SIZE(a)
653 bsize = SIZE(b)
654 IF (asize /= bsize) THEN
655 differ = 0 /= pure_print('warning: comparing arrays of different size')
656 ELSE IF (asize > 0) THEN
657 asize_c = int(asize, c_int)
658 differ = cmp_dbl_arrays(asize_c, a, b)
659 ELSE
660 differ = .false.
661 END IF
662 END FUNCTION cmp_dbl_arrays_a1d_a1d
663
664 PURE FUNCTION cmp_dbl_arrays_a2d_a2d(a, b) RESULT(differ)
665 DOUBLE PRECISION, INTENT(in) :: a(:,:), b(:,:)
666 LOGICAL :: differ
667 INTEGER :: asize, bsize
668 INTEGER(c_int) :: asize_c
669
670 asize = SIZE(a)
671 bsize = SIZE(b)
672 IF (asize /= bsize) THEN
673 differ = 0 /= pure_print('warning: comparing arrays of different size')
674 ELSE IF (asize > 0) THEN
675 asize_c = int(asize, c_int)
676 differ = cmp_dbl_arrays(asize_c, a, b)
677 ELSE
678 differ = .false.
679 END IF
680 END FUNCTION cmp_dbl_arrays_a2d_a2d
681
682 PURE FUNCTION cmp_dbl_arrays_a3d_a3d(a, b) RESULT(differ)
683 DOUBLE PRECISION, INTENT(in) :: a(:,:,:), b(:,:,:)
684 LOGICAL :: differ
685 INTEGER :: asize, bsize
686 INTEGER(c_int) :: asize_c
687
688 asize = SIZE(a)
689 bsize = SIZE(b)
690 IF (asize /= bsize) THEN
691 differ = 0 /= pure_print('warning: comparing arrays of different size')
692 ELSE IF (asize > 0) THEN
693 asize_c = int(asize, c_int)
694 differ = cmp_dbl_arrays(asize_c, a, b)
695 ELSE
696 differ = .false.
697 END IF
698 END FUNCTION cmp_dbl_arrays_a3d_a3d
699
700 PURE FUNCTION cmp_i2_arrays_a1d_a1d(a, b) RESULT(differ)
701 INTEGER(i2), INTENT(in) :: a(:), b(:)
702 LOGICAL :: differ
703 INTEGER :: asize, bsize
704 INTEGER(c_int) :: asize_c
705
706 asize = SIZE(a)
707 bsize = SIZE(b)
708 IF (asize /= bsize) THEN
709 differ = 0 /= pure_print('warning: comparing arrays of different size')
710 ELSE IF (asize > 0) THEN
711 asize_c = int(asize, c_int)
712 differ = cmp_int16_arrays(asize_c, a, b)
713 ELSE
714 differ = .false.
715 END IF
716 END FUNCTION cmp_i2_arrays_a1d_a1d
717
718 PURE FUNCTION cmp_i4_arrays_a1d_a1d(a, b) RESULT(differ)
719 INTEGER(i4), INTENT(in) :: a(:), b(:)
720 LOGICAL :: differ
721 INTEGER :: asize, bsize
722 INTEGER(c_int) :: asize_c
723
724 asize = SIZE(a)
725 bsize = SIZE(b)
726 IF (asize /= bsize) THEN
727 differ = 0 /= pure_print('warning: comparing arrays of different size')
728 ELSE IF (asize > 0) THEN
729 asize_c = int(asize, c_int)
730 differ = cmp_int32_arrays(asize_c, a, b)
731 ELSE
732 differ = .false.
733 END IF
734 END FUNCTION cmp_i4_arrays_a1d_a1d
735
736 PURE FUNCTION cmp_i8_arrays_a1d_a1d(a, b) RESULT(differ)
737 INTEGER(i8), INTENT(in) :: a(:), b(:)
738 LOGICAL :: differ
739 INTEGER :: asize, bsize
740 INTEGER(c_int) :: asize_c
741
742 asize = SIZE(a)
743 bsize = SIZE(b)
744 IF (asize /= bsize) THEN
745 differ = 0 /= pure_print('warning: comparing arrays of different size')
746 ELSE IF (asize > 0) THEN
747 asize_c = int(asize, c_int)
748 differ = cmp_int64_arrays(asize_c, a, b)
749 ELSE
750 differ = .false.
751 END IF
752 END FUNCTION cmp_i8_arrays_a1d_a1d
753
754 PURE FUNCTION cmp_i2_arrays_a2d_a2d(a, b) RESULT(differ)
755 INTEGER(i2), INTENT(in) :: a(:,:), b(:,:)
756 LOGICAL :: differ
757 INTEGER :: asize, bsize
758 INTEGER(c_int) :: asize_c
759
760 asize = SIZE(a)
761 bsize = SIZE(b)
762 IF (asize /= bsize) THEN
763 differ = 0 /= pure_print('warning: comparing arrays of different size')
764 ELSE IF (asize > 0) THEN
765 asize_c = int(asize, c_int)
766 differ = cmp_int16_arrays(asize_c, a, b)
767 ELSE
768 differ = .false.
769 END IF
770 END FUNCTION cmp_i2_arrays_a2d_a2d
771
772 PURE FUNCTION cmp_i4_arrays_a2d_a2d(a, b) RESULT(differ)
773 INTEGER(i4), INTENT(in) :: a(:,:), b(:,:)
774 LOGICAL :: differ
775 INTEGER :: asize, bsize
776 INTEGER(c_int) :: asize_c
777
778 asize = SIZE(a)
779 bsize = SIZE(b)
780 IF (asize /= bsize) THEN
781 differ = 0 /= pure_print('warning: comparing arrays of different size')
782 ELSE IF (asize > 0) THEN
783 asize_c = int(asize, c_int)
784 differ = cmp_int32_arrays(asize_c, a, b)
785 ELSE
786 differ = .false.
787 END IF
788 END FUNCTION cmp_i4_arrays_a2d_a2d
789
790 PURE FUNCTION cmp_i8_arrays_a2d_a2d(a, b) RESULT(differ)
791 INTEGER(i8), INTENT(in) :: a(:,:), b(:,:)
792 LOGICAL :: differ
793 INTEGER :: asize, bsize
794 INTEGER(c_int) :: asize_c
795
796 asize = SIZE(a)
797 bsize = SIZE(b)
798 IF (asize /= bsize) THEN
799 differ = 0 /= pure_print('warning: comparing arrays of different size')
800 ELSE IF (asize > 0) THEN
801 asize_c = int(asize, c_int)
802 differ = cmp_int64_arrays(asize_c, a, b)
803 ELSE
804 differ = .false.
805 END IF
806 END FUNCTION cmp_i8_arrays_a2d_a2d
807
808 PURE FUNCTION cmp_i2_arrays_a3d_a3d(a, b) RESULT(differ)
809 INTEGER(i2), INTENT(in) :: a(:,:,:), b(:,:,:)
810 LOGICAL :: differ
811 INTEGER :: asize, bsize
812 INTEGER(c_int) :: asize_c
813
814 asize = SIZE(a)
815 bsize = SIZE(b)
816 IF (asize /= bsize) THEN
817 differ = 0 /= pure_print('warning: comparing arrays of different size')
818 ELSE IF (asize > 0) THEN
819 asize_c = int(asize, c_int)
820 differ = cmp_int16_arrays(asize_c, a, b)
821 ELSE
822 differ = .false.
823 END IF
824 END FUNCTION cmp_i2_arrays_a3d_a3d
825
826 PURE FUNCTION cmp_i4_arrays_a3d_a3d(a, b) RESULT(differ)
827 INTEGER(i4), INTENT(in) :: a(:,:,:), b(:,:,:)
828 LOGICAL :: differ
829 INTEGER :: asize, bsize
830 INTEGER(c_int) :: asize_c
831
832 asize = SIZE(a)
833 bsize = SIZE(b)
834 IF (asize /= bsize) THEN
835 differ = 0 /= pure_print('warning: comparing arrays of different size')
836 ELSE IF (asize > 0) THEN
837 asize_c = int(asize, c_int)
838 differ = cmp_int32_arrays(asize_c, a, b)
839 ELSE
840 differ = .false.
841 END IF
842 END FUNCTION cmp_i4_arrays_a3d_a3d
843
844 PURE FUNCTION cmp_i4_i2_arrays_a2d_a2d(a, b) RESULT(differ)
845 INTEGER(i4), INTENT(in) :: a(:,:)
846 INTEGER(i2), INTENT(in) :: b(:,:)
847 LOGICAL :: differ
848 INTEGER :: i, j, m, n
849
850 m = SIZE(a, 1)
851 n = SIZE(a, 2)
852 IF (m /= SIZE(b, 1) .OR. n /= SIZE(b, 2)) THEN
853 differ = 0 /= pure_print('warning: comparing arrays of different size')
854 ELSE IF (SIZE(a) > 0) THEN
855 differ = .false.
856 DO j = 1, n
857 DO i = 1, m
858 differ = differ .OR. a(i, j) /= int(b(i, j), i4)
859 END DO
860 END DO
861 ELSE
862 differ = .false.
863 END IF
864 END FUNCTION cmp_i4_i2_arrays_a2d_a2d
865
866 PURE FUNCTION cmp_i4_i8_arrays_a2d_a2d(a, b) RESULT(differ)
867 INTEGER(i4), INTENT(in) :: a(:,:)
868 INTEGER(i8), INTENT(in) :: b(:,:)
869 LOGICAL :: differ
870 INTEGER :: i, j, m, n
871
872 m = SIZE(a, 1)
873 n = SIZE(a, 2)
874 IF (m /= SIZE(b, 1) .OR. n /= SIZE(b, 2)) THEN
875 differ = 0 /= pure_print('warning: comparing arrays of different size')
876 ELSE IF (SIZE(a) > 0) THEN
877 differ = .false.
878 DO j = 1, n
879 DO i = 1, m
880 differ = differ .OR. int(a(i, j), i8) /= b(i, j)
881 END DO
882 END DO
883 ELSE
884 differ = .false.
885 END IF
886 END FUNCTION cmp_i4_i8_arrays_a2d_a2d
887
888 PURE FUNCTION cmp_i4_i2_arrays_a3d_a3d(a, b) RESULT(differ)
889 INTEGER(i4), INTENT(in) :: a(:,:,:)
890 INTEGER(i2), INTENT(in) :: b(:,:,:)
891 LOGICAL :: differ
892 INTEGER :: i, j, k, m, n, o
893
894 m = SIZE(a, 1)
895 n = SIZE(a, 2)
896 o = SIZE(a, 3)
897 IF (m /= SIZE(b, 1) .OR. n /= SIZE(b, 2) .OR. o /= SIZE(b, 3)) THEN
898 differ = 0 /= pure_print('warning: comparing arrays of different size')
899 ELSE IF (SIZE(a) > 0) THEN
900 differ = .false.
901 DO k = 1, o
902 DO j = 1, n
903 DO i = 1, m
904 differ = differ .OR. a(i, j, k) /= int(b(i, j, k), i4)
905 END DO
906 END DO
907 END DO
908 ELSE
909 differ = .false.
910 END IF
911 END FUNCTION cmp_i4_i2_arrays_a3d_a3d
912
913 PURE FUNCTION cmp_i4_i8_arrays_a3d_a3d(a, b) RESULT(differ)
914 INTEGER(i4), INTENT(in) :: a(:,:,:)
915 INTEGER(i8), INTENT(in) :: b(:,:,:)
916 LOGICAL :: differ
917 INTEGER :: i, j, k, m, n, o
918
919 m = SIZE(a, 1)
920 n = SIZE(a, 2)
921 o = SIZE(a, 3)
922 IF (m /= SIZE(b, 1) .OR. n /= SIZE(b, 2) .OR. o /= SIZE(b, 3)) THEN
923 differ = 0 /= pure_print('warning: comparing arrays of different size')
924 ELSE IF (SIZE(a) > 0) THEN
925 differ = .false.
926 DO k = 1, o
927 DO j = 1, n
928 DO i = 1, m
929 differ = differ .OR. int(a(i, j, k), i8) /= b(i, j, k)
930 END DO
931 END DO
932 END DO
933 ELSE
934 differ = .false.
935 END IF
936 END FUNCTION cmp_i4_i8_arrays_a3d_a3d
937
938 PURE FUNCTION cmp_i8_arrays_a3d_a3d(a, b) RESULT(differ)
939 INTEGER(i8), INTENT(in) :: a(:,:,:), b(:,:,:)
940 LOGICAL :: differ
941 INTEGER :: asize, bsize
942 INTEGER(c_int) :: asize_c
943
944 asize = SIZE(a)
945 bsize = SIZE(b)
946 IF (asize /= bsize) THEN
947 differ = 0 /= pure_print('warning: comparing arrays of different size')
948 ELSE IF (asize > 0) THEN
949 asize_c = int(asize, c_int)
950 differ = cmp_int64_arrays(asize_c, a, b)
951 ELSE
952 differ = .false.
953 END IF
954 END FUNCTION cmp_i8_arrays_a3d_a3d
955
956 SUBROUTINE factorize(c, a, b)
957 INTEGER, INTENT(in) :: c
958 INTEGER, INTENT(out) :: a, b ! c = a*b
959
960 INTEGER :: x0, i
961
962 IF (c<1) CALL test_abort('factorize: invalid process space', &
963 filename, __line__)
964 IF (c <= 3 .OR. c == 5 .OR. c == 7) THEN
965 a = c
966 b = 1
967 RETURN
968 ENDIF
969
970 ! simple approach, we try to be near c = (2*x) * x
971 x0 = int(sqrt(0.5 * real(c)) + 0.5)
972 a = 2*x0
973 f_loop: DO i = a, 1, -1
974 IF (mod(c,i) == 0) THEN
975 a = i
976 b = c/i
977 EXIT f_loop
978 ENDIF
979 ENDDO f_loop
980
981 END SUBROUTINE factorize
982
983 !> computes n**(1/3)
984 FUNCTION icbrt_i2(n) RESULT(icbrt)
985 INTEGER(i2), INTENT(in) :: n
986 INTEGER(i2) :: icbrt
987 INTEGER(i2), PARAMETER :: nbits = bit_size(n)-1_i2
988 INTEGER(i2) :: s
989 INTEGER(i2) :: b, x
990
991 x = abs(n)
992 icbrt = 0_i2
993 DO s = nbits, 0_i2, -3_i2
994 icbrt = icbrt + icbrt
995 b = 3_i2 * icbrt * (icbrt + 1_i2) + 1_i2
996 IF (ishft(x, -s) >= b) THEN
997 x = x - ishft(b, s)
998 icbrt = icbrt + 1_i2
999 END IF
1000 END DO
1001 icbrt = sign(icbrt, n)
1002 END FUNCTION icbrt_i2
1003
1004 !> computes n**(1/3)
1005 FUNCTION icbrt_i4(n) RESULT(icbrt)
1006 INTEGER(i4), INTENT(in) :: n
1007 INTEGER(i4) :: icbrt
1008 INTEGER(i4), PARAMETER :: nbits = bit_size(n)-1_i4
1009 INTEGER(i4) :: s
1010 INTEGER(i4) :: b, x
1011
1012 x = abs(n)
1013 icbrt = 0_i4
1014 DO s = nbits-1, 0_i4, -3_i4
1015 icbrt = icbrt + icbrt
1016 b = 3_i4 * icbrt * (icbrt + 1_i4) + 1_i4
1017 IF (ishft(x, -s) >= b) THEN
1018 x = x - ishft(b, s)
1019 icbrt = icbrt + 1_i4
1020 END IF
1021 END DO
1022 icbrt = sign(icbrt, n)
1023 END FUNCTION icbrt_i4
1024
1025 !> computes n**(1/3)
1026 FUNCTION icbrt_i8(n) RESULT(icbrt)
1027 INTEGER(i8), INTENT(in) :: n
1028 INTEGER(i8) :: icbrt
1029 INTEGER(i8), PARAMETER :: nbits = bit_size(n)-1_i8
1030 INTEGER(i8) :: s
1031 INTEGER(i8) :: b, x
1032
1033 x = abs(n)
1034 icbrt = 0_i8
1035 DO s = nbits, 0_i8, -3_i8
1036 icbrt = icbrt + icbrt
1037 b = 3_i8 * icbrt * (icbrt + 1_i8) + 1_i8
1038 IF (ishft(x, -s) >= b) THEN
1039 x = x - ishft(b, s)
1040 icbrt = icbrt + 1_i8
1041 END IF
1042 END DO
1043 icbrt = sign(icbrt, n)
1044 END FUNCTION icbrt_i8
1045
1046
1047 SUBROUTINE regular_deco(g_cn, c0, cn)
1048 INTEGER, INTENT(in) :: g_cn
1049 INTEGER, INTENT(out) :: c0(0:), cn(0:)
1050
1051 ! convention: process space coords start at 0, grid point coords start at 1
1052
1053 integer :: tn
1054 INTEGER :: d, m
1055 INTEGER :: it
1056
1057 tn = SIZE(c0)
1058 IF (tn<0) CALL test_abort('(tn<0)', filename, __line__)
1059 IF (tn>g_cn) CALL test_abort('regular_deco: too many task for such a core&
1060 & region', &
1061 filename, __line__)
1062
1063 d = g_cn/tn
1064 m = mod(g_cn, tn)
1065
1066 DO it = 0, m-1
1067 cn(it) = d + 1
1068 ENDDO
1069 DO it = m, tn-1
1070 cn(it) = d
1071 ENDDO
1072
1073 c0(0)=0
1074 DO it = 1, tn-1
1075 c0(it) = c0(it-1) + cn(it-1)
1076 ENDDO
1077 IF (c0(tn-1)+cn(tn-1) /= g_cn) &
1078 CALL test_abort('regular_deco: internal error 1', filename, __line__)
1079 END SUBROUTINE regular_deco
1080
1081 FUNCTION run_randomized_tests() RESULT(fully_random_tests)
1082 LOGICAL :: fully_random_tests
1083 CHARACTER(len=32) :: envval
1084 INTEGER :: envlen, envstat
1085 CALL get_environment_variable("YAXT_FULLY_RANDOM_TESTS", envval, envlen, &
1086 status=envstat)
1087 IF (envstat == 0 .AND. (envlen == 1 .OR. envlen == 3)) THEN
1088 IF (envlen == 1 .AND. (envval(1:1) == 'y' .OR. envval(1:1) == 'Y' &
1089 & .OR. envval(1:1) == '1')) THEN
1090 fully_random_tests = .true.
1091 ELSE IF (str2lower(envval(1:3)) == 'yes') THEN
1092 fully_random_tests = .true.
1093 ELSE
1094 fully_random_tests = .false.
1095 END IF
1096 ELSE
1097 fully_random_tests = .false.
1098 END IF
1099 END FUNCTION run_randomized_tests
1100
1101 FUNCTION str2lower(s) RESULT(t)
1102 CHARACTER(len=*), INTENT(in) :: s
1103 CHARACTER(len=LEN(s)) :: t
1104 INTEGER, PARAMETER :: idel = ichar('a')-ichar('A')
1105 INTEGER :: i
1106 DO i = 1, len_trim(s)
1107 t(i:i) = char( ichar(s(i:i)) &
1108 + merge(idel, 0, ichar(s(i:i)) >= ichar('A') &
1109 & .AND. ichar(s(i:i)) <= ichar('Z')))
1110 ENDDO
1111 END FUNCTION str2lower
1112
1113 SUBROUTINE init_fortran_random(full_random)
1114 LOGICAL, INTENT(in) :: full_random
1115 INTEGER, ALLOCATABLE :: rseed(:)
1116
1117 INTEGER :: rseed_size, i
1118 CHARACTER(len=32) :: fmt
1119 INTEGER :: tparts(8), timeseed
1120 INTEGER :: days_per_month(12), days_prefix
1121 INTEGER, PARAMETER :: tparts_mult(7) = (/ &
1122 365 * 24 * 60 * 60, & ! year
1123 0, & ! sum over days_per_month added to day
1124 24 * 60 * 60, & ! day
1125 0, & ! ignore timezone offset
1126 60 * 60, & ! hour of day
1127 60, & ! minute of hour
1128 1 /) ! seconnd
1129 CHARACTER(len=32) :: envval
1130 INTEGER :: envlen, envstat
1131
1132 CALL random_seed(size=rseed_size)
1133 ALLOCATE(rseed(rseed_size))
1134 DO i = 1, rseed_size
1135 rseed(i) = 4711
1136 END DO
1137 IF (full_random) THEN
1138
1139 CALL date_and_time(values=tparts)
1140 days_per_month( 1) = 31
1141 days_per_month( 2) = merge(28, 29, &
1142 mod(tparts(1), 4) == 0 .AND. ( mod(tparts(1), 100) /= 0 &
1143 & .OR. mod(tparts(1), 400) == 0))
1144 days_per_month( 3) = 31
1145 days_per_month( 4) = 30
1146 days_per_month( 5) = 31
1147 days_per_month( 6) = 30
1148 days_per_month( 7) = 31
1149 days_per_month( 8) = 31
1150 days_per_month( 9) = 30
1151 days_per_month(10) = 31
1152 days_per_month(11) = 30
1153 days_per_month(12) = 31
1154 tparts(1) = tparts(1) - 1970
1155 days_prefix = sum(days_per_month(1:tparts(2)-1))
1156 tparts(3) = tparts(3) + days_prefix - 1
1157 tparts(2) = 0
1158 timeseed = sum(tparts(1:7) * tparts_mult)
1159 timeseed = ieor(tparts(8), timeseed) ! mix in microseconds
1160 rseed(1) = timeseed
1161 CALL get_environment_variable("YAXT_RANDOM_SEED", envval, envlen, &
1162 status=envstat)
1163 IF (envstat == 0) THEN
1164 WRITE (fmt, '(a,i0,a)') '(i', digits(rseed), ')'
1165 READ(envval(1:envlen), fmt) rseed(1)
1166 END IF
1167 WRITE(0, '(a,i0)') 'used extra seed=', rseed(1)
1168 FLUSH(0)
1169 END IF
1170 CALL random_seed(put=rseed)
1171 END SUBROUTINE init_fortran_random
1172
1173 SUBROUTINE random_fill_i2(a)
1174 INTEGER(i2), ALLOCATABLE :: a(:)
1175 INTEGER :: n, nb, i, j
1176 INTEGER, PARAMETER :: block_len=8
1177 REAL(c_double) :: rand_nums(block_len)
1178 REAL(c_double) :: sc
1179 n = SIZE(a)
1180 nb = n/block_len
1181 sc = real(huge(a), c_double)
1182 DO j = 1, nb
1183 CALL random_number(rand_nums)
1184 DO i = 1, block_len
1185 a(i+(j-1)*block_len) = nint(rand_nums(i) * sc, i2)
1186 END DO
1187 END DO
1188 CALL random_number(rand_nums)
1189 DO i = 1, mod(n, block_len)
1190 a(i+nb*block_len) = nint(rand_nums(i) * sc, i2)
1191 END DO
1192 END SUBROUTINE random_fill_i2
1193
1194 SUBROUTINE random_fill_i4(a)
1195 INTEGER(i4), ALLOCATABLE :: a(:)
1196 INTEGER :: n, nb, i, j
1197 INTEGER, PARAMETER :: block_len=8
1198 REAL(c_double) :: rand_nums(block_len)
1199 REAL(c_double) :: sc
1200 n = SIZE(a)
1201 nb = n/block_len
1202 sc = real(huge(a), c_double)
1203 DO j = 1, nb
1204 CALL random_number(rand_nums)
1205 DO i = 1, block_len
1206 a(i+(j-1)*block_len) = nint(rand_nums(i) * sc, i4)
1207 END DO
1208 END DO
1209 CALL random_number(rand_nums)
1210 DO i = 1, mod(n, block_len)
1211 a(i+nb*block_len) = nint(rand_nums(i) * sc, i4)
1212 END DO
1213 END SUBROUTINE random_fill_i4
1214
1215 SUBROUTINE random_fill_i8(a)
1216 INTEGER(i8), ALLOCATABLE :: a(:)
1217 INTEGER :: n, nb, i, j
1218 INTEGER, PARAMETER :: block_len=8
1219 ! the least signifant digits of huge(a) do not fit into an ieee754
1220 ! double precision real and gfortran warns about that. For that
1221 ! reason the below expression limits scale_val to those 1-bits
1222 ! which can be represented in a double precision constant
1223#if ! defined __PGI || __PGIC__ > 24
1224 INTEGER(i8), PARAMETER :: scale_val &
1225 = iand(huge(a), not(ishft(1_i8, bit_size(1_i8) &
1226 & - digits(0.0_c_double)) - 1 ) )
1227#else
1228 INTEGER(i8) :: scale_val
1229#endif
1230 REAL(c_double) :: rand_nums(block_len), sc
1231 n = SIZE(a)
1232 nb = n/block_len
1233#if defined __PGI && __PGIC__ <= 24
1234 scale_val &
1235 = iand(huge(a), not(ishft(1_i8, bit_size(1_i8) &
1236 & - digits(0.0_c_double)) - 1 ) )
1237#endif
1238
1239 sc = real(scale_val, c_double)
1240 DO j = 1, nb
1241 CALL random_number(rand_nums)
1242 DO i = 1, block_len
1243 a(i+(j-1)*block_len) = nint(rand_nums(i) * sc, i8)
1244 END DO
1245 END DO
1246 CALL random_number(rand_nums)
1247 DO i = 1, mod(n, block_len)
1248 a(i+nb*block_len) = nint(rand_nums(i) * sc, i8)
1249 END DO
1250 END SUBROUTINE random_fill_i8
1251
1252 FUNCTION memcrc_i2(a) RESULT(crcval)
1253 INTEGER(i2), TARGET, INTENT(in) :: a(:)
1254 INTEGER(c_int32_t) :: crcval
1255 TYPE(c_ptr) :: p
1256 INTEGER(c_size_t) :: n
1257 xt_slice_c_loc(a, p)
1258 n = SIZE(a) * 2_c_size_t
1259 crcval = memcrc(p, n)
1260 END FUNCTION memcrc_i2
1261
1262 FUNCTION memcrc_i4(a) RESULT(crcval)
1263 INTEGER(i4), TARGET, INTENT(in) :: a(:)
1264 INTEGER(c_int32_t) :: crcval
1265 TYPE(c_ptr) :: p
1266 INTEGER(c_size_t) :: n
1267 xt_slice_c_loc(a, p)
1268 n = SIZE(a) * 4_c_size_t
1269 crcval = memcrc(p, n)
1270 END FUNCTION memcrc_i4
1271
1272 FUNCTION memcrc_i8(a) RESULT(crcval)
1273 INTEGER(i8), TARGET, INTENT(in) :: a(:)
1274 INTEGER(c_int32_t) :: crcval
1275 TYPE(c_ptr) :: p
1276 INTEGER(c_size_t) :: n
1277 xt_slice_c_loc(a, p)
1278 n = SIZE(a) * 8_c_size_t
1279 crcval = memcrc(p, n)
1280 END FUNCTION memcrc_i8
1281
1282 SUBROUTINE permute_i2(a, p)
1283 INTEGER(i2), INTENT(inout) :: a(:)
1284 INTEGER, INTENT(inout) :: p(:)
1285
1286 INTEGER :: i, n, next, temp
1287 INTEGER(i2) :: t
1288
1289 n = SIZE(a)
1290 DO i = 1, n
1291 next = i
1292 ! Check if it is already
1293 ! considered in cycle
1294 DO WHILE (p(next) >= 0)
1295 ! Swap the current element according
1296 ! to the permutation in P
1297 t = a(i)
1298 temp = p(next)
1299 a(i) = a(temp)
1300 a(temp) = t
1301 p(next) = -temp
1302 next = temp
1303 END DO
1304 END DO
1305 END SUBROUTINE permute_i2
1306
1307 SUBROUTINE permute_i4(a, p)
1308 INTEGER(i4), INTENT(inout) :: a(:)
1309 INTEGER, INTENT(inout) :: p(:)
1310
1311 INTEGER :: i, n, next, temp
1312 INTEGER(i4) :: t
1313
1314 n = SIZE(a)
1315 DO i = 1, n
1316 next = i
1317 ! Check if it is already
1318 ! considered in cycle
1319 DO WHILE (p(next) >= 0)
1320 ! Swap the current element according
1321 ! to the permutation in P
1322 t = a(i)
1323 temp = p(next)
1324 a(i) = a(temp)
1325 a(temp) = t
1326 p(next) = -temp
1327 next = temp
1328 END DO
1329 END DO
1330 END SUBROUTINE permute_i4
1331
1332 SUBROUTINE permute_i8(a, p)
1333 INTEGER(i8), INTENT(inout) :: a(:)
1334 INTEGER, INTENT(inout) :: p(:)
1335
1336 INTEGER :: i, n, next, temp
1337 INTEGER(i8) :: t
1338
1339 n = SIZE(a)
1340 DO i = 1, n
1341 next = i
1342 ! Check if it is already
1343 ! considered in cycle
1344 DO WHILE (p(next) >= 0)
1345 ! Swap the current element according
1346 ! to the permutation in P
1347 t = a(i)
1348 temp = p(next)
1349 a(i) = a(temp)
1350 a(temp) = t
1351 p(next) = -temp
1352 next = temp
1353 END DO
1354 END DO
1355 END SUBROUTINE permute_i8
1356
1357END MODULE ftest_common
1358!
1359! Local Variables:
1360! f90-continuation-indent: 5
1361! coding: utf-8
1362! indent-tabs-mode: nil
1363! show-trailing-whitespace: t
1364! require-trailing-newline: t
1365! End:
1366!
uint32_t SymPrefix memcrc(const unsigned char *b, size_t n)
Definition cksum.c:115