FMS 2025.01.02-dev
Flexible Modeling System
All Classes Namespaces Files Functions Variables Modules Pages
compressed_write.inc
1!***********************************************************************
2!* GNU Lesser General Public License
3!*
4!* This file is part of the GFDL Flexible Modeling System (FMS).
5!*
6!* FMS is free software: you can redistribute it and/or modify it under
7!* the terms of the GNU Lesser General Public License as published by
8!* the Free Software Foundation, either version 3 of the License, or (at
9!* your option) any later version.
10!*
11!* FMS is distributed in the hope that it will be useful, but WITHOUT
12!* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13!* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14!* for more details.
15!*
16!* You should have received a copy of the GNU Lesser General Public
17!* License along with FMS. If not, see <http://www.gnu.org/licenses/>.
18!***********************************************************************
19!> @file
20!> @brief Compressed write routines for @ref write_data interface
21
22!> @addtogroup fms2_io_mod
23!> @{
24
25!> @brief For variables without a compressed dimension, this routine simply wraps
26!! netcdf_write data. If the variable does have a compressed axis, the I/O root
27!! gathers the data from the rest of the ranks and then writes the combined data to
28!! the netcdf file.
29subroutine compressed_write_0d(fileobj, variable_name, cdata, unlim_dim_level, &
30 corner)
31
32 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
33 character(len=*), intent(in) :: variable_name !< Variable name.
34 class(*), intent(in) :: cdata !< Compressed data that
35 !! will be gathered and
36 !! written to the
37 !! netcdf file.
38 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
39 !! dimension.
40 integer, intent(in), optional :: corner !< Array of starting
41 !! indices describing
42 !! where the data
43 !! will be written to.
44
45 integer, dimension(2) :: compressed_dim_index
46
47 compressed_dim_index = get_variable_compressed_dimension_index(fileobj, variable_name)
48 if (compressed_dim_index(1) .eq. dimension_not_found) then
49 call netcdf_write_data(fileobj, variable_name, cdata, &
50 unlim_dim_level=unlim_dim_level, corner=corner)
51 return
52 endif
53end subroutine compressed_write_0d
54
55
56!> @brief Wrapper to distinguish interfaces.
57subroutine compressed_write_0d_wrap(fileobj, variable_name, cdata, unlim_dim_level, &
58 corner)
59
60 type(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
61 character(len=*), intent(in) :: variable_name !< Variable name.
62 class(*), intent(in) :: cdata !< Compressed data that
63 !! will be gathered and
64 !! written to the
65 !! netcdf file.
66 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
67 !! dimension.
68 integer, intent(in), optional :: corner !< Array of starting
69 !! indices describing
70 !! where the data
71 !! will be written to.
72
73 call compressed_write(fileobj, variable_name, cdata, unlim_dim_level, corner)
74end subroutine compressed_write_0d_wrap
75
76
77!> @brief For variables without a compressed dimension, this routine simply wraps
78!! netcdf_write data. If the variable does have a compressed axis, the I/O root
79!! gathers the data from the rest of the ranks and then writes the combined data to
80!! the netcdf file.
81subroutine compressed_write_1d(fileobj, variable_name, cdata, unlim_dim_level, &
82 corner, edge_lengths)
83
84 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
85 character(len=*), intent(in) :: variable_name !< Variable name.
86 class(*), dimension(:), intent(in) :: cdata !< Compressed data that
87 !! will be gathered and
88 !! written to the
89 !! netcdf file.
90 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
91 !! dimension.
92 integer, dimension(1), intent(in), optional :: corner !< Array of starting
93 !! indices describing
94 !! where the data
95 !! will be written to.
96 integer, dimension(1), intent(in), optional :: edge_lengths !< The number of
97 !! elements that
98 !! will be written
99 !! in each dimension.
100
101 integer, dimension(2) :: compressed_dim_index !< index of the compressed dimension
102 !! compressed_dim_index(1) relative to cdata
103 !! compressed_dim_index(2) relative to the fileobj
104 integer, dimension(1) :: e !< "edges" number of points to read
105
106 integer(kind=i4_kind), dimension(:), allocatable :: buf_i4_kind !< Global buffer of data
107 integer(kind=i8_kind), dimension(:), allocatable :: buf_i8_kind !< Global buffer of data
108 real (kind=r4_kind), dimension(:), allocatable :: buf_r4_kind !< Global buffer of data
109 real (kind=r8_kind), dimension(:), allocatable :: buf_r8_kind !< Global buffer of data
110
111 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
112
113 append_error_msg = "compressed_write_1d: file:"//trim(fileobj%path)//" variable:"//trim(variable_name)
114
115 compressed_dim_index = get_variable_compressed_dimension_index(fileobj, variable_name)
116
117 if (compressed_dim_index(1) .eq. dimension_not_found) then
118 call netcdf_write_data(fileobj, variable_name, cdata, &
119 unlim_dim_level=unlim_dim_level, corner=corner, &
120 edge_lengths=edge_lengths)
121 return
122 endif
123
124 e(:) = shape(cdata)
125 !The root pe creates a buffer big enough to store the data:
126 if (fileobj%is_root) then
127 e(compressed_dim_index(1)) = sum(fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems)
128 select type(cdata)
129 type is (integer(kind=i4_kind))
130 call allocate_array(buf_i4_kind, e)
131 type is (integer(kind=i8_kind))
132 call allocate_array(buf_i8_kind, e)
133 type is (real(kind=r4_kind))
134 call allocate_array(buf_r4_kind, e)
135 type is (real(kind=r8_kind))
136 call allocate_array(buf_r8_kind, e)
137 class default
138 call error("unsupported variable type: "//trim(append_error_msg))
139 end select
140 endif
141
142 !Gather the data onto the I/O root and write it out.
143 select type(cdata)
144 type is (integer(kind=i4_kind))
145 call mpp_gather(cdata, size(cdata), buf_i4_kind, fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems, &
146 fileobj%pelist)
147 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_i4_kind, &
148 unlim_dim_level=unlim_dim_level)
149 type is (integer(kind=i8_kind))
150 call mpp_gather(cdata, size(cdata), buf_i8_kind, fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems, &
151 fileobj%pelist)
152 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_i8_kind, &
153 unlim_dim_level=unlim_dim_level)
154 type is (real(kind=r4_kind))
155 call mpp_gather(cdata, size(cdata), buf_r4_kind, fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems, &
156 fileobj%pelist)
157 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_r4_kind, &
158 unlim_dim_level=unlim_dim_level)
159 type is (real(kind=r8_kind))
160 call mpp_gather(cdata, size(cdata), buf_r8_kind, fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems, &
161 fileobj%pelist)
162 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_r8_kind, &
163 unlim_dim_level=unlim_dim_level)
164 class default
165 call error("unsupported variable type: "//trim(append_error_msg))
166 end select
167 if (fileobj%is_root) then
168 if (allocated(buf_i4_kind)) deallocate(buf_i4_kind)
169 if (allocated(buf_i8_kind)) deallocate(buf_i8_kind)
170 if (allocated(buf_r4_kind)) deallocate(buf_r4_kind)
171 if (allocated(buf_r8_kind)) deallocate(buf_r8_kind)
172 endif
173end subroutine compressed_write_1d
174
175
176!> @brief For variables without a compressed dimension, this routine simply wraps
177!! netcdf_write data. If the variable does have a compressed axis, the I/O root
178!! gathers the data from the rest of the ranks and then writes the combined data to
179!! the netcdf file.
180subroutine compressed_write_2d(fileobj, variable_name, cdata, unlim_dim_level, &
181 corner, edge_lengths)
182
183 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
184 character(len=*), intent(in) :: variable_name !< Variable name.
185 class(*), dimension(:,:), intent(in) :: cdata !< Compressed data that
186 !! will be gathered and
187 !! written to the
188 !! netcdf file.
189 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
190 !! dimension.
191 integer, dimension(2), intent(in), optional :: corner !< Array of starting
192 !! indices describing
193 !! where the data
194 !! will be written to.
195 integer, dimension(2), intent(in), optional :: edge_lengths !< The number of
196 !! elements that
197 !! will be written
198 !! in each dimension.
199
200 integer, dimension(2) :: compressed_dim_index !< index of the compressed dimension
201 !! compressed_dim_index(1) relative to cdata
202 !! compressed_dim_index(2) relative to the fileobj
203 integer, dimension(2) :: c !! corners of the data to read
204 integer, dimension(2) :: e !< "edges" number of points to read
205
206 integer(kind=i4_kind), dimension(:,:), allocatable :: buf_i4_kind !< Global buffer of data
207 integer(kind=i8_kind), dimension(:,:), allocatable :: buf_i8_kind !< Global buffer of data
208 real (kind=r4_kind), dimension(:,:), allocatable :: buf_r4_kind !< Global buffer of data
209 real (kind=r8_kind), dimension(:,:), allocatable :: buf_r8_kind !< Global buffer of data
210
211 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
212
213 integer :: index(1) !< index of the PE in the pelist
214
215 integer :: is !< Starting index of the first dimension
216 integer :: ie !< Ending index of the first dimension
217 integer :: js !< Starting index of the second dimension
218 integer :: je !< Ending index of the second dimension
219
220 append_error_msg = "compressed_write_2d: file:"//trim(fileobj%path)//" variable:"//trim(variable_name)
221
222 compressed_dim_index = get_variable_compressed_dimension_index(fileobj, variable_name)
223 if (compressed_dim_index(1) .eq. dimension_not_found) then
224 call netcdf_write_data(fileobj, variable_name, cdata, &
225 unlim_dim_level=unlim_dim_level, corner=corner, &
226 edge_lengths=edge_lengths)
227 return
228 endif
229
230 e(:) = shape(cdata)
231 !The root pe creates a buffer big enough to store the data:
232 if (fileobj%is_root) then
233 e(compressed_dim_index(1)) = sum(fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems)
234 select type(cdata)
235 type is (integer(kind=i4_kind))
236 call allocate_array(buf_i4_kind, e)
237 type is (integer(kind=i8_kind))
238 call allocate_array(buf_i8_kind, e)
239 type is (real(kind=r4_kind))
240 call allocate_array(buf_r4_kind, e)
241 type is (real(kind=r8_kind))
242 call allocate_array(buf_r8_kind, e)
243 class default
244 call error("unsupported variable type: "//trim(append_error_msg))
245 end select
246 endif
247
248 c(:) = 1
249 index = findloc(fileobj%pelist, mpp_pe())
250
251 c(compressed_dim_index(1)) = fileobj%compressed_dims(compressed_dim_index(2))%npes_corner(index(1))
252 e(compressed_dim_index(1)) = fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems(index(1))
253
254 if (compressed_dim_index(1) .eq. 1) then
255 is = c(compressed_dim_index(1))
256 ie = c(compressed_dim_index(1)) + e(compressed_dim_index(1)) - 1
257 js = 1
258 je = size(cdata,2)
259 else
260 js = c(compressed_dim_index(1))
261 je = c(compressed_dim_index(1)) + e(compressed_dim_index(1)) - 1
262 is = 1
263 ie = size(cdata,2)
264 endif
265
266 select type(cdata)
267 type is (integer(kind=i4_kind))
268 call mpp_gather(is, ie, js, je, fileobj%pelist, cdata, buf_i4_kind, fileobj%is_root)
269 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_i4_kind, &
270 unlim_dim_level=unlim_dim_level)
271 type is (integer(kind=i8_kind))
272 call mpp_gather(is, ie, js, je, fileobj%pelist, cdata, buf_i8_kind, fileobj%is_root)
273 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_i8_kind, &
274 unlim_dim_level=unlim_dim_level)
275 type is (real(kind=r4_kind))
276 call mpp_gather(is, ie, js, je, fileobj%pelist, cdata, buf_r4_kind, fileobj%is_root)
277 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_r4_kind, &
278 unlim_dim_level=unlim_dim_level)
279 type is (real(kind=r8_kind))
280 call mpp_gather(is, ie, js, je, fileobj%pelist, cdata, buf_r8_kind, fileobj%is_root)
281 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_r8_kind, &
282 unlim_dim_level=unlim_dim_level)
283 class default
284 call error("unsupported variable type: "//trim(append_error_msg))
285 end select
286
287 if (fileobj%is_root) then
288 if (allocated(buf_i4_kind)) deallocate(buf_i4_kind)
289 if (allocated(buf_i8_kind)) deallocate(buf_i8_kind)
290 if (allocated(buf_r4_kind)) deallocate(buf_r4_kind)
291 if (allocated(buf_r8_kind)) deallocate(buf_r8_kind)
292 endif
293end subroutine compressed_write_2d
294
295
296!> @brief For variables without a compressed dimension, this routine simply wraps
297!! netcdf_write data. If the variable does have a compressed axis, the I/O root
298!! gathers the data from the rest of the ranks and then writes the combined data to
299!! the netcdf file.
300subroutine compressed_write_3d(fileobj, variable_name, cdata, unlim_dim_level, &
301 corner, edge_lengths)
302
303 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
304 character(len=*), intent(in) :: variable_name !< Variable name.
305 class(*), contiguous, intent(in), target :: cdata(:,:,:) !< Compressed data that
306 !! will be gathered and
307 !! written to the
308 !! netcdf file.
309 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
310 !! dimension.
311 integer, dimension(3), intent(in), optional :: corner !< Array of starting
312 !! indices describing
313 !! where the data
314 !! will be written to.
315 integer, dimension(3), intent(in), optional :: edge_lengths !< The number of
316 !! elements that
317 !! will be written
318 !! in each dimension.
319
320 integer, dimension(2) :: compressed_dim_index !< index of the compressed dimension
321 !! compressed_dim_index(1) relative to cdata
322 !! compressed_dim_index(2) relative to the fileobj
323 integer, dimension(3) :: c !! corners of the data to read
324 integer, dimension(3) :: e !< "edges" number of points to read
325
326 integer(kind=i4_kind), dimension(:,:,:), allocatable :: buf_i4_kind !< Global buffer of data
327 integer(kind=i8_kind), dimension(:,:,:), allocatable :: buf_i8_kind !< Global buffer of data
328 real (kind=r4_kind), dimension(:,:,:), allocatable :: buf_r4_kind !< Global buffer of data
329 real (kind=r8_kind), dimension(:,:,:), allocatable :: buf_r8_kind !< Global buffer of data
330
331 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
332 class(*), pointer :: cdata_dummy(:,:,:,:)
333
334 integer :: index(1) !< index of the PE in the pelist
335
336 integer :: is !< Starting index of the first dimension
337 integer :: ie !< Ending index of the first dimension
338 integer :: js !< Starting index of the second dimension
339 integer :: je !< Ending index of the second dimension
340
341 append_error_msg = "compressed_write_3d: file:"//trim(fileobj%path)//" variable:"//trim(variable_name)
342
343 compressed_dim_index = get_variable_compressed_dimension_index(fileobj, variable_name)
344 if (compressed_dim_index(1) .eq. dimension_not_found) then
345 call netcdf_write_data(fileobj, variable_name, cdata, &
346 unlim_dim_level=unlim_dim_level, corner=corner, &
347 edge_lengths=edge_lengths)
348 return
349 else if (compressed_dim_index(1) .eq. 3) then
350 cdata_dummy(1:size(cdata,1), 1:size(cdata,2), 1:size(cdata,3), 1:1) => cdata(:,:,:)
351 call compressed_write_4d(fileobj, variable_name, cdata_dummy, unlim_dim_level)
352 endif
353
354 e(:) = shape(cdata)
355 !The root pe creates a buffer big enough to store the data:
356 if (fileobj%is_root) then
357 e(compressed_dim_index(1)) = sum(fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems)
358 select type(cdata)
359 type is (integer(kind=i4_kind))
360 call allocate_array(buf_i4_kind, e)
361 type is (integer(kind=i8_kind))
362 call allocate_array(buf_i8_kind, e)
363 type is (real(kind=r4_kind))
364 call allocate_array(buf_r4_kind, e)
365 type is (real(kind=r8_kind))
366 call allocate_array(buf_r8_kind, e)
367 class default
368 call error("unsupported variable type: "//trim(append_error_msg))
369 end select
370 endif
371
372 c(:) = 1
373 index = findloc(fileobj%pelist, mpp_pe())
374
375 c(compressed_dim_index(1)) = fileobj%compressed_dims(compressed_dim_index(2))%npes_corner(index(1))
376 e(compressed_dim_index(1)) = fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems(index(1))
377
378 if (compressed_dim_index(1) .eq. 1) then
379 is = c(compressed_dim_index(1))
380 ie = c(compressed_dim_index(1)) + e(compressed_dim_index(1)) - 1
381 js = 1
382 je = size(cdata,2)
383 else
384 js = c(compressed_dim_index(1))
385 je = c(compressed_dim_index(1)) + e(compressed_dim_index(1)) - 1
386 is = 1
387 ie = size(cdata,2)
388 endif
389
390 select type(cdata)
391 type is (integer(kind=i4_kind))
392 call mpp_gather(is, ie, js, je, size(cdata,3), &
393 fileobj%pelist, cdata, buf_i4_kind, fileobj%is_root)
394 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_i4_kind, &
395 unlim_dim_level=unlim_dim_level)
396 type is (integer(kind=i8_kind))
397 call mpp_gather(is, ie, js, je, size(cdata,3), &
398 fileobj%pelist, cdata, buf_i8_kind, fileobj%is_root)
399 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_i8_kind, &
400 unlim_dim_level=unlim_dim_level)
401 type is (real(kind=r4_kind))
402 call mpp_gather(is, ie, js, je, size(cdata,3), &
403 fileobj%pelist, cdata, buf_r4_kind, fileobj%is_root)
404 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_r4_kind, &
405 unlim_dim_level=unlim_dim_level)
406 type is (real(kind=r8_kind))
407 call mpp_gather(is, ie, js, je, size(cdata,3), &
408 fileobj%pelist, cdata, buf_r8_kind, fileobj%is_root)
409 if (fileobj%is_root) call netcdf_write_data(fileobj, variable_name, buf_r8_kind, &
410 unlim_dim_level=unlim_dim_level)
411 class default
412 call error("unsupported variable type: "//trim(append_error_msg))
413 end select
414
415 if (fileobj%is_root) then
416 if (allocated(buf_i4_kind)) deallocate(buf_i4_kind)
417 if (allocated(buf_i8_kind)) deallocate(buf_i8_kind)
418 if (allocated(buf_r4_kind)) deallocate(buf_r4_kind)
419 if (allocated(buf_r8_kind)) deallocate(buf_r8_kind)
420 endif
421end subroutine compressed_write_3d
422
423
424!> @brief For variables without a compressed dimension, this routine simply wraps
425!! netcdf_write data. If the variable does have a compressed axis, the I/O root
426!! gathers the data from the rest of the ranks and then writes the combined data to
427!! the netcdf file.
428subroutine compressed_write_4d(fileobj, variable_name, cdata, unlim_dim_level, &
429 corner, edge_lengths)
430
431 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
432 character(len=*), intent(in) :: variable_name !< Variable name.
433 class(*), dimension(:,:,:,:), intent(in) :: cdata !< Compressed data that
434 !! will be gathered and
435 !! written to the
436 !! netcdf file.
437 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
438 !! dimension.
439 integer, dimension(4), intent(in), optional :: corner !< Array of starting
440 !! indices describing
441 !! where the data
442 !! will be written to.
443 integer, dimension(4), intent(in), optional :: edge_lengths !< The number of
444 !! elements that
445 !! will be written
446 !! in each dimension.
447
448 integer, dimension(2) :: compressed_dim_index
449 integer, dimension(4) :: c
450 integer, dimension(4) :: e
451 integer :: i
452 integer(kind=i4_kind), dimension(:,:,:,:), allocatable :: buf_i4_kind
453 integer(kind=i8_kind), dimension(:,:,:,:), allocatable :: buf_i8_kind
454 real(kind=r4_kind), dimension(:,:,:,:), allocatable :: buf_r4_kind
455 real(kind=r8_kind), dimension(:,:,:,:), allocatable :: buf_r8_kind
456 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
457
458 append_error_msg = "compressed_write_4d: file:"//trim(fileobj%path)//" variable:"//trim(variable_name)
459
460 compressed_dim_index = get_variable_compressed_dimension_index(fileobj, variable_name)
461 if (compressed_dim_index(1) .eq. dimension_not_found) then
462 call netcdf_write_data(fileobj, variable_name, cdata, &
463 unlim_dim_level=unlim_dim_level, corner=corner, &
464 edge_lengths=edge_lengths)
465 return
466 endif
467
468 !Gather the data onto the I/O root and write it out.
469 if (fileobj%is_root) then
470 c(:) = 1
471 e(:) = shape(cdata)
472 do i = 1, size(fileobj%pelist)
473 c(compressed_dim_index(1)) = fileobj%compressed_dims(compressed_dim_index(2))%npes_corner(i)
474 e(compressed_dim_index(1)) = fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems(i)
475 if (i .eq. 1) then
476 call netcdf_write_data(fileobj, variable_name, cdata, &
477 unlim_dim_level=unlim_dim_level, corner=c, &
478 edge_lengths=e)
479 else
480 select type(cdata)
481 type is (integer(kind=i4_kind))
482 call allocate_array(buf_i4_kind, e)
483 call mpp_recv(buf_i4_kind, size(buf_i4_kind), fileobj%pelist(i), block=.true.)
484 call netcdf_write_data(fileobj, variable_name, buf_i4_kind, &
485 unlim_dim_level=unlim_dim_level, corner=c, &
486 edge_lengths=e)
487 deallocate(buf_i4_kind)
488 type is (integer(kind=i8_kind))
489 call allocate_array(buf_i8_kind, e)
490 call mpp_recv(buf_i8_kind, size(buf_i8_kind), fileobj%pelist(i), block=.true.)
491 call netcdf_write_data(fileobj, variable_name, buf_i8_kind, &
492 unlim_dim_level=unlim_dim_level, corner=c, &
493 edge_lengths=e)
494 deallocate(buf_i8_kind)
495 type is (real(kind=r4_kind))
496 call allocate_array(buf_r4_kind, e)
497 call mpp_recv(buf_r4_kind, size(buf_r4_kind), fileobj%pelist(i), block=.true.)
498 call netcdf_write_data(fileobj, variable_name, buf_r4_kind, &
499 unlim_dim_level=unlim_dim_level, corner=c, &
500 edge_lengths=e)
501 deallocate(buf_r4_kind)
502 type is (real(kind=r8_kind))
503 call allocate_array(buf_r8_kind, e)
504 call mpp_recv(buf_r8_kind, size(buf_r8_kind), fileobj%pelist(i), block=.true.)
505 call netcdf_write_data(fileobj, variable_name, buf_r8_kind, &
506 unlim_dim_level=unlim_dim_level, corner=c, &
507 edge_lengths=e)
508 deallocate(buf_r8_kind)
509 class default
510 call error("unsupported variable type: "//trim(append_error_msg))
511 end select
512 endif
513 enddo
514 else
515 select type(cdata)
516 type is (integer(kind=i4_kind))
517 call mpp_send(cdata, size(cdata), fileobj%io_root)
518 type is (integer(kind=i8_kind))
519 call mpp_send(cdata, size(cdata), fileobj%io_root)
520 type is (real(kind=r4_kind))
521 call mpp_send(cdata, size(cdata), fileobj%io_root)
522 type is (real(kind=r8_kind))
523 call mpp_send(cdata, size(cdata), fileobj%io_root)
524 class default
525 call error("unsupported variable type: "//trim(append_error_msg))
526 end select
527 call mpp_sync_self(check=event_send)
528 endif
529end subroutine compressed_write_4d
530
531
532!> @brief For variables without a compressed dimension, this routine simply wraps
533!! netcdf_write data. If the variable does have a compressed axis, the I/O root
534!! gathers the data from the rest of the ranks and then writes the combined data to
535!! the netcdf file.
536subroutine compressed_write_5d(fileobj, variable_name, cdata, unlim_dim_level, &
537 corner, edge_lengths)
538
539 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
540 character(len=*), intent(in) :: variable_name !< Variable name.
541 class(*), dimension(:,:,:,:,:), intent(in) :: cdata !< Compressed data that
542 !! will be gathered and
543 !! written to the
544 !! netcdf file.
545 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
546 !! dimension.
547 integer, dimension(5), intent(in), optional :: corner !< Array of starting
548 !! indices describing
549 !! where the data
550 !! will be written to.
551 integer, dimension(5), intent(in), optional :: edge_lengths !< The number of
552 !! elements that
553 !! will be written
554 !! in each dimension.
555
556 integer, dimension(2) :: compressed_dim_index
557 integer, dimension(5) :: c
558 integer, dimension(5) :: e
559 integer :: i
560 integer(kind=i4_kind), dimension(:,:,:,:,:), allocatable :: buf_i4_kind
561 integer(kind=i8_kind), dimension(:,:,:,:,:), allocatable :: buf_i8_kind
562 real(kind=r4_kind), dimension(:,:,:,:,:), allocatable :: buf_r4_kind
563 real(kind=r8_kind), dimension(:,:,:,:,:), allocatable :: buf_r8_kind
564 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
565
566 append_error_msg = "compressed_write_5d: file:"//trim(fileobj%path)//" variable:"//trim(variable_name)
567
568 compressed_dim_index = get_variable_compressed_dimension_index(fileobj, variable_name)
569 if (compressed_dim_index(1) .eq. dimension_not_found) then
570 call netcdf_write_data(fileobj, variable_name, cdata, &
571 unlim_dim_level=unlim_dim_level, corner=corner, &
572 edge_lengths=edge_lengths)
573 return
574 endif
575
576 !Gather the data onto the I/O root and write it out.
577 if (fileobj%is_root) then
578 c(:) = 1
579 e(:) = shape(cdata)
580 do i = 1, size(fileobj%pelist)
581 c(compressed_dim_index(1)) = fileobj%compressed_dims(compressed_dim_index(2))%npes_corner(i)
582 e(compressed_dim_index(1)) = fileobj%compressed_dims(compressed_dim_index(2))%npes_nelems(i)
583 if (i .eq. 1) then
584 call netcdf_write_data(fileobj, variable_name, cdata, &
585 unlim_dim_level=unlim_dim_level, corner=c, &
586 edge_lengths=e)
587 else
588 select type(cdata)
589 type is (integer(kind=i4_kind))
590 call allocate_array(buf_i4_kind, e)
591 call mpp_recv(buf_i4_kind, size(buf_i4_kind), fileobj%pelist(i), block=.true.)
592 call netcdf_write_data(fileobj, variable_name, buf_i4_kind, &
593 unlim_dim_level=unlim_dim_level, corner=c, &
594 edge_lengths=e)
595 deallocate(buf_i4_kind)
596 type is (integer(kind=i8_kind))
597 call allocate_array(buf_i8_kind, e)
598 call mpp_recv(buf_i8_kind, size(buf_i8_kind), fileobj%pelist(i), block=.true.)
599 call netcdf_write_data(fileobj, variable_name, buf_i8_kind, &
600 unlim_dim_level=unlim_dim_level, corner=c, &
601 edge_lengths=e)
602 deallocate(buf_i8_kind)
603 type is (real(kind=r4_kind))
604 call allocate_array(buf_r4_kind, e)
605 call mpp_recv(buf_r4_kind, size(buf_r4_kind), fileobj%pelist(i), block=.true.)
606 call netcdf_write_data(fileobj, variable_name, buf_r4_kind, &
607 unlim_dim_level=unlim_dim_level, corner=c, &
608 edge_lengths=e)
609 deallocate(buf_r4_kind)
610 type is (real(kind=r8_kind))
611 call allocate_array(buf_r8_kind, e)
612 call mpp_recv(buf_r8_kind, size(buf_r8_kind), fileobj%pelist(i), block=.true.)
613 call netcdf_write_data(fileobj, variable_name, buf_r8_kind, &
614 unlim_dim_level=unlim_dim_level, corner=c, &
615 edge_lengths=e)
616 deallocate(buf_r8_kind)
617 class default
618 call error("unsupported variable type: "//trim(append_error_msg))
619 end select
620 endif
621 enddo
622 else
623 select type(cdata)
624 type is (integer(kind=i4_kind))
625 call mpp_send(cdata, size(cdata), fileobj%io_root)
626 type is (integer(kind=i8_kind))
627 call mpp_send(cdata, size(cdata), fileobj%io_root)
628 type is (real(kind=r4_kind))
629 call mpp_send(cdata, size(cdata), fileobj%io_root)
630 type is (real(kind=r8_kind))
631 call mpp_send(cdata, size(cdata), fileobj%io_root)
632 class default
633 call error("unsupported variable type: "//trim(append_error_msg))
634 end select
635 call mpp_sync_self(check=event_send)
636 endif
637end subroutine compressed_write_5d
638
639
640!> @brief Wrapper to distinguish interfaces.
641subroutine compressed_write_1d_wrap(fileobj, variable_name, cdata, unlim_dim_level, &
642 corner, edge_lengths)
643
644 type(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
645 character(len=*), intent(in) :: variable_name !< Variable name.
646 class(*), dimension(:), intent(in) :: cdata !< Compressed data that
647 !! will be gathered and
648 !! written to the
649 !! netcdf file.
650 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
651 !! dimension.
652 integer, dimension(1), intent(in), optional :: corner !< Array of starting
653 !! indices describing
654 !! where the data
655 !! will be written to.
656 integer, dimension(1), intent(in), optional :: edge_lengths !< The number of
657 !! elements that
658 !! will be written
659 !! in each dimension.
660
661 call compressed_write(fileobj, variable_name, cdata, unlim_dim_level, corner, &
662 edge_lengths=edge_lengths)
663
664end subroutine compressed_write_1d_wrap
665
666
667!> @brief Wrapper to distinguish interfaces.
668subroutine compressed_write_2d_wrap(fileobj, variable_name, cdata, unlim_dim_level, &
669 corner, edge_lengths)
670
671 type(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
672 character(len=*), intent(in) :: variable_name !< Variable name.
673 class(*), dimension(:,:), intent(in) :: cdata !< Compressed data that
674 !! will be gathered and
675 !! written to the
676 !! netcdf file.
677 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
678 !! dimension.
679 integer, dimension(2), intent(in), optional :: corner !< Array of starting
680 !! indices describing
681 !! where the data
682 !! will be written to.
683 integer, dimension(2), intent(in), optional :: edge_lengths !< The number of
684 !! elements that
685 !! will be written
686 !! in each dimension.
687
688 call compressed_write(fileobj, variable_name, cdata, unlim_dim_level, corner, &
689 edge_lengths=edge_lengths)
690end subroutine compressed_write_2d_wrap
691
692
693!> @brief Wrapper to distinguish interfaces.
694subroutine compressed_write_3d_wrap(fileobj, variable_name, cdata, unlim_dim_level, &
695 corner, edge_lengths)
696
697 type(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
698 character(len=*), intent(in) :: variable_name !< Variable name.
699 class(*), dimension(:,:,:), intent(in) :: cdata !< Compressed data that
700 !! will be gathered and
701 !! written to the
702 !! netcdf file.
703 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
704 !! dimension.
705 integer, dimension(3), intent(in), optional :: corner !< Array of starting
706 !! indices describing
707 !! where the data
708 !! will be written to.
709 integer, dimension(3), intent(in), optional :: edge_lengths !< The number of
710 !! elements that
711 !! will be written
712 !! in each dimension.
713
714 call compressed_write(fileobj, variable_name, cdata, unlim_dim_level, corner, &
715 edge_lengths=edge_lengths)
716end subroutine compressed_write_3d_wrap
717
718
719!> @brief Wrapper to distinguish interfaces.
720subroutine compressed_write_4d_wrap(fileobj, variable_name, cdata, unlim_dim_level, &
721 corner, edge_lengths)
722
723 type(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
724 character(len=*), intent(in) :: variable_name !< Variable name.
725 class(*), dimension(:,:,:,:), intent(in) :: cdata !< Compressed data that
726 !! will be gathered and
727 !! written to the
728 !! netcdf file.
729 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
730 !! dimension.
731 integer, dimension(4), intent(in), optional :: corner !< Array of starting
732 !! indices describing
733 !! where the data
734 !! will be written to.
735 integer, dimension(4), intent(in), optional :: edge_lengths !< The number of
736 !! elements that
737 !! will be written
738 !! in each dimension.
739
740 call compressed_write(fileobj, variable_name, cdata, unlim_dim_level, corner, &
741 edge_lengths=edge_lengths)
742end subroutine compressed_write_4d_wrap
743
744
745!> @brief Wrapper to distinguish interfaces.
746subroutine compressed_write_5d_wrap(fileobj, variable_name, cdata, unlim_dim_level, &
747 corner, edge_lengths)
748
749 type(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
750 character(len=*), intent(in) :: variable_name !< Variable name.
751 class(*), dimension(:,:,:,:,:), intent(in) :: cdata !< Compressed data that
752 !! will be gathered and
753 !! written to the
754 !! netcdf file.
755 integer, intent(in), optional :: unlim_dim_level !< Level for the unlimited
756 !! dimension.
757 integer, dimension(5), intent(in), optional :: corner !< Array of starting
758 !! indices describing
759 !! where the data
760 !! will be written to.
761 integer, dimension(5), intent(in), optional :: edge_lengths !< The number of
762 !! elements that
763 !! will be written
764 !! in each dimension.
765
766 call compressed_write(fileobj, variable_name, cdata, unlim_dim_level, corner, &
767 edge_lengths=edge_lengths)
768end subroutine compressed_write_5d_wrap
769!> @}
subroutine compressed_write_1d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine compressed_write_3d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_5d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine compressed_write_0d(fileobj, variable_name, cdata, unlim_dim_level, corner)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_5d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_0d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner)
Wrapper to distinguish interfaces.
subroutine compressed_write_4d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_2d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine compressed_write_1d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_4d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine compressed_write_2d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_3d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine mpp_sync_self(pelist, check, request, msg_size, msg_type)
This is to check if current PE's outstanding puts are complete but we can't use shmem_fence because w...