New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 9977 for NEMO/branches/UKMO/dev_r9888_proto_GO8_package/src/ICE/icethd_dh.F90 – NEMO

Ignore:
Timestamp:
2018-07-20T10:24:45+02:00 (6 years ago)
Author:
davestorkey
Message:

UKMO/dev_r9888_proto_GO8_package branch: merge in changes from trunk to rev 9922.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/UKMO/dev_r9888_proto_GO8_package/src/ICE/icethd_dh.F90

    r9892 r9977  
    8585 
    8686      REAL(wp), DIMENSION(jpij) ::   zqprec      ! energy of fallen snow                       (J.m-3) 
    87       REAL(wp), DIMENSION(jpij) ::   zq_su       ! heat for surface ablation                   (J.m-2) 
    88       REAL(wp), DIMENSION(jpij) ::   zq_bo       ! heat for bottom ablation                    (J.m-2) 
     87      REAL(wp), DIMENSION(jpij) ::   zq_top      ! heat for surface ablation                   (J.m-2) 
     88      REAL(wp), DIMENSION(jpij) ::   zq_bot      ! heat for bottom ablation                    (J.m-2) 
    8989      REAL(wp), DIMENSION(jpij) ::   zq_rema     ! remaining heat at the end of the routine    (J.m-2) 
    9090      REAL(wp), DIMENSION(jpij) ::   zf_tt       ! Heat budget to determine melting or freezing(W.m-2) 
     
    131131         ! 
    132132         DO ji = 1, npti 
    133             zq_su(ji)      = MAX( 0._wp, qml_ice_1d(ji) * rdt_ice ) 
     133            zq_top(ji)     = MAX( 0._wp, qml_ice_1d(ji) * rdt_ice ) 
    134134         END DO 
    135135         ! 
    136       CASE( np_jules_EMULE ) 
     136      CASE( np_jules_OFF , np_jules_EMULE ) 
    137137         ! 
    138138         DO ji = 1, npti 
    139             zdum           = qns_ice_1d(ji) + qsr_ice_1d(ji) - qsr_ice_tr_1d(ji) - fc_su(ji) 
     139            zdum           = qns_ice_1d(ji) + qsr_ice_1d(ji) - qtr_ice_top_1d(ji) - qcn_ice_top_1d(ji) 
    140140            qml_ice_1d(ji) = zdum * MAX( 0._wp , SIGN( 1._wp, t_su_1d(ji) - rt0 ) ) 
    141             zq_su(ji)      = MAX( 0._wp, qml_ice_1d(ji) * rdt_ice ) 
    142          END DO 
    143          ! 
    144       CASE( np_jules_OFF )  
    145          ! 
    146          DO ji = 1, npti 
    147             zdum           = qns_ice_1d(ji) + qsr_ice_1d(ji) - qsr_ice_tr_1d(ji) - fc_su(ji)  
    148             qml_ice_1d(ji) = zdum * MAX( 0._wp , SIGN( 1._wp, t_su_1d(ji) - rt0 ) ) 
    149             zq_su(ji)      = MAX( 0._wp, qml_ice_1d(ji) * rdt_ice ) 
     141            zq_top(ji)     = MAX( 0._wp, qml_ice_1d(ji) * rdt_ice ) 
    150142         END DO 
    151143         ! 
     
    153145      ! 
    154146      DO ji = 1, npti 
    155          zf_tt(ji)         = fc_bo_i(ji) + fhtur_1d(ji) + fhld_1d(ji)  
    156          zq_bo(ji)         = MAX( 0._wp, zf_tt(ji) * rdt_ice ) 
     147         zf_tt(ji)         = qcn_ice_bot_1d(ji) + qsb_ice_bot_1d(ji) + fhld_1d(ji)  
     148         zq_bot(ji)        = MAX( 0._wp, zf_tt(ji) * rdt_ice ) 
    157149      END DO 
    158150 
     
    210202            ! --- melt of falling snow --- 
    211203            rswitch              = MAX( 0._wp , SIGN( 1._wp , zqprec(ji) - epsi20 ) ) 
    212             zdeltah       (ji,1) = - rswitch * zq_su(ji) / MAX( zqprec(ji) , epsi20 )   ! thickness change 
    213             zdeltah       (ji,1) = MAX( - zdh_s_pre(ji), zdeltah(ji,1) )                ! bound melting  
     204            zdeltah       (ji,1) = - rswitch * zq_top(ji) / MAX( zqprec(ji) , epsi20 )   ! thickness change 
     205            zdeltah       (ji,1) = MAX( - zdh_s_pre(ji), zdeltah(ji,1) )                 ! bound melting  
    214206            hfx_snw_1d    (ji)   = hfx_snw_1d    (ji) - zdeltah(ji,1) * a_i_1d(ji) * zqprec(ji)    * r1_rdtice   ! heat used to melt snow (W.m-2, >0) 
    215207            wfx_snw_sum_1d(ji)   = wfx_snw_sum_1d(ji) - rhosn         * a_i_1d(ji) * zdeltah(ji,1) * r1_rdtice   ! snow melting only = water into the ocean (then without snow precip), >0 
     
    217209            ! updates available heat + precipitations after melting 
    218210            dh_s_mlt (ji) = dh_s_mlt(ji) + zdeltah(ji,1) 
    219             zq_su    (ji) = MAX( 0._wp , zq_su (ji) + zdeltah(ji,1) * zqprec(ji) )       
     211            zq_top   (ji) = MAX( 0._wp , zq_top (ji) + zdeltah(ji,1) * zqprec(ji) )       
    220212            zdh_s_pre(ji) = zdh_s_pre(ji) + zdeltah(ji,1) 
    221213             
     
    240232      ! Snow melting 
    241233      ! ------------ 
    242       ! If heat still available (zq_su > 0), then melt more snow 
     234      ! If heat still available (zq_top > 0), then melt more snow 
    243235      zdeltah(1:npti,:) = 0._wp 
    244236      zdh_s_mel(1:npti) = 0._wp 
    245237      DO jk = 1, nlay_s 
    246238         DO ji = 1, npti 
    247             IF( zh_s(ji,jk) > 0._wp .AND. zq_su(ji) > 0._wp ) THEN 
     239            IF( zh_s(ji,jk) > 0._wp .AND. zq_top(ji) > 0._wp ) THEN 
    248240               ! 
    249241               rswitch          = MAX( 0._wp, SIGN( 1._wp, e_s_1d(ji,jk) - epsi20 ) ) 
    250                zdeltah  (ji,jk) = - rswitch * zq_su(ji) / MAX( e_s_1d(ji,jk), epsi20 )   ! thickness change 
    251                zdeltah  (ji,jk) = MAX( zdeltah(ji,jk) , - zh_s(ji,jk) )                  ! bound melting 
     242               zdeltah  (ji,jk) = - rswitch * zq_top(ji) / MAX( e_s_1d(ji,jk), epsi20 )   ! thickness change 
     243               zdeltah  (ji,jk) = MAX( zdeltah(ji,jk) , - zh_s(ji,jk) )                   ! bound melting 
    252244               zdh_s_mel(ji)    = zdh_s_mel(ji) + zdeltah(ji,jk) 
    253245                
     
    257249               ! updates available heat + thickness 
    258250               dh_s_mlt(ji)    = dh_s_mlt(ji) + zdeltah(ji,jk) 
    259                zq_su   (ji)    = MAX( 0._wp , zq_su (ji) + zdeltah(ji,jk) * e_s_1d(ji,jk) ) 
     251               zq_top  (ji)    = MAX( 0._wp , zq_top(ji) + zdeltah(ji,jk) * e_s_1d(ji,jk) ) 
    260252               h_s_1d  (ji)    = MAX( 0._wp , h_s_1d(ji) + zdeltah(ji,jk) ) 
    261253               zh_s    (ji,jk) = MAX( 0._wp , zh_s(ji,jk) + zdeltah(ji,jk) ) 
     
    349341               zdE            =    zEi - zEw                          ! Specific enthalpy difference < 0 
    350342                
    351                zfmdt          = - zq_su(ji) / zdE                     ! Mass flux to the ocean [kg/m2, >0] 
     343               zfmdt          = - zq_top(ji) / zdE                    ! Mass flux to the ocean [kg/m2, >0] 
    352344                
    353345               zdeltah(ji,jk) = - zfmdt * r1_rhoic                    ! Melt of layer jk [m, <0] 
     
    355347               zdeltah(ji,jk) = MIN( 0._wp , MAX( zdeltah(ji,jk) , - zh_i(ji,jk) ) )    ! Melt of layer jk cannot exceed the layer thickness [m, <0] 
    356348                
    357                zq_su(ji)      = MAX( 0._wp , zq_su(ji) - zdeltah(ji,jk) * rhoic * zdE ) ! update available heat 
     349               zq_top(ji)      = MAX( 0._wp , zq_top(ji) - zdeltah(ji,jk) * rhoic * zdE ) ! update available heat 
    358350                
    359351               dh_i_sum(ji)   = dh_i_sum(ji) + zdeltah(ji,jk)         ! Cumulate surface melt 
     
    416408      !------------------ 
    417409      ! Basal growth is driven by heat imbalance at the ice-ocean interface, 
    418       ! between the inner conductive flux  (fc_bo_i), from the open water heat flux  
    419       ! (fhld) and the turbulent ocean flux (fhtur).  
    420       ! fc_bo_i is positive downwards. fhtur and fhld are positive to the ice  
     410      ! between the inner conductive flux  (qcn_ice_bot), from the open water heat flux  
     411      ! (fhld) and the sensible ice-ocean flux (qsb_ice_bot).  
     412      ! qcn_ice_bot is positive downwards. qsb_ice_bot and fhld are positive to the ice  
    421413 
    422414      ! If salinity varies in time, an iterative procedure is required, because 
     
    515507                  zdE             = zEi - zEw                                                 ! Specific enthalpy difference   (J/kg, <0) 
    516508 
    517                   zfmdt           = - zq_bo(ji) / zdE                                         ! Mass flux x time step (kg/m2, >0) 
     509                  zfmdt           = - zq_bot(ji) / zdE                                        ! Mass flux x time step (kg/m2, >0) 
    518510 
    519511                  zdeltah(ji,jk)  = - zfmdt * r1_rhoic                                        ! Gross thickness change 
     
    521513                  zdeltah(ji,jk)  = MIN( 0._wp , MAX( zdeltah(ji,jk), - zh_i(ji,jk) ) )       ! bound thickness change 
    522514                   
    523                   zq_bo(ji)       = MAX( 0._wp , zq_bo(ji) - zdeltah(ji,jk) * rhoic * zdE )   ! update available heat. MAX is necessary for roundup errors 
    524  
    525                   dh_i_bom(ji)    = dh_i_bom(ji) + zdeltah(ji,jk)                            ! Update basal melt 
     515                  zq_bot(ji)      = MAX( 0._wp , zq_bot(ji) - zdeltah(ji,jk) * rhoic * zdE )  ! update available heat. MAX is necessary for roundup errors 
     516 
     517                  dh_i_bom(ji)    = dh_i_bom(ji) + zdeltah(ji,jk)                             ! Update basal melt 
    526518 
    527519                  zfmdt           = - zdeltah(ji,jk) * rhoic                                  ! Mass flux x time step > 0 
     
    556548      zdeltah(1:npti,:) = 0._wp ! important 
    557549      DO ji = 1, npti 
    558          zq_rema (ji)   = zq_su(ji) + zq_bo(ji)  
     550         zq_rema (ji)   = zq_top(ji) + zq_bot(ji)  
    559551         rswitch        = 1._wp - MAX( 0._wp, SIGN( 1._wp, - h_s_1d(ji) ) )   ! =1 if snow 
    560552         rswitch        = rswitch * MAX( 0._wp, SIGN( 1._wp, e_s_1d(ji,1) - epsi20 ) ) 
     
    570562         !     
    571563         ! Remaining heat flux (W.m-2) is sent to the ocean heat budget 
    572          hfx_out_1d(ji) = hfx_out_1d(ji) + ( zq_rema(ji) * a_i_1d(ji) ) * r1_rdtice 
     564         qt_oce_ai_1d(ji) = qt_oce_ai_1d(ji) + ( zq_rema(ji) * a_i_1d(ji) ) * r1_rdtice 
    573565 
    574566         IF( ln_icectl .AND. zq_rema(ji) < 0. .AND. lwp ) WRITE(numout,*) 'ALERTE zq_rema <0 = ', zq_rema(ji) 
Note: See TracChangeset for help on using the changeset viewer.