Skip to content
Draft
Show file tree
Hide file tree
Changes from 17 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
44 changes: 22 additions & 22 deletions src/framework/mpas_io_streams.F
Original file line number Diff line number Diff line change
Expand Up @@ -3335,14 +3335,14 @@ subroutine MPAS_writeStream(stream, frame, ierr)

if (field_cursor % field_type == FIELD_0D_INT) then

!call mpas_log_write('Writing out field '//trim(field_cursor % int0dField % fieldName))
!call mpas_log_write(' > is the field decomposed? $l', logicArgs=(/field_cursor % isDecomposed/))
!call mpas_log_write(' > outer dimension size $i', intArgs=(/field_cursor % totalDimSize/))
! call mpas_log_write('Writing out field '//trim(field_cursor % int0dField % fieldName))
! call mpas_log_write(' > is the field decomposed? $l', logicArgs=(/field_cursor % isDecomposed/))
! call mpas_log_write(' > outer dimension size $i', intArgs=(/field_cursor % totalDimSize/))

!call mpas_log_write('Copying field from first block')
!call mpas_log_write('Copying field from first block')
int0d_temp = field_cursor % int0dField % scalar

!call mpas_log_write('MGD calling MPAS_io_put_var now...')
!call mpas_log_write('MGD calling MPAS_io_put_var now...')
call MPAS_io_put_var(stream % fileHandle, field_cursor % int0dField % fieldName, int0d_temp, io_err)
call MPAS_io_err_mesg(stream % fileHandle % ioContext, io_err, .false.)
if (io_err /= MPAS_IO_NOERR .and. present(ierr)) ierr = MPAS_IO_ERR
Expand Down Expand Up @@ -3376,7 +3376,7 @@ subroutine MPAS_writeStream(stream, frame, ierr)
end if

if (field_cursor % int1dField % isVarArray) then
! I suspect we will never hit this code, as it doesn't make sense, really
! I suspect we will never hit this code, as it doesn't make sense, really
int0d_temp = field_1dint_ptr % array(j)
else
int1d_temp(i:i+ownedSize-1) = field_1dint_ptr % array(1:ownedSize)
Expand Down Expand Up @@ -3543,14 +3543,14 @@ subroutine MPAS_writeStream(stream, frame, ierr)

else if (field_cursor % field_type == FIELD_0D_REAL) then

!call mpas_log_write('Writing out field '//trim(field_cursor % real0dField % fieldName))
!call mpas_log_write(' > is the field decomposed? $l', logicArgs=(/field_cursor % isDecomposed/))
!call mpas_log_write(' > outer dimension size $i', intArgs=(/field_cursor % totalDimSize/))
!call mpas_log_write('Writing out field '//trim(field_cursor % real0dField % fieldName))
!call mpas_log_write(' > is the field decomposed? $l', logicArgs=(/field_cursor % isDecomposed/))
!call mpas_log_write(' > outer dimension size $i', intArgs=(/field_cursor % totalDimSize/))

!call mpas_log_write('Copying field from first block')
!call mpas_log_write('Copying field from first block')
real0d_temp = field_cursor % real0dField % scalar

!call mpas_log_write('MGD calling MPAS_io_put_var now...')
!call mpas_log_write('MGD calling MPAS_io_put_var now...')
call MPAS_io_put_var(stream % fileHandle, field_cursor % real0dField % fieldName, real0d_temp, io_err)
call MPAS_io_err_mesg(stream % fileHandle % ioContext, io_err, .false.)
if (io_err /= MPAS_IO_NOERR .and. present(ierr)) ierr = MPAS_IO_ERR
Expand Down Expand Up @@ -3584,7 +3584,7 @@ subroutine MPAS_writeStream(stream, frame, ierr)
end if

if (field_cursor % real1dField % isVarArray) then
! I suspect we will never hit this code, as it doesn't make sense, really
! I suspect we will never hit this code, as it doesn't make sense, really
real0d_temp = field_1dreal_ptr % array(j)
else
real1d_temp(i:i+ownedSize-1) = field_1dreal_ptr % array(1:ownedSize)
Expand Down Expand Up @@ -3891,24 +3891,24 @@ subroutine MPAS_writeStream(stream, frame, ierr)

else if (field_cursor % field_type == FIELD_0D_CHAR) then

!call mpas_log_write('Writing out field '//trim(field_cursor % char0dField % fieldName))
!call mpas_log_write(' > is the field decomposed? $l', logicArgs=(/field_cursor % isDecomposed/))
!call mpas_log_write(' > outer dimension size $i', intArgs=(/field_cursor % totalDimSize/))
!call mpas_log_write('Writing out field '//trim(field_cursor % char0dField % fieldName))
!call mpas_log_write(' > is the field decomposed? $l', logicArgs=(/field_cursor % isDecomposed/))
!call mpas_log_write(' > outer dimension size $i', intArgs=(/field_cursor % totalDimSize/))

!call mpas_log_write('Copying field from first block')
!call mpas_log_write('MGD calling MPAS_io_put_var now...')
!call mpas_log_write('Copying field from first block')
!call mpas_log_write('MGD calling MPAS_io_put_var now...')
call MPAS_io_put_var(stream % fileHandle, field_cursor % char0dField % fieldName, field_cursor % char0dField % scalar, io_err)
call MPAS_io_err_mesg(stream % fileHandle % ioContext, io_err, .false.)
if (io_err /= MPAS_IO_NOERR .and. present(ierr)) ierr = MPAS_IO_ERR

else if (field_cursor % field_type == FIELD_1D_CHAR) then

!call mpas_log_write('Writing out field '//trim(field_cursor % char1dField % fieldName))
!call mpas_log_write(' > is the field decomposed? $l', logicArgs=(/field_cursor % isDecomposed/))
!call mpas_log_write(' > outer dimension size $i', intArgs=(/field_cursor % totalDimSize/))
!call mpas_log_write('Writing out field '//trim(field_cursor % char1dField % fieldName))
!call mpas_log_write(' > is the field decomposed? $l', logicArgs=(/field_cursor % isDecomposed/))
!call mpas_log_write(' > outer dimension size $i', intArgs=(/field_cursor % totalDimSize/))

!call mpas_log_write('Copying field from first block')
!call mpas_log_write('MGD calling MPAS_io_put_var now...')
!call mpas_log_write('Copying field from first block')
!call mpas_log_write('MGD calling MPAS_io_put_var now...')
call MPAS_io_put_var(stream % fileHandle, field_cursor % char1dField % fieldName, field_cursor % char1dField % array, io_err)
call MPAS_io_err_mesg(stream % fileHandle % ioContext, io_err, .false.)
if (io_err /= MPAS_IO_NOERR .and. present(ierr)) ierr = MPAS_IO_ERR
Expand Down
166 changes: 165 additions & 1 deletion src/framework/mpas_stream_manager.F
Original file line number Diff line number Diff line change
Expand Up @@ -3084,7 +3084,7 @@ subroutine write_stream(manager, stream, blockID, timeLevel, mgLevel, forceWrite
character (len=StrKIND) :: now_string, time_string
character (len=StrKIND) :: temp_filename, actualWhen
character (len=StrKIND) :: err_string
logical :: ringing_alarm, recordSeek, swapRecords
logical :: ringing_alarm, recordSeek, swapRecords, use_cdf2
logical :: clobberRecords, clobberFiles, truncateFiles
integer :: maxRecords, tempRecord
integer :: local_ierr, threadNum
Expand Down Expand Up @@ -3182,6 +3182,7 @@ subroutine write_stream(manager, stream, blockID, timeLevel, mgLevel, forceWrite
!
! Build stream from pools of fields and attributes
!
call check_max_var_size(stream % field_pool, manager % allFields, timeLevel, mgLevel, use_cdf2)
allocate(stream % stream)
call MPAS_createStream(stream % stream, manager % ioContext, stream % filename, stream % io_type, MPAS_IO_WRITE, &
precision=stream % precision, clobberRecords=clobberRecords, &
Expand Down Expand Up @@ -4325,6 +4326,7 @@ end subroutine gen_random
timeLevel = 1
end if


select case (info % fieldType)
case (MPAS_POOL_REAL)
select case (info % nDims)
Expand Down Expand Up @@ -4495,6 +4497,123 @@ subroutine update_stream(stream, allFields, timeLevelIn, mgLevelIn, ierr) !{{{
end subroutine update_stream !}}}


subroutine check_max_var_size(field_pool, allFields, timeLevelIn, mgLevelIn, use_cdf2) !{{{
use iso_c_binding, only : c_sizeof

implicit none

type (mpas_pool_type), intent(inout) :: field_pool
type (MPAS_Pool_type), intent(in) :: allFields
integer, intent(in) :: timeLevelIn
integer, intent(in) :: mgLevelIn
logical, intent(out) :: use_cdf2

type (MPAS_Pool_iterator_type) :: itr
type (mpas_pool_field_info_type) :: info
integer :: timeLevel

type (field5DReal), pointer :: real5d
type (field4DReal), pointer :: real4d
type (field3DReal), pointer :: real3d
type (field2DReal), pointer :: real2d
type (field1DReal), pointer :: real1d
type (field0DReal), pointer :: real0d

type (field3DInteger), pointer :: int3d
type (field2DInteger), pointer :: int2d
type (field1DInteger), pointer :: int1d
type (field0DInteger), pointer :: int0d

type (field1DChar), pointer :: char1d
type (field0DChar), pointer :: char0d

integer :: max_field_size, field_size
integer :: tmp_int
real(kind=RKIND) :: tmp_real
character :: tmp_char
integer(kind=I8KIND) :: max_field_bytes, field_bytes
integer(kind=I8KIND), parameter :: int_size = c_sizeof(tmp_int)
integer(kind=I8KIND), parameter :: real_size = c_sizeof(tmp_real)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It may be possible to avoid creating temporary variables only to check their size. I think it might work, for example, to write

integer(kind=I8KIND), parameter :: real_size = c_sizeof(1.0_RKIND)

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This seems to work for RKIND, but I wasn't sure about how it works for character and int. Got some build errors.

integer(kind=I8KIND), parameter :: char_size = c_sizeof(tmp_char)
integer(kind=I8KIND), parameter :: cdf2_max_bytes = 2**32 - 4
character(len=StrKIND):: message

call mpas_pool_begin_iteration(field_pool)
max_field_bytes = 0
field_bytes = 0
do while ( mpas_pool_get_next_member(field_pool, itr) )

if (itr % memberType == MPAS_POOL_CONFIG) then

! To avoid accidentally matching in case statements below...
info % fieldType = -1
call mpas_pool_get_field_info(allFields, itr % memberName, info)
! Set time level to read
if (info % nTimeLevels >= timeLevelIn) then
timeLevel = timeLevelIn
else
timeLevel = 1
end if
call mpas_log_write('In check_max_var_size, field '//trim(itr % memberName)//' ndims: $i',intArgs=(/info % nDims/))
select case (info % fieldType)
case (MPAS_POOL_REAL)
select case (info % nDims)
case (0)
field_size = 1
case (1)
call mpas_pool_get_field(allFields, itr % memberName, real1d, timeLevel)
field_size = global_dim_size(real1d % block, real1d % dimNames, real1d % dimSizes, info % nDims)
case (2)
call mpas_pool_get_field(allFields, itr % memberName, real2d, timeLevel)
field_size = global_dim_size(real2d % block, real2d % dimNames, real2d % dimSizes, info % nDims)
case (3)
call mpas_pool_get_field(allFields, itr % memberName, real3d, timeLevel)
field_size = global_dim_size(real3d % block, real3d % dimNames, real3d % dimSizes, info % nDims)
case (4)
call mpas_pool_get_field(allFields, itr % memberName, real4d, timeLevel)
field_size = global_dim_size(real4d % block, real4d % dimNames, real4d % dimSizes, info % nDims)
case (5)
call mpas_pool_get_field(allFields, itr % memberName, real5d, timeLevel)
field_size = global_dim_size(real5d % block, real5d % dimNames, real5d % dimSizes, info % nDims)
end select
field_bytes = int(field_size, kind=I8KIND) * real_size
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Perhaps field_size should be declared as an I8KIND, and this type conversion should be eliminated.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Isn't there a type conversion involved whether it's in this function or inside global_dim_size . Not sure if mpas_pool_get_dimension can work with I8KIND int pointers?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I've now changed this. The only type conversion happens in global_dim_size , due to mpas_pool_get_dimension forcing a regular int. And I've changed to mpas_dmpar_sum_int8 now

case (MPAS_POOL_INTEGER)
select case (info % nDims)
case (0)
field_size = 1
case (1)
call mpas_pool_get_field(allFields, itr % memberName, int1d, timeLevel)
field_size = global_dim_size(int1d % block, int1d % dimNames, int1d % dimSizes, info % nDims)
case (2)
call mpas_pool_get_field(allFields, itr % memberName, int2d, timeLevel)
field_size = global_dim_size(int2d % block, int2d % dimNames, int2d % dimSizes, info % nDims)
case (3)
call mpas_pool_get_field(allFields, itr % memberName, int3d, timeLevel)
field_size = global_dim_size(int3d % block, int3d % dimNames, int3d % dimSizes, info % nDims)
end select
field_bytes = int(field_size, kind=I8KIND) * int_size
case (MPAS_POOL_CHARACTER)
select case (info % nDims)
case (0)
field_size = 1
case (1)
! call mpas_pool_get_field(allFields, itr % memberName, char1d, timeLevel)
call mpas_log_write('In check_max_var_size, unsupported type field1DChar.', MPAS_LOG_ERR)
end select
field_bytes = int(field_size, kind=I8KIND) * char_size
end select
max_field_bytes = merge(field_bytes, max_field_bytes, field_bytes > max_field_bytes)
write(message,fmt='(A,i14,A,i18,A,i18)') 'check_max_var_size.. field_bytes =',field_bytes,' max_field_bytes =',max_field_bytes, &
& ' cdf2_max_bytes =',cdf2_max_bytes
call mpas_log_write(message)
end if
end do
use_cdf2 = cdf2_max_bytes > max_field_bytes
call mpas_log_write('In check_max_var_size use_cdf2: $l',logicArgs=(/use_cdf2/))

end subroutine check_max_var_size !}}}


!-----------------------------------------------------------------------
! routine stream_active_pkg_check
!
Expand Down Expand Up @@ -4845,6 +4964,51 @@ logical function is_decomposed_dim(dimName) !{{{

end function is_decomposed_dim !}}}



integer function global_dim_size(block, dimNames, dimSizes, nDim) !{{{
use iso_c_binding, only : c_sizeof

implicit none

character(len=*), intent(in), dimension(:) :: dimNames
integer, intent(in), dimension(:) :: dimSizes
type(block_type), intent(in) :: block
integer, intent(in) :: nDim
integer, pointer :: block_dim_size
integer :: iDim, sum_block_dim_size
character(len=StrKIND):: message

global_dim_size = 1
do iDim = 1, nDim
if ( is_decomposed_dim(dimNames(iDim))) then
if (trim(dimNames(iDim)) == 'nCells') then
call mpas_pool_get_dimension(block % dimensions, 'nCellsSolve', block_dim_size)
else if (trim(dimNames(iDim)) == 'nEdges') then
call mpas_pool_get_dimension(block % dimensions, 'nEdgesSolve', block_dim_size)
else if (trim(dimNames(iDim)) == 'nVertices') then
call mpas_pool_get_dimension(block % dimensions, 'nVerticesSolve', block_dim_size)
else
global_dim_size = -1
end if

call mpas_dmpar_sum_int(block % domain % dminfo, block_dim_size, sum_block_dim_size)
call mpas_log_write('----- from global_dim_size ----- Dimname is decomposed '//trim(dimNames(iDim))//' local size $i global size $i ',intArgs=(/block_dim_size, sum_block_dim_size/))

global_dim_size = global_dim_size * sum_block_dim_size
else
global_dim_size = global_dim_size * dimSizes(iDim)
call mpas_log_write('----- from global_dim_size ----- Dimname is not decomposed '//trim(dimNames(iDim))//' Size $i',intArgs=(/dimSizes(iDim)/))
end if
end do

call mpas_log_write('----- from global_dim_size: cumulative global_dim_size $i',intArgs=(/global_dim_size/))


end function global_dim_size !}}}




!-----------------------------------------------------------------------
! routine prewrite_reindex
Expand Down