subroutine ESMF_FieldRedistStoreNF(srcField, dstField, &
routehandle, keywordEnforcer, srcToDstTransposeMap, &
ignoreUnmatchedIndices, rc)
!
! !ARGUMENTS:
type(ESMF_Field), intent(in) :: srcField
type(ESMF_Field), intent(inout) :: dstField
type(ESMF_RouteHandle), intent(inout) :: routehandle
type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below
integer, intent(in), optional :: srcToDstTransposeMap(:)
logical, intent(in), optional :: ignoreUnmatchedIndices
integer, intent(out), optional :: rc
!
! !STATUS:
! \begin{itemize}
! \item\apiStatusCompatibleVersion{5.2.0r}
! \end{itemize}
!
! !DESCRIPTION:
!
! \label{FieldRedistStoreNF}
! {\tt ESMF\_FieldRedistStore()} is a collective method across all PETs of the
! current Component. The interface of the method is overloaded, allowing
! -- in principle -- each PET to call into {\tt ESMF\_FieldRedistStore()}
! through a different entry point. Restrictions apply as to which combinations
! are sensible. All other combinations result in ESMF run time errors. The
! complete semantics of the {\tt ESMF\_FieldRedistStore()} method, as provided
! through the separate entry points shown in \ref{FieldRedistStoreTK} and
! \ref{FieldRedistStoreNF}, is described in the following paragraphs as a whole.
!
! Store a Field redistribution operation from {\tt srcField} to {\tt dstField}.
! Interface \ref{FieldRedistStoreTK} allows PETs to specify a {\tt factor}
! argument. PETs not specifying a {\tt factor} argument call into interface
! \ref{FieldRedistStoreNF}. If multiple PETs specify the {\tt factor} argument,
! its type and kind, as well as its value must match across all PETs. If none
! of the PETs specify a {\tt factor} argument the default will be a factor of
! 1. The resulting factor is applied to all of the source data during
! redistribution, allowing scaling of the data, e.g. for unit transformation.
!
! Both {\tt srcField} and {\tt dstField} are interpreted as sequentialized
! vectors. The sequence is defined by the order of DistGrid dimensions and the
! order of tiles within the DistGrid or by user-supplied arbitrary sequence
! indices. See section \ref{Array:SparseMatMul} for details on the definition
! of {\em sequence indices}.
!
! Source Field, destination Field, and the factor may be of different
! <type><kind>. Further, source and destination Fields may differ in shape,
! however, the number of elements must match.
!
! The default redistribution operation, when {\tt srcToDstTransposeMap} is not
! specified, corresponds to the identity mapping: each element of the
! sequentialized source Field is copied to the sequentialized
! destination Field element in order.
!
! \begin{sloppypar}
! If the {\tt srcToDstTransposeMap} argument is provided it must be identical
! across all PETs. The {\tt srcToDstTransposeMap} allows source and destination
! Field dimensions to be transposed during the redistribution. To support this
! option, the number of source and destination Field dimensions must be equal
! and the size of the associated dimensions must match.
! See section \ref{Array:Redist:TransposeMode} for more details about the
! use of the {\tt srcToDstTransposeMap} argument.
! \end{sloppypar}
!
! It is erroneous to specify the identical Field object for {\tt srcField} and
! {\tt dstField} arguments.
!
! The routine returns an {\tt ESMF\_RouteHandle} that can be used to call
! {\tt ESMF\_FieldRedist()} on any pair of Fields that matches
! {\tt srcField} and {\tt dstField} in {\em type}, {\em kind}, and
! memory layout of the {\em gridded} dimensions. However, the size, number,
! and index order of {\em ungridded} dimensions may be different. See section
! \ref{RH:Reusability} for a more detailed discussion of RouteHandle
! reusability.
!
! This call is {\em collective} across the current VM.
!
! For examples and associated documentation regarding this method see Section
! \ref{sec:field:usage:redist_1dptr}.
!
! The arguments are:
! \begin{description}
! \item [srcField]
! {\tt ESMF\_Field} with source data.
! \item [dstField]
! {\tt ESMF\_Field} with destination data. The data in this Field may be
! destroyed by this call.
! \item [routehandle]
! Handle to the precomputed Route.
! \item [{[srcToDstTransposeMap]}]
! A list with as many entries as there are dimensions in {\tt srcField}, or
! {\tt tileCount} times this many entries.
! Each entry maps the corresponding {\tt srcField} dimension against the
! specified {\tt dstField} dimension. Mixing distributed and
! undistributed dimensions is supported.
! Negative entries reverse the order of elements along the specified
! dimension when going from source to destination.
! When providing $rank \times tileCount$ elements in
! {\tt srcToDstTransposeMap}, each block of size {\tt rank} is associated
! with the corresponding tile (in order), and interpreted as the
! tile-specific transpose map.
! \item [{[ignoreUnmatchedIndices]}]
! A logical flag that affects the behavior for when not all elements match
! between the {\tt srcField} and {\tt dstField} side. The default setting
! is {\tt .false.}, indicating that it is an error when such a situation is
! encountered. Setting {\tt ignoreUnmatchedIndices} to {\tt .true.} ignores
! unmatched indices.
! \item [{[rc]}]
! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors.
! \end{description}
!
!EOP
!----------------------------------------------------------------------------
! local variables as temporary input/output arguments
! internal local variables
integer :: localrc
type(ESMF_Array) :: srcArray, dstArray
! Initialize return code; assume routine not implemented
localrc = ESMF_RC_NOT_IMPL
if(present(rc)) rc = ESMF_RC_NOT_IMPL
! check variable: focus on field and farray
! rely on ArrayRedist to check the sanity of other variables
ESMF_INIT_CHECK_DEEP(ESMF_FieldGetInit, srcField, rc)
ESMF_INIT_CHECK_DEEP(ESMF_FieldGetInit, dstField, rc)
#define PROGRESSLOG_on
#ifdef PROGRESSLOG_on
call ESMF_LogWrite("ESMF_FieldRedistStoreNF(): Just entered routine.", ESMF_LOGMSG_INFO)
#endif
! Retrieve source and destination arrays.
call ESMF_FieldGet(srcField, array=srcArray, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_FieldGet(dstField, array=dstArray, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! perform redist through array
! For performance consideration:
! Rely on ArrayRedist to perform sanity checking of the other parameters
call ESMF_ArrayRedistStore(srcArray, dstArray, routehandle, &
srcToDstTransposeMap=srcToDstTransposeMap, &
ignoreUnmatchedIndices=ignoreUnmatchedIndices, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
#ifdef PROGRESSLOG_on
call ESMF_LogWrite("ESMF_FieldRedistStoreNF(): Leaving routine.", ESMF_LOGMSG_INFO)
#endif
if (present(rc)) rc = ESMF_SUCCESS
end subroutine ESMF_FieldRedistStoreNF