subroutine MAPL_CFIOWriteBundleWait( MCFIO, CLOCK, RC )
type(MAPL_CFIO ), intent(INOUT) :: MCFIO
type(ESMF_CLOCK), intent(INOUT) :: CLOCK
integer, optional, intent( OUT) :: RC
! Locals
!-------
integer :: status
integer :: L, K, NN
logical :: AmRoot, MyGlobal
real, pointer :: Gptr2Out(:,: )
real, pointer :: PtrGlob (:,: )
integer :: counts(5)
integer :: IM0,JM0,I,IP
logical :: FixPole
integer :: lm, nv
logical :: transAlreadyDone
type(Ptr2Arr), allocatable :: globPtrArr(:)
type(Ptr2Arr) :: PtrTypeIn(2)
type(Ptr2Arr) :: PtrTypeOut(2)
_UNUSED_DUMMY(CLOCK)
! Space for global arrays is allocated everywhere, even if not used.
!------------------------------------------------------------------
_ASSERT(MCFIO%CREATED, 'MCFIO%CREATED is false')
! Allocate global 2d and 3d arrays at the writing resolution
! Note that everybody allocated these.
!-----------------------------------------------------------
call MAPL_GridGet( MCFIO%GRID, globalCellCountPerDim=COUNTS, RC=STATUS)
_VERIFY(STATUS)
IM0 = COUNTS(1)
JM0 = COUNTS(2)
!if(any(mCFIO%myPE==mCFIO%Krank)) then
!allocate(Gptr3Out(Mcfio%IM, Mcfio%JM,1), stat=STATUS)
!_VERIFY(STATUS)
!Gptr2Out => Gptr3Out(:,:,1)
!Gptr2Out(:,:) = 0.0
!end if
nn = 0
AmRoot = mCFIO%myPE==MCFIO%RootRank
allocate(globPtrArr(size(mCFIO%reqs)), stat=status)
_VERIFY(STATUS)
COLCTVWAIT: do nn=1,size(mCFIO%reqs)
! Wait on request for slice nn
!-----------------------------
call MAPL_CollectiveWait(MCFIO%reqs(nn), DstArray=PtrGlob, rc=status)
_VERIFY(STATUS)
globPtrArr(nn)%ptr => PtrGlob ! this is valid only if myGlobal is .true.
end do COLCTVWAIT
nn = 0
VARIABLES: do L=1,size(MCFIO%VarDims)
FixPole = (MCFIO%VarType(L) == MAPL_VectorField) .and. &
(JM0 == 6*IM0) .and. &
(Mcfio%JM /= 6*mcfio%IM)
RANK: if (MCFIO%VarDims(L)==2) then
LM = 1
else if (MCFIO%VarDims(L)==3) then
LM = MCFIO%lm
else
LM = 0
end if RANK
LEVELS: do k=1,LM
nn = nn + 1
MyGlobal = mCFIO%Krank(nn) == MCFIO%MYPE
PtrGlob => globPtrArr(nn)%ptr
! Horizontal Interpolation and Shaving on PEs with global data
! ------------------------------------------------------------
if( MyGlobal ) then
nv = mCFIO%pairList(nn)
VECTORTEST: if (nv == 0) then
! scalar
allocate( MCFIO%reqs(nn)%Trans_Array(Mcfio%IM, Mcfio%JM, 1), stat=STATUS )
_VERIFY(STATUS)
Gptr2Out => MCFIO%reqs(nn)%Trans_Array(:,:,1)
PtrTypeIn (1)%ptr => globPtrArr(nn)%ptr
PtrTypeOut(1)%ptr => Gptr2Out
call TransShaveAndSend(PtrTypeIn(1:1),PtrTypeOut(1:1),MCFIO%reqs(nn)%s_rqst,doTrans=.true.,IdxOut=1)
_VERIFY(status)
else if (nv > 0) then
! I am U part of vector
if (associated(MCFIO%reqs(nn)%Trans_Array)) then
_ASSERT(associated(MCFIO%reqs(nv)%Trans_Array), 'Trans_Array not associated')
TransAlreadyDone = .true.
else
TransAlreadyDone = .false.
allocate( MCFIO%reqs(nn)%Trans_Array(Mcfio%IM, Mcfio%JM, 1), stat=STATUS )
_VERIFY(STATUS)
allocate( MCFIO%reqs(nv)%Trans_Array(Mcfio%IM, Mcfio%JM, 1), stat=STATUS )
_VERIFY(STATUS)
endif
PtrTypeIn (1)%ptr => globPtrArr(nn)%ptr
PtrTypeIn (2)%ptr => globPtrArr(nv)%ptr
PtrTypeOut(1)%ptr => MCFIO%reqs(nn)%Trans_Array(:,:,1)
PtrTypeOut(2)%ptr => MCFIO%reqs(nv)%Trans_Array(:,:,1)
call TransShaveAndSend(PtrTypeIn(1:2),PtrTypeOut(1:2),MCFIO%reqs(nn)%s_rqst,doTrans=.not.TransAlreadyDone,IdxOut=1)
_VERIFY(status)
else
! I am V part of vector
nv = abs(nv)
if (associated(MCFIO%reqs(nn)%Trans_Array)) then
_ASSERT(associated(MCFIO%reqs(nv)%Trans_Array), 'Trans_Array not associated')
TransAlreadyDone = .true.
else
TransAlreadyDone = .false.
allocate( MCFIO%reqs(nn)%Trans_Array(Mcfio%IM, Mcfio%JM, 1), stat=STATUS )
_VERIFY(STATUS)
allocate( MCFIO%reqs(nv)%Trans_Array(Mcfio%IM, Mcfio%JM, 1), stat=STATUS )
_VERIFY(STATUS)
endif
PtrTypeIn (1)%ptr => globPtrArr(nv)%ptr
PtrTypeIn (2)%ptr => globPtrArr(nn)%ptr
PtrTypeOut(1)%ptr => MCFIO%reqs(nv)%Trans_Array(:,:,1)
PtrTypeOut(2)%ptr => MCFIO%reqs(nn)%Trans_Array(:,:,1)
call TransShaveAndSend(PtrTypeIn(1:2),PtrTypeOut(1:2),MCFIO%reqs(nn)%s_rqst,doTrans=.not.TransAlreadyDone,IdxOut=2)
_VERIFY(status)
end if VECTORTEST
endif
end do LEVELS
end do VARIABLES
! do nn=1,size(mCFIO%reqs)
! MyGlobal = MCFIO%Krank(nn) == MCFIO%MYPE
! if (myGlobal) then
! deallocate(globPtrArr(nn)%ptr)
! NULLIFY(globPtrArr(nn)%ptr)
! end if
! end do
deallocate(globPtrArr)
!if(AmRoot) then
! write(6,'(1X,"TransShaveAndSend: ",i6," Slices (",i3," Nodes, ",i2," CoresPerNode) to File: ",a)') &
! size(MCFIO%reqs),mCFIO%partsize/mCFIO%numcores,mCFIO%numcores,trim(mCFIO%fName)
!endif
!if (any(mCFIO%myPE==mCFIO%Krank)) then
!deallocate(Gptr3Out, stat=STATUS)
!_VERIFY(STATUS)
!end if
_RETURN(ESMF_SUCCESS)
contains
subroutine TransShaveAndSend(PtrIn,PtrOut,request,doTrans,idxOut)
type(Ptr2Arr) :: PtrIn(:)
type(Ptr2Arr) :: PtrOut(:)
integer :: request
logical :: doTrans
integer :: idxOut
real, pointer :: Gin (:,:)
real, pointer :: Gout(:,:)
real, dimension(:,:,:), pointer :: uin, uout, vin, vout
integer :: im, jm
type(c_ptr) :: cptr
if (size(PtrIn) == 1) then
_ASSERT(idxOut ==1, 'idxOut /= 1')
Gin => PtrIn(1)%ptr
Gout => PtrOut(1)%ptr
if (associated(mCFIO%regridder)) then
if (mCFIO%regridConservative) then
call mCFIO%regridder%regrid(Gin, Gout, rc=status)
_VERIFY(STATUS)
else
call mCFIO%regridder%set_undef_value(MAPL_undef)
call mCFIO%regridder%regrid(Gin, Gout, rc=status)
_VERIFY(status)
end if
if (mcfio%vinterp .and. (lm .ne. 1) ) then
if (mcfio%ascending) then
where(mcfio%surfaceLayer<mcfio%unmodifiedLevs(k)) gout=MAPL_UNDEF
else
where(mcfio%surfaceLayer>mcfio%unmodifiedLevs(k)) gout=MAPL_UNDEF
endif
end if
else
_ASSERT( all(shape(gout)==shape(gin)), 'in-out shape mismatch')
gout=gin
end if
! if going from CS to LAT-LON pole winds are wrong, approximate fix below
if (FixPole) then
do i=1,mcfio%im
ip = i+(mcfio%im/2)
if (ip > mcfio%im) ip = ip - mcfio%im
if ( (gout(i,mcfio%jm-1) == MAPL_UNDEF) .or. (gout(ip,mcfio%jm-1) == MAPL_UNDEF)) then
gout(i,mcfio%jm) = MAPL_UNDEF
else
gout(i,mcfio%jm)=(gout(i,mcfio%jm-1)-gout(ip,mcfio%jm-1))/2.0
end if
if ( (gout(i,2) == MAPL_UNDEF) .or. (gout(ip,2) == MAPL_UNDEF)) then
gout(i,1) = MAPL_UNDEF
else
gout(i,1)=(gout(i,2)-gout(ip,2))/2.0
endif
enddo
endif
deallocate(Gin)
nullify (Gin)
else
_ASSERT(size(PtrIn) == 2, 'if not scalar, ptrIn must be 2-vector')
_ASSERT(size(PtrOut) == 2, 'if not scalar, ptrOut must be 2-vector')
Gout => PtrOut(idxOut)%ptr
! TLC: Probably do not need this conditional now that there are identity regridders
if (doTrans) then
_ASSERT(associated(mcfio%regridder), 'mcfio%regridder not associated')
im = size(PtrIn(1)%ptr,1)
jm = size(PtrIn(1)%ptr,2)
! MAT PGI cannot handle C_LOC call inside C_F_POINTER
cptr = C_loc(PtrIn(1)%ptr(1,1))
call C_F_pointer (cptr, uin,[im,jm,1])
cptr = C_loc(PtrIn(2)%ptr(1,1))
call C_F_pointer (cptr, vin,[im,jm,1])
!@# allocate(uin(im,jm,1), vin(im,jm,1))
!@# uin(:,:,1) = PtrIn(1)%ptr
!@# vin(:,:,1) = PtrIn(2)%ptr
im = size(PtrOut(1)%ptr,1)
jm = size(PtrOut(1)%ptr,2)
cptr = C_loc(PtrOut(1)%ptr(1,1))
call C_F_pointer (cptr, uout,[im,jm,1])
cptr = C_loc(PtrOut(2)%ptr(1,1))
call C_F_pointer (cptr, vout,[im,jm,1])
!@# allocate(uout(im,jm,1), vout(im,jm,1))
call mCFIO%regridder%set_undef_value(MAPL_undef)
call mCFIO%regridder%regrid(uin, vin, uout, vout, rc=status)
_VERIFY(status)
if (mcfio%vinterp .and. (lm .ne. 1) ) then
if (mcfio%ascending) then
where(mcfio%surfaceLayer<mcfio%unmodifiedLevs(k)) uout(:,:,1)=MAPL_UNDEF
where(mcfio%surfaceLayer<mcfio%unmodifiedLevs(k)) vout(:,:,1)=MAPL_UNDEF
else
where(mcfio%surfaceLayer>mcfio%unmodifiedLevs(k)) uout(:,:,1)=MAPL_UNDEF
where(mcfio%surfaceLayer>mcfio%unmodifiedLevs(k)) vout(:,:,1)=MAPL_UNDEF
endif
end if
deallocate(PtrIn(1)%ptr)
nullify(PtrIn(1)%ptr)
deallocate(PtrIn(2)%ptr)
nullify(PtrIn(2)%ptr)
end if
end if
if(mCFIO%NBITS < 24) then
call ESMF_CFIODownBit ( Gout, Gout, mCFIO%NBITS, undef=MAPL_undef, rc=STATUS )
_VERIFY(STATUS)
end if
call MPI_ISend(Gout, size(Gout), MPI_REAL, MCFIO%RootRank, &
trans_tag, mCFIO%comm, request, STATUS)
_VERIFY(STATUS)
return
end subroutine TransShaveAndSend
end subroutine MAPL_CFIOWriteBundleWait