C $Header: /u/gcmpack/MITgcm/pkg/seaice/seaice_itd_redist.F,v 1.4 2014/10/20 03:20:57 gforget Exp $ C $Name: $ #include "SEAICE_OPTIONS.h" #ifdef ALLOW_AUTODIFF # include "AUTODIFF_OPTIONS.h" #endif C !ROUTINE: SEAICE_ITD_REDIST C !INTERFACE: ========================================================== SUBROUTINE SEAICE_ITD_REDIST( I bi, bj, myTime, myIter, myThid ) C !DESCRIPTION: \bv C *===========================================================* C | SUBROUTINE SEAICE_ITD_REDIST C | o checks if absolute ice thickness in any category C | exceeds its category limits C | o redistributes sea ice area and volume C | and associated ice properties in thickness space C | C | Torge Martin, Feb. 2012, torge@mit.edu C *===========================================================* C \ev C !USES: =============================================================== IMPLICIT NONE C === Global variables to be checked and redistributed === C AREAITD :: sea ice area by category C HEFFITD :: sea ice thickness by category C C === Global variables to be redistributed === C HSNOWITD :: snow thickness by category C enthalpy ? C temperature ? C salinity ? C age ? C #include "SIZE.h" #include "EEPARAMS.h" #include "PARAMS.h" #include "GRID.h" #include "SEAICE_SIZE.h" #include "SEAICE_PARAMS.h" #include "SEAICE.h" #ifdef ALLOW_AUTODIFF_TAMC # include "tamc.h" #endif C !INPUT PARAMETERS: =================================================== C === Routine arguments === C bi, bj :: outer loop counters C myTime :: current time C myIter :: iteration number C myThid :: Thread no. that called this routine. _RL myTime INTEGER bi,bj INTEGER myIter INTEGER myThid CEndOfInterface #ifdef SEAICE_ITD C !LOCAL VARIABLES: ==================================================== C === Local variables === C i,j,k :: inner loop counters C nITD :: number of sea ice thickness categories C openwater :: open water area fraction C INTEGER i, j, k #ifdef ALLOW_AUTODIFF_TAMC INTEGER itmpkey #endif /* ALLOW_AUTODIFF_TAMC */ #ifdef SEAICE_AGE INTEGER iTracer #endif _RL openWater (1-OLx:sNx+OLx,1-OLy:sNy+OLy) C---+-|--1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| c DO bj=myByLo(myThid),myByHi(myThid) c DO bi=myBxLo(myThid),myBxHi(myThid) C must now be called within bi,bj loop C calculate area of open water DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx openWater(i,j) = ONE ENDDO ENDDO DO k=1,nITD DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx openWater(i,j) = openWater(i,j) - AREAITD(i,j,k,bi,bj) ENDDO ENDDO ENDDO C ---------------------------------------------------- C | redistribute/"advect" sea ice in thickness space | C | as described in Bitz et al. (2001) | C ---------------------------------------------------- C---+-|--1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| C-- Hibler-type "ridging", i.e. cut back excessive ice area fraction --- C in case ice concentration exceeds 100% assume that C convergence of floe field has eliminated all open water C and eventual rafting occured in thinnest category: DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx IF (openWater(i,j) .lt. 0.0) & AREAITD(i,j,1,bi,bj) = openWater(i,j) + AREAITD(i,j,1,bi,bj) ENDDO ENDDO C C the following steps only make sense if there are actually C multi-categories IF (nITD .gt. 1) THEN C C-- check if more thicker ice needs to be rafted to accomodate area excess: DO k=1,nITD-1 DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx IF (AREAITD(i,j,k,bi,bj) .lt. 0.0) THEN C-- pass concentration deficit up to next thicker category C-- since all quantities are extensive, we add instead of average AREAITD (i,j,k+1,bi,bj) = AREAITD (i,j,k+1,bi,bj) & + AREAITD (i,j,k,bi,bj) AREAITD (i,j,k ,bi,bj) = ZERO HEFFITD (i,j,k+1,bi,bj) = HEFFITD (i,j,k+1,bi,bj) & + HEFFITD (i,j,k,bi,bj) HEFFITD (i,j,k ,bi,bj) = ZERO HSNOWITD(i,j,k+1,bi,bj) = HSNOWITD(i,j,k+1,bi,bj) & + HSNOWITD(i,j,k,bi,bj) HSNOWITD(i,j,k ,bi,bj) = ZERO C t1(k+1) = t1(k+1)+t1(k); t1(k) = ZERO C t2(k+1) = t2(k+1)+t2(k); t2(k) = ZERO C age(k+1)=age(k+1)+age(k);age(k)= ZERO C this is for ridged sea ice volume fraction C IF (PRESENT(rdg)) THEN C rdg(k+1)=rdg(k+1)+rdg(k); rdg(k)= ZERO C ENDIF ENDIF ENDDO ENDDO ENDDO C --- ice thickness redistribution --- C now check that ice thickness stays within category limits DO k=1,nITD-1 DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx IF (HEFFITD(i,j,k,bi,bj) .gt. & Hlimit(k)*AREAITD(i,j,k,bi,bj)) THEN C-- the upper thickness limit is exceeded: move ice up to next C thicker category AREAITD (i,j,k+1,bi,bj) = AREAITD (i,j,k+1,bi,bj) & + AREAITD (i,j,k,bi,bj) AREAITD (i,j,k ,bi,bj) = ZERO HEFFITD (i,j,k+1,bi,bj) = HEFFITD (i,j,k+1,bi,bj) & + HEFFITD (i,j,k,bi,bj) HEFFITD (i,j,k ,bi,bj) = ZERO HSNOWITD(i,j,k+1,bi,bj) = HSNOWITD(i,j,k+1,bi,bj) & + HSNOWITD(i,j,k,bi,bj) HSNOWITD(i,j,k ,bi,bj) = ZERO C t1(k+1) = t1(k+1)+t1(k); t1(k) = ZERO C t2(k+1) = t2(k+1)+t2(k); t2(k) = ZERO C age(k+1)=age(k+1)+age(k);age(k)= ZERO C IF (PRESENT(rdg)) THEN C rdg(k+1)=rdg(k+1)+rdg(k);rdg(k)= ZERO C ENDIF ENDIF ENDDO ENDDO ENDDO C DO k=nITD,2,-1 DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx IF (HEFFITD(i,j,k,bi,bj) .lt. & Hlimit(k-1)*AREAITD(i,j,k,bi,bj)) THEN C-- the lower thickness limit is exceeded: move ice down to next thinner C category AREAITD (i,j,k-1,bi,bj) = AREAITD (i,j,k-1,bi,bj) & + AREAITD (i,j,k,bi,bj) AREAITD (i,j,k ,bi,bj) = ZERO HEFFITD (i,j,k-1,bi,bj) = HEFFITD (i,j,k-1,bi,bj) & + HEFFITD (i,j,k,bi,bj) HEFFITD (i,j,k ,bi,bj) = ZERO HSNOWITD(i,j,k-1,bi,bj) = HSNOWITD(i,j,k-1,bi,bj) & + HSNOWITD(i,j,k,bi,bj) HSNOWITD(i,j,k ,bi,bj) = ZERO c snow(k-1) = snow(k-1)+snow(k); snow(k) = ZERO C t1(k-1) = t1(k-1)+t1(k); t1(k) = ZERO C t2(k-1) = t2(k-1)+t2(k); t2(k) = ZERO C age(k-1)=age(k-1)+age(k);age(k)= ZERO C IF (PRESENT(rdg)) THEN C rdg(k-1)=rdg(k-1)+rdg(k);rdg(k)= ZERO C ENDIF ENDIF ENDDO ENDDO ENDDO C C end nITD>1 constraint ENDIF C end bi,bj loop c ENDDO c ENDDO C---+-|--1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----| #endif /* SEAICE_ITD */ RETURN END