pkg/seaice/seaice_itd_redist.F

86b84a92fc Patr*0001 #include "SEAICE_OPTIONS.h"
                0002
e7af59f6fd Jean*0003 CBOP
86b84a92fc Patr*0004 C !ROUTINE: SEAICE_ITD_REDIST
                0005
                0006 C !INTERFACE: ==========================================================
                0007       SUBROUTINE SEAICE_ITD_REDIST(
                0008      I     bi, bj, myTime, myIter, myThid )
                0009
                0010 C !DESCRIPTION: \bv
                0011 C     *===========================================================*
                0012 C     | SUBROUTINE SEAICE_ITD_REDIST
b8f04b8c26 Jean*0013 C     | o checks if absolute ice thickness in any category
86b84a92fc Patr*0014 C     |   exceeds its category limits
b8f04b8c26 Jean*0015 C     | o redistributes sea ice area and volume
86b84a92fc Patr*0016 C     |   and associated ice properties in thickness space
                0017 C     |
                0018 C     | Torge Martin, Feb. 2012, torge@mit.edu
                0019 C     *===========================================================*
                0020 C \ev
                0021
                0022 C !USES: ===============================================================
                0023       IMPLICIT NONE
                0024
                0025 C     === Global variables to be checked and redistributed ===
                0026 C     AREAITD   :: sea ice area      by category
                0027 C     HEFFITD   :: sea ice thickness by category
                0028 C
                0029 C     === Global variables to be redistributed ===
                0030 C     HSNOWITD  :: snow thickness    by category
                0031 C     enthalpy ?
                0032 C     temperature ?
                0033 C     salinity ?
                0034 C     age ?
                0035 C
                0036 #include "SIZE.h"
                0037 #include "EEPARAMS.h"
e7af59f6fd Jean*0038 c#include "PARAMS.h"
86b84a92fc Patr*0039 #include "SEAICE_SIZE.h"
                0040 #include "SEAICE_PARAMS.h"
                0041 #include "SEAICE.h"
                0042
                0043 C !INPUT PARAMETERS: ===================================================
                0044 C     === Routine arguments ===
                0045 C     bi, bj    :: outer loop counters
                0046 C     myTime    :: current time
                0047 C     myIter    :: iteration number
                0048 C     myThid    :: Thread no. that called this routine.
                0049       _RL myTime
                0050       INTEGER bi,bj
                0051       INTEGER myIter
                0052       INTEGER myThid
                0053
                0054 #ifdef SEAICE_ITD
                0055
                0056 C !LOCAL VARIABLES: ====================================================
                0057 C     === Local variables ===
                0058 C     i,j,k       :: inner loop counters
                0059 C     nITD        :: number of sea ice thickness categories
                0060 C     openwater   :: open water area fraction
                0061 C
                0062       INTEGER i, j, k
524979f562 Mart*0063       _RL openWater  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
e7af59f6fd Jean*0064 CEOP
86b84a92fc Patr*0065 C---+-|--1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0066
                0067 c      DO bj=myByLo(myThid),myByHi(myThid)
                0068 c       DO bi=myBxLo(myThid),myBxHi(myThid)
                0069 C must now be called within bi,bj loop
                0070
                0071 C       calculate area of open water
524979f562 Mart*0072       DO j=1-OLy,sNy+OLy
                0073        DO i=1-OLx,sNx+OLx
                0074         openWater(i,j) = ONE
                0075        ENDDO
                0076       ENDDO
                0077       DO k=1,nITD
                0078        DO j=1-OLy,sNy+OLy
                0079         DO i=1-OLx,sNx+OLx
                0080          openWater(i,j) = openWater(i,j) - AREAITD(i,j,k,bi,bj)
86b84a92fc Patr*0081         ENDDO
524979f562 Mart*0082        ENDDO
                0083       ENDDO
e7af59f6fd Jean*0084
86b84a92fc Patr*0085 C     ----------------------------------------------------
b8f04b8c26 Jean*0086 C     | redistribute/"advect" sea ice in thickness space |
86b84a92fc Patr*0087 C     | as described in Bitz et al. (2001)               |
                0088 C     ----------------------------------------------------
                0089
                0090 C---+-|--1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0091
524979f562 Mart*0092 C--   Hibler-type "ridging", i.e. cut back excessive ice area fraction ---
                0093 C     in case ice concentration exceeds 100% assume that
                0094 C     convergence of floe field has eliminated all open water
                0095 C     and eventual rafting occured in thinnest category:
                0096       DO j=1-OLy,sNy+OLy
                0097        DO i=1-OLx,sNx+OLx
                0098         IF (openWater(i,j) .lt. 0.0)
                0099      &     AREAITD(i,j,1,bi,bj) = openWater(i,j) + AREAITD(i,j,1,bi,bj)
                0100        ENDDO
                0101       ENDDO
86b84a92fc Patr*0102 C
e7af59f6fd Jean*0103 C     the following steps only make sense if there are actually
                0104 C     multi-categories
524979f562 Mart*0105       IF (nITD .gt. 1) THEN
86b84a92fc Patr*0106 C
524979f562 Mart*0107 C--   check if more thicker ice needs to be rafted to accomodate area excess:
                0108        DO k=1,nITD-1
                0109         DO j=1-OLy,sNy+OLy
                0110          DO i=1-OLx,sNx+OLx
                0111           IF (AREAITD(i,j,k,bi,bj) .lt. 0.0) THEN
                0112 C--   pass concentration deficit up to next thicker category
                0113 C--   since all quantities are extensive, we add instead of average
                0114            AREAITD (i,j,k+1,bi,bj) = AREAITD (i,j,k+1,bi,bj)
                0115      &                             + AREAITD (i,j,k,bi,bj)
                0116            AREAITD (i,j,k  ,bi,bj) = ZERO
                0117            HEFFITD (i,j,k+1,bi,bj) = HEFFITD (i,j,k+1,bi,bj)
                0118      &                             + HEFFITD (i,j,k,bi,bj)
                0119            HEFFITD (i,j,k  ,bi,bj) = ZERO
                0120            HSNOWITD(i,j,k+1,bi,bj) = HSNOWITD(i,j,k+1,bi,bj)
                0121      &                             + HSNOWITD(i,j,k,bi,bj)
                0122            HSNOWITD(i,j,k  ,bi,bj) = ZERO
86b84a92fc Patr*0123 C            t1(k+1) = t1(k+1)+t1(k); t1(k) = ZERO
                0124 C            t2(k+1) = t2(k+1)+t2(k); t2(k) = ZERO
                0125 C            age(k+1)=age(k+1)+age(k);age(k)= ZERO
                0126 C this is for ridged sea ice volume fraction
                0127 C            IF (PRESENT(rdg)) THEN
b8f04b8c26 Jean*0128 C             rdg(k+1)=rdg(k+1)+rdg(k); rdg(k)= ZERO
86b84a92fc Patr*0129 C            ENDIF
524979f562 Mart*0130           ENDIF
86b84a92fc Patr*0131          ENDDO
                0132         ENDDO
524979f562 Mart*0133        ENDDO
86b84a92fc Patr*0134
                0135 C     --- ice thickness redistribution ---
                0136 C         now check that ice thickness stays within category limits
524979f562 Mart*0137        DO k=1,nITD-1
                0138         DO j=1-OLy,sNy+OLy
                0139          DO i=1-OLx,sNx+OLx
                0140           IF (HEFFITD(i,j,k,bi,bj) .gt.
86b84a92fc Patr*0141      &         Hlimit(k)*AREAITD(i,j,k,bi,bj)) THEN
e7af59f6fd Jean*0142 C--   the upper thickness limit is exceeded: move ice up to next
524979f562 Mart*0143 C     thicker category
                0144            AREAITD (i,j,k+1,bi,bj) = AREAITD (i,j,k+1,bi,bj)
                0145      &                             + AREAITD (i,j,k,bi,bj)
                0146            AREAITD (i,j,k  ,bi,bj) = ZERO
                0147            HEFFITD (i,j,k+1,bi,bj) = HEFFITD (i,j,k+1,bi,bj)
                0148      &                             + HEFFITD (i,j,k,bi,bj)
                0149            HEFFITD (i,j,k  ,bi,bj) = ZERO
                0150            HSNOWITD(i,j,k+1,bi,bj) = HSNOWITD(i,j,k+1,bi,bj)
                0151      &                             + HSNOWITD(i,j,k,bi,bj)
                0152            HSNOWITD(i,j,k  ,bi,bj) = ZERO
86b84a92fc Patr*0153 C            t1(k+1) = t1(k+1)+t1(k); t1(k) = ZERO
                0154 C            t2(k+1) = t2(k+1)+t2(k); t2(k) = ZERO
                0155 C            age(k+1)=age(k+1)+age(k);age(k)= ZERO
                0156 C            IF (PRESENT(rdg)) THEN
b8f04b8c26 Jean*0157 C             rdg(k+1)=rdg(k+1)+rdg(k);rdg(k)= ZERO
86b84a92fc Patr*0158 C            ENDIF
524979f562 Mart*0159           ENDIF
86b84a92fc Patr*0160          ENDDO
                0161         ENDDO
524979f562 Mart*0162        ENDDO
e7af59f6fd Jean*0163 C
524979f562 Mart*0164        DO k=nITD,2,-1
                0165         DO j=1-OLy,sNy+OLy
                0166          DO i=1-OLx,sNx+OLx
                0167           IF (HEFFITD(i,j,k,bi,bj) .lt.
86b84a92fc Patr*0168      &         Hlimit(k-1)*AREAITD(i,j,k,bi,bj)) THEN
e7af59f6fd Jean*0169 C--   the lower thickness limit is exceeded: move ice down to next thinner
524979f562 Mart*0170 C     category
                0171            AREAITD (i,j,k-1,bi,bj) = AREAITD (i,j,k-1,bi,bj)
                0172      &                             + AREAITD (i,j,k,bi,bj)
                0173            AREAITD (i,j,k  ,bi,bj) = ZERO
                0174            HEFFITD (i,j,k-1,bi,bj) = HEFFITD (i,j,k-1,bi,bj)
                0175      &                             + HEFFITD (i,j,k,bi,bj)
                0176            HEFFITD (i,j,k  ,bi,bj) = ZERO
                0177            HSNOWITD(i,j,k-1,bi,bj) = HSNOWITD(i,j,k-1,bi,bj)
                0178      &                             + HSNOWITD(i,j,k,bi,bj)
                0179            HSNOWITD(i,j,k  ,bi,bj) = ZERO
86b84a92fc Patr*0180 c            snow(k-1) = snow(k-1)+snow(k); snow(k) = ZERO
                0181 C            t1(k-1) = t1(k-1)+t1(k); t1(k) = ZERO
                0182 C            t2(k-1) = t2(k-1)+t2(k); t2(k) = ZERO
                0183 C            age(k-1)=age(k-1)+age(k);age(k)= ZERO
                0184 C            IF (PRESENT(rdg)) THEN
b8f04b8c26 Jean*0185 C            rdg(k-1)=rdg(k-1)+rdg(k);rdg(k)= ZERO
86b84a92fc Patr*0186 C            ENDIF
524979f562 Mart*0187           ENDIF
86b84a92fc Patr*0188          ENDDO
                0189         ENDDO
524979f562 Mart*0190        ENDDO
86b84a92fc Patr*0191 C
524979f562 Mart*0192 C     end nITD>1 constraint
                0193       ENDIF
86b84a92fc Patr*0194
                0195 C     end bi,bj loop
                0196 c       ENDDO
                0197 c      ENDDO
                0198
                0199 C---+-|--1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0200 #endif /* SEAICE_ITD */
                0201       RETURN
b8f04b8c26 Jean*0202       END