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 12949 for NEMO/branches/2020/r4.0-HEAD_r12713_clem_dan_fixcpl/src/ICE/icedyn_adv_umx.F90 – NEMO

Ignore:
Timestamp:
2020-05-18T23:58:40+02:00 (4 years ago)
Author:
clem
Message:

make sure ice and snow temperatures do not overshoot after advection (both with Prather or UMx), as well as ice salinity. Include the salinity bounds (rn_simin, rn_simax) inside the icethd_sal routine. There is no more problems with super cold temperatures now. Sette tests passed.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/r4.0-HEAD_r12713_clem_dan_fixcpl/src/ICE/icedyn_adv_umx.F90

    r12832 r12949  
    9393      REAL(wp) ::   zamsk                   ! 1 if advection of concentration, 0 if advection of other tracers 
    9494      REAL(wp) ::   zdt, zvi_cen 
    95       REAL(wp), DIMENSION(1)           ::   zcflprv, zcflnow   ! for global communication 
    96       REAL(wp), DIMENSION(jpi,jpj)     ::   zudy, zvdx, zcu_box, zcv_box 
    97       REAL(wp), DIMENSION(jpi,jpj)     ::   zati1, zati2 
    98       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zu_cat, zv_cat 
    99       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zua_ho, zva_ho, zua_ups, zva_ups 
    100       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   z1_ai , z1_aip, zhvar 
    101       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zhi_max, zhs_max, zhip_max 
     95      REAL(wp), DIMENSION(1)                  ::   zcflprv, zcflnow   ! for global communication 
     96      REAL(wp), DIMENSION(jpi,jpj)            ::   zudy, zvdx, zcu_box, zcv_box 
     97      REAL(wp), DIMENSION(jpi,jpj)            ::   zati1, zati2 
     98      REAL(wp), DIMENSION(jpi,jpj,jpl)        ::   zu_cat, zv_cat 
     99      REAL(wp), DIMENSION(jpi,jpj,jpl)        ::   zua_ho, zva_ho, zua_ups, zva_ups 
     100      REAL(wp), DIMENSION(jpi,jpj,jpl)        ::   z1_ai , z1_aip, zhvar 
     101      REAL(wp), DIMENSION(jpi,jpj,jpl)        ::   zhi_max, zhs_max, zhip_max, zs_i, zsi_max 
     102      REAL(wp), DIMENSION(jpi,jpj,nlay_i,jpl) ::   ze_i, zei_max 
     103      REAL(wp), DIMENSION(jpi,jpj,nlay_s,jpl) ::   ze_s, zes_max 
    102104      ! 
    103105      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zuv_ho, zvv_ho, zuv_ups, zvv_ups, z1_vi, z1_vs  
     
    106108      IF( kt == nit000 .AND. lwp )   WRITE(numout,*) '-- ice_dyn_adv_umx: Ultimate-Macho advection scheme' 
    107109      ! 
    108       ! --- Record max of the surrounding 9-pts ice thick. (for call Hbig) --- ! 
     110      ! --- Record max of the surrounding 9-pts (for call Hbig) --- ! 
     111      ! thickness and salinity 
     112      WHERE( pv_i(:,:,:) >= epsi10 ) ; zs_i(:,:,:) = psv_i(:,:,:) / pv_i(:,:,:) 
     113      ELSEWHERE                      ; zs_i(:,:,:) = 0._wp 
     114      END WHERE 
    109115      DO jl = 1, jpl 
    110116         DO jj = 2, jpjm1 
     
    122128                  &                                               ph_s (ji+1,jj+1,jl), ph_s (ji-1,jj-1,jl), & 
    123129                  &                                               ph_s (ji+1,jj-1,jl), ph_s (ji-1,jj+1,jl) ) 
    124             END DO 
    125          END DO 
    126       END DO 
    127       CALL lbc_lnk_multi( 'icedyn_adv_umx', zhi_max, 'T', 1., zhs_max, 'T', 1., zhip_max, 'T', 1. ) 
     130               zsi_max (ji,jj,jl) = MAX( epsi20, zs_i (ji,jj,jl), zs_i (ji+1,jj  ,jl), zs_i (ji  ,jj+1,jl), & 
     131                  &                                               zs_i (ji-1,jj  ,jl), zs_i (ji  ,jj-1,jl), & 
     132                  &                                               zs_i (ji+1,jj+1,jl), zs_i (ji-1,jj-1,jl), & 
     133                  &                                               zs_i (ji+1,jj-1,jl), zs_i (ji-1,jj+1,jl) ) 
     134            END DO 
     135         END DO 
     136      END DO 
     137      CALL lbc_lnk_multi( 'icedyn_adv_umx', zhi_max, 'T', 1., zhs_max, 'T', 1., zhip_max, 'T', 1., zsi_max, 'T', 1. ) 
     138      ! 
     139      ! enthalpies 
     140      DO jk = 1, nlay_i 
     141         WHERE( pv_i(:,:,:) >= epsi10 ) ; ze_i(:,:,jk,:) = pe_i(:,:,jk,:) / pv_i(:,:,:) 
     142         ELSEWHERE                      ; ze_i(:,:,jk,:) = 0._wp 
     143         END WHERE 
     144      END DO 
     145      DO jk = 1, nlay_s 
     146         WHERE( pv_s(:,:,:) >= epsi10 ) ; ze_s(:,:,jk,:) = pe_s(:,:,jk,:) / pv_s(:,:,:) 
     147         ELSEWHERE                      ; ze_s(:,:,jk,:) = 0._wp 
     148         END WHERE 
     149      END DO 
     150      DO jl = 1, jpl 
     151         DO jk = 1, nlay_i 
     152            DO jj = 2, jpjm1 
     153               DO ji = fs_2, fs_jpim1 
     154                  zei_max(ji,jj,jk,jl) = MAX( epsi20, ze_i(ji,jj,jk,jl), ze_i(ji+1,jj  ,jk,jl), ze_i(ji  ,jj+1,jk,jl), & 
     155                     &                                                   ze_i(ji-1,jj  ,jk,jl), ze_i(ji  ,jj-1,jk,jl), & 
     156                     &                                                   ze_i(ji+1,jj+1,jk,jl), ze_i(ji-1,jj-1,jk,jl), & 
     157                     &                                                   ze_i(ji+1,jj-1,jk,jl), ze_i(ji-1,jj+1,jk,jl) ) 
     158               END DO 
     159            END DO 
     160         END DO 
     161      END DO 
     162      DO jl = 1, jpl 
     163         DO jk = 1, nlay_s 
     164            DO jj = 2, jpjm1 
     165               DO ji = fs_2, fs_jpim1 
     166                  zes_max(ji,jj,jk,jl) = MAX( epsi20, ze_s(ji,jj,jk,jl), ze_s(ji+1,jj  ,jk,jl), ze_s(ji  ,jj+1,jk,jl), & 
     167                     &                                                   ze_s(ji-1,jj  ,jk,jl), ze_s(ji  ,jj-1,jk,jl), & 
     168                     &                                                   ze_s(ji+1,jj+1,jk,jl), ze_s(ji-1,jj-1,jk,jl), & 
     169                     &                                                   ze_s(ji+1,jj-1,jk,jl), ze_s(ji-1,jj+1,jk,jl) ) 
     170               END DO 
     171            END DO 
     172         END DO 
     173      END DO 
     174      CALL lbc_lnk( 'icedyn_adv_pra', zei_max, 'T', 1. ) 
     175      CALL lbc_lnk( 'icedyn_adv_pra', zes_max, 'T', 1. ) 
    128176      ! 
    129177      ! 
     
    362410         ! --- Make sure ice thickness is not too big --- ! 
    363411         !     (because ice thickness can be too large where ice concentration is very small) 
    364          CALL Hbig( zdt, zhi_max, zhs_max, zhip_max, pv_i, pv_s, pa_i, pa_ip, pv_ip, pe_s ) 
     412         CALL Hbig( zdt, zhi_max, zhs_max, zhip_max, zsi_max, zes_max, zei_max, & 
     413            &            pv_i, pv_s, pa_i, pa_ip, pv_ip, psv_i, pe_s, pe_i ) 
    365414         ! 
    366415         ! --- Ensure snow load is not too big --- ! 
     
    15251574 
    15261575 
    1527    SUBROUTINE Hbig( pdt, phi_max, phs_max, phip_max, pv_i, pv_s, pa_i, pa_ip, pv_ip, pe_s ) 
     1576   SUBROUTINE Hbig( pdt, phi_max, phs_max, phip_max, psi_max, pes_max, pei_max, & 
     1577      &                  pv_i, pv_s, pa_i, pa_ip, pv_ip, psv_i, pe_s, pe_i ) 
    15281578      !!------------------------------------------------------------------- 
    15291579      !!                  ***  ROUTINE Hbig  *** 
     
    15391589      !! ** input   : Max thickness of the surrounding 9-points 
    15401590      !!------------------------------------------------------------------- 
    1541       REAL(wp)                    , INTENT(in   ) ::   pdt                          ! tracer time-step 
    1542       REAL(wp), DIMENSION(:,:,:)  , INTENT(in   ) ::   phi_max, phs_max, phip_max   ! max ice thick from surrounding 9-pts 
    1543       REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pv_i, pv_s, pa_i, pa_ip, pv_ip 
     1591      REAL(wp)                    , INTENT(in   ) ::   pdt                                   ! tracer time-step 
     1592      REAL(wp), DIMENSION(:,:,:)  , INTENT(in   ) ::   phi_max, phs_max, phip_max, psi_max   ! max ice thick from surrounding 9-pts 
     1593      REAL(wp), DIMENSION(:,:,:,:), INTENT(in   ) ::   pes_max 
     1594      REAL(wp), DIMENSION(:,:,:,:), INTENT(in   ) ::   pei_max 
     1595      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pv_i, pv_s, pa_i, pa_ip, pv_ip, psv_i 
    15441596      REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) ::   pe_s 
    1545       ! 
    1546       INTEGER  ::   ji, jj, jl         ! dummy loop indices 
    1547       REAL(wp) ::   z1_dt, zhip, zhi, zhs, zfra 
     1597      REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) ::   pe_i 
     1598      ! 
     1599      INTEGER  ::   ji, jj, jk, jl         ! dummy loop indices 
     1600      REAL(wp) ::   z1_dt, zhip, zhi, zhs, zsi, zes, zei, zfra 
    15481601      !!------------------------------------------------------------------- 
    15491602      ! 
     
    15511604      ! 
    15521605      DO jl = 1, jpl 
    1553  
    15541606         DO jj = 1, jpj 
    15551607            DO ji = 1, jpi 
     
    15851637                  ENDIF            
    15861638                  !                   
     1639                  !                               ! -- check s_i -- ! 
     1640                  ! if s_i is larger than the surrounding 9 pts => put salt excess in the ocean 
     1641                  zsi = psv_i(ji,jj,jl) / pv_i(ji,jj,jl) 
     1642                  IF( zsi > psi_max(ji,jj,jl) .AND. pa_i(ji,jj,jl) < 0.15 ) THEN 
     1643                     zfra = psi_max(ji,jj,jl) / zsi 
     1644                     sfx_res(ji,jj) = sfx_res(ji,jj) + psv_i(ji,jj,jl) * ( 1._wp - zfra ) * rhoi * z1_dt 
     1645                     psv_i(ji,jj,jl) = psv_i(ji,jj,jl) * zfra 
     1646                  ENDIF 
     1647                  ! 
    15871648               ENDIF 
    15881649            END DO 
    15891650         END DO 
    15901651      END DO  
     1652      ! 
     1653      !                                           ! -- check e_i/v_i -- ! 
     1654      DO jl = 1, jpl 
     1655         DO jk = 1, nlay_i 
     1656            DO jj = 1, jpj 
     1657               DO ji = 1, jpi 
     1658                  IF ( pv_i(ji,jj,jl) > 0._wp ) THEN 
     1659                     ! if e_i/v_i is larger than the surrounding 9 pts => put the heat excess in the ocean 
     1660                     zei = pe_i(ji,jj,jk,jl) / pv_i(ji,jj,jl) 
     1661                     IF( zei > pei_max(ji,jj,jk,jl) .AND. pa_i(ji,jj,jl) < 0.15 ) THEN 
     1662                        zfra = pei_max(ji,jj,jk,jl) / zei 
     1663                        hfx_res(ji,jj) = hfx_res(ji,jj) - pe_i(ji,jj,jk,jl) * ( 1._wp - zfra ) * z1_dt ! W.m-2 <0 
     1664                        pe_i(ji,jj,jk,jl) = pe_i(ji,jj,jk,jl) * zfra 
     1665                     ENDIF 
     1666                  ENDIF 
     1667               END DO 
     1668            END DO 
     1669         END DO 
     1670      END DO 
     1671      !                                           ! -- check e_s/v_s -- ! 
     1672      DO jl = 1, jpl 
     1673         DO jk = 1, nlay_s 
     1674            DO jj = 1, jpj 
     1675               DO ji = 1, jpi 
     1676                  IF ( pv_s(ji,jj,jl) > 0._wp ) THEN 
     1677                     ! if e_s/v_s is larger than the surrounding 9 pts => put the heat excess in the ocean 
     1678                     zes = pe_s(ji,jj,jk,jl) / pv_s(ji,jj,jl) 
     1679                     IF( zes > pes_max(ji,jj,jk,jl) .AND. pa_i(ji,jj,jl) < 0.15 ) THEN 
     1680                        zfra = pes_max(ji,jj,jk,jl) / zes 
     1681                        hfx_res(ji,jj) = hfx_res(ji,jj) - pe_s(ji,jj,jk,jl) * ( 1._wp - zfra ) * z1_dt ! W.m-2 <0 
     1682                        pe_s(ji,jj,jk,jl) = pe_s(ji,jj,jk,jl) * zfra 
     1683                     ENDIF 
     1684                  ENDIF 
     1685               END DO 
     1686            END DO 
     1687         END DO 
     1688      END DO 
    15911689      ! 
    15921690   END SUBROUTINE Hbig 
Note: See TracChangeset for help on using the changeset viewer.