1! @file benchmark_hdf5_f90.F90
2! @author M. Scot Breitenfeld <brtnfld@hdfgroup.org>
3! @version 0.1
4!
5! @section LICENSE
6! BSD style license
7!
8! @section DESCRIPTION
9! F90 benchmarking program for pcgns library
10
11MODULE testing_functions
12  !
13  ! Contains functions to verify values
14  !
15  INTERFACE check_eq
16     MODULE PROCEDURE c_float_eq, c_double_eq, c_long_eq, c_long_long_eq
17  END INTERFACE
18
19CONTAINS
20
21  LOGICAL FUNCTION c_float_eq(a,b)
22    IMPLICIT NONE
23    ! Check if two C_FLOAT reals are equivalent
24    REAL*4, INTENT(IN):: a,b
25    REAL*4, PARAMETER :: eps = 1.e-8
26    c_float_eq = ABS(a-b) .LT. eps
27  END FUNCTION c_float_eq
28
29  LOGICAL FUNCTION c_double_eq(a,b)
30    IMPLICIT NONE
31    ! Check if two C_DOUBLE reals are equivalent
32    REAL*8, INTENT(IN):: a,b
33    REAL*8, PARAMETER :: eps = 1.e-8
34    c_double_eq = ABS(a-b) .LT. eps
35  END FUNCTION c_double_eq
36
37  LOGICAL FUNCTION c_long_eq(a,b)
38    IMPLICIT NONE
39    ! Check if two C_LONG integers are equivalent
40    INTEGER*4, INTENT(IN):: a,b
41    c_long_eq = a-b .EQ. 0
42  END FUNCTION c_long_eq
43
44  LOGICAL FUNCTION c_long_long_eq(a,b)
45    IMPLICIT NONE
46    ! Check if two C_LONG_LONG integers are equivalent
47    INTEGER*8, INTENT(IN):: a,b
48    c_long_long_eq = a-b .EQ. 0
49  END FUNCTION c_long_long_eq
50
51END MODULE testing_functions
52
53PROGRAM benchmark_hdf5_f90
54
55  USE mpi
56  USE ISO_C_BINDING
57  USE testing_functions
58  USE cgns
59  IMPLICIT NONE
60
61#ifdef WINNT
62  INCLUDE 'cgnswin_f.h'
63#endif
64
65#include "cgnstypes_f03.h"
66
67  INTEGER, PARAMETER :: dp = KIND(1.d0)
68  ! Use powers of 2
69  INTEGER(CGSIZE_T), PARAMETER :: Nelem = 65536 ! 33554432 ! Use multiples of number of cores per node
70  INTEGER(CGSIZE_T), PARAMETER :: NodePerElem = 6
71
72  INTEGER(CGSIZE_T) :: Nnodes
73  INTEGER(C_INT) :: comm_size, comm_rank, comm_info, mpi_err
74  INTEGER :: err
75  INTEGER :: fn
76  INTEGER :: B
77  INTEGER :: Z
78  INTEGER :: S
79  INTEGER :: Cx,Cy,Cz, Fx, Fy, Fz, Ar, Ai
80  INTEGER :: cell_dim = 3
81  INTEGER :: phys_dim = 3
82  INTEGER :: r_cell_dim = 0
83  INTEGER :: r_phys_dim = 0
84  INTEGER(CGSIZE_T), DIMENSION(1:3) :: nijk, sizes
85  INTEGER(CGSIZE_T) :: min, max
86  INTEGER(CGSIZE_T) :: k, count
87  ! For writing and reading data
88  REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: Coor_x
89  REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: Coor_y
90  REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: Coor_z
91  REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: Data_Fx
92  REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: Data_Fy
93  REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: Data_Fz
94  REAL(KIND=dp), DIMENSION(:), ALLOCATABLE :: Array_r
95  INTEGER(CGSIZE_T), DIMENSION(:), ALLOCATABLE :: Array_i
96  INTEGER(CGSIZE_T) :: start, iend, emin, emax
97  INTEGER(CGSIZE_T), DIMENSION(:), ALLOCATABLE :: elements
98  CHARACTER(LEN=32) :: fname, name
99  CHARACTER(LEN=180) :: bname, zname
100  INTEGER :: indx_null
101  LOGICAL :: queue, debug
102  REAL(KIND=dp) t0, t1, t2
103
104  ! Timing storage convention:
105  ! timing(0) = Total program time
106  ! timing(1) = Time to write nodal coordinates
107  ! timing(2) = Time to write connectivity table
108  ! timing(3) = Time to write solution data (field data)
109  ! timing(4) = Time to write array data
110  ! timing(5) = Time to read nodal coordinates
111  ! timing(6) = Time to read connectivity table
112  ! timing(7) = Time to read solution data (field data)
113  ! timing(8) = Time to read array data
114  ! timing(9) = Time for cgp_open, CG_MODE_WRITE
115  ! timing(10) = Time for cg_base_write
116  ! timing(11) = Time for cg_zone_write
117  ! timing(12) = Time for cgp_open, CG_MODE_READ
118  ! timing(13) = Time for cg_base_read
119  ! timing(14) = Time for cg_zone_read
120
121  REAL(KIND=dp), DIMENSION(1:15) :: xtiming, timing, timingMin, timingMax
122
123  INTEGER :: ierr
124  INTEGER :: istat
125  CHARACTER(LEN=6) :: ichr6
126  INTEGER, DIMENSION(1:3) :: Cvec, Fvec
127  INTEGER, DIMENSION(1:2) :: Avec
128
129  CALL MPI_Init(mpi_err)
130  CALL MPI_Comm_size(MPI_COMM_WORLD,comm_size,mpi_err)
131  CALL MPI_Comm_rank(MPI_COMM_WORLD,comm_rank,mpi_err)
132
133  CALL MPI_Info_Create(comm_info,mpi_err)
134  ! set this to what your GPFS block size actually is
135  CALL MPI_Info_set(comm_info, "striping_unit", "8388608",mpi_err)
136
137  WRITE(ichr6,'(I6.6)') comm_size
138
139  ! parameters
140  queue = .FALSE.
141  debug = .FALSE.
142
143  t0 = MPI_Wtime()
144
145  CALL cgp_mpi_info_f(comm_info, ierr)
146  CALL cgp_pio_mode_f(CGP_COLLECTIVE, ierr)
147
148  Nnodes = Nelem*NodePerElem
149
150  nijk(1) = Nnodes ! Number of vertices
151  nijk(2) = Nelem  ! Number of cells
152  nijk(3) = 0      ! Number of boundary vertices
153
154  ! ======================================
155  ! ==    **WRITE THE CGNS FILE **      ==
156  ! ======================================
157
158  t1 = MPI_Wtime()
159  CALL cgp_open_f("benchmark_"//ichr6//".cgns", CG_MODE_WRITE, fn, err)
160  IF(err.NE.CG_OK)THEN
161     PRINT*,'*FAILED* cgp_open_f'
162     CALL cgp_error_exit_f()
163  ENDIF
164  t2 = MPI_Wtime()
165  xtiming(10) = t2-t1
166
167  t1 = MPI_Wtime()
168  CALL cg_base_write_f(fn, "Base 1", cell_dim, phys_dim, B, err)
169  IF(err.NE.CG_OK)THEN
170     PRINT*,'*FAILED* cg_base_write_f'
171     CALL cgp_error_exit_f()
172  ENDIF
173  t2 = MPI_Wtime()
174  xtiming(11) = t2-t1
175
176  t1 = MPI_Wtime()
177  CALL cg_zone_write_f(fn, B, "Zone 1", nijk, Unstructured, Z, err)
178  IF(err.NE.CG_OK)THEN
179     PRINT*,'*FAILED* cg_zone_write_f'
180     CALL cgp_error_exit_f()
181  ENDIF
182  t2 = MPI_Wtime()
183  xtiming(12) = t2-t1
184
185  ! ======================================
186  ! == (A) WRITE THE NODAL COORDINATES  ==
187  ! ======================================
188
189  count = nijk(1)/comm_size
190
191  ALLOCATE(Coor_x(1:count), STAT = istat)
192  IF (istat.NE.0)THEN
193     PRINT*, '*FAILED* allocation of Coor_x'
194     CALL cgp_error_exit_f()
195  ENDIF
196  ALLOCATE(Coor_y(1:count), STAT = istat)
197  IF (istat.NE.0)THEN
198     PRINT*, '*FAILED* allocation of Coor_y'
199     CALL cgp_error_exit_f()
200  ENDIF
201  ALLOCATE(Coor_z(1:count), STAT = istat)
202  IF (istat.NE.0)THEN
203     PRINT*, '*FAILED* allocation of Coor_z'
204     CALL cgp_error_exit_f()
205  ENDIF
206
207  min = count*comm_rank+1
208  max = count*(comm_rank+1)
209
210  DO k = 1, count
211     Coor_x(k) = REAL(comm_rank*count + k, KIND=dp) + 0.1_dp
212     Coor_y(k) = Coor_x(k) + 0.1_dp
213     Coor_z(k) = Coor_y(k) + 0.1_dp
214  ENDDO
215
216  CALL cgp_coord_write_f(fn,B,Z,RealDouble,"CoordinateX",Cx,err)
217  IF(err.NE.CG_OK)THEN
218     PRINT*,'*FAILED* cgp_coord_write_f (Coor_X)'
219     CALL cgp_error_exit_f()
220  ENDIF
221  CALL cgp_coord_write_f(fn,B,Z,RealDouble,"CoordinateY",Cy,err)
222  IF(err.NE.CG_OK)THEN
223     PRINT*,'*FAILED* cgp_coord_write_f (Coor_Y)'
224     CALL cgp_error_exit_f()
225  ENDIF
226  CALL cgp_coord_write_f(fn,B,Z,RealDouble,"CoordinateZ",Cz,err)
227  IF(err.NE.CG_OK)THEN
228     PRINT*,'*FAILED* cgp_coord_write_f (Coor_Z)'
229     CALL cgp_error_exit_f()
230  ENDIF
231
232  t1 = MPI_Wtime()
233#if HDF5_HAVE_MULTI_DATASETS
234  Cvec(1:3) = (/Cx,Cy,Cz/)
235  CALL cgp_coord_multi_write_data_f(fn,B,Z,Cvec,min,max,Coor_x,Coor_y,Coor_z,err)
236  IF(err.NE.CG_OK)THEN
237     PRINT*,'*FAILED* cgp_coord_multi_write_data_f'
238     CALL cgp_error_exit_f()
239  ENDIF
240#else
241  CALL cgp_coord_write_data_f(fn,B,Z,Cx,min,max,Coor_x,err)
242  IF(err.NE.CG_OK)THEN
243     PRINT*,'*FAILED* cgp_coord_write_data_f (Coor_x)'
244     CALL cgp_error_exit_f()
245  ENDIF
246  CALL cgp_coord_write_data_f(fn,B,Z,Cy,min,max,Coor_y,err)
247  IF(err.NE.CG_OK)THEN
248     PRINT*,'*FAILED* cgp_coord_write_data_f (Coor_y)'
249     CALL cgp_error_exit_f()
250  ENDIF
251  CALL cgp_coord_write_data_f(fn,B,Z,Cz,min,max,Coor_z,err)
252  IF(err.NE.CG_OK)THEN
253     PRINT*,'*FAILED* cgp_coord_write_data_f (Coor_z)'
254     CALL cgp_error_exit_f()
255  ENDIF
256#endif
257  t2 = MPI_Wtime()
258  xtiming(2) = t2-t1
259
260  DEALLOCATE(Coor_x, Coor_y, Coor_z)
261
262  ! ======================================
263  ! == (B) WRITE THE CONNECTIVITY TABLE ==
264  ! ======================================
265
266  start = 1
267  iend = nijk(2)
268  CALL cgp_section_write_f(fn,B,Z,"Elements",PENTA_6,start,iend,0,S,err)
269  IF(err.NE.CG_OK)THEN
270     PRINT*,'*FAILED* cgp_section_write_f'
271     CALL cgp_error_exit_f()
272  ENDIF
273
274  count = nijk(2)/comm_size
275  ALLOCATE(elements(1:count*NodePerElem), STAT = istat)
276  IF (istat.NE.0)THEN
277     PRINT*, '*FAILED* allocation of elements'
278     CALL cgp_error_exit_f()
279  ENDIF
280
281  ! Create ridiculous connectivity table ...
282  DO k = 1, count*NodePerElem
283     elements(k) = comm_rank*count*NodePerElem + k
284  ENDDO
285
286  emin = count*comm_rank+1
287  emax = count*(comm_rank+1)
288
289  t1 = MPI_Wtime()
290  CALL cgp_elements_write_data_f(fn, B, Z, S, emin, emax, elements,err)
291  IF(err.NE.CG_OK)THEN
292     PRINT*,'*FAILED* cgp_elements_write_data_f (elements)'
293     CALL cgp_error_exit_f()
294  ENDIF
295
296  t2 = MPI_Wtime()
297  xtiming(3) = t2-t1
298
299  DEALLOCATE(elements)
300
301  ! ======================================
302  ! == (C) WRITE THE FIELD DATA         ==
303  ! ======================================
304
305  count = nijk(1)/comm_size
306
307  ALLOCATE(Data_Fx(1:count), STAT = istat)
308  IF (istat.NE.0)THEN
309     PRINT*, '*FAILED* allocation of Data_Fx'
310     CALL cgp_error_exit_f()
311  ENDIF
312  ALLOCATE(Data_Fy(1:count), STAT = istat)
313  IF (istat.NE.0)THEN
314     PRINT*, '*FAILED* allocation of Data_Fy'
315     CALL cgp_error_exit_f()
316  ENDIF
317  ALLOCATE(Data_Fz(1:count), STAT = istat)
318  IF (istat.NE.0)THEN
319     PRINT*, '*FAILED* allocation of Data_Fz'
320     CALL cgp_error_exit_f()
321  ENDIF
322
323  DO k = 1, count
324     Data_Fx(k) = REAL(comm_rank*count+k, KIND=dp) + .01_dp
325     Data_Fy(k) = REAL(comm_rank*count+k, KIND=dp) + .02_dp
326     Data_Fz(k) = REAL(comm_rank*count+k, KIND=dp) + .03_dp
327  ENDDO
328
329  CALL cg_sol_write_f(fn, B, Z, "Solution", Vertex, S, err)
330  IF(err.NE.CG_OK)THEN
331     PRINT*,'*FAILED* cg_sol_write_f'
332     CALL cgp_error_exit_f()
333  ENDIF
334
335  CALL cgp_field_write_f(fn,B,Z,S,RealDouble,"MomentumX",Fx, err)
336  IF(err.NE.CG_OK)THEN
337     PRINT*,'*FAILED* cgp_field_write  (MomentumX)'
338     CALL cgp_error_exit_f()
339  ENDIF
340  CALL cgp_field_write_f(fn,B,Z,S,RealDouble,"MomentumY",Fy, err)
341  IF(err.NE.CG_OK)THEN
342     PRINT*,'*FAILED* cgp_field_write (MomentumY)'
343     CALL cgp_error_exit_f()
344  ENDIF
345  CALL cgp_field_write_f(fn,B,Z,S,RealDouble,"MomentumZ",Fz, err)
346  IF(err.NE.CG_OK)THEN
347     PRINT*,'*FAILED* cgp_field_write (MomentumZ)'
348     CALL cgp_error_exit_f()
349  ENDIF
350
351  t1 = MPI_Wtime()
352#if HDF5_HAVE_MULTI_DATASETS
353  Fvec(1:3) = (/Fx,Fy,Fz/)
354  CALL cgp_field_multi_write_data_f(fn,B,Z,S,Fvec,min,max,err,3,Data_Fx,Data_Fy,Data_Fz)
355  IF(err.NE.CG_OK)THEN
356     PRINT*,'*FAILED* cgp_field_multi_write_data_f'
357     CALL cgp_error_exit_f()
358  ENDIF
359#else
360  call cgp_field_write_data_f(fn,B,Z,S,Fx,min,max,Data_Fx, err)
361  IF(err.NE.CG_OK)THEN
362     PRINT*,'*FAILED* cgp_field_write_data (Data_Fx)'
363     CALL cgp_error_exit_f()
364  ENDIF
365  call cgp_field_write_data_f(fn,B,Z,S,Fy,min,max,Data_Fy, err)
366  IF(err.NE.CG_OK)THEN
367     PRINT*,'*FAILED* cgp_field_write_data (Data_Fy)'
368     CALL cgp_error_exit_f()
369  ENDIF
370  call cgp_field_write_data_f(fn,B,Z,S,Fz,min,max,Data_Fz, err)
371  IF(err.NE.CG_OK)THEN
372     PRINT*,'*FAILED* cgp_field_write_data (Data_Fz)'
373     CALL cgp_error_exit_f()
374  ENDIF
375#endif
376  t2 = MPI_Wtime()
377  xtiming(4) = t2-t1
378
379  DEALLOCATE(Data_Fx,Data_Fy,Data_Fz)
380
381  ! ======================================
382  ! == (D) WRITE THE ARRAY DATA         ==
383  ! ======================================
384
385  count = nijk(1)/comm_size
386
387  ALLOCATE(Array_r(1:count), STAT = istat)
388  IF (istat.NE.0)THEN
389     PRINT*, '*FAILED* allocation of Array_r'
390     CALL cgp_error_exit_f()
391  ENDIF
392  ALLOCATE(Array_i(1:count), STAT = istat)
393  IF (istat.NE.0)THEN
394     PRINT*, '*FAILED* allocation of Array_i'
395     CALL cgp_error_exit_f()
396  ENDIF
397
398  min = count*comm_rank+1
399  max = count*(comm_rank+1)
400
401  DO k = 1, count
402     Array_r(k) = REAL(comm_rank*count + k, KIND=dp) + .001_dp
403     Array_i(k) = comm_rank*count*NodePerElem + k
404  ENDDO
405
406  CALL cg_goto_f(fn, B, err, 'Zone 1', 0, 'end')
407  IF(err.NE.CG_OK)THEN
408     PRINT*,'*FAILED* cg_goto_f'
409     CALL cgp_error_exit_f()
410  ENDIF
411  CALL cg_user_data_write_f("User Data", err)
412  IF(err.NE.CG_OK)THEN
413     PRINT*,'*FAILED* cg_user_data_write_f'
414     CALL cgp_error_exit_f()
415  ENDIF
416
417  CALL cg_gorel_f(fn,err,'User Data',0,'end')
418  IF(err.NE.CG_OK)THEN
419     PRINT*,'*FAILED* cg_gorel_f'
420     CALL cgp_error_exit_f()
421  ENDIF
422
423  CALL cgp_array_write_f("ArrayR",cg_get_type(Array_r(1)),1,INT(nijk(1),cgsize_t),Ar, err)
424  IF(err.NE.CG_OK)THEN
425     PRINT*,'*FAILED* cgp_array_write_f (Array_Ar)'
426     CALL cgp_error_exit_f()
427  ENDIF
428
429  CALL cgp_array_write_f("ArrayI",cg_get_type(Array_i(1)),1,INT(nijk(1),cgsize_t),Ai, err)
430  IF(err.NE.CG_OK)THEN
431     PRINT*,'*FAILED* cgp_array_write_f   (Array_Ai)'
432     CALL cgp_error_exit_f()
433  ENDIF
434
435  t1 = MPI_Wtime()
436
437#if HDF5_HAVE_MULTI_DATASETS
438  Avec = (/Ai,Ar/)
439  CALL cgp_array_multi_write_data_f(fn,Avec,min,max,err,2,Array_i,Array_r)
440  IF(err.NE.CG_OK)THEN
441     PRINT*,'*FAILED* cgp_array_multi_write_data_f'
442     CALL cgp_error_exit_f()
443  ENDIF
444#else
445  CALL cgp_array_write_data_f(Ai,min,max,Array_i, err)
446  IF(err.NE.CG_OK)THEN
447     PRINT*,'*FAILED* cgp_array_write_data_f  (Array_Ai)'
448     CALL cgp_error_exit_f()
449  ENDIF
450  call cgp_array_write_data_f(Ar,min,max,Array_r, err)
451  IF(err.NE.CG_OK)THEN
452     PRINT*,'*FAILED* cgp_array_write_data_f (Array_Ar)'
453     CALL cgp_error_exit_f()
454  ENDIF
455#endif
456  t2 = MPI_Wtime()
457  xtiming(5) = t2-t1
458
459  DEALLOCATE(Array_r,Array_i)
460
461  CALL cgp_close_f(fn, err)
462  IF(err.NE.CG_OK)THEN
463     PRINT*,'*FAILED* cgp_close_f'
464     CALL cgp_error_exit_f()
465  ENDIF
466
467  ! ======================================
468  ! ==    **  READ THE CGNS FILE **     ==
469  ! ======================================
470
471  CALL MPI_Barrier(MPI_COMM_WORLD, mpi_err)
472
473  t1 = MPI_Wtime()
474  ! Open the cgns file for reading
475  CALL cgp_open_f("benchmark_"//ichr6//".cgns", CG_MODE_MODIFY, fn, err)
476  IF(err.NE.CG_OK)THEN
477     PRINT*,'*FAILED* cgp_open_f'
478     CALL cgp_error_exit_f()
479  ENDIF
480  t2 = MPI_Wtime()
481  xtiming(13) = t2-t1
482
483  ! Read the base information
484  t1 = MPI_Wtime()
485  CALL cg_base_read_f(fn, B, bname, r_cell_dim, r_phys_dim, err)
486  IF(err.NE.CG_OK)THEN
487     PRINT*,'*FAILED* cg_base_read_f'
488     CALL cgp_error_exit_f()
489  ENDIF
490  t2 = MPI_Wtime()
491  xtiming(14) = t2-t1
492  IF(r_cell_dim.NE.cell_dim.OR.r_phys_dim.NE.phys_dim)THEN
493     WRITE(*,'(2(A,I0))') '*FAILED* bad cell dim = ',r_cell_dim,'or phy dim =',r_phys_dim
494     CALL cgp_error_exit_f()
495  ENDIF
496
497  IF(TRIM(bname).NE.'Base 1')THEN
498     WRITE(*,'(A,A)') '*FAILED* bad base name = ',TRIM(bname)
499     CALL cgp_error_exit_f()
500  ENDIF
501
502  ! Read the zone information
503  t1 = MPI_Wtime()
504  CALL cg_zone_read_f(fn, B, Z, zname, sizes, err)
505  IF(err.NE.CG_OK) PRINT*,'*FAILED* cgp_zone_read_f'
506  t2 = MPI_Wtime()
507  xtiming(15) = t2-t1
508
509  ! Check the read zone information is correct
510  IF(sizes(1).NE.Nnodes)THEN
511     WRITE(*,'(A,I0)') '*FAILED* bad num points = ',sizes(1)
512     CALL cgp_error_exit_f()
513  ENDIF
514
515  IF(sizes(2).NE.Nelem)THEN
516     WRITE(*,'(A,I0)') '*FAILED* bad num elements = ',sizes(2)
517     CALL cgp_error_exit_f()
518  ENDIF
519  IF(sizes(3).NE.0)THEN
520     WRITE(*,'(A,I0)') '*FAILED* bad num elements = ',sizes(3)
521     CALL cgp_error_exit_f()
522  ENDIF
523
524  IF(TRIM(zname).NE.'Zone 1')THEN
525     WRITE(*,'(A,A)') '*FAILED* bad zone name = ',TRIM(zname)
526     CALL cgp_error_exit_f()
527  ENDIF
528
529
530  ! ======================================
531  ! ==  (A) READ THE NODAL COORDINATES  ==
532  ! ======================================
533
534  count = nijk(1)/comm_size
535  ALLOCATE(Coor_x(1:count), STAT = istat)
536  IF (istat.NE.0)THEN
537     PRINT*, '*FAILED* allocation of Coor_x'
538     CALL cgp_error_exit_f()
539  ENDIF
540  ALLOCATE(Coor_y(1:count), STAT = istat)
541  IF (istat.NE.0)THEN
542     PRINT*, '*FAILED* allocation of Coor_y'
543     CALL cgp_error_exit_f()
544  ENDIF
545  ALLOCATE(Coor_z(1:count), STAT = istat)
546  IF (istat.NE.0)THEN
547     PRINT*, '*FAILED* allocation of Coor_z'
548     CALL cgp_error_exit_f()
549  ENDIF
550
551  min = count*comm_rank+1
552  max = count*(comm_rank+1)
553
554  t1 = MPI_Wtime()
555#if HDF5_HAVE_MULTI_DATASETS
556  Cvec(1:3) = (/Cx,Cy,Cz/)
557  CALL cgp_coord_multi_read_data_f(fn,B,Z,Cvec,min,max,Coor_x,Coor_y,Coor_z,err)
558  IF(err.NE.CG_OK)THEN
559     PRINT*,'*FAILED* cgp_coord_multi_read_data_f'
560     CALL cgp_error_exit_f()
561  ENDIF
562#else
563  CALL cgp_coord_read_data_f(fn,B,Z,Cx,min,max,Coor_x,err)
564  IF(err.NE.CG_OK)THEN
565     PRINT*,'*FAILED* cgp_coord_read_data_f (Reading Coor_x)'
566     CALL cgp_error_exit_f()
567  ENDIF
568  CALL cgp_coord_read_data_f(fn,B,Z,Cy,min,max,Coor_y,err)
569  IF(err.NE.CG_OK)THEN
570     PRINT*,'*FAILED* cgp_coord_read_data_f (Reading Coor_y)'
571     CALL cgp_error_exit_f()
572  ENDIF
573  CALL cgp_coord_read_data_f(fn,B,Z,Cz,min,max,Coor_z,err)
574  IF(err.NE.CG_OK)THEN
575     PRINT*,'*FAILED* cgp_coord_read_data_f (Reading Coor_z)'
576     CALL cgp_error_exit_f()
577  ENDIF
578#endif
579  t2 = MPI_Wtime()
580  xtiming(6) = t2-t1
581
582  ! Check if read the data back correctly
583  IF(debug)THEN
584     DO k = 1, count
585        IF(.NOT.check_eq(Coor_x(k), REAL(comm_rank*count + k, KIND=DP) + 0.1_DP).OR. &
586             .NOT.check_eq(Coor_y(k), REAL(comm_rank*count + k, KIND=DP) + 0.2_DP).OR. &
587             .NOT.check_eq(Coor_z(k), REAL(comm_rank*count + k, KIND=DP) + 0.3_DP)) THEN
588           PRINT*,'*FAILED* cgp_coord_read_data values are incorrect'
589           CALL cgp_error_exit_f()
590        ENDIF
591     ENDDO
592  ENDIF
593
594  DEALLOCATE(Coor_x, Coor_y, Coor_z)
595
596  ! ======================================
597  ! == (B) READ THE CONNECTIVITY TABLE  ==
598  ! ======================================
599
600  count = nijk(2)/comm_size
601  ALLOCATE(elements(1:count*NodePerElem), STAT = istat)
602  IF (istat.NE.0)THEN
603     PRINT*, '*FAILED* allocation of elements'
604     CALL cgp_error_exit_f()
605  ENDIF
606
607  emin = count*comm_rank+1
608  emax = count*(comm_rank+1)
609
610  t1 = MPI_Wtime()
611  CALL cgp_elements_read_data_f(fn, B, Z, S, emin, emax, elements, err)
612  IF(err.NE.CG_OK)THEN
613     PRINT*,'*FAILED* cgp_elements_read_data_f ( Reading elements)'
614     CALL cgp_error_exit_f()
615  ENDIF
616  t2 = MPI_Wtime()
617  xtiming(7) = t2-t1
618  IF(debug)THEN
619     DO k = 1, count
620        IF(.NOT.check_eq(elements(k), comm_rank*count*NodePerElem + k)) THEN
621           PRINT*,'*FAILED* cgp_elements_read_data values are incorrect'
622           CALL cgp_error_exit_f()
623        ENDIF
624     ENDDO
625  ENDIF
626
627  DEALLOCATE(elements)
628
629  ! ======================================
630  ! == (C) READ THE FIELD DATA          ==
631  ! ======================================
632  count = nijk(1)/comm_size
633
634  ALLOCATE(Data_Fx(1:count), STAT = istat)
635  IF (istat.NE.0)THEN
636     PRINT*, '*FAILED* allocation of Data_Fx'
637     CALL cgp_error_exit_f()
638  ENDIF
639  ALLOCATE(Data_Fy(1:count), STAT = istat)
640  IF (istat.NE.0)THEN
641     PRINT*, '*FAILED* allocation of Data_Fy'
642     CALL cgp_error_exit_f()
643  ENDIF
644  ALLOCATE(Data_Fz(1:count), STAT = istat)
645  IF (istat.NE.0)THEN
646     PRINT*, '*FAILED* allocation of Data_Fz'
647     CALL cgp_error_exit_f()
648  ENDIF
649
650  t1 = MPI_Wtime()
651#if HDF5_HAVE_MULTI_DATASETS
652  Fvec(1:3) = (/Fx,Fy,Fz/)
653  CALL cgp_field_multi_read_data_f(fn,B,Z,S,Fvec,min,max,err,3,Data_Fx,Data_Fy,Data_Fz)
654  IF(err.NE.CG_OK)THEN
655     PRINT*,'*FAILED* cgp_field_multi_read_data_f'
656     CALL cgp_error_exit_f()
657  ENDIF
658#else
659  CALL cgp_field_read_data_f(fn,B,Z,S,Fx,min,max,Data_Fx,err)
660  IF(err.NE.CG_OK)THEN
661     PRINT*,'*FAILED* cgp_field_read_data (Data_Fx)'
662     CALL cgp_error_exit_f()
663  ENDIF
664  CALL cgp_field_read_data_f(fn,B,Z,S,Fy,min,max,Data_Fy,err)
665  IF(err.NE.CG_OK)THEN
666     PRINT*,'*FAILED* cgp_field_read_data (Data_Fy)'
667     CALL cgp_error_exit_f()
668  ENDIF
669  CALL cgp_field_read_data_f(fn,B,Z,S,Fz,min,max,Data_Fz,err)
670  IF(err.NE.CG_OK)THEN
671     PRINT*,'*FAILED* cgp_field_read_data (Data_Fz)'
672     CALL cgp_error_exit_f()
673  ENDIF
674#endif
675  t2 = MPI_Wtime()
676  xtiming(8) = t2-t1
677
678  ! Check if read the data back correctly
679  IF(debug)THEN
680     DO k = 1, count
681        IF(.NOT.check_eq(Data_Fx(k), REAL(comm_rank*count + k, KIND=dp) + 0.01_dp).OR. &
682             .NOT.check_eq(Data_Fy(k), REAL(comm_rank*count + k, KIND=dp) + 0.02_dp).OR. &
683             .NOT.check_eq(Data_Fz(k), REAL(comm_rank*count + k, KIND=dp) + 0.03_dp)) THEN
684           PRINT*,'*FAILED* cgp_field_read_data_f values are incorrect'
685           CALL cgp_error_exit_f()
686        ENDIF
687     ENDDO
688  ENDIF
689
690  DEALLOCATE(Data_Fx,Data_Fy,Data_Fz)
691
692  ! ======================================
693  ! == (D) READ THE ARRAY DATA          ==
694  ! ======================================
695
696  count = nijk(1)/comm_size
697
698  ALLOCATE(Array_r(1:count), STAT = istat)
699  IF (istat.NE.0)THEN
700     PRINT*, '*FAILED* allocation of Array_r'
701     CALL cgp_error_exit_f()
702  ENDIF
703  ALLOCATE(Array_i(1:count), STAT = istat)
704  IF (istat.NE.0)THEN
705     PRINT*, '*FAILED* allocation of Array_i'
706     CALL cgp_error_exit_f()
707  ENDIF
708
709  min = count*comm_rank+1
710  max = count*(comm_rank+1)
711
712  CALL cg_goto_f(fn, B, err, "Zone_t",Z,"UserDefinedData_t",1,'end')
713  IF(err.NE.CG_OK)THEN
714     PRINT*,'*FAILED* cg_goto (User Defined Data)'
715     CALL cgp_error_exit_f()
716  ENDIF
717
718  t1 = MPI_Wtime()
719#if HDF5_HAVE_MULTI_DATASETS
720  Avec = (/Ai,Ar/)
721  CALL cgp_array_multi_read_data_f(fn,Avec,min,max,err,2,Array_i,Array_r)
722  IF(err.NE.CG_OK)THEN
723     PRINT*,'*FAILED* cgp_array_multi_read_data_f'
724     CALL cgp_error_exit_f()
725  ENDIF
726#else
727  CALL cgp_array_read_data_f(Ar, min, max, Array_r, err)
728  IF(err.NE.CG_OK)THEN
729     PRINT*,'*FAILED* cgp_field_read_data (Array_r)'
730     CALL cgp_error_exit_f()
731  ENDIF
732  CALL cgp_array_read_data_f(Ai, min, max, Array_i, err)
733  IF(err.NE.CG_OK)THEN
734     PRINT*,'*FAILED* cgp_field_read_data (Array_i)'
735     CALL cgp_error_exit_f()
736  ENDIF
737#endif
738  t2 = MPI_Wtime()
739  xtiming(9) = t2-t1
740
741  ! Check if read the data back correctly
742  IF(debug)THEN
743     DO k = 1, count
744        IF(.NOT.check_eq(Array_r(k), REAL(comm_rank*count + k, KIND=dp) + 0.001_dp).OR. &
745             .NOT.check_eq(Array_i(k), comm_rank*count*NodePerElem + k)) THEN
746           PRINT*,'*FAILED* cgp_array_read_data values are incorrect'
747           CALL cgp_error_exit_f()
748        ENDIF
749     ENDDO
750  ENDIF
751
752  DEALLOCATE(Array_r, Array_i)
753
754  CALL cgp_close_f(fn, err)
755  IF(err.NE.CG_OK)THEN
756     PRINT*,'*FAILED* cgp_close'
757     CALL cgp_error_exit_f()
758  ENDIF
759  t2 = MPI_Wtime()
760
761  xtiming(1) = t2-t0
762
763  CALL MPI_Reduce(xtiming, timing, 15, MPI_DOUBLE, &
764       MPI_SUM, 0, MPI_COMM_WORLD, mpi_err)
765  CALL MPI_Reduce(xtiming, timingMin, 15, MPI_DOUBLE, &
766       MPI_MIN, 0, MPI_COMM_WORLD, mpi_err)
767  CALL MPI_Reduce(xtiming, timingMax, 15, MPI_DOUBLE, &
768       MPI_MAX, 0, MPI_COMM_WORLD, mpi_err)
769
770  IF(comm_rank.EQ.0)THEN
771     OPEN(10,FILE="timing_"//ichr6//".dat", FORM='formatted')
772     WRITE(10,'(A)')"#nprocs, wcoord, welem, wfield, warray, rcoord, relem, rfield, rarray"
773     WRITE(10,'(i0,100(x,3(f14.7)))') comm_size, &
774          timing(1)/DBLE(comm_size), timingMin(1), timingMax(1), &
775          timing(2)/DBLE(comm_size), timingMin(2), timingMax(2), &
776          timing(3)/DBLE(comm_size), timingMin(3), timingMax(3), &
777          timing(4)/DBLE(comm_size), timingMin(4), timingMax(4), &
778          timing(5)/DBLE(comm_size), timingMin(5), timingMax(5), &
779          timing(6)/DBLE(comm_size), timingMin(6), timingMax(6), &
780          timing(7)/DBLE(comm_size), timingMin(7), timingMax(7), &
781          timing(8)/DBLE(comm_size), timingMin(8), timingMax(8), &
782          timing(9)/DBLE(comm_size), timingMin(9), timingMax(9), &
783          timing(10)/DBLE(comm_size), timingMin(10), timingMax(10), &
784          timing(11)/DBLE(comm_size), timingMin(11), timingMax(11), &
785          timing(12)/DBLE(comm_size), timingMin(12), timingMax(12), &
786          timing(13)/DBLE(comm_size), timingMin(13), timingMax(13), &
787          timing(14)/DBLE(comm_size), timingMin(14), timingMax(14), &
788          timing(15)/DBLE(comm_size), timingMin(15), timingMax(15)
789  ENDIF
790
791  CALL MPI_FINALIZE(mpi_err)
792
793END PROGRAM benchmark_hdf5_f90
794
795