1program TestBPWriteReadHeatMap3D
2  use mpi
3  use adios2
4
5  implicit none
6
7  integer(kind=8) :: sum_i1, sum_i2
8  type(adios2_adios) :: adios
9  type(adios2_io) :: ioPut, ioGet
10  type(adios2_engine) :: bpWriter, bpReader
11  type(adios2_variable), dimension(6) :: var_temperatures, var_temperaturesIn
12
13  integer(kind=1), dimension(:, :, :), allocatable :: temperatures_i1, &
14                                                      sel_temperatures_i1
15
16  integer(kind=2), dimension(:, :, :), allocatable :: temperatures_i2, &
17                                                      sel_temperatures_i2
18
19  integer(kind=4), dimension(:, :, :), allocatable :: temperatures_i4, &
20                                                      sel_temperatures_i4
21
22  integer(kind=8), dimension(:, :, :), allocatable :: temperatures_i8, &
23                                                      sel_temperatures_i8
24
25  real(kind=4), dimension(:, :, :), allocatable :: temperatures_r4, &
26                                                   sel_temperatures_r4
27
28  real(kind=8), dimension(:, :, :), allocatable :: temperatures_r8, &
29                                                   sel_temperatures_r8
30
31  integer(kind=8), dimension(3) :: ishape, istart, icount
32  integer(kind=8), dimension(3) :: sel_start, sel_count
33  integer :: ierr, irank, isize
34  integer :: in1, in2, in3
35  integer :: i1, i2, i3
36
37  call MPI_INIT(ierr)
38  call MPI_COMM_RANK(MPI_COMM_WORLD, irank, ierr)
39  call MPI_COMM_SIZE(MPI_COMM_WORLD, isize, ierr)
40
41  in1 = 5
42  in2 = 10
43  in3 = 10
44
45  icount = (/in1, in2, in3/)
46  istart = (/0, 0, in3*irank/)
47  ishape = (/in1, in2, in3*isize/)
48
49  allocate (temperatures_i1(in1, in2, in3))
50  allocate (temperatures_i2(in1, in2, in3))
51  allocate (temperatures_i4(in1, in2, in3))
52  allocate (temperatures_i8(in1, in2, in3))
53  allocate (temperatures_r4(in1, in2, in3))
54  allocate (temperatures_r8(in1, in2, in3))
55
56  temperatures_i1 = 1
57  temperatures_i2 = 1
58  temperatures_i4 = 1
59  temperatures_i8 = 1_8
60  temperatures_r4 = 1.0
61  temperatures_r8 = 1.0_8
62
63  ! Start adios2 Writer
64  call adios2_init(adios, MPI_COMM_WORLD, adios2_debug_mode_on, ierr)
65  call adios2_declare_io(ioPut, adios, 'HeatMapWrite', ierr)
66
67  call adios2_define_variable(var_temperatures(1), ioPut, &
68                              'temperatures_i1', adios2_type_integer1, &
69                              3, ishape, istart, icount, &
70                              adios2_constant_dims, ierr)
71
72  call adios2_define_variable(var_temperatures(2), ioPut, &
73                              'temperatures_i2', adios2_type_integer2, &
74                              3, ishape, istart, icount, &
75                              adios2_constant_dims, ierr)
76
77  call adios2_define_variable(var_temperatures(3), ioPut, &
78                              'temperatures_i4', adios2_type_integer4, &
79                              3, ishape, istart, icount, &
80                              adios2_constant_dims, ierr)
81
82  call adios2_define_variable(var_temperatures(4), ioPut, &
83                              'temperatures_i8', adios2_type_integer8, &
84                              3, ishape, istart, icount, &
85                              adios2_constant_dims, ierr)
86
87  call adios2_define_variable(var_temperatures(5), ioPut, &
88                              'temperatures_r4', adios2_type_real, &
89                              3, ishape, istart, icount, &
90                              adios2_constant_dims, ierr)
91
92  call adios2_define_variable(var_temperatures(6), ioPut, &
93                              'temperatures_r8', adios2_type_dp, &
94                              3, ishape, istart, icount, &
95                              adios2_constant_dims, ierr)
96
97  call adios2_open(bpWriter, ioPut, 'HeatMap3D_f.bp', adios2_mode_write, &
98                   ierr)
99
100  call adios2_put(bpWriter, var_temperatures(1), temperatures_i1, ierr)
101  call adios2_put(bpWriter, var_temperatures(2), temperatures_i2, ierr)
102  call adios2_put(bpWriter, var_temperatures(3), temperatures_i4, ierr)
103  call adios2_put(bpWriter, var_temperatures(4), temperatures_i8, ierr)
104  call adios2_put(bpWriter, var_temperatures(5), temperatures_r4, ierr)
105  call adios2_put(bpWriter, var_temperatures(6), temperatures_r8, ierr)
106
107  call adios2_close(bpWriter, ierr)
108
109  if (allocated(temperatures_i1)) deallocate (temperatures_i1)
110  if (allocated(temperatures_i2)) deallocate (temperatures_i2)
111  if (allocated(temperatures_i4)) deallocate (temperatures_i4)
112  if (allocated(temperatures_i8)) deallocate (temperatures_i8)
113  if (allocated(temperatures_r4)) deallocate (temperatures_r4)
114  if (allocated(temperatures_r8)) deallocate (temperatures_r8)
115
116  ! Start adios2 Reader in rank 0
117  if (irank == 0) then
118
119    call adios2_declare_io(ioGet, adios, 'HeatMapRead', ierr)
120
121    call adios2_open(bpReader, ioGet, 'HeatMap3D_f.bp', &
122                     adios2_mode_read, MPI_COMM_SELF, ierr)
123
124    call adios2_inquire_variable(var_temperaturesIn(1), ioGet, &
125                                 'temperatures_i1', ierr)
126    call adios2_inquire_variable(var_temperaturesIn(2), ioGet, &
127                                 'temperatures_i2', ierr)
128    call adios2_inquire_variable(var_temperaturesIn(3), ioGet, &
129                                 'temperatures_i4', ierr)
130    call adios2_inquire_variable(var_temperaturesIn(4), ioGet, &
131                                 'temperatures_i8', ierr)
132    call adios2_inquire_variable(var_temperaturesIn(5), ioGet, &
133                                 'temperatures_r4', ierr)
134    call adios2_inquire_variable(var_temperaturesIn(6), ioGet, &
135                                 'temperatures_r8', ierr)
136
137    sel_start = (/0, 0, 0/)
138    sel_count = (/ishape(1), ishape(2), ishape(3)/)
139
140    allocate (sel_temperatures_i1(ishape(1), ishape(2), ishape(3)))
141    allocate (sel_temperatures_i2(ishape(1), ishape(2), ishape(3)))
142    allocate (sel_temperatures_i4(ishape(1), ishape(2), ishape(3)))
143    allocate (sel_temperatures_i8(ishape(1), ishape(2), ishape(3)))
144    allocate (sel_temperatures_r4(ishape(1), ishape(2), ishape(3)))
145    allocate (sel_temperatures_r8(ishape(1), ishape(2), ishape(3)))
146
147    sel_temperatures_i1 = 0
148    sel_temperatures_i2 = 0
149    sel_temperatures_i4 = 0
150    sel_temperatures_i8 = 0_8
151    sel_temperatures_r4 = 0.0_4
152    sel_temperatures_r8 = 0.0_8
153
154    call adios2_set_selection(var_temperaturesIn(1), 3, sel_start, sel_count, &
155                              ierr)
156    call adios2_set_selection(var_temperaturesIn(2), 3, sel_start, sel_count, &
157                              ierr)
158    call adios2_set_selection(var_temperaturesIn(3), 3, sel_start, sel_count, &
159                              ierr)
160    call adios2_set_selection(var_temperaturesIn(4), 3, sel_start, sel_count, &
161                              ierr)
162    call adios2_set_selection(var_temperaturesIn(5), 3, sel_start, sel_count, &
163                              ierr)
164    call adios2_set_selection(var_temperaturesIn(6), 3, sel_start, sel_count, &
165                              ierr)
166
167    call adios2_get(bpReader, var_temperaturesIn(1), sel_temperatures_i1, ierr)
168    call adios2_get(bpReader, var_temperaturesIn(2), sel_temperatures_i2, ierr)
169    call adios2_get(bpReader, var_temperaturesIn(3), sel_temperatures_i4, ierr)
170    call adios2_get(bpReader, var_temperaturesIn(4), sel_temperatures_i8, ierr)
171    call adios2_get(bpReader, var_temperaturesIn(5), sel_temperatures_r4, ierr)
172    call adios2_get(bpReader, var_temperaturesIn(6), sel_temperatures_r8, ierr)
173
174    call adios2_close(bpReader, ierr)
175
176    sum_i1 = 0
177    sum_i2 = 0
178
179    do i3 = 1, sel_count(3)
180      do i2 = 1, sel_count(2)
181        do i1 = 1, sel_count(1)
182          sum_i1 = sum_i1 + sel_temperatures_i1(i1, i2, i3)
183          sum_i2 = sum_i2 + sel_temperatures_i2(i1, i2, i3)
184        end do
185      end do
186    end do
187
188    if (sum_i1 /= 500*isize) stop 'Test failed integer*1'
189    if (sum_i2 /= 500*isize) stop 'Test failed integer*2'
190    if (sum(sel_temperatures_i4) /= 500*isize) stop 'Test failed integer*4'
191    if (sum(sel_temperatures_i8) /= 500*isize) stop 'Test failed integer*8'
192    if (sum(sel_temperatures_r4) /= 500*isize) stop 'Test failed real*4'
193    if (sum(sel_temperatures_r8) /= 500*isize) stop 'Test failed real*8'
194
195    if (allocated(sel_temperatures_i1)) deallocate (sel_temperatures_i1)
196    if (allocated(sel_temperatures_i2)) deallocate (sel_temperatures_i2)
197    if (allocated(sel_temperatures_i4)) deallocate (sel_temperatures_i4)
198    if (allocated(sel_temperatures_i8)) deallocate (sel_temperatures_i8)
199    if (allocated(sel_temperatures_r4)) deallocate (sel_temperatures_r4)
200    if (allocated(sel_temperatures_r8)) deallocate (sel_temperatures_r8)
201
202  end if
203
204  call adios2_finalize(adios, ierr)
205  call MPI_Finalize(ierr)
206
207end program TestBPWriteReadHeatMap3D
208