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