Back to home page

MITgcm

 
 

    


File indexing completed on 2020-07-29 05:11:13 UTC

view on githubraw file Latest commit b9dadda2 on 2020-07-28 16:49:33 UTC
ca64f811ee Jean*0001 #include "CPP_EEOPTIONS.h"
                0002 #include "W2_OPTIONS.h"
                0003 
                0004 CBOP 0
                0005 C !ROUTINE: W2_SET_MAP_CUMSUM
                0006 
                0007 C !INTERFACE:
                0008       SUBROUTINE W2_SET_MAP_CUMSUM( myThid )
                0009 
                0010 C     !DESCRIPTION:
                0011 C     Set-up mapping for global cumulated sum.
                0012 
                0013 C     !USES:
                0014       IMPLICIT NONE
                0015 
                0016 C      Tile topology settings data structures
                0017 #include "SIZE.h"
                0018 #include "EEPARAMS.h"
                0019 #include "W2_EXCH2_SIZE.h"
                0020 #include "W2_EXCH2_PARAMS.h"
                0021 #include "W2_EXCH2_TOPOLOGY.h"
                0022 
                0023 C     !INPUT PARAMETERS:
                0024 C     myThid  :: my Thread Id number
                0025 C               (Note: not relevant since threading has not yet started)
                0026       INTEGER myThid
                0027 
                0028 C     !LOCAL VARIABLES:
                0029 C     === Local variables ===
                0030 C     facetXYSum   :: Tile to Facet Matrix for facet-Increment in X & Y dir.
                0031 C     facet_CSum   :: Tile to Facet Matrix for CumSum @ facet-origin
                0032 C     msgBuf       :: Informational/error message buffer
                0033       INTEGER fNx, fNy, nbTx, nbTy
                0034       INTEGER nActiveFacets
                0035       INTEGER fCnt, prev_fCnt
                0036       INTEGER npass, nType
                0037       LOGICAL prtFlag
                0038       LOGICAL fIsSet(0:W2_maxNbFacets)
575f97366a Jean*0039       INTEGER tN, i, j, k, ii, jj
                0040 #ifdef W2_CUMSUM_USE_MATRIX
                0041       INTEGER tS, bi, bj, l, is, ie
ca64f811ee Jean*0042       INTEGER facetXYSum(2,W2_maxNbTiles,W2_maxNbFacets)
                0043       INTEGER facet_CSum(2,W2_maxNbTiles,W2_maxNbFacets)
575f97366a Jean*0044 #endif
ca64f811ee Jean*0045       CHARACTER*(MAX_LEN_MBUF) msgBuf
                0046 CEOP
                0047 c      DATA edge / 'N' , 'S' , 'E' , 'W' /
                0048 
                0049       WRITE(msgBuf,'(2A)') 'W2_SET_MAP_CUMSUM: ',
                0050      &       'setting Facet Matrix for CUMUL-SUM'
                0051       CALL PRINT_MESSAGE( msgBuf, W2_oUnit, SQUEEZE_RIGHT, myThid )
                0052       prtFlag = ABS(W2_printMsg).GE.2
                0053      &       .OR. ( W2_printMsg .NE.0 .AND. myProcId.EQ.0 )
                0054 
                0055 C--   Initialise Common-block
                0056       W2_tMC1 = 0
                0057       W2_tMC2 = 0
575f97366a Jean*0058       DO j=1,W2_maxNbFacets
                0059        DO i=1,W2_maxNbFacets
                0060          W2_cumSum_facet(1,i,j) = 0
                0061          W2_cumSum_facet(2,i,j) = 0
                0062        ENDDO
                0063       ENDDO
                0064 #ifdef W2_CUMSUM_USE_MATRIX
ca64f811ee Jean*0065       DO j=1,W2_maxNbTiles
                0066        DO i=1,W2_maxNbTiles
                0067          W2_cumSum_tiles(1,i,j) = 0
                0068          W2_cumSum_tiles(2,i,j) = 0
                0069        ENDDO
                0070       ENDDO
575f97366a Jean*0071 #endif /* W2_CUMSUM_USE_MATRIX */
ca64f811ee Jean*0072 
                0073 C--   Start setting cumul-sum Face mapping
                0074       fCnt = 0
                0075       fIsSet(0) = .TRUE.
                0076       DO j=1,nFacets
                0077         fIsSet(j) = .FALSE.
                0078       ENDDO
                0079 
                0080 C--   Start with first non-empty face:
                0081       nActiveFacets = 0
                0082       DO j=1,nFacets
                0083         IF ( facet_dims(2*j-1)*facet_dims(2*j).GE.1 ) THEN
                0084           nActiveFacets = nActiveFacets + 1
                0085           IF ( fCnt.EQ.0 ) THEN
                0086             fIsSet(j) = .TRUE.
                0087             fCnt = 1
                0088             IF ( ABS(W2_printMsg).GE.2 ) WRITE(W2_oUnit,'(A,I4)')
                0089      &        ' CumSum starts @ SW.corner of facet #', j
                0090           ENDIF
                0091         ENDIF
                0092       ENDDO
                0093 
                0094 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0095 
                0096 C--   Go through list of connections:
                0097       prev_fCnt = 0
                0098       npass = 0
                0099       DO WHILE ( fCnt.GT.prev_fCnt )
                0100        npass = npass + 1
                0101        prev_fCnt = fCnt
                0102        DO nType=1,3
                0103         IF ( fCnt.EQ.prev_fCnt ) THEN
                0104          DO j=1,nFacets
                0105           IF ( fIsSet(j) ) THEN
                0106            DO i=1,4
                0107 C-    connected to:
                0108             jj = INT(facet_link(i,j))
                0109             ii = MOD( NINT(facet_link(i,j)*10.), 10 )
                0110 C--   1) with same orientation ( N <-> S or E <-> W ), forward progression
                0111             IF ( .NOT.fIsSet(jj) .AND. nType.EQ.1 ) THEN
                0112 C-    N <-- S ;
                0113              IF ( i.EQ.W2_NORTH .AND. ii.EQ.W2_SOUTH ) THEN
                0114               DO k=1,nFacets
575f97366a Jean*0115                W2_cumSum_facet(1,k,jj) = W2_cumSum_facet(1,k,j)
                0116                W2_cumSum_facet(2,k,jj) = W2_cumSum_facet(2,k,j)
ca64f811ee Jean*0117               ENDDO
575f97366a Jean*0118               W2_cumSum_facet(2,j,jj) = W2_cumSum_facet(2,j,jj) + 1
ca64f811ee Jean*0119               fCnt = fCnt + 1
                0120               fIsSet(jj) = .TRUE.
                0121               IF ( ABS(W2_printMsg).GE.2 ) WRITE(W2_oUnit,'(5(A,I4))')
                0122      &         ' CumSum SW.corner of facet #', jj,' set from facet',j,
                0123      &         ' (pass,type=', npass,',', nType,')'
                0124              ENDIF
                0125 C-    E <-- W ;
                0126              IF ( i.EQ.W2_EAST .AND. ii.EQ.W2_WEST ) THEN
                0127               DO k=1,nFacets
575f97366a Jean*0128                W2_cumSum_facet(1,k,jj) = W2_cumSum_facet(1,k,j)
                0129                W2_cumSum_facet(2,k,jj) = W2_cumSum_facet(2,k,j)
ca64f811ee Jean*0130               ENDDO
575f97366a Jean*0131               W2_cumSum_facet(1,j,jj) = W2_cumSum_facet(1,j,jj) + 1
ca64f811ee Jean*0132               fCnt = fCnt + 1
                0133               fIsSet(jj) = .TRUE.
                0134               IF ( ABS(W2_printMsg).GE.2 ) WRITE(W2_oUnit,'(5(A,I4))')
                0135      &         ' CumSum SW.corner of facet #', jj,' set from facet',j,
                0136      &         ' (pass,type=', npass,',', nType,')'
                0137              ENDIF
                0138             ENDIF
                0139 C--   2) with same orientation ( N <-> S or E <-> W ), backward progression
                0140             IF ( .NOT.fIsSet(jj) .AND. nType.EQ.2 ) THEN
                0141 C-    S <-- N ;
                0142              IF ( i.EQ.W2_SOUTH .AND. ii.EQ.W2_NORTH ) THEN
                0143               DO k=1,nFacets
575f97366a Jean*0144                W2_cumSum_facet(1,k,jj) = W2_cumSum_facet(1,k,j)
                0145                W2_cumSum_facet(2,k,jj) = W2_cumSum_facet(2,k,j)
ca64f811ee Jean*0146               ENDDO
575f97366a Jean*0147               W2_cumSum_facet(2,jj,jj) = W2_cumSum_facet(2,jj,jj) - 1
ca64f811ee Jean*0148               fCnt = fCnt + 1
                0149               fIsSet(jj) = .TRUE.
                0150               IF ( ABS(W2_printMsg).GE.2 ) WRITE(W2_oUnit,'(5(A,I4))')
                0151      &         ' CumSum SW.corner of facet #', jj,' set from facet',j,
                0152      &         ' (pass,type=', npass,',', nType,')'
                0153              ENDIF
                0154 C-    W <-- E ;
                0155              IF ( i.EQ.W2_WEST .AND. ii.EQ.W2_EAST ) THEN
                0156               DO k=1,nFacets
575f97366a Jean*0157                W2_cumSum_facet(1,k,jj) = W2_cumSum_facet(1,k,j)
                0158                W2_cumSum_facet(2,k,jj) = W2_cumSum_facet(2,k,j)
ca64f811ee Jean*0159               ENDDO
575f97366a Jean*0160               W2_cumSum_facet(1,jj,jj) = W2_cumSum_facet(1,jj,jj) - 1
ca64f811ee Jean*0161               fCnt = fCnt + 1
                0162               fIsSet(jj) = .TRUE.
                0163               IF ( ABS(W2_printMsg).GE.2 ) WRITE(W2_oUnit,'(5(A,I4))')
                0164      &         ' CumSum SW.corner of facet #', jj,' set from facet',j,
                0165      &         ' (pass,type=', npass,',', nType,')'
                0166              ENDIF
                0167             ENDIF
                0168 C--   3) with different orientation ( N <-> W or S <-> E )
                0169 C- Note: cannot rely on these connections for Cumul-Sum @ grid-cell center
                0170             IF ( .NOT.fIsSet(jj) .AND. nType.EQ.3
                0171      &                           .AND. fCnt.EQ.prev_fCnt ) THEN
                0172 C-    N <-- W or W <-- N ;
                0173              IF (  ( i.EQ.W2_NORTH .AND. ii.EQ.W2_WEST  ) .OR.
                0174      &             ( i.EQ.W2_WEST  .AND. ii.EQ.W2_NORTH ) ) THEN
                0175               DO k=1,nFacets
575f97366a Jean*0176                W2_cumSum_facet(1,k,jj) = W2_cumSum_facet(1,k,j)
                0177                W2_cumSum_facet(2,k,jj) = W2_cumSum_facet(2,k,j)
ca64f811ee Jean*0178               ENDDO
575f97366a Jean*0179               W2_cumSum_facet(2,j, jj) = W2_cumSum_facet(2,j, jj) + 1
                0180               W2_cumSum_facet(2,jj,jj) = W2_cumSum_facet(2,jj,jj) - 1
ca64f811ee Jean*0181               fCnt = fCnt + 1
                0182               fIsSet(jj) = .TRUE.
                0183               IF ( ABS(W2_printMsg).GE.2 ) WRITE(W2_oUnit,'(5(A,I4))')
                0184      &         ' CumSum SW.corner of facet #', jj,' set from facet',j,
                0185      &         ' (pass,type=', npass,',', nType,')'
                0186              ENDIF
                0187 C-    E <-- S or S <-- E ;
                0188              IF (  ( i.EQ.W2_EAST  .AND. ii.EQ.W2_SOUTH ) .OR.
                0189      &             ( i.EQ.W2_SOUTH .AND. ii.EQ.W2_EAST  ) ) THEN
                0190               DO k=1,nFacets
575f97366a Jean*0191                W2_cumSum_facet(1,k,jj) = W2_cumSum_facet(1,k,j)
                0192                W2_cumSum_facet(2,k,jj) = W2_cumSum_facet(2,k,j)
ca64f811ee Jean*0193               ENDDO
575f97366a Jean*0194               W2_cumSum_facet(1,j, jj) = W2_cumSum_facet(1,j, jj) + 1
                0195               W2_cumSum_facet(1,jj,jj) = W2_cumSum_facet(1,jj,jj) - 1
ca64f811ee Jean*0196               fCnt = fCnt + 1
                0197               fIsSet(jj) = .TRUE.
                0198               IF ( ABS(W2_printMsg).GE.2 ) WRITE(W2_oUnit,'(5(A,I4))')
                0199      &         ' CumSum SW.corner of facet #', jj,' set from facet',j,
                0200      &         ' (pass,type=', npass,',', nType,')'
                0201              ENDIF
                0202             ENDIF
                0203 C          end i loop
                0204            ENDDO
                0205           ENDIF
                0206 C        end facets loop
                0207          ENDDO
                0208          IF ( fCnt.GT.prev_fCnt ) THEN
                0209           WRITE(msgBuf,'(2A,3(I4,A),I2,A)') 'W2_SET_MAP_CUMSUM: ',
                0210      &     'set ', fCnt - prev_fCnt, ' /', nActiveFacets,
                0211      &     ' active facets (pass,type=', npass, ',', nType, ')'
                0212           CALL PRINT_MESSAGE( msgBuf, W2_oUnit, SQUEEZE_RIGHT, myThid )
                0213          ENDIF
                0214         ENDIF
                0215 C      end nType loop
                0216        ENDDO
                0217 C-    end do while (npass count)
                0218       ENDDO
                0219       IF ( fCnt.LT.nActiveFacets ) THEN
                0220         WRITE(msgBuf,'(2A,2(I4,A))') 'W2_SET_MAP_CUMSUM: ',
                0221      &  'Only get', fCnt, ' /', nActiveFacets,' active facets done'
                0222         CALL PRINT_MESSAGE( msgBuf, W2_oUnit, SQUEEZE_RIGHT, myThid )
                0223         WRITE(msgBuf,'(2A)') '** WARNING ** W2_SET_MAP_CUMSUM: ',
                0224      &    ' missing connections in Cumulated Sum'
                0225         CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
                0226      &                      SQUEEZE_RIGHT, myThid )
                0227       ENDIF
                0228       IF ( myProcId.EQ.0 ) THEN
                0229         WRITE(W2_oUnit,'(2A,2(I4,A))')' Facet Matrix for CUMUL-SUM (',
                0230      &         'nFacets=',nFacets, ', nActive=', nActiveFacets, ' ):'
                0231         DO j=1,nFacets
                0232           WRITE(W2_oUnit,'(A,I3,A,30(2I3,A))') '- facet', j, ' :',
575f97366a Jean*0233      &         (W2_cumSum_facet(1,i,j),W2_cumSum_facet(2,i,j),' ,',
                0234      &                                 i=1,nFacets)
ca64f811ee Jean*0235         ENDDO
                0236       ENDIF
                0237 
                0238 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
575f97366a Jean*0239 C--   record "missing corner" tile
                0240 
                0241       IF ( useCubedSphereExchange ) THEN
                0242        DO j=1,nFacets
                0243         fNx = facet_dims(2*j-1)
                0244         fNy = facet_dims( 2*j )
                0245         IF ( fNx*fNy .GE. 1 ) THEN
                0246          nbTx = fNx/sNx
                0247          nbTy = fNy/sNy
                0248          tN = facet_owns(1,j) - 1 + nbTx
                0249          IF ( W2_tMC2.EQ.0 .AND. MOD(j,2).EQ.0
                0250      &                     .AND. exch2_myFace(tN).NE.0 ) W2_tMC2 = tN
                0251          tN = facet_owns(1,j) + (nbTy-1)*nbTx
                0252          IF ( W2_tMC1.EQ.0 .AND. MOD(j,2).EQ.1
                0253      &                     .AND. exch2_myFace(tN).NE.0 ) W2_tMC1 = tN
                0254         ENDIF
                0255        ENDDO
b9dadda204 Mart*0256        IF ( myProcId.EQ.0 ) WRITE(W2_oUnit,'(3(A,I8))')
575f97366a Jean*0257      &   ' missing-corner Tile for CUMUL-SUM (nTiles=', exch2_nTiles,
                0258      &   ' ): W2_tMC1=', W2_tMC1, ' , W2_tMC2=', W2_tMC2
                0259       ENDIF
                0260 
                0261 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
ca64f811ee Jean*0262 C--   Now Set cumul-sum Tile mapping
                0263 
575f97366a Jean*0264 #ifdef W2_CUMSUM_USE_MATRIX
ca64f811ee Jean*0265       DO j=1,W2_maxNbFacets
                0266        DO i=1,W2_maxNbTiles
                0267         facetXYSum(1,i,j) = 0
                0268         facetXYSum(2,i,j) = 0
                0269         facet_CSum(1,i,j) = 0
                0270         facet_CSum(2,i,j) = 0
                0271        ENDDO
                0272       ENDDO
                0273 
                0274 C-    First within each face:
                0275       DO j=1,nFacets
                0276        fNx = facet_dims(2*j-1)
                0277        fNy = facet_dims( 2*j )
                0278        IF ( fNx*fNy .GE. 1 ) THEN
                0279         nbTx = fNx/sNx
                0280         nbTy = fNy/sNy
                0281 
                0282         DO bj=1,nbTy
                0283          DO bi=1,nbTx-1
                0284            tS = facet_owns(1,j) - 1 + bi
                0285            tN = tS + 1 + (bj-1)*nbTx
                0286            DO k=facet_owns(1,j),tS
                0287             W2_cumSum_tiles(1,k,tN) = 1
                0288            ENDDO
                0289          ENDDO
                0290         ENDDO
                0291         tN = facet_owns(1,j) - 1 + nbTx
                0292         facetXYSum(1,tN,j) = 1
                0293         DO k=facet_owns(1,j),tN-1
                0294           facetXYSum(1,k,j) = W2_cumSum_tiles(1,k,tN)
                0295         ENDDO
                0296 
                0297         DO bj=1,nbTy-1
                0298          DO bi=1,nbTx
                0299           tS = facet_owns(1,j) - 1 + bi
                0300           tN = tS + bj*nbTx
                0301           DO k=1,bj
                0302            l = tS + (k-1)*nbTx
                0303            W2_cumSum_tiles(2,l,tN) = 1
                0304           ENDDO
                0305          ENDDO
                0306         ENDDO
                0307         tN = facet_owns(1,j) + (nbTy-1)*nbTx
                0308         facetXYSum(2,tN,j) = 1
                0309         DO k=facet_owns(1,j),tN-1
                0310           facetXYSum(2,k,j) = W2_cumSum_tiles(2,k,tN)
                0311         ENDDO
                0312 
                0313        ENDIF
                0314       ENDDO
                0315 
                0316 C-    Then across facet:
                0317       DO j=1,nFacets
                0318        DO k=1,exch2_nTiles
                0319         DO i=1,nFacets
                0320         facet_CSum(1,k,j) = facet_CSum(1,k,j)
575f97366a Jean*0321      &                    + W2_cumSum_facet(1,i,j)*facetXYSum(1,k,i)
ca64f811ee Jean*0322         facet_CSum(2,k,j) = facet_CSum(2,k,j)
575f97366a Jean*0323      &                    + W2_cumSum_facet(2,i,j)*facetXYSum(2,k,i)
ca64f811ee Jean*0324         ENDDO
                0325        ENDDO
                0326       ENDDO
                0327 
                0328 C-    Finally, account for cumulated sum at facet origin:
                0329       DO j=1,nFacets
                0330        DO tN=facet_owns(1,j),facet_owns(2,j)
                0331         DO k=1,exch2_nTiles
                0332          W2_cumSum_tiles(1,k,tN) = W2_cumSum_tiles(1,k,tN)
                0333      &                           + facet_CSum(1,k,j)
                0334          W2_cumSum_tiles(2,k,tN) = W2_cumSum_tiles(2,k,tN)
                0335      &                           + facet_CSum(2,k,j)
                0336         ENDDO
                0337        ENDDO
                0338       ENDDO
                0339 
                0340       IF (prtFlag) THEN
b9dadda204 Mart*0341        WRITE(W2_oUnit,'(A,I8,A)')
575f97366a Jean*0342      &    ' Tile Matrix for CUMUL-SUM (nTiles=', exch2_nTiles, ' ):'
ca64f811ee Jean*0343        DO j=1,exch2_nTiles
                0344         DO is=1,exch2_nTiles,10
                0345          ie = MIN(is+9,exch2_nTiles)
                0346          IF ( is.EQ.1 ) THEN
b9dadda204 Mart*0347           WRITE(W2_oUnit,'(3(I8,A),10(2I3,A))') j,' ,',is,' ->',ie,' :',
ca64f811ee Jean*0348      &     (W2_cumSum_tiles(1,i,j),W2_cumSum_tiles(2,i,j),' ,',i=is,ie)
                0349          ELSE
b9dadda204 Mart*0350           WRITE(W2_oUnit,'(8X,2(I8,A),10(2I3,A))') is,' ->',ie,' :',
ca64f811ee Jean*0351      &     (W2_cumSum_tiles(1,i,j),W2_cumSum_tiles(2,i,j),' ,',i=is,ie)
                0352          ENDIF
                0353         ENDDO
                0354        ENDDO
                0355       ENDIF
575f97366a Jean*0356 
ca64f811ee Jean*0357       WRITE(msgBuf,'(2A)') 'W2_SET_MAP_CUMSUM: ',
                0358      &  'setting Tile Matrix for CUMUL-SUM : done'
                0359       CALL PRINT_MESSAGE( msgBuf, W2_oUnit, SQUEEZE_RIGHT, myThid )
                0360 
575f97366a Jean*0361 #else /* W2_CUMSUM_USE_MATRIX */
                0362       WRITE(msgBuf,'(2A)') 'W2_SET_MAP_CUMSUM: ',
                0363      &  'done (skip Tile Matrix setting)'
                0364       CALL PRINT_MESSAGE( msgBuf, W2_oUnit, SQUEEZE_RIGHT, myThid )
                0365 #endif /* W2_CUMSUM_USE_MATRIX */
                0366 
ca64f811ee Jean*0367       RETURN
                0368       END