/[MITgcm]/MITgcm_contrib/torge/itd/code/seaice_growth.F
ViewVC logotype

Diff of /MITgcm_contrib/torge/itd/code/seaice_growth.F

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

revision 1.4 by torge, Fri Sep 28 19:01:17 2012 UTC revision 1.6 by torge, Tue Oct 2 17:06:04 2012 UTC
# Line 141  c     The change of mean ice thickness d Line 141  c     The change of mean ice thickness d
141        _RL d_HEFFbyRLX         (1:sNx,1:sNy)        _RL d_HEFFbyRLX         (1:sNx,1:sNy)
142  #endif  #endif
143    
 #ifdef SEAICE_ITD  
 c     The change of mean ice area due to out-of-bounds values following  
 c     sea ice dynamics  
       _RL d_AREAbyNEG         (1:sNx,1:sNy)  
 #endif  
144  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
145  c     sea ice dynamics  c     sea ice dynamics
146        _RL d_HEFFbyNEG         (1:sNx,1:sNy)        _RL d_HEFFbyNEG         (1:sNx,1:sNy)
# Line 179  C                     as EVAP (positive Line 174  C                     as EVAP (positive
174        _RL d_HEFFbySublim      (1:sNx,1:sNy)        _RL d_HEFFbySublim      (1:sNx,1:sNy)
175        _RL d_HSNWbySublim      (1:sNx,1:sNy)        _RL d_HSNWbySublim      (1:sNx,1:sNy)
176    
177  #if (defined(SEAICE_CAP_SUBLIM) || defined(SEAICE_ITD))  #ifdef SEAICE_CAP_SUBLIM
178  C     The latent heat flux which will sublimate all snow and ice  C     The latent heat flux which will sublimate all snow and ice
179  C     over one time step  C     over one time step
180        _RL latentHeatFluxMax   (1:sNx,1:sNy)        _RL latentHeatFluxMax   (1:sNx,1:sNy)
# Line 233  C temporary variables available for the Line 228  C temporary variables available for the
228    
229        INTEGER ilockey        INTEGER ilockey
230        INTEGER it        INTEGER it
 #ifdef SEAICE_ITD  
       INTEGER K  
 #endif  
231        _RL pFac        _RL pFac
232        _RL ticeInMult          (1:sNx,1:sNy,MULTDIM)        _RL ticeInMult          (1:sNx,1:sNy,MULTDIM)
233        _RL ticeOutMult         (1:sNx,1:sNy,MULTDIM)        _RL ticeOutMult         (1:sNx,1:sNy,MULTDIM)
# Line 367  C ===================== Line 359  C =====================
359            d_HEFFbyRLX(I,J)           = 0.0 _d 0            d_HEFFbyRLX(I,J)           = 0.0 _d 0
360  #endif  #endif
361    
 #ifdef SEAICE_ITD  
           d_AREAbyNEG(I,J)           = 0.0 _d 0  
 #endif  
362            d_HEFFbyNEG(I,J)           = 0.0 _d 0            d_HEFFbyNEG(I,J)           = 0.0 _d 0
363            d_HEFFbyOCNonICE(I,J)      = 0.0 _d 0            d_HEFFbyOCNonICE(I,J)      = 0.0 _d 0
364            d_HEFFbyATMonOCN(I,J)      = 0.0 _d 0            d_HEFFbyATMonOCN(I,J)      = 0.0 _d 0
# Line 403  c Line 392  c
392              a_QbyATMmult_cover(I,J,IT)    = 0.0 _d 0              a_QbyATMmult_cover(I,J,IT)    = 0.0 _d 0
393              a_QSWbyATMmult_cover(I,J,IT)  = 0.0 _d 0              a_QSWbyATMmult_cover(I,J,IT)  = 0.0 _d 0
394              a_FWbySublimMult(I,J,IT)      = 0.0 _d 0              a_FWbySublimMult(I,J,IT)      = 0.0 _d 0
395    #ifdef SEAICE_CAP_SUBLIM
396                latentHeatFluxMaxMult(I,J,IT) = 0.0 _d 0
397    #endif
398  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
399              r_QbyATMmult_cover (I,J,IT)   = 0.0 _d 0              r_QbyATMmult_cover (I,J,IT)   = 0.0 _d 0
400              r_FWbySublimMult(I,J,IT)      = 0.0 _d 0              r_FWbySublimMult(I,J,IT)      = 0.0 _d 0
401  #endif  #endif
 #if (defined(SEAICE_CAP_SUBLIM) || defined(SEAICE_ITD))  
             latentHeatFluxMaxMult(I,J,IT) = 0.0 _d 0  
 #endif  
402            ENDDO            ENDDO
403           ENDDO           ENDDO
404          ENDDO          ENDDO
# Line 444  C ====================================== Line 433  C ======================================
433           ENDDO           ENDDO
434          ENDDO          ENDDO
435  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
436          DO K=1,nITD          DO IT=1,nITD
437           DO J=1,sNy           DO J=1,sNy
438            DO I=1,sNx            DO I=1,sNx
439             HEFFITDpreTH(I,J,K)=HEFFITD(I,J,K,bi,bj)             HEFFITDpreTH(I,J,IT)=HEFFITD(I,J,IT,bi,bj)
440             HSNWITDpreTH(I,J,K)=HSNOWITD(I,J,K,bi,bj)             HSNWITDpreTH(I,J,IT)=HSNOWITD(I,J,IT,bi,bj)
441             AREAITDpreTH(I,J,K)=AREAITD(I,J,K,bi,bj)             AREAITDpreTH(I,J,IT)=AREAITD(I,J,IT,bi,bj)
442            ENDDO            ENDDO
443           ENDDO           ENDDO
444          ENDDO          ENDDO
# Line 504  CADJ STORE area(:,:,bi,bj) = comlev1_bib Line 493  CADJ STORE area(:,:,bi,bj) = comlev1_bib
493          DO J=1,sNy          DO J=1,sNy
494           DO I=1,sNx           DO I=1,sNx
495  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
496            DO K=1,nITD            DO IT=1,nITD
            tmpscal1=0. _d 0  
497             tmpscal2=0. _d 0             tmpscal2=0. _d 0
498             tmpscal3=0. _d 0             tmpscal3=0. _d 0
499             tmpscal2=MAX(-HEFFITD(I,J,K,bi,bj),0. _d 0)             tmpscal2=MAX(-HEFFITD(I,J,IT,bi,bj),0. _d 0)
500             HEFFITD(I,J,K,bi,bj)=HEFFITD(I,J,K,bi,bj)+tmpscal2             HEFFITD(I,J,IT,bi,bj)=HEFFITD(I,J,IT,bi,bj)+tmpscal2
501             d_HEFFbyNEG(I,J)=d_HEFFbyNEG(I,J)+tmpscal2             d_HEFFbyNEG(I,J)=d_HEFFbyNEG(I,J)+tmpscal2
502             tmpscal3=MAX(-HSNOWITD(I,J,K,bi,bj),0. _d 0)             tmpscal3=MAX(-HSNOWITD(I,J,IT,bi,bj),0. _d 0)
503             HSNOWITD(I,J,K,bi,bj)=HSNOWITD(I,J,K,bi,bj)+tmpscal3             HSNOWITD(I,J,IT,bi,bj)=HSNOWITD(I,J,IT,bi,bj)+tmpscal3
504             d_HSNWbyNEG(I,J)=d_HSNWbyNEG(I,J)+tmpscal3             d_HSNWbyNEG(I,J)=d_HSNWbyNEG(I,J)+tmpscal3
505             tmpscal1=MAX(-AREAITD(I,J,K,bi,bj),0. _d 0)             AREAITD(I,J,IT,bi,bj)=MAX(-AREAITD(I,J,IT,bi,bj),0. _d 0)
            AREAITD(I,J,K,bi,bj)=AREAITD(I,J,K,bi,bj)+tmpscal1  
            d_AREAbyNEG(I,J)=d_AREAbyNEG(I,J)+tmpscal1  
506            ENDDO            ENDDO
507  CToM      AREA, HEFF, and HSNOW will be updated at end of PART 1  CToM      AREA, HEFF, and HSNOW will be updated at end of PART 1
508  C         by calling SEAICE_ITD_SUM  C         by calling SEAICE_ITD_SUM
509  #else  #else
510            d_HEFFbyNEG(I,J)=MAX(-HEFF(I,J,bi,bj),0. _d 0)            d_HEFFbyNEG(I,J)=MAX(-HEFF(I,J,bi,bj),0. _d 0)
           d_HSNWbyNEG(I,J)=MAX(-HSNOW(I,J,bi,bj),0. _d 0)  
           AREA(I,J,bi,bj)=MAX(AREA(I,J,bi,bj),0. _d 0)  
511            HEFF(I,J,bi,bj)=HEFF(I,J,bi,bj)+d_HEFFbyNEG(I,J)            HEFF(I,J,bi,bj)=HEFF(I,J,bi,bj)+d_HEFFbyNEG(I,J)
512              d_HSNWbyNEG(I,J)=MAX(-HSNOW(I,J,bi,bj),0. _d 0)
513            HSNOW(I,J,bi,bj)=HSNOW(I,J,bi,bj)+d_HSNWbyNEG(I,J)            HSNOW(I,J,bi,bj)=HSNOW(I,J,bi,bj)+d_HSNWbyNEG(I,J)
514              AREA(I,J,bi,bj)=MAX(AREA(I,J,bi,bj),0. _d 0)
515  #endif  #endif
516           ENDDO           ENDDO
517          ENDDO          ENDDO
# Line 538  CADJ STORE heff(:,:,bi,bj)  = comlev1_bi Line 524  CADJ STORE heff(:,:,bi,bj)  = comlev1_bi
524          DO J=1,sNy          DO J=1,sNy
525           DO I=1,sNx           DO I=1,sNx
526  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
527            DO K=1,nITD            DO IT=1,nITD
528  #endif  #endif
529            tmpscal2=0. _d 0            tmpscal2=0. _d 0
530            tmpscal3=0. _d 0            tmpscal3=0. _d 0
531  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
532             IF (HEFFITD(I,J,K,bi,bj).LE.siEps) THEN             IF (HEFFITD(I,J,IT,bi,bj).LE.siEps) THEN
533              tmpscal2=-HEFFITD(I,J,K,bi,bj)              tmpscal2=-HEFFITD(I,J,IT,bi,bj)
534              tmpscal3=-HSNOWITD(I,J,K,bi,bj)              tmpscal3=-HSNOWITD(I,J,IT,bi,bj)
535              TICES(I,J,K,bi,bj)=celsius2K              TICES(I,J,IT,bi,bj)=celsius2K
536  CToM        TICE will be updated at end of Part 1 together with AREA and HEFF  CToM        TICE will be updated at end of Part 1 together with AREA and HEFF
537             ENDIF             ENDIF
538             HEFFITD(I,J,K,bi,bj) =HEFFITD(I,J,K,bi,bj) +tmpscal2             HEFFITD(I,J,IT,bi,bj) =HEFFITD(I,J,IT,bi,bj) +tmpscal2
539             HSNOWITD(I,J,K,bi,bj)=HSNOWITD(I,J,K,bi,bj)+tmpscal3             HSNOWITD(I,J,IT,bi,bj)=HSNOWITD(I,J,IT,bi,bj)+tmpscal3
540  #else  #else
541            IF (HEFF(I,J,bi,bj).LE.siEps) THEN            IF (HEFF(I,J,bi,bj).LE.siEps) THEN
542             tmpscal2=-HEFF(I,J,bi,bj)             tmpscal2=-HEFF(I,J,bi,bj)
# Line 580  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi Line 566  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi
566          DO J=1,sNy          DO J=1,sNy
567           DO I=1,sNx           DO I=1,sNx
568  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
569            DO K=1,nITD            DO IT=1,nITD
570             IF ((HEFFITD(i,j,k,bi,bj).EQ.0. _d 0).AND.             IF ((HEFFITD(I,J,IT,bi,bj).EQ.0. _d 0).AND.
571       &         (HSNOWITD(i,j,k,bi,bj).EQ.0. _d 0))       &         (HSNOWITD(I,J,IT,bi,bj).EQ.0. _d 0))
572       &      AREAITD(I,J,K,bi,bj)=0. _d 0       &      AREAITD(I,J,IT,bi,bj)=0. _d 0
573            ENDDO            ENDDO
574  #else  #else
575            IF ((HEFF(i,j,bi,bj).EQ.0. _d 0).AND.            IF ((HEFF(i,j,bi,bj).EQ.0. _d 0).AND.
# Line 601  CADJ STORE area(:,:,bi,bj)  = comlev1_bi Line 587  CADJ STORE area(:,:,bi,bj)  = comlev1_bi
587          DO J=1,sNy          DO J=1,sNy
588           DO I=1,sNx           DO I=1,sNx
589  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
590            DO K=1,nITD            DO IT=1,nITD
591             IF ((HEFFITD(i,j,k,bi,bj).GT.0).OR.             IF ((HEFFITD(I,J,IT,bi,bj).GT.0).OR.
592       &         (HSNOWITD(i,j,k,bi,bj).GT.0)) THEN       &         (HSNOWITD(I,J,IT,bi,bj).GT.0)) THEN
593  CToM        SEAICE_area_floor*nITD cannot be allowed to exceed 1  CToM        SEAICE_area_floor*nITD cannot be allowed to exceed 1
594  C           hence use SEAICE_area_floor devided by nITD  C           hence use SEAICE_area_floor devided by nITD
595  C           (or install a warning in e.g. seaice_readparms.F)  C           (or install a warning in e.g. seaice_readparms.F)
596              AREAITD(I,J,K,bi,bj)=              AREAITD(I,J,IT,bi,bj)=
597       &         MAX(AREAITD(I,J,K,bi,bj),SEAICE_area_floor/float(nITD))       &         MAX(AREAITD(I,J,IT,bi,bj),SEAICE_area_floor/float(nITD))
598             ENDIF             ENDIF
599            ENDDO            ENDDO
600  #else  #else
# Line 639  CADJ STORE area(:,:,bi,bj)  = comlev1_bi Line 625  CADJ STORE area(:,:,bi,bj)  = comlev1_bi
625            AREA(I,J,bi,bj)=MIN(AREA(I,J,bi,bj),SEAICE_area_max)            AREA(I,J,bi,bj)=MIN(AREA(I,J,bi,bj),SEAICE_area_max)
626           ENDDO           ENDDO
627          ENDDO          ENDDO
628  #endif /* SEAICE_ITD */  #endif /* notSEAICE_ITD */
629    
630  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
631  CToM catch up with items 1.25 and 2.5 involving category sums AREA and HEFF  CToM catch up with items 1.25 and 2.5 involving category sums AREA and HEFF
 C    first, update AREA and HEFF:  
         CALL SEAICE_ITD_SUM(bi, bj, myTime, myIter, myThid)  
 C  
632          DO J=1,sNy          DO J=1,sNy
633           DO I=1,sNx           DO I=1,sNx
634  C    TICES was changed above (item 1.25), now update TICE as ice volume  C    TICES was changed above (item 1.25), now update TICE as ice volume
635  C     weighted average of TICES  C     weighted average of TICES
636    C    also compute total of AREAITD (needed for finishing item 2.5, see below)
637            tmpscal1 = 0. _d 0            tmpscal1 = 0. _d 0
638            tmpscal2 = 0. _d 0            tmpscal2 = 0. _d 0
639            DO K=1,nITD            tmpscal3 = 0. _d 0
640             tmpscal1=tmpscal1 + TICES(I,J,K,bi,bj)*HEFFITD(I,J,K,bi,bj)            DO IT=1,nITD
641             tmpscal2=tmpscal2 + HEFFITD(I,J,K,bi,bj)             tmpscal1=tmpscal1 + TICES(I,J,IT,bi,bj)*HEFFITD(I,J,IT,bi,bj)
642               tmpscal2=tmpscal2 + HEFFITD(I,J,IT,bi,bj)
643               tmpscal3=tmpscal3 + AREAITD(I,J,IT,bi,bj)
644            ENDDO            ENDDO
645            TICE(I,J,bi,bj)=tmpscal1/tmpscal2            TICE(I,J,bi,bj)=tmpscal1/tmpscal2
646  C    lines of item 2.5 that were omitted:  C    lines of item 2.5 that were omitted:
# Line 662  C    in 2.5 these lines are executed bef Line 648  C    in 2.5 these lines are executed bef
648  C    hence we execute them here before SEAICE_ITD_REDIST is called  C    hence we execute them here before SEAICE_ITD_REDIST is called
649  C    although this means that AREA has not been completely regularized  C    although this means that AREA has not been completely regularized
650  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
651            DIAGarrayA(I,J) = AREA(I,J,bi,bj)            DIAGarrayA(I,J) = tmpscal3
652  #endif  #endif
653  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
654            SItrAREA(I,J,bi,bj,1)=AREA(I,J,bi,bj)            SItrAREA(I,J,bi,bj,1)=tmpscal3
655  #endif  #endif
656           ENDDO           ENDDO
657          ENDDO          ENDDO
# Line 678  C    and update AREA, HEFF, and HSNOW Line 664  C    and update AREA, HEFF, and HSNOW
664    
665  #endif  #endif
666  #if (defined ALLOW_AUTODIFF_TAMC && defined SEAICE_MODIFY_GROWTH_ADJ)  #if (defined ALLOW_AUTODIFF_TAMC && defined SEAICE_MODIFY_GROWTH_ADJ)
667  C        ENDIF SEAICEadjMODE.EQ.0  C        end SEAICEadjMODE.EQ.0 statement:
668          ENDIF          ENDIF
669  #endif  #endif
670    
# Line 700  C 3) store regularized values of heff, h Line 686  C 3) store regularized values of heff, h
686           ENDDO           ENDDO
687          ENDDO          ENDDO
688  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
689          DO K=1,nITD          DO IT=1,nITD
690           DO J=1,sNy           DO J=1,sNy
691            DO I=1,sNx            DO I=1,sNx
692             HEFFITDpreTH(I,J,K)=HEFFITD(I,J,K,bi,bj)             HEFFITDpreTH(I,J,IT)=HEFFITD(I,J,IT,bi,bj)
693             HSNWITDpreTH(I,J,K)=HSNOWITD(I,J,K,bi,bj)             HSNWITDpreTH(I,J,IT)=HSNOWITD(I,J,IT,bi,bj)
694             AREAITDpreTH(I,J,K)=AREAITD(I,J,K,bi,bj)             AREAITDpreTH(I,J,IT)=AREAITD(I,J,IT,bi,bj)
695    
696  C memorize areal and volume fraction of each ITD category  C memorize areal and volume fraction of each ITD category
697             IF (AREA(I,J,bi,bj).GT.0) THEN             IF (AREA(I,J,bi,bj).GT.0) THEN
698              areaFracFactor(I,J,K)=AREAITD(I,J,K,bi,bj)/AREA(I,J,bi,bj)              areaFracFactor(I,J,IT)=AREAITD(I,J,IT,bi,bj)/AREA(I,J,bi,bj)
699             ELSE             ELSE
700              areaFracFactor(I,J,K)=ZERO              areaFracFactor(I,J,IT)=ZERO
701             ENDIF             ENDIF
702             IF (HEFF(I,J,bi,bj).GT.0) THEN             IF (HEFF(I,J,bi,bj).GT.0) THEN
703              heffFracFactor(I,J,K)=HEFFITD(I,J,K,bi,bj)/HEFF(I,J,bi,bj)              heffFracFactor(I,J,IT)=HEFFITD(I,J,IT,bi,bj)/HEFF(I,J,bi,bj)
704             ELSE             ELSE
705              heffFracFactor(I,J,K)=ZERO              heffFracFactor(I,J,IT)=ZERO
706             ENDIF             ENDIF
707            ENDDO            ENDDO
708           ENDDO           ENDDO
709          ENDDO          ENDDO
710  C prepare SItrHEFF to be computed as cumulative sum  C prepare SItrHEFF to be computed as cumulative sum
711          DO K=2,5          DO iTr=2,5
712           DO J=1,sNy           DO J=1,sNy
713            DO I=1,sNx            DO I=1,sNx
714             SItrHEFF(I,J,bi,bj,K)=ZERO             SItrHEFF(I,J,bi,bj,iTr)=ZERO
715            ENDDO            ENDDO
716           ENDDO           ENDDO
717          ENDDO          ENDDO
# Line 791  Cgf no additional dependency of air-sea Line 777  Cgf no additional dependency of air-sea
777            ENDDO            ENDDO
778           ENDDO           ENDDO
779  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
780           DO K=1,nITD           DO IT=1,nITD
781            DO J=1,sNy            DO J=1,sNy
782             DO I=1,sNx             DO I=1,sNx
783              HEFFITDpreTH(I,J,K) = 0. _d 0              HEFFITDpreTH(I,J,IT) = 0. _d 0
784              HSNWITDpreTH(I,J,K) = 0. _d 0              HSNWITDpreTH(I,J,IT) = 0. _d 0
785              AREAITDpreTH(I,J,K) = 0. _d 0              AREAITDpreTH(I,J,IT) = 0. _d 0
786             ENDDO             ENDDO
787            ENDDO            ENDDO
788           ENDDO           ENDDO
# Line 821  CADJ STORE HEFFpreTH = comlev1_bibj, key Line 807  CADJ STORE HEFFpreTH = comlev1_bibj, key
807  CADJ STORE HSNWpreTH = comlev1_bibj, key = iicekey, byte = isbyte  CADJ STORE HSNWpreTH = comlev1_bibj, key = iicekey, byte = isbyte
808  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
809  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
810          DO K=1,nITD          DO IT=1,nITD
811  #endif  #endif
812          DO J=1,sNy          DO J=1,sNy
813           DO I=1,sNx           DO I=1,sNx
814    
815  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
816            IF (HEFFITDpreTH(I,J,K) .GT. ZERO) THEN            IF (HEFFITDpreTH(I,J,IT) .GT. ZERO) THEN
817  #ifdef SEAICE_GROWTH_LEGACY  #ifdef SEAICE_GROWTH_LEGACY
818              tmpscal1          = MAX(SEAICE_area_reg/float(nITD),              tmpscal1          = MAX(SEAICE_area_reg/float(nITD),
819       &                              AREAITDpreTH(I,J,K))       &                              AREAITDpreTH(I,J,IT))
820              hsnowActualMult(I,J,K) = HSNWITDpreTH(I,J,K)/tmpscal1              hsnowActualMult(I,J,IT) = HSNWITDpreTH(I,J,IT)/tmpscal1
821              tmpscal2               = HEFFITDpreTH(I,J,K)/tmpscal1              tmpscal2               = HEFFITDpreTH(I,J,IT)/tmpscal1
822              heffActualMult(I,J,K)  = MAX(tmpscal2,SEAICE_hice_reg)              heffActualMult(I,J,IT)  = MAX(tmpscal2,SEAICE_hice_reg)
823  #else /* SEAICE_GROWTH_LEGACY */  #else /* SEAICE_GROWTH_LEGACY */
824  cif        regularize AREA with SEAICE_area_reg  cif        regularize AREA with SEAICE_area_reg
825             tmpscal1 = SQRT(AREAITDpreTH(I,J,K) * AREAITDpreTH(I,J,K)             tmpscal1 = SQRT(AREAITDpreTH(I,J,IT) * AREAITDpreTH(I,J,IT)
826       &                     + area_reg_sq)       &                     + area_reg_sq)
827  cif        heffActual calculated with the regularized AREA  cif        heffActual calculated with the regularized AREA
828             tmpscal2 = HEFFITDpreTH(I,J,K) / tmpscal1             tmpscal2 = HEFFITDpreTH(I,J,IT) / tmpscal1
829  cif        regularize heffActual with SEAICE_hice_reg (add lower bound)  cif        regularize heffActual with SEAICE_hice_reg (add lower bound)
830             heffActualMult(I,J,K) = SQRT(tmpscal2 * tmpscal2             heffActualMult(I,J,IT) = SQRT(tmpscal2 * tmpscal2
831       &                                  + hice_reg_sq)       &                                  + hice_reg_sq)
832  cif        hsnowActual calculated with the regularized AREA  cif        hsnowActual calculated with the regularized AREA
833             hsnowActualMult(I,J,K) = HSNWITDpreTH(I,J,K) / tmpscal1             hsnowActualMult(I,J,IT) = HSNWITDpreTH(I,J,IT) / tmpscal1
834  #endif /* SEAICE_GROWTH_LEGACY */  #endif /* SEAICE_GROWTH_LEGACY */
835  cif        regularize the inverse of heffActual by hice_reg  cif        regularize the inverse of heffActual by hice_reg
836             recip_heffActualMult(I,J,K)  = AREAITDpreTH(I,J,K) /             recip_heffActualMult(I,J,IT)  = AREAITDpreTH(I,J,IT) /
837       &                 sqrt(HEFFITDpreTH(I,J,K) * HEFFITDpreTH(I,J,K)       &                 sqrt(HEFFITDpreTH(I,J,IT) * HEFFITDpreTH(I,J,IT)
838       &                      + hice_reg_sq)       &                      + hice_reg_sq)
839  cif       Do not regularize when HEFFpreTH = 0  cif       Do not regularize when HEFFpreTH = 0
840            ELSE            ELSE
841              heffActualMult(I,J,K) = ZERO              heffActualMult(I,J,IT) = ZERO
842              hsnowActualMult(I,J,K) = ZERO              hsnowActualMult(I,J,IT) = ZERO
843              recip_heffActualMult(I,J,K)  = ZERO              recip_heffActualMult(I,J,IT)  = ZERO
844            ENDIF            ENDIF
845  #else /* SEAICE_ITD */  #else /* SEAICE_ITD */
846            IF (HEFFpreTH(I,J) .GT. ZERO) THEN            IF (HEFFpreTH(I,J) .GT. ZERO) THEN
# Line 900  cif       Do not regularize when HEFFpre Line 886  cif       Do not regularize when HEFFpre
886  C 5) COMPUTE MAXIMUM LATENT HEAT FLUXES FOR THE CURRENT ICE  C 5) COMPUTE MAXIMUM LATENT HEAT FLUXES FOR THE CURRENT ICE
887  C    AND SNOW THICKNESS  C    AND SNOW THICKNESS
888  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
889          DO K=1,nITD          DO IT=1,nITD
890  #endif  #endif
891          DO J=1,sNy          DO J=1,sNy
892           DO I=1,sNx           DO I=1,sNx
# Line 908  c        The latent heat flux over the s Line 894  c        The latent heat flux over the s
894  c        will sublimate all of the snow and ice over one time  c        will sublimate all of the snow and ice over one time
895  c        step (W/m^2)  c        step (W/m^2)
896  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
897            IF (HEFFITDpreTH(I,J,K) .GT. ZERO) THEN            IF (HEFFITDpreTH(I,J,IT) .GT. ZERO) THEN
898              latentHeatFluxMaxMult(I,J,K) = lhSublim*recip_deltaTtherm *              latentHeatFluxMaxMult(I,J,IT) = lhSublim*recip_deltaTtherm *
899       &         (HEFFITDpreTH(I,J,K)*SEAICE_rhoIce +       &         (HEFFITDpreTH(I,J,IT)*SEAICE_rhoIce +
900       &          HSNWITDpreTH(I,J,K)*SEAICE_rhoSnow)/AREAITDpreTH(I,J,K)       &          HSNWITDpreTH(I,J,IT)*SEAICE_rhoSnow)
901         &         /AREAITDpreTH(I,J,IT)
902            ELSE            ELSE
903              latentHeatFluxMaxMult(I,J,K) = ZERO              latentHeatFluxMaxMult(I,J,IT) = ZERO
904            ENDIF            ENDIF
905  #else  #else
906            IF (HEFFpreTH(I,J) .GT. ZERO) THEN            IF (HEFFpreTH(I,J) .GT. ZERO) THEN
# Line 1084  C      the temperature here is a result Line 1071  C      the temperature here is a result
1071  C      computed individually for each single category in SEAICE_SOLVE4TEMP  C      computed individually for each single category in SEAICE_SOLVE4TEMP
1072  C      and hence is averaged area weighted [areaFracFactor])  C      and hence is averaged area weighted [areaFracFactor])
1073             TICE(I,J,bi,bj) = TICE(I,J,bi,bj)             TICE(I,J,bi,bj) = TICE(I,J,bi,bj)
1074       &          +  ticeOutMult(I,J,IT)*areaFracFactor(I,J,K)       &          +  ticeOutMult(I,J,IT)*areaFracFactor(I,J,IT)
1075  #else  #else
1076             TICE(I,J,bi,bj) = TICE(I,J,bi,bj)             TICE(I,J,bi,bj) = TICE(I,J,bi,bj)
1077       &          +  ticeOutMult(I,J,IT)*recip_multDim       &          +  ticeOutMult(I,J,IT)*recip_multDim
# Line 1095  C     average over categories Line 1082  C     average over categories
1082  C     calculate area weighted mean  C     calculate area weighted mean
1083  C     (fluxes are per unit (ice surface) area and are thus area weighted)  C     (fluxes are per unit (ice surface) area and are thus area weighted)
1084             a_QbyATM_cover   (I,J) = a_QbyATM_cover(I,J)             a_QbyATM_cover   (I,J) = a_QbyATM_cover(I,J)
1085       &          + a_QbyATMmult_cover(I,J,IT)*areaFracFactor(I,J,K)       &          + a_QbyATMmult_cover(I,J,IT)*areaFracFactor(I,J,IT)
1086             a_QSWbyATM_cover (I,J) = a_QSWbyATM_cover(I,J)             a_QSWbyATM_cover (I,J) = a_QSWbyATM_cover(I,J)
1087       &          + a_QSWbyATMmult_cover(I,J,IT)*areaFracFactor(I,J,K)       &          + a_QSWbyATMmult_cover(I,J,IT)*areaFracFactor(I,J,IT)
1088             a_FWbySublim     (I,J) = a_FWbySublim(I,J)             a_FWbySublim     (I,J) = a_FWbySublim(I,J)
1089       &          + a_FWbySublimMult(I,J,IT)*areaFracFactor(I,J,K)       &          + a_FWbySublimMult(I,J,IT)*areaFracFactor(I,J,IT)
1090  #else  #else
1091             a_QbyATM_cover   (I,J) = a_QbyATM_cover(I,J)             a_QbyATM_cover   (I,J) = a_QbyATM_cover(I,J)
1092       &          + a_QbyATMmult_cover(I,J,IT)*recip_multDim       &          + a_QbyATMmult_cover(I,J,IT)*recip_multDim
# Line 1141  CADJ STORE a_FWbySublim    = comlev1_bib Line 1128  CADJ STORE a_FWbySublim    = comlev1_bib
1128    
1129  C switch heat fluxes from W/m2 to 'effective' ice meters  C switch heat fluxes from W/m2 to 'effective' ice meters
1130  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1131          DO K=1,nITD          DO IT=1,nITD
1132           DO J=1,sNy           DO J=1,sNy
1133            DO I=1,sNx            DO I=1,sNx
1134             a_QbyATMmult_cover(I,J,K)   = a_QbyATMmult_cover(I,J,K)             a_QbyATMmult_cover(I,J,IT)   = a_QbyATMmult_cover(I,J,IT)
1135       &          * convertQ2HI * AREAITDpreTH(I,J,K)       &          * convertQ2HI * AREAITDpreTH(I,J,IT)
1136             a_QSWbyATMmult_cover(I,J,K) = a_QSWbyATMmult_cover(I,J,K)             a_QSWbyATMmult_cover(I,J,IT) = a_QSWbyATMmult_cover(I,J,IT)
1137       &          * convertQ2HI * AREAITDpreTH(I,J,K)       &          * convertQ2HI * AREAITDpreTH(I,J,IT)
1138  C and initialize r_QbyATM_cover  C and initialize r_QbyATM_cover
1139             r_QbyATMmult_cover(I,J,K)=a_QbyATMmult_cover(I,J,K)             r_QbyATMmult_cover(I,J,IT)=a_QbyATMmult_cover(I,J,IT)
1140  C     Convert fresh water flux by sublimation to 'effective' ice meters.  C     Convert fresh water flux by sublimation to 'effective' ice meters.
1141  C     Negative sublimation is resublimation and will be added as snow.  C     Negative sublimation is resublimation and will be added as snow.
1142  #ifdef SEAICE_DISABLE_SUBLIM  #ifdef SEAICE_DISABLE_SUBLIM
1143             a_FWbySublimMult(I,J,K) = ZERO             a_FWbySublimMult(I,J,IT) = ZERO
1144  #endif  #endif
1145             a_FWbySublimMult(I,J,K) = SEAICE_deltaTtherm*recip_rhoIce             a_FWbySublimMult(I,J,IT) = SEAICE_deltaTtherm*recip_rhoIce
1146       &            * a_FWbySublimMult(I,J,K)*AREAITDpreTH(I,J,K)       &            * a_FWbySublimMult(I,J,IT)*AREAITDpreTH(I,J,IT)
1147             r_FWbySublimMult(I,J,K)=a_FWbySublimMult(I,J,K)             r_FWbySublimMult(I,J,IT)=a_FWbySublimMult(I,J,IT)
1148            ENDDO            ENDDO
1149           ENDDO           ENDDO
1150          ENDDO          ENDDO
# Line 1214  CADJ STORE r_FWbySublim    = comlev1_bib Line 1201  CADJ STORE r_FWbySublim    = comlev1_bib
1201  Cgf no additional dependency through ice cover  Cgf no additional dependency through ice cover
1202          IF ( SEAICEadjMODE.GE.3 ) THEN          IF ( SEAICEadjMODE.GE.3 ) THEN
1203  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1204           DO K=1,nITD           DO IT=1,nITD
1205            DO J=1,sNy            DO J=1,sNy
1206             DO I=1,sNx             DO I=1,sNx
1207              a_QbyATMmult_cover(I,J,K)   = 0. _d 0              a_QbyATMmult_cover(I,J,IT)   = 0. _d 0
1208              r_QbyATMmult_cover(I,J,K)   = 0. _d 0              r_QbyATMmult_cover(I,J,IT)   = 0. _d 0
1209              a_QSWbyATMmult_cover(I,J,K) = 0. _d 0              a_QSWbyATMmult_cover(I,J,IT) = 0. _d 0
1210             ENDDO             ENDDO
1211            ENDDO            ENDDO
1212           ENDDO           ENDDO
# Line 1296  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi Line 1283  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi
1283  CADJ STORE r_FWbySublim     = comlev1_bibj,key=iicekey,byte=isbyte  CADJ STORE r_FWbySublim     = comlev1_bibj,key=iicekey,byte=isbyte
1284  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1285  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1286          DO K=1,nITD          DO IT=1,nITD
1287  #endif  #endif
1288          DO J=1,sNy          DO J=1,sNy
1289           DO I=1,sNx           DO I=1,sNx
1290  C     First sublimate/deposite snow  C     First sublimate/deposite snow
1291            tmpscal2 =            tmpscal2 =
1292  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1293       &     MAX(MIN(r_FWbySublimMult(I,J,K),HSNOWITD(I,J,K,bi,bj)       &     MAX(MIN(r_FWbySublimMult(I,J,IT),HSNOWITD(I,J,IT,bi,bj)
1294       &             *SNOW2ICE),ZERO)       &             *SNOW2ICE),ZERO)
1295  C         accumulate change over ITD categories  C         accumulate change over ITD categories
1296            d_HSNWbySublim(I,J)     = d_HSNWbySublim(I,J)     - tmpscal2            d_HSNWbySublim(I,J)     = d_HSNWbySublim(I,J)      - tmpscal2
1297       &                                                       *ICE2SNOW       &                                                       *ICE2SNOW
1298            HSNOWITD(I,J,K,bi,bj)   = HSNOWITD(I,J,K,bi,bj)   - tmpscal2            HSNOWITD(I,J,IT,bi,bj)  = HSNOWITD(I,J,IT,bi,bj)   - tmpscal2
1299       &                                                       *ICE2SNOW       &                                                       *ICE2SNOW
1300            r_FWbySublimMult(I,J,K) = r_FWbySublimMult(I,J,K) - tmpscal2            r_FWbySublimMult(I,J,IT)= r_FWbySublimMult(I,J,IT) - tmpscal2
1301  C         keep total up to date, too  C         keep total up to date, too
1302            r_FWbySublim(I,J)       = r_FWbySublim(I,J)       - tmpscal2            r_FWbySublim(I,J)       = r_FWbySublim(I,J)        - tmpscal2
1303  #else  #else
1304       &     MAX(MIN(r_FWbySublim(I,J),HSNOW(I,J,bi,bj)*SNOW2ICE),ZERO)       &     MAX(MIN(r_FWbySublim(I,J),HSNOW(I,J,bi,bj)*SNOW2ICE),ZERO)
1305            d_HSNWbySublim(I,J) = - tmpscal2 * ICE2SNOW            d_HSNWbySublim(I,J) = - tmpscal2 * ICE2SNOW
# Line 1330  CADJ STORE r_FWbySublim    = comlev1_bib Line 1317  CADJ STORE r_FWbySublim    = comlev1_bib
1317  C     If anything is left, sublimate ice  C     If anything is left, sublimate ice
1318            tmpscal2 =            tmpscal2 =
1319  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1320       &     MAX(MIN(r_FWbySublimMult(I,J,K),HEFFITD(I,J,K,bi,bj)),ZERO)       &     MAX(MIN(r_FWbySublimMult(I,J,IT),HEFFITD(I,J,IT,bi,bj)),ZERO)
1321  C         accumulate change over ITD categories  C         accumulate change over ITD categories
1322            d_HSNWbySublim(I,J)     = d_HSNWbySublim(I,J)     - tmpscal2            d_HSNWbySublim(I,J)     = d_HSNWbySublim(I,J)     - tmpscal2
1323            HEFFITD(I,J,K,bi,bj)    = HEFFITD(I,J,K,bi,bj)    - tmpscal2            HEFFITD(I,J,IT,bi,bj)    = HEFFITD(I,J,IT,bi,bj)    - tmpscal2
1324            r_FWbySublimMult(I,J,K) = r_FWbySublimMult(I,J,K) - tmpscal2            r_FWbySublimMult(I,J,IT) = r_FWbySublimMult(I,J,IT) - tmpscal2
1325  C         keep total up to date, too  C         keep total up to date, too
1326            r_FWbySublim(I,J)       = r_FWbySublim(I,J)       - tmpscal2            r_FWbySublim(I,J)       = r_FWbySublim(I,J)       - tmpscal2
1327  #else  #else
# Line 1351  C     If anything is left, it will be ev Line 1338  C     If anything is left, it will be ev
1338  C     Since a_QbyATM_cover was computed for sublimation, not simple evaporation, we need to  C     Since a_QbyATM_cover was computed for sublimation, not simple evaporation, we need to
1339  C     remove the fusion part for the residual (that happens to be precisely r_FWbySublim).  C     remove the fusion part for the residual (that happens to be precisely r_FWbySublim).
1340  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1341            a_QbyATMmult_cover(I,J,K) = a_QbyATMmult_cover(I,J,K)            a_QbyATMmult_cover(I,J,IT) = a_QbyATMmult_cover(I,J,IT)
1342       &                              - r_FWbySublimMult(I,J,K)       &                               - r_FWbySublimMult(I,J,IT)
1343            r_QbyATMmult_cover(I,J,K) = r_QbyATMmult_cover(I,J,K)            r_QbyATMmult_cover(I,J,IT) = r_QbyATMmult_cover(I,J,IT)
1344       &                              - r_FWbySublimMult(I,J,K)       &                               - r_FWbySublimMult(I,J,IT)
1345           ENDDO           ENDDO
1346          ENDDO          ENDDO
1347  C       end K loop  C       end IT loop
1348          ENDDO          ENDDO
1349  C       then update totals        C       then update totals      
1350          DO J=1,sNy          DO J=1,sNy
# Line 1389  CADJ STORE r_QbyOCN = comlev1_bibj,key=i Line 1376  CADJ STORE r_QbyOCN = comlev1_bibj,key=i
1376  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1377    
1378  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1379          DO K=1,nITD          DO IT=1,nITD
1380           DO J=1,sNy           DO J=1,sNy
1381            DO I=1,sNx            DO I=1,sNx
1382  C          ice growth/melt due to ocean heat is equally distributed under the ice  C          ice growth/melt due to ocean heat is equally distributed under the ice
1383  C          and hence weighted by fractional area of each thickness category  C          and hence weighted by fractional area of each thickness category
1384             tmpscal1=MAX(r_QbyOCN(i,j)*areaFracFactor(I,J,K),             tmpscal1=MAX(r_QbyOCN(i,j)*areaFracFactor(I,J,IT),
1385       &                               -HEFFITD(I,J,K,bi,bj))       &                               -HEFFITD(I,J,IT,bi,bj))
1386             d_HEFFbyOCNonICE(I,J)= d_HEFFbyOCNonICE(I,J) + tmpscal1             d_HEFFbyOCNonICE(I,J)= d_HEFFbyOCNonICE(I,J) + tmpscal1
1387             r_QbyOCN(I,J)        = r_QbyOCN(I,J)         - tmpscal1             r_QbyOCN(I,J)        = r_QbyOCN(I,J)         - tmpscal1
1388             HEFFITD(I,J,K,bi,bj) = HEFFITD(I,J,K,bi,bj)  + tmpscal1             HEFFITD(I,J,IT,bi,bj) = HEFFITD(I,J,IT,bi,bj)  + tmpscal1
1389  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
1390             SItrHEFF(I,J,bi,bj,2) = SItrHEFF(I,J,bi,bj,2)             SItrHEFF(I,J,bi,bj,2) = SItrHEFF(I,J,bi,bj,2)
1391       &                           + HEFFITD(I,J,K,bi,bj)       &                           + HEFFITD(I,J,IT,bi,bj)
1392  #endif  #endif
1393            ENDDO            ENDDO
1394           ENDDO           ENDDO
# Line 1440  CADJ STORE r_QbyATM_cover = comlev1_bibj Line 1427  CADJ STORE r_QbyATM_cover = comlev1_bibj
1427  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1428    
1429  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1430          DO K=1,nITD          DO IT=1,nITD
1431           DO J=1,sNy           DO J=1,sNy
1432            DO I=1,sNx            DO I=1,sNx
1433  C     Convert to standard units (meters of ice) rather than to meters  C     Convert to standard units (meters of ice) rather than to meters
1434  C     of snow. This appears to be more robust.  C     of snow. This appears to be more robust.
1435            tmpscal1=MAX(r_QbyATMmult_cover(I,J,K),-HSNOWITD(I,J,K,bi,bj)             tmpscal1=MAX(r_QbyATMmult_cover(I,J,IT),
1436       &                                           *SNOW2ICE)       &                  -HSNOWITD(I,J,IT,bi,bj)*SNOW2ICE)
1437            tmpscal2=MIN(tmpscal1,0. _d 0)             tmpscal2=MIN(tmpscal1,0. _d 0)
1438  #ifdef SEAICE_MODIFY_GROWTH_ADJ  #ifdef SEAICE_MODIFY_GROWTH_ADJ
1439  Cgf no additional dependency through snow  Cgf no additional dependency through snow
1440            IF ( SEAICEadjMODE.GE.2 ) tmpscal2 = 0. _d 0             IF ( SEAICEadjMODE.GE.2 ) tmpscal2 = 0. _d 0
1441  #endif  #endif
1442            d_HSNWbyATMonSNW(I,J)=d_HSNWbyATMonSNW(I,J)+tmpscal2*ICE2SNOW             d_HSNWbyATMonSNW(I,J) = d_HSNWbyATMonSNW(I,J)
1443            HSNOWITD(I,J,K,bi,bj)=HSNOWITD(I,J,K,bi,bj)+tmpscal2*ICE2SNOW       &                           + tmpscal2*ICE2SNOW
1444            r_QbyATMmult_cover(I,J,K)=r_QbyATMmult_cover(I,J,K) - tmpscal2             HSNOWITD(I,J,IT,bi,bj)= HSNOWITD(I,J,IT,bi,bj)
1445  C         keep the total up to date, too       &                           + tmpscal2*ICE2SNOW
1446            r_QbyATM_cover(I,J)=r_QbyATM_cover(I,J) - tmpscal2             r_QbyATMmult_cover(I,J,IT)=r_QbyATMmult_cover(I,J,IT)
1447         &                           - tmpscal2
1448    C          keep the total up to date, too
1449               r_QbyATM_cover(I,J)=r_QbyATM_cover(I,J) - tmpscal2
1450            ENDDO            ENDDO
1451           ENDDO           ENDDO
1452          ENDDO          ENDDO
# Line 1504  Cgf the v1.81=>v1.82 revision would chan Line 1494  Cgf the v1.81=>v1.82 revision would chan
1494  Cgf warming conditions, the lab_sea results were not changed.  Cgf warming conditions, the lab_sea results were not changed.
1495    
1496  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1497          DO K=1,nITD          DO IT=1,nITD
1498           DO J=1,sNy           DO J=1,sNy
1499            DO I=1,sNx            DO I=1,sNx
1500  #ifdef SEAICE_GROWTH_LEGACY  #ifdef SEAICE_GROWTH_LEGACY
1501             tmpscal2 =MAX(-HEFFITD(I,J,K,bi,bj),r_QbyATMmult_cover(I,J,K))             tmpscal2 = MAX(-HEFFITD(I,J,IT,bi,bj),
1502         &                     r_QbyATMmult_cover(I,J,IT))
1503  #else  #else
1504             tmpscal2 =MAX(-HEFFITD(I,J,K,bi,bj),r_QbyATMmult_cover(I,J,K)             tmpscal2 = MAX(-HEFFITD(I,J,IT,bi,bj),
1505         &                     r_QbyATMmult_cover(I,J,IT)
1506  c         Limit ice growth by potential melt by ocean  c         Limit ice growth by potential melt by ocean
1507       &     + AREAITDpreTH(I,J,K) * r_QbyOCN(I,J))       &              + AREAITDpreTH(I,J,IT) * r_QbyOCN(I,J))
1508  #endif /* SEAICE_GROWTH_LEGACY */  #endif /* SEAICE_GROWTH_LEGACY */
1509             d_HEFFbyATMonOCN_cover(I,J) = d_HEFFbyATMonOCN_cover(I,J)             d_HEFFbyATMonOCN_cover(I,J) = d_HEFFbyATMonOCN_cover(I,J)
1510       &                                 + tmpscal2       &                                 + tmpscal2
# Line 1520  c         Limit ice growth by potential Line 1512  c         Limit ice growth by potential
1512       &                                 + tmpscal2       &                                 + tmpscal2
1513             r_QbyATM_cover(I,J)         = r_QbyATM_cover(I,J)             r_QbyATM_cover(I,J)         = r_QbyATM_cover(I,J)
1514       &                                 - tmpscal2       &                                 - tmpscal2
1515             HEFFITD(I,J,K,bi,bj) = HEFFITD(I,J,K,bi,bj) + tmpscal2             HEFFITD(I,J,IT,bi,bj) = HEFFITD(I,J,IT,bi,bj) + tmpscal2
1516    
1517  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
1518             SItrHEFF(I,J,bi,bj,3) = SItrHEFF(I,J,bi,bj,3)             SItrHEFF(I,J,bi,bj,3) = SItrHEFF(I,J,bi,bj,3)
1519       &                           + HEFFITD(I,J,K,bi,bj)       &                           + HEFFITD(I,J,IT,bi,bj)
1520  #endif  #endif
1521            ENDDO            ENDDO
1522           ENDDO           ENDDO
# Line 1593  C           add precip to the fresh wate Line 1585  C           add precip to the fresh wate
1585           ENDDO           ENDDO
1586          ENDDO          ENDDO
1587  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1588          DO K=1,nITD          DO IT=1,nITD
1589           DO J=1,sNy           DO J=1,sNy
1590            DO I=1,sNx            DO I=1,sNx
1591             HSNOWITD(I,J,K,bi,bj) = HSNOWITD(I,J,K,bi,bj)             HSNOWITD(I,J,IT,bi,bj) = HSNOWITD(I,J,IT,bi,bj)
1592       &     + d_HSNWbyRAIN(I,J)*areaFracFactor(I,J,K)       &     + d_HSNWbyRAIN(I,J)*areaFracFactor(I,J,IT)
1593            ENDDO            ENDDO
1594           ENDDO           ENDDO
1595          ENDDO          ENDDO
# Line 1637  CADJ STORE r_QbyOCN = comlev1_bibj,key=i Line 1629  CADJ STORE r_QbyOCN = comlev1_bibj,key=i
1629  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1630    
1631  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1632          DO K=1,nITD          DO IT=1,nITD
1633           DO J=1,sNy           DO J=1,sNy
1634            DO I=1,sNx            DO I=1,sNx
1635             tmpscal1=MAX(r_QbyOCN(i,j)*ICE2SNOW*areaFracFactor(I,J,K),             tmpscal1=MAX(r_QbyOCN(i,j)*ICE2SNOW*areaFracFactor(I,J,IT),
1636       &                  -HSNOWITD(I,J,K,bi,bj))       &                  -HSNOWITD(I,J,IT,bi,bj))
1637             tmpscal2=MIN(tmpscal1,0. _d 0)             tmpscal2=MIN(tmpscal1,0. _d 0)
1638  #ifdef SEAICE_MODIFY_GROWTH_ADJ  #ifdef SEAICE_MODIFY_GROWTH_ADJ
1639  Cgf no additional dependency through snow  Cgf no additional dependency through snow
# Line 1649  Cgf no additional dependency through sno Line 1641  Cgf no additional dependency through sno
1641  #endif  #endif
1642             d_HSNWbyOCNonSNW(I,J) = d_HSNWbyOCNonSNW(I,J) + tmpscal2             d_HSNWbyOCNonSNW(I,J) = d_HSNWbyOCNonSNW(I,J) + tmpscal2
1643             r_QbyOCN(I,J)=r_QbyOCN(I,J) - tmpscal2*SNOW2ICE             r_QbyOCN(I,J)=r_QbyOCN(I,J) - tmpscal2*SNOW2ICE
1644             HSNOWITD(I,J,K,bi,bj) = HSNOWITD(I,J,K,bi,bj) + tmpscal2             HSNOWITD(I,J,IT,bi,bj) = HSNOWITD(I,J,IT,bi,bj) + tmpscal2
1645            ENDDO            ENDDO
1646           ENDDO           ENDDO
1647          ENDDO          ENDDO
# Line 1717  C         determine thickness of new ice Line 1709  C         determine thickness of new ice
1709  C         considering the entire open water area to refreeze  C         considering the entire open water area to refreeze
1710            tmpscal1 = tmpscal3/tmpscal0            tmpscal1 = tmpscal3/tmpscal0
1711  C         then add new ice volume to appropriate thickness category  C         then add new ice volume to appropriate thickness category
1712            DO K=1,nITD            DO IT=1,nITD
1713             IF (tmpscal1.LT.Hlimit(K)) THEN             IF (tmpscal1.LT.Hlimit(IT)) THEN
1714              HEFFITD(I,J,K,bi,bj) = HEFFITD(I,J,K,bi,bj) + tmpscal3              HEFFITD(I,J,IT,bi,bj) = HEFFITD(I,J,IT,bi,bj) + tmpscal3
1715              tmpscal3=ZERO              tmpscal3=ZERO
1716  C not sure if AREAITD should be changed here since AREA is incremented  C not sure if AREAITD should be changed here since AREA is incremented
1717  C   in PART 4 below in non-itd code  C   in PART 4 below in non-itd code
1718  C in this scenario all open water is covered by new ice instantaneously,  C in this scenario all open water is covered by new ice instantaneously,
1719  C   i.e. no delayed lead closing is concidered such as is associated with  C   i.e. no delayed lead closing is concidered such as is associated with
1720  C   Hibler's h_0 parameter  C   Hibler's h_0 parameter
1721              AREAITD(I,J,K,bi,bj) = AREAITD(I,J,K,bi,bj)              AREAITD(I,J,IT,bi,bj) = AREAITD(I,J,IT,bi,bj)
1722       &                           + tmpscal0       &                           + tmpscal0
1723              tmpscal0=ZERO              tmpscal0=ZERO
1724             ENDIF             ENDIF
# Line 1739  C   Hibler's h_0 parameter Line 1731  C   Hibler's h_0 parameter
1731    
1732  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
1733  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1734          DO K=1,nITD          DO IT=1,nITD
1735           DO J=1,sNy           DO J=1,sNy
1736            DO I=1,sNx            DO I=1,sNx
1737  c needs to be here to allow use also with LEGACY branch  c needs to be here to allow use also with LEGACY branch
1738             SItrHEFF(I,J,bi,bj,4) = SItrHEFF(I,J,bi,bj,4)             SItrHEFF(I,J,bi,bj,4) = SItrHEFF(I,J,bi,bj,4)
1739       &                           + HEFFITD(I,J,K,bi,bj)       &                           + HEFFITD(I,J,IT,bi,bj)
1740            ENDDO            ENDDO
1741           ENDDO           ENDDO
1742          ENDDO          ENDDO
# Line 1781  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi Line 1773  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi
1773  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1774          IF ( SEAICEuseFlooding ) THEN          IF ( SEAICEuseFlooding ) THEN
1775  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1776           DO K=1,nITD           DO IT=1,nITD
1777            DO J=1,sNy            DO J=1,sNy
1778             DO I=1,sNx             DO I=1,sNx
1779              tmpscal0 = (HSNOWITD(I,J,K,bi,bj)*SEAICE_rhoSnow              tmpscal0 = (HSNOWITD(I,J,IT,bi,bj)*SEAICE_rhoSnow
1780       &               +HEFFITD(I,J,K,bi,bj)*SEAICE_rhoIce)*recip_rhoConst       &               +  HEFFITD(I,J,IT,bi,bj) *SEAICE_rhoIce)
1781              tmpscal1 = MAX( 0. _d 0, tmpscal0 - HEFFITD(I,J,K,bi,bj))       &                                        *recip_rhoConst
1782              d_HEFFbyFLOODING(I,J) = d_HEFFbyFLOODING(I,J) + tmpscal1              tmpscal1 = MAX( 0. _d 0, tmpscal0 - HEFFITD(I,J,IT,bi,bj))
1783              HEFFITD(I,J,K,bi,bj)  = HEFFITD(I,J,K,bi,bj)  + tmpscal1              d_HEFFbyFLOODING(I,J) = d_HEFFbyFLOODING(I,J)  + tmpscal1
1784              HSNOWITD(I,J,K,bi,bj) = HSNOWITD(I,J,K,bi,bj) - tmpscal1              HEFFITD(I,J,IT,bi,bj) = HEFFITD(I,J,IT,bi,bj)  + tmpscal1
1785                HSNOWITD(I,J,IT,bi,bj)= HSNOWITD(I,J,IT,bi,bj) - tmpscal1
1786       &                            * ICE2SNOW       &                            * ICE2SNOW
1787             ENDDO             ENDDO
1788            ENDDO            ENDDO
# Line 1859  C       does not account for lateral mel Line 1852  C       does not account for lateral mel
1852  C  C
1853  C        AREAITD is incremented in section "gain of new ice over open water" above  C        AREAITD is incremented in section "gain of new ice over open water" above
1854  C  C
1855          DO K=1,nITD          DO IT=1,nITD
1856           DO J=1,sNy           DO J=1,sNy
1857            DO I=1,sNx            DO I=1,sNx
1858             IF (HEFFITD(I,J,K,bi,bj).LE.ZERO) THEN             IF (HEFFITD(I,J,IT,bi,bj).LE.ZERO) THEN
1859              AREAITD(I,J,K,bi,bj)=ZERO              AREAITD(I,J,IT,bi,bj)=ZERO
1860             ENDIF             ENDIF
1861  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
1862             SItrAREA(I,J,bi,bj,3) = SItrAREA(I,J,bi,bj,3)             SItrAREA(I,J,bi,bj,3) = SItrAREA(I,J,bi,bj,3)
1863       &                           + AREAITD(I,J,K,bi,bj)       &                           + AREAITD(I,J,IT,bi,bj)
1864  #endif /* ALLOW_SITRACER */  #endif /* ALLOW_SITRACER */
1865            ENDDO            ENDDO
1866           ENDDO           ENDDO
# Line 1948  C apply tendency Line 1941  C apply tendency
1941  Cgf 'bulk' linearization of area=f(HEFF)  Cgf 'bulk' linearization of area=f(HEFF)
1942          IF ( SEAICEadjMODE.GE.1 ) THEN          IF ( SEAICEadjMODE.GE.1 ) THEN
1943  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1944           DO K=1,nITD           DO IT=1,nITD
1945            DO J=1,sNy            DO J=1,sNy
1946             DO I=1,sNx             DO I=1,sNx
1947              AREAITD(I,J,K,bi,bj) = AREAITDpreTH(I,J,K) + 0.1 _d 0 *              AREAITD(I,J,IT,bi,bj) = AREAITDpreTH(I,J,IT) + 0.1 _d 0 *
1948       &               ( HEFFITD(I,J,K,bi,bj) - HEFFITDpreTH(I,J,K) )       &               ( HEFFITD(I,J,IT,bi,bj) - HEFFITDpreTH(I,J,IT) )
1949             ENDDO             ENDDO
1950            ENDDO            ENDDO
1951           ENDDO           ENDDO

Legend:
Removed from v.1.4  
changed lines
  Added in v.1.6

  ViewVC Help
Powered by ViewVC 1.1.22