ESMF_SparseMatca Subroutine

private subroutine ESMF_SparseMatca(sparseMats, sparseMatd, ngrid, tag, rc)

Arguments

Type IntentOptional Attributes Name
type(ESMF_XGridSpec), intent(in) :: sparseMats(:)
type(ESMF_XGridSpec), pointer :: sparseMatd(:)
integer, intent(in) :: ngrid
character(len=*), intent(in) :: tag
integer, intent(out), optional :: rc

Source Code

subroutine ESMF_SparseMatca(sparseMats, sparseMatd, ngrid, tag, rc)

!
! !ARGUMENTS:
    type(ESMF_XGridSpec), intent(in)    :: sparseMats(:)
    type(ESMF_XGridSpec),     pointer   :: sparseMatd(:)
    integer, intent(in)                 :: ngrid
    character(len=*), intent(in)        :: tag
    integer, intent(out), optional      :: rc

!
! !DESCRIPTION:
!      Allocate internal SMM parameters and copy from src.
!
!     The arguments are:
!     \begin{description}
!     \item [sparseMats]
!           the source {\tt ESMF\_XGridSpec} object.
!     \item [sparseMatd]
!           the destination {\tt ESMF\_XGridSpec} object.
!     \item [ngrid]
!           number of grid, redundency check.
!     \item [tag]
!           A string to indicate which one of the 4 SMM parameters is used.
!     \item [{[rc]}]
!           Return code; equals {\tt ESMF\_SUCCESS} only if successful.
!     \end{description}
!
!EOPI
    integer                             :: i, localrc

    ! Initialize
    localrc = ESMF_RC_NOT_IMPL

    ! Initialize return code   
    if(present(rc)) rc = ESMF_RC_NOT_IMPL

    if(size(sparseMats,1) /= ngrid) then
        call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, & 
           msg="- number of Grids different from size of sparseMat for "//tag, &
           ESMF_CONTEXT, rcToReturn=rc) 
        return
    endif

    do i = 1, ngrid
        if(.not. associated(sparseMats(i)%factorIndexList) .or. &
           .not. associated(sparseMats(i)%factorList)) then
            call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, & 
               msg="- sparseMat not initiailzed properly for "//tag, &
               ESMF_CONTEXT, rcToReturn=rc) 
            return
        endif

        if(size(sparseMats(i)%factorIndexList, 2) /= size(sparseMats(i)%factorList, 1)) then
            call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, & 
               msg="- sparseMat factorIndexList and factorList sizes not consistent "//tag, &
               ESMF_CONTEXT, rcToReturn=rc) 
            return
        endif
    enddo
        
    allocate(sparseMatd(ngrid), stat=localrc)
    if (ESMF_LogFoundAllocError(localrc, &
        msg="Allocating xgtype%"//tag, &
        ESMF_CONTEXT, rcToReturn=rc)) return
    do i = 1, ngrid
      sparseMatd(i) = sparseMats(i)
    enddo

    if(present(rc)) rc = ESMF_SUCCESS

end subroutine ESMF_SparseMatca