Back to home page

MITgcm
 
 

    

File indexing completed on 2020-02-28 06:11:23 UTC

view on githubraw file Latest commit 3b867959 on 2020-02-11 01:31:16 UTC
9952f046d7 dngo*0001 #include "SHELFICE_OPTIONS.h"
                0002 
                0003 CBOP
                0004 C     !ROUTINE: SHELFICE_REMESH_C_MASK
                0005 C     !INTERFACE:
                0006       SUBROUTINE SHELFICE_REMESH_C_MASK(
7155c3b3df Jean*0007      O                    k1SurfC, mrgFacC,
9952f046d7 dngo*0008      I                    myTime, myIter, myThid )
                0009 C     !DESCRIPTION: \bv
                0010 C     *==========================================================*
                0011 C     | SUBROUTINE SHELFICE_REMESH_C_MASK
                0012 C     | o Loops through top level cells and determines those where
                0013 C     |   hFac is too large and hence splits into two cells,
                0014 C     |   and cells where hFac is too small, and merges cell with
                0015 C     |   below.
7155c3b3df Jean*0016 C     | o Update surface ref position Ro_surf and anomaly Eta
3b86795949 Jean*0017 C     |   as well as R_shelfIce.
7155c3b3df Jean*0018 C     | o Recomputes h0FacC using formula from initialisation
3b86795949 Jean*0019 C     |   based on new R_shelfIce
9952f046d7 dngo*0020 C     *==========================================================*
                0021 C     \ev
                0022 
                0023 C     !USES:
                0024       IMPLICIT NONE
                0025 C     === Global variables ===
                0026 #include "SIZE.h"
                0027 #include "EEPARAMS.h"
                0028 #include "PARAMS.h"
                0029 #include "GRID.h"
                0030 #include "SURFACE.h"
                0031 #include "DYNVARS.h"
                0032 #include "SHELFICE.h"
                0033 
                0034 C     !INPUT/OUTPUT PARAMETERS:
7155c3b3df Jean*0035 C     k1SurfC   :: surface level index (at cell center) before remeshing
                0036 C     mrgFacC   :: merging weight for tracers (cell center)
9952f046d7 dngo*0037 C     myTime    :: Current time in simulation
                0038 C     myIter    :: Current iteration number
                0039 C     myThid    :: my Thread Id number
7155c3b3df Jean*0040       INTEGER k1SurfC(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0041       _RL mrgFacC(1-OLx:sNx+OLx,1-OLy:sNy+OLy,2,nSx,nSy)
9952f046d7 dngo*0042       _RL myTime
                0043       INTEGER myIter
                0044       INTEGER myThid
                0045 CEOP
                0046 
                0047 #ifdef ALLOW_SHELFICE_REMESHING
7155c3b3df Jean*0048 #ifdef NONLIN_FRSURF
9952f046d7 dngo*0049 C     !LOCAL VARIABLES:
                0050 C     bi,bj   :: tile indices
                0051 C     i,j,k   :: Loop counters
7155c3b3df Jean*0052 C     k1, k2  :: previous and new surface level
9952f046d7 dngo*0053 C     ioUnit  :: temp for writing msg unit
                0054 C     msgBuf  :: Informational/error message buffer
                0055       INTEGER bi, bj
7155c3b3df Jean*0056       INTEGER i, j
                0057       INTEGER ks, k1, k2
9952f046d7 dngo*0058 #ifdef SHELFICE_REMESH_PRINT
7155c3b3df Jean*0059       LOGICAL prtFirst, prtPoint
9952f046d7 dngo*0060       INTEGER ioUnit
                0061       CHARACTER*(MAX_LEN_MBUF) msgBuf
                0062 #endif
                0063       _RL sync_fac, stag_fac
7155c3b3df Jean*0064       _RL adjust, r_newDz, tmpVar
9952f046d7 dngo*0065 
                0066 #ifdef SHELFICE_REMESH_PRINT
                0067       prtFirst = .TRUE.
                0068       ioUnit = -1
                0069       IF ( debugLevel.GE.debLevB ) ioUnit = standardMessageUnit
                0070 #endif
                0071       IF ( staggerTimestep ) THEN
                0072         sync_fac = 0.0
                0073         stag_fac = 1.0
                0074       ELSE
                0075         sync_fac = 1.0
                0076         stag_fac = 0.0
                0077       ENDIF
                0078 
                0079       DO bj = myByLo(myThid), myByHi(myThid)
                0080        DO bi = myBxLo(myThid), myBxHi(myThid)
                0081 
7155c3b3df Jean*0082         DO j = 1-OLy,sNy+OLy
                0083          DO i = 1-OLx,sNx+OLx
                0084           k1SurfC(i,j,bi,bj) = kSurfC(i,j,bi,bj)
                0085           mrgFacC(i,j,1,bi,bj) = 0.
                0086           mrgFacC(i,j,2,bi,bj) = 0.
                0087           ks = kTopC(i,j,bi,bj)
                0088 
                0089           IF ( ks.NE.0 ) THEN
9952f046d7 dngo*0090 C-- SPLIT CELLS
7155c3b3df Jean*0091            IF ( ks.GT.1 .AND.
9952f046d7 dngo*0092      &          hFac_surfC(i,j,bi,bj).GT.SHELFICEsplitThreshold ) THEN
7155c3b3df Jean*0093             IF ( (hFac_surfC(i,j,bi,bj)-1)*drF(ks)*recip_drF(ks-1)
                0094      &           .GT.SHELFICEmergeThreshold ) THEN
9952f046d7 dngo*0095 
7155c3b3df Jean*0096              k1 = ks
                0097              k2 = ks-1
                0098              kSurfC(i,j,bi,bj) = k2
3b86795949 Jean*0099              adjust = rF(k2) - R_shelfIce(i,j,bi,bj)
9952f046d7 dngo*0100 #ifdef SHELFICE_REMESH_PRINT
7155c3b3df Jean*0101              prtPoint = ( ioUnit.GE.0 ) .AND.
                0102      &          ( i.GE.1 .AND. i.LE.sNx .AND. j.GE.1 .AND. j.LE.sNy )
                0103              IF ( prtPoint ) THEN
9952f046d7 dngo*0104               IF ( prtFirst ) THEN
                0105                WRITE(msgBuf,'(A,I10)') 'SHI_REMESH at it=', myIter
                0106                CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0107                prtFirst = .FALSE.
                0108               ENDIF
                0109               WRITE(msgBuf,'(A,2I5,2I4,A,1P2E12.3)') '--> REMESH in:',
                0110      &            i, j, bi, bj, ' , x,y=', XC(i,j,bi,bj),YC(i,j,bi,bj)
                0111               CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0112               WRITE(msgBuf,'(2A,I4,3(A,1P1E12.3))') ' before:',
7155c3b3df Jean*0113      &        '  ks=', k1, ' Ro_s=', Ro_surf(i,j,bi,bj),
9952f046d7 dngo*0114      &        ' eta=', etaH(i,j,bi,bj), ' hFac=', hFac_surfC(i,j,bi,bj)
                0115               CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0116              ENDIF
7155c3b3df Jean*0117 #endif /* SHELFICE_REMESH_PRINT */
3b86795949 Jean*0118 C-    decrement eta and increment R_shelfIce and Ro_surf
9952f046d7 dngo*0119              etaN(i,j,bi,bj) = etaN(i,j,bi,bj)- adjust
                0120              etaH(i,j,bi,bj) = etaH(i,j,bi,bj)- adjust
                0121              etaHnm1(i,j,bi,bj) = etaHnm1(i,j,bi,bj)- adjust
                0122              R_shelfIce(i,j,bi,bj) = R_shelfIce(i,j,bi,bj)+adjust
                0123              Ro_surf(i,j,bi,bj) = Ro_surf(i,j,bi,bj)+adjust
7155c3b3df Jean*0124 C     update cell-centered grid factors and mask: maskC, h0FacC, recip_hFacC
                0125              maskC(i,j,k2,bi,bj)  = oneRS
                0126              h0FacC(i,j,k1,bi,bj) = oneRS
                0127              h0FacC(i,j,k2,bi,bj) = oneRS
                0128              recip_hFacC(i,j,k1,bi,bj) = oneRS
                0129              recip_hFacC(i,j,k2,bi,bj) = oneRS
                0130              IF ( k1 .EQ. kLowC(i,j,bi,bj) ) THEN
                0131                h0FacC(i,j,k1,bi,bj) = ( rF(k1) - R_low(i,j,bi,bj) )
                0132      &                              *recip_drF(k1)
                0133                recip_hFacC(i,j,k1,bi,bj) = oneRS / h0FacC(i,j,k1,bi,bj)
9952f046d7 dngo*0134              ENDIF
                0135 #ifdef SHELFICE_REMESH_PRINT
7155c3b3df Jean*0136              IF ( prtPoint ) THEN
9952f046d7 dngo*0137               WRITE(msgBuf,'(2A,I4,3(A,1P1E12.3))') ' after :',
7155c3b3df Jean*0138      &        '  ks=', k2, ' Ro_s=', Ro_surf(i,j,bi,bj),
9952f046d7 dngo*0139      &        ' eta=', etaH(i,j,bi,bj), ' hFac=',
7155c3b3df Jean*0140      &           h0FacC(i,j,k2,bi,bj)+etaH(i,j,bi,bj)*recip_drF(k2)
9952f046d7 dngo*0141               CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0142              ENDIF
7155c3b3df Jean*0143 #endif /* SHELFICE_REMESH_PRINT */
9952f046d7 dngo*0144             ENDIF
                0145            ENDIF
                0146           ENDIF
                0147 
                0148 C-- MERGE CELLS
7155c3b3df Jean*0149           IF ( ks.NE.0 .AND. ks.LT.kLowC (i,j,bi,bj) ) THEN
                0150            IF ( hFac_surfC(i,j,bi,bj).LT.SHELFICEmergeThreshold ) THEN
                0151             IF ( (hFac_surfC(i,j,bi,bj)*drF(ks)*recip_drF(ks+1)+1)
                0152      &           .LT.SHELFICEsplitThreshold ) THEN
                0153 
                0154              k1 = ks
                0155              k2 = ks+1
                0156              kSurfC(i,j,bi,bj) = k2
3b86795949 Jean*0157              adjust = R_shelfIce(i,j,bi,bj)-rF(k2)
9952f046d7 dngo*0158 #ifdef SHELFICE_REMESH_PRINT
7155c3b3df Jean*0159              prtPoint = ( ioUnit.GE.0 ) .AND.
                0160      &          ( i.GE.1 .AND. i.LE.sNx .AND. j.GE.1 .AND. j.LE.sNy )
                0161              IF ( prtPoint ) THEN
9952f046d7 dngo*0162               IF ( prtFirst ) THEN
                0163                WRITE(msgBuf,'(A,I10)') 'SHI_REMESH at it=', myIter
                0164                CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0165                prtFirst = .FALSE.
                0166               ENDIF
                0167               WRITE(msgBuf,'(A,2I5,2I4,A,1P2E12.3)') '--> REMESH in:',
                0168      &            i, j, bi, bj, ' , x,y=', XC(i,j,bi,bj),YC(i,j,bi,bj)
                0169               CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0170               WRITE(msgBuf,'(2A,I4,3(A,1P1E12.3))') ' before:',
7155c3b3df Jean*0171      &        '  ks=', k1, ' Ro_s=', Ro_surf(i,j,bi,bj),
9952f046d7 dngo*0172      &        ' eta=', etaH(i,j,bi,bj), ' hFac=', hFac_surfC(i,j,bi,bj)
                0173               CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0174              ENDIF
7155c3b3df Jean*0175 #endif /* SHELFICE_REMESH_PRINT */
3b86795949 Jean*0176 C-    increment eta and decrement R_shelfIce and Ro_surf
9952f046d7 dngo*0177              etaN(i,j,bi,bj) = etaN(i,j,bi,bj) + adjust
                0178              etaH(i,j,bi,bj) = etaH(i,j,bi,bj) + adjust
                0179              etaHnm1(i,j,bi,bj) = etaHnm1(i,j,bi,bj) + adjust
3b86795949 Jean*0180              R_shelfIce(i,j,bi,bj) = R_shelfIce(i,j,bi,bj)-adjust
9952f046d7 dngo*0181              Ro_surf(i,j,bi,bj) = Ro_surf(i,j,bi,bj) - adjust
                0182 
7155c3b3df Jean*0183 C-    Compute merging weights for tracer in new (=k2) top cell, accounting
                0184 C     for former (=k1) top cell content:
                0185 C      Tr(k2) <-- mrgFacC(1)*Tr(k1) + mrgFacC(2)*Tr(k2) ; Tr(k1) <-- 0.
                0186 C     first find the inverse thickness of the new top cell:
                0187              IF ( k2 .LT. kLowC(i,j,bi,bj) ) THEN
                0188               r_newDz = oneRL / ( drF(k2) + etaH(i,j,bi,bj) )
9952f046d7 dngo*0189              ELSE
                0190               r_newDz = oneRL
7155c3b3df Jean*0191      &                / ( rF(k2)-r_low(i,j,bi,bj) + etaH(i,j,bi,bj) )
9952f046d7 dngo*0192              ENDIF
7155c3b3df Jean*0193              mrgFacC(i,j,1,bi,bj) = ( stag_Fac*hFac_surfC(i,j,bi,bj)
                0194      &                              + sync_fac*hFacC(i,j,k1,bi,bj)
                0195      &                              )*drF(k1)*r_newDz
                0196              mrgFacC(i,j,2,bi,bj) = hFacC(i,j,k2,bi,bj)*drF(k2)*r_newDz
                0197 C     update cell-centered grid factors and mask: maskC, h0FacC, recip_hFacC
                0198              maskC(i,j,k1,bi,bj)  = zeroRS
                0199              h0FacC(i,j,k1,bi,bj) = zeroRS
                0200              h0FacC(i,j,k2,bi,bj) = oneRS
                0201              recip_hFacC(i,j,k1,bi,bj) = zeroRS
                0202              recip_hFacC(i,j,k2,bi,bj) = oneRS
                0203              IF ( k2 .EQ. kLowC(i,j,bi,bj) ) THEN
                0204                h0FacC(i,j,k2,bi,bj) = ( rF(k2) - R_low(i,j,bi,bj) )
                0205      &                              * recip_drF(k2)
                0206                recip_hFacC(i,j,k2,bi,bj) = oneRS / h0FacC(i,j,k2,bi,bj)
9952f046d7 dngo*0207              ENDIF
                0208 #ifdef SHELFICE_REMESH_PRINT
7155c3b3df Jean*0209              IF ( prtPoint ) THEN
9952f046d7 dngo*0210               WRITE(msgBuf,'(2A,I4,3(A,1P1E12.3))') ' after :',
7155c3b3df Jean*0211      &        '  ks=', k2, ' Ro_s=', Ro_surf(i,j,bi,bj),
9952f046d7 dngo*0212      &        ' eta=', etaH(i,j,bi,bj), ' hFac=',
7155c3b3df Jean*0213      &           h0FacC(i,j,k2,bi,bj)+etaH(i,j,bi,bj)*recip_drF(k2)
9952f046d7 dngo*0214               CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0215              ENDIF
7155c3b3df Jean*0216 #endif /* SHELFICE_REMESH_PRINT */
9952f046d7 dngo*0217             ENDIF
                0218            ENDIF
                0219           ENDIF
                0220 
7155c3b3df Jean*0221 C-  Update remaining cell-centered grid fields
                0222 C   Note: CALC_SURF_DR & UPDATE_SURF_DR will update hFac_surfC and surface hFacC
                0223           IF ( kSurfC(i,j,bi,bj).NE.k1SurfC(i,j,bi,bj) ) THEN
                0224             k1 = k1SurfC(i,j,bi,bj)
                0225             k2 =  kSurfC(i,j,bi,bj)
                0226             hFacC(i,j,k1,bi,bj) = h0FacC(i,j,k1,bi,bj)
                0227             hFacC(i,j,k2,bi,bj) = h0FacC(i,j,k2,bi,bj)
                0228             tmpVar = Ro_surf(i,j,bi,bj) - R_low(i,j,bi,bj)
                0229             recip_Rcol(i,j,bi,bj) = 1. _d 0 / tmpVar
                0230             kTopC(i,j,bi,bj)  = k2
                0231           ENDIF
                0232 
9952f046d7 dngo*0233          ENDDO
                0234         ENDDO
7155c3b3df Jean*0235 C---
9952f046d7 dngo*0236        ENDDO
                0237       ENDDO
                0238 #ifdef SHELFICE_REMESH_PRINT
                0239       IF ( ioUnit.GE.0 .AND. .NOT. prtFirst ) THEN
                0240         WRITE(msgBuf,'(A,A)') 'SHI_REMESH :', ' end of report'
                0241         CALL PRINT_MESSAGE( msgBuf,ioUnit,SQUEEZE_RIGHT,myThid )
                0242       ENDIF
7155c3b3df Jean*0243 #endif /* SHELFICE_REMESH_PRINT */
9952f046d7 dngo*0244 
                0245 #endif /* NONLIN_FRSURF */
7155c3b3df Jean*0246 #endif /* ALLOW_SHELFICE_REMESHING */
9952f046d7 dngo*0247       RETURN
                0248       END
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated from https://github.com/MITgcm/MITgcm by the 2.2.1-MITgcm-0.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!