/[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.3 by torge, Wed Sep 26 17:50:17 2012 UTC revision 1.5 by torge, Tue Oct 2 16:47:41 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)  #ifdef SEAICE_ITD
181        _RL latentHeatFluxMaxMult  (1:sNx,1:sNy,MULTDIM)        _RL latentHeatFluxMaxMult  (1:sNx,1:sNy,MULTDIM)
182    #else
183          _RL latentHeatFluxMax   (1:sNx,1:sNy)
184    #endif
185  #endif  #endif
186    
187  C     actual ice thickness (with upper and lower limit)  C     actual ice thickness (with upper and lower limit)
# Line 233  C temporary variables available for the Line 231  C temporary variables available for the
231    
232        INTEGER ilockey        INTEGER ilockey
233        INTEGER it        INTEGER it
 #ifdef SEAICE_ITD  
       INTEGER K  
 #endif  
234        _RL pFac        _RL pFac
235        _RL ticeInMult          (1:sNx,1:sNy,MULTDIM)        _RL ticeInMult          (1:sNx,1:sNy,MULTDIM)
236        _RL ticeOutMult         (1:sNx,1:sNy,MULTDIM)        _RL ticeOutMult         (1:sNx,1:sNy,MULTDIM)
# Line 367  C ===================== Line 362  C =====================
362            d_HEFFbyRLX(I,J)           = 0.0 _d 0            d_HEFFbyRLX(I,J)           = 0.0 _d 0
363  #endif  #endif
364    
 #ifdef SEAICE_ITD  
           d_AREAbyNEG(I,J)           = 0.0 _d 0  
 #endif  
365            d_HEFFbyNEG(I,J)           = 0.0 _d 0            d_HEFFbyNEG(I,J)           = 0.0 _d 0
366            d_HEFFbyOCNonICE(I,J)      = 0.0 _d 0            d_HEFFbyOCNonICE(I,J)      = 0.0 _d 0
367            d_HEFFbyATMonOCN(I,J)      = 0.0 _d 0            d_HEFFbyATMonOCN(I,J)      = 0.0 _d 0
# Line 387  C ===================== Line 379  C =====================
379            d_HEFFbySublim(I,J)        = 0.0 _d 0            d_HEFFbySublim(I,J)        = 0.0 _d 0
380            d_HSNWbySublim(I,J)        = 0.0 _d 0            d_HSNWbySublim(I,J)        = 0.0 _d 0
381  #ifdef SEAICE_CAP_SUBLIM  #ifdef SEAICE_CAP_SUBLIM
382    #ifdef SEAICE_ITD
383              DO IT=1,SEAICE_multDim
384               latentHeatFluxMaxMult(I,J,IT) = 0.0 _d 0
385              ENDDO
386    #else
387            latentHeatFluxMax(I,J)     = 0.0 _d 0            latentHeatFluxMax(I,J)     = 0.0 _d 0
388  #endif  #endif
389    #endif
390  c  c
391            d_HFRWbyRAIN(I,J)          = 0.0 _d 0            d_HFRWbyRAIN(I,J)          = 0.0 _d 0
392    
# Line 407  c Line 405  c
405              r_QbyATMmult_cover (I,J,IT)   = 0.0 _d 0              r_QbyATMmult_cover (I,J,IT)   = 0.0 _d 0
406              r_FWbySublimMult(I,J,IT)      = 0.0 _d 0              r_FWbySublimMult(I,J,IT)      = 0.0 _d 0
407  #endif  #endif
 #if (defined(SEAICE_CAP_SUBLIM) || defined(SEAICE_ITD))  
             latentHeatFluxMaxMult(I,J,IT) = 0.0 _d 0  
 #endif  
408            ENDDO            ENDDO
409           ENDDO           ENDDO
410          ENDDO          ENDDO
# Line 444  C ====================================== Line 439  C ======================================
439           ENDDO           ENDDO
440          ENDDO          ENDDO
441  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
442          DO K=1,nITD          DO IT=1,nITD
443           DO J=1,sNy           DO J=1,sNy
444            DO I=1,sNx            DO I=1,sNx
445             HEFFITDpreTH(I,J,K)=HEFFITD(I,J,K,bi,bj)             HEFFITDpreTH(I,J,IT)=HEFFITD(I,J,IT,bi,bj)
446             HSNWITDpreTH(I,J,K)=HSNOWITD(I,J,K,bi,bj)             HSNWITDpreTH(I,J,IT)=HSNOWITD(I,J,IT,bi,bj)
447             AREAITDpreTH(I,J,K)=AREAITD(I,J,K,bi,bj)             AREAITDpreTH(I,J,IT)=AREAITD(I,J,IT,bi,bj)
448            ENDDO            ENDDO
449           ENDDO           ENDDO
450          ENDDO          ENDDO
# Line 504  CADJ STORE area(:,:,bi,bj) = comlev1_bib Line 499  CADJ STORE area(:,:,bi,bj) = comlev1_bib
499          DO J=1,sNy          DO J=1,sNy
500           DO I=1,sNx           DO I=1,sNx
501  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
502            DO K=1,nITD            DO IT=1,nITD
            tmpscal1=0. _d 0  
503             tmpscal2=0. _d 0             tmpscal2=0. _d 0
504             tmpscal3=0. _d 0             tmpscal3=0. _d 0
505             tmpscal2=MAX(-HEFFITD(I,J,K,bi,bj),0. _d 0)             tmpscal2=MAX(-HEFFITD(I,J,IT,bi,bj),0. _d 0)
506             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
507             d_HEFFbyNEG(I,J)=d_HEFFbyNEG(I,J)+tmpscal2             d_HEFFbyNEG(I,J)=d_HEFFbyNEG(I,J)+tmpscal2
508             tmpscal3=MAX(-HSNOWITD(I,J,K,bi,bj),0. _d 0)             tmpscal3=MAX(-HSNOWITD(I,J,IT,bi,bj),0. _d 0)
509             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
510             d_HSNWbyNEG(I,J)=d_HSNWbyNEG(I,J)+tmpscal3             d_HSNWbyNEG(I,J)=d_HSNWbyNEG(I,J)+tmpscal3
511             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  
512            ENDDO            ENDDO
513  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
514  C         by calling SEAICE_ITD_SUM  C         by calling SEAICE_ITD_SUM
515  #else  #else
516            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)  
517            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)
518              d_HSNWbyNEG(I,J)=MAX(-HSNOW(I,J,bi,bj),0. _d 0)
519            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)
520              AREA(I,J,bi,bj)=MAX(AREA(I,J,bi,bj),0. _d 0)
521  #endif  #endif
522           ENDDO           ENDDO
523          ENDDO          ENDDO
# Line 538  CADJ STORE heff(:,:,bi,bj)  = comlev1_bi Line 530  CADJ STORE heff(:,:,bi,bj)  = comlev1_bi
530          DO J=1,sNy          DO J=1,sNy
531           DO I=1,sNx           DO I=1,sNx
532  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
533            DO K=1,nITD            DO IT=1,nITD
534  #endif  #endif
535            tmpscal2=0. _d 0            tmpscal2=0. _d 0
536            tmpscal3=0. _d 0            tmpscal3=0. _d 0
537  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
538             IF (HEFFITD(I,J,K,bi,bj).LE.siEps) THEN             IF (HEFFITD(I,J,IT,bi,bj).LE.siEps) THEN
539              tmpscal2=-HEFFITD(I,J,K,bi,bj)              tmpscal2=-HEFFITD(I,J,IT,bi,bj)
540              tmpscal3=-HSNOWITD(I,J,K,bi,bj)              tmpscal3=-HSNOWITD(I,J,IT,bi,bj)
541              TICES(I,J,K,bi,bj)=celsius2K              TICES(I,J,IT,bi,bj)=celsius2K
542  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
543             ENDIF             ENDIF
544             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
545             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
546  #else  #else
547            IF (HEFF(I,J,bi,bj).LE.siEps) THEN            IF (HEFF(I,J,bi,bj).LE.siEps) THEN
548             tmpscal2=-HEFF(I,J,bi,bj)             tmpscal2=-HEFF(I,J,bi,bj)
# Line 580  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi Line 572  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi
572          DO J=1,sNy          DO J=1,sNy
573           DO I=1,sNx           DO I=1,sNx
574  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
575            DO K=1,nITD            DO IT=1,nITD
576             IF ((HEFFITD(i,j,k,bi,bj).EQ.0. _d 0).AND.             IF ((HEFFITD(I,J,IT,bi,bj).EQ.0. _d 0).AND.
577       &         (HSNOWITD(i,j,k,bi,bj).EQ.0. _d 0))       &         (HSNOWITD(I,J,IT,bi,bj).EQ.0. _d 0))
578       &      AREAITD(I,J,K,bi,bj)=0. _d 0       &      AREAITD(I,J,IT,bi,bj)=0. _d 0
579            ENDDO            ENDDO
580  #else  #else
581            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 593  CADJ STORE area(:,:,bi,bj)  = comlev1_bi
593          DO J=1,sNy          DO J=1,sNy
594           DO I=1,sNx           DO I=1,sNx
595  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
596            DO K=1,nITD            DO IT=1,nITD
597             IF ((HEFFITD(i,j,k,bi,bj).GT.0).OR.             IF ((HEFFITD(I,J,IT,bi,bj).GT.0).OR.
598       &         (HSNOWITD(i,j,k,bi,bj).GT.0)) THEN       &         (HSNOWITD(I,J,IT,bi,bj).GT.0)) THEN
599  CToM        SEAICE_area_floor*nITD cannot be allowed to exceed 1  CToM        SEAICE_area_floor*nITD cannot be allowed to exceed 1
600  C           hence use SEAICE_area_floor devided by nITD  C           hence use SEAICE_area_floor devided by nITD
601  C           (or install a warning in e.g. seaice_readparms.F)  C           (or install a warning in e.g. seaice_readparms.F)
602              AREAITD(I,J,K,bi,bj)=              AREAITD(I,J,IT,bi,bj)=
603       &         MAX(AREAITD(I,J,K,bi,bj),SEAICE_area_floor/float(nITD))       &         MAX(AREAITD(I,J,IT,bi,bj),SEAICE_area_floor/float(nITD))
604             ENDIF             ENDIF
605            ENDDO            ENDDO
606  #else  #else
# Line 639  CADJ STORE area(:,:,bi,bj)  = comlev1_bi Line 631  CADJ STORE area(:,:,bi,bj)  = comlev1_bi
631            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)
632           ENDDO           ENDDO
633          ENDDO          ENDDO
634  #endif /* SEAICE_ITD */  #endif /* notSEAICE_ITD */
635    
636  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
637  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  
638          DO J=1,sNy          DO J=1,sNy
639           DO I=1,sNx           DO I=1,sNx
640  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
641  C     weighted average of TICES  C     weighted average of TICES
642    C    also compute total of AREAITD (needed for finishing item 2.5, see below)
643            tmpscal1 = 0. _d 0            tmpscal1 = 0. _d 0
644            tmpscal2 = 0. _d 0            tmpscal2 = 0. _d 0
645            DO K=1,nITD            tmpscal3 = 0. _d 0
646             tmpscal1=tmpscal1 + TICES(I,J,K,bi,bj)*HEFFITD(I,J,K,bi,bj)            DO IT=1,nITD
647             tmpscal2=tmpscal2 + HEFFITD(I,J,K,bi,bj)             tmpscal1=tmpscal1 + TICES(I,J,IT,bi,bj)*HEFFITD(I,J,IT,bi,bj)
648               tmpscal2=tmpscal2 + HEFFITD(I,J,IT,bi,bj)
649               tmpscal3=tmpscal3 + AREAITD(I,J,IT,bi,bj)
650            ENDDO            ENDDO
651            TICE(I,J,bi,bj)=tmpscal1/tmpscal2            TICE(I,J,bi,bj)=tmpscal1/tmpscal2
652  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 654  C    in 2.5 these lines are executed bef
654  C    hence we execute them here before SEAICE_ITD_REDIST is called  C    hence we execute them here before SEAICE_ITD_REDIST is called
655  C    although this means that AREA has not been completely regularized  C    although this means that AREA has not been completely regularized
656  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
657            DIAGarrayA(I,J) = AREA(I,J,bi,bj)            DIAGarrayA(I,J) = tmpscal3
658  #endif  #endif
659  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
660            SItrAREA(I,J,bi,bj,1)=AREA(I,J,bi,bj)            SItrAREA(I,J,bi,bj,1)=tmpscal3
661  #endif  #endif
662           ENDDO           ENDDO
663          ENDDO          ENDDO
# Line 678  C    and update AREA, HEFF, and HSNOW Line 670  C    and update AREA, HEFF, and HSNOW
670    
671  #endif  #endif
672  #if (defined ALLOW_AUTODIFF_TAMC && defined SEAICE_MODIFY_GROWTH_ADJ)  #if (defined ALLOW_AUTODIFF_TAMC && defined SEAICE_MODIFY_GROWTH_ADJ)
673  C        ENDIF SEAICEadjMODE.EQ.0  C        end SEAICEadjMODE.EQ.0 statement:
674          ENDIF          ENDIF
675  #endif  #endif
676    
# Line 700  C 3) store regularized values of heff, h Line 692  C 3) store regularized values of heff, h
692           ENDDO           ENDDO
693          ENDDO          ENDDO
694  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
695          DO K=1,nITD          DO IT=1,nITD
696           DO J=1,sNy           DO J=1,sNy
697            DO I=1,sNx            DO I=1,sNx
698             HEFFITDpreTH(I,J,K)=HEFFITD(I,J,K,bi,bj)             HEFFITDpreTH(I,J,IT)=HEFFITD(I,J,IT,bi,bj)
699             HSNWITDpreTH(I,J,K)=HSNOWITD(I,J,K,bi,bj)             HSNWITDpreTH(I,J,IT)=HSNOWITD(I,J,IT,bi,bj)
700             AREAITDpreTH(I,J,K)=AREAITD(I,J,K,bi,bj)             AREAITDpreTH(I,J,IT)=AREAITD(I,J,IT,bi,bj)
701    
702  C memorize areal and volume fraction of each ITD category  C memorize areal and volume fraction of each ITD category
703             IF (AREA(I,J,bi,bj).GT.0) THEN             IF (AREA(I,J,bi,bj).GT.0) THEN
704              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)
705             ELSE             ELSE
706              areaFracFactor(I,J,K)=ZERO              areaFracFactor(I,J,IT)=ZERO
707             ENDIF             ENDIF
708             IF (HEFF(I,J,bi,bj).GT.0) THEN             IF (HEFF(I,J,bi,bj).GT.0) THEN
709              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)
710             ELSE             ELSE
711              heffFracFactor(I,J,K)=ZERO              heffFracFactor(I,J,IT)=ZERO
712             ENDIF             ENDIF
713            ENDDO            ENDDO
714           ENDDO           ENDDO
715          ENDDO          ENDDO
716  C prepare SItrHEFF to be computed as cumulative sum  C prepare SItrHEFF to be computed as cumulative sum
717          DO K=2,5          DO iTr=2,5
718           DO J=1,sNy           DO J=1,sNy
719            DO I=1,sNx            DO I=1,sNx
720             SItrHEFF(I,J,bi,bj,K)=ZERO             SItrHEFF(I,J,bi,bj,iTr)=ZERO
721            ENDDO            ENDDO
722           ENDDO           ENDDO
723          ENDDO          ENDDO
# Line 791  Cgf no additional dependency of air-sea Line 783  Cgf no additional dependency of air-sea
783            ENDDO            ENDDO
784           ENDDO           ENDDO
785  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
786           DO K=1,nITD           DO IT=1,nITD
787            DO J=1,sNy            DO J=1,sNy
788             DO I=1,sNx             DO I=1,sNx
789              HEFFITDpreTH(I,J,K) = 0. _d 0              HEFFITDpreTH(I,J,IT) = 0. _d 0
790              HSNWITDpreTH(I,J,K) = 0. _d 0              HSNWITDpreTH(I,J,IT) = 0. _d 0
791              AREAITDpreTH(I,J,K) = 0. _d 0              AREAITDpreTH(I,J,IT) = 0. _d 0
792             ENDDO             ENDDO
793            ENDDO            ENDDO
794           ENDDO           ENDDO
# Line 821  CADJ STORE HEFFpreTH = comlev1_bibj, key Line 813  CADJ STORE HEFFpreTH = comlev1_bibj, key
813  CADJ STORE HSNWpreTH = comlev1_bibj, key = iicekey, byte = isbyte  CADJ STORE HSNWpreTH = comlev1_bibj, key = iicekey, byte = isbyte
814  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
815  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
816          DO K=1,nITD          DO IT=1,nITD
817  #endif  #endif
818          DO J=1,sNy          DO J=1,sNy
819           DO I=1,sNx           DO I=1,sNx
820    
821  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
822            IF (HEFFITDpreTH(I,J,K) .GT. ZERO) THEN            IF (HEFFITDpreTH(I,J,IT) .GT. ZERO) THEN
823  #ifdef SEAICE_GROWTH_LEGACY  #ifdef SEAICE_GROWTH_LEGACY
824              tmpscal1          = MAX(SEAICE_area_reg/float(nITD),              tmpscal1          = MAX(SEAICE_area_reg/float(nITD),
825       &                              AREAITDpreTH(I,J,K))       &                              AREAITDpreTH(I,J,IT))
826              hsnowActualMult(I,J,K) = HSNWITDpreTH(I,J,K)/tmpscal1              hsnowActualMult(I,J,IT) = HSNWITDpreTH(I,J,IT)/tmpscal1
827              tmpscal2               = HEFFITDpreTH(I,J,K)/tmpscal1              tmpscal2               = HEFFITDpreTH(I,J,IT)/tmpscal1
828              heffActualMult(I,J,K)  = MAX(tmpscal2,SEAICE_hice_reg)              heffActualMult(I,J,IT)  = MAX(tmpscal2,SEAICE_hice_reg)
829  #else /* SEAICE_GROWTH_LEGACY */  #else /* SEAICE_GROWTH_LEGACY */
830  cif        regularize AREA with SEAICE_area_reg  cif        regularize AREA with SEAICE_area_reg
831             tmpscal1 = SQRT(AREAITDpreTH(I,J,K) * AREAITDpreTH(I,J,K)             tmpscal1 = SQRT(AREAITDpreTH(I,J,IT) * AREAITDpreTH(I,J,IT)
832       &                     + area_reg_sq)       &                     + area_reg_sq)
833  cif        heffActual calculated with the regularized AREA  cif        heffActual calculated with the regularized AREA
834             tmpscal2 = HEFFITDpreTH(I,J,K) / tmpscal1             tmpscal2 = HEFFITDpreTH(I,J,IT) / tmpscal1
835  cif        regularize heffActual with SEAICE_hice_reg (add lower bound)  cif        regularize heffActual with SEAICE_hice_reg (add lower bound)
836             heffActualMult(I,J,K) = SQRT(tmpscal2 * tmpscal2             heffActualMult(I,J,IT) = SQRT(tmpscal2 * tmpscal2
837       &                                  + hice_reg_sq)       &                                  + hice_reg_sq)
838  cif        hsnowActual calculated with the regularized AREA  cif        hsnowActual calculated with the regularized AREA
839             hsnowActualMult(I,J,K) = HSNWITDpreTH(I,J,K) / tmpscal1             hsnowActualMult(I,J,IT) = HSNWITDpreTH(I,J,IT) / tmpscal1
840  #endif /* SEAICE_GROWTH_LEGACY */  #endif /* SEAICE_GROWTH_LEGACY */
841  cif        regularize the inverse of heffActual by hice_reg  cif        regularize the inverse of heffActual by hice_reg
842             recip_heffActualMult(I,J,K)  = AREAITDpreTH(I,J,K) /             recip_heffActualMult(I,J,IT)  = AREAITDpreTH(I,J,IT) /
843       &                 sqrt(HEFFITDpreTH(I,J,K) * HEFFITDpreTH(I,J,K)       &                 sqrt(HEFFITDpreTH(I,J,IT) * HEFFITDpreTH(I,J,IT)
844       &                      + hice_reg_sq)       &                      + hice_reg_sq)
845  cif       Do not regularize when HEFFpreTH = 0  cif       Do not regularize when HEFFpreTH = 0
846            ELSE            ELSE
847              heffActualMult(I,J,K) = ZERO              heffActualMult(I,J,IT) = ZERO
848              hsnowActualMult(I,J,K) = ZERO              hsnowActualMult(I,J,IT) = ZERO
849              recip_heffActualMult(I,J,K)  = ZERO              recip_heffActualMult(I,J,IT)  = ZERO
850            ENDIF            ENDIF
851  #else /* SEAICE_ITD */  #else /* SEAICE_ITD */
852            IF (HEFFpreTH(I,J) .GT. ZERO) THEN            IF (HEFFpreTH(I,J) .GT. ZERO) THEN
# Line 900  cif       Do not regularize when HEFFpre Line 892  cif       Do not regularize when HEFFpre
892  C 5) COMPUTE MAXIMUM LATENT HEAT FLUXES FOR THE CURRENT ICE  C 5) COMPUTE MAXIMUM LATENT HEAT FLUXES FOR THE CURRENT ICE
893  C    AND SNOW THICKNESS  C    AND SNOW THICKNESS
894  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
895          DO K=1,nITD          DO IT=1,nITD
896  #endif  #endif
897          DO J=1,sNy          DO J=1,sNy
898           DO I=1,sNx           DO I=1,sNx
# Line 908  c        The latent heat flux over the s Line 900  c        The latent heat flux over the s
900  c        will sublimate all of the snow and ice over one time  c        will sublimate all of the snow and ice over one time
901  c        step (W/m^2)  c        step (W/m^2)
902  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
903            IF (HEFFITDpreTH(I,J,K) .GT. ZERO) THEN            IF (HEFFITDpreTH(I,J,IT) .GT. ZERO) THEN
904              latentHeatFluxMaxMult(I,J,K) = lhSublim*recip_deltaTtherm *              latentHeatFluxMaxMult(I,J,IT) = lhSublim*recip_deltaTtherm *
905       &         (HEFFITDpreTH(I,J,K)*SEAICE_rhoIce +       &         (HEFFITDpreTH(I,J,IT)*SEAICE_rhoIce +
906       &          HSNWITDpreTH(I,J,K)*SEAICE_rhoSnow)/AREAITDpreTH(I,J,K)       &          HSNWITDpreTH(I,J,IT)*SEAICE_rhoSnow)
907         &         /AREAITDpreTH(I,J,IT)
908            ELSE            ELSE
909              latentHeatFluxMaxMult(I,J,K) = ZERO              latentHeatFluxMaxMult(I,J,IT) = ZERO
910            ENDIF            ENDIF
911  #else  #else
912            IF (HEFFpreTH(I,J) .GT. ZERO) THEN            IF (HEFFpreTH(I,J) .GT. ZERO) THEN
# Line 1084  C      the temperature here is a result Line 1077  C      the temperature here is a result
1077  C      computed individually for each single category in SEAICE_SOLVE4TEMP  C      computed individually for each single category in SEAICE_SOLVE4TEMP
1078  C      and hence is averaged area weighted [areaFracFactor])  C      and hence is averaged area weighted [areaFracFactor])
1079             TICE(I,J,bi,bj) = TICE(I,J,bi,bj)             TICE(I,J,bi,bj) = TICE(I,J,bi,bj)
1080       &          +  ticeOutMult(I,J,IT)*areaFracFactor(I,J,K)       &          +  ticeOutMult(I,J,IT)*areaFracFactor(I,J,IT)
1081  #else  #else
1082             TICE(I,J,bi,bj) = TICE(I,J,bi,bj)             TICE(I,J,bi,bj) = TICE(I,J,bi,bj)
1083       &          +  ticeOutMult(I,J,IT)*recip_multDim       &          +  ticeOutMult(I,J,IT)*recip_multDim
# Line 1095  C     average over categories Line 1088  C     average over categories
1088  C     calculate area weighted mean  C     calculate area weighted mean
1089  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)
1090             a_QbyATM_cover   (I,J) = a_QbyATM_cover(I,J)             a_QbyATM_cover   (I,J) = a_QbyATM_cover(I,J)
1091       &          + a_QbyATMmult_cover(I,J,IT)*areaFracFactor(I,J,K)       &          + a_QbyATMmult_cover(I,J,IT)*areaFracFactor(I,J,IT)
1092             a_QSWbyATM_cover (I,J) = a_QSWbyATM_cover(I,J)             a_QSWbyATM_cover (I,J) = a_QSWbyATM_cover(I,J)
1093       &          + a_QSWbyATMmult_cover(I,J,IT)*areaFracFactor(I,J,K)       &          + a_QSWbyATMmult_cover(I,J,IT)*areaFracFactor(I,J,IT)
1094             a_FWbySublim     (I,J) = a_FWbySublim(I,J)             a_FWbySublim     (I,J) = a_FWbySublim(I,J)
1095       &          + a_FWbySublimMult(I,J,IT)*areaFracFactor(I,J,K)       &          + a_FWbySublimMult(I,J,IT)*areaFracFactor(I,J,IT)
1096  #else  #else
1097             a_QbyATM_cover   (I,J) = a_QbyATM_cover(I,J)             a_QbyATM_cover   (I,J) = a_QbyATM_cover(I,J)
1098       &          + 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 1134  CADJ STORE a_FWbySublim    = comlev1_bib
1134    
1135  C switch heat fluxes from W/m2 to 'effective' ice meters  C switch heat fluxes from W/m2 to 'effective' ice meters
1136  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1137          DO K=1,nITD          DO IT=1,nITD
1138           DO J=1,sNy           DO J=1,sNy
1139            DO I=1,sNx            DO I=1,sNx
1140             a_QbyATMmult_cover(I,J,K)   = a_QbyATMmult_cover(I,J,K)             a_QbyATMmult_cover(I,J,IT)   = a_QbyATMmult_cover(I,J,IT)
1141       &          * convertQ2HI * AREAITDpreTH(I,J,K)       &          * convertQ2HI * AREAITDpreTH(I,J,IT)
1142             a_QSWbyATMmult_cover(I,J,K) = a_QSWbyATMmult_cover(I,J,K)             a_QSWbyATMmult_cover(I,J,IT) = a_QSWbyATMmult_cover(I,J,IT)
1143       &          * convertQ2HI * AREAITDpreTH(I,J,K)       &          * convertQ2HI * AREAITDpreTH(I,J,IT)
1144  C and initialize r_QbyATM_cover  C and initialize r_QbyATM_cover
1145             r_QbyATMmult_cover(I,J,K)=a_QbyATMmult_cover(I,J,K)             r_QbyATMmult_cover(I,J,IT)=a_QbyATMmult_cover(I,J,IT)
1146  C     Convert fresh water flux by sublimation to 'effective' ice meters.  C     Convert fresh water flux by sublimation to 'effective' ice meters.
1147  C     Negative sublimation is resublimation and will be added as snow.  C     Negative sublimation is resublimation and will be added as snow.
1148  #ifdef SEAICE_DISABLE_SUBLIM  #ifdef SEAICE_DISABLE_SUBLIM
1149             a_FWbySublimMult(I,J,K) = ZERO             a_FWbySublimMult(I,J,IT) = ZERO
1150  #endif  #endif
1151             a_FWbySublimMult(I,J,K) = SEAICE_deltaTtherm*recip_rhoIce             a_FWbySublimMult(I,J,IT) = SEAICE_deltaTtherm*recip_rhoIce
1152       &            * a_FWbySublimMult(I,J,K)*AREAITDpreTH(I,J,K)       &            * a_FWbySublimMult(I,J,IT)*AREAITDpreTH(I,J,IT)
1153             r_FWbySublimMult(I,J,K)=a_FWbySublimMult(I,J,K)             r_FWbySublimMult(I,J,IT)=a_FWbySublimMult(I,J,IT)
1154            ENDDO            ENDDO
1155           ENDDO           ENDDO
1156          ENDDO          ENDDO
1157          DO J=1,sNy          DO J=1,sNy
1158           DO I=1,sNx           DO I=1,sNx
1159              a_QbyATM_open(I,J)    = a_QbyATM_open(I,J)
1160         &         * convertQ2HI * ( ONE - AREApreTH(I,J) )
1161              a_QSWbyATM_open(I,J)  = a_QSWbyATM_open(I,J)
1162         &         * convertQ2HI * ( ONE - AREApreTH(I,J) )
1163  C and initialize r_QbyATM_open  C and initialize r_QbyATM_open
1164            r_QbyATM_open(I,J)=a_QbyATM_open(I,J)            r_QbyATM_open(I,J)=a_QbyATM_open(I,J)
1165           ENDDO           ENDDO
# Line 1210  CADJ STORE r_FWbySublim    = comlev1_bib Line 1207  CADJ STORE r_FWbySublim    = comlev1_bib
1207  Cgf no additional dependency through ice cover  Cgf no additional dependency through ice cover
1208          IF ( SEAICEadjMODE.GE.3 ) THEN          IF ( SEAICEadjMODE.GE.3 ) THEN
1209  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1210           DO K=1,nITD           DO IT=1,nITD
1211            DO J=1,sNy            DO J=1,sNy
1212             DO I=1,sNx             DO I=1,sNx
1213              a_QbyATMmult_cover(I,J,K)   = 0. _d 0              a_QbyATMmult_cover(I,J,IT)   = 0. _d 0
1214              r_QbyATMmult_cover(I,J,K)   = 0. _d 0              r_QbyATMmult_cover(I,J,IT)   = 0. _d 0
1215              a_QSWbyATMmult_cover(I,J,K) = 0. _d 0              a_QSWbyATMmult_cover(I,J,IT) = 0. _d 0
1216             ENDDO             ENDDO
1217            ENDDO            ENDDO
1218           ENDDO           ENDDO
# Line 1272  c available turbulent flux Line 1269  c available turbulent flux
1269            a_QbyOCN(i,j) =            a_QbyOCN(i,j) =
1270       &         tmpscal1 * tmpscal2 * MixedLayerTurbulenceFactor       &         tmpscal1 * tmpscal2 * MixedLayerTurbulenceFactor
1271            r_QbyOCN(i,j) = a_QbyOCN(i,j)            r_QbyOCN(i,j) = a_QbyOCN(i,j)
 ctm  
           if (i.eq.20 .and. j.eq.20) then  
             print *, HeatCapacity_Cp  
             print *, rhoConst  
             print *, recip_QI  
             print *, theta(20,20,kSurface,bi,bj)  
             print *, tempFrz  
             print *, SEAICE_deltaTtherm  
             print *, maskC(20,20,kSurface,bi,bj)  
             print *, tmpscal2  
             print *, a_QbyOCN(20,20)  
           endif  
 ctm  
1272           ENDDO           ENDDO
1273          ENDDO          ENDDO
1274    
# Line 1305  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi Line 1289  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi
1289  CADJ STORE r_FWbySublim     = comlev1_bibj,key=iicekey,byte=isbyte  CADJ STORE r_FWbySublim     = comlev1_bibj,key=iicekey,byte=isbyte
1290  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1291  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1292          DO K=1,nITD          DO IT=1,nITD
1293  #endif  #endif
1294          DO J=1,sNy          DO J=1,sNy
1295           DO I=1,sNx           DO I=1,sNx
1296  C     First sublimate/deposite snow  C     First sublimate/deposite snow
1297            tmpscal2 =            tmpscal2 =
1298  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1299       &     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)
1300       &             *SNOW2ICE),ZERO)       &             *SNOW2ICE),ZERO)
1301  C         accumulate change over ITD categories  C         accumulate change over ITD categories
1302            d_HSNWbySublim(I,J)     = d_HSNWbySublim(I,J)     - tmpscal2            d_HSNWbySublim(I,J)     = d_HSNWbySublim(I,J)      - tmpscal2
1303       &                                                       *ICE2SNOW       &                                                       *ICE2SNOW
1304            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
1305       &                                                       *ICE2SNOW       &                                                       *ICE2SNOW
1306            r_FWbySublimMult(I,J,K) = r_FWbySublimMult(I,J,K) - tmpscal2            r_FWbySublimMult(I,J,IT)= r_FWbySublimMult(I,J,IT) - tmpscal2
1307  C         keep total up to date, too  C         keep total up to date, too
1308            r_FWbySublim(I,J)       = r_FWbySublim(I,J)       - tmpscal2            r_FWbySublim(I,J)       = r_FWbySublim(I,J)        - tmpscal2
1309  #else  #else
1310       &     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)
1311            d_HSNWbySublim(I,J) = - tmpscal2 * ICE2SNOW            d_HSNWbySublim(I,J) = - tmpscal2 * ICE2SNOW
# Line 1339  CADJ STORE r_FWbySublim    = comlev1_bib Line 1323  CADJ STORE r_FWbySublim    = comlev1_bib
1323  C     If anything is left, sublimate ice  C     If anything is left, sublimate ice
1324            tmpscal2 =            tmpscal2 =
1325  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1326       &     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)
1327  C         accumulate change over ITD categories  C         accumulate change over ITD categories
1328            d_HSNWbySublim(I,J)     = d_HSNWbySublim(I,J)     - tmpscal2            d_HSNWbySublim(I,J)     = d_HSNWbySublim(I,J)     - tmpscal2
1329            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
1330            r_FWbySublimMult(I,J,K) = r_FWbySublimMult(I,J,K) - tmpscal2            r_FWbySublimMult(I,J,IT) = r_FWbySublimMult(I,J,IT) - tmpscal2
1331  C         keep total up to date, too  C         keep total up to date, too
1332            r_FWbySublim(I,J)       = r_FWbySublim(I,J)       - tmpscal2            r_FWbySublim(I,J)       = r_FWbySublim(I,J)       - tmpscal2
1333  #else  #else
# Line 1360  C     If anything is left, it will be ev Line 1344  C     If anything is left, it will be ev
1344  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
1345  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).
1346  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1347            a_QbyATMmult_cover(I,J,K) = a_QbyATMmult_cover(I,J,K)            a_QbyATMmult_cover(I,J,IT) = a_QbyATMmult_cover(I,J,IT)
1348       &                              - r_FWbySublimMult(I,J,K)       &                               - r_FWbySublimMult(I,J,IT)
1349            r_QbyATMmult_cover(I,J,K) = r_QbyATMmult_cover(I,J,K)            r_QbyATMmult_cover(I,J,IT) = r_QbyATMmult_cover(I,J,IT)
1350       &                              - r_FWbySublimMult(I,J,K)       &                               - r_FWbySublimMult(I,J,IT)
1351           ENDDO           ENDDO
1352          ENDDO          ENDDO
1353  C       end K loop  C       end IT loop
1354          ENDDO          ENDDO
1355  C       then update totals        C       then update totals      
1356          DO J=1,sNy          DO J=1,sNy
# Line 1378  C       then update totals Line 1362  C       then update totals
1362          ENDDO          ENDDO
1363  c ToM<<< debug seaice_growth  c ToM<<< debug seaice_growth
1364          WRITE(msgBuf,'(A,7F6.2)')          WRITE(msgBuf,'(A,7F6.2)')
1365    #ifdef SEAICE_ITD
1366       &    ' SEAICE_GROWTH: Heff increments 1, HEFFITD = ',       &    ' SEAICE_GROWTH: Heff increments 1, HEFFITD = ',
1367       &     HEFFITD(20,20,:,bi,bj)       &     HEFFITD(20,20,:,bi,bj)
1368    #else
1369         &    ' SEAICE_GROWTH: Heff increments 1, HEFF = ',
1370         &     HEFF(20,20,bi,bj)
1371    #endif
1372          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1373       &    SQUEEZE_RIGHT , myThid)       &    SQUEEZE_RIGHT , myThid)
1374  c ToM>>>  c ToM>>>
# Line 1393  CADJ STORE r_QbyOCN = comlev1_bibj,key=i Line 1382  CADJ STORE r_QbyOCN = comlev1_bibj,key=i
1382  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1383    
1384  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1385          DO K=1,nITD          DO IT=1,nITD
1386           DO J=1,sNy           DO J=1,sNy
1387            DO I=1,sNx            DO I=1,sNx
1388  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
1389  C          and hence weighted by fractional area of each thickness category  C          and hence weighted by fractional area of each thickness category
1390             tmpscal1=MAX(r_QbyOCN(i,j)*areaFracFactor(I,J,K),             tmpscal1=MAX(r_QbyOCN(i,j)*areaFracFactor(I,J,IT),
1391       &                               -HEFFITD(I,J,K,bi,bj))       &                               -HEFFITD(I,J,IT,bi,bj))
1392             d_HEFFbyOCNonICE(I,J)= d_HEFFbyOCNonICE(I,J) + tmpscal1             d_HEFFbyOCNonICE(I,J)= d_HEFFbyOCNonICE(I,J) + tmpscal1
1393             r_QbyOCN(I,J)        = r_QbyOCN(I,J)         - tmpscal1             r_QbyOCN(I,J)        = r_QbyOCN(I,J)         - tmpscal1
1394             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
1395  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
1396             SItrHEFF(I,J,bi,bj,2) = SItrHEFF(I,J,bi,bj,2)             SItrHEFF(I,J,bi,bj,2) = SItrHEFF(I,J,bi,bj,2)
1397       &                           + HEFFITD(I,J,K,bi,bj)       &                           + HEFFITD(I,J,IT,bi,bj)
1398  #endif  #endif
1399            ENDDO            ENDDO
1400           ENDDO           ENDDO
1401          ENDDO          ENDDO
 c ToM<<< debug seaice_growth  
         WRITE(msgBuf,'(A,7F9.6)')  
      &    ' SEAICE_GROWTH: d_HEFFbyOCNonICE w/ITD: ',  
      &     d_HEFFbyOCNonICE(20,20)  
         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,  
      &    SQUEEZE_RIGHT , myThid)  
 c ToM>>>  
1402  #else /* SEAICE_ITD */  #else /* SEAICE_ITD */
1403          DO J=1,sNy          DO J=1,sNy
1404           DO I=1,sNx           DO I=1,sNx
# Line 1428  c ToM>>> Line 1410  c ToM>>>
1410  #endif  #endif
1411           ENDDO           ENDDO
1412          ENDDO          ENDDO
 c ToM<<< debug seaice_growth  
         WRITE(msgBuf,'(A,7F9.6)')  
      &    ' SEAICE_GROWTH: d_HEFFbyOCNonICE w/o ITD: ',  
      &     d_HEFFbyOCNonICE(20,20)  
         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,  
      &    SQUEEZE_RIGHT , myThid)  
 c ToM>>>  
1413  #endif /* SEAICE_ITD */  #endif /* SEAICE_ITD */
1414  c ToM<<< debug seaice_growth  c ToM<<< debug seaice_growth
1415          WRITE(msgBuf,'(A,7F6.2)')          WRITE(msgBuf,'(A,7F6.2)')
1416    #ifdef SEAICE_ITD
1417       &    ' SEAICE_GROWTH: Heff increments 2, HEFFITD = ',       &    ' SEAICE_GROWTH: Heff increments 2, HEFFITD = ',
1418       &     HEFFITD(20,20,:,bi,bj)       &     HEFFITD(20,20,:,bi,bj)
1419    #else
1420         &    ' SEAICE_GROWTH: Heff increments 2, HEFF = ',
1421         &     HEFF(20,20,bi,bj)
1422    #endif
1423          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1424       &    SQUEEZE_RIGHT , myThid)       &    SQUEEZE_RIGHT , myThid)
1425  c ToM>>>  c ToM>>>
# Line 1453  CADJ STORE r_QbyATM_cover = comlev1_bibj Line 1433  CADJ STORE r_QbyATM_cover = comlev1_bibj
1433  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1434    
1435  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1436          DO K=1,nITD          DO IT=1,nITD
1437           DO J=1,sNy           DO J=1,sNy
1438            DO I=1,sNx            DO I=1,sNx
1439  C     Convert to standard units (meters of ice) rather than to meters  C     Convert to standard units (meters of ice) rather than to meters
1440  C     of snow. This appears to be more robust.  C     of snow. This appears to be more robust.
1441            tmpscal1=MAX(r_QbyATMmult_cover(I,J,K),-HSNOWITD(I,J,K,bi,bj)             tmpscal1=MAX(r_QbyATMmult_cover(I,J,IT),
1442       &                                           *SNOW2ICE)       &                  -HSNOWITD(I,J,IT,bi,bj)*SNOW2ICE)
1443            tmpscal2=MIN(tmpscal1,0. _d 0)             tmpscal2=MIN(tmpscal1,0. _d 0)
1444  #ifdef SEAICE_MODIFY_GROWTH_ADJ  #ifdef SEAICE_MODIFY_GROWTH_ADJ
1445  Cgf no additional dependency through snow  Cgf no additional dependency through snow
1446            IF ( SEAICEadjMODE.GE.2 ) tmpscal2 = 0. _d 0             IF ( SEAICEadjMODE.GE.2 ) tmpscal2 = 0. _d 0
1447  #endif  #endif
1448            d_HSNWbyATMonSNW(I,J)=d_HSNWbyATMonSNW(I,J)+tmpscal2*ICE2SNOW             d_HSNWbyATMonSNW(I,J) = d_HSNWbyATMonSNW(I,J)
1449            HSNOWITD(I,J,K,bi,bj)=HSNOWITD(I,J,K,bi,bj)+tmpscal2*ICE2SNOW       &                           + tmpscal2*ICE2SNOW
1450            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)
1451  C         keep the total up to date, too       &                           + tmpscal2*ICE2SNOW
1452            r_QbyATM_cover(I,J)=r_QbyATM_cover(I,J) - tmpscal2             r_QbyATMmult_cover(I,J,IT)=r_QbyATMmult_cover(I,J,IT)
1453         &                           - tmpscal2
1454    C          keep the total up to date, too
1455               r_QbyATM_cover(I,J)=r_QbyATM_cover(I,J) - tmpscal2
1456            ENDDO            ENDDO
1457           ENDDO           ENDDO
1458          ENDDO          ENDDO
# Line 1492  Cgf no additional dependency through sno Line 1475  Cgf no additional dependency through sno
1475  #endif /* SEAICE_ITD */  #endif /* SEAICE_ITD */
1476  c ToM<<< debug seaice_growth  c ToM<<< debug seaice_growth
1477          WRITE(msgBuf,'(A,7F6.2)')          WRITE(msgBuf,'(A,7F6.2)')
1478    #ifdef SEAICE_ITD
1479       &    ' SEAICE_GROWTH: Heff increments 3, HEFFITD = ',       &    ' SEAICE_GROWTH: Heff increments 3, HEFFITD = ',
1480       &     HEFFITD(20,20,:,bi,bj)       &     HEFFITD(20,20,:,bi,bj)
1481    #else
1482         &    ' SEAICE_GROWTH: Heff increments 3, HEFF = ',
1483         &     HEFF(20,20,bi,bj)
1484    #endif
1485          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1486       &    SQUEEZE_RIGHT , myThid)       &    SQUEEZE_RIGHT , myThid)
1487  c ToM>>>  c ToM>>>
# Line 1512  Cgf the v1.81=>v1.82 revision would chan Line 1500  Cgf the v1.81=>v1.82 revision would chan
1500  Cgf warming conditions, the lab_sea results were not changed.  Cgf warming conditions, the lab_sea results were not changed.
1501    
1502  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1503          DO K=1,nITD          DO IT=1,nITD
1504           DO J=1,sNy           DO J=1,sNy
1505            DO I=1,sNx            DO I=1,sNx
1506  #ifdef SEAICE_GROWTH_LEGACY  #ifdef SEAICE_GROWTH_LEGACY
1507             tmpscal2 =MAX(-HEFFITD(I,J,K,bi,bj),r_QbyATMmult_cover(I,J,K))             tmpscal2 = MAX(-HEFFITD(I,J,IT,bi,bj),
1508         &                     r_QbyATMmult_cover(I,J,IT))
1509  #else  #else
1510             tmpscal2 =MAX(-HEFFITD(I,J,K,bi,bj),r_QbyATMmult_cover(I,J,K)             tmpscal2 = MAX(-HEFFITD(I,J,IT,bi,bj),
1511         &                     r_QbyATMmult_cover(I,J,IT)
1512  c         Limit ice growth by potential melt by ocean  c         Limit ice growth by potential melt by ocean
1513       &     + AREAITDpreTH(I,J,K) * r_QbyOCN(I,J))       &              + AREAITDpreTH(I,J,IT) * r_QbyOCN(I,J))
1514  #endif /* SEAICE_GROWTH_LEGACY */  #endif /* SEAICE_GROWTH_LEGACY */
1515             d_HEFFbyATMonOCN_cover(I,J) = d_HEFFbyATMonOCN_cover(I,J)             d_HEFFbyATMonOCN_cover(I,J) = d_HEFFbyATMonOCN_cover(I,J)
1516       &                                 + tmpscal2       &                                 + tmpscal2
# Line 1528  c         Limit ice growth by potential Line 1518  c         Limit ice growth by potential
1518       &                                 + tmpscal2       &                                 + tmpscal2
1519             r_QbyATM_cover(I,J)         = r_QbyATM_cover(I,J)             r_QbyATM_cover(I,J)         = r_QbyATM_cover(I,J)
1520       &                                 - tmpscal2       &                                 - tmpscal2
1521             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
1522    
1523  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
1524             SItrHEFF(I,J,bi,bj,3) = SItrHEFF(I,J,bi,bj,3)             SItrHEFF(I,J,bi,bj,3) = SItrHEFF(I,J,bi,bj,3)
1525       &                           + HEFFITD(I,J,K,bi,bj)       &                           + HEFFITD(I,J,IT,bi,bj)
1526  #endif  #endif
1527            ENDDO            ENDDO
1528           ENDDO           ENDDO
# Line 1562  c         Limit ice growth by potential Line 1552  c         Limit ice growth by potential
1552  #endif /* SEAICE_ITD */  #endif /* SEAICE_ITD */
1553  c ToM<<< debug seaice_growth  c ToM<<< debug seaice_growth
1554          WRITE(msgBuf,'(A,7F6.2)')          WRITE(msgBuf,'(A,7F6.2)')
1555    #ifdef SEAICE_ITD
1556       &    ' SEAICE_GROWTH: Heff increments 4, HEFFITD = ',       &    ' SEAICE_GROWTH: Heff increments 4, HEFFITD = ',
1557       &     HEFFITD(20,20,:,bi,bj)       &     HEFFITD(20,20,:,bi,bj)
1558    #else
1559         &    ' SEAICE_GROWTH: Heff increments 4, HEFF = ',
1560         &     HEFF(20,20,bi,bj)
1561    #endif
1562          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1563       &    SQUEEZE_RIGHT , myThid)       &    SQUEEZE_RIGHT , myThid)
1564  c ToM>>>  c ToM>>>
# Line 1596  C           add precip to the fresh wate Line 1591  C           add precip to the fresh wate
1591           ENDDO           ENDDO
1592          ENDDO          ENDDO
1593  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1594          DO K=1,nITD          DO IT=1,nITD
1595           DO J=1,sNy           DO J=1,sNy
1596            DO I=1,sNx            DO I=1,sNx
1597             HSNOWITD(I,J,K,bi,bj) = HSNOWITD(I,J,K,bi,bj)             HSNOWITD(I,J,IT,bi,bj) = HSNOWITD(I,J,IT,bi,bj)
1598       &     + d_HSNWbyRAIN(I,J)*areaFracFactor(I,J,K)       &     + d_HSNWbyRAIN(I,J)*areaFracFactor(I,J,IT)
1599            ENDDO            ENDDO
1600           ENDDO           ENDDO
1601          ENDDO          ENDDO
# Line 1617  Cgf rain to snow is not a surface proces Line 1612  Cgf rain to snow is not a surface proces
1612  #endif /* ALLOW_ATM_TEMP */  #endif /* ALLOW_ATM_TEMP */
1613  c ToM<<< debug seaice_growth  c ToM<<< debug seaice_growth
1614          WRITE(msgBuf,'(A,7F6.2)')          WRITE(msgBuf,'(A,7F6.2)')
1615    #ifdef SEAICE_ITD
1616       &    ' SEAICE_GROWTH: Heff increments 5, HEFFITD = ',       &    ' SEAICE_GROWTH: Heff increments 5, HEFFITD = ',
1617       &     HEFFITD(20,20,:,bi,bj)       &     HEFFITD(20,20,:,bi,bj)
1618    #else
1619         &    ' SEAICE_GROWTH: Heff increments 5, HEFF = ',
1620         &     HEFF(20,20,bi,bj)
1621    #endif
1622          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1623       &    SQUEEZE_RIGHT , myThid)       &    SQUEEZE_RIGHT , myThid)
1624  c ToM>>>  c ToM>>>
# Line 1635  CADJ STORE r_QbyOCN = comlev1_bibj,key=i Line 1635  CADJ STORE r_QbyOCN = comlev1_bibj,key=i
1635  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1636    
1637  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1638          DO K=1,nITD          DO IT=1,nITD
1639           DO J=1,sNy           DO J=1,sNy
1640            DO I=1,sNx            DO I=1,sNx
1641             tmpscal1=MAX(r_QbyOCN(i,j)*ICE2SNOW*areaFracFactor(I,J,K),             tmpscal1=MAX(r_QbyOCN(i,j)*ICE2SNOW*areaFracFactor(I,J,IT),
1642       &                  -HSNOWITD(I,J,K,bi,bj))       &                  -HSNOWITD(I,J,IT,bi,bj))
1643             tmpscal2=MIN(tmpscal1,0. _d 0)             tmpscal2=MIN(tmpscal1,0. _d 0)
1644  #ifdef SEAICE_MODIFY_GROWTH_ADJ  #ifdef SEAICE_MODIFY_GROWTH_ADJ
1645  Cgf no additional dependency through snow  Cgf no additional dependency through snow
# Line 1647  Cgf no additional dependency through sno Line 1647  Cgf no additional dependency through sno
1647  #endif  #endif
1648             d_HSNWbyOCNonSNW(I,J) = d_HSNWbyOCNonSNW(I,J) + tmpscal2             d_HSNWbyOCNonSNW(I,J) = d_HSNWbyOCNonSNW(I,J) + tmpscal2
1649             r_QbyOCN(I,J)=r_QbyOCN(I,J) - tmpscal2*SNOW2ICE             r_QbyOCN(I,J)=r_QbyOCN(I,J) - tmpscal2*SNOW2ICE
1650             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
1651            ENDDO            ENDDO
1652           ENDDO           ENDDO
1653          ENDDO          ENDDO
# Line 1671  Cgf no additional dependency through sno Line 1671  Cgf no additional dependency through sno
1671  Cph)  Cph)
1672  c ToM<<< debug seaice_growth  c ToM<<< debug seaice_growth
1673          WRITE(msgBuf,'(A,7F6.2)')          WRITE(msgBuf,'(A,7F6.2)')
1674    #ifdef SEAICE_ITD
1675       &    ' SEAICE_GROWTH: Heff increments 6, HEFFITD = ',       &    ' SEAICE_GROWTH: Heff increments 6, HEFFITD = ',
1676       &     HEFFITD(20,20,:,bi,bj)       &     HEFFITD(20,20,:,bi,bj)
1677    #else
1678         &    ' SEAICE_GROWTH: Heff increments 6, HEFF = ',
1679         &     HEFF(20,20,bi,bj)
1680    #endif
1681          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1682       &    SQUEEZE_RIGHT , myThid)       &    SQUEEZE_RIGHT , myThid)
1683  c ToM>>>  c ToM>>>
# Line 1710  C         determine thickness of new ice Line 1715  C         determine thickness of new ice
1715  C         considering the entire open water area to refreeze  C         considering the entire open water area to refreeze
1716            tmpscal1 = tmpscal3/tmpscal0            tmpscal1 = tmpscal3/tmpscal0
1717  C         then add new ice volume to appropriate thickness category  C         then add new ice volume to appropriate thickness category
1718            DO K=1,nITD            DO IT=1,nITD
1719             IF (tmpscal1.LT.Hlimit(K)) THEN             IF (tmpscal1.LT.Hlimit(IT)) THEN
1720              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
1721              tmpscal3=ZERO              tmpscal3=ZERO
1722  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
1723  C   in PART 4 below in non-itd code  C   in PART 4 below in non-itd code
1724  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,
1725  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
1726  C   Hibler's h_0 parameter  C   Hibler's h_0 parameter
1727              AREAITD(I,J,K,bi,bj) = AREAITD(I,J,K,bi,bj)              AREAITD(I,J,IT,bi,bj) = AREAITD(I,J,IT,bi,bj)
1728       &                           + tmpscal0       &                           + tmpscal0
1729              tmpscal0=ZERO              tmpscal0=ZERO
1730             ENDIF             ENDIF
# Line 1732  C   Hibler's h_0 parameter Line 1737  C   Hibler's h_0 parameter
1737    
1738  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
1739  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1740          DO K=1,nITD          DO IT=1,nITD
1741           DO J=1,sNy           DO J=1,sNy
1742            DO I=1,sNx            DO I=1,sNx
1743  c needs to be here to allow use also with LEGACY branch  c needs to be here to allow use also with LEGACY branch
1744             SItrHEFF(I,J,bi,bj,4) = SItrHEFF(I,J,bi,bj,4)             SItrHEFF(I,J,bi,bj,4) = SItrHEFF(I,J,bi,bj,4)
1745       &                           + HEFFITD(I,J,K,bi,bj)       &                           + HEFFITD(I,J,IT,bi,bj)
1746            ENDDO            ENDDO
1747           ENDDO           ENDDO
1748          ENDDO          ENDDO
# Line 1752  c needs to be here to allow use also wit Line 1757  c needs to be here to allow use also wit
1757  #endif /* ALLOW_SITRACER */  #endif /* ALLOW_SITRACER */
1758  c ToM<<< debug seaice_growth  c ToM<<< debug seaice_growth
1759          WRITE(msgBuf,'(A,7F6.2)')          WRITE(msgBuf,'(A,7F6.2)')
1760    #ifdef SEAICE_ITD
1761       &    ' SEAICE_GROWTH: Heff increments 7, HEFFITD = ',       &    ' SEAICE_GROWTH: Heff increments 7, HEFFITD = ',
1762       &     HEFFITD(20,20,:,bi,bj)       &     HEFFITD(20,20,:,bi,bj)
1763    #else
1764         &    ' SEAICE_GROWTH: Heff increments 7, HEFF = ',
1765         &     HEFF(20,20,bi,bj)
1766    #endif
1767          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1768       &    SQUEEZE_RIGHT , myThid)       &    SQUEEZE_RIGHT , myThid)
1769  c ToM>>>  c ToM>>>
# Line 1769  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi Line 1779  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi
1779  #endif /* ALLOW_AUTODIFF_TAMC */  #endif /* ALLOW_AUTODIFF_TAMC */
1780          IF ( SEAICEuseFlooding ) THEN          IF ( SEAICEuseFlooding ) THEN
1781  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1782           DO K=1,nITD           DO IT=1,nITD
1783            DO J=1,sNy            DO J=1,sNy
1784             DO I=1,sNx             DO I=1,sNx
1785              tmpscal0 = (HSNOWITD(I,J,K,bi,bj)*SEAICE_rhoSnow              tmpscal0 = (HSNOWITD(I,J,IT,bi,bj)*SEAICE_rhoSnow
1786       &               +HEFFITD(I,J,K,bi,bj)*SEAICE_rhoIce)*recip_rhoConst       &               +  HEFFITD(I,J,IT,bi,bj) *SEAICE_rhoIce)
1787              tmpscal1 = MAX( 0. _d 0, tmpscal0 - HEFFITD(I,J,K,bi,bj))       &                                        *recip_rhoConst
1788              d_HEFFbyFLOODING(I,J) = d_HEFFbyFLOODING(I,J) + tmpscal1              tmpscal1 = MAX( 0. _d 0, tmpscal0 - HEFFITD(I,J,IT,bi,bj))
1789              HEFFITD(I,J,K,bi,bj)  = HEFFITD(I,J,K,bi,bj)  + tmpscal1              d_HEFFbyFLOODING(I,J) = d_HEFFbyFLOODING(I,J)  + tmpscal1
1790              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
1791                HSNOWITD(I,J,IT,bi,bj)= HSNOWITD(I,J,IT,bi,bj) - tmpscal1
1792       &                            * ICE2SNOW       &                            * ICE2SNOW
1793             ENDDO             ENDDO
1794            ENDDO            ENDDO
# Line 1799  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi Line 1810  CADJ STORE hsnow(:,:,bi,bj) = comlev1_bi
1810  #endif /* SEAICE_GROWTH_LEGACY */  #endif /* SEAICE_GROWTH_LEGACY */
1811  c ToM<<< debug seaice_growth  c ToM<<< debug seaice_growth
1812          WRITE(msgBuf,'(A,7F6.2)')          WRITE(msgBuf,'(A,7F6.2)')
1813    #ifdef SEAICE_ITD
1814       &    ' SEAICE_GROWTH: Heff increments 8, HEFFITD = ',       &    ' SEAICE_GROWTH: Heff increments 8, HEFFITD = ',
1815       &     HEFFITD(20,20,:,bi,bj)       &     HEFFITD(20,20,:,bi,bj)
1816    #else
1817         &    ' SEAICE_GROWTH: Heff increments 8, HEFF = ',
1818         &     HEFF(20,20,bi,bj)
1819    #endif
1820          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
1821       &    SQUEEZE_RIGHT , myThid)       &    SQUEEZE_RIGHT , myThid)
1822  c ToM>>>  c ToM>>>
# Line 1842  C       does not account for lateral mel Line 1858  C       does not account for lateral mel
1858  C  C
1859  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
1860  C  C
1861          DO K=1,nITD          DO IT=1,nITD
1862           DO J=1,sNy           DO J=1,sNy
1863            DO I=1,sNx            DO I=1,sNx
1864             IF (HEFFITD(I,J,K,bi,bj).LE.ZERO) THEN             IF (HEFFITD(I,J,IT,bi,bj).LE.ZERO) THEN
1865              AREAITD(I,J,K,bi,bj)=ZERO              AREAITD(I,J,IT,bi,bj)=ZERO
1866             ENDIF             ENDIF
1867  #ifdef ALLOW_SITRACER  #ifdef ALLOW_SITRACER
1868             SItrAREA(I,J,bi,bj,3) = SItrAREA(I,J,bi,bj,3)             SItrAREA(I,J,bi,bj,3) = SItrAREA(I,J,bi,bj,3)
1869       &                           + AREAITD(I,J,K,bi,bj)       &                           + AREAITD(I,J,IT,bi,bj)
1870  #endif /* ALLOW_SITRACER */  #endif /* ALLOW_SITRACER */
1871            ENDDO            ENDDO
1872           ENDDO           ENDDO
# Line 1931  C apply tendency Line 1947  C apply tendency
1947  Cgf 'bulk' linearization of area=f(HEFF)  Cgf 'bulk' linearization of area=f(HEFF)
1948          IF ( SEAICEadjMODE.GE.1 ) THEN          IF ( SEAICEadjMODE.GE.1 ) THEN
1949  #ifdef SEAICE_ITD  #ifdef SEAICE_ITD
1950           DO K=1,nITD           DO IT=1,nITD
1951            DO J=1,sNy            DO J=1,sNy
1952             DO I=1,sNx             DO I=1,sNx
1953              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 *
1954       &               ( HEFFITD(I,J,K,bi,bj) - HEFFITDpreTH(I,J,K) )       &               ( HEFFITD(I,J,IT,bi,bj) - HEFFITDpreTH(I,J,IT) )
1955             ENDDO             ENDDO
1956            ENDDO            ENDDO
1957           ENDDO           ENDDO
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.5
  ViewVC Help
Powered by ViewVC 1.1.22