Skip to content
Merged
Changes from all 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
77 changes: 75 additions & 2 deletions streams/dshr_strdata_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -344,7 +344,7 @@ subroutine shr_strdata_init_from_inline(sdat, my_task, logunit, compname, &
sdat%logunit, trim(compname), sdat%mainproc, src_mask, dst_mask)

! Now finish initializing sdat
call shr_strdata_init(sdat, model_clock, stream_name, rc)
call shr_strdata_init(sdat, model_clock, stream_name, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

end subroutine shr_strdata_init_from_inline
Expand Down Expand Up @@ -500,7 +500,23 @@ subroutine shr_strdata_init(sdat, model_clock, stream_name, rc)
endif

! We do not yet have mask information, but we are required to set it here and change it later.
if (filename /= 'none') then
if (trim(sdat%stream(ns)%mapalgo) == 'redist' .and. filename /= 'none') then
! A 'redist' stream is on the same grid as the model, so reuse the already-built model mesh here
! instead of reading the stream mesh file and constructing a duplicate full ESMF mesh -- that
! duplicate is a large init-time memory and time cost at high resolution.
!
! Why this is safe (subtle): we call the stream->model mapping a 'redist', but in fact the stream
! decomposition is the same as the model decomposition, so the redist ends up simply being a copy.
! That is what makes it safe to set the stream mesh equal to the model mesh -- both the underlying
! grids/meshes and their decompositions are the same, so the result is identical to building a
! separate stream mesh from the file. (Only stream_mesh is read below -- no MeshSet/MeshDestroy --
! so sharing the model-mesh handle is also safe.)
!
! First verify the grids really match (redist requires identical global sizes) and error clearly if not.
call shr_strdata_check_redist_size(sdat, ns, trim(filename), rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
stream_mesh = sdat%model_mesh
Comment on lines +503 to +518

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

@mvertens - I especially would welcome your input on whether the assignment of stream_mesh = sdat%model_mesh seems general and safe here.

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

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

This looks okay to me. There is already a copy being done in the inline call - so this is should be fine as well.

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

@sjsprecious - I talked with @mvertens about this. She has convinced me that it's fine. However, the reason this is safe is potentially subtle, and I think warrants an additional comment - which could add to the comment above (or possibly replace some of the text in the existing comment). Can you add something like this, if this sounds good / correct to you?

We call this a redist, but in fact the model decomposition is the same as the stream decomposition, so the redist ends up simply being a copy. This makes it safe to set the stream mesh equal to the model mesh - since both the underlying grids/meshes and their decompositions are the same.

else if (filename /= 'none') then
stream_mesh = ESMF_MeshCreate(trim(filename), fileformat=ESMF_FILEFORMAT_ESMFMESH, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
endif
Expand Down Expand Up @@ -754,6 +770,63 @@ subroutine shr_strdata_init(sdat, model_clock, stream_name, rc)
endif
end subroutine shr_strdata_init

!===============================================================================
subroutine shr_strdata_check_redist_size(sdat, ns, meshfile, rc)

! Verify that a 'redist' stream really is on the model grid before reusing the
! model mesh in place of the stream mesh file. 'redist' is a 1-to-1 index copy,
! so the stream mesh file and the model mesh must have the same global element
! count; otherwise the redistribution is invalid. Only the elementCount
! dimension of the (ESMF-format) mesh file is read here -- cheap metadata, not
! the full mesh, so the memory/time saving of reusing the model mesh is kept.

! input/output variables
type(shr_strdata_type) , intent(in) :: sdat
integer , intent(in) :: ns ! stream index
character(len=*) , intent(in) :: meshfile ! stream mesh filename
integer , intent(out) :: rc

! local variables
type(file_desc_t) :: pioid
integer :: rcode
integer :: dimid
integer :: stream_gsize
integer :: old_handle ! previous setting of pio error handling
character(len=*), parameter :: subname = '(shr_strdata_check_redist_size) '
! ----------------------------------------------

rc = ESMF_SUCCESS

! Nothing to verify against if there is no mesh file.
if (len_trim(meshfile) == 0 .or. trim(meshfile) == 'none') return

! Read the global element count from the mesh file (collective over the pio
! subsystem; all ranks already hold the same meshfile and model_gsize, so the
! comparison below is identical on every rank).
rcode = pio_openfile(sdat%pio_subsystem, pioid, sdat%io_type, trim(meshfile), pio_nowrite)
call pio_seterrorhandling(pioid, PIO_BCAST_ERROR, old_handle)
rcode = pio_inq_dimid(pioid, 'elementCount', dimid)
if (rcode == PIO_NOERR) then
rcode = pio_inq_dimlen(pioid, dimid, stream_gsize)
end if
call pio_seterrorhandling(pioid, old_handle)
call pio_closefile(pioid)

! If the mesh file is not in ESMF format (no elementCount dimension) we cannot
! verify the size here; leave any problem for the downstream data read to catch.
if (rcode /= PIO_NOERR) return

if (stream_gsize /= sdat%model_gsize) then
call shr_log_error(subname//"ERROR: mapalgo='redist' for stream "//toString(ns)// &
" requires the model grid and the stream mesh file to be the same size, but the model"// &
" grid has "//toString(sdat%model_gsize)//" elements while the stream mesh file '"// &
trim(meshfile)//"' has elementCount = "//toString(stream_gsize)//". Use a non-redist"// &
" mapalgo (e.g. consf, consd, nn or bilinear) when the stream is on a different grid.", rc=rc)
return
end if

end subroutine shr_strdata_check_redist_size

!===============================================================================
subroutine shr_strdata_get_stream_nlev(sdat, stream_index, rc)

Expand Down
Loading