25 subroutine scatter_data_bc_2d(fileobj, varname, vdata, bc_info, unlim_dim_level, ignore_checksum)
26 class(fmsnetcdffile_t),
intent(inout) :: fileobj
27 character(len=*),
intent(in) :: varname
28 class(*),
dimension(:,:),
intent(inout) :: vdata
29 type(bc_information),
intent(inout) :: bc_info
30 integer,
intent(in),
optional :: unlim_dim_level
31 logical,
intent(in),
optional :: ignore_checksum
33 real(kind=r8_kind),
dimension(:,:),
allocatable :: global_buf_r8_kind
34 real(kind=r8_kind),
dimension(:,:),
allocatable :: local_buf_r8_kind
35 real(kind=r4_kind),
dimension(:,:),
allocatable :: global_buf_r4_kind
36 real(kind=r4_kind),
dimension(:,:),
allocatable :: local_buf_r4_kind
38 integer(kind=i8_kind) :: chksum_val
39 character(len=32) :: chksum
40 character(len=32) :: chksum_read
41 logical :: chksum_ignore = .false.
44 integer :: isc, iec, jsc, jec
45 integer :: i1, i2, j1, j2
46 integer :: i_add, j_add
51 if (
PRESENT(ignore_checksum)) chksum_ignore = ignore_checksum
54 isc = bc_info%indices(1)
55 iec = bc_info%indices(2)
56 jsc = bc_info%indices(3)
57 jec = bc_info%indices(4)
60 i1 = 1 + bc_info%x_halo
62 j1 = 1 + bc_info%y_halo
66 i_add = bc_info%ishift
67 j_add = bc_info%jshift
69 if (fileobj%is_root)
then
70 i_glob = bc_info%global_size(1)
71 j_glob = bc_info%global_size(2)
75 type is (real(kind=r4_kind))
76 if (fileobj%is_root)
then
78 allocate(global_buf_r4_kind(i_glob, j_glob))
80 call netcdf_read_data(fileobj, varname, &
82 unlim_dim_level=unlim_dim_level, &
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()/))
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:"// &
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
105 allocate(local_buf_r4_kind(i1:i2,j1:j2))
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)
113 if (fileobj%is_root .and. .not. bc_info%data_on_file_root)
then
114 deallocate(local_buf_r4_kind)
119 vdata(i1:i2,j1:j2) = local_buf_r4_kind(i1:i2,j1:j2)
121 deallocate(local_buf_r4_kind)
122 if (fileobj%is_root)
deallocate(global_buf_r4_kind)
124 type is (real(kind=r8_kind))
125 if (fileobj%is_root)
then
127 allocate(global_buf_r8_kind(i_glob, j_glob))
129 call netcdf_read_data(fileobj, varname, &
130 global_buf_r8_kind, &
131 unlim_dim_level=unlim_dim_level, &
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()/))
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)
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
153 allocate(local_buf_r8_kind(i1:i2,j1:j2))
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)
161 if (fileobj%is_root .and. .not. bc_info%data_on_file_root)
then
162 deallocate(local_buf_r8_kind)
167 vdata(i1:i2,j1:j2) = local_buf_r8_kind(i1:i2,j1:j2)
169 deallocate(local_buf_r8_kind)
170 if (fileobj%is_root)
deallocate(global_buf_r8_kind)
173 call error(
"scatter_data_bc_2d: unsupported type. Currently only r8_kind and r4_kinds are supported")
179 class(fmsnetcdffile_t),
intent(inout) :: fileobj
180 character(len=*),
intent(in) :: varname
181 class(*),
dimension(:,:,:),
intent(inout) :: vdata
182 type(bc_information),
intent(inout) :: bc_info
183 integer,
intent(in),
optional :: unlim_dim_level
184 logical,
intent(in),
optional :: ignore_checksum
186 real(kind=r8_kind),
dimension(:,:,:),
allocatable :: global_buf_r8_kind
187 real(kind=r8_kind),
dimension(:,:,:),
allocatable :: local_buf_r8_kind
188 real(kind=r4_kind),
dimension(:,:,:),
allocatable :: global_buf_r4_kind
189 real(kind=r4_kind),
dimension(:,:,:),
allocatable :: local_buf_r4_kind
191 integer(kind=i8_kind) :: chksum_val
192 character(len=32) :: chksum
193 character(len=32) :: chksum_read
194 logical :: chksum_ignore = .false.
197 integer :: isc, iec, jsc, jec
198 integer :: i1, i2, j1, j2
199 integer :: i_add, j_add
205 if (
PRESENT(ignore_checksum)) chksum_ignore = ignore_checksum
208 isc = bc_info%indices(1)
209 iec = bc_info%indices(2)
210 jsc = bc_info%indices(3)
211 jec = bc_info%indices(4)
214 i1 = 1 + bc_info%x_halo
216 j1 = 1 + bc_info%y_halo
220 i_add = bc_info%ishift
221 j_add = bc_info%jshift
223 if (fileobj%is_root)
then
224 i_glob = bc_info%global_size(1)
225 j_glob = bc_info%global_size(2)
228 k_glob = bc_info%global_size(3)
231 type is (real(kind=r4_kind))
232 if (fileobj%is_root)
then
234 allocate(global_buf_r4_kind(i_glob, j_glob, k_glob))
236 call netcdf_read_data(fileobj, varname, &
237 global_buf_r4_kind, &
238 unlim_dim_level=unlim_dim_level, &
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()/))
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:"// &
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
261 allocate(local_buf_r4_kind(i1:i2,j1:j2,k_glob))
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)
269 if (fileobj%is_root .and. .not. bc_info%data_on_file_root)
then
270 deallocate(local_buf_r4_kind)
275 vdata(i1:i2,j1:j2,:) = local_buf_r4_kind(i1:i2,j1:j2,:)
277 deallocate(local_buf_r4_kind)
278 if (fileobj%is_root)
deallocate(global_buf_r4_kind)
280 type is (real(kind=r8_kind))
281 if (fileobj%is_root)
then
283 allocate(global_buf_r8_kind(i_glob, j_glob, k_glob))
285 call netcdf_read_data(fileobj, varname, &
286 global_buf_r8_kind, &
287 unlim_dim_level=unlim_dim_level, &
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()/))
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)
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;
309 allocate(local_buf_r8_kind(i1:i2,j1:j2,k_glob))
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)
317 if (fileobj%is_root .and. .not. bc_info%data_on_file_root)
then
318 deallocate(local_buf_r8_kind)
323 vdata(i1:i2,j1:j2,:) = local_buf_r8_kind(i1:i2,j1:j2,:)
325 deallocate(local_buf_r8_kind)
326 if (fileobj%is_root)
deallocate(global_buf_r8_kind)
329 call error(
"scatter_data_bc_3d: unsupported type. Currently only r8_kind and r4_kinds are supported")
integer function mpp_pe()
Returns processor ID.
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)