Ignore:
Timestamp:
2015-06-19T17:18:00+02:00 (5 years ago)
Author:
davestorkey
Message:

Update 2015/dev_r5021_UKMO1_CICE_coupling branch to revision 5442 of the trunk.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5021_UKMO1_CICE_coupling/NEMOGCM/NEMO/LIM_SRC_3/limupdate2.F90

    • Property svn:keywords set to Id
    r5234 r5443  
    55   !!====================================================================== 
    66   !! History :  3.0  !  2006-04  (M. Vancoppenolle) Original code 
    7    !!            3.6  !  2014-06  (C. Rousset)       Complete rewriting/cleaning 
     7   !!            3.5  !  2014-06  (C. Rousset)       Complete rewriting/cleaning 
    88   !!---------------------------------------------------------------------- 
    99#if defined key_lim3 
     
    1313   !!    lim_update2   : computes update of sea-ice global variables from trend terms 
    1414   !!---------------------------------------------------------------------- 
    15    USE limrhg          ! ice rheology 
    16  
    17    USE dom_oce 
    18    USE oce             ! dynamics and tracers variables 
    19    USE in_out_manager 
    2015   USE sbc_oce         ! Surface boundary condition: ocean fields 
    2116   USE sbc_ice         ! Surface boundary condition: ice fields 
    2217   USE dom_ice 
     18   USE dom_oce 
    2319   USE phycst          ! physical constants 
    2420   USE ice 
    25    USE limdyn 
    26    USE limtrp 
    27    USE limthd 
    28    USE limsbc 
    29    USE limdiahsb 
    30    USE limwri 
    31    USE limrst 
    3221   USE thd_ice         ! LIM thermodynamic sea-ice variables 
    33    USE par_ice 
    3422   USE limitd_th 
    35    USE limitd_me 
    3623   USE limvar 
    37    USE prtctl           ! Print control 
    38    USE lbclnk           ! lateral boundary condition - MPP exchanges 
    39    USE wrk_nemo         ! work arrays 
    40    USE lib_fortran     ! glob_sum 
     24   USE prtctl          ! Print control 
     25   USE lbclnk          ! lateral boundary condition - MPP exchanges 
     26   USE wrk_nemo        ! work arrays 
    4127   USE timing          ! Timing 
    42    USE limcons        ! conservation tests 
     28   USE limcons         ! conservation tests 
     29   USE limctl 
     30   USE lib_mpp         ! MPP library 
     31   USE lib_fortran     ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
     32   USE in_out_manager 
    4333 
    4434   IMPLICIT NONE 
     
    5646CONTAINS 
    5747 
    58    SUBROUTINE lim_update2 
     48   SUBROUTINE lim_update2( kt ) 
    5949      !!------------------------------------------------------------------- 
    6050      !!               ***  ROUTINE lim_update2  *** 
     
    6454      !! 
    6555      !!--------------------------------------------------------------------- 
    66       INTEGER  ::   ji, jj, jk, jl    ! dummy loop indices 
    67       INTEGER  ::   i_ice_switch 
    68       REAL(wp) ::   zh, zsal 
    69       ! 
    70       REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b  
     56      INTEGER, INTENT(in) ::   kt    ! number of iteration 
     57      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices 
     58      REAL(wp) ::   zsal 
     59      REAL(wp) ::   zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b  
    7160      !!------------------------------------------------------------------- 
    7261      IF( nn_timing == 1 )  CALL timing_start('limupdate2') 
    7362 
     63      IF( kt == nit000 .AND. lwp ) THEN 
     64         WRITE(numout,*) ' lim_update2 ' 
     65         WRITE(numout,*) ' ~~~~~~~~~~~ ' 
     66      ENDIF 
     67 
    7468      ! conservation test 
    7569      IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limupdate2', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    7670 
    77       !----------------- 
    78       ! zap small values 
    79       !----------------- 
    80       CALL lim_itd_me_zapsmall 
    81  
    82       CALL lim_var_glo2eqv 
    83  
    84       !---------------------------------------------------- 
    85       ! Rebin categories with thickness out of bounds 
    86       !---------------------------------------------------- 
    87       IF ( jpl > 1 )   CALL lim_itd_th_reb(1, jpl) 
    88  
    8971      !---------------------------------------------------------------------- 
    90       ! Constrain the thickness of the smallest category above hiclim 
     72      ! Constrain the thickness of the smallest category above himin 
    9173      !---------------------------------------------------------------------- 
    9274      DO jj = 1, jpj  
    9375         DO ji = 1, jpi 
    94             IF( v_i(ji,jj,1) > 0._wp .AND. ht_i(ji,jj,1) < hiclim ) THEN 
    95                zh             = hiclim / ht_i(ji,jj,1) 
    96                ht_s(ji,jj,1) = ht_s(ji,jj,1) * zh 
    97                ht_i(ji,jj,1) = ht_i(ji,jj,1) * zh 
    98                a_i (ji,jj,1) = a_i(ji,jj,1)  / zh 
     76            rswitch = MAX( 0._wp , SIGN( 1._wp, a_i(ji,jj,1) - epsi20 ) )   !0 if no ice and 1 if yes 
     77            ht_i(ji,jj,1) = v_i (ji,jj,1) / MAX( a_i(ji,jj,1) , epsi20 ) * rswitch 
     78            IF( v_i(ji,jj,1) > 0._wp .AND. ht_i(ji,jj,1) < rn_himin ) THEN 
     79               a_i (ji,jj,1) = a_i (ji,jj,1) * ht_i(ji,jj,1) / rn_himin 
     80               oa_i(ji,jj,1) = oa_i(ji,jj,1) * ht_i(ji,jj,1) / rn_himin 
    9981            ENDIF 
    10082         END DO 
     
    11294         DO jj = 1, jpj 
    11395            DO ji = 1, jpi 
    114                IF( at_i(ji,jj) > amax .AND. a_i(ji,jj,jl) > 0._wp ) THEN 
    115                   a_i(ji,jj,jl)  = a_i(ji,jj,jl) * ( 1._wp - ( 1._wp - amax / at_i(ji,jj) ) ) 
    116                   ht_i(ji,jj,jl) = v_i(ji,jj,jl) / a_i(ji,jj,jl) 
     96               IF( at_i(ji,jj) > rn_amax .AND. a_i(ji,jj,jl) > 0._wp ) THEN 
     97                  a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax / at_i(ji,jj) ) ) 
     98                  oa_i(ji,jj,jl) = oa_i(ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax / at_i(ji,jj) ) ) 
    11799               ENDIF 
    118100            END DO 
     
    120102      END DO 
    121103 
    122       at_i(:,:) = 0.0 
    123       DO jl = 1, jpl 
    124          at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 
    125       END DO 
    126  
    127       ! -------------------------------------- 
    128       ! Final thickness distribution rebinning 
    129       ! -------------------------------------- 
    130       IF ( jpl > 1 ) CALL lim_itd_th_reb( 1, jpl ) 
    131  
    132       !----------------- 
    133       ! zap small values 
    134       !----------------- 
    135       CALL lim_itd_me_zapsmall 
    136  
    137104      !--------------------- 
    138       ! 2.11) Ice salinity 
     105      ! Ice salinity 
    139106      !--------------------- 
    140       IF (  num_sal == 2  ) THEN  
     107      IF (  nn_icesal == 2  ) THEN  
    141108         DO jl = 1, jpl 
    142109            DO jj = 1, jpj  
    143110               DO ji = 1, jpi 
    144111                  zsal            = smv_i(ji,jj,jl) 
    145                   smv_i(ji,jj,jl) = sm_i(ji,jj,jl) * v_i(ji,jj,jl) 
    146112                  ! salinity stays in bounds 
    147                   i_ice_switch    = 1._wp - MAX( 0._wp, SIGN( 1._wp, - v_i(ji,jj,jl) ) ) 
    148                   smv_i(ji,jj,jl) = i_ice_switch * MAX( MIN( s_i_max * v_i(ji,jj,jl), smv_i(ji,jj,jl) ), s_i_min * v_i(ji,jj,jl) ) !+ s_i_min * ( 1._wp - i_ice_switch ) * v_i(ji,jj,jl) 
     113                  rswitch         = 1._wp - MAX( 0._wp, SIGN( 1._wp, - v_i(ji,jj,jl) ) ) 
     114                  smv_i(ji,jj,jl) = rswitch * MAX( MIN( rn_simax * v_i(ji,jj,jl), smv_i(ji,jj,jl) ), rn_simin * v_i(ji,jj,jl) ) 
    149115                  ! associated salt flux 
    150116                  sfx_res(ji,jj) = sfx_res(ji,jj) - ( smv_i(ji,jj,jl) - zsal ) * rhoic * r1_rdtice 
    151                END DO ! ji 
    152             END DO ! jj 
    153          END DO !jl 
     117               END DO 
     118            END DO 
     119         END DO 
    154120      ENDIF 
    155121 
     122      !---------------------------------------------------- 
     123      ! Rebin categories with thickness out of bounds 
     124      !---------------------------------------------------- 
     125      IF ( jpl > 1 ) CALL lim_itd_th_reb( 1, jpl ) 
     126 
     127      !----------------- 
     128      ! zap small values 
     129      !----------------- 
     130      CALL lim_var_zapsmall 
     131 
    156132      !------------------------------------------------------------------------------ 
    157       ! 2) Corrections to avoid wrong values                                        | 
     133      ! Corrections to avoid wrong values                                        | 
    158134      !------------------------------------------------------------------------------ 
    159135      ! Ice drift 
     
    173149      CALL lbc_lnk( v_ice(:,:), 'V', -1. ) 
    174150      !mask velocities 
    175       u_ice(:,:) = u_ice(:,:) * tmu(:,:) 
    176       v_ice(:,:) = v_ice(:,:) * tmv(:,:) 
     151      u_ice(:,:) = u_ice(:,:) * umask(:,:,1) 
     152      v_ice(:,:) = v_ice(:,:) * vmask(:,:,1) 
    177153  
    178154      ! ------------------------------------------------- 
    179155      ! Diagnostics 
    180156      ! ------------------------------------------------- 
    181       d_a_i_thd(:,:,:)   = a_i(:,:,:)   - a_i_b(:,:,:)  
    182       d_v_s_thd(:,:,:)   = v_s(:,:,:)   - v_s_b(:,:,:) 
    183       d_v_i_thd(:,:,:)   = v_i(:,:,:)   - v_i_b(:,:,:)   
    184       d_e_s_thd(:,:,:,:) = e_s(:,:,:,:) - e_s_b(:,:,:,:)  
    185       d_e_i_thd(:,:,1:nlay_i,:) = e_i(:,:,1:nlay_i,:) - e_i_b(:,:,1:nlay_i,:) 
    186       !?? d_oa_i_thd(:,:,:)  = oa_i (:,:,:) - oa_i_b (:,:,:) 
    187       d_smv_i_thd(:,:,:) = 0._wp 
    188       IF( num_sal == 2 )   d_smv_i_thd(:,:,:) = smv_i(:,:,:) - smv_i_b(:,:,:) 
    189       ! diag only (clem) 
    190       dv_dt_thd(:,:,:) = d_v_i_thd(:,:,:) * r1_rdtice * rday 
    191  
    192       ! heat content variation (W.m-2) 
     157      DO jl  = 1, jpl 
     158         oa_i(:,:,jl) = oa_i(:,:,jl) + a_i(:,:,jl) * rdt_ice / rday   ! ice natural aging 
     159         afx_thd(:,:) = afx_thd(:,:) + ( a_i(:,:,jl) - a_i_b(:,:,jl) ) * r1_rdtice 
     160      END DO 
     161      afx_tot = afx_thd + afx_dyn 
     162 
    193163      DO jj = 1, jpj 
    194164         DO ji = 1, jpi             
    195             diag_heat_dhc(ji,jj) = ( SUM( d_e_i_trp(ji,jj,1:nlay_i,:) + d_e_i_thd(ji,jj,1:nlay_i,:) ) +  &  
    196                &                     SUM( d_e_s_trp(ji,jj,1:nlay_s,:) + d_e_s_thd(ji,jj,1:nlay_s,:) )    & 
    197                &                   ) * unit_fac * r1_rdtice / area(ji,jj)    
     165            ! heat content variation (W.m-2) 
     166            diag_heat(ji,jj) = diag_heat(ji,jj) -  & 
     167               &               ( SUM( e_i(ji,jj,1:nlay_i,:) - e_i_b(ji,jj,1:nlay_i,:) ) +  &  
     168               &                 SUM( e_s(ji,jj,1:nlay_s,:) - e_s_b(ji,jj,1:nlay_s,:) )    & 
     169               &               ) * r1_rdtice    
     170            ! salt, volume 
     171            diag_smvi(ji,jj) = diag_smvi(ji,jj) + SUM( smv_i(ji,jj,:) - smv_i_b(ji,jj,:) ) * rhoic * r1_rdtice 
     172            diag_vice(ji,jj) = diag_vice(ji,jj) + SUM( v_i  (ji,jj,:) - v_i_b  (ji,jj,:) ) * rhoic * r1_rdtice 
     173            diag_vsnw(ji,jj) = diag_vsnw(ji,jj) + SUM( v_s  (ji,jj,:) - v_s_b  (ji,jj,:) ) * rhosn * r1_rdtice 
    198174         END DO 
    199175      END DO 
    200176 
    201177      ! conservation test 
    202       IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limupdate2', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
     178      IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limupdate2', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
     179 
     180      ! necessary calls (at least for coupling) 
     181      CALL lim_var_glo2eqv 
     182      CALL lim_var_agg(2) 
     183 
     184      ! ------------------------------------------------- 
     185      ! control prints 
     186      ! ------------------------------------------------- 
     187      IF( ln_icectl )   CALL lim_prt( kt, iiceprt, jiceprt, 2, ' - Final state - ' )   ! control print 
    203188 
    204189      IF(ln_ctl) THEN   ! Control print 
     
    206191         CALL prt_ctl_info(' - Cell values : ') 
    207192         CALL prt_ctl_info('   ~~~~~~~~~~~~~ ') 
    208          CALL prt_ctl(tab2d_1=area       , clinfo1=' lim_update2  : cell area   :') 
     193         CALL prt_ctl(tab2d_1=e12t       , clinfo1=' lim_update2  : cell area   :') 
    209194         CALL prt_ctl(tab2d_1=at_i       , clinfo1=' lim_update2  : at_i        :') 
    210195         CALL prt_ctl(tab2d_1=vt_i       , clinfo1=' lim_update2  : vt_i        :') 
     
    226211            CALL prt_ctl(tab2d_1=a_i        (:,:,jl)        , clinfo1= ' lim_update2  : a_i         : ') 
    227212            CALL prt_ctl(tab2d_1=a_i_b      (:,:,jl)        , clinfo1= ' lim_update2  : a_i_b       : ') 
    228             CALL prt_ctl(tab2d_1=d_a_i_thd  (:,:,jl)        , clinfo1= ' lim_update2  : d_a_i_thd   : ') 
    229213            CALL prt_ctl(tab2d_1=v_i        (:,:,jl)        , clinfo1= ' lim_update2  : v_i         : ') 
    230214            CALL prt_ctl(tab2d_1=v_i_b      (:,:,jl)        , clinfo1= ' lim_update2  : v_i_b       : ') 
    231             CALL prt_ctl(tab2d_1=d_v_i_thd  (:,:,jl)        , clinfo1= ' lim_update2  : d_v_i_thd   : ') 
    232215            CALL prt_ctl(tab2d_1=v_s        (:,:,jl)        , clinfo1= ' lim_update2  : v_s         : ') 
    233216            CALL prt_ctl(tab2d_1=v_s_b      (:,:,jl)        , clinfo1= ' lim_update2  : v_s_b       : ') 
    234             CALL prt_ctl(tab2d_1=d_v_s_thd  (:,:,jl)        , clinfo1= ' lim_update2  : d_v_s_thd   : ') 
    235             CALL prt_ctl(tab2d_1=e_i        (:,:,1,jl)/1.0e9, clinfo1= ' lim_update2  : e_i1        : ') 
    236             CALL prt_ctl(tab2d_1=e_i_b      (:,:,1,jl)/1.0e9, clinfo1= ' lim_update2  : e_i1_b      : ') 
    237             CALL prt_ctl(tab2d_1=d_e_i_thd  (:,:,1,jl)/1.0e9, clinfo1= ' lim_update2  : de_i1_thd   : ') 
    238             CALL prt_ctl(tab2d_1=e_i        (:,:,2,jl)/1.0e9, clinfo1= ' lim_update2  : e_i2        : ') 
    239             CALL prt_ctl(tab2d_1=e_i_b      (:,:,2,jl)/1.0e9, clinfo1= ' lim_update2  : e_i2_b      : ') 
    240             CALL prt_ctl(tab2d_1=d_e_i_thd  (:,:,2,jl)/1.0e9, clinfo1= ' lim_update2  : de_i2_thd   : ') 
     217            CALL prt_ctl(tab2d_1=e_i        (:,:,1,jl)      , clinfo1= ' lim_update2  : e_i1        : ') 
     218            CALL prt_ctl(tab2d_1=e_i_b      (:,:,1,jl)      , clinfo1= ' lim_update2  : e_i1_b      : ') 
     219            CALL prt_ctl(tab2d_1=e_i        (:,:,2,jl)      , clinfo1= ' lim_update2  : e_i2        : ') 
     220            CALL prt_ctl(tab2d_1=e_i_b      (:,:,2,jl)      , clinfo1= ' lim_update2  : e_i2_b      : ') 
    241221            CALL prt_ctl(tab2d_1=e_s        (:,:,1,jl)      , clinfo1= ' lim_update2  : e_snow      : ') 
    242222            CALL prt_ctl(tab2d_1=e_s_b      (:,:,1,jl)      , clinfo1= ' lim_update2  : e_snow_b    : ') 
    243             CALL prt_ctl(tab2d_1=d_e_s_thd  (:,:,1,jl)/1.0e9, clinfo1= ' lim_update2  : d_e_s_thd   : ') 
    244223            CALL prt_ctl(tab2d_1=smv_i      (:,:,jl)        , clinfo1= ' lim_update2  : smv_i       : ') 
    245224            CALL prt_ctl(tab2d_1=smv_i_b    (:,:,jl)        , clinfo1= ' lim_update2  : smv_i_b     : ') 
    246             CALL prt_ctl(tab2d_1=d_smv_i_thd(:,:,jl)        , clinfo1= ' lim_update2  : d_smv_i_thd : ') 
    247225            CALL prt_ctl(tab2d_1=oa_i       (:,:,jl)        , clinfo1= ' lim_update2  : oa_i        : ') 
    248226            CALL prt_ctl(tab2d_1=oa_i_b     (:,:,jl)        , clinfo1= ' lim_update2  : oa_i_b      : ') 
    249             CALL prt_ctl(tab2d_1=d_oa_i_thd (:,:,jl)        , clinfo1= ' lim_update2  : d_oa_i_thd  : ') 
    250227 
    251228            DO jk = 1, nlay_i 
Note: See TracChangeset for help on using the changeset viewer.