!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| program POP_HaloTest !---------------------------------------------------------------------- ! ! this program tests POP global reduction operations, including ! global sums, minval, maxval, etc. This test code relies on the ! redistribution mod working correctly for gather/scatter operations. ! !---------------------------------------------------------------------- use POP_KindsMod use POP_ErrorMod use POP_CommMod use POP_DomainSizeMod use POP_BlocksMod use POP_DistributionMod use POP_RedistributeMod use POP_ReductionsMod use POP_FieldMod use POP_GridHorzMod use POP_HaloMod implicit none !---------------------------------------------------------------------- ! ! local variables ! !---------------------------------------------------------------------- integer (POP_i4) :: & errorCode, &! error flag istat, &! allocate error flag icount, &! counter for validation nProcsClinic, &! number of processors in clinic distribution nBlocksClinic, &! number of local blocks in clinic distribution blockID, &! global block id ib,ie,jb,je, &! beg,end indices for physical domain on block ig, jg, &! global indices for array i,j,k,l,n ! loop indices type (POP_block) :: & thisBlock ! block info for local block type (POP_distrb) :: & distrbClinic ! baroclinic distribution type (POP_halo) :: & haloTripole ! halo data type for tripole updates integer (POP_i4), dimension(:), allocatable :: & workPerBlock ! test values integer (POP_i4), dimension(:), pointer :: & iGlobal ! global index for this block integer (POP_i4), dimension(:), pointer :: & jGlobal ! global index for this block integer (POP_i4), & dimension(POP_nxBlock, POP_nyBlock, POP_maxBlocksClinic) :: & i4Test, i4Expect, i4Orig real (POP_r4), & dimension(POP_nxBlock, POP_nyBlock, POP_maxBlocksClinic) :: & r4Test, r4Expect, r4Orig ! test values real (POP_r8), & dimension(POP_nxBlock, POP_nyBlock, POP_maxBlocksClinic) :: & r8Test, r8Expect, r8Orig ! test values logical (POP_logical), & dimension(POP_nxBlock, POP_nyBlock, POP_maxBlocksClinic) :: & mask integer (POP_i4), dimension(POP_nxBlock, POP_nyBlock, & POP_km, POP_maxBlocksClinic) :: & i4Test3D real (POP_r4), dimension(POP_nxBlock, POP_nyBlock, & POP_km, POP_maxBlocksClinic) :: & r4Test3D real (POP_r8), dimension(POP_nxBlock, POP_nyBlock, & POP_km, POP_maxBlocksClinic) :: & r8Test3D integer (POP_i4), dimension(POP_nxBlock, POP_nyBlock, & POP_km, POP_nt, & POP_maxBlocksClinic) :: & i4Test4D real (POP_r4), dimension(POP_nxBlock, POP_nyBlock, & POP_km, POP_nt, & POP_maxBlocksClinic) :: & r4Test4D real (POP_r8), dimension(POP_nxBlock, POP_nyBlock, & POP_km, POP_nt, & POP_maxBlocksClinic) :: & r8Test4D integer (POP_i4), dimension(POP_nxGlobal, POP_nyGlobal) :: & i4ArrayG !---------------------------------------------------------------------- ! ! initialize communcation environment ! !---------------------------------------------------------------------- call POP_CommInitMessageEnvironment call POP_CommInit errorCode = POP_Success !---------------------------------------------------------------------- ! ! create block decomposition and distributions ! !---------------------------------------------------------------------- !*** create block decomposition call POP_BlocksCreate(POP_nxGlobal, POP_nyGlobal, & 'cyclic', 'tripole', errorCode) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'HaloTest: error creating block decomposition') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif !*** initialize a global array with land in the middle do j=1,POP_nyGlobal do i=1,POP_nxGlobal if (i > POP_nxGlobal/2 + POP_nxBlock .or. & i < POP_nxGlobal/2 - POP_nxBlock .or. & j > POP_nyGlobal/2 + POP_nyBlock .or. & j < POP_nyGlobal/2 - POP_nyBlock) then i4ArrayG(i,j) = (i+j)*100_POP_i4 else i4ArrayG(i,j) = 0_POP_i4 endif end do end do !*** create artificial work per block with varying work and !*** at least one empty block allocate(workPerBlock(POP_numBlocks), stat=istat) if (istat > 0) then call POP_ErrorSet(errorCode, & 'HaloTest: error allocating work per block') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif do blockID=1,POP_numBlocks call POP_BlocksGetBlockInfo(blockID, errorCode, & ib=ib, ie=ie, jb=jb, je=je, & iGlobal=iGlobal, jGlobal=jGlobal) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'HaloTest: error getting block info') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif workPerBlock(blockID) = 0 do j=jb,je do i=ib,ie ig = iGlobal(i) jg = jGlobal(j) if (ig > 0 .and. jg > 0) then if (i4ArrayG(iGlobal(i),jGlobal(j)) /= 0) then workPerBlock(blockID) = workPerBlock(blockID) + 1 endif endif end do end do end do !*** create baroclinic distribution nprocsClinic = POP_CommGetNumProcs(POP_Communicator) distrbClinic = POP_DistributionCreate(POP_distribMethodRake, & nProcsClinic, & workPerBlock, errorCode) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'HaloTest: error creating clinic distribution') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif !---------------------------------------------------------------------- ! ! initialize halo updates ! !---------------------------------------------------------------------- haloTripole = POP_HaloCreate(distrbClinic, 'tripole', 'cyclic', & POP_nxGlobal, errorCode) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'HaloTest: error creating tripole halo') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif !---------------------------------------------------------------------- ! ! initialize test arrays by scattering global arrays to distribution ! !---------------------------------------------------------------------- i4Test = -999 r4Test = -999.0_POP_r4 r8Test = -999.0_POP_r8 call POP_RedistributeScatter(i4Test, i4ArrayG, & POP_masterTask, distrbClinic, errorCode) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'HaloTest: error in i4 scatter') endif mask = (i4Test /= 0) r4Test = i4Test*10.0_POP_r4 r8Test = i4Test*100.0_POP_r8 !---------------------------------------------------------------------- ! ! set up expected values for tripole case ! !---------------------------------------------------------------------- call POP_DistributionGet(distrbClinic, errorCode, & numLocalBlocks=nBlocksClinic) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'HaloTest: error getting local block count') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif if (nBlocksClinic > POP_maxBlocksClinic) then call POP_ErrorSet(errorCode, & 'HaloTest: num local blocks exceeds maxBlocks parameter') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'HaloTest: error getting block ID') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'HaloTest: error getting block') call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment stop endif do j=1,POP_nyBlock do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) jg = thisBlock%jGlobal(j) if (ig <= 0 .or. jg <= 0) then ! outside physical domain i4Orig(i,j,n) = 0_POP_i4 else i4Orig(i,j,n) = i4ArrayG(ig,jg) endif end do end do end do i4Expect = i4Orig r4Orig = i4Orig*10.0_POP_r4 r4Expect = i4Orig*10.0_POP_r4 r8Orig = i4Orig*100.0_POP_r8 r8Expect = i4Orig*100.0_POP_r8 !---------------------------------------------------------------------- ! ! test halo updates at cell centers ! !---------------------------------------------------------------------- !*** set expected values in northern halo region do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary do j=1,POP_haloWidth jg = POP_nyGlobal + 1 - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig /= 0) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal-ig+1,jg) end if end do end do endif end do r4Expect = i4Expect*10.0_POP_r4 r8Expect = i4Expect*100.0_POP_r8 do n=1,nBlocksClinic do k=1,POP_km r8Test3d(:,:,k,n) = r8Test(:,:,n)*k r4Test3d(:,:,k,n) = r4Test(:,:,n)*k i4Test3d(:,:,k,n) = i4Test(:,:,n)*k end do end do do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km r8Test4d(:,:,k,l,n) = r8Test(:,:,n)*k*l r4Test4d(:,:,k,l,n) = r4Test(:,:,n)*k*l i4Test4d(:,:,k,l,n) = i4Test(:,:,n)*k*l end do end do end do !*** test 2d arrays call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole halo at cell center') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole halo - centers') endif call POP_HaloUpdate(r4Test, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r4 tripole halo at cell center') endif if (count(r4Test /= r4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r4 tripole halo - centers') endif call POP_HaloUpdate(i4Test, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d i4 tripole halo at cell center') endif if (count(i4Test /= i4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d i4 tripole halo - centers') endif !*** test 3d arrays call POP_HaloUpdate(r8Test3D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r8 tripole halo at cell center') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r8Test3D(:,:,k,n) /= r8Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r8 tripole halo - centers') endif call POP_HaloUpdate(r4Test3D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r4 tripole halo at cell center') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r4Test3D(:,:,k,n) /= r4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r4 tripole halo - centers') endif call POP_HaloUpdate(i4Test3D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d i4 tripole halo at cell center') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(i4Test3D(:,:,k,n) /= i4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d i4 tripole halo - centers') endif !*** test 4d arrays call POP_HaloUpdate(r8Test4D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r8 tripole halo at cell center') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r8Test4D(:,:,k,l,n) /= & r8Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r8 tripole halo - centers') endif call POP_HaloUpdate(r4Test4D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r4 tripole halo at cell center') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r4Test4D(:,:,k,l,n) /= & r4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r4 tripole halo - centers') endif call POP_HaloUpdate(i4Test4D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d i4 tripole halo at cell center') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(i4Test4D(:,:,k,l,n) /= & i4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d i4 tripole halo - centers') endif !---------------------------------------------------------------------- ! ! test halo updates at east faces ! !---------------------------------------------------------------------- !*** set expected values in northern halo region i4Expect = i4Orig i4Test = i4Orig r4Test = r4Orig r8Test = r8Orig do n=1,nBlocksClinic do k=1,POP_km r8Test3d(:,:,k,n) = r8Test(:,:,n)*k r4Test3d(:,:,k,n) = r4Test(:,:,n)*k i4Test3d(:,:,k,n) = i4Test(:,:,n)*k end do end do do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km r8Test4d(:,:,k,l,n) = r8Test(:,:,n)*k*l r4Test4d(:,:,k,l,n) = r4Test(:,:,n)*k*l i4Test4d(:,:,k,l,n) = i4Test(:,:,n)*k*l end do end do end do do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary do j=1,POP_haloWidth jg = POP_nyGlobal + 1 - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig < POP_nxGlobal .and. & ig /= POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal-ig,jg) else if (ig == POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal/2,jg) else if (ig == POP_nxGlobal) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal,jg) else i4Expect(i,thisBlock%je+j,n) = 0_POP_i4 end if end do end do endif end do r4Expect = i4Expect*10.0_POP_r4 r8Expect = i4Expect*100.0_POP_r8 call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole halo on east faces') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole halo - east face') endif call POP_HaloUpdate(r4Test, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r4 tripole halo on east faces') endif if (count(r4Test /= r4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r4 tripole halo - east face') endif call POP_HaloUpdate(i4Test, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d i4 tripole halo on east faces') endif if (count(i4Test /= i4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d i4 tripole halo - east face') endif !*** test 3d arrays call POP_HaloUpdate(r8Test3D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r8 tripole halo at east face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r8Test3D(:,:,k,n) /= r8Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r8 tripole halo - east face') endif call POP_HaloUpdate(r4Test3D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r4 tripole halo at east face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r4Test3D(:,:,k,n) /= r4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r4 tripole halo - east face') endif call POP_HaloUpdate(i4Test3D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d i4 tripole halo at east face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(i4Test3D(:,:,k,n) /= i4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d i4 tripole halo - east face') endif !*** test 4d arrays call POP_HaloUpdate(r8Test4D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r8 tripole halo at east face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r8Test4D(:,:,k,l,n) /= & r8Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r8 tripole halo - east face') endif call POP_HaloUpdate(r4Test4D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r4 tripole halo at east face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r4Test4D(:,:,k,l,n) /= & r4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r4 tripole halo - east face') endif call POP_HaloUpdate(i4Test4D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d i4 tripole halo at cell east face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(i4Test4D(:,:,k,l,n) /= & i4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d i4 tripole halo - east face') endif !---------------------------------------------------------------------- ! ! test halo updates at north faces ! !---------------------------------------------------------------------- !*** set expected values in northern halo region i4Expect = i4Orig i4Test = i4Orig r4Test = r4Orig r8Test = r8Orig do n=1,nBlocksClinic do k=1,POP_km r8Test3d(:,:,k,n) = r8Test(:,:,n)*k r4Test3d(:,:,k,n) = r4Test(:,:,n)*k i4Test3d(:,:,k,n) = i4Test(:,:,n)*k end do end do do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km r8Test4d(:,:,k,l,n) = r8Test(:,:,n)*k*l r4Test4d(:,:,k,l,n) = r4Test(:,:,n)*k*l i4Test4d(:,:,k,l,n) = i4Test(:,:,n)*k*l end do end do end do do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary !*** set ghost cells do j=1,POP_haloWidth jg = POP_nyGlobal - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig /= 0) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal-ig+1,jg) end if end do end do !*** enforce symmetry at top physical row jg = POP_nyGlobal do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig /= 0) then i4Expect(i,thisBlock%je,n) = 0.5_POP_r8*( & i4ArrayG( ig ,jg) + & i4ArrayG(POP_nxGlobal-ig+1,jg)) end if end do endif end do r4Expect = i4Expect*10.0_POP_r4 r8Expect = i4Expect*100.0_POP_r8 call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole halo on north faces') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole halo - north faces') endif call POP_HaloUpdate(r4Test, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r4 tripole halo on north faces') endif if (count(r4Test /= r4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r4 tripole halo - north faces') endif call POP_HaloUpdate(i4Test, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d i4 tripole halo on north faces') endif if (count(i4Test /= i4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d i4 tripole halo - north faces') endif !*** test 3d arrays call POP_HaloUpdate(r8Test3D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r8 tripole halo at north face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r8Test3D(:,:,k,n) /= r8Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r8 tripole halo - north face') endif call POP_HaloUpdate(r4Test3D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r4 tripole halo at north face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r4Test3D(:,:,k,n) /= r4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r4 tripole halo - north face') endif call POP_HaloUpdate(i4Test3D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d i4 tripole halo at north face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(i4Test3D(:,:,k,n) /= i4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d i4 tripole halo - north face') endif !*** test 4d arrays call POP_HaloUpdate(r8Test4D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r8 tripole halo at north face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r8Test4D(:,:,k,l,n) /= & r8Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r8 tripole halo - north face') endif call POP_HaloUpdate(r4Test4D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r4 tripole halo at north face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r4Test4D(:,:,k,l,n) /= & r4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r4 tripole halo - north face') endif call POP_HaloUpdate(i4Test4D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d i4 tripole halo at north face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(i4Test4D(:,:,k,l,n) /= & i4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d i4 tripole halo - north face') endif !---------------------------------------------------------------------- ! ! test halo updates at northeast corners ! !---------------------------------------------------------------------- !*** set expected values in northern halo region and reset land !*** mask i4Expect = i4Orig i4Test = i4Orig r4Test = r4Orig r8Test = r8Orig do n=1,nBlocksClinic do k=1,POP_km r8Test3d(:,:,k,n) = r8Test(:,:,n)*k r4Test3d(:,:,k,n) = r4Test(:,:,n)*k i4Test3d(:,:,k,n) = i4Test(:,:,n)*k end do end do do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km r8Test4d(:,:,k,l,n) = r8Test(:,:,n)*k*l r4Test4d(:,:,k,l,n) = r4Test(:,:,n)*k*l i4Test4d(:,:,k,l,n) = i4Test(:,:,n)*k*l end do end do end do do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary !*** set halo points do j=1,POP_haloWidth jg = POP_nyGlobal - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig /= POP_nxGlobal .and. & ig /= POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal-ig,jg) else if (ig == POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal/2,jg) else if (ig == POP_nxGlobal) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal,jg) else i4Expect(i,thisBlock%je+j,n) = 0_POP_i4 end if end do end do !*** enforce symmetry jg = POP_nyGlobal do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig < POP_nxGlobal .and. & ig /= POP_nxGlobal/2) then i4Expect(i,thisBlock%je,n) = 0.5_POP_r8*( & i4ArrayG( ig,jg) + & i4ArrayG(POP_nxGlobal-ig,jg)) else if (ig == POP_nxGlobal .or. & ig == POP_nxGlobal/2) then i4Expect(i,thisBlock%je,n) = i4ArrayG(ig,jg) else i4Expect(i,thisBlock%je,n) = 0_POP_i4 end if end do endif end do r4Expect = i4Expect*10.0_POP_r4 r8Expect = i4Expect*100.0_POP_r8 call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole halo') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole halo - NE Corner') endif call POP_HaloUpdate(r4Test, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r4 tripole halo') endif if (count(r4Test /= r4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r4 tripole halo - NE Corner') endif call POP_HaloUpdate(i4Test, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d i4 tripole halo') endif if (count(i4Test /= i4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d i4 tripole halo - NE Corner') endif !*** test 3d arrays call POP_HaloUpdate(r8Test3D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r8 tripole halo at NE corner') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r8Test3D(:,:,k,n) /= r8Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r8 tripole halo - NE corner') endif call POP_HaloUpdate(r4Test3D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r4 tripole halo at NE corner') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r4Test3D(:,:,k,n) /= r4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r4 tripole halo - NE corner') endif call POP_HaloUpdate(i4Test3D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d i4 tripole halo at NE corner') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(i4Test3D(:,:,k,n) /= i4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d i4 tripole halo - NE corner') endif !*** test 4d arrays call POP_HaloUpdate(r8Test4D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r8 tripole halo at NE corner') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r8Test4D(:,:,k,l,n) /= & r8Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r8 tripole halo - NE corner') endif call POP_HaloUpdate(r4Test4D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r4 tripole halo at NE corner') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r4Test4D(:,:,k,l,n) /= & r4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r4 tripole halo - NE corner') endif call POP_HaloUpdate(i4Test4D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d i4 tripole halo at NE corner') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(i4Test4D(:,:,k,l,n) /= & i4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d i4 tripole halo - NE corner') endif !---------------------------------------------------------------------- ! ! test roundoff in symmetry enforcement by doing repeated halo ! updates and checking for equivalence ! use halo updates at northeast corners ! !---------------------------------------------------------------------- !*** set expected values in northern halo region and reset land !*** mask i4Expect = i4Orig i4Test = i4Orig r8Test = r8Orig do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary !*** set halo points do j=1,POP_haloWidth jg = POP_nyGlobal - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig /= POP_nxGlobal .and. & ig /= POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal-ig,jg) else if (ig == POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal/2,jg) else if (ig == POP_nxGlobal) then i4Expect(i,thisBlock%je+j,n) = i4ArrayG(POP_nxGlobal,jg) else i4Expect(i,thisBlock%je+j,n) = 0_POP_i4 end if end do end do !*** enforce symmetry jg = POP_nyGlobal do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig < POP_nxGlobal .and. & ig /= POP_nxGlobal/2) then i4Expect(i,thisBlock%je,n) = 0.5_POP_r8*( & i4ArrayG( ig,jg) + & i4ArrayG(POP_nxGlobal-ig,jg)) else if (ig == POP_nxGlobal .or. & ig == POP_nxGlobal/2) then i4Expect(i,thisBlock%je,n) = i4ArrayG(ig,jg) else i4Expect(i,thisBlock%je,n) = 0_POP_i4 end if end do endif end do r8Expect = i4Expect*100.0_POP_r8 call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole halo') endif call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindScalar, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole halo') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole halo - roundoff') endif !---------------------------------------------------------------------- ! ! test halo updates for vectors at cell centers ! !---------------------------------------------------------------------- !*** set expected values in northern halo region and reset land !*** mask i4Expect = i4Orig i4Test = i4Orig r4Test = r4Orig r8Test = r8Orig do n=1,nBlocksClinic do k=1,POP_km r8Test3d(:,:,k,n) = r8Test(:,:,n)*k r4Test3d(:,:,k,n) = r4Test(:,:,n)*k i4Test3d(:,:,k,n) = i4Test(:,:,n)*k end do end do do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km r8Test4d(:,:,k,l,n) = r8Test(:,:,n)*k*l r4Test4d(:,:,k,l,n) = r4Test(:,:,n)*k*l i4Test4d(:,:,k,l,n) = i4Test(:,:,n)*k*l end do end do end do do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary !*** set halo points do j=1,POP_haloWidth jg = POP_nyGlobal + 1 - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig <= POP_nxGlobal) then i4Expect(i,thisBlock%je+j,n) = -i4ArrayG(POP_nxGlobal-ig+1,jg) else i4Expect(i,thisBlock%je+j,n) = 0_POP_i4 end if end do end do endif end do r4Expect = i4Expect*10.0_POP_r4 r8Expect = i4Expect*100.0_POP_r8 call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole vector') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole vector - center') endif call POP_HaloUpdate(r4Test, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r4 tripole vector') endif if (count(r4Test /= r4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r4 tripole vector - center') endif call POP_HaloUpdate(i4Test, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d i4 tripole vector') endif if (count(i4Test /= i4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d i4 tripole vector - center') endif !*** test 3d arrays call POP_HaloUpdate(r8Test3D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r8 tripole vector at center') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r8Test3D(:,:,k,n) /= r8Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r8 tripole vector - center') endif call POP_HaloUpdate(r4Test3D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r4 tripole vector at NE corner') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r4Test3D(:,:,k,n) /= r4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r4 tripole vector - center') endif call POP_HaloUpdate(i4Test3D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d i4 tripole vector at center') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(i4Test3D(:,:,k,n) /= i4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d i4 tripole vector - center') endif !*** test 4d arrays call POP_HaloUpdate(r8Test4D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r8 tripole vector at center') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r8Test4D(:,:,k,l,n) /= & r8Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r8 tripole vector - center') endif call POP_HaloUpdate(r4Test4D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r4 tripole vector at center') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r4Test4D(:,:,k,l,n) /= & r4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r4 tripole vector - center') endif call POP_HaloUpdate(i4Test4D, haloTripole, & POP_gridHorzLocCenter, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d i4 tripole vector at center') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(i4Test4D(:,:,k,l,n) /= & i4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d i4 tripole vector - center') endif !---------------------------------------------------------------------- ! ! test halo updates for vectors at east faces ! !---------------------------------------------------------------------- !*** set expected values in northern halo region and reset land !*** mask i4Expect = i4Orig i4Test = i4Orig r4Test = r4Orig r8Test = r8Orig do n=1,nBlocksClinic do k=1,POP_km r8Test3d(:,:,k,n) = r8Test(:,:,n)*k r4Test3d(:,:,k,n) = r4Test(:,:,n)*k i4Test3d(:,:,k,n) = i4Test(:,:,n)*k end do end do do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km r8Test4d(:,:,k,l,n) = r8Test(:,:,n)*k*l r4Test4d(:,:,k,l,n) = r4Test(:,:,n)*k*l i4Test4d(:,:,k,l,n) = i4Test(:,:,n)*k*l end do end do end do do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary !*** set halo points do j=1,POP_haloWidth jg = POP_nyGlobal + 1 - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig < POP_nxGlobal .and. & ig /= POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = -i4ArrayG(POP_nxGlobal-ig,jg) else if (ig == POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = -i4ArrayG(POP_nxGlobal/2,jg) else if (ig == POP_nxGlobal) then i4Expect(i,thisBlock%je+j,n) = -i4ArrayG(POP_nxGlobal,jg) else i4Expect(i,thisBlock%je+j,n) = 0_POP_i4 end if end do end do endif end do r4Expect = i4Expect*10.0_POP_r4 r8Expect = i4Expect*100.0_POP_r8 call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole vector') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole vector - E Face') endif call POP_HaloUpdate(r4Test, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r4 tripole vector') endif if (count(r4Test /= r4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r4 tripole vector - E Face') endif call POP_HaloUpdate(i4Test, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d i4 tripole vector') endif if (count(i4Test /= i4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d i4 tripole vector - E Face') endif !*** test 3d arrays call POP_HaloUpdate(r8Test3D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r8 tripole vector at E Face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r8Test3D(:,:,k,n) /= r8Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r8 tripole vector - E Face') endif call POP_HaloUpdate(r4Test3D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r4 tripole vector at E Face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r4Test3D(:,:,k,n) /= r4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r4 tripole vector - E Face') endif call POP_HaloUpdate(i4Test3D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d i4 tripole vector at E Face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(i4Test3D(:,:,k,n) /= i4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d i4 tripole vector - E Face') endif !*** test 4d arrays call POP_HaloUpdate(r8Test4D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r8 tripole vector at E Face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r8Test4D(:,:,k,l,n) /= & r8Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r8 tripole vector - E Face') endif call POP_HaloUpdate(r4Test4D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r4 tripole vector at E Face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r4Test4D(:,:,k,l,n) /= & r4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r4 tripole vector - E Face') endif call POP_HaloUpdate(i4Test4D, haloTripole, & POP_gridHorzLocEFace, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d i4 tripole vector at E Face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(i4Test4D(:,:,k,l,n) /= & i4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d i4 tripole vector - E Face') endif !---------------------------------------------------------------------- ! ! test halo updates for vectors at north faces ! !---------------------------------------------------------------------- !*** set expected values in northern halo region and reset land !*** mask i4Expect = i4Orig i4Test = i4Orig r4Test = r4Orig r8Test = r8Orig do n=1,nBlocksClinic do k=1,POP_km r8Test3d(:,:,k,n) = r8Test(:,:,n)*k r4Test3d(:,:,k,n) = r4Test(:,:,n)*k i4Test3d(:,:,k,n) = i4Test(:,:,n)*k end do end do do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km r8Test4d(:,:,k,l,n) = r8Test(:,:,n)*k*l r4Test4d(:,:,k,l,n) = r4Test(:,:,n)*k*l i4Test4d(:,:,k,l,n) = i4Test(:,:,n)*k*l end do end do end do do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary !*** set halo points do j=1,POP_haloWidth jg = POP_nyGlobal - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig <= POP_nxGlobal) then i4Expect(i,thisBlock%je+j,n) = -i4ArrayG(POP_nxGlobal+1-ig,jg) else i4Expect(i,thisBlock%je+j,n) = 0_POP_i4 end if end do end do !*** enforce symmetry jg = POP_nyGlobal do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig <= POP_nxGlobal) then i4Expect(i,thisBlock%je,n) = 0.5_POP_r8*( & i4ArrayG( ig,jg) + & i4ArrayG(POP_nxGlobal+1-ig,jg)) else i4Expect(i,thisBlock%je,n) = 0_POP_i4 end if end do endif end do r4Expect = i4Expect*10.0_POP_r4 r8Expect = i4Expect*100.0_POP_r8 call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole vector') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole vector - N Face') endif call POP_HaloUpdate(r4Test, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r4 tripole vector') endif if (count(r4Test /= r4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r4 tripole vector - N Face') endif call POP_HaloUpdate(i4Test, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d i4 tripole vector') endif if (count(i4Test /= i4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d i4 tripole vector - N Face') endif !*** test 3d arrays call POP_HaloUpdate(r8Test3D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r8 tripole vector at N Face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r8Test3D(:,:,k,n) /= r8Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r8 tripole vector - N Face') endif call POP_HaloUpdate(r4Test3D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r4 tripole vector at N Face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r4Test3D(:,:,k,n) /= r4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r4 tripole vector - N Face') endif call POP_HaloUpdate(i4Test3D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d i4 tripole vector at N Face') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(i4Test3D(:,:,k,n) /= i4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d i4 tripole vector - N Face') endif !*** test 4d arrays call POP_HaloUpdate(r8Test4D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r8 tripole vector at N Face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r8Test4D(:,:,k,l,n) /= & r8Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r8 tripole vector - N Face') endif call POP_HaloUpdate(r4Test4D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r4 tripole vector at N Face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r4Test4D(:,:,k,l,n) /= & r4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r4 tripole vector - N Face') endif call POP_HaloUpdate(i4Test4D, haloTripole, & POP_gridHorzLocNFace, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d i4 tripole vector at N Face') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(i4Test4D(:,:,k,l,n) /= & i4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d i4 tripole vector - N face') endif !---------------------------------------------------------------------- ! ! test halo updates for vectors at northeast corners ! !---------------------------------------------------------------------- !*** set expected values in northern halo region and reset land !*** mask i4Expect = i4Orig i4Test = i4Orig r4Test = r4Orig r8Test = r8Orig do n=1,nBlocksClinic do k=1,POP_km r8Test3d(:,:,k,n) = r8Test(:,:,n)*k r4Test3d(:,:,k,n) = r4Test(:,:,n)*k i4Test3d(:,:,k,n) = i4Test(:,:,n)*k end do end do do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km r8Test4d(:,:,k,l,n) = r8Test(:,:,n)*k*l r4Test4d(:,:,k,l,n) = r4Test(:,:,n)*k*l i4Test4d(:,:,k,l,n) = i4Test(:,:,n)*k*l end do end do end do do n=1,nBlocksClinic call POP_DistributionGetBlockID(distrbClinic, n, blockID, & errorCode) thisBlock = POP_BlocksGetBlock(blockID, errorCode) if (thisBlock%jGlobal(thisBlock%je+1) < 0) then !*** this block contains the northern boundary !*** set halo points do j=1,POP_haloWidth jg = POP_nyGlobal - j do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig < POP_nxGlobal .and. & ig /= POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = -i4ArrayG(POP_nxGlobal-ig,jg) else if (ig == POP_nxGlobal/2) then i4Expect(i,thisBlock%je+j,n) = -i4ArrayG(POP_nxGlobal/2,jg) else if (ig == POP_nxGlobal) then i4Expect(i,thisBlock%je+j,n) = -i4ArrayG(POP_nxGlobal,jg) else i4Expect(i,thisBlock%je+j,n) = 0_POP_i4 end if end do end do !*** enforce symmetry jg = POP_nyGlobal do i=1,POP_nxBlock ig = thisBlock%iGlobal(i) if (ig > 0 .and. ig < POP_nxGlobal .and. & ig /= POP_nxGlobal/2) then i4Expect(i,thisBlock%je,n) = 0.5_POP_r8*( & i4ArrayG( ig,jg) + & i4ArrayG(POP_nxGlobal-ig,jg)) else if (ig == POP_nxGlobal .or. & ig == POP_nxGlobal/2) then i4Expect(i,thisBlock%je,n) = i4ArrayG(ig,jg) else i4Expect(i,thisBlock%je,n) = 0_POP_i4 end if end do endif end do r4Expect = i4Expect*10.0_POP_r4 r8Expect = i4Expect*100.0_POP_r8 call POP_HaloUpdate(r8Test, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r8 tripole vector') endif if (count(r8Test /= r8Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r8 tripole vector - NE Corner') endif call POP_HaloUpdate(r4Test, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d r4 tripole vector') endif if (count(r4Test /= r4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d r4 tripole vector - NE Corner') endif call POP_HaloUpdate(i4Test, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 2d i4 tripole vector') endif if (count(i4Test /= i4Expect) > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 2d i4 tripole vector - NE Corner') endif !*** test 3d arrays call POP_HaloUpdate(r8Test3D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r8 tripole vector at NE corner') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r8Test3D(:,:,k,n) /= r8Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r8 tripole vector - NE corner') endif call POP_HaloUpdate(r4Test3D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d r4 tripole vector at NE corner') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(r4Test3D(:,:,k,n) /= r4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d r4 tripole vector - NE corner') endif call POP_HaloUpdate(i4Test3D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 3d i4 tripole vector at NE corner') endif icount = 0 do n=1,nBlocksClinic do k=1,POP_km icount = icount + count(i4Test3D(:,:,k,n) /= i4Expect(:,:,n)*k) end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 3d i4 tripole vector - NE corner') endif !*** test 4d arrays call POP_HaloUpdate(r8Test4D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r8) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r8 tripole vector at NE corner') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r8Test4D(:,:,k,l,n) /= & r8Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r8 tripole vector - NE corner') endif call POP_HaloUpdate(r4Test4D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0.0_POP_r4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d r4 tripole vector at NE corner') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(r4Test4D(:,:,k,l,n) /= & r4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d r4 tripole vector - NE corner') endif call POP_HaloUpdate(i4Test4D, haloTripole, & POP_gridHorzLocNECorner, & POP_fieldKindVector, errorCode, & fillValue=0_POP_i4) if (errorCode /= POP_Success) then call POP_ErrorSet(errorCode, & 'Halo test: error in 4d i4 tripole vector at NE corner') endif icount = 0 do n=1,nBlocksClinic do l=1,POP_nt do k=1,POP_km icount = icount + count(i4Test4D(:,:,k,l,n) /= & i4Expect(:,:,n)*k*l) end do end do end do if (icount > 0) then call POP_ErrorSet(errorCode, & 'Halo test: bad results from 4d i4 tripole vector - NE corner') endif !---------------------------------------------------------------------- ! ! clean up ! !---------------------------------------------------------------------- call POP_HaloDestroy(haloTripole, errorCode) call POP_ErrorPrint(errorCode) call POP_CommExitMessageEnvironment !---------------------------------------------------------------------- end program POP_HaloTest !||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||