FMS  2025.04
Flexible Modeling System
scatter_data_bc.inc
1 !***********************************************************************
2 !* Apache License 2.0
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* Licensed under the Apache License, Version 2.0 (the "License");
7 !* you may not use this file except in compliance with the License.
8 !* You may obtain a copy of the License at
9 !*
10 !* http://www.apache.org/licenses/LICENSE-2.0
11 !*
12 !* FMS is distributed in the hope that it will be useful, but WITHOUT
13 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
14 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
15 !* PARTICULAR PURPOSE. See the License for the specific language
16 !* governing permissions and limitations under the License.
17 !***********************************************************************
18 !> @file
19 !> @brief Routines for the @ref gather_data_bc interface
20 
21 !> @addtogroup netcdf_io_mod
22 !> @{
23 
24 !< Root pe reads the data in the netcdf file and scatters it to the other pes
25 subroutine scatter_data_bc_2d(fileobj, varname, vdata, bc_info, unlim_dim_level, ignore_checksum)
26  class(fmsnetcdffile_t), intent(inout) :: fileobj !< Fms2io netcdf fileobj
27  character(len=*), intent(in) :: varname !< Variable name.
28  class(*), dimension(:,:), intent(inout) :: vdata !< scattered data
29  type(bc_information), intent(inout) :: bc_info !< information about the boundary condition variable
30  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension level.
31  logical, intent(in), optional :: ignore_checksum !< Checksum data integrity flag.
32 
33  real(kind=r8_kind), dimension(:,:), allocatable :: global_buf_r8_kind !< buffer with a data read in from file
34  real(kind=r8_kind), dimension(:,:), allocatable :: local_buf_r8_kind !< current pe's portion of the data
35  real(kind=r4_kind), dimension(:,:), allocatable :: global_buf_r4_kind !< buffer with a data read in from file
36  real(kind=r4_kind), dimension(:,:), allocatable :: local_buf_r4_kind !< current pe's portion of the data
37 
38  integer(kind=i8_kind) :: chksum_val !< Checksum value calculated from the read data
39  character(len=32) :: chksum !< Cheksum value converted to a string
40  character(len=32) :: chksum_read !< String checksum read in from file
41  logical :: chksum_ignore = .false. !< local variable for data integrity checks
42  !! default: .FALSE. - checks enabled
43 
44  integer :: isc, iec, jsc, jec !< current PE's indices
45  integer :: i1, i2, j1, j2 !< current PE's indices relative to the global domain
46  integer :: i_add, j_add !< current PE's shift
47 
48  integer :: i_glob !< Size of the global domain in x
49  integer :: j_glob !< Size of the global domain in y
50 
51  if (PRESENT(ignore_checksum)) chksum_ignore = ignore_checksum
52 
53  !> Set the indices
54  isc = bc_info%indices(1)
55  iec = bc_info%indices(2)
56  jsc = bc_info%indices(3)
57  jec = bc_info%indices(4)
58 
59  !> This is the section of the PE that will actually be read from the global_buffer
60  i1 = 1 + bc_info%x_halo
61  i2 = i1 + (iec-isc)
62  j1 = 1 + bc_info%y_halo
63  j2 = j1 + (jec-jsc)
64 
65  !> Set up index shifts for global array
66  i_add = bc_info%ishift
67  j_add = bc_info%jshift
68 
69  if (fileobj%is_root) then
70  i_glob = bc_info%global_size(1)
71  j_glob = bc_info%global_size(2)
72  endif
73 
74  select type(vdata)
75  type is (real(kind=r4_kind))
76  if (fileobj%is_root) then
77  !> If you are the file root, read in the data
78  allocate(global_buf_r4_kind(i_glob, j_glob))
79 
80  call netcdf_read_data(fileobj, varname, &
81  global_buf_r4_kind, &
82  unlim_dim_level=unlim_dim_level, &
83  broadcast=.false.)
84  !> If the checksum exists read it and compare it with the calculated from the data that was read
85  if (.not.chksum_ignore .and. variable_att_exists(fileobj, varname, "checksum", broadcast=.false.)) then
86  call get_variable_attribute(fileobj, varname, "checksum", chksum_read, broadcast=.false.)
87  chksum_val = mpp_chksum(global_buf_r4_kind, (/mpp_pe()/))
88  chksum = ""
89  write(chksum, "(Z16)") chksum_val
90  if (.not. string_compare(adjustl(chksum), adjustl(chksum_read))) then
91  call error("scatter_data_bc_2d: "//trim(varname)//" chksum_in:"//(chksum)//" chksum_out:"// &
92  & (chksum_read))
93  endif
94  endif
95  endif
96 
97  !> If the fileobj's root is not the same as the variable's root, then no data will be read for the
98  !! file root
99  if (fileobj%is_root .and. .not. bc_info%data_on_file_root) then
100  allocate(local_buf_r4_kind(1,1))
101  local_buf_r4_kind = 0.
102  isc = 1; i_add=0; iec=1; jsc=1; j_add=0; jec=1
103  i1=1; i2=1; j1=1; j2=1
104  else
105  allocate(local_buf_r4_kind(i1:i2,j1:j2))
106  endif
107 
108  !> Scatter the data to the other PE
109  call mpp_scatter(isc+i_add, iec+i_add, jsc+j_add, jec+j_add, bc_info%pelist, &
110  local_buf_r4_kind(i1:i2,j1:j2), global_buf_r4_kind, fileobj%is_root)
111 
112  !> Return if fileobj's root is not the same as the variable's root
113  if (fileobj%is_root .and. .not. bc_info%data_on_file_root) then
114  deallocate(local_buf_r4_kind)
115  return
116  endif
117 
118  !> Set vdata to equal the scattered data
119  vdata(i1:i2,j1:j2) = local_buf_r4_kind(i1:i2,j1:j2)
120 
121  deallocate(local_buf_r4_kind)
122  if (fileobj%is_root) deallocate(global_buf_r4_kind)
123 
124  type is (real(kind=r8_kind))
125  if (fileobj%is_root) then
126  !> If you are the file root, read in the data
127  allocate(global_buf_r8_kind(i_glob, j_glob))
128 
129  call netcdf_read_data(fileobj, varname, &
130  global_buf_r8_kind, &
131  unlim_dim_level=unlim_dim_level, &
132  broadcast=.false.)
133  !> If the checksum exists read it and compare it with the calculated from the data that was read
134  if (.not.chksum_ignore .and. variable_att_exists(fileobj, varname, "checksum", broadcast=.false.)) then
135  call get_variable_attribute(fileobj, varname, "checksum", chksum_read, broadcast=.false.)
136  chksum_val = mpp_chksum(global_buf_r8_kind, (/mpp_pe()/))
137  chksum = ""
138  write(chksum, "(Z16)") chksum_val
139  if (.not. string_compare(adjustl(chksum), adjustl(chksum_read))) then
140  call error("scatter_data_bc_2d: "//trim(varname)//" chksum_in:"//chksum//" chksum_out:"//chksum_read)
141  endif
142  endif
143  endif
144 
145  !> If the fileobj's root is not the same as the variable's root, then no data will be read for the
146  !! file root
147  if (fileobj%is_root .and. .not. bc_info%data_on_file_root) then
148  allocate(local_buf_r8_kind(1,1))
149  local_buf_r8_kind = 0.
150  isc = 1; i_add=0; iec=1; jsc=1; j_add=0; jec=1
151  i1=1; i2=1; j1=1; j2=1
152  else
153  allocate(local_buf_r8_kind(i1:i2,j1:j2))
154  endif
155 
156  !> Scatter the data to the other PE
157  call mpp_scatter(isc+i_add, iec+i_add, jsc+j_add, jec+j_add, bc_info%pelist, &
158  local_buf_r8_kind(i1:i2,j1:j2), global_buf_r8_kind, fileobj%is_root)
159 
160  !> Return if fileobj's root is not the same as the variable's root
161  if (fileobj%is_root .and. .not. bc_info%data_on_file_root) then
162  deallocate(local_buf_r8_kind)
163  return
164  endif
165 
166  !> Set vdata to equal the scattered data
167  vdata(i1:i2,j1:j2) = local_buf_r8_kind(i1:i2,j1:j2)
168 
169  deallocate(local_buf_r8_kind)
170  if (fileobj%is_root) deallocate(global_buf_r8_kind)
171 
172  class default
173  call error("scatter_data_bc_2d: unsupported type. Currently only r8_kind and r4_kinds are supported")
174  end select
175 end subroutine scatter_data_bc_2d
176 
177 !< Root pe reads the data in the netcdf file and scatters it to the other pes
178 subroutine scatter_data_bc_3d(fileobj, varname, vdata, bc_info, unlim_dim_level, ignore_checksum)
179  class(fmsnetcdffile_t), intent(inout) :: fileobj !< Fms2io netcdf fileobj
180  character(len=*), intent(in) :: varname !< Variable name.
181  class(*), dimension(:,:,:), intent(inout) :: vdata !< scattered data
182  type(bc_information), intent(inout) :: bc_info !< information about the boundary condition variable
183  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension level.
184  logical, intent(in), optional :: ignore_checksum !< Checksum data integrity flag.
185 
186  real(kind=r8_kind), dimension(:,:,:), allocatable :: global_buf_r8_kind !< buffer with a data read in from file
187  real(kind=r8_kind), dimension(:,:,:), allocatable :: local_buf_r8_kind !< current pe's portion of the data
188  real(kind=r4_kind), dimension(:,:,:), allocatable :: global_buf_r4_kind !< buffer with a data read in from file
189  real(kind=r4_kind), dimension(:,:,:), allocatable :: local_buf_r4_kind !< current pe's portion of the data
190 
191  integer(kind=i8_kind) :: chksum_val !< Checksum value calculated from the read data
192  character(len=32) :: chksum !< Cheksum value converted to a string
193  character(len=32) :: chksum_read !< String checksum read in from file
194  logical :: chksum_ignore = .false. !< local variable for data integrity checks
195  !! default: .FALSE. - checks enabled
196 
197  integer :: isc, iec, jsc, jec !< current PE's indices
198  integer :: i1, i2, j1, j2 !< current PE's indices relative to the global domain
199  integer :: i_add, j_add !< current PE's shift
200 
201  integer :: i_glob !< Size of the global domain in x
202  integer :: j_glob !< Size of the global domain in y
203  integer :: k_glob !< Size of the z axis
204 
205  if (PRESENT(ignore_checksum)) chksum_ignore = ignore_checksum
206 
207  !> Set the indices
208  isc = bc_info%indices(1)
209  iec = bc_info%indices(2)
210  jsc = bc_info%indices(3)
211  jec = bc_info%indices(4)
212 
213  !> This is the section of the PE that will actually be read from the global_buffer
214  i1 = 1 + bc_info%x_halo
215  i2 = i1 + (iec-isc)
216  j1 = 1 + bc_info%y_halo
217  j2 = j1 + (jec-jsc)
218 
219  !> Set up index shifts for global array
220  i_add = bc_info%ishift
221  j_add = bc_info%jshift
222 
223  if (fileobj%is_root) then
224  i_glob = bc_info%global_size(1)
225  j_glob = bc_info%global_size(2)
226  endif
227 
228  k_glob = bc_info%global_size(3)
229 
230  select type(vdata)
231  type is (real(kind=r4_kind))
232  if (fileobj%is_root) then
233  !> If you are the file root, read in the data
234  allocate(global_buf_r4_kind(i_glob, j_glob, k_glob))
235 
236  call netcdf_read_data(fileobj, varname, &
237  global_buf_r4_kind, &
238  unlim_dim_level=unlim_dim_level, &
239  broadcast=.false.)
240  !> If the checksum exists read it and compare it with the calculated from the data that was read
241  if (.not.chksum_ignore .and. variable_att_exists(fileobj, varname, "checksum", broadcast=.false.)) then
242  call get_variable_attribute(fileobj, varname, "checksum", chksum_read, broadcast=.false.)
243  chksum_val = mpp_chksum(global_buf_r4_kind, (/mpp_pe()/))
244  chksum = ""
245  write(chksum, "(Z16)") chksum_val
246  if (.not. string_compare(adjustl(chksum), adjustl(chksum_read))) then
247  call error("scatter_data_bc_2d: "//trim(varname)//" chksum_in:"//(chksum)//" chksum_out:"// &
248  & (chksum_read))
249  endif
250  endif
251  endif
252 
253  !> If the fileobj's root is not the same as the variable's root, then no data will be read for the
254  !! file root
255  if (fileobj%is_root .and. .not. bc_info%data_on_file_root) then
256  allocate(local_buf_r4_kind(1,1,k_glob))
257  local_buf_r4_kind = 0.
258  isc = 1; i_add=0; iec=1; jsc=1; j_add=0; jec=1
259  i1=1; i2=1; j1=1; j2=1
260  else
261  allocate(local_buf_r4_kind(i1:i2,j1:j2,k_glob))
262  endif
263 
264  !> Scatter the data to the other PE
265  call mpp_scatter(isc+i_add, iec+i_add, jsc+j_add, jec+j_add, k_glob, bc_info%pelist, &
266  local_buf_r4_kind(i1:i2,j1:j2,:), global_buf_r4_kind, fileobj%is_root)
267 
268  !> Return if fileobj's root is not the same as the variable's root
269  if (fileobj%is_root .and. .not. bc_info%data_on_file_root) then
270  deallocate(local_buf_r4_kind)
271  return
272  endif
273 
274  !> Set vdata to equal the scattered data
275  vdata(i1:i2,j1:j2,:) = local_buf_r4_kind(i1:i2,j1:j2,:)
276 
277  deallocate(local_buf_r4_kind)
278  if (fileobj%is_root) deallocate(global_buf_r4_kind)
279 
280  type is (real(kind=r8_kind))
281  if (fileobj%is_root) then
282  !> If you are the file root, read in the data
283  allocate(global_buf_r8_kind(i_glob, j_glob, k_glob))
284 
285  call netcdf_read_data(fileobj, varname, &
286  global_buf_r8_kind, &
287  unlim_dim_level=unlim_dim_level, &
288  broadcast=.false.)
289  !> If the checksum exists read it and compare it with the calculated from the data that was read
290  if (.not.chksum_ignore .and. variable_att_exists(fileobj, varname, "checksum", broadcast=.false.)) then
291  call get_variable_attribute(fileobj, varname, "checksum", chksum_read, broadcast=.false.)
292  chksum_val = mpp_chksum(global_buf_r8_kind, (/mpp_pe()/))
293  chksum = ""
294  write(chksum, "(Z16)") chksum_val
295  if (.not. string_compare(adjustl(chksum), adjustl(chksum_read))) then
296  call error("scatter_data_bc_2d: "//trim(varname)//" chksum_in:"//chksum//" chksum_out:"//chksum_read)
297  endif
298  endif
299  endif
300 
301  !> If the fileobj's root is not the same as the variable's root, then no data will be read for the
302  !! file root
303  if (fileobj%is_root .and. .not. bc_info%data_on_file_root) then
304  allocate(local_buf_r8_kind(1,1,k_glob))
305  local_buf_r8_kind = 0.
306  isc = 1; i_add=0; iec=1; jsc=1; j_add=0; jec=1
307  i1=1; i2=1; j1=1; j2=1;
308  else
309  allocate(local_buf_r8_kind(i1:i2,j1:j2,k_glob))
310  endif
311 
312  !> Scatter the data to the other PE
313  call mpp_scatter(isc+i_add, iec+i_add, jsc+j_add, jec+j_add, k_glob, bc_info%pelist, &
314  local_buf_r8_kind(i1:i2,j1:j2,:), global_buf_r8_kind, fileobj%is_root)
315 
316  !> Return if fileobj's root is not the same as the variable's root
317  if (fileobj%is_root .and. .not. bc_info%data_on_file_root) then
318  deallocate(local_buf_r8_kind)
319  return
320  endif
321 
322  !> Set vdata to equal the scattered data
323  vdata(i1:i2,j1:j2,:) = local_buf_r8_kind(i1:i2,j1:j2,:)
324 
325  deallocate(local_buf_r8_kind)
326  if (fileobj%is_root) deallocate(global_buf_r8_kind)
327 
328  class default
329  call error("scatter_data_bc_3d: unsupported type. Currently only r8_kind and r4_kinds are supported")
330  end select
331 end subroutine scatter_data_bc_3d
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:406
subroutine scatter_data_bc_2d(fileobj, varname, vdata, bc_info, unlim_dim_level, ignore_checksum)
subroutine scatter_data_bc_3d(fileobj, varname, vdata, bc_info, unlim_dim_level, ignore_checksum)