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

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

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

revision 1.29 by mlosch, Tue Feb 7 14:01:10 2012 UTC revision 1.30 by gforget, Sat Feb 11 03:35:01 2012 UTC
# Line 11  C     !ROUTINE: SEAICE_SOLVE4TEMP Line 11  C     !ROUTINE: SEAICE_SOLVE4TEMP
11  C     !INTERFACE:  C     !INTERFACE:
12        SUBROUTINE SEAICE_SOLVE4TEMP(        SUBROUTINE SEAICE_SOLVE4TEMP(
13       I   UG, HICE_ACTUAL, HSNOW_ACTUAL,       I   UG, HICE_ACTUAL, HSNOW_ACTUAL,
14  #ifdef SEAICE_ADD_SUBLIMATION_TO_FWBUDGET  #ifdef SEAICE_CAP_SUBLIM
15       I   F_lh_max,       I   F_lh_max,
16  #endif  #endif
17       U   TSURF,       U   TSURF,
# Line 39  C     === Global variables === Line 39  C     === Global variables ===
39  #include "SEAICE_SIZE.h"  #include "SEAICE_SIZE.h"
40  #include "SEAICE_PARAMS.h"  #include "SEAICE_PARAMS.h"
41  #include "SEAICE.h"  #include "SEAICE.h"
 #ifdef SEAICE_VARIABLE_FREEZING_POINT  
42  #include "DYNVARS.h"  #include "DYNVARS.h"
 #endif /* SEAICE_VARIABLE_FREEZING_POINT */  
43  #ifdef ALLOW_EXF  #ifdef ALLOW_EXF
44  # include "EXF_FIELDS.h"  # include "EXF_FIELDS.h"
45  #endif  #endif
# Line 72  C---------- Line 70  C----------
70        _RL UG          (1:sNx,1:sNy)        _RL UG          (1:sNx,1:sNy)
71        _RL HICE_ACTUAL (1:sNx,1:sNy)        _RL HICE_ACTUAL (1:sNx,1:sNy)
72        _RL HSNOW_ACTUAL(1:sNx,1:sNy)        _RL HSNOW_ACTUAL(1:sNx,1:sNy)
73  #ifdef SEAICE_ADD_SUBLIMATION_TO_FWBUDGET  #ifdef SEAICE_CAP_SUBLIM
74        _RL F_lh_max    (1:sNx,1:sNy)        _RL F_lh_max    (1:sNx,1:sNy)
75  #endif  #endif
76        _RL TSURF       (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)        _RL TSURF       (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
# Line 88  CEOP Line 86  CEOP
86  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
87  C     === Local variables ===  C     === Local variables ===
88  C     i, j  :: Loop counters  C     i, j  :: Loop counters
89  C     kSrf  :: vertical index of surface layer  C     kSurface  :: vertical index of surface layer
90        INTEGER i, j        INTEGER i, j
91  #ifdef SEAICE_VARIABLE_FREEZING_POINT        INTEGER kSurface
       INTEGER kSrf  
 #endif /* SEAICE_VARIABLE_FREEZING_POINT */  
92        INTEGER ITER        INTEGER ITER
93  C     TB :: ocean temperature in contact with ice (=seawater freezing point) (K)  C     tempFrz :: ocean temperature in contact with ice (=seawater freezing point) (K)
94        _RL TB         (1:sNx,1:sNy)        _RL tempFrz         (1:sNx,1:sNy)
95        _RL D1, D1I        _RL D1, D1I
96        _RL D3(1:sNx,1:sNy)        _RL D3(1:sNx,1:sNy)
97        _RL TMELT, XKI, XKS, HCUT, XIO        _RL TMELT, XKI, XKS, HCUT, XIO
# Line 178  C     cc0 = TEN ** aa2 Line 174  C     cc0 = TEN ** aa2
174        cc1 = cc0*aa1*bb1*Ppascals*lnTEN        cc1 = cc0*aa1*bb1*Ppascals*lnTEN
175        cc2 = cc0*bb2        cc2 = cc0*bb2
176    
177  #ifdef SEAICE_VARIABLE_FREEZING_POINT        IF ( buoyancyRelation .EQ. 'OCEANICP' ) THEN
178        kSrf = 1         kSurface        = Nr
179  #endif /* SEAICE_VARIABLE_FREEZING_POINT */        ELSE
180           kSurface        = 1
181          ENDIF
182    
183  C     SENSIBLE HEAT CONSTANT  C     SENSIBLE HEAT CONSTANT
184        D1=SEAICE_dalton*SEAICE_cpAir*SEAICE_rhoAir        D1=SEAICE_dalton*SEAICE_cpAir*SEAICE_rhoAir
# Line 233  c       tsurfLoc (I,J) = MIN( celsius2K+ Line 231  c       tsurfLoc (I,J) = MIN( celsius2K+
231          lwdownLoc(I,J) = MAX( MIN_LWDOWN, LWDOWN(I,J,bi,bj) )          lwdownLoc(I,J) = MAX( MIN_LWDOWN, LWDOWN(I,J,bi,bj) )
232          atempLoc (I,J) = MAX( celsius2K+MIN_ATEMP, ATEMP(I,J,bi,bj) )          atempLoc (I,J) = MAX( celsius2K+MIN_ATEMP, ATEMP(I,J,bi,bj) )
233    
234  C     FREEZING TEMPERATURE OF SEAWATER  c     FREEZING TEMP. OF SEA WATER (K)
235  #ifdef SEAICE_VARIABLE_FREEZING_POINT          tempFrz(I,J) = SEAICE_dTempFrz_dS *salt(I,J,kSurface,bi,bj)
236  C     Use a variable seawater freezing point       &     + SEAICE_tempFrz0 + celsius2K
237          TB(I,J) = -0.0575 _d 0*salt(I,J,kSrf,bi,bj) + 0.0901 _d 0      
      &       + celsius2K  
 #else  
 C     Use a constant freezing temperature (SEAICE_VARIABLE_FREEZING_POINT undef)  
 C     old SOLVE4TEMP_LEGACY setting (not consistent with seaice_growth value)  
 c       TB(I,J) = 271.2 _d 0  
         TB(I,J) = celsius2K + SEAICE_freeze  
 #endif /* SEAICE_VARIABLE_FREEZING_POINT */  
   
238  C     Now determine fixed (relative to tsurf) forcing term in heat budget  C     Now determine fixed (relative to tsurf) forcing term in heat budget
239    
240          IF(HSNOW_ACTUAL(I,J).GT.0.0) THEN          IF(HSNOW_ACTUAL(I,J).GT.0.0) THEN
# Line 396  C     d(qh)/d(TICE) Line 386  C     d(qh)/d(TICE)
386            ENDIF            ENDIF
387    
388  C     Calculate the flux terms based on the updated tsurfLoc  C     Calculate the flux terms based on the updated tsurfLoc
389            F_c(I,J)    = effConduct(I,J)*(TB(I,J)-tsurfLoc(I,J))            F_c(I,J)    = effConduct(I,J)*(tempFrz(I,J)-tsurfLoc(I,J))
390            F_lh(I,J)   = D1I*UG(I,J)*(qhice(I,J)-AQH(I,J,bi,bj))            F_lh(I,J)   = D1I*UG(I,J)*(qhice(I,J)-AQH(I,J,bi,bj))
391  #ifdef SEAICE_ADD_SUBLIMATION_TO_FWBUDGET  #ifdef SEAICE_CAP_SUBLIM
392  C     if the latent heat flux implied by tsurfLoc exceeds  C     if the latent heat flux implied by tsurfLoc exceeds
393  C     F_lh_max, cap F_lh and decouple the flux magnitude from TICE  C     F_lh_max, cap F_lh and decouple the flux magnitude from TICE
394            IF (F_lh(I,J) .GT. F_lh_max(I,J)) THEN            IF (F_lh(I,J) .GT. F_lh_max(I,J)) THEN
395               F_lh(I,J)  = F_lh_max(I,J)               F_lh(I,J)  = F_lh_max(I,J)
396               dqh_dTs(I,J) = ZERO               dqh_dTs(I,J) = ZERO
397            ENDIF            ENDIF
398  #endif /* SEAICE_ADD_SUBLIMATION_TO_FWBUDGET */  #endif /* SEAICE_CAP_SUBLIM */
399    
400            F_lwu(I,J) = t4 * D3(I,J)            F_lwu(I,J) = t4 * D3(I,J)
401            F_sens(I,J)= D1 * UG(I,J) * (t1 - atempLoc(I,J))            F_sens(I,J)= D1 * UG(I,J) * (t1 - atempLoc(I,J))
# Line 495  C     The following form does the same, Line 485  C     The following form does the same,
485  C     over ice specific humidity  C     over ice specific humidity
486             qhice(I,J) = bb1*mm_pi/( Ppascals -(1.0 _d 0 - bb1)*mm_pi )             qhice(I,J) = bb1*mm_pi/( Ppascals -(1.0 _d 0 - bb1)*mm_pi )
487            ENDIF            ENDIF
488            F_c(I,J)  = effConduct(I,J) * (TB(I,J) - t1)            F_c(I,J)  = effConduct(I,J) * (tempFrz(I,J) - t1)
489            F_lh(I,J) = D1I * UG(I,J)*(qhice(I,J)-AQH(I,J,bi,bj))            F_lh(I,J) = D1I * UG(I,J)*(qhice(I,J)-AQH(I,J,bi,bj))
490  #ifdef SEAICE_ADD_SUBLIMATION_TO_FWBUDGET  #ifdef SEAICE_CAP_SUBLIM
491            IF (F_lh(I,J) .GT. F_lh_max(I,J)) THEN            IF (F_lh(I,J) .GT. F_lh_max(I,J)) THEN
492               F_lh(I,J)  = F_lh_max(I,J)               F_lh(I,J)  = F_lh_max(I,J)
493            ENDIF            ENDIF
494  #endif /* SEAICE_ADD_SUBLIMATION_TO_FWBUDGET */  #endif /* SEAICE_CAP_SUBLIM */
495            F_lwu(I,J)  = t4 * D3(I,J)            F_lwu(I,J)  = t4 * D3(I,J)
496            F_sens(I,J) = D1 * UG(I,J) * (t1 - atempLoc(I,J))            F_sens(I,J) = D1 * UG(I,J) * (t1 - atempLoc(I,J))
497  C     The flux between the ice/snow surface and the atmosphere.  C     The flux between the ice/snow surface and the atmosphere.
# Line 511  C     The flux between the ice/snow surf Line 501  C     The flux between the ice/snow surf
501           ELSEIF ( postSolvTempIter.EQ.1 ) THEN           ELSEIF ( postSolvTempIter.EQ.1 ) THEN
502  C     Update fluxes (consistent with the linearized formulation)  C     Update fluxes (consistent with the linearized formulation)
503            delTsurf  = tsurfLoc(I,J)-tsurfPrev(I,J)            delTsurf  = tsurfLoc(I,J)-tsurfPrev(I,J)
504            F_c(I,J)  = effConduct(I,J)*(TB(I,J)-tsurfLoc(I,J))            F_c(I,J)  = effConduct(I,J)*(tempFrz(I,J)-tsurfLoc(I,J))
505            F_ia(I,J) = F_ia(I,J) + dFia_dTs(I,J)*delTsurf            F_ia(I,J) = F_ia(I,J) + dFia_dTs(I,J)*delTsurf
506            F_lh(I,J) = F_lh(I,J)            F_lh(I,J) = F_lh(I,J)
507       &              + D1I*UG(I,J)*dqh_dTs(I,J)*delTsurf       &              + D1I*UG(I,J)*dqh_dTs(I,J)*delTsurf
# Line 557  C     Calculate the net ice-ocean and ic Line 547  C     Calculate the net ice-ocean and ic
547           print '(A,4(1x,D24.15))',           print '(A,4(1x,D24.15))',
548       &         'ibi F(c,lh,sens)           ',       &         'ibi F(c,lh,sens)           ',
549       &         F_c(I,J), F_lh(I,J), F_sens(I,J)       &         F_c(I,J), F_lh(I,J), F_sens(I,J)
550  #ifdef SEAICE_ADD_SUBLIMATION_TO_FWBUDGET  #ifdef SEAICE_CAP_SUBLIM
551           IF (F_lh_max(I,J) .GT. ZERO) THEN           IF (F_lh_max(I,J) .GT. ZERO) THEN
552               print '(A,4(1x,D24.15))',               print '(A,4(1x,D24.15))',
553       &         'ibi F_lh_max,  F_lh/lhmax) ',       &         'ibi F_lh_max,  F_lh/lhmax) ',
# Line 570  C     Calculate the net ice-ocean and ic Line 560  C     Calculate the net ice-ocean and ic
560       &         'ibi FWsub, FWsubm*dT/rhoI  ',       &         'ibi FWsub, FWsubm*dT/rhoI  ',
561       &          FWsublim(I,J),       &          FWsublim(I,J),
562       &          FWsublim(I,J)*SEAICE_deltaTtherm/SEAICE_rhoICE       &          FWsublim(I,J)*SEAICE_deltaTtherm/SEAICE_rhoICE
563  #endif /* SEAICE_ADD_SUBLIMATION_TO_FWBUDGET */  #endif /* SEAICE_CAP_SUBLIM */
564            print '(A,4(1x,D24.15))',            print '(A,4(1x,D24.15))',
565       &         'ibi F_ia, F_ia_net, F_c    ',       &         'ibi F_ia, F_ia_net, F_c    ',
566       &         F_ia(I,J), F_ia_net, F_c(I,J)       &         F_ia(I,J), F_ia_net, F_c(I,J)
Legend:
Removed from v.1.29  
changed lines
  Added in v.1.30
  ViewVC Help
Powered by ViewVC 1.1.22