/[MITgcm]/MITgcm/pkg/seaice/seaice_growth.F
ViewVC logotype

Diff of /MITgcm/pkg/seaice/seaice_growth.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.122 by gforget, Fri May 27 23:27:15 2011 UTC revision 1.123 by gforget, Sat May 28 22:39:20 2011 UTC
# Line 120  C conversion factors to go from precip ( Line 120  C conversion factors to go from precip (
120    
121  C ICE/SNOW stocks tendencies associated with the various melt/freeze processes  C ICE/SNOW stocks tendencies associated with the various melt/freeze processes
122        _RL d_AREAbyATM         (1:sNx,1:sNy)        _RL d_AREAbyATM         (1:sNx,1:sNy)
 #ifdef FENTY_AREA_EXPANSION_CONTRACTION  
123        _RL d_AREAbyOCN         (1:sNx,1:sNy)        _RL d_AREAbyOCN         (1:sNx,1:sNy)
124        _RL d_AREAbyICE         (1:sNx,1:sNy)        _RL d_AREAbyICE         (1:sNx,1:sNy)
 #endif  
125    
126  c     The change of mean ice thickness due to out-of-bounds values following  c     The change of mean ice thickness due to out-of-bounds values following
127  c     sea ice dyhnamics  c     sea ice dyhnamics
# Line 166  C     actual snow thickness Line 164  C     actual snow thickness
164        _RL hsnowActual         (1:sNx,1:sNy)        _RL hsnowActual         (1:sNx,1:sNy)
165  C     actual ice thickness (with lower limit only) Reciprocal  C     actual ice thickness (with lower limit only) Reciprocal
166        _RL recip_heffActual    (1:sNx,1:sNy)        _RL recip_heffActual    (1:sNx,1:sNy)
167    C     local value (=HO or HO_south)
168          _RL HO_loc
169    
170  C     AREA_PRE :: hold sea-ice fraction field before any seaice-thermo update  C     AREA_PRE :: hold sea-ice fraction field before any seaice-thermo update
171        _RL AREApreTH           (1:sNx,1:sNy)        _RL AREApreTH           (1:sNx,1:sNy)
# Line 311  C ===================== Line 311  C =====================
311            r_QbyOCN (I,J)             = 0.0 _d 0            r_QbyOCN (I,J)             = 0.0 _d 0
312    
313            d_AREAbyATM(I,J)           = 0.0 _d 0            d_AREAbyATM(I,J)           = 0.0 _d 0
 #ifdef FENTY_AREA_EXPANSION_CONTRACTION  
314            d_AREAbyICE(I,J)           = 0.0 _d 0            d_AREAbyICE(I,J)           = 0.0 _d 0
315            d_AREAbyOCN(I,J)           = 0.0 _d 0            d_AREAbyOCN(I,J)           = 0.0 _d 0
 #endif  
316    
317              d_HEFFbyNEG(I,J)           = 0.0 _d 0
318            d_HEFFbyOCNonICE(I,J)      = 0.0 _d 0            d_HEFFbyOCNonICE(I,J)      = 0.0 _d 0
319            d_HEFFbyATMonOCN(I,J)      = 0.0 _d 0            d_HEFFbyATMonOCN(I,J)      = 0.0 _d 0
320            d_HEFFbyFLOODING(I,J)      = 0.0 _d 0            d_HEFFbyFLOODING(I,J)      = 0.0 _d 0
# Line 323  C ===================== Line 322  C =====================
322            d_HEFFbyATMonOCN_open(I,J) = 0.0 _d 0            d_HEFFbyATMonOCN_open(I,J) = 0.0 _d 0
323            d_HEFFbyATMonOCN_cover(I,J) = 0.0 _d 0            d_HEFFbyATMonOCN_cover(I,J) = 0.0 _d 0
324    
325              d_HSNWbyNEG(I,J)           = 0.0 _d 0
326            d_HSNWbyATMonSNW(I,J)      = 0.0 _d 0            d_HSNWbyATMonSNW(I,J)      = 0.0 _d 0
327            d_HSNWbyOCNonSNW(I,J)      = 0.0 _d 0            d_HSNWbyOCNonSNW(I,J)      = 0.0 _d 0
328            d_HSNWbyRAIN(I,J)          = 0.0 _d 0            d_HSNWbyRAIN(I,J)          = 0.0 _d 0
# Line 843  c          The turbulent ocean-ice heat Line 843  c          The turbulent ocean-ice heat
843  c          of potential ice melt  c          of potential ice melt
844             a_QbyOCN(I,J) = a_QbyOCN(I,J) * convertQ2HI             a_QbyOCN(I,J) = a_QbyOCN(I,J) * convertQ2HI
845    
846  c          by design r_QbyOCN .LE. 0. so that initial ice growth cannot  c          by design a_QbyOCN .LE. 0. so that initial ice growth cannot
847  c          be triggered by this term, which Ian says is better for adjoint  c          be triggered by this term, which Ian says is better for adjoint
848  #else  #else
849    
# Line 1218  CADJ STORE AREA(:,:,bi,bj) = comlev1_bib Line 1218  CADJ STORE AREA(:,:,bi,bj) = comlev1_bib
1218          DO J=1,sNy          DO J=1,sNy
1219           DO I=1,sNx           DO I=1,sNx
1220    
1221  C         compute ice melt due to ATM (and OCN) heat stocks            IF ( YC(I,J,bi,bj) .LT. ZERO ) THEN
1222  #ifdef FENTY_AREA_EXPANSION_CONTRACTION             HO_loc=HO_south
1223              ELSE
1224  #ifdef ALLOW_DIAGNOSTICS             HO_loc=HO
           DIAGarrayA(I,J) = ZERO  
           DIAGarrayB(I,J) = ZERO  
           DIAGarrayC(I,J) = ZERO  
           DIAGarray(I,J) = ZERO  
 #endif  
   
 c         the various thickness tendency terms  
           tmpscal1 = d_HEFFbyATMOnOCN_open(I,J)  
           tmpscal2 = d_HEFFbyATMOnOCN_cover(I,J)  
           tmpscal3 = d_HEFFbyOCNonICE(I,J)  
   
 c         Part 1: Treat expansion/contraction of ice cover in open-water areas  
   
 C         All new ice cover is created from divergent air-sea heat fluxes.  Divergent  
 c         air-sea heat fluxes must exceed the potential convergent ocean-ice heat fluxes  
 c         for ice to form.    
           IF (tmpscal1 .GT. ZERO) then  
             IF ( YC(I,J,bi,bj) .LT. ZERO ) THEN  
              d_AREAbyATM(I,J)=tmpscal1/HO_south  
             ELSE  
              d_AREAbyATM(I,J)=tmpscal1/HO  
             ENDIF  
   
 c         If there are convergent air-sea heat fluxes and convergent air-sea  
 c         heat fluxes are permitted to melt ice, tmpscal1 < 0.  
           ELSEIF (AREApreTH(I,J) .GT. ZERO) THEN  
              d_AREAbyATM(i,J) = HALF*tmpscal1*recip_heffActual(I,J)  
1225            ENDIF            ENDIF
1226    
1227  c         Part 2: Reduce ice concentration from ice thinning from above  C         compute contraction/expansion from melt/growth
   
 c         Ice concentrations are reduced whenever existing ice thins from surface  
 c         heat flux convergence.  
           if (tmpscal2 .LE. ZERO) then  
              d_AREAbyICE(I,J) = HALF * tmpscal2 *recip_heffActual(I,J)  
           endif  
   
 c         Part 3: Reduce ice concentration from ice thinning from below  
   
 c         Sensible heat transfer from the ocean to the sea ice thins it and  
 c         reduces concentrations  
           if (tmpscal3 .LE.ZERO) then  
              d_AREAbyOCN(I,J) = HALF * tmpscal3 *recip_heffActual(I,J)  
           endif  
   
 #ifdef ALLOW_DIAGNOSTICS  
           DIAGarrayA(I,J) = d_AREAbyATM(I,J)  
           DIAGarrayB(I,J) = d_AREAbyICE(I,J)  
           DIAGarrayC(I,J) = d_AREAbyOCN(I,J)  
           DIAGarray(I,J) = d_AREAbyICE(I,J) + d_AREAbyATM(I,J)  
      &                     + d_AREAbyOCN(I,J)  
 #endif  
   
 #else /* FENTY_AREA_EXPANSION_CONTRACTION */  
   
1228  #ifdef SEAICE_GROWTH_LEGACY  #ifdef SEAICE_GROWTH_LEGACY
1229    
1230  C compute heff after ice melt by ocn:  C compute heff after ice melt by ocn:
1231            tmpscal0=HEFF(I,J,bi,bj)            tmpscal0=HEFF(I,J,bi,bj) - d_HEFFbyATMonOCN(I,J)
      &            - d_HEFFbyATMonOCN(I,J) - d_HEFFbyFLOODING(I,J)  
1232  C compute available heat left after snow melt by atm:  C compute available heat left after snow melt by atm:
1233            tmpscal1= a_QbyATM_open(I,J)+a_QbyATM_cover(I,J)            tmpscal1= a_QbyATM_open(I,J)+a_QbyATM_cover(I,J)
1234       &            - d_HSNWbyATMonSNW(I,J)/ICE2SNOW       &            - d_HSNWbyATMonSNW(I,J)/ICE2SNOW
# Line 1292  C (cannot melt more than all the ice) Line 1239  C (cannot melt more than all the ice)
1239            DIAGarray(I,J) = tmpscal2            DIAGarray(I,J) = tmpscal2
1240  #endif  #endif
1241  C gain of new ice over open water  C gain of new ice over open water
1242            tmpscal4 = MAX(ZERO,a_QbyATM_open(I,J))            tmpscal4=MAX(ZERO,a_QbyATM_open(I,J))
1243    C compute cover fraction tendency
1244              d_AREAbyATM(I,J)=tmpscal4/HO_loc+HALF*tmpscal3
1245         &         *AREApreTH(I,J) /(tmpscal0+.00001 _d 0)
1246    
1247  #else /* SEAICE_GROWTH_LEGACY */  #else /* SEAICE_GROWTH_LEGACY */
1248    
1249  # ifdef SEAICE_OCN_MELT_ACT_ON_AREA  #  ifdef FENTY_AREA_EXPANSION_CONTRACTION
1250    
1251    C         the various volume tendency terms
1252              tmpscal1 = d_HEFFbyATMOnOCN_open(I,J)
1253              tmpscal2 = d_HEFFbyATMOnOCN_cover(I,J)
1254              tmpscal3 = d_HEFFbyOCNonICE(I,J)
1255    
1256    C         Part 1: expansion/contraction of ice cover in open-water areas.
1257              IF (tmpscal1 .GT. ZERO) then
1258    C         Ice cover is all created from open ocean-water air-sea heat fluxes
1259                 d_AREAbyATM(I,J)=tmpscal1/HO_loc
1260              ELSEIF (AREApreTH(I,J) .GT. ZERO) THEN
1261    C         that can also act to remove ice cover.
1262                 d_AREAbyATM(i,J) = HALF*tmpscal1*recip_heffActual(I,J)
1263              ENDIF
1264    
1265    C         Part 2: contraction from ice thinning from above
1266              if (tmpscal2 .LE. ZERO) d_AREAbyICE(I,J) =
1267         &            HALF * tmpscal2 * recip_heffActual(I,J)
1268    
1269    C         Part 3: contraction from ice thinning from below
1270              if (tmpscal3 .LE.ZERO) d_AREAbyOCN(I,J) =
1271         &         HALF * tmpscal3* recip_heffActual(I,J)
1272    
1273    # else /* FENTY_AREA_EXPANSION_CONTRACTION */
1274    
1275    #   ifdef SEAICE_OCN_MELT_ACT_ON_AREA
1276  C ice cover reduction by joint OCN+ATM melt  C ice cover reduction by joint OCN+ATM melt
1277            tmpscal3 = MIN( 0. _d 0 ,            tmpscal3 = MIN( 0. _d 0 ,
1278       &              d_HEFFbyATMonOCN(I,J)+d_HEFFbyOCNonICE(I,J) )       &              d_HEFFbyATMonOCN(I,J)+d_HEFFbyOCNonICE(I,J) )
1279  # else  #   else
1280  C ice cover reduction by ATM melt only -- as in legacy code  C ice cover reduction by ATM melt only -- as in legacy code
1281            tmpscal3 = MIN( 0. _d 0 , d_HEFFbyATMonOCN(I,J) )            tmpscal3 = MIN( 0. _d 0 , d_HEFFbyATMonOCN(I,J) )
1282  # endif  #   endif
1283  C gain of new ice over open water  C gain of new ice over open water
1284    
1285  # ifdef SEAICE_DO_OPEN_WATER_GROWTH  #   ifdef SEAICE_DO_OPEN_WATER_GROWTH
1286  C the one effectively used to increment HEFF  C the one effectively used to increment HEFF
1287            tmpscal4 = d_HEFFbyATMonOCN_open(I,J)            tmpscal4 = MAX(ZERO,d_HEFFbyATMonOCN_open(I,J))
1288  # else  #   else
1289  C the virtual one -- as in legcy code  C the virtual one -- as in legcy code
1290            tmpscal4 = MAX(ZERO,a_QbyATM_open(I,J))            tmpscal4 = MAX(ZERO,a_QbyATM_open(I,J))
1291  # endif  #   endif
 #endif /* SEAICE_GROWTH_LEGACY */  
1292    
1293  C compute cover fraction tendency  C compute cover fraction tendency
1294            IF ( YC(I,J,bi,bj) .LT. ZERO ) THEN            d_AREAbyATM(I,J)=tmpscal4/HO_loc+
1295             d_AREAbyATM(I,J)=tmpscal4/HO_south       &                     HALF*tmpscal3*recip_heffActual(I,J)
1296            ELSE  
1297             d_AREAbyATM(I,J)=tmpscal4/HO  #   endif /* FENTY_AREA_EXPANSION_CONTRACTION */
           ENDIF  
           d_AREAbyATM(I,J)=d_AREAbyATM(I,J)  
 #ifdef SEAICE_GROWTH_LEGACY  
      &         +HALF*tmpscal3*AREApreTH(I,J)  
      &         /(tmpscal0+.00001 _d 0)  
 #else  
      &         +HALF*tmpscal3*recip_heffActual(I,J)  
 #endif  
1298    
1299  #endif /* FENTY_AREA_EXPANSION_CONTRACTION */  #endif /* SEAICE_GROWTH_LEGACY */
1300    
1301  C apply tendency  C apply tendency
1302            IF ( (HEFF(i,j,bi,bj).GT.0. _d 0).OR.            IF ( (HEFF(i,j,bi,bj).GT.0. _d 0).OR.
1303       &        (HSNOW(i,j,bi,bj).GT.0. _d 0) ) THEN       &        (HSNOW(i,j,bi,bj).GT.0. _d 0) ) THEN
1304             AREA(I,J,bi,bj)=max(0. _d 0 , min( 1. _d 0,             AREA(I,J,bi,bj)=max(0. _d 0, min( 1. _d 0, AREA(I,J,bi,bj)
1305       &                     AREA(I,J,bi,bj)+d_AREAbyATM(I,J)       &     +d_AREAbyATM(I,J) + d_AREAbyOCN(I,J) + d_AREAbyICE(I,J) ))
 #ifdef FENTY_AREA_EXPANSION_CONTRACTION  
      &                   + d_AREAbyOCN(I,J) + d_AREAbyICE(I,J)  
 #endif  
      &                                       ))  
1306            ELSE            ELSE
1307             AREA(I,J,bi,bj)=0. _d 0             AREA(I,J,bi,bj)=0. _d 0
1308            ENDIF            ENDIF
1309    
1310    #ifdef ALLOW_DIAGNOSTICS
1311              DIAGarray(I,J) = d_AREAbyICE(I,J) + d_AREAbyATM(I,J)
1312         &                     + d_AREAbyOCN(I,J)
1313    #endif
1314    
1315           ENDDO           ENDDO
1316          ENDDO          ENDDO
1317    
# Line 1353  C apply tendency Line 1322  C apply tendency
1322       &      'SIdA    ',0,1,3,bi,bj,myThid)       &      'SIdA    ',0,1,3,bi,bj,myThid)
1323           ENDIF           ENDIF
1324           IF ( DIAGNOSTICS_IS_ON('SIdAbATO',myThid) ) THEN           IF ( DIAGNOSTICS_IS_ON('SIdAbATO',myThid) ) THEN
1325            CALL DIAGNOSTICS_FILL(DIAGarrayA,            CALL DIAGNOSTICS_FILL(d_AREAbyATM,
1326       &      'SIdAbATO',0,1,3,bi,bj,myThid)       &      'SIdAbATO',0,1,3,bi,bj,myThid)
1327           ENDIF           ENDIF
1328           IF ( DIAGNOSTICS_IS_ON('SIdAbATC',myThid) ) THEN           IF ( DIAGNOSTICS_IS_ON('SIdAbATC',myThid) ) THEN
1329            CALL DIAGNOSTICS_FILL(DIAGarrayB,            CALL DIAGNOSTICS_FILL(d_AREAbyICE,
1330       &      'SIdAbATC',0,1,3,bi,bj,myThid)       &      'SIdAbATC',0,1,3,bi,bj,myThid)
1331           ENDIF           ENDIF
1332           IF ( DIAGNOSTICS_IS_ON('SIdAbOCN',myThid) ) THEN           IF ( DIAGNOSTICS_IS_ON('SIdAbOCN',myThid) ) THEN
1333            CALL DIAGNOSTICS_FILL(DIAGarrayC,            CALL DIAGNOSTICS_FILL(d_AREAbyOCN,
1334       &      'SIdAbOCN',0,1,3,bi,bj,myThid)       &      'SIdAbOCN',0,1,3,bi,bj,myThid)
1335           ENDIF           ENDIF
1336          ENDIF          ENDIF
Legend:
Removed from v.1.122  
changed lines
  Added in v.1.123
  ViewVC Help
Powered by ViewVC 1.1.22