! $Id: ESMF_ArrayCreate.cppF90,v 1.84 2011/06/30 18:45:48 w6ws Exp $ ! ! Earth System Modeling Framework ! Copyright 2002-2011, University Corporation for Atmospheric Research, ! Massachusetts Institute of Technology, Geophysical Fluid Dynamics ! Laboratory, University of Michigan, National Centers for Environmental ! Prediction, Los Alamos National Laboratory, Argonne National Laboratory, ! NASA Goddard Space Flight Center. ! Licensed under the University of Illinois-NCSA License. ! !============================================================================== ^define ESMF_FILENAME "ESMF_ArrayCreate.F90" !============================================================================== #if 0 !============================================================================== ! TKR overloading macros #endif #include "ESMF_TypeKindRankMacros.hcppF90" !============================================================================== ! ESMF ArrayCreate module module ESMF_ArrayCreateMod ! !============================================================================== ! ! This file contains the ArrayCreate() methods. ! !------------------------------------------------------------------------------ ! INCLUDES ^include "ESMF.h" !------------------------------------------------------------------------------ !BOPI ! !MODULE: ESMF_ArrayCreateMod - Provide TKR overloading for ESMF_ArrayCreate() ! ! !DESCRIPTION: ! ! The code in this file is part of the {\tt ESMF\_Array} class Fortran API. ! ! !------------------------------------------------------------------------------ ! !USES: use ESMF_UtilTypesMod ! ESMF utility types use ESMF_InitMacrosMod ! ESMF initializer macros use ESMF_BaseMod ! ESMF base class use ESMF_LogErrMod ! ESMF error handling use ESMF_LocalArrayMod use ESMF_ArraySpecMod use ESMF_VMMod use ESMF_DELayoutMod use ESMF_DistGridMod use ESMF_RHandleMod use ESMF_F90InterfaceMod ! ESMF Fortran-C++ interface helper implicit none !------------------------------------------------------------------------------ ! !PRIVATE TYPES: private !------------------------------------------------------------------------------ ! ! ESMF_Array ! !------------------------------------------------------------------------------ ! Fortran class type to hold pointer to C++ object type ESMF_Array sequence private type(ESMF_Pointer) :: this ESMF_INIT_DECLARE end type !------------------------------------------------------------------------------ ! !PUBLIC TYPES: public ESMF_Array !------------------------------------------------------------------------------ ! ! !PUBLIC MEMBER FUNCTIONS: ! - ESMF-public methods: public operator(==) public operator(/=) public ESMF_ArrayCreate public ESMF_ArrayDestroy ! - ESMF-internal methods: public ESMF_ArrayGetInit public ESMF_ArraySetInitCreated public ESMF_ArrayGetThis public ESMF_ArraySetThis public ESMF_ArraySetThisNull public ESMF_ArraySetPLocalDeInternal public ESMF_ArrayCopyThis ! ESMF_ArrayCreateFromGrid - defined in ESMF_Grid.F90 due to circularity issues !EOPI !------------------------------------------------------------------------------ !------------------------------------------------------------------------------ ! The following line turns the CVS identifier string into a printable variable. character(*), parameter, private :: version = & '$Id: ESMF_ArrayCreate.cppF90,v 1.84 2011/06/30 18:45:48 w6ws Exp $' !============================================================================== ! ! INTERFACE BLOCKS ! !============================================================================== ! -------------------------- ESMF-public method ------------------------------- !BOPI ! !IROUTINE: ESMF_ArrayCreate -- Generic interface ! !INTERFACE: interface ESMF_ArrayCreate ! !PRIVATE MEMBER FUNCTIONS: ! TypeKindRankInterfaceMacro(ArrayCreateFromPtr) TypeKindRankInterfaceMacro(ArrayCreateFromPtrArb) TypeKindRankInterfaceMacro(ArrayCreateAssmdShape) module procedure ESMF_ArrayCreateLocalArray module procedure ESMF_ArrayCreateLocalArrayArb module procedure ESMF_ArrayCreateAllocate module procedure ESMF_ArrayCreateAllocateArb module procedure ESMF_ArrayCreateAllocateAS module procedure ESMF_ArrayCreateAllocateASArb module procedure ESMF_ArrayCreateCopy module procedure ESMF_ArrayCreateFromLA !1st prototype ! !DESCRIPTION: ! This interface provides a single entry point for the various ! types of {\tt ESMF\_ArrayCreate} functions. !EOPI end interface !=============================================================================== ! ArrayOperator() interfaces !=============================================================================== ! -------------------------- ESMF-public method ------------------------------- !BOP ! !IROUTINE: ESMF_ArrayAssignment(=) - Array assignment ! ! !INTERFACE: ! interface assignment(=) ! array1 = array2 ! ! !ARGUMENTS: ! type(ESMF_Array) :: array1 ! type(ESMF_Array) :: array2 ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Assign array1 as an alias to the same ESMF Array object in memory ! as array2. If array2 is invalid, then array1 will be equally invalid after ! the assignment. ! ! The arguments are: ! \begin{description} ! \item[array1] ! The {\tt ESMF\_Array} object on the left hand side of the assignment. ! \item[array2] ! The {\tt ESMF\_Array} object on the right hand side of the assignment. ! \end{description} ! !EOP !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- !BOP ! !IROUTINE: ESMF_ArrayOperator(==) - Array equality operator ! ! !INTERFACE: interface operator(==) ! if (array1 == array2) then ... endif ! OR ! result = (array1 == array2) ! !RETURN VALUE: ! logical :: result ! ! !ARGUMENTS: ! type(ESMF_Array), intent(in) :: array1 ! type(ESMF_Array), intent(in) :: array2 ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Test whether array1 and array2 are valid aliases to the same ESMF ! Array object in memory. For a more general comparison of two ESMF Arrays, ! going beyond the simple alias test, the ESMF\_ArrayMatch() function (not yet ! implemented) must be used. ! ! The arguments are: ! \begin{description} ! \item[array1] ! The {\tt ESMF\_Array} object on the left hand side of the equality ! operation. ! \item[array2] ! The {\tt ESMF\_Array} object on the right hand side of the equality ! operation. ! \end{description} ! !EOP module procedure ESMF_ArrayEQ end interface !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- !BOP ! !IROUTINE: ESMF_ArrayOperator(/=) - Array not equal operator ! ! !INTERFACE: interface operator(/=) ! if (array1 /= array2) then ... endif ! OR ! result = (array1 /= array2) ! !RETURN VALUE: ! logical :: result ! ! !ARGUMENTS: ! type(ESMF_Array), intent(in) :: array1 ! type(ESMF_Array), intent(in) :: array2 ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Test whether array1 and array2 are {\it not} valid aliases to the ! same ESMF Array object in memory. For a more general comparison of two ESMF ! Arrays, going beyond the simple alias test, the ESMF\_ArrayMatch() function ! (not yet implemented) must be used. ! ! The arguments are: ! \begin{description} ! \item[array1] ! The {\tt ESMF\_Array} object on the left hand side of the non-equality ! operation. ! \item[array2] ! The {\tt ESMF\_Array} object on the right hand side of the non-equality ! operation. ! \end{description} ! !EOP module procedure ESMF_ArrayNE end interface !------------------------------------------------------------------------------ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !------------------------------------------------------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayEQ()" !BOPI ! !IROUTINE: ESMF_ArrayEQ - Compare two Calendars for equality ! ! !INTERFACE: function ESMF_ArrayEQ(array1, array2) ! ! !RETURN VALUE: logical :: ESMF_ArrayEQ ! !ARGUMENTS: type(ESMF_Array), intent(in) :: array1 type(ESMF_Array), intent(in) :: array2 ! !DESCRIPTION: ! Test if both {\tt array1} and {\tt array2} alias the same ESMF Array ! object. ! !EOPI !------------------------------------------------------------------------------- ESMF_INIT_TYPE ainit1, ainit2 integer :: localrc1, localrc2 logical :: lval1, lval2 ! Use the following logic, rather than "ESMF-INIT-CHECK-DEEP", to gain ! init checks on both args, and in the case where both are uninitialized, ! to distinguish equality based on uninitialized type (uncreated, ! deleted). ! TODO: Consider moving this logic to C++: use Base class? status? ! Or replicate logic for C interface also. ! check inputs ainit1 = ESMF_ArrayGetInit(array1) ainit2 = ESMF_ArrayGetInit(array2) ! TODO: this line must remain split in two for SunOS f90 8.3 127000-03 if (ainit1 .eq. ESMF_INIT_CREATED .and. & ainit2 .eq. ESMF_INIT_CREATED) then ESMF_ArrayEQ = array1%this .eq. array2%this else ESMF_ArrayEQ = .false. endif end function ESMF_ArrayEQ !------------------------------------------------------------------------------- !------------------------------------------------------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayNE()" !BOPI ! !IROUTINE: ESMF_ArrayNE - Compare two Calendars for non-equality ! ! !INTERFACE: function ESMF_ArrayNE(array1, array2) ! ! !RETURN VALUE: logical :: ESMF_ArrayNE ! !ARGUMENTS: type(ESMF_Array), intent(in) :: array1 type(ESMF_Array), intent(in) :: array2 ! !DESCRIPTION: ! Test if both {\tt array1} and {\tt array2} alias the same ESMF Array ! object. ! !EOPI !------------------------------------------------------------------------------- ESMF_INIT_TYPE ainit1, ainit2 integer :: localrc1, localrc2 logical :: lval1, lval2 ! Use the following logic, rather than "ESMF-INIT-CHECK-DEEP", to gain ! init checks on both args, and in the case where both are uninitialized, ! to distinguish equality based on uninitialized type (uncreated, ! deleted). ESMF_ArrayNE = .not.ESMF_ArrayEQ(array1, array2) end function ESMF_ArrayNE !------------------------------------------------------------------------------- !=============================================================================== ! ArrayCreate() interfaces !=============================================================================== #define ArrayCreateFromPtrDoc() \ ! -------------------------- ESMF-public method ----------------------------- @\ !BOP @\ ! !IROUTINE: ESMF_ArrayCreate - Create Array object from Fortran array pointer @\ ! @\ ! !INTERFACE: @\ ! ! Private name; call using ESMF_ArrayCreate() @\ ! function ESMF_ArrayCreateFromPtr(distgrid, farrayPtr, & @\ ! keywordEnforcer, datacopyflag, distgridToArrayMap, computationalEdgeLWidth, & @\ ! computationalEdgeUWidth, computationalLWidth, & @\ ! computationalUWidth, totalLWidth, & @\ ! totalUWidth, name, rc) @\ ! @\ ! !ARGUMENTS: @\ ! type(ESMF_DistGrid), intent(in) :: distgrid @\ ! (ESMF_KIND_), pointer :: farrayPtr() @\ ! type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below @\ ! type(ESMF_DataCopy_Flag), intent(in), optional :: datacopyflag @\ ! integer, intent(in), optional :: distgridToArrayMap(:) @\ ! integer, intent(in), optional :: computationalEdgeLWidth(:) @\ ! integer, intent(in), optional :: computationalEdgeUWidth(:) @\ ! integer, intent(in), optional :: computationalLWidth(:) @\ ! integer, intent(in), optional :: computationalUWidth(:) @\ ! integer, intent(in), optional :: totalLWidth(:) @\ ! integer, intent(in), optional :: totalUWidth(:) @\ ! character (len=*), intent(in), optional :: name @\ ! integer, intent(out), optional :: rc @\ ! @\ ! !RETURN VALUE: @\ ! type(ESMF_Array) :: ESMF_ArrayCreateDataPtr @\ ! @\ ! !STATUS: @\ ! \apiStatusCompatible @\ ! @\ ! !DESCRIPTION: @\ ! Create an {\tt ESMF\_Array} object from existing local native Fortran @\ ! arrays with pointer attribute, according to distgrid. Besides @\ ! {\tt farrayPtr} each PET must issue this call with identical arguments in @\ ! order to create a consistent Array object. The bounds of the local arrays @\ ! are preserved by this call and determine the bounds of the total region of @\ ! the resulting Array object. Bounds of the DE-local exclusive regions are @\ ! set to be consistent with the total regions and the specified distgrid @\ ! argument. Bounds for Array dimensions that are not distributed are @\ ! automatically set to the bounds provided by {\tt farrayPtr}. @\ ! @\ ! This interface requires a 1 DE per PET decomposition. The Array object will @\ ! not be created and an error will be returned if this condition is not met. @\ ! @\ ! The not distributed Array dimensions form a tensor of rank = array.rank - @\ ! distgrid.dimCount. By default all tensor elements are associated with @\ ! stagger location 0. The widths of the computational region are set to @\ ! the provided value, or zero by default, for all tensor elements. Use @\ ! {\tt ESMF\_ArraySet()} to change these default settings after the @\ ! Array object has been created. @\ ! @\ ! The return value is the newly created {\tt ESMF\_Array} object. @\ ! @\ ! The arguments are: @\ ! \begin{description} @\ ! \item[distgrid] @\ ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and @\ ! distributed over DEs. The dimCount of distgrid must be smaller or equal @\ ! to the rank of {\tt farrayPtr}. @\ ! \item[farrayPtr] @\ ! Valid native Fortran array with pointer attribute. Memory must be @\ ! associated with the actual argument. The type/kind/rank information of @\ ! {\tt farrayPtr} will be used to set {\tt Array}|s properties @\ ! accordingly. The shape of {\tt farrayPtr} will be checked against the @\ ! information contained in the {\tt distgrid}. The bounds of @\ ! {\tt farrayPtr} will be preserved by this call and the bounds of the @\ ! resulting Array object are set accordingly. @\ ! \item[{[datacopyflag]}] @\ ! Specifies whether the Array object will reference the memory allocation @\ ! provided by {\tt farrayPtr} directly or will copy the data from @\ ! {\tt farrayPtr} into a new memory allocation. Valid options are @\ ! {\tt ESMF\_DATACOPY\_REFERENCE} (default) or {\tt ESMF\_DATACOPY\_VALUE}. @\ ! Depending on the specific situation the {\tt ESMF\_DATACOPY\_REFERENCE} option @\ ! may be unsafe when specifying an array slice for {\tt farrayPtr}. @\ ! \item[{[distgridToArrayMap]}] @\ ! List that contains as many elements as is indicated by @\ ! {\tt distgrids}|s dimCount. The list elements map each dimension of @\ ! the DistGrid object to a dimension in {\tt farrayPtr} by specifying the @\ ! appropriate Array dimension index. The default is to map all of @\ ! {\tt distgrid}|s dimensions against the lower dimensions of the @\ ! {\tt farrayPtr} argument in sequence, i.e. {\tt distgridToArrayMap = @\ ! (/1, 2, .../)}. @\ ! Unmapped {\tt farrayPtr} dimensions are not decomposed dimensions and @\ ! form a tensor of rank = Array.rank - DistGrid.dimCount. @\ ! All {\tt distgridToArrayMap} entries must be greater than or equal @\ ! to zero and smaller than or equal to the Array rank. It is erroneous @\ ! to specify the same entry multiple times unless it is zero. @\ ! If the Array rank is less than the DistGrid dimCount then the default @\ ! distgridToArrayMap will contain zeros for the dimCount - rank @\ ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} @\ ! indicates that the particular DistGrid dimension will be replicating @\ ! the Array across the DEs along this direction. @\ ! \item[{[computationalEdgeLWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the lower corner of the computational @\ ! region with respect to the lower corner of the exclusive region for DEs @\ ! that are located on the edge of a tile. @\ ! The default is a zero vector. @\ ! \item[{[computationalEdgeUWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the upper corner of the computational @\ ! region with respect to the upper corner of the exclusive region for DEs @\ ! that are located on the edge of a tile. @\ ! The default is a zero vector. @\ ! \item[{[computationalLWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the lower corner of the computational @\ ! region with respect to the lower corner of the exclusive region. @\ ! The default is a zero vector. @\ ! \item[{[computationalUWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the upper corner of the computational @\ ! region with respect to the upper corner of the exclusive region. @\ ! The default is a zero vector. @\ ! \item[{[totalLWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the lower corner of the total memory @\ ! region with respect to the lower corner of the computational region. @\ ! The default is to accommodate the union of exclusive and computational @\ ! region exactly. @\ ! \item[{[totalUWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the upper corner of the total memory @\ ! region with respect to the upper corner of the computational region. @\ ! The default is a vector that contains the remaining number of elements @\ ! in each direction as to fit the union of exclusive and computational @\ ! region into the memory region provided by the {\tt farrayPtr} argument. @\ ! \item[{[name]}] @\ ! Name of the Array object. @\ ! \item[{[rc]}] @\ ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. @\ ! \end{description} @\ ! @\ !EOP @\ !---------------------------------------------------------------------------- @\ #define ArrayCreateFromPtrMacro(mtype, mtypekind, mrank, mdim, mlen, mrng, mloc) \ ! -------------------------- ESMF-public method ----------------------------- @\ ^undef ESMF_METHOD @\ ^define ESMF_METHOD "ESMF_ArrayCreateFromPtr" @\ function ESMF_ArrayCreateFromPtr##mrank##D##mtypekind(distgrid, farrayPtr, & @\ keywordEnforcer, datacopyflag, distgridToArrayMap, & @\ computationalEdgeLWidth, computationalEdgeUWidth, & @\ computationalLWidth, computationalUWidth, totalLWidth, & @\ totalUWidth, name, rc) @\ @\ type(ESMF_Array) :: ESMF_ArrayCreateFromPtr##mrank##D##mtypekind @\ @\ type(ESMF_DistGrid), intent(in) :: distgrid @\ mtype (ESMF_KIND_##mtypekind),dimension(mdim),pointer :: farrayPtr @\ type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below @\ type(ESMF_DataCopy_Flag), intent(in), optional :: datacopyflag @\ integer, intent(in), optional :: distgridToArrayMap(:) @\ integer, intent(in), optional :: computationalEdgeLWidth(:) @\ integer, intent(in), optional :: computationalEdgeUWidth(:) @\ integer, intent(in), optional :: computationalLWidth(:) @\ integer, intent(in), optional :: computationalUWidth(:) @\ integer, intent(in), optional :: totalLWidth(:) @\ integer, intent(in), optional :: totalUWidth(:) @\ character (len=*), intent(in), optional :: name @\ integer, intent(out), optional :: rc @\ @\ ! Local variables @\ integer :: localrc ! local return code @\ type(ESMF_Array) :: array ! new Array @\ type(ESMF_LocalArray), allocatable :: localarrayList(:) ! helper variable @\ type(ESMF_DataCopy_Flag) :: datacopyflag_opt ! helper variable @\ @\ ! Initialize return code; assume failure until success is certain @\ localrc = ESMF_RC_NOT_IMPL @\ if (present(rc)) rc = ESMF_RC_NOT_IMPL @\ @\ ! Check init status of arguments @\ ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) @\ @\ ! Mark this Array object as invalid @\ array%this = ESMF_NULL_POINTER @\ @\ ! Set copy/ref behavior @\ datacopyflag_opt = ESMF_DATACOPY_REFERENCE ! default @\ if (present(datacopyflag)) datacopyflag_opt = datacopyflag @\ @\ ! Prepare the LocalArray list to be used in the ArrayCreate() call. @\ ! Must pass the datacopyflag_opt flag into LocalArrayCreate() to ensure @\ ! that the LocalArray object contained in the localarrayList, passed to @\ ! ArrayCreate() below, is suitable to be placed inside of the Array @\ ! object directly. ArrayCreate() does _not_ create LocalArray objects @\ ! internally for ESMF_INDEX_USER. @\ allocate(localarrayList(1)) @\ localarrayList(1) = ESMF_LocalArrayCreate(farrayPtr, datacopyflag=datacopyflag_opt, rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ @\ ! Create the Array object by calling into overloaded ArrayCreate() interf. @\ array = ESMF_ArrayCreate(distgrid, localarrayList, indexflag=ESMF_INDEX_USER, & @\ datacopyflag=datacopyflag_opt, distgridToArrayMap=distgridToArrayMap, & @\ computationalEdgeLWidth=computationalEdgeLWidth, & @\ computationalEdgeUWidth=computationalEdgeUWidth, & @\ computationalLWidth=computationalLWidth, & @\ computationalUWidth=computationalUWidth, & @\ totalLWidth=totalLWidth, totalUWidth=totalUWidth, & @\ name=name, rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ @\ ! Garbage collection. The ArrayCreate() call will have placed the @\ ! LocalArray object created above into the Array object. It will be @\ ! destroyed when the Array object is destroyed. Only the localarrayList @\ ! must be cleaned up here. @\ deallocate(localarrayList) @\ @\ ! Set return value @\ ESMF_ArrayCreateFromPtr##mrank##D##mtypekind = array @\ @\ ! Set init code @\ ESMF_INIT_SET_CREATED(ESMF_ArrayCreateFromPtr##mrank##D##mtypekind) @\ @\ ! Return successfully @\ if (present(rc)) rc = ESMF_SUCCESS @\ @\ end function ESMF_ArrayCreateFromPtr##mrank##D##mtypekind @\ !---------------------------------------------------------------------------- @\ TypeKindRankDeclarationMacro(ArrayCreateFromPtr) #define ArrayCreateFromPtrArbDoc() \ ! -------------------------- ESMF-public method ----------------------------- @\ !BOP @\ ! !IROUTINE: ESMF_ArrayCreate - Create Array object from Fortran array pointer w/ arbitrary seqIndices for halo@\ ! @\ ! !INTERFACE: @\ ! ! Private name; call using ESMF_ArrayCreate() @\ ! function ESMF_ArrayCreateFromPtrArb(distgrid, farrayPtr, & @\ ! haloSeqIndexList, keywordEnforcer datacopyflag, distgridToArrayMap, & @\ ! name, rc) @\ ! @\ ! !ARGUMENTS: @\ ! type(ESMF_DistGrid), intent(in) :: distgrid @\ ! (ESMF_KIND_), pointer :: farrayPtr() @\ ! integer, intent(in) :: haloSeqIndexList(:) @\ ! type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below @\ ! type(ESMF_DataCopy_Flag), intent(in), optional :: datacopyflag @\ ! integer, intent(in), optional :: distgridToArrayMap(:) @\ ! character (len=*), intent(in), optional :: name @\ ! integer, intent(out), optional :: rc @\ ! @\ ! !RETURN VALUE: @\ ! type(ESMF_Array) :: ESMF_ArrayCreateDataPtrArb @\ ! @\ ! !STATUS: @\ ! \apiStatusCompatible @\ ! @\ ! !DESCRIPTION: @\ ! Create an {\tt ESMF\_Array} object from existing local native Fortran @\ ! arrays with pointer attribute, according to distgrid. Besides @\ ! {\tt farrayPtr} each PET must issue this call with identical arguments in @\ ! order to create a consistent Array object. The bounds of the local arrays @\ ! are preserved by this call and determine the bounds of the total region of @\ ! the resulting Array object. Bounds of the DE-local exclusive regions are @\ ! set to be consistent with the total regions and the specified distgrid @\ ! argument. Bounds for Array dimensions that are not distributed are @\ ! automatically set to the bounds provided by {\tt farrayPtr}. @\ ! @\ ! This interface requires a 1 DE per PET decomposition. The Array object will @\ ! not be created and an error will be returned if this condition is not met. @\ ! @\ ! The not distributed Array dimensions form a tensor of rank = array.rank - @\ ! distgrid.dimCount. By default all tensor elements are associated with @\ ! stagger location 0. The widths of the computational region are set to @\ ! the provided value, or zero by default, for all tensor elements. Use @\ ! {\tt ESMF\_ArraySet()} to change these default settings after the @\ ! Array object has been created. @\ ! @\ ! The return value is the newly created {\tt ESMF\_Array} object. @\ ! @\ ! The arguments are: @\ ! \begin{description} @\ ! \item[distgrid] @\ ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and @\ ! distributed over DEs. The dimCount of distgrid must be smaller or equal @\ ! to the rank of {\tt farrayPtr}. @\ ! \item[farrayPtr] @\ ! Valid native Fortran array with pointer attribute. Memory must be @\ ! associated with the actual argument. The type/kind/rank information of @\ ! {\tt farrayPtr} will be used to set {\tt Array}|s properties @\ ! accordingly. The shape of {\tt farrayPtr} will be checked against the @\ ! information contained in the {\tt distgrid}. The bounds of @\ ! {\tt farrayPtr} will be preserved by this call and the bounds of the @\ ! resulting Array object are set accordingly. @\ ! \item[haloSeqIndexList] @\ ! One dimensional array containing sequence indices of local halo region. @\ ! The size (and content) of {\tt haloSeqIndexList} can (and typically will) @\ ! be different on each PET. @\ ! \item[{[datacopyflag]}] @\ ! Specifies whether the Array object will reference the memory allocation @\ ! provided by {\tt farrayPtr} directly or will copy the data from @\ ! {\tt farrayPtr} into a new memory allocation. Valid options are @\ ! {\tt ESMF\_DATACOPY\_REFERENCE} (default) or {\tt ESMF\_DATACOPY\_VALUE}. @\ ! Depending on the specific situation the {\tt ESMF\_DATACOPY\_REFERENCE} option @\ ! may be unsafe when specifying an array slice for {\tt farrayPtr}. @\ ! \item[{[distgridToArrayMap]}] @\ ! List that contains as many elements as is indicated by @\ ! {\tt distgrids}|s dimCount. The list elements map each dimension of @\ ! the DistGrid object to a dimension in {\tt farrayPtr} by specifying the @\ ! appropriate Array dimension index. The default is to map all of @\ ! {\tt distgrid}|s dimensions against the lower dimensions of the @\ ! {\tt farrayPtr} argument in sequence, i.e. {\tt distgridToArrayMap = @\ ! (/1, 2, .../)}. @\ ! Unmapped {\tt farrayPtr} dimensions are not decomposed dimensions and @\ ! form a tensor of rank = Array.rank - DistGrid.dimCount. @\ ! All {\tt distgridToArrayMap} entries must be greater than or equal @\ ! to zero and smaller than or equal to the Array rank. It is erroneous @\ ! to specify the same entry multiple times unless it is zero. @\ ! If the Array rank is less than the DistGrid dimCount then the default @\ ! distgridToArrayMap will contain zeros for the dimCount - rank @\ ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} @\ ! indicates that the particular DistGrid dimension will be replicating @\ ! the Array across the DEs along this direction. @\ ! \item[{[name]}] @\ ! Name of the Array object. @\ ! \item[{[rc]}] @\ ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. @\ ! \end{description} @\ ! @\ !EOP @\ !---------------------------------------------------------------------------- @\ #define ArrayCreateFromPtrArbMacro(mtype, mtypekind, mrank, mdim, mlen, mrng, mloc) \ ! -------------------------- ESMF-public method ----------------------------- @\ ^undef ESMF_METHOD @\ ^define ESMF_METHOD "ESMF_ArrayCreateFromPtrArb" @\ function ESMF_ArrayCreateFromPtrArb##mrank##D##mtypekind(distgrid, farrayPtr, & @\ haloSeqIndexList, keywordEnforcer, datacopyflag, distgridToArrayMap, & @\ name, rc) @\ @\ type(ESMF_Array) :: ESMF_ArrayCreateFromPtrArb##mrank##D##mtypekind @\ @\ type(ESMF_DistGrid), intent(in) :: distgrid @\ mtype (ESMF_KIND_##mtypekind),dimension(mdim),pointer :: farrayPtr @\ integer, intent(in) :: haloSeqIndexList(:) @\ type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below @\ type(ESMF_DataCopy_Flag), intent(in), optional :: datacopyflag @\ integer, intent(in), optional :: distgridToArrayMap(:) @\ character (len=*), intent(in), optional :: name @\ integer, intent(out), optional :: rc @\ @\ ! Local variables @\ integer :: localrc ! local return code @\ integer :: localHaloCount @\ integer :: dimCount, localDeCount @\ type(ESMF_DELayout) :: delayout @\ @\ ! Initialize return code; assume failure until success is certain @\ localrc = ESMF_RC_NOT_IMPL @\ if (present(rc)) rc = ESMF_RC_NOT_IMPL @\ @\ ! Check init status of arguments used directly @\ ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) @\ @\ ! Ensure incoming DistGrid is 1D @\ call ESMF_DistGridGet(distgrid, dimCount=dimCount, delayout=delayout, & @\ rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ @\ if (dimCount /= 1) then @\ call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & @\ msg="- this ArrayCreate() interface requires 1D DistGrid", & @\ ESMF_CONTEXT, rcToReturn=rc) @\ return @\ endif @\ @\ ! Ensure there is only 1 DE/PET @\ call ESMF_DELayoutGet(delayout, localDeCount=localDeCount, rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ @\ if (localDeCount /= 1) then @\ call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & @\ msg="- this ArrayCreate() interface only supports 1DE/PETs", & @\ ESMF_CONTEXT, rcToReturn=rc) @\ return @\ endif @\ @\ ! Prepare for actual ArrayCreate() call @\ localHaloCount = size(haloSeqIndexList) @\ @\ ! Call into conventional interface to create the Array object @\ ESMF_ArrayCreateFromPtrArb##mrank##D##mtypekind = & @\ ESMF_ArrayCreate(distgrid, farrayPtr, datacopyflag=datacopyflag, & @\ distgridToArrayMap=distgridToArrayMap, totalUWidth=(/localHaloCount/), & @\ name=name, rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ @\ ! Set the arbitrary sequence indices in the Array halo rim @\ call ESMF_ArraySetPLocalDeInternal( & @\ ESMF_ArrayCreateFromPtrArb##mrank##D##mtypekind, & @\ localDe=0, rimSeqIndex=haloSeqIndexList, rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ @\ ! Return successfully @\ if (present(rc)) rc = ESMF_SUCCESS @\ @\ end function ESMF_ArrayCreateFromPtrArb##mrank##D##mtypekind @\ !---------------------------------------------------------------------------- @\ TypeKindRankDeclarationMacro(ArrayCreateFromPtrArb) #define ArrayCreateAssmdShapeDoc() \ ! -------------------------- ESMF-public method ----------------------------- @\ !BOP @\ ! !IROUTINE: ESMF_ArrayCreate - Create Array object from Fortran array @\ ! @\ ! !INTERFACE: @\ ! ! Private name; call using ESMF_ArrayCreate() @\ ! function ESMF_ArrayCreateAssmdShape(distgrid, farray, & @\ ! indexflag, keywordEnforcer, datacopyflag, distgridToArrayMap, & @\ ! computationalEdgeLWidth, computationalEdgeUWidth, computationalLWidth, & @\ ! computationalUWidth, totalLWidth, & @\ ! totalUWidth, undistLBound, undistUBound, name, rc) @\ ! @\ ! !ARGUMENTS: @\ ! type(ESMF_DistGrid), intent(in) :: distgrid @\ ! (ESMF_KIND_), intent(in), target :: farray() @\ ! type(ESMF_Index_Flag), intent(in) :: indexflag @\ ! type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below @\ ! type(ESMF_DataCopy_Flag), intent(in), optional :: datacopyflag @\ ! integer, intent(in), optional :: distgridToArrayMap(:) @\ ! integer, intent(in), optional :: computationalEdgeLWidth(:) @\ ! integer, intent(in), optional :: computationalEdgeUWidth(:) @\ ! integer, intent(in), optional :: computationalLWidth(:) @\ ! integer, intent(in), optional :: computationalUWidth(:) @\ ! integer, intent(in), optional :: totalLWidth(:) @\ ! integer, intent(in), optional :: totalUWidth(:) @\ ! integer, intent(in), optional :: undistLBound(:) @\ ! integer, intent(in), optional :: undistUBound(:) @\ ! character (len=*), intent(in), optional :: name @\ ! integer, intent(out), optional :: rc @\ ! @\ ! !RETURN VALUE: @\ ! type(ESMF_Array) :: ESMF_ArrayCreateDataAssmdShape @\ ! @\ ! !STATUS: @\ ! \apiStatusCompatible @\ ! @\ ! !DESCRIPTION: @\ ! Create an {\tt ESMF\_Array} object from an existing local native Fortran @\ ! array according to distgrid. Besides {\tt farray} each PET must issue this @\ ! call with identical arguments in order to create a consistent Array object. @\ ! The local arrays provided must be dimensioned according to the DE-local @\ ! total region. Bounds of the exclusive regions are set as specified in the @\ ! distgrid argument. Bounds for Array dimensions that are not distributed can @\ ! be chosen freely using the {\tt undistLBound} and {\tt undistUBound} @\ ! arguments. @\ ! @\ ! This interface requires a 1 DE per PET decomposition. The Array object will @\ ! not be created and an error will be returned if this condition is not met. @\ ! @\ ! The not distributed Array dimensions form a tensor of rank = array.rank - @\ ! distgrid.dimCount. By default all tensor elements are associated with @\ ! stagger location 0. The widths of the computational region are set to @\ ! the provided value, or zero by default, for all tensor elements. Use @\ ! {\tt ESMF\_ArraySet()} to change these default settings after the @\ ! Array object has been created. @\ ! @\ ! The return value is the newly created {\tt ESMF\_Array} object. @\ ! @\ ! The arguments are: @\ ! \begin{description} @\ ! \item[distgrid] @\ ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and @\ ! distributed over DEs. The dimCount of distgrid must be smaller or equal @\ ! to the rank of farray. @\ ! \item[farray] @\ ! Valid native Fortran array, i.e. memory must be associated with the @\ ! actual argument. The type/kind/rank information of {\tt farray} will be @\ ! used to set {\tt Array}|s properties accordingly. The shape of @\ ! {\tt farray} will be checked against the information contained in the @\ ! {\tt distgrid}. @\ ! \item[indexflag] @\ ! Indicate how DE-local indices are defined. See section @\ ! \ref{const:indexflag} for a list of valid indexflag options. @\ ! \item[{[datacopyflag]}] @\ ! Specifies whether the Array object will reference the memory allocation @\ ! provided by {\tt farray} directly or will copy the data from @\ ! {\tt farray} into a new memory allocation. Valid options are @\ ! {\tt ESMF\_DATACOPY\_REFERENCE} (default) or {\tt ESMF\_DATACOPY\_VALUE}. @\ ! Depending on the specific situation the {\tt ESMF\_DATACOPY\_REFERENCE} option @\ ! may be unsafe when specifying an array slice for {\tt farray}. @\ ! \item[{[distgridToArrayMap]}] @\ ! List that contains as many elements as is indicated by @\ ! {\tt distgrids}|s dimCount. The list elements map each dimension of @\ ! the DistGrid object to a dimension in {\tt farray} by specifying the @\ ! appropriate Array dimension index. The default is to map all of @\ ! {\tt distgrid}|s dimensions against the lower dimensions of the @\ ! {\tt farray} argument in sequence, i.e. {\tt distgridToArrayMap = @\ ! (/1, 2, .../)}. @\ ! Unmapped {\tt farray} dimensions are not decomposed dimensions and @\ ! form a tensor of rank = Array.rank - DistGrid.dimCount. @\ ! All {\tt distgridToArrayMap} entries must be greater than or equal @\ ! to zero and smaller than or equal to the Array rank. It is erroneous @\ ! to specify the same entry multiple times unless it is zero. @\ ! If the Array rank is less than the DistGrid dimCount then the default @\ ! distgridToArrayMap will contain zeros for the dimCount - rank @\ ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} @\ ! indicates that the particular DistGrid dimension will be replicating @\ ! the Array across the DEs along this direction. @\ ! \item[{[computationalEdgeLWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the lower corner of the computational @\ ! region with respect to the lower corner of the exclusive region for DEs @\ ! that are located on the edge of a tile. @\ ! The default is a zero vector. @\ ! \item[{[computationalEdgeUWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the upper corner of the computational @\ ! region with respect to the upper corner of the exclusive region for DEs @\ ! that are located on the edge of a tile. @\ ! The default is a zero vector. @\ ! \item[{[computationalLWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the lower corner of the computational @\ ! region with respect to the lower corner of the exclusive region. @\ ! The default is a zero vector. @\ ! \item[{[computationalUWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the upper corner of the computational @\ ! region with respect to the upper corner of the exclusive region. @\ ! The default is a zero vector. @\ ! \item[{[totalLWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the lower corner of the total memory @\ ! region with respect to the lower corner of the computational region. @\ ! The default is to accommodate the union of exclusive and computational @\ ! region exactly. @\ ! \item[{[totalUWidth]}] @\ ! This vector argument must have dimCount elements, where dimCount is @\ ! specified in distgrid. It specifies the upper corner of the total memory @\ ! region with respect to the upper corner of the computational region. @\ ! The default is a vector that contains the remaining number of elements @\ ! in each direction as to fit the union of exclusive and computational @\ ! region into the memory region provided by the {\tt farray} argument. @\ ! \item[{[undistLBound]}] @\ ! Lower bounds for the array dimensions that are not distributed. @\ ! By default lbound is 1. @\ ! \item[{[undistUBound]}] @\ ! Upper bounds for the array dimensions that are not distributed. @\ ! By default ubound is equal to the extent of the corresponding @\ ! dimension in {\tt farray}. @\ ! \item[{[name]}] @\ ! Name of the Array object. @\ ! \item[{[rc]}] @\ ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. @\ ! \end{description} @\ ! @\ !EOP @\ !---------------------------------------------------------------------------- @\ #define ArrayCreateAssmdShapeMacro(mtype, mtypekind, mrank, mdim, mlen, mrng, mloc) \ ! -------------------------- ESMF-public method ----------------------------- @\ ^undef ESMF_METHOD @\ ^define ESMF_METHOD "ESMF_ArrayCreateAssmdShape" @\ function ESMF_ArrayCreateAssmdShape##mrank##D##mtypekind(distgrid, farray, & @\ indexflag, keywordEnforcer, datacopyflag, distgridToArrayMap, & @\ computationalEdgeLWidth, computationalEdgeUWidth, & @\ computationalLWidth, computationalUWidth, totalLWidth, & @\ totalUWidth, undistLBound, undistUBound, name, rc) @\ @\ type(ESMF_Array) :: ESMF_ArrayCreateAssmdShape##mrank##D##mtypekind @\ @\ type(ESMF_DistGrid), intent(in) :: distgrid @\ mtype (ESMF_KIND_##mtypekind),dimension(mdim),intent(in),target :: farray @\ type(ESMF_Index_Flag), intent(in) :: indexflag @\ type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below @\ type(ESMF_DataCopy_Flag), intent(in), optional :: datacopyflag @\ integer, intent(in), optional :: distgridToArrayMap(:) @\ integer, intent(in), optional :: computationalEdgeLWidth(:) @\ integer, intent(in), optional :: computationalEdgeUWidth(:) @\ integer, intent(in), optional :: computationalLWidth(:) @\ integer, intent(in), optional :: computationalUWidth(:) @\ integer, intent(in), optional :: totalLWidth(:) @\ integer, intent(in), optional :: totalUWidth(:) @\ integer, intent(in), optional :: undistLBound(:) @\ integer, intent(in), optional :: undistUBound(:) @\ character (len=*), intent(in), optional :: name @\ integer, intent(out), optional :: rc @\ @\ ! Local variables @\ integer :: localrc ! local return code @\ type(ESMF_Array) :: array ! new Array @\ type(ESMF_LocalArray), allocatable :: localarrayList(:) ! helper variable @\ mtype (ESMF_KIND_##mtypekind),dimension(mdim),pointer :: fptr ! helper variable @\ type(ESMF_DataCopy_Flag) :: datacopyflag_opt ! helper variable @\ @\ ! Initialize return code; assume failure until success is certain @\ localrc = ESMF_RC_NOT_IMPL @\ if (present(rc)) rc = ESMF_RC_NOT_IMPL @\ @\ ! Check init status of arguments @\ ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) @\ @\ ! Mark this Array object as invalid @\ array%this = ESMF_NULL_POINTER @\ @\ ! Set copy/ref behavior @\ datacopyflag_opt = ESMF_DATACOPY_REFERENCE ! default @\ if (present(datacopyflag)) datacopyflag_opt = datacopyflag @\ @\ ! Check that indexflag does not specify mode that is not supported here @\ if (indexflag==ESMF_INDEX_USER) then @\ call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & @\ msg="- ESMF_INDEX_USER not supported as input value", & @\ ESMF_CONTEXT, rcToReturn=rc) @\ return @\ endif @\ @\ ! Prepare the LocalArray list to be used in the ArrayCreate() call. @\ ! Do not pass the datacopyflag_opt flag into LocalArrayCreate() to prevent @\ ! unncessary data movements at this point. The ArrayCreate() below will @\ ! create a new LocalArray internally (to be used inside of the Array @\ ! object), and that is where the datacopyflag_opt will need to be considered. @\ allocate(localarrayList(1)) @\ fptr => farray @\ localarrayList(1) = ESMF_LocalArrayCreate(fptr, rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ @\ ! Create the Array object by calling into overloaded ArrayCreate() interf. @\ array = ESMF_ArrayCreate(distgrid, localarrayList, indexflag=indexflag, & @\ datacopyflag=datacopyflag_opt, distgridToArrayMap=distgridToArrayMap, & @\ computationalEdgeLWidth=computationalEdgeLWidth, & @\ computationalEdgeUWidth=computationalEdgeUWidth, & @\ computationalLWidth=computationalLWidth, & @\ computationalUWidth=computationalUWidth, & @\ totalLWidth=totalLWidth, totalUWidth=totalUWidth, & @\ undistLBound=undistLBound, undistUBound=undistUBound, & @\ name=name, rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ @\ ! Garbage collection. The ArrayCreate() call will have created a new @\ ! LocalArray object internally, so the temporary object created above @\ ! must be deleted here. @\ call ESMF_LocalArrayDestroy(localarrayList(1), rc=localrc) @\ if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & @\ ESMF_CONTEXT, rcToReturn=rc)) return @\ deallocate(localarrayList) @\ @\ ! Set return value @\ ESMF_ArrayCreateAssmdShape##mrank##D##mtypekind = array @\ @\ ! Set init code @\ ESMF_INIT_SET_CREATED(ESMF_ArrayCreateAssmdShape##mrank##D##mtypekind) @\ @\ ! Return successfully @\ if (present(rc)) rc = ESMF_SUCCESS @\ @\ end function ESMF_ArrayCreateAssmdShape##mrank##D##mtypekind @\ !---------------------------------------------------------------------------- @\ TypeKindRankDeclarationMacro(ArrayCreateAssmdShape) ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCreateLocalArray()" !BOP ! !IROUTINE: ESMF_ArrayCreate - Create Array object from a list of LocalArray objects ! !INTERFACE: ! Private name; call using ESMF_ArrayCreate() function ESMF_ArrayCreateLocalArray(distgrid, localarrayList, keywordEnforcer, & indexflag, datacopyflag, distgridToArrayMap, computationalEdgeLWidth, & computationalEdgeUWidth, computationalLWidth, computationalUWidth, & totalLWidth, totalUWidth, undistLBound, undistUBound, name, rc) ! ! !ARGUMENTS: type(ESMF_DistGrid), intent(in) :: distgrid type(ESMF_LocalArray), intent(in) :: localarrayList(:) type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below type(ESMF_Index_Flag), intent(in), optional :: indexflag type(ESMF_DataCopy_Flag), intent(in), optional :: datacopyflag integer, intent(in), optional :: distgridToArrayMap(:) integer, intent(in), optional :: computationalEdgeLWidth(:) integer, intent(in), optional :: computationalEdgeUWidth(:) integer, intent(in), optional :: computationalLWidth(:) integer, intent(in), optional :: computationalUWidth(:) integer, intent(in), optional :: totalLWidth(:) integer, intent(in), optional :: totalUWidth(:) integer, intent(in), optional :: undistLBound(:) integer, intent(in), optional :: undistUBound(:) character (len=*), intent(in), optional :: name integer, intent(out), optional :: rc ! ! !RETURN VALUE: type(ESMF_Array) :: ESMF_ArrayCreateLocalArray ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! \begin{sloppypar} ! Create an {\tt ESMF\_Array} object from existing {\tt ESMF\_LocalArray} ! objects according to distgrid. Besides {\tt localarrayList} each PET must issue ! this call with identical arguments in order to create a consistent Array ! object. The local arrays provided must be dimensioned according to the ! DE-local total region. Bounds of the exclusive regions are set as specified ! in the distgrid argument. Bounds for array dimensions that are not distributed ! can be chosen freely using the {\tt undistLBound} and {\tt undistUBound} ! arguments. ! \end{sloppypar} ! ! This interface is able to handle multiple DEs per PET. ! ! The not distributed Array dimensions form a tensor of rank = array.rank - ! distgrid.dimCount. By default all tensor elements are associated with ! stagger location 0. The widths of the computational region are set to ! the provided value, or zero by default, for all tensor elements. Use ! {\tt ESMF\_ArraySet()} to change these default settings after the ! Array object has been created. ! ! The return value is the newly created {\tt ESMF\_Array} object. ! ! The arguments are: ! \begin{description} ! \item[distgrid] ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and ! distributed over DEs. The dimCount of distgrid must be smaller or equal ! to the rank specified in arrayspec, otherwise a runtime ESMF error will be ! raised. ! \item[localarrayList] ! List of valid {\tt ESMF\_LocalArray} objects, i.e. memory must be ! associated with the actual arguments. The type/kind/rank information of ! all {\tt localarrayList} elements must be identical and will ! be used to set {\tt Array}|s properties accordingly. The shape of each ! {\tt localarrayList} element will be checked against the information ! contained in the {\tt distgrid}. ! \item[{[indexflag]}] ! Indicate how DE-local indices are defined. By default, the exclusive ! region of each DE is placed to start at the local index space origin, ! i.e. (1, 1, ..., 1). Alternatively the DE-local index space can be ! aligned with the global index space, if a global index space is well ! defined by the associated DistGrid. See section \ref{const:indexflag} ! for a list of valid indexflag options. ! \item[{[datacopyflag]}] ! Specifies whether the Array object will reference the memory allocation ! provided by {\tt farray} directly or will copy the data from ! {\tt farray} into a new memory allocation. Valid options are ! {\tt ESMF\_DATACOPY\_REFERENCE} (default) or {\tt ESMF\_DATACOPY\_VALUE}. ! Depending on the specific situation the {\tt ESMF\_DATACOPY\_REFERENCE} option ! may be unsafe when specifying an array slice for {\tt farray}. ! \item[{[distgridToArrayMap]}] ! List that contains as many elements as is indicated by ! {\tt distgrids}|s dimCount. The list elements map each dimension of ! the DistGrid object to a dimension in the {\tt localarrayList} elements ! by specifying the appropriate Array dimension index. The default is ! to map all of {\tt distgrid}|s dimensions against the lower dimensions ! of the {\tt localarrayList} elements in sequence, i.e. ! {\tt distgridToArrayMap = (/1, 2, .../)}. ! Unmapped dimensions in the {\tt localarrayList} elements are not ! decomposed dimensions and form a tensor of ! rank = Array.rank - DistGrid.dimCount. ! All {\tt distgridToArrayMap} entries must be greater than or equal ! to zero and smaller than or equal to the Array rank. It is erroneous ! to specify the same entry multiple times unless it is zero. ! If the Array rank is less than the DistGrid dimCount then the default ! distgridToArrayMap will contain zeros for the dimCount - rank ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} ! indicates that the particular DistGrid dimension will be replicating ! the Array across the DEs along this direction. ! \item[{[computationalEdgeLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the computational ! region with respect to the lower corner of the exclusive region for DEs ! that are located on the edge of a tile. ! \item[{[computationalEdgeUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the computational ! region with respect to the upper corner of the exclusive region for DEs ! that are located on the edge of a tile. ! \item[{[computationalLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the computational ! region with respect to the lower corner of the exclusive region. ! The default is a zero vector. ! \item[{[computationalUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the computational ! region with respect to the upper corner of the exclusive region. ! The default is a zero vector. ! \item[{[totalLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the total memory ! region with respect to the lower corner of the computational region. ! The default is to accommodate the union of exclusive and computational ! region exactly. ! \item[{[totalUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the total memory ! region with respect to the upper corner of the exclusive region. ! The default is a vector that contains the remaining number of elements ! in each direction as to fit the union of exclusive and computational ! region into the memory region provided by the {\tt localarrayList} argument. ! \item[{[undistLBound]}] ! Lower bounds for the array dimensions that are not distributed. ! By default lbound is 1. ! \item[{[undistUBound]}] ! Upper bounds for the array dimensions that are not distributed. ! By default ubound is equal to the extent of the corresponding ! dimension in {\tt localarrayList}. ! \item[{[name]}] ! Name of the Array object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOP !------------------------------------------------------------------------------ integer :: localrc ! local return code type(ESMF_Array) :: array ! opaque pointer to new C++ Array integer :: localarrayCount, i ! helper variable type(ESMF_InterfaceInt) :: distgridToArrayMapArg ! helper variable type(ESMF_InterfaceInt) :: computationalEdgeLWidthArg ! helper variable type(ESMF_InterfaceInt) :: computationalEdgeUWidthArg ! helper variable type(ESMF_InterfaceInt) :: computationalLWidthArg ! helper variable type(ESMF_InterfaceInt) :: computationalUWidthArg ! helper variable type(ESMF_InterfaceInt) :: totalLWidthArg ! helper variable type(ESMF_InterfaceInt) :: totalUWidthArg ! helper variable type(ESMF_InterfaceInt) :: undistLBoundArg ! helper variable type(ESMF_InterfaceInt) :: undistUBoundArg ! helper variable type(ESMF_Pointer), allocatable :: localarrayPointerList(:) ! helper variable integer :: len_name ! helper variable type(ESMF_DataCopy_Flag) :: datacopyflag_opt ! helper variable ! Initialize return code; assume failure until success is certain localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Determine the number of LocalArray elements in list localarrayCount = size(localarrayList) ! Check init status of arguments ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) do i=1, localarrayCount ESMF_INIT_CHECK_DEEP_SHORT(ESMF_LocalArrayGetInit, localarrayList(i), rc) enddo ! Copy C++ pointers of deep objects into a simple ESMF_Pointer array ! This is necessary in order to strip off the Fortran init check members ! when passing into C++ allocate(localarrayPointerList(localarrayCount)) do i=1, localarrayCount call ESMF_LocalArrayGetThis(localarrayList(i), localarrayPointerList(i), localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return enddo ! Set copy/ref behavior datacopyflag_opt = ESMF_DATACOPY_REFERENCE ! default if (present(datacopyflag)) datacopyflag_opt = datacopyflag ! Check that indexflag does not specify mode that is not supported here !TODO: cannot check for this here because ArrayCreateFromPtr() calls into here! ! if (present(indexflag)) then ! if (indexflag==ESMF_INDEX_USER) then ! call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & ! "- ESMF_INDEX_USER not supported as input value", & ! ESMF_CONTEXT, rcToReturn=rc) ! return ! endif ! endif ! Deal with (optional) array arguments distgridToArrayMapArg = ESMF_InterfaceIntCreate(distgridToArrayMap, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return computationalEdgeLWidthArg = & ESMF_InterfaceIntCreate(computationalEdgeLWidth, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return computationalEdgeUWidthArg = & ESMF_InterfaceIntCreate(computationalEdgeUWidth, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return computationalLWidthArg = ESMF_InterfaceIntCreate(computationalLWidth, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return computationalUWidthArg = ESMF_InterfaceIntCreate(computationalUWidth, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return totalLWidthArg = ESMF_InterfaceIntCreate(totalLWidth, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return totalUWidthArg = ESMF_InterfaceIntCreate(totalUWidth, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return undistLBoundArg = ESMF_InterfaceIntCreate(undistLBound, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return undistUBoundArg = ESMF_InterfaceIntCreate(undistUBound, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Mark this Array object as invalid array%this = ESMF_NULL_POINTER ! Call into the C++ interface, which will sort out optional arguments ! Optional name argument requires separate calls into C++ if (present(name)) then len_name = len(name) call c_ESMC_ArrayCreateLocalArray(array, localarrayPointerList, localarrayCount, & distgrid, datacopyflag_opt, distgridToArrayMapArg, & computationalEdgeLWidthArg, computationalEdgeUWidthArg, & computationalLWidthArg, computationalUWidthArg, totalLWidthArg, & totalUWidthArg, indexflag, & undistLBoundArg, undistUBoundArg, name, len_name, localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return else len_name = 0 call c_ESMC_ArrayCreateLocalArray(array, localarrayPointerList, localarrayCount, & distgrid, datacopyflag_opt, distgridToArrayMapArg, & computationalEdgeLWidthArg, computationalEdgeUWidthArg, & computationalLWidthArg, computationalUWidthArg, totalLWidthArg, & totalUWidthArg, indexflag, & undistLBoundArg, undistUBoundArg, "", len_name, localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return endif ! Garbage collection deallocate(localarrayPointerList) call ESMF_InterfaceIntDestroy(distgridToArrayMapArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(computationalEdgeLWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(computationalEdgeUWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(computationalLWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(computationalUWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(totalLWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(totalUWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(undistLBoundArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(undistUBoundArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Set return value ESMF_ArrayCreateLocalArray = array ! Set init code ESMF_INIT_SET_CREATED(ESMF_ArrayCreateLocalArray) ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end function ESMF_ArrayCreateLocalArray !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCreateLocalArrayArb()" !BOP ! !IROUTINE: ESMF_ArrayCreate - Create Array object from a list of LocalArray objects w/ arbitrary seqIndices for halo ! !INTERFACE: ! Private name; call using ESMF_ArrayCreate() function ESMF_ArrayCreateLocalArrayArb(distgrid, localarrayList, & haloSeqIndexList, keywordEnforcer, datacopyflag, distgridToArrayMap, undistLBound, & undistUBound, name, rc) ! ! !ARGUMENTS: type(ESMF_DistGrid), intent(in) :: distgrid type(ESMF_LocalArray), intent(in) :: localarrayList(:) integer, intent(in) :: haloSeqIndexList(:) type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below type(ESMF_DataCopy_Flag), intent(in), optional :: datacopyflag integer, intent(in), optional :: distgridToArrayMap(:) integer, intent(in), optional :: undistLBound(:) integer, intent(in), optional :: undistUBound(:) character (len=*), intent(in), optional :: name integer, intent(out), optional :: rc ! ! !RETURN VALUE: type(ESMF_Array) :: ESMF_ArrayCreateLocalArrayArb ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Create an {\tt ESMF\_Array} object from existing {\tt ESMF\_LocalArray} ! objects according to distgrid. Each PET must issue this call in unison ! in order to create a consistent Array object. The local arrays provided must ! be dimensioned according to the DE-local total region. Bounds of the ! exclusive regions are set as specified in the distgrid argument. Bounds ! for array dimensions that are not distributed can be chosen freely using ! the {\tt undistLBound} and {\tt undistUBound} arguments. ! ! The return value is the newly created {\tt ESMF\_Array} object. ! ! The arguments are: ! \begin{description} ! \item[distgrid] ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and ! distributed over DEs. The dimCount of distgrid must be smaller or equal ! to the rank specified in arrayspec, otherwise a runtime ESMF error will be ! raised. ! \item[localarrayList] ! List of valid {\tt ESMF\_LocalArray} objects, i.e. memory must be ! associated with the actual arguments. The type/kind/rank information of ! all {\tt localarrayList} elements must be identical and will ! be used to set {\tt Array}|s properties accordingly. The shape of each ! {\tt localarrayList} element will be checked against the information ! contained in the {\tt distgrid}. ! \item[haloSeqIndexList] ! One dimensional array containing sequence indices of local halo region. ! The size (and content) of {\tt haloSeqIndexList} can (and typically will) ! be different on each PET. ! \item[{[datacopyflag]}] ! Specifies whether the Array object will reference the memory allocation ! provided by {\tt farray} directly or will copy the data from ! {\tt farray} into a new memory allocation. Valid options are ! {\tt ESMF\_DATACOPY\_REFERENCE} (default) or {\tt ESMF\_DATACOPY\_VALUE}. ! Depending on the specific situation the {\tt ESMF\_DATACOPY\_REFERENCE} option ! may be unsafe when specifying an array slice for {\tt farray}. ! \item[{[distgridToArrayMap]}] ! List that contains as many elements as is indicated by ! {\tt distgrids}|s dimCount. The list elements map each dimension of ! the DistGrid object to a dimension in the {\tt localarrayList} elements ! by specifying the appropriate Array dimension index. The default is ! to map all of {\tt distgrid}|s dimensions against the lower dimensions ! of the {\tt localarrayList} elements in sequence, i.e. ! {\tt distgridToArrayMap = (/1, 2, .../)}. ! Unmapped dimensions in the {\tt localarrayList} elements are not ! decomposed dimensions and form a tensor of ! rank = Array.rank - DistGrid.dimCount. ! All {\tt distgridToArrayMap} entries must be greater than or equal ! to zero and smaller than or equal to the Array rank. It is erroneous ! to specify the same entry multiple times unless it is zero. ! If the Array rank is less than the DistGrid dimCount then the default ! distgridToArrayMap will contain zeros for the dimCount - rank ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} ! indicates that the particular DistGrid dimension will be replicating ! the Array across the DEs along this direction. ! \item[{[undistLBound]}] ! Lower bounds for the array dimensions that are not distributed. ! By default lbound is 1. ! \item[{[undistUBound]}] ! Upper bounds for the array dimensions that are not distributed. ! By default ubound is equal to the extent of the corresponding ! dimension in {\tt localarrayList}. ! \item[{[name]}] ! Name of the Array object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOP !------------------------------------------------------------------------------ integer :: localrc ! local return code integer :: localHaloCount integer :: dimCount, localDeCount type(ESMF_DELayout) :: delayout ! Initialize return code; assume failure until success is certain localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Check init status of arguments used directly ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) ! Ensure incoming DistGrid is 1D call ESMF_DistGridGet(distgrid, dimCount=dimCount, delayout=delayout, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return if (dimCount /= 1) then call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & msg="- this ArrayCreate() interface requires 1D DistGrid", & ESMF_CONTEXT, rcToReturn=rc) return endif ! Ensure there is only 1 DE/PET call ESMF_DELayoutGet(delayout, localDeCount=localDeCount, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return if (localDeCount /= 1) then call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & msg="- this ArrayCreate() interface only supports 1DE/PETs", & ESMF_CONTEXT, rcToReturn=rc) return endif ! Prepare for actual ArrayCreate() call localHaloCount = size(haloSeqIndexList) ! Call into conventional interface to create the Array object ESMF_ArrayCreateLocalArrayArb = ESMF_ArrayCreate(distgrid, localarrayList, & datacopyflag=datacopyflag, distgridToArrayMap=distgridToArrayMap, & totalUWidth=(/localHaloCount/), undistLBound=undistLBound, & undistUBound=undistUBound, name=name, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Set the arbitrary sequence indices in the Array halo rim call ESMF_ArraySetPLocalDeInternal(ESMF_ArrayCreateLocalArrayArb, & localDe=0, rimSeqIndex=haloSeqIndexList, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end function ESMF_ArrayCreateLocalArrayArb !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCreateAllocate()" !BOP ! !IROUTINE: ESMF_ArrayCreate - Create Array object from typekind (allocate memory) ! !INTERFACE: ! Private name; call using ESMF_ArrayCreate() function ESMF_ArrayCreateAllocate(distgrid, typekind, keywordEnforcer, indexflag, & distgridToArrayMap, computationalEdgeLWidth, computationalEdgeUWidth, & computationalLWidth, computationalUWidth, totalLWidth, totalUWidth, & undistLBound, undistUBound, name, rc) ! ! !ARGUMENTS: type(ESMF_DistGrid), intent(in) :: distgrid type(ESMF_TypeKind_Flag), intent(in) :: typekind type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below type(ESMF_Index_Flag), intent(in), optional :: indexflag integer, intent(in), optional :: distgridToArrayMap(:) integer, intent(in), optional :: computationalEdgeLWidth(:) integer, intent(in), optional :: computationalEdgeUWidth(:) integer, intent(in), optional :: computationalLWidth(:) integer, intent(in), optional :: computationalUWidth(:) integer, intent(in), optional :: totalLWidth(:) integer, intent(in), optional :: totalUWidth(:) integer, intent(in), optional :: undistLBound(:) integer, intent(in), optional :: undistUBound(:) character (len=*), intent(in), optional :: name integer, intent(out), optional :: rc ! ! !RETURN VALUE: type(ESMF_Array) :: ESMF_ArrayCreateAllocate ! ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Create an {\tt ESMF\_Array} object and allocate uninitialized data space ! according to typekind and distgrid. The Array rank is indirectly determined ! by the incoming information. Each PET must issue this call in unison in order ! to create a consistent Array object. DE-local allocations are made according ! to the total region defined by the {\tt distgrid} and the optional {\tt Width} ! arguments. ! ! The return value is the newly created {\tt ESMF\_Array} object. ! ! The arguments are: ! \begin{description} ! \item[distgrid] ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and ! distributed over DEs. The dimCount of distgrid must be smaller or equal ! to the rank specified in arrayspec, otherwise a runtime ESMF error will be ! raised. ! \item[typekind] ! The typekind of the Array. ! \item[{[indexflag]}] ! Indicate how DE-local indices are defined. By default, the exclusive ! region of each DE is placed to start at the local index space origin, ! i.e. (1, 1, ..., 1). Alternatively the DE-local index space can be ! aligned with the global index space, if a global index space is well ! defined by the associated DistGrid. See section \ref{const:indexflag} ! for a list of valid indexflag options. ! \item[{[distgridToArrayMap]}] ! List that contains as many elements as is indicated by ! {\tt distgrids}|s dimCount. The list elements map each dimension of ! the DistGrid object to a dimension in the newly allocated Array object ! by specifying the appropriate Array dimension index. The default is ! to map all of {\tt distgrid}|s dimensions against the lower dimensions ! of the Array object in sequence, i.e. {\tt distgridToArrayMap = ! (/1, 2, .../)}. ! Unmapped dimensions in the Array object are not decomposed dimensions ! and form a tensor of rank = Array.rank - DistGrid.dimCount. ! All {\tt distgridToArrayMap} entries must be greater than or equal ! to zero and smaller than or equal to the Array rank. It is erroneous ! to specify the same entry multiple times unless it is zero. ! If the Array rank is less than the DistGrid dimCount then the default ! distgridToArrayMap will contain zeros for the dimCount - rank ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} ! indicates that the particular DistGrid dimension will be replicating ! the Array across the DEs along this direction. ! \item[{[computationalEdgeLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the computational ! region with respect to the lower corner of the exclusive region for DEs ! that are located on the edge of a tile. ! \item[{[computationalEdgeUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the computational ! region with respect to the upper corner of the exclusive region for DEs ! that are located on the edge of a tile. ! \item[{[computationalLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the computational ! region with respect to the lower corner of the exclusive region. ! The default is a zero vector. ! \item[{[computationalUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the computational ! region with respect to the upper corner of the exclusive region. ! The default is a zero vector. ! \item[{[totalLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the total memory ! region with respect to the lower corner of the computational region. ! The default is to accommodate the union of exclusive and computational ! region. ! \item[{[totalUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the total memory ! region with respect to the upper corner of the computational region. ! The default is to accommodate the union of exclusive and computational ! region. ! \item[{[undistLBound]}] ! Lower bounds for the array dimensions that are not distributed. ! \item[{[undistUBound]}] ! Upper bounds for the array dimensions that are not distributed. ! \item[{[name]}] ! Name of the Array object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOP !------------------------------------------------------------------------------ integer :: localrc ! local return code integer :: rank, dimCount, i type(ESMF_ArraySpec) :: arrayspec ! Initialize return code; assume failure until success is certain localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Check init status of arguments used directly ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) ! Determine rank from input information call ESMF_DistGridGet(distgrid, dimCount=dimCount, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return rank = dimCount ! default value if (present(distgridToArrayMap)) then rank = 0 ! reset do i=1, size(distgridToArrayMap) if (distgridToArrayMap(i)>0) rank = rank+1 ! mapping against Array dim enddo endif if (present(undistLBound)) then rank = rank + size(undistLBound) else if (present(undistUBound)) then rank = rank + size(undistUBound) endif ! Set ArraySpec from rank and typekind input call ESMF_ArraySpecSet(arrayspec, rank=rank, typekind=typekind, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Call into ArraySpec overloaded interface ESMF_ArrayCreateAllocate = ESMF_ArrayCreate(distgrid, arrayspec, & indexflag=indexflag, distgridToArrayMap=distgridToArrayMap, & computationalEdgeLWidth=computationalEdgeLWidth, & computationalEdgeUWidth=computationalEdgeUWidth, & computationalLWidth=computationalLWidth, & computationalUWidth=computationalUWidth, & totalLWidth=totalLWidth, totalUWidth=totalUWidth, & undistLBound=undistLBound, undistUBound=undistUBound, & name=name, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end function ESMF_ArrayCreateAllocate !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCreateAllocateArb()" !BOP ! !IROUTINE: ESMF_ArrayCreate - Create Array object from typekind (allocate memory) w/ arbitrary seqIndices for halo ! !INTERFACE: ! Private name; call using ESMF_ArrayCreate() function ESMF_ArrayCreateAllocateArb(distgrid, typekind, & haloSeqIndexList, keywordEnforcer, distgridToArrayMap, & undistLBound, undistUBound, name, rc) ! ! !ARGUMENTS: type(ESMF_DistGrid), intent(in) :: distgrid type(ESMF_TypeKind_Flag), intent(in) :: typekind integer, intent(in) :: haloSeqIndexList(:) type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below integer, intent(in), optional :: distgridToArrayMap(:) integer, intent(in), optional :: undistLBound(:) integer, intent(in), optional :: undistUBound(:) character (len=*), intent(in), optional :: name integer, intent(out), optional :: rc ! ! !RETURN VALUE: type(ESMF_Array) :: ESMF_ArrayCreateAllocateArb ! ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Create an {\tt ESMF\_Array} object and allocate uninitialized data space ! according to typekind and distgrid. The Array rank is indirectly determined ! by the incoming information. Each PET must issue this call in unison in order ! to create a consistent Array object. DE-local allocations are made according ! to the total region defined by the {\tt distgrid} and {\tt haloSeqIndexList} ! arguments. ! ! The return value is the newly created {\tt ESMF\_Array} object. ! ! The arguments are: ! \begin{description} ! \item[distgrid] ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and ! distributed over DEs. The dimCount of distgrid must be smaller or equal ! to the rank specified in arrayspec, otherwise a runtime ESMF error will be ! raised. ! \item[typekind] ! The typekind of the Array. ! \item[haloSeqIndexList] ! One dimensional array containing sequence indices of local halo region. ! The size (and content) of {\tt haloSeqIndexList} can (and typically will) ! be different on each PET. ! \item[{[distgridToArrayMap]}] ! List that contains as many elements as is indicated by ! {\tt distgrids}|s dimCount. The list elements map each dimension of ! the DistGrid object to a dimension in the newly allocated Array object ! by specifying the appropriate Array dimension index. The default is ! to map all of {\tt distgrid}|s dimensions against the lower dimensions ! of the Array object in sequence, i.e. {\tt distgridToArrayMap = ! (/1, 2, .../)}. ! Unmapped dimensions in the Array object are not decomposed dimensions ! and form a tensor of rank = Array.rank - DistGrid.dimCount. ! All {\tt distgridToArrayMap} entries must be greater than or equal ! to zero and smaller than or equal to the Array rank. It is erroneous ! to specify the same entry multiple times unless it is zero. ! If the Array rank is less than the DistGrid dimCount then the default ! distgridToArrayMap will contain zeros for the dimCount - rank ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} ! indicates that the particular DistGrid dimension will be replicating ! the Array across the DEs along this direction. ! \item[{[undistLBound]}] ! Lower bounds for the array dimensions that are not distributed. ! \item[{[undistUBound]}] ! Upper bounds for the array dimensions that are not distributed. ! \item[{[name]}] ! Name of the Array object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOP !------------------------------------------------------------------------------ integer :: localrc ! local return code integer :: rank, dimCount, i type(ESMF_ArraySpec) :: arrayspec ! Initialize return code; assume failure until success is certain localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Check init status of arguments used directly ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) ! Determine rank from input information call ESMF_DistGridGet(distgrid, dimCount=dimCount, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return rank = dimCount ! default value if (present(distgridToArrayMap)) then rank = 0 ! reset do i=1, size(distgridToArrayMap) if (distgridToArrayMap(i)>0) rank = rank+1 ! mapping against Array dim enddo endif if (present(undistLBound)) then rank = rank + size(undistLBound) else if (present(undistUBound)) then rank = rank + size(undistUBound) endif ! Set ArraySpec from rank and typekind input call ESMF_ArraySpecSet(arrayspec, rank=rank, typekind=typekind, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Call into ArraySpec overloaded interface ESMF_ArrayCreateAllocateArb = ESMF_ArrayCreate(distgrid, arrayspec, & haloSeqIndexList, distgridToArrayMap=distgridToArrayMap, & undistLBound=undistLBound, undistUBound=undistUBound, & name=name, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end function ESMF_ArrayCreateAllocateArb !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCreateAllocateAS()" !BOP ! !IROUTINE: ESMF_ArrayCreate - Create Array object from ArraySpec (allocate memory) ! !INTERFACE: ! Private name; call using ESMF_ArrayCreate() function ESMF_ArrayCreateAllocateAS(distgrid, arrayspec, keywordEnforcer, indexflag, & distgridToArrayMap, computationalEdgeLWidth, computationalEdgeUWidth, & computationalLWidth, computationalUWidth, totalLWidth, totalUWidth, & undistLBound, undistUBound, name, rc) ! ! !ARGUMENTS: type(ESMF_DistGrid), intent(in) :: distgrid type(ESMF_ArraySpec), intent(in) :: arrayspec type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below type(ESMF_Index_Flag), intent(in), optional :: indexflag integer, intent(in), optional :: distgridToArrayMap(:) integer, intent(in), optional :: computationalEdgeLWidth(:) integer, intent(in), optional :: computationalEdgeUWidth(:) integer, intent(in), optional :: computationalLWidth(:) integer, intent(in), optional :: computationalUWidth(:) integer, intent(in), optional :: totalLWidth(:) integer, intent(in), optional :: totalUWidth(:) integer, intent(in), optional :: undistLBound(:) integer, intent(in), optional :: undistUBound(:) character (len=*), intent(in), optional :: name integer, intent(out), optional :: rc ! ! !RETURN VALUE: type(ESMF_Array) :: ESMF_ArrayCreateAllocateAS ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Create an {\tt ESMF\_Array} object and allocate uninitialized data space ! according to arrayspec and distgrid. Each PET must issue ! this call with identical arguments in order to create a consistent Array ! object. DE-local allocations are made according to the total region defined ! by the arguments to this call: {\tt distgrid} and the optional {\tt Width} ! arguments. ! ! The return value is the newly created {\tt ESMF\_Array} object. ! ! The arguments are: ! \begin{description} ! \item[distgrid] ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and ! distributed over DEs. The dimCount of distgrid must be smaller or equal ! to the rank specified in arrayspec, otherwise a runtime ESMF error will be ! raised. ! \item[arrayspec] ! {\tt ESMF\_ArraySpec} object containing the type/kind/rank information. ! \item[{[indexflag]}] ! Indicate how DE-local indices are defined. By default, the exclusive ! region of each DE is placed to start at the local index space origin, ! i.e. (1, 1, ..., 1). Alternatively the DE-local index space can be ! aligned with the global index space, if a global index space is well ! defined by the associated DistGrid. See section \ref{const:indexflag} ! for a list of valid indexflag options. ! \item[{[distgridToArrayMap]}] ! List that contains as many elements as is indicated by ! {\tt distgrids}|s dimCount. The list elements map each dimension of ! the DistGrid object to a dimension in the newly allocated Array object ! by specifying the appropriate Array dimension index. The default is ! to map all of {\tt distgrid}|s dimensions against the lower dimensions ! of the Array object in sequence, i.e. {\tt distgridToArrayMap = ! (/1, 2, .../)}. ! Unmapped dimensions in the Array object are not decomposed dimensions ! and form a tensor of rank = Array.rank - DistGrid.dimCount. ! All {\tt distgridToArrayMap} entries must be greater than or equal ! to zero and smaller than or equal to the Array rank. It is erroneous ! to specify the same entry multiple times unless it is zero. ! If the Array rank is less than the DistGrid dimCount then the default ! distgridToArrayMap will contain zeros for the dimCount - rank ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} ! indicates that the particular DistGrid dimension will be replicating ! the Array across the DEs along this direction. ! \item[{[computationalEdgeLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the computational ! region with respect to the lower corner of the exclusive region for DEs ! that are located on the edge of a tile. ! \item[{[computationalEdgeUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the computational ! region with respect to the upper corner of the exclusive region for DEs ! that are located on the edge of a tile. ! \item[{[computationalLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the computational ! region with respect to the lower corner of the exclusive region. ! The default is a zero vector. ! \item[{[computationalUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the computational ! region with respect to the upper corner of the exclusive region. ! The default is a zero vector. ! \item[{[totalLWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the lower corner of the total memory ! region with respect to the lower corner of the computational region. ! The default is to accommodate the union of exclusive and computational ! region. ! \item[{[totalUWidth]}] ! This vector argument must have dimCount elements, where dimCount is ! specified in distgrid. It specifies the upper corner of the total memory ! region with respect to the upper corner of the computational region. ! The default is to accommodate the union of exclusive and computational ! region. ! \item[{[undistLBound]}] ! Lower bounds for the array dimensions that are not distributed. ! \item[{[undistUBound]}] ! Upper bounds for the array dimensions that are not distributed. ! \item[{[name]}] ! Name of the Array object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOP !------------------------------------------------------------------------------ integer :: localrc ! local return code type(ESMF_Array) :: array ! opaque pointer to new C++ Array type(ESMF_InterfaceInt) :: distgridToArrayMapArg ! helper variable type(ESMF_InterfaceInt) :: computationalEdgeLWidthArg ! helper variable type(ESMF_InterfaceInt) :: computationalEdgeUWidthArg ! helper variable type(ESMF_InterfaceInt) :: computationalLWidthArg ! helper variable type(ESMF_InterfaceInt) :: computationalUWidthArg ! helper variable type(ESMF_InterfaceInt) :: totalLWidthArg ! helper variable type(ESMF_InterfaceInt) :: totalUWidthArg ! helper variable type(ESMF_InterfaceInt) :: undistLBoundArg ! helper variable type(ESMF_InterfaceInt) :: undistUBoundArg ! helper variable integer :: len_name ! helper variable ! Initialize return code; assume failure until success is certain localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Check init status of arguments ESMF_INIT_CHECK_SHALLOW(ESMF_ArraySpecGetInit, arrayspec, rc) ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) ! Check that indexflag does not specify mode that is not supported here if (present(indexflag)) then if (indexflag==ESMF_INDEX_USER) then call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & msg="- ESMF_INDEX_USER not supported as input value", & ESMF_CONTEXT, rcToReturn=rc) return endif endif ! Deal with (optional) array arguments distgridToArrayMapArg = ESMF_InterfaceIntCreate(distgridToArrayMap, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return computationalEdgeLWidthArg = & ESMF_InterfaceIntCreate(computationalEdgeLWidth, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return computationalEdgeUWidthArg = & ESMF_InterfaceIntCreate(computationalEdgeUWidth, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return computationalLWidthArg = ESMF_InterfaceIntCreate(computationalLWidth, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return computationalUWidthArg = ESMF_InterfaceIntCreate(computationalUWidth, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return totalLWidthArg = ESMF_InterfaceIntCreate(totalLWidth, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return totalUWidthArg = ESMF_InterfaceIntCreate(totalUWidth, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return undistLBoundArg = ESMF_InterfaceIntCreate(undistLBound, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return undistUBoundArg = ESMF_InterfaceIntCreate(undistUBound, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Mark this Array object as invalid array%this = ESMF_NULL_POINTER ! Call into the C++ interface, which will sort out optional arguments ! Optional name argument requires separate calls into C++ if (present(name)) then len_name = len(name) call c_ESMC_ArrayCreateAllocate(array, arrayspec, distgrid, & distgridToArrayMapArg, & computationalEdgeLWidthArg, computationalEdgeUWidthArg, & computationalLWidthArg, computationalUWidthArg, totalLWidthArg, & totalUWidthArg, indexflag, & undistLBoundArg, undistUBoundArg, name, len_name, localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return else len_name = 0 call c_ESMC_ArrayCreateAllocate(array, arrayspec, distgrid, & distgridToArrayMapArg, & computationalEdgeLWidthArg, computationalEdgeUWidthArg, & computationalLWidthArg, computationalUWidthArg, totalLWidthArg, & totalUWidthArg, indexflag, & undistLBoundArg, undistUBoundArg, "", len_name, localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return endif ! Garbage collection call ESMF_InterfaceIntDestroy(distgridToArrayMapArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(computationalEdgeLWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(computationalEdgeUWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(computationalLWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(computationalUWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(totalLWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(totalUWidthArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(undistLBoundArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return call ESMF_InterfaceIntDestroy(undistUBoundArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Set return value ESMF_ArrayCreateAllocateAS = array ! Set init code ESMF_INIT_SET_CREATED(ESMF_ArrayCreateAllocateAS) ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end function ESMF_ArrayCreateAllocateAS !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCreateAllocateASArb()" !BOP ! !IROUTINE: ESMF_ArrayCreate - Create Array object from ArraySpec (allocate memory) w/ arbitrary seqIndices for halo ! !INTERFACE: ! Private name; call using ESMF_ArrayCreate() function ESMF_ArrayCreateAllocateASArb(distgrid, arrayspec, & haloSeqIndexList, keywordEnforcer, distgridToArrayMap, & undistLBound, undistUBound, name, rc) ! ! !ARGUMENTS: type(ESMF_DistGrid), intent(in) :: distgrid type(ESMF_ArraySpec), intent(in) :: arrayspec integer, intent(in) :: haloSeqIndexList(:) type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below integer, intent(in), optional :: distgridToArrayMap(:) integer, intent(in), optional :: undistLBound(:) integer, intent(in), optional :: undistUBound(:) character (len=*), intent(in), optional :: name integer, intent(out), optional :: rc ! ! !RETURN VALUE: type(ESMF_Array) :: ESMF_ArrayCreateAllocateASArb ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Create an {\tt ESMF\_Array} object and allocate uninitialized data space ! according to arrayspec and distgrid. Each PET must issue this call in unison ! in order to create a consistent Array object. DE-local allocations are made ! according to the total region defined by the arguments to this call: ! {\tt distgrid} and {\tt haloSeqIndexList} arguments. ! ! The return value is the newly created {\tt ESMF\_Array} object. ! ! The arguments are: ! \begin{description} ! \item[distgrid] ! {\tt ESMF\_DistGrid} object that describes how the array is decomposed and ! distributed over DEs. The dimCount of distgrid must be smaller or equal ! to the rank specified in arrayspec, otherwise a runtime ESMF error will be ! raised. ! \item[arrayspec] ! {\tt ESMF\_ArraySpec} object containing the type/kind/rank information. ! \item[haloSeqIndexList] ! One dimensional array containing sequence indices of local halo region. ! The size (and content) of {\tt haloSeqIndexList} can (and typically will) ! be different on each PET. ! \item[{[distgridToArrayMap]}] ! List that contains as many elements as is indicated by ! {\tt distgrids}|s dimCount. The list elements map each dimension of ! the DistGrid object to a dimension in the newly allocated Array object ! by specifying the appropriate Array dimension index. The default is ! to map all of {\tt distgrid}|s dimensions against the lower dimensions ! of the Array object in sequence, i.e. {\tt distgridToArrayMap = ! (/1, 2, .../)}. ! Unmapped dimensions in the Array object are not decomposed dimensions ! and form a tensor of rank = Array.rank - DistGrid.dimCount. ! All {\tt distgridToArrayMap} entries must be greater than or equal ! to zero and smaller than or equal to the Array rank. It is erroneous ! to specify the same entry multiple times unless it is zero. ! If the Array rank is less than the DistGrid dimCount then the default ! distgridToArrayMap will contain zeros for the dimCount - rank ! rightmost entries. A zero entry in the {\tt distgridToArrayMap} ! indicates that the particular DistGrid dimension will be replicating ! the Array across the DEs along this direction. ! \item[{[undistLBound]}] ! Lower bounds for the array dimensions that are not distributed. ! \item[{[undistUBound]}] ! Upper bounds for the array dimensions that are not distributed. ! \item[{[name]}] ! Name of the Array object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOP !------------------------------------------------------------------------------ integer :: localrc ! local return code integer :: localHaloCount integer :: dimCount, localDeCount type(ESMF_DELayout) :: delayout ! Initialize return code; assume failure until success is certain localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Check init status of arguments used directly ESMF_INIT_CHECK_DEEP_SHORT(ESMF_DistGridGetInit, distgrid, rc) ! Ensure incoming DistGrid is 1D call ESMF_DistGridGet(distgrid, dimCount=dimCount, delayout=delayout, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return if (dimCount /= 1) then call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & msg="- this ArrayCreate() interface requires 1D DistGrid", & ESMF_CONTEXT, rcToReturn=rc) return endif ! Ensure there is only 1 DE/PET call ESMF_DELayoutGet(delayout, localDeCount=localDeCount, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return if (localDeCount /= 1) then call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_VALUE, & msg="- this ArrayCreate() interface only supports 1DE/PETs", & ESMF_CONTEXT, rcToReturn=rc) return endif ! Prepare for actual ArrayCreate() call localHaloCount = size(haloSeqIndexList) ! Call into conventional interface to create the Array object ESMF_ArrayCreateAllocateASArb = ESMF_ArrayCreate(distgrid, arrayspec, & distgridToArrayMap=distgridToArrayMap, totalUWidth=(/localHaloCount/), & undistLBound=undistLBound, undistUBound=undistUBound, name=name, & rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Set the arbitrary sequence indices in the Array halo rim call ESMF_ArraySetPLocalDeInternal(ESMF_ArrayCreateAllocateASArb, & localDe=0, rimSeqIndex=haloSeqIndexList, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end function ESMF_ArrayCreateAllocateASArb !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCreateCopy()" !BOP ! !IROUTINE: ESMF_ArrayCreate - Create Array object as copy of existing Array object ! !INTERFACE: ! Private name; call using ESMF_ArrayCreate() function ESMF_ArrayCreateCopy(array, keywordEnforcer, rc) ! ! !ARGUMENTS: type(ESMF_Array), intent(in) :: array type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below integer, intent(out), optional :: rc ! ! !RETURN VALUE: type(ESMF_Array) :: ESMF_ArrayCreateCopy ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Create an {\tt ESMF\_Array} object as the copy of an existing Array. ! ! The return value is the newly created {\tt ESMF\_Array} object. ! ! The arguments are: ! \begin{description} ! \item[array] ! {\tt ESMF\_Array} object to be copied. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOP !------------------------------------------------------------------------------ integer :: localrc ! local return code type(ESMF_Array) :: arrayOut ! opaque pointer to new C++ Array ! Initialize return code; assume failure until success is certain localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Mark this Array object as invalid arrayOut%this = ESMF_NULL_POINTER ! Call into the C++ interface call c_ESMC_ArrayCreateCopy(array, arrayOut, localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Set return value ESMF_ArrayCreateCopy = arrayOut ! Set init code ESMF_INIT_SET_CREATED(ESMF_ArrayCreateCopy) ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end function ESMF_ArrayCreateCopy !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayDestroy()" !BOP ! !IROUTINE: ESMF_ArrayDestroy - Release resources associated with an Array object ! !INTERFACE: subroutine ESMF_ArrayDestroy(array, keywordEnforcer, rc) ! ! !ARGUMENTS: type(ESMF_Array), intent(inout) :: array type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below integer, intent(out), optional :: rc ! ! !STATUS: ! \apiStatusCompatible ! ! !DESCRIPTION: ! Destroy an {\tt ESMF\_Array}, releasing the resources associated with ! the object. ! ! The arguments are: ! \begin{description} ! \item[array] ! {\tt ESMF\_Array} object to be destroyed. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOP !------------------------------------------------------------------------------ integer :: localrc ! local return code ! Initialize return code; assume failure until success is certain localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Check init status of arguments ESMF_INIT_CHECK_DEEP_SHORT(ESMF_ArrayGetInit, array, rc) ! Call into the C++ interface, which will sort out optional arguments. call c_ESMC_ArrayDestroy(array, localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Mark this Array object as invalid array%this = ESMF_NULL_POINTER ! Set init code ESMF_INIT_SET_DELETED(array) ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end subroutine ESMF_ArrayDestroy !------------------------------------------------------------------------------ ! -------------------------- ESMF-internal method ----------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayGetInit" !BOPI ! !IROUTINE: ESMF_ArrayGetInit - Internal access routine for init code ! ! !INTERFACE: function ESMF_ArrayGetInit(array) ! ! !RETURN VALUE: ESMF_INIT_TYPE :: ESMF_ArrayGetInit ! ! !ARGUMENTS: type(ESMF_Array), intent(in), optional :: array ! ! !DESCRIPTION: ! Access deep object init code. ! ! The arguments are: ! \begin{description} ! \item [array] ! Array object. ! \end{description} ! !EOPI if (present(array)) then ESMF_ArrayGetInit = ESMF_INIT_GET(array) else ESMF_ArrayGetInit = ESMF_INIT_CREATED endif end function ESMF_ArrayGetInit !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArraySetInitCreated()" !BOPI ! !IROUTINE: ESMF_ArraySetInitCreated - Set Array init code to "CREATED" ! !INTERFACE: subroutine ESMF_ArraySetInitCreated(array, rc) ! ! !ARGUMENTS: type(ESMF_Array), intent(inout) :: array integer, intent(out), optional :: rc ! ! ! !DESCRIPTION: ! Set init code in Array object to "CREATED". ! ! The arguments are: ! \begin{description} ! \item[array] ! Specified {\tt ESMF\_Array} object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOPI !------------------------------------------------------------------------------ integer :: localrc ! local return code ! Assume failure until success if (present(rc)) rc = ESMF_RC_NOT_IMPL localrc = ESMF_RC_NOT_IMPL ! Set init code ESMF_INIT_SET_CREATED(array) ! Return success if (present(rc)) rc = ESMF_SUCCESS end subroutine ESMF_ArraySetInitCreated !------------------------------------------------------------------------------ ! -------------------------- ESMF-internal method ----------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayGetThis()" !BOPI ! !IROUTINE: ESMF_ArrayGetThis - Internal access routine for C++ pointer ! !INTERFACE: subroutine ESMF_ArrayGetThis(array, this, rc) ! ! !ARGUMENTS: type(ESMF_Array), intent(in) :: array type(ESMF_Pointer), intent(out) :: this integer, intent(out), optional :: rc ! ! ! !DESCRIPTION: ! Internal access routine for C++ pointer. ! ! The arguments are: ! \begin{description} ! \item[array] ! Specified {\tt ESMF\_Array} object. ! \item[this] ! C++ pointer. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOPI !------------------------------------------------------------------------------ ! initialize return code; assume routine not implemented if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Copy C++ pointer this = array%this ! return successfully if (present(rc)) rc = ESMF_SUCCESS end subroutine ESMF_ArrayGetThis !------------------------------------------------------------------------------ ! -------------------------- ESMF-internal method ----------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArraySetThis()" !BOPI ! !IROUTINE: ESMF_ArraySetThis - Set C++ pointer in VM ! !INTERFACE: subroutine ESMF_ArraySetThis(array, this, rc) ! ! !ARGUMENTS: type(ESMF_Array), intent(inout) :: array type(ESMF_Pointer), intent(in) :: this integer, intent(out), optional :: rc ! ! ! !DESCRIPTION: ! Set C++ pointer in VM. ! ! The arguments are: ! \begin{description} ! \item[array] ! Specified {\tt ESMF\_Array} object. ! \item[this] ! C++ pointer. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOPI !------------------------------------------------------------------------------ ! initialize return code; assume routine not implemented if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Copy C++ pointer array%this = this ! return successfully if (present(rc)) rc = ESMF_SUCCESS end subroutine ESMF_ArraySetThis !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArraySetThisNull()" !BOPI ! !IROUTINE: ESMF_ArraySetThisNull - Set Array this member to ESMF_NULL_POINTER ! !INTERFACE: subroutine ESMF_ArraySetThisNull(array, rc) ! ! !ARGUMENTS: type(ESMF_Array), intent(inout) :: array integer, intent(out), optional :: rc ! ! ! !DESCRIPTION: ! Set Array this member to ESMF_NULL_POINTER. ! ! The arguments are: ! \begin{description} ! \item[array] ! Specified {\tt ESMF\_Array} object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOPI !------------------------------------------------------------------------------ integer :: localrc ! local return code ! Assume failure until success if (present(rc)) rc = ESMF_RC_NOT_IMPL localrc = ESMF_RC_NOT_IMPL ! Set this member to NULL array%this = ESMF_NULL_POINTER ! Return success if (present(rc)) rc = ESMF_SUCCESS end subroutine ESMF_ArraySetThisNull !------------------------------------------------------------------------------ ! -------------------------- ESMF-internal method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArraySetPLocalDeInternal()" !BOPI ! !IROUTINE: ESMF_ArraySetPLocalDeInternal - Set Array properties ! ! !INTERFACE: subroutine ESMF_ArraySetPLocalDeInternal(array, localDe, rimSeqIndex, rc) ! ! !ARGUMENTS: type(ESMF_Array), intent(inout) :: array integer, intent(in) :: localDe integer, intent(in), optional :: rimSeqIndex(:) integer, intent(out), optional :: rc ! ! !DESCRIPTION: ! Sets adjustable settings in an {\tt ESMF\_Array} object for a specific ! localDe. ! ! The arguments are: ! \begin{description} ! \item [array] ! {\tt ESMF\_Array} object for which to set properties. ! \item [localDe] ! Local DE for which to set values. ! \item[{[rimSeqIndex]}] ! Sequence indices in the halo rim of localDe. ! \item [{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOPI !------------------------------------------------------------------------------ integer :: localrc ! local return code type(ESMF_InterfaceInt) :: rimSeqIndexArg ! helper variable ! initialize return code; assume routine not implemented localrc = ESMF_RC_NOT_IMPL if (present(rc)) rc = ESMF_RC_NOT_IMPL ! Check init status of arguments ESMF_INIT_CHECK_DEEP(ESMF_ArrayGetInit, array, rc) ! Deal with (optional) array arguments rimSeqIndexArg = ESMF_InterfaceIntCreate(farray1D=rimSeqIndex, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Call into the C++ interface, which will sort out optional arguments call c_esmc_arraysetplocalde(array, localDe, rimSeqIndexArg, localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Garbage collection call ESMF_InterfaceIntDestroy(rimSeqIndexArg, rc=localrc) if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! return successfully if (present(rc)) rc = ESMF_SUCCESS end subroutine ESMF_ArraySetPLocalDeInternal !------------------------------------------------------------------------------ ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCopyThis()" !BOPI ! !IROUTINE: ESMF_ArrayCopyThis - Copy Array this member ! !INTERFACE: subroutine ESMF_ArrayCopyThis(arrayIn, arrayOut, rc) ! ! !ARGUMENTS: type(ESMF_Array), intent(in) :: arrayIn type(ESMF_Array), intent(inout) :: arrayOut integer, intent(out), optional :: rc ! ! ! !DESCRIPTION: ! Copy Array this member. Do not set init code. ! ! The arguments are: ! \begin{description} ! \item[arrayIn] ! Input {\tt ESMF\_Array} object. ! \item[arrayOut] ! Output {\tt ESMF\_Array} object. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOPI !------------------------------------------------------------------------------ integer :: localrc ! local return code ! Assume failure until success if (present(rc)) rc = ESMF_RC_NOT_IMPL localrc = ESMF_RC_NOT_IMPL ! Copy this member arrayOut%this = arrayIn%this ! Return success if (present(rc)) rc = ESMF_SUCCESS end subroutine ESMF_ArrayCopyThis !------------------------------------------------------------------------------ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!! old-style newArray calls of 1st prototype calls !!!!!!!!!!!!!!!!! ! -------------------------- ESMF-public method ------------------------------- ^undef ESMF_METHOD ^define ESMF_METHOD "ESMF_ArrayCreateFromLA()" !BOPI ! !IROUTINE: ESMF_ArrayCreate - Create from LocalArray ! !INTERFACE: ! Private name; call using ESMF_ArrayCreate() function ESMF_ArrayCreateFromLA(localarray, haloWidth, deCount, rootPET, rc) ! ! !ARGUMENTS: type(ESMF_LocalArray), intent(in), optional :: localarray integer, target, intent(in), optional :: haloWidth(:) integer, intent(in), optional :: deCount integer, intent(in) :: rootPET integer, intent(out), optional :: rc ! ! !RETURN VALUE: type(ESMF_Array) :: ESMF_ArrayCreateFromLA ! ! !DESCRIPTION: ! Create a new {\tt ESMF\_Array} from an {\tt ESMF\_LocalArray} on ! {\tt rootPET}. ! ! The arguments are: ! \begin{description} ! \item[{[localarray]}] ! The {\tt ESMF\_LocalArray} object that will be used to determine ! the index space of the {\tt ESMF\_Array}. The data in ! {\tt rootPET}|s {\tt localarray} will be scattered across the newly ! created {\tt ESMF\_Array} object before returning. Only ! {\tt rootPET} must provide a (valid) {\tt localarray} argument. ! \item[{[haloWidth]}] ! Array holding the halo width for each array dimension. Default is ! a halo width of 0 in all dimensions. A negative halo width indicates ! that the respective array dimension is not to be decomposed. Only ! {\tt rootPET} must provide a (valid) {\tt haloWidth} argument. ! \item[{[deCount]}] ! Number of DEs into which this Array is to be decomposed. By default ! deCount will equal petCount of the current VM context. Only ! {\tt rootPET} must provide a (valid) {\tt deCount} argument. ! \item[rootPET] ! PET that holds the valid data in {\tt localarray}. ! \item[{[rc]}] ! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors. ! \end{description} ! !EOPI !------------------------------------------------------------------------------ integer :: localrc ! local return code integer :: laRank integer, pointer :: opt_haloWidth(:) integer :: len_haloWidth integer, target :: dummy(0) ! used to satisfy the C interface... type(ESMF_VM):: vm integer:: localPET type(ESMF_Array):: array ! opaque pointer to new C++ Array object ! Assume failure until success if (present(rc)) rc = ESMF_RC_NOT_IMPL localrc = ESMF_RC_NOT_IMPL ! Check init status of arguments ESMF_INIT_CHECK_DEEP_SHORT(ESMF_LocalArrayGetInit, localarray, rc) ! Initialize the pointer to NULL array%this = ESMF_NULL_POINTER if (present(haloWidth)) then opt_haloWidth => haloWidth len_haloWidth = size(haloWidth) ! For rootPET check if size(haloWidth) matches the rank of the LocalArray call ESMF_VMGetCurrent(vm, rc=localrc) call ESMF_VMGet(vm, localPET=localPET, rc=localrc) if (localPET==rootPET) then call ESMF_LocalArrayGet(localarray, rank=laRank, rc=localrc) if (len_haloWidth/=laRank) then if (ESMF_LogFoundError(ESMF_RC_ARG_VALUE, & msg="size of haloWidth array does not match rank of LocalArray", & ESMF_CONTEXT, rcToReturn=rc)) return endif endif else opt_haloWidth => dummy len_haloWidth = 0 endif !DUMMY TEST TO QUIET DOWN COMPILER WARNINGS !TODO: Remove the following dummy test when dummy argument actually used if (size(opt_haloWidth) == size(opt_haloWidth)) continue ! Call into the C++ interface, which will sort out optional arguments. ! call c_ESMC_ArrayCreate(array, localarray, opt_haloWidth, len_haloWidth, & ! deCount, rootPET, localrc) ! Use LogErr to handle return code if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, & ESMF_CONTEXT, rcToReturn=rc)) return ! Set return value ESMF_ArrayCreateFromLA = array ! Set init code ESMF_INIT_SET_CREATED(ESMF_ArrayCreateFromLA) ! Return successfully if (present(rc)) rc = ESMF_SUCCESS end function ESMF_ArrayCreateFromLA !------------------------------------------------------------------------------ !------------------------------------------------------------------------------ end module ESMF_ArrayCreateMod