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 8327 for branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd_da.F90 – NEMO

Ignore:
Timestamp:
2017-07-13T11:29:29+02:00 (7 years ago)
Author:
clem
Message:

STEP4 (3): put all thermodynamics in 1D (limthd_da OK)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd_da.F90

    r8325 r8327  
    1212   !!   lim_thd_da   : sea ice lateral melting 
    1313   !!---------------------------------------------------------------------- 
    14    USE par_oce                ! ocean parameters 
    15    USE phycst                 ! physical constants (ocean directory) 
    16    USE sbc_oce, ONLY: sst_m   ! Surface boundary condition: ocean fields 
    17    USE ice                    ! LIM variables 
    18    USE lib_mpp                ! MPP library 
    19    USE wrk_nemo               ! work arrays 
    20    USE lib_fortran            ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
     14   USE par_oce        ! ocean parameters 
     15   USE phycst         ! physical constants (ocean directory) 
     16   USE sbc_oce , ONLY : sst_m 
     17   USE ice            ! LIM variables 
     18   USE thd_ice        ! thermodynamic sea-ice variables 
     19   USE limtab         ! 1D <==> 2D transformation 
     20   ! 
     21   USE lib_mpp        ! MPP library 
     22   USE wrk_nemo       ! work arrays 
     23   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
    2124 
    2225   IMPLICIT NONE 
     
    97100      !!              Phil. Trans. R. Soc. A, 373(2052), 20140167. 
    98101      !!--------------------------------------------------------------------- 
    99       INTEGER             ::   ji, jj, jl      ! dummy loop indices 
     102      INTEGER             ::   ji, jj, jk, jl     ! dummy loop indices 
     103      INTEGER             ::   nidx 
    100104      REAL(wp)            ::   zastar, zdfloe, zperi, zwlat, zda 
    101105      REAL(wp), PARAMETER ::   zdmax = 300._wp 
     
    104108      REAL(wp), PARAMETER ::   zm2   = 1.36_wp 
    105109      ! 
    106       REAL(wp), POINTER, DIMENSION(:,:) ::   zda_tot 
     110      REAL(wp), DIMENSION(jpij) ::   zda_tot 
    107111      !!--------------------------------------------------------------------- 
    108       CALL wrk_alloc( jpi,jpj, zda_tot ) 
     112 
     113      ! select ice covered grid points 
     114      nidx = 0 ; idxice(:) = 0 
     115      DO jj = 2, jpjm1 
     116         DO ji = 2, jpim1 
     117            IF ( at_i(ji,jj) > epsi10 ) THEN 
     118               nidx         = nidx  + 1 
     119               idxice(nidx) = (jj - 1) * jpi + ji 
     120            ENDIF 
     121         END DO 
     122      END DO 
    109123 
    110124      !------------------------------------------------------------! 
     
    112126      !------------------------------------------------------------! 
    113127      zastar = 1._wp / ( 1._wp - (rn_dmin / zdmax)**(1._wp/rn_beta) ) 
    114        
    115       DO jj = 1, jpj 
    116          DO ji = 1, jpi 
    117              
    118             ! Mean floe caliper diameter [m] 
    119             zdfloe = rn_dmin * ( zastar / ( zastar - at_i(ji,jj) ) )**rn_beta 
    120              
    121             ! Mean perimeter of the floe = N*pi*D = (A/cs*D^2)*pi*D [m.m-2] 
    122             zperi = at_i(ji,jj) * rpi / ( zcs * zdfloe ) 
    123              
    124             ! Melt speed rate [m/s] 
    125             zwlat = zm1 * ( MAX( 0._wp, sst_m(ji,jj) - ( t_bo(ji,jj) - rt0 ) ) )**zm2 
    126              
    127             ! sea ice concentration decrease 
    128             zda_tot(ji,jj) = - MIN( zwlat * zperi * rdt_ice, at_i(ji,jj) ) 
    129              
    130          END DO 
     128 
     129      CALL tab_2d_1d( nidx, at_i_1d(1:nidx), at_i , jpi, jpj, idxice(1:nidx) ) 
     130      CALL tab_2d_1d( nidx, t_bo_1d(1:nidx), t_bo , jpi, jpj, idxice(1:nidx) ) 
     131      CALL tab_2d_1d( nidx, sst_1d (1:nidx), sst_m, jpi, jpj, idxice(1:nidx) ) 
     132      DO ji = 1, nidx    
     133         zdfloe = rn_dmin * ( zastar / ( zastar - at_i_1d(ji) ) )**rn_beta         ! Mean floe caliper diameter [m] 
     134         zperi  = at_i_1d(ji) * rpi / ( zcs * zdfloe )                             ! Mean perimeter of the floe = N*pi*D = (A/cs*D^2)*pi*D [m.m-2] 
     135         zwlat  = zm1 * ( MAX( 0._wp, sst_1d(ji) - ( t_bo_1d(ji) - rt0 ) ) )**zm2  ! Melt speed rate [m/s] 
     136          
     137         zda_tot(ji) = - MIN( zwlat * zperi * rdt_ice, at_i_1d(ji) )               ! sea ice concentration decrease 
    131138      END DO 
    132139       
     
    134141      ! --- Distribute reduction among ice categories and calculate associated ice-ocean fluxes --- ! 
    135142      !---------------------------------------------------------------------------------------------! 
    136       DO jl = jpl, 1, -1 
    137          DO jj = 1, jpj 
    138             DO ji = 1, jpi 
    139                 
    140                ! decrease of concentration for the category jl 
    141                !    1st option: each category contributes to melting in proportion to its concentration 
    142                rswitch = MAX( 0._wp , SIGN( 1._wp, at_i(ji,jj) - epsi10 ) ) 
    143                zda     = rswitch * zda_tot(ji,jj) * a_i(ji,jj,jl) / MAX( at_i(ji,jj), epsi10 ) 
    144                !    2d option: melting of the upper cat first 
    145                !!zda = MAX( zda_tot(ji,jj), - a_i(ji,jj,jl) ) 
    146                !!zda_tot(ji,jj) = zda_tot(ji,jj) + zda 
    147                 
    148                ! Contribution to salt flux 
    149                sfx_lam(ji,jj) = sfx_lam(ji,jj) - rhoic *  ht_i(ji,jj,jl) * zda * sm_i(ji,jj,jl) * r1_rdtice 
    150                 
    151                ! Contribution to heat flux into the ocean [W.m-2], <0   
    152 !clemX               hfx_thd(ji,jj) = hfx_thd(ji,jj) + zda * r1_rdtice * ( ht_i(ji,jj,jl) * SUM( e_i(ji,jj,:,jl) ) * r1_nlay_i  & 
    153 !                  &                                                + ht_s(ji,jj,jl) *      e_s(ji,jj,1,jl)   * r1_nlay_s ) 
    154                hfx_thd(ji,jj) = hfx_thd(ji,jj) + rswitch * zda_tot(ji,jj) / MAX( at_i(ji,jj), epsi10 )  & 
    155                   &                                      * r1_rdtice * ( SUM( e_i(ji,jj,:,jl) ) + e_s(ji,jj,1,jl) ) 
    156                 
    157                ! Contribution to mass flux 
    158                wfx_lam(ji,jj) =  wfx_lam(ji,jj) - zda * r1_rdtice * ( rhoic * ht_i(ji,jj,jl) + rhosn * ht_s(ji,jj,jl) ) 
    159                 
    160                ! new concentration 
    161                a_i(ji,jj,jl) = a_i(ji,jj,jl) + zda 
    162             END DO 
    163          END DO 
     143      DO jl = 1, jpl 
     144         CALL tab_2d_1d( nidx, a_i_1d    (1:nidx), a_i(:,:,jl) , jpi, jpj, idxice(1:nidx) ) 
     145         CALL tab_2d_1d( nidx, ht_i_1d   (1:nidx), ht_i(:,:,jl), jpi, jpj, idxice(1:nidx) ) 
     146         CALL tab_2d_1d( nidx, sm_i_1d   (1:nidx), sm_i(:,:,jl), jpi, jpj, idxice(1:nidx) ) 
     147         CALL tab_2d_1d( nidx, sfx_lam_1d(1:nidx), sfx_lam     , jpi, jpj, idxice(1:nidx) ) 
     148         CALL tab_2d_1d( nidx, hfx_thd_1d(1:nidx), hfx_thd     , jpi, jpj, idxice(1:nidx) ) 
     149         CALL tab_2d_1d( nidx, wfx_lam_1d(1:nidx), wfx_lam     , jpi, jpj, idxice(1:nidx) ) 
     150         DO jk = 1, nlay_i 
     151            CALL tab_2d_1d( nidx, e_i_1d(1:nidx,jk), e_i(:,:,jk,jl)  , jpi, jpj, idxice(1:nidx) ) 
     152         END DO 
     153         DO jk = 1, nlay_s 
     154            CALL tab_2d_1d( nidx, e_s_1d(1:nidx,jk), e_s(:,:,jk,jl)  , jpi, jpj, idxice(1:nidx) ) 
     155         END DO 
     156 
     157         DO ji = 1, nidx 
     158            ! decrease of concentration for the category jl 
     159            !    each category contributes to melting in proportion to its concentration 
     160            zda     = zda_tot(ji) * a_i_1d(ji) / at_i_1d(ji) 
     161             
     162            ! Contribution to salt flux 
     163            sfx_lam_1d(ji) = sfx_lam_1d(ji) - rhoic *  ht_i_1d(ji) * zda * sm_i_1d(ji) * r1_rdtice 
     164             
     165            ! Contribution to heat flux into the ocean [W.m-2], <0   
     166            !clemX               hfx_thd_1d(ji) = hfx_thd_1d(ji) + zda * r1_rdtice * ( ht_i_1d(ji) * SUM( e_i_1d(ji,:) ) * r1_nlay_i  & 
     167            !                  &                                                     + ht_s_1d(ji) *      e_s_1d(ji,1)   * r1_nlay_s ) 
     168            hfx_thd_1d(ji) = hfx_thd_1d(ji) + zda_tot(ji) / at_i_1d(ji) * r1_rdtice * ( SUM( e_i_1d(ji,:) ) + e_s_1d(ji,1) ) 
     169             
     170            ! Contribution to mass flux 
     171            wfx_lam_1d(ji) =  wfx_lam_1d(ji) - zda * r1_rdtice * ( rhoic * ht_i_1d(ji) + rhosn * ht_s_1d(ji) ) 
     172             
     173            ! new concentration 
     174            a_i_1d(ji) = a_i_1d(ji) + zda 
     175 
     176            ! ensure that ht_i = 0 where a_i = 0 
     177            IF( a_i_1d(ji) == 0._wp )   ht_i_1d(ji) = 0._wp   
     178         END DO 
     179 
     180         CALL tab_1d_2d( nidx, a_i (:,:,jl), idxice, a_i_1d     (1:nidx), jpi, jpj ) 
     181         CALL tab_1d_2d( nidx, ht_i(:,:,jl), idxice, ht_i_1d    (1:nidx), jpi, jpj ) 
     182         CALL tab_1d_2d( nidx, sfx_lam     , idxice, sfx_lam_1d(1:nidx) , jpi, jpj ) 
     183         CALL tab_1d_2d( nidx, hfx_thd     , idxice, hfx_thd_1d(1:nidx) , jpi, jpj ) 
     184         CALL tab_1d_2d( nidx, wfx_lam     , idxice, wfx_lam_1d(1:nidx) , jpi, jpj ) 
     185 
    164186      END DO 
    165        
    166       ! total concentration 
    167       at_i(:,:) = SUM( a_i(:,:,:), dim=3 ) 
    168        
    169       ! --- ensure that ht_i = 0 where a_i = 0 --- 
    170       WHERE( a_i == 0._wp ) ht_i = 0._wp 
    171       ! 
    172       CALL wrk_dealloc( jpi,jpj, zda_tot ) 
     187             
    173188      ! 
    174189   END SUBROUTINE lim_thd_da 
Note: See TracChangeset for help on using the changeset viewer.